1.   引言

信号源的无源定位，具有作用距离远、隐蔽性高、生存能力强的特点，受到国内外广泛的研究^[1-9]。传统的单站无源定位技术有到达方向(DOA)估计^[10]、到达时间(TOA)估计^[11]和TOA/DOA联合估计^[12]。早期的DOA估计主要通过波束形成法来完成，受瑞利限的限制，后来发展的许多高分辨DOA估计算法，如Capon算法、Pisarenko算法和MUSIC算法等极大地提高了估计的分辨力和精度^[13]，但是DOA估计的分辨能力始终受到天线阵列孔径的限制，空间源数的估计也影响部分DOA估计算法的性能。TOA估计算法测量信号到达观测平台的时间估计信号源的距离，检测门限的设置是TOA的关键，不合适的检测门限会导致TOA估计误差的加大^[14]。TOA/DOA联合估计通过测量到的角度和距离进行目标的定位，由于DOA分辨能力的固有限制，信号源距离越远，估计性能越差^[15]。

为了减小测量误差，提高定位精度，国际上众多学者又发展了多个观测站测量的方法：到达时间差(TDOA)和到达频率差(FDOA)^[16]。多站TDOA定位通过测量信号到达不同观测平台的时间差，进行双曲面定位^[17]，多站TDOA标准算法需要观测平台间时钟同步，时间不同步产生的误差，会导致定位的误差；TDOA/FDOA联合定位技术利用时差曲面和频差曲面进行定位^[18]，只需要两个平台。并且在星载应用中，卫星移动速度快，多普勒频差大，定位精度。目前已经提出了很多TDOA/FDOA联合定位的算法，要达到这些算法的定位精度，需要很高的时差频差估计精度^[19]。目前宽带信号的FDOA估计的研究主要基于CAF^[20]算法，CAF算法是一种双星联合估计算法，在信号带宽不满足窄带假设的情况下，估计效果变差。

针对这些问题，本文提出一种基于短合成孔径的双星干涉的时频差精确估计方法，该方法结合合成孔径的思想，将SAR成像中的相干积累的原理引入到时差/频差估计中。在短孔径的条件下，接收信号的时差线性变化，通过最小二乘拟合提高时差估计精度；信号的多普勒频率表现为方位向的线性相位，通过多个脉冲相干积累，提高了频差的测量精度。对于窄带信号，两个卫星可以单独处理数据，单独处理得到的多普勒频率做差得到频差；对于宽带信号，可以通过双通道联合处理获得高的估计精度，有效解决CAF算法的不足。对于窄带信号，本文算法可以获得单个通道的多普勒频率，因此，本文算法也适用于窄带信号的多普勒频率估计，进行信号的到达角估计。本文的内容安排如下：第2节介绍了短合成孔径的双星信号模型，第3节介绍了时频差精确估计算法及定位方法，第4节分析了算法的性能并进行了算法的仿真，第5节通过STK仿真数据验证了算法的有效性。

2.   短合成孔径的双星信号模型

合成孔径双星的信号模型如图1所示，为了表述的方便，这里假设两个卫星的轨道高度相同，速度均为$v$。在短合成孔径中心时刻，两个卫星运动方向与信号源P所在方向的夹角分别为${\theta _1}$, ${\theta _2}$。卫星的轨迹如图1中粗实线所示。利用斜视SAR中的斜距模型，点目标到卫星1的斜距历程可以写成，${R_1}({t_{\rm a}}) \approx$${R_0} - v\cos {\theta _1}({t_{\rm a}} - {t_{\rm pc1}})+ \dfrac{{{v^2}{{\sin }^2}{\theta _1}}}{{2{R_0}}}{({t_{\rm a}} - {t_{\rm pc1}})^2}$，其中${t_{\rm pc1}} = ({X_{\rm p}} - {R_{\rm b}}\cot {\theta _1})/v$为合成孔径中心时刻，${R_0} = {R_{\rm b}}\sqrt {1 + {{\cot }^2}{\theta _1}}$为合成孔径中心时刻的斜距，${R_{\rm b}}$为信号源到航迹的最近距离。短合成孔径的条件为斜距历程中的2次项引起的相位变化小于$\pi /4$，即$\dfrac{{2\pi {f_{\rm c}}}}{{\rm{c}}}\dfrac{{{v^2}{{\sin }^2}{\theta _1}}}{{2{R_0}}}{({t_{\rm a}} - {t_{\rm pc1}})^2} < \dfrac{\pi }{4}$，其中${\rm{c}}$为光速。由此可以得到短合成孔径的时间限制：${t_{\rm syn}} < \sqrt {\dfrac{{\lambda {R_0}}}{{{v^2}{{\sin }^2}{\theta _1}}}}$。

