CN101609147B

CN101609147B - 一种基于时频域预滤波的中心频率估计方法

Info

Publication number: CN101609147B
Application number: CN2008101151802A
Authority: CN
Inventors: 胡东辉; 仇晓兰; 丁赤飚
Original assignee: Institute of Electronics of CAS
Current assignee: Institute of Electronics of CAS
Priority date: 2008-06-18
Filing date: 2008-06-18
Publication date: 2012-02-08
Anticipated expiration: 2028-06-18
Also published as: CN101609147A

Abstract

本发明是一种基于时频域预滤波的中心频率估计方法。所述方法首先采用时频分析方法将合成孔径雷达的方位信号变换到时频域，然后根据中心频率和多普勒调频率初值对时频域信号进行预滤波以去除不完整孔径信号，接着将预滤波后的信号变回时域，并采用基于频谱对称性的中心频率估计方法进行估计。本发明有效的解决了中心频率估计精度对于场景的依赖性，大大提高了不均匀场景尤其是两方位边缘具有强对比度的场景下的估计精度。

Description

一种基于时频域预滤波的中心频率估计方法

技术领域

本发明涉及合成孔径雷达信号处理领域，特别涉及一种基于时频域预滤波的中心频率估计方法。

背景技术

中心频率是合成孔径雷达(SAR)信号处理中一个至关重要的参数。首先它是SAR成像处理的关键参数，中心频率误差的存在，会导致图像信噪比的降低、方位模糊度的升高以及分辨率的下降；其次它是图像定位的重要参数，中心频率误差将直接导致定位精度的下降；第三正确的多普勒参数估计值也能为姿轨控制误差提供有效的检验手段。在星载SAR的信号处理中，通过卫星平台姿轨信息并利用星地几何关系，可以计算得到多普勒中心频率和调频率的初值。然而由于姿轨控制存在误差，该计算值可能存在较大的偏差。为此，通回波数据估计多普勒参数往往是必不可少的。

在现有技术中，多普勒中心频率估计的算法可大致分为两类，一类是利用特显点，从特显点的距离徙动曲线以及相位历程获取多普勒参数的信息(参见文献[1]：王宏艳，吴彦鸿，贾鑫，“星载SAR多普勒参数估计的改进方法”，《现代雷达》，Vol.28，No.12，pp52-54，2006年12月)；另一类则是基于天线方向图(从而回波数据的方位谱)的对称性进行估计。前一种方法需要场景中有比较明显的特显点，故对场景有一定的要求，并且特显点的自动选取也是一个非常困难的问题。后一种方法较前一种方法简单易行，但不管是利用频谱幅度(如能量均衡法，参见文献[2]：Curlander，J.C.，1991.Synthetic Aperture Radar Systems and SignalProcessing，New York，John Wiley & Sons，Inc)实现，还是利用相关函数的相位信息(如CDE算法(参见文献[3]：Madsen，S.N.，1989.Estimating the Doppler centroid of SAR Data.IEEETransaction on Aerospace and Electronic Systems，25(2)134-145)实现，都会受到位于场景边缘处目标的不完整孔径信号的影响。尤其当场景两端对比度较大(如海陆交界场景)时，估计值会存在很大的偏差。随着星载SAR分辨率的进一步提高和偏航控制等技术的引入，对多普勒中心频率估计的精确性和稳健性要求都将随之提高，因此需要一种更高精度更具稳健性的中心频率估计方法。

发明内容

鉴于上述问题的存在，本发明的目的是提供一种基于时频域预滤波的中心频率估计方法，以克服现有算法受不完整孔径信号影响大的缺点，减小估计精度对于场景内容的依赖性，从而提高估计的精度和算法的稳健性。所述方法的步骤包括：

1)通过计算或估计的方法得到多普勒中心频率和调频率的初值；

2)对估计区间的每个方位信号进行时频分析；

