基于傅里叶梅林变换的宽带非线性调频信号时差/尺度差估计方法

张 亮， 于洪波， 孙殿星， 程树强， 张翔宇， 王国宏

(1.海军航空大学，山东 烟台 264000； 2.中国人民解放军94326部队，济南 250000；3.中国人民解放军94201部队，济南 250000)

0 引言

电子战是信息化战争的重要组成部分，对目标的精确定位是电子战的一项重要内容。目标定位系统通常可分为有源定位系统和无源定位系统，典型的有源定位系统为雷达，由于需要对外辐射电磁信号以获取目标回波，战场上容易被敌方侦察而遭受反辐射武器打击。相较而言，无源定位系统自身不辐射电磁信号，隐蔽性好，典型体制为双星无源定位[1]。宽带条件下，双星无源定位系统根据两路接收信号估计时差(TDOA)和尺度差(SDOA)，再利用卫星、目标和地球的几何关系对辐射源进行精确定位。为估计时差/尺度差，传统方法[2-3]是计算两路接收信号宽带互模糊函数(Wideband Cross Ambiguity Function,WBCAF)，该方法能够同时估计2个参数，且适用于多种源信号形式，而存在的问题是计算量大。WBCAF计算过程包含两步，首先计算一路信号的尺度版本，然后将其与另一路信号进行相关处理，由于第2步可利用FFT，IFFT在频域快速实现，因此WBCAF快速实现的难点为第1步，针对该问题，文献[2]提出多速率采样重构滤波方法，但是当尺度因子为多位小数时，该方法无法获得精确的尺度版本。部分研究考虑仅在变换域计算WBCAF，以规避接收信号尺度版本获取难题，其中，文献[3]利用小波变换内积定理，提出基于交互小波变换的WBCAF实现方法，但母小波的选取对计算结果影响较大，文献[4]进一步提出基于梅林变换的WBCAF实现方法，但梅林变换的数值计算方法有多种[5-9]，所用方法计算精度不高。另外，如果源信号形式已知，时差/尺度差估计问题会变得相对容易，例如针对线性调频(LFM)信号，文献[10-11]利用分数阶傅里叶变换(Fractional Fourier Transform,FRFT)测量两路接收信号调频斜率，根据两路信号调频斜率比例关系估计尺度差，再对第1路信号“拉伸”或者“压缩”，与第2路信号进行相关处理估计时差。该改进方法通过两步级联处理降低了运算整体复杂度，低信噪比(SNR)条件下2个参数估计均方误差逼近克拉美罗下界(Cramer-Rao Lower Bound,CRLB)，遗憾的是当源信号为非线性调频(NLFM)信号时，该方法不适用。

针对上述问题，本文引入时间尺度(TS)、尺度估计(SE)基本概念，利用两路信号频谱包络是彼此尺度版本的特点，提出基于傅里叶梅林变换的宽带非线性调频信号时差/尺度差估计方法。

1 传统与改进方法

宽带条件下，利用双星无源定位系统对地面一个固定辐射源进行定位，以第1路信号为参考基准，则2个独立接收机接收到的信号可表示为[10-11]

(1)

式中：x1(t)，x2(t)分别为第1路和第2路接收信号；A1，A2为源信号s(t)到达2个接收机的路径增益；Δt为两路信号时差(TDOA)；Δα为尺度差(SDOA)，与卫星相对辐射源的径向速度及光速有关；n1(t)，n2(t)分别为2个接收机基底复高斯白噪声。辐射源相对卫星的径向速度体现为接收信号尺度的变化，传统估计方法[2-3]是利用宽带互模糊函数(WBCAF)实现。连续信号z(t)与x(t)的WBCAF定义为

(2)

式中：γ为WBCAF尺度因子；τ为时移。令z(t)=x2(t)，x(t)=x1(t)，易知，当γ=1+Δα，τ=Δt时，|Xzx(γ,τ)|

取得最大值，峰值搜索可得到时差和尺度差(需减1)。WBCAF可表示成卷积形式频域快速实现，即