图 1  双星信号模型

下载: 全尺寸图片 幻灯片

由上述分析可以得到，短合成孔径时间的条件下，${R_1}({t_{\rm a}})$线性变化，假设信号源以${T_{\rm r}}$为脉冲重复间隔发送信号，每个脉冲的提前量为$m{T_{\rm r}}v\cos {\theta _1}/{\rm{c}}$。考虑到信号的多普勒频率${f_{\rm dc1}} = v\cos {\theta _1}{f_{\rm c}}/{\rm{c}}$，每个脉冲的提前量可以表示为$m{T_{\rm r}}{f_{\rm dc1}}/{f_{\rm c}}$。实际测量的脉冲重复周期为${T_{\rm r1}}$，一般情况下${T_{\rm r1}} \ne {T_{\rm r}}$，接收的信号写成2维的形式，方位慢时间${t_m} = m{T_{\rm r1}}$，信号提前量为${t_m}\dfrac{{{T_{\rm r}}}}{{{T_{\rm r1}}}}\dfrac{{{f_{\rm dc}}}}{{{f_{\rm c}}}}$，由于脉冲重复周期估计误差带来的信号时延量为${t_m}(1 - {T_{\rm r}}/{T_{\rm r1}})$，定义$\alpha = \left(1 - \dfrac{{{T_{\rm r}}}}{{{T_{\rm r1}}}}\right) - \dfrac{{{T_{\rm r}}{f_{\rm dc}}}}{{{T_{\rm r1}}{f_{\rm c}}}}$为方位延时因子，延时量可以表示为$\alpha {t_m}$。假设发射信号为线性调频信号，卫星1接收的信号为

其中，${t_1}$表示信号出现的初始时刻，${f_{\rm c1}}$为估计的载频，$-{\rm j}2\pi {f_{\rm c1}}t$项表示去载频，$\gamma$为调频率。两个卫星同时开始录取数据，卫星2接收的的信号表达式为

$\Delta t$表示两个卫星由于斜距不同产生的时间差。在较短的时间内，假设$\Delta t$保持不变。双星干涉要求系统是同步的，并且时钟是短期稳定的。双星系统间的不同步会导致$\Delta t$的变化，从而导致时差估计产生误差。在短合成孔径的时间内，通常小于1 s，时钟的稳定度优于$8 \times {10^{ - 14}}$，时间同步能达到0.15 ns^[21]。

3.   短孔径双星定位算法

3.1   双星干涉测频差

短合成孔径条件下，信号的包络移动量很小，所以可以忽略包络的影响，对于线性调频信号来说，那么其信号可以简化表示为

式(3)中第1项的相位是${t_m}$的线性函数，取两列信号相关后，可以得到其斜率${\alpha _1}$, ${\alpha _1}$中包含信号的真实多普勒频率。对信号进行方位相关处理得到

为了利用全部数据并且兼顾减小计算量，可以取${t_0} = \dfrac{M}{2}{T_{\rm r1}}$, $M$为总的脉冲个数。相关处理后，式(4)的第1项为常数项，第2项为线性调频信号方位相关后的残余相位。将式(4)其重新写为

其中，${f_{\rm c}}{\alpha _1}$为信号多普勒的估计值${\hat f_{\rm dc1}}$，从式(5)可以看出对于窄带而言，可以得到$\gamma (t - {t_1}) \ll {f_{\rm c}}$和$\gamma ({t_0} + 2{t_m}) \ll {f_{\rm c}}$，并且信号多普勒频率和载频的比值${f_{\rm c}}{\alpha _1}/{f_{\rm c}}$一般也很小，所以后一项的相位很小，可以忽略，对式(5)取相位，并通过求快时间和慢时间的平均值，于是得到

双星处理的数据相减即可得到真实的频差。