3)确定时频域中完整信号所在区域，并将区域外信号置零，完成时频域预滤波；

4)对每个子段进行逆傅立叶变换，得到滤波后的时域信号，并进行中心频率估计，得到基带中心频率估计值。

进一步，所述对估计区间的每个方位信号进行时频分析具体为：先根据

及

尽量接近

的准则选择短时傅立叶变换的子段长度

，然后对每个子段进行傅立叶变换从而得到时频域信号。

进一步，所述确定时频域中完整信号所在区域是根据时频域关系

0 \leq \frac{f - f_{DC}}{f_{r}} + t_{a} < T_{a},

&ForAll; f &Element; [f_{DC} - f_{P} / 2, f_{DC} - f_{P} / 2]

得到的。

进一步，所述中心频率估计是采用基于方位谱对称特性的中心频率估计方法进行的。

进一步，首先在所述时频域预滤波去除不完整孔径信号。

进一步，步骤4)中所述的对时频域滤波后的中心频率估计方法采用基于方位谱对称特性的估计方法，后续采用CDE方法、或能量均衡法、或频谱拟合法、或Madsen方法。

附图说明

图1为方位回波信号时频域分布示意图；

图2为子段长度过大时的时频域示意图；

图3为子段长度过小时的时频域示意图；

图4为回波仿真场景；

图5(a)至图5(f)为仿真数据方位信号经时频域滤波前后，在时域、频域及时频域的幅度图；

图6为仿真数据中心频率估计误差；

图7为估计精度与予段长度的关系；

图8为RadarSat-1实测数据场景内容及方位区间划分；

图9(a)至图9(h)为RadarSat-1实测数据8个方位区间在采用传统CDE方法、时频域滤波的CDE方法以及两次迭代的时频域滤波CDE方法时的中心频率估计值。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚明白，以下结合具体实施例，并参照附图，对本发明进一步详细说明。

本发明提供了一种基于时频域预滤波的中心频率估计方法，所述方法的具体步骤包括：

1)根据卫星平台姿轨信息或星历数据计算多普勒中心频率和调频率的初值。

中心频率：

f_{DC} = - \frac{2 V_{st} \cdot R}{λR}

调频率：

f_{r} = - \frac{2}{λ} [\frac{V_{st} \cdot V_{st}}{R} + \frac{A_{st} \cdot R}{R} - \frac{{(V_{st} \cdot R)}^{3}}{R^{3}}]

其中

粗黑体表示矢量，细斜体表示标量。

为多普勒中心频率，

为多普勒调频率，λ为波长。表示

卫星的位置矢量，

为目标的位置矢量，

为卫星速度矢量，

为目标速度矢量。

为卫星加速度矢量，

为目标加速度矢量。上述矢量均位于同一坐标系(如惯性坐标系)下。为目标至卫星的位置矢量，为目标与卫星的相对速度矢量，为目标与卫星的相对加速度矢量。下标

表示卫星，下标

表示目标。

2)将回波数据(距离压缩前)沿距离向划分成不同的区间，方位向长度大于一个(最好大于2个)合成孔径长度，取其中部分或全部区间作为估计区间。

3)对每个估计区间内的每条方位信号进行短时傅立叶分析(STFT)，即将方位信号分成长度为的时间子段，并对每一个子段进行傅立叶变换，从而构成了一个轴为子段序号(时间)，另一个轴为频率的二维时频域。子段长度

应这样选择：首先保证子段信号时间带宽积

，其中

为脉冲重复频率；然后使

尽量接近其中为一个合成孔径时间内的脉冲个数。

4)根据多普勒中心频率和调频率的初值，计算时频域中完整信号所在区间，该区间的时间长度由满足下式的时间

的取值区间决定。

0 \leq \frac{f - f_{DC}}{f_{r}} + t_{a} < T_{a},

&ForAll; f &Element; [f_{DC} - f_{P} / 2, f_{DC} - f_{P} / 2]