Xzx(γ,τ)=z(t)⊗x*(-γt)=F-1{F[z(t)]F*[x(γt)]}

(3)

式中：⊗为卷积符号；F[·]，F-1[·]分别为傅里叶变换与逆变换表示符号。由式(3)可知，为计算两路接收信号WBCAF，需要设定尺度因子搜索区间和间隔，然后根据不同的尺度因子对第1路信号进行“拉伸”或“压缩”，搜索区间的扩大会导致运算量成倍增加。由式(1)可知，如果源信号为单频信号，则尺度上的差异反映为接收信号频率上的不同，如果源信号为LFM信号，则反映为调频斜率上的不同，通过估计两路信号频率或者调频斜率可得SDOA，然后根据估计的SDOA对x1(t)进行“拉伸”或者“压缩”，再与x2(t)进行时域相关处理可得TDOA，该改进方法虽然能够显著减少运算量，但当源信号为NLFM信号时则不适用。雷达、通信系统常用NLFM信号包括二次调频(QFM)、指数调频(EFM)、双曲调频(HFM)和正弦调频(SFM)等[12]，复数形式依次为

(4)

(5)

(6)

(7)

式中：T为时宽，f1,f2分别为4种信号的起始频率和截止频率;T0为HFM延时时间。对式(4)～(7)相位求导可得QFM,EFM,HFM和SFM瞬时频率分别为

(8)

(9)

(10)

(11)

2 时间尺度与尺度估计

为介绍所提方法，本章引入2个基本概念，即时间尺度(TS)和尺度估计(SE)，然后利用梅林变换(MT)的特例尺度变换(ST)给出TS和SE快速实现方案。

2.1 基本概念

所谓时间尺度，即连续信号x(t)到y(t)的映射过程[13]，可表示为

(12)

(13)

由式(12)可知，当0<α<1时，信号时域扩张，当α>1时，信号时域收缩。根据傅里叶变换尺度性质可知

(14)

式中，Y(f)，X(f)分别为y(t)，x(t)的傅里叶变换，f为频率自变量。可以看出，对信号进行TS还会导致带宽成比例减小(0<α<1)或增大(α>1)，但信号时宽-带宽积不变[14]。设连续信号z(t)=TSα1[x(t)]，α1为尺度因子。所谓尺度估计，即z(t)和x(t)已知时对α1的估计，可表示为

(15)

(16)

式中：◇为尺度互相关符号；*表示共轭转置。由式(2)可知，当τ=0时，WBCAF与尺度互相关函数近似，区别在于不同的积分区间。

2.2 实现方法

尺度变换(ST)是梅林变换(MT)复自变量实部取0.5的特例，表达式[13]为

(17)

式中：S{·}为ST表示符号；Df(c)为信号f(t)的ST，c为尺度自变量。逆尺度变换(Inverse Scale Transform,IST)为

(18)

式中,S-1{·}为IST表示符号。令t=eu，代入式(17)，得到

(19)

(20)

IST同样可利用IFFT快速实现，得到快速逆尺度变换(Inverse Fast Scale Transform,IFST)，计算过程与FST相反，先计算Df(c)的逆傅里叶变换，再进行对数采样。为解决TS快速实现问题，引用ST尺度不变性，即

Dz(c)=eicln α1Dx(c)

(21)

式中,Dx(c),Dz(c)分别为x(t),z(t)的ST，两者包络相同、相位不同。计算式(21)的IST，易得

z(t)=S-1{eicln α1Dx(c);t} 。

(22)

利用式(22)可得x(t)在某一尺度因子下的TS后信号。SE同样利用ST快速实现[15]，具体为

(23)

3 基于傅里叶梅林变换的时差/尺度差估计

3.1 算法原理

为介绍本文所提方法，暂不考虑接收机噪声，重写式(1)为

(24)

利用傅里叶变换时移不变性消除时差影响，容易得到

(25)

式中，X1(f)，X2(f)分别为x1(t)，x2(t)的傅里叶变换。易知

(26)