对于宽带线性调频信号，上述近似不再成立。主要是因为大的调频率$\gamma$导致式(5)的第2项相位不能忽略，考虑到宽带线性调频信号在频域的表现形式相同，可以通过两颗卫星的数据做干涉消除调频率$\gamma$的影响。通过驻定相位原理，将接收到的每个脉冲变换到频域，可以得到

可以看到，距离频域中信号的包络相同，第1项相位为线性调频项，此项在两个卫星接收的信号中相同，第2项为载频决定的相位，第3项为随距离频率变化的相位，这两项均为方位慢时间的线性函数。根据式(2)，第2个卫星接收信号的频域形式会在${t_1}$的基础上额外引入一个固定时间$\Delta t$。由于线性调频项在不同卫星中相同，可以通过双星数据干涉去掉第1项相位，对两个卫星的数据进行干涉可以得到

式(8)中的第1项为频差有关的相位，为方位慢时间的1次函数，并随距离频率线性变化，第2项为两个卫星的时差带来的额外相位，在不同的方位时刻其为定值，可以通过方位相关去掉此项。取双星数据干涉后的两列信号进行方位相关可以得到

其中，${t_0}$为已知量，可以取${t_0} = \dfrac{M}{2}{T_{\rm r1}}$。可以看到两颗卫星的数据干涉后的结果只包含频差项，其随距离频率线性变化。取相位操作可以得到精确的频差。对式(9)取相位得到

当载频的估计误差较小时，通过2维平均可以得到频差$\Delta {f_{\rm dc}} \approx \dfrac{1}{{2\pi {t_0}}}\angle {C_{F12}}(f,{t_0})$。算法流程如图2所示。

图 2  频差测量算法流程

下载: 全尺寸图片 幻灯片

3.2   相关测时差

时间差可以通过双星数据包络相关进行测量，两个矩形窗进行相关得到的结果是一个三角窗函数，通过搜索峰值，可以得到两个卫星的时差。两个卫星的数据进行包络相关处理，得到

当两个信号重合的时候，取得最大值，通过峰值搜索，即可得到峰值的位置

峰值出现的位置包括两部分，一部分是两个卫星由于距离差产生的时差$\Delta t$，另一部分是两个卫星多普勒频率不同产生的每个脉冲的提前量的差。时差为

时差是慢时间${t_m}$的线性函数，通过最小二乘法即可求得斜线的斜率，可以得到${t_m}$时刻的时差。

3.3   双星时差频差定位算法

时差定位技术通过时差参数，获得定位曲面，3星时差定位通过两个时差参数，获得两个定位曲面，与地球表面一起确定信号源的位置。频差定位技术与时差定位相同，通过频差曲面定位信号源的位置。双星时差频差定位通过两颗卫星确定一个时差曲面和一个频差曲面，与地球表面一起确定信号源的位置，与单纯的时差定位技术和频差定位技术相比，减少了对卫星个数的要求。

如图3所示，采用空间直角坐标系，两个卫星的位置为${S_1}$和${S_2}$，坐标为${({x_1},{y_1},{z_1})^{\rm{T}}}$和${({x_2},{y_2},{z_2})^{\rm{T}}}$，速度矢量为${{{v}}_1} = {({v_{x1}},{v_{y1}},{v_{z1}})^{\rm{T}}}$和${{{v}}_2} = {({v_{x2}},{v_{y2}},{v_{z2}})^{\rm{T}}}$，目标的位置为${{R}} = {(x,y,z)^{\rm{T}}}$，两个卫星编队飞行，它们的速度相同${{{v}}_1} \approx {{{v}}_2} \approx {{v}}\underline{\underline \Delta } {({v_x},{v_y},{v_z})^{\rm{T}}}$。利用时差参数、频差参数可得方程组

图 3  双星定位模型

下载: 全尺寸图片 幻灯片

其中，$\Delta {v_r}$为两颗卫星相对于信号源的径向速度差，$\Delta r$为两颗卫星到信号源的距离差，$R$为地球半径。这是一个3元高次非线性的方程组，文献[22]将该方程组近似简化为解4次1元方程组，得到方程的解析解，解方程组后，可得得到4组解，去除虚根和模糊根即可得到信号源的位置。

4.   仿真及性能分析

4.1   频差误差分析