其中

为方位信号时间长度，为方位频率，

为多普勒中心频率的初值，

为多普勒调频率的初值，

为脉冲重复频率。

5)将上述区间外的时频域区间置零，然后将每个子段信号变回到时域，得到时频域滤波后的方位信号。

6)对滤波后的方位信号进行CDE估计(参见文献[3])，得到基带中心频率的估计值

估计方法为

其中

通常取1，

为方位信号自相关函数。

7)如估计精度仍然较差，则将上述估计值作为初值，返回步骤3)再次进行估计。

对各估计区间的估计值进行多项式拟合，得到回波数据每个距离门的中心频率估计值。

图1形象的示意了不完整孔径信号对于方位频谱的影响。由于SAR方位回波信号可以近似为线性调频信号，故其在时频域的分布图如图1所示。位于图中阴影四边形内的信号具有完整的多普勒历程，而该区间外的信号的多普勒历程则是不完整的。如果场景比较均匀，则场景两端不完整信号强度相当，故信号的方位谱能够很好的体现天线方向图的形状，在中心频率处出现峰值，且两边具有很好的对称性(如图中右边虚线所示)，从而CDE等依赖频谱对称性的估计方法能够得到很准确的估计值。如场景两端信号强度悬殊(图中用线条粗细表示信号强弱)，则方位信号功率谱会偏向一端(如图右边实线所示)，从而使CDE估计结果存在很大的偏差。本发明针对这一现象，在时频域滤除不完整信号，仅保留阴影四边形内的完整信号，从而能改善频谱偏移现象，提高中心频率估计精度。

另外需要说明的是，由于完整孔径信号区间需要通过多普勒中心频率和调频率的初值确定，故初值的误差会导致不能完全滤除不完整信号，从而影响估计精度。然而由于时间子段长度一般较小，子段间隔时间较脉冲重复间隔大许多，故时频域的时间分辨率和频率分辨率都比较低，因此对多普勒参数初值的误差并不敏感，如当调频率初值误差为5％，中心频率初值误差10％个

时，估计精度仍能在要求以内(后面的仿真性能分析将证实这一点)。

时频域滤波性能的另一个影响因素是子段长度

，

太大则子段间的时间间隔太长(如图2)，影响时域分辨率，从而导致时频域滤波时不完整信号区域锯齿过大，从而影响滤波精度。同样，当子段长度

太小时(如图3)，频域分辨率将很低，同时由于时间带宽积太小，还会引起频谱展宽，从而影响滤波精度。因此应该选取合适的

值，使得滤波精度最高，从而使中心频率估计精度达到最佳。根据实验，

的取值首先应该保证时间带宽积

，然后应该尽量使的时间分辨率和频率分辨率能够折衷，也即

尽量接近合成孔径长度的开根号值。

1)仿真性能分析

为了观察算法在场景两端对比度极大情况下的结果，仿真采用了如图4的海陆交界场景，并在考虑方位向第一旁瓣(即考虑方位模糊)的前提下仿真了星载X波段的回波数据。仿真数据合成孔径长度

为1509个脉冲，估计时采用方位信号长度为4096个脉冲，距离区间深度为256个距离门，脉冲重复频率为3460Hz。

图5给出了滤波前后方位信号在时域、频域以及时频域的幅度形式。其中图5(a)和图5(b)分别为时频域预滤波前后，方位时域信号的归一化幅度均值；图5(c)和图5(d)分别为时频域预滤波前后，方位信号频谱的归一化幅度均值；图5(e)和图5(f)则分别为时频域预滤波前后，方位信号在时频域的幅度图。由滤波前信号的时域幅度可见场景两端信号能量的悬殊，同时由滤波前的方位谱可见其频谱中心偏离中心频率真值。而经过时频域滤波后，方位频谱中心则非常接近中心频率真值。图6给出了16个距离估计区间的估计误差，时频分析采用的子段长度