可以看出，|X2(f)|实际为|X1(f)|的尺度版本，尺度因子为1/(1+Δα)，当|X1(f)|和|X2(f)|均已知时，对1/(1+Δα)的估计属SE概念，利用式(23)可估计1/(1+Δα)，进而得到尺度差Δα，然后对x1(t)进行时间尺度(尺度因子取1+Δα)，与x2(t)进行时域相关处理可进一步估计TDOA。为直观显示，图1给出了本文所提方法实现流程。

图1 时差/尺度差估计流程

3.2 计算复杂度

首先，分析时间尺度的计算复杂度。TS实现过程包含“FST+相位修正+IFST”3个环节，计算复杂度为O[(NlnN)lb(NlnN)]，总共需要进行3NlnN+NlnN·lb(NlnN)次复数乘法。其次，分析SE计算复杂度。SE实现的关键在于计算尺度互相关函数，具体需要进行2次FST，1次复数乘积和1次IFST，计算复杂度同样为O[(NlnN)lb(NlnN)]，总共需要进行4NlnN+1.5·NlnNlb(NlnN)次复数乘法。再次，分析传统基于宽带互模糊函数的时差/尺度差估计方法计算复杂度。对于1个尺度因子，为计算式(3)需要进行1次TS，2次FFT，1次 IFFT和1次点乘，总共需要进行MN+1.5MN·lbN+3MNlnN+MNlnNlb(NlnN)次复乘运算，计算复杂度为O[(MNlnN)lb(NlnN)]，M为尺度因子取值个数。最后，分析所提方法计算复杂度。为估计尺度差需要进行2次FFT和1次SE，而估计时差需要进行1次TS和1次时域相关处理(2次FFT，1次IFFT和1次点乘)，因此，总共需要进行N+2.5NlbN+7NlnN+2.5NlnN·lb(NlnN)次复数乘法运算，计算复杂度为O[(NlnN)·lb(NlnN)]。为直观显示，设M分别为10，100和1000，图2给出了复乘次数随信号采样点数变化曲线。所提方法计算量明显低于传统方法。上述分析未考虑峰值搜索环节，由于传统方法需要进行二维搜索，实际运算量会更大。

图2 复乘次数随信号采样点数变化曲线

4 仿真结果与分析

4.1 参数设置

辐射源发射NLFM信号，时宽2.5 μs，带宽250 MHz，采样频率1 GHz，双星无源定位系统接收的两路NLFM信号时差/尺度差分别为3.75 μs和0.1，QFM，EFM，HFM和SFM起始频率和截止频率分别为10 MHz和260 MHz，HFM延迟时间为1个时域采样间隔。4种信号时频分布如图3所示。

图3 NLFM信号时频分布

4.2 可行性分析

以SFM信号时差/尺度差估计为例，利用2个仿真试验依次对第1章所提WBCAF实现方法和第3章所提基于傅里叶梅林变换的时差/尺度差估计方法的可行性进行验证。为展示细节，接收信号中均不添加噪声(4.3节将进行抗噪效能分析)。

4.2.1 仿真试验1

设接收信号时长为10 μs，两路接收信号时域波形和频谱如图4所示。第2路信号中的SFM信号相比第1路信号中的SFM信号存在时域收缩和频谱扩张现象，与2.1节中的TS概念一致。

图4 接收信号时域波形和频谱

设WBCAF尺度因子取值0.5～1.5，间隔0.01，利用第1章方法，计算两路接收信号宽带互模糊函数，结果如图5所示。

图5 接收信号宽带互模糊函数

由图5可以看出，两路信号WBCAF呈明显冲击特性，峰值搜索可得时差和尺度差分别为3.745 μs，0.096，与仿真所设3.75 μs，0.1基本一致。

需要注意的是，宽带互模糊函数X轴实际为尺度因子，与尺度差存在减1的关系，图5尺度差坐标已进行减1处理。

4.2.2 仿真试验2