载频和脉冲重复周期的估计误差会引起多普勒频差估计的误差，从式(10)可以得到双星干涉得到的多普勒频差${\hat f_{\rm dc}}$为$\dfrac{{{f_{\rm c1}}}}{{{f_{\rm c}}}}\dfrac{{{T_{\rm r}}}}{{{T_{\rm r1}}}}\Delta {f_{\rm dc}}$，当载频和脉冲重复周期的估计有误差时，${\hat f_{\rm dc}}$可以写为

通常${T_{\rm r1}} - {T_{\rm r}}$约为1～5个采样点，而地对空雷达的${T_{\rm r1}}$为20000～200000个采样点，由重复周期估计误差引起的频差估计误差最多为$\Delta {f_{\rm dc}}/4000$，${f_{\rm c}}$的估计误差一般在100 Hz以内，当载频在GHz时，误差小于$\Delta {f_{\rm dc}}/{10^7}$，几乎可以忽略不计。载频越高，载频估计误差对频差估计的影响越小，但是考虑到短孔径假设条件，载频高波长短，短孔径的合成时间${t_{\rm syn}}$变短。因此在实际应用中，信号的载频未知，为了适应不同载频信号的频差测量，数据录取的时间应小于极限合成孔径时间${t_{\rm syn}}$，并留有余量。同时对于窄带信号，高的载频能更充分的满足3.1节信号带宽远小于载频的假设，使窄带频差估计的算法更有效。

仿真采样率为100 MHz，脉宽为50 μs，带宽90 MHz，重复周期为100 μs，载频为1.5 GHz，多普勒频率为115853.4 Hz和123526.3 Hz，合成孔径时间为0.2 s。图4为带有噪声的情况下的干涉相位，噪声为高斯噪声，SNR为20 dB，图5为某个脉冲的距离向剖面图，可以看出干涉相位随着频率线性变化。

图 4  有噪声的情况下干涉的相位图

下载: 全尺寸图片 幻灯片

图 5  某个脉冲距离向剖面图

下载: 全尺寸图片 幻灯片

噪声的情况下干涉相位随着频率变化，通过线性拟合可以得到干涉相位随频率变化的斜线，用所得的拟合斜线对图4的干涉相位进行补偿，经过多次相位补偿和拟合最终可以得到相位随距离变化的结果。图4有噪声情况下的相位还可以采用最大似然估计法得到中心的相位，从而得到估计的差频。

脉冲积累数为100时，通过CAF算法和本文算法，得到频差的均方根误差(RMSE)随信噪比的变化情况，如图6所示。从图6中可以看出，本文算法有较好的宽带信号处理能力，误差可以达到10 Hz以内，并且信噪比越高，均方根误差越小。

图 6  LFM信号频差均方误差值随信噪比变化

下载: 全尺寸图片 幻灯片

本文所提算法同样适用于单频信号的测量，采样率为100 MHz，脉宽为50 μs，重复周期为100 μs，载频为1.5 GHz，多普勒频率为115853.4 Hz和123526.3 Hz，通过1000次蒙特卡洛实验，得到基于FFT的Rife算法和本文算法的均方根误差随信噪比变化的情况如图7所示。

图 7  单频信号频差均方误差值随信噪比变化

下载: 全尺寸图片 幻灯片

4.2   时差误差分析

由式(13)可知时差是慢时间${t_m}$的线性函数，斜率为${T_{\rm r}}\Delta {f_{\rm dc}}/({T_{\rm r1}}{f_{\rm c}})$，时差主要用过包络的相关获得，影响信号包络的主要因素是信噪比，信噪比在10 dB时，时由于噪声引起的时差误差不超过30 ns。双星TDOA的估计方法为检测每个卫星接收的脉冲信号的TOA，双星得到的结果相减获得双星的TDOA。常用方法有^[23]：基于单点滑动乘法累加的TOA估计，基于分段DFT的TOA估计，基于自相关的TOA估计。

仿真采样率为100 MHz时，通过1000次蒙特卡洛实验，得到常规算法和本文提出算法的均方根误差(RMSE)随信噪比变化的情况如图8所示。

图 8  时差均值随信噪比的变化

下载: 全尺寸图片 幻灯片

4.3   定位误差分析