＝146。由图，不采用时频域预滤波时，估计值严重偏离真值，估计误差甚至达到近300Hz

。而进行时频域预滤波后，估计误差则在10Hz以内，即使调频率误差为5％时，中心频率估计误差仍在30Hz以内。表1给出了16个估计区间的估计误差统计值，再次证明了算法的优越性。

表1 仿真数据估计结果

图7则给出了估计精度与子段长度的关系，可见当

值过大或过小时，估计精度都较差，而在中间比较长的区间内，估计精度都比较理想，而本发明采用的M值选取方法基本上能够确保估计精度处于相对理想的状态。

2)实际数据验证

为了验证算法应用于实际数据时的稳健性，将该方法应用于RadarSat-1在温哥华地区(如图8)的实测回波数据。方位合成孔径长度为704个脉冲，每个估计区间的方位信号长度为2048，距离向深度为256个距离门。共对方位向8个连续区间进行估计，每个区间沿距离向的子带个数为24。第一次时频域滤波采用常规CDE算法的估计值为初值，并进行2次迭代估计，估计结果用二次多项式进行拟和，并统计估计值与拟和曲线的误差均方值(见表2)。

表2 中心频率估计值与二次拟和曲线的误差均方值(Hz)

图9(a)至图9(h)分别给出了图8中所示8个估计区间在1)不进行时频域预滤波而直接采用CDE方法进行中心频率估计，2)采用1)的估计值作为初值进行时频域预滤波，再采用CDE方法进行中心频率估计，3)采用2)的估计值作为初值，再次进行时频域预滤波，然后采用CDE方法进行中心频率估计，这三种情况下的估计结果和拟和曲线。由图9可见，如不进行时频域滤波直接采用CDE方法进行中心频率估计，估计结果在两端对比强烈的场景(如区间1，2，5)中存在剧烈的波动，而进行时频滤波后，波动明显减小。另外，当场景比较均匀时，滤波后的估计结果则与直接估计的结果非常一致(如区间6)。在8个区间中，除了信噪比较低的区间1，2(大部分是海面)，其他场景区间的误差均方值都在10Hz以内，说明本发明受场景的影响很小，具有很好的稳健性。此外，第二次迭代的估计结果较第一次迭代只是稍有改善，说明该方法对于中心频率初值的误差不太敏感，因此一般情况下迭代一次已经足够。

Claims

1.一种基于时频域预滤波的合成孔径雷达中心频率估计方法，包括如下步骤：

2)对估计区间的每个方位信号进行时频分析；

4)对每个子段进行逆傅立叶变换，得到滤波后的时域信号，并进行中心频率估计，得到基带中心频率估计值；

所述对估计区间的每个方位信号进行时频分析，具体为：先根据D＝|f_rM²/f_p ²|＞10及M尽量接近的准则选择短时傅立叶变换的子段长度M，然后对每个子段进行傅立叶变换从而得到时频域信号；其中：D为子段信号时间带宽积；f_r为多普勒调频率的初值；f_p为脉冲重复频率；N_syn为一个合成孔径时间内的脉冲个数。

2.根据权利要求1所述的方法，其特征在于，所述确定时频域中完整信号所在区域是根据时频域关系

&ForAll; f &Element; [f_{DC} - f_{P} / 2, f_{DC} - f_{P} / 2]

得到的；其中：

f为方位频率；f_DC为多普勒中心频率的初值；T_a为方位信号时间长度；t_a为时间。

3.根据权利要求1所述的方法，其特征在于，所述中心频率估计是采用基于方位谱对称特性的中心频率估计方法进行的。

4.根据权利要求1所述的方法，其特征在于，首先在所述时频域预滤波去除不完整孔径信号。

5.根据权利要求1所述的方法，其特征在于，所述中心频率估计是采用CDE方法、或能量均衡法、或频谱拟合法。