按照图1流程，对所提基于傅里叶梅林变换的时差/尺度差估计方法可行性进行验证。首先，利用2.2节方法计算图4(b)两路信号频谱尺度互相关函数(SCF)结果,如图6所示。可以看出，对频谱取包络时，估计得到的频域尺度因子为0.909 2，对应时域尺度因子为1.099 9，尺度差为0.099 9，与仿真所设尺度差一致，而不对频谱取包络，估计结果完全错误。

图6 接收信号频谱尺度互相关函数

根据估计的时域尺度因子，对第1路信号进行时间尺度，消除尺度差对时差估计的影响，接收信号时域波形和频谱如图7所示。

图7 时间尺度后的接收信号时域波形和频谱

两路接收信号中的SFM时宽、带宽基本一致，不再存在收缩、扩展现象。以第1路信号为匹配信号，对第2路信号进行相关处理，结果如图8所示，峰值搜索估计得到时差为3.75 μs，与仿真所设参数一致。

图8 相关处理结果

4.3 抗噪效能分析

对传统基于宽带互模糊函数时差/尺度差估计方法[2-3]和所提基于傅里叶梅林变换的时差/尺度差估计方法的抗噪效能进行验证，使用的评估指标为平均相对误差(MRE)[16-17]，同时,为验证所提算法的适用范围，选取LFM信号(脉宽、初始频率、截止频率同4.1节),4种NLFM信号(QFM,EFM,HFM,SFM)和零均值高斯白噪声信号(脉宽同4.1节)作为源信号。SNR取值-10～10 dB，间隔1 dB，运行蒙特卡罗仿真500次，得到尺度差估计平均相对误差随信噪比变化曲线如图9所示，图10为时差估计平均相对误差随信噪比变化曲线。

图9 尺度差估计平均相对误差随信噪比变化曲线

图10 时差估计平均相对误差随信噪比变化曲线

传统方法与所提方法对不同信号形式时差/尺度差估计效能不尽相同，总体上看，传统方法更适合低SNR条件，当SNR大于-4 dB时，对LFM，QFM，EFM和SFM这4种信号形式尺度差估计MRE小于5%，时差估计MRE小于1%，而所提方法则要求SNR大于5 dB，需要注意的是，高SNR条件下所提方法估计精度要高于传统方法，原因是计算WBCAF需人为设定尺度因子取值范围和间隔(本节仿真中，尺度因子取0.5～1.5，间隔0.01)，当取值范围一定时，降低间隔能够提高估计精度，但同时也意味着更大的计算量，所提方法不存在此问题。

同时可以看出，传统方法对HFM尺度差估计效果极不理想，不同SNR条件下尺度差估计MRE大于13%，估计误差明显高于所提方法，分析原因是HFM具备“多普勒不变性”，其宽带模糊函数脊线近似为线性，脊线幅度较为接近[18-19]，峰值搜索无法准确估计时差/尺度差，而所提方法是对信号频谱包络进行尺度估计，实现原理不同。传统方法与所提方法对噪声信号尺度差估计MRE小于5%的临界SNR分别为2 dB和3 dB，时差估计MRE小于1%的临界SNR分别为1 dB和2 dB，两种方法估计效能相差不大。

5 结束语

针对双星无源定位系统中的时差/尺度差估计问题，提出基于傅里叶梅林变换的宽带非线性调频信号时差/尺度差估计方法，核心思想是利用傅里叶变换“时移不变性”、梅林尺度变换“尺度不变性”对时差/尺度差进行级联估计。仿真表明，当SNR大于5 dB时，所提方法能够准确估计4种宽带非线性调频信号时差和尺度差，时差估计平均相对误差小于1%，尺度差估计平均相对误差小于5%；相比传统基于宽带互模糊函数的时差/尺度差估计方法，所提方法计算复杂度更低；相比基于FRFT的时差/尺度差估计方法，所提方法可适用LFM，QFM，EFM，HFM和SFM，甚至噪声脉冲，理论上，所提方法可适用任意信号形式，限于篇幅不做进一步仿真。所提方法要求SNR大于5 dB，低SNR条件下的时差/尺度差估计问题是下一步研究工作的重点。