双站TDOA/FDOA定位方程中，影响定位精度的测量误差有：TDOA/FDOA测量误差，卫星的三轴位置误差和速度误差。在TDOA/FDOA定位系统中，定位误差对FDOA的测量精度更敏感，要达到1 km的定位精度，需要Hz级的频差测量精度。当信噪比为5 dB时，在1$\sigma$原则下，时差误差为17.78 ns，频差误差为3.134 Hz，卫星的3维位置误差为0.5 m，卫星的3维速度误差为0.1 m/s。图9是定位的几何精度因子(Geometric Dilution Of Precision, GDOP)图，在卫星飞行方向两侧较大的范围内，定位误差在1 km以内。

图 9  双星时差频差定位的GDOP (km)

下载: 全尺寸图片 幻灯片

5.   STK仿真数据处理

5.1   STK仿真数据处理

通过STK软件进行卫星轨道和信号源3维位置的仿真，将仿真的轨道数据和信号源位置用于信号的生成，仿真数据中加入高斯噪声，信噪比为10 dB。

表1为通过本文的方法处理STK仿真数据的结果，第1组数据的载频为1.59 GHz，脉宽10 μs，脉冲重复周期100 μs，信号的带宽10 MHz，合成孔径时间为0.2 s；第2组数据的载频为1.78 GHz，其他参数与第1组信号参数相同。从测量的结果可以看出，单个卫星的基频测量误差为30～40 Hz，误差较大。通过双星干涉测量的频差，误差在10 Hz以内。STK仿真数据处理结果表明，双星干涉测量频差的结果更加精确。

表 1  STK仿真数据处理结果

载频(GHz)	脉宽(μs)	脉冲重复周期(μs)	带宽(MHz)	基频测量结果(Hz)	基频测量误差(Hz)	相对误差(Hz)
1.59	10	100	10	–4781.5	33.9742	2.16
1.59	10	100	10	1710.5	31.8189
1.78	10	100	10	–1120.3	40.6852	4.20
1.78	10	100	10	2486.2	36.4891

下载: 导出CSV 

| 显示表格

通过CAF对两组数据进行频差估计，第1组数据估计的频差为6531.8 Hz，误差为41.93 Hz；第2组数据估计的频差为6366.4 Hz，误差为31.29 Hz。

5.2   定位结果

利用上述频差估计的方法，进行定位误差的仿真。地球模型采用球体模型，地球半径为6378.14 km，卫星轨道高度为1000 km。两颗卫星相距100 km，两个卫星的坐标为(单位m)：卫星1(7377830.58, 0, 50000)，卫星2(7377830.58, 0, –50000)，卫星速度(m/s): (0, 7000, 0)。辐射源目标位于地球表面，位置参数如表2所示。通过测量时差和频差，代入定位方程，求解目标的位置，目标的定位结果如表2所示，定位误差小于1 km。

表 2  时差频差定位结果

载频(GHz)	脉宽(μs)	脉冲重复周期(μs)	带宽(MHz)	目标真实位置(m)	目标测量位置(m)	相对误差(m)
1.59	10	100	10	(6303534, –92215, 968286)	(6303532, –92876, 968234)	662
1.59	10	100	10	(6348796, –27959, 610439)	(6348797, –27448, 610446)	510
1.78	10	100	10	(6339522, –16218, 700588)	(6339520, –16829, 700590)	611
1.78	10	100	10	(6374947, 30978, 199294)	(6374945, 30530, 199409)	462

下载: 导出CSV 

| 显示表格

通过CAF的方法进行频差估计的定位结果如表3所示。

表 3  常规方法定位结果

目标真实位置(m)	目标测量位置(m)	相对误差(m)
(6303534, –92215, 968286)	(6303808, –91936, 970069)	1825
(6348796, –27959, 610439)	(6348976, –28173, 612353)	1935
(6339522, –16218, 700588)	(6339615, –16496, 701470)	930
(6374947, 30978, 199294)	(6374905, 30788, 198028)	1280

下载: 导出CSV 

| 显示表格

6.   结论

本文从当前定位模式出发，分析了传统定位方法的局限性和双星干涉测量时差频差的稳定性，在此基础上提出了基于合成孔径的双星时差频差定位算法，并分析了算法的性能。STK仿真数据的处理结果证实了本文提出的新方法在大范围内可以实现1 km以内的定位精度。