CN104076353B - 一种面目标回波波束中心速度测量方法 - Google Patents

Abstract

本发明公开了一种面目标回波波束中心速度测量方法，首先对回波进行加速度补偿，然后判断速度滤波值是处于近区测量模式还是远区测量模式，近区测量模式经过FFT运算得到回波的频谱之后，利用Chirp‑Z方法对多普勒速度附近的谱线进行细化操作，得到回波频谱，远区测量模式经过FFT运算得到回波频谱，最后将回波频谱进行平滑处理和包络截取以得到波束中心对应的频率值，从而计算得到速度值，最后采用Kalman滤波的方法对速度进行跟踪，并对下一时刻的速度，加速度等信息进行预测。该方法兼顾了大动态范围和高精度的要求，提高了测量精度。

Description

一种面目标回波波束中心速度测量方法

技术领域

本发明涉及一种回波波束中心速度测量方法，尤其涉及一种连续波多普勒雷达高精度面目标回波波束中心速度测量方法，属于空间微波遥感技术领域。

背景技术

连续波测速雷达是月球深空探测着陆工程中的关键单机，它为卫星制导、导航和控制系统提供精确的卫星速度信息，从而保证着陆的安全。连续波测速雷达的测量精度较高，采用单一利用高精度面目标回波多普勒波束中心求取方法，对面目标多普勒频谱计算出波束中心的多普勒频率。另外，连续波测速雷达在导弹制导以及舰载飞机、无人机起降方面都有广阔的应用前景。

天线发射的电磁波是具有一定宽度的波束，波束内各个分量方向上速度的径向分量不同，也就对应不同的多普勒频移，因此天线收到的回波信号是由多个具有不同频率的分量叠加而成，由于着陆器存在加速度，会导致多普勒带宽展宽，同时由于散射等因素的影响，会使回波信号中峰值分量的位置发生偏移。

为了提高多普勒信号的测量精度，必须消除着陆器运动加速度引起的多普勒频谱展宽的影响，所以需要进行加速度的估计和补偿，目前现有文献都是对点目标机动运动条件下的运动补偿，可采用最大似然方法、时频分析方法、基于加速度模板的解线性调频方法估计目标加速度，这些方法所需要的运算量都过大,不利于实时处理；二阶多项相位变换(DPT)是一种运算量较小的估计加速度的方法，但对于面目标的加速度估计，它的精度将大大降低。另外，对于面目标多普勒频率的检测，现有的文献中采用前沿检测、能量重心法或几何中心法，但是这些方法对于求取天线波束中心对应的多普勒频率，估计误差较差。

发明内容

本发明解决的技术问题是：克服现有技术的不足，提供一种面目标回波波束中心速度测量方法，兼顾了大动态范围和高精度的要求，提高了测量精度。

本发明的技术解决方案是：一种面目标回波波束中心速度测量方法，步骤如下：

(1)天线接收到的射频回波信号经过微波通道变频解调后得到基带模拟信号，该基带模拟信号经过AD采集模块转换为雷达数字回波信号s(n)；

(2)利用卡尔曼滤波对前一个测量周期面目标回波波束中心的速度进行跟踪滤波，并对本周期面目标回波波束中心的加速度进行估计，得到前一个周期面目标回波波束中心的速度滤波值v和本周期面目标回波波束中心的加速度估计值a；

(3)利用步骤(2)得到的本周期面目标回波波束中心的加速度估计值a，对步骤(1)得到的雷达数字回波信号s(n)进行加速度补偿；

(4)判断步骤(2)得到的前一个周期面目标回波波束中心的速度滤波值v是处于近区测量模式还是远区测量模式，如果是处于近区测量模式，执行步骤(5)；否则执行步骤(6)；

(5)对加速度补偿后的雷达数字回波信号进行FFT谱分析得到雷达数字回波多普勒频谱，根据雷达数字回波多普勒频谱对加速度补偿后的雷达数字回波信号进行Chirp-Z变换，得到频谱细化后的雷达数字回波多普勒频谱，执行步骤(7)；

(6)对加速度补偿后的雷达数字回波信号进行FFT谱分析，得到雷达数字回波多普勒频谱，执行步骤(7)；

(7)对雷达数字回波多普勒频谱进行平滑处理，对平滑处理后的雷达数字回波多普勒频谱进行包络截取，在此基础上计算得到面目标回波波束中心速度值；

(8)重复执行步骤(1)—(7)，得到面目标回波波束中心的实时速度值。

所述步骤(2)的实现方式为：

设速度滤波值的三阶矩阵为x_filter[]，速度估计值的三阶矩阵为x_estimate[]，速度测量值为v_m，则x_filter＝x_estimate+K*z；

其中

z＝v_m-H*x_estimate

$K=\left(\frac{1}{H*P*H^T+R}\right)*P*H^T$

P＝φ*P_上一周期*φ^T+Q

x_estimate＝φ*x_filter

P_上一周期为上一个周期更新后的误差协方差矩阵；

每个周期计算完成后，按照下式更新误差协方差矩阵：

P＝(I-K*H)*P，其中I为单位矩阵；

在上述基础上，得到速度滤波值v＝x_filter[0]；加速度估计值a＝x_estimate[1]；

其中，H为观测矩阵， $H=\left[\begin{array}{c}1\\ 0\\ 0\end{array}\right];$ φ为状态转移矩阵， $\mathrm{\φ}=\left[\begin{array}{ccc}1& T& \frac{T^2}{2}\\ 0& 1& T\\ 0& 0& 1\end{array}\right];$ Q为动态噪声协方差矩阵， $Q=\left[\begin{array}{ccc}\frac{T^4}{4}& \frac{T^3}{2}& \frac{T^2}{2}\\ \frac{T^3}{2}& T^2& T\\ \frac{T^2}{2}& T& 1\end{array}\right]{\mathrm{\σ}}_w^2;$ P为误差协方差矩阵， $P=\left[\begin{array}{ccc}1& \frac{1}{T}& \frac{1}{T^2}\\ \frac{1}{T}& \frac{2}{T^2}& \frac{3}{T^3}\\ \frac{2}{T^2}& \frac{3}{T^3}& \frac{6}{T^4}\end{array}\right]q_v^2;$ R为观测噪声协方差矩阵初始值，T为周期间隔，为目标速度动态噪声方差，为观测速度方差。

所述步骤(3)的实现方式为：

加速度补偿后的雷达数字回波信号为

$s^{\text{'}}\left(n\right)=s\left(n\right)\×\exp (-j\frac{4\mathrm{\πa}}{\mathrm{\λ}}{n}^2)=A_s\exp \left(j\frac{4\mathrm{\πv}}{\mathrm{\λ}}n\right)$

其中AS为数字回波信号幅度，λ为雷达发射连续波波长，n为离散数字采样序列号， $s\left(n\right)=A_s\exp (j\frac{4\mathrm{\πv}}{\mathrm{\λ}}n+j\frac{4\mathrm{\πa}}{\mathrm{\λ}}{n}^2).$

所述步骤(4)的实现方式为：

当前一个周期波束中心的速度滤波值v∈[80m/s～500m/s]时，为远区测量模式，当v∈[-36m/s～80m/s]时，为近区测量模式。

所述步骤(5)中根据雷达数字回波多普勒频谱对加速度补偿后的雷达数字回波信号进行Chirp-Z变换，得到频谱细化后的雷达数字回波多普勒频谱的实现方式为：

设加速度补偿后的雷达数字回波信号中的时域序列x(n)长度为N，则利用Chirp-Z变换算法得到Z平面上M点频谱采样值X(z_k)的方式为： $X\left(z_k\right)=W^{k^2/2}\underset{n=0}{\overset{N-1}{\mathrm{\Σ}}}y\left(n\right)h(k-n),0\≤k\≤M-1;$

其中

$y\left(n\right)=x\left(n\right)A^{-n}{W}^{n^2/2}$

$h(k-n)=W^{-{(k-n)}^2/2}$

A和W为极坐标下的参数，A₀和W₀为1；θ₀为频谱序列x(n)起始取样点的相角，通过从FFT谱分析得到的信息提取得到；为取样间隔，M为取样点数。

所述步骤(7)的实现方式为：

对雷达数字回波多普勒频谱进行平滑处理的方式为：设频谱序列f(n')共N₁个点，n'为频谱序列号，平滑点数为2L+1，那么平滑处理后序列F(n')表示为：

$F\left(n^{\text{'}}\right)=\left\{\begin{array}{cc}f\left(n^{\text{'}}\right)& n^{\text{'}=\mathrm{0,1},...,L-1}\\ \frac{1}{2L+1}\underset{m=n^{\text{'}}-L}{\overset{m=n^{=\text{'}}+L}{\mathrm{\Σ}}}& n^{\text{'}}=L,...,N_1-L-1\\ f\left(n^{\text{'}}\right)& n^{\text{'}}=N_1-L,...,N_1-1\end{array}\right\}$

对平滑处理后的雷达数字回波多普勒频谱进行包络截取的方式为：在平滑处理后的雷达数字回波多普勒频谱上选择该段频谱幅度的峰值以及该幅度峰值对应的频率坐标，并根据信噪比和天线波束的主瓣宽度确定幅度门限，在频谱上从幅度峰值对应的频率坐标处向两边进行搜索，当搜索到的频谱幅度小于设定的门限时，得到两个频率坐标值，分别为最小多普勒频率f_min和最大多普勒频率f_max；

计算得到波束中心频率估计值的方式为：其中f_beam为多普勒中心频率，θ表示天线波束主瓣宽度对应角度；从而得到波束中心速度值为其中λ为雷达发射连续波波长

本发明与现有技术相比具有如下有益效果：

(1)本发明为了提高速度测量精度，采用加速度补偿技术，补偿掉回波信号中加速度引起的多普勒分量，消除了回波多普勒频谱因为加速度的存在而产生展宽的情况；使用卡尔曼滤波同时提供了速度测量滤波值和加速度的估计值；

(2)本发明为了兼顾大动态范围和高精度的要求，将速度测量区域划分为远区处理模式和近区处理模式，利用不同的采样率覆盖不同的多普勒频率范围，并且在精度要求较高的近区处理模式中，采用Chirp-Z变换的方法，提高频谱分辨率；

(3)求取多普勒中心时，先进行平滑处理，其优点是能够提高多普勒中心估计的精度，同时能够克服地形变化引起的多普勒频谱不连续的情况；在多普勒功率谱上进行自顶向下的包络截取方法，减小了天线旁瓣引入强散射目标的可能。

附图说明

图1为面目标回波波束中心速度测量流程图；

图2为Chirp-z变换在z平面上的螺旋采样示意图；

图3为Chirp-z变换算法原理框图；

图4为面目标回波波束中心频率求取示意图。

具体实施方式

如图1所示，本发明提出了一种高精度面目标回波波束中心速度测量方法，步骤如下：

(1)天线接收到的射频回波信号经过微波通道变频解调后得到基带模拟信号，经过AD采集模块转换为雷达数字回波信号s(n)；

(2)利用卡尔曼滤波对前一个测量周期波束中心的速度进行跟踪滤波，并对本周期波束中心的加速度进行估计，得到前一个周期波束中心的速度滤波值v和本周期波束中心的加速度估计值a；

(3)利用步骤(2)得到的本周期波束中心的加速度估计值a，对步骤(1)得到的雷达数字回波信号s(n)进行加速度补偿；

当不考虑加速度时，连续波多普勒回波信号为：

s(n)＝A_Sexp(j2πf_dn)

$f_d=\frac{2v}{\mathrm{\λ}}$

当考虑加速度时，连续波多普勒回波信号为：

s(n)＝A_Sexp(j2πf_dn)

$f_d=\frac{2(v+\mathrm{an})}{\mathrm{\λ}}$

$s\left(n\right)=A_s\exp (j\frac{4\mathrm{\πv}}{\mathrm{\λ}}n+j\frac{4\mathrm{\πa}}{\mathrm{\λ}}{n}^2)$

所以，考虑均加速度运动模型的基础上进行加速度补偿，需要给回波信号乘以因子 $\exp (-j\frac{4\mathrm{\πa}}{\mathrm{\λ}}{n}^2);$

加速度补偿后的雷达数字回波信号

$s^{\text{'}}\left(n\right)=s\left(n\right)\×\exp (-j\frac{4\mathrm{\πa}}{\mathrm{\λ}}{n}^2)=A_s\exp \left(j\frac{4\mathrm{\πv}}{\mathrm{\λ}}n\right)$

其中A_S为数字回波信号幅度，λ为雷达发射连续波波长，n为离散数字采样序列号。

(4)根据步骤(2)得到的前一个周期波束中心的速度滤波值v判断该速度滤波值是处于近区测量模式还是远区测量模式，如果是处于近区测量模式，执行步骤(5)；否则执行步骤(6)；

根据速度测量范围以及所要求的精度，将速度测量范围划分为远区测量模式和近区测量模式，远区测量模式速度范围80～500m/s，测量要求精度1σ为0.33％V，近区测量模式速度范围为-36～80m/s，测量精度最高1σ为0.05m/s；

根据上一个周期速度测量结果对下一个周期测量模式进行切换，当前一个周期波束中心的速度滤波值v大于500m/s或小于-36m/s时，测量模式不变；按照Ka频段8mm波长计算，远区模式多普勒范围为10KHz～125KHz，在满足采样定理的前提下，设置采用率为600KHz，近区模式多普勒范围为-4.5KHz～10KHz，设置采样率为300KHz。

(5)对加速度补偿后的雷达数字回波信号进行FFT谱分析，得到雷达数字回波多普勒频谱，根据雷达数字回波多普勒频谱对加速度补偿后的雷达数字回波信号进行Chirp-Z变换，得到频谱细化后的雷达数字回波多普勒频谱，执行步骤(7)；

设加速度补偿后的雷达数字回波时域序列x(n)的长度为N，傅里叶变换(FFT)写成如下形式：

$X\left(k^{\text{'}}\right)=\underset{n=0}{\overset{N-1}{\mathrm{\Σ}}}x\left(n\right)e^{-j2\mathrm{\πn}{k}^{\text{'}}/N},k^{\text{'}}\∈[0,N-1]$

设时域序列x(n)的长度为N，要获得Z平面上M点频谱采样值X(z_k)，0≤k≤M-1，则其采用的Chirp-Z变换算法具体如下。

设z_k＝AW^-k,0≤k≤M-1

式中A和W均为复数，用极坐标形式表示为

式中A₀和W₀为实数，当k＝0时，有

$z_0=A_0{e}^{j{\mathrm{\θ}}_0}$

由此可见，在上式中，A决定谱分析起始点z₀的位置；W₀的值决定分析路径的盘旋趋势，则表示两个相邻分析点之间的夹角。如果W₀<1，则随着k增大，分析点z_k以为步长向外旋转；W₀>1时则向内旋转。如图2所示，称为螺旋采样。特别的当A₀＝1，并且W₀＝1时则沿着单位圆旋转。

将z_k代入x(n)的Z变换公式得到

$\begin{array}{c}X\left(z_k\right)=\underset{n=0}{\overset{N-1}{\mathrm{\&Sigma;}}}x\left(n\right){\left[{\mathrm{AW}}^{-k}\right]}^{-n}\\ =\underset{n=0}{\overset{N-1}{\mathrm{\&Sigma;}}}x\left(n\right)A^{-n}{W}^{\mathrm{kn}},0\&le;k\&le;M-1\end{array}$

利用如下的关系式

$\mathrm{nk}=\frac{1}{2}[n^2+k^2-{(k-n)}^2]$

可以得到

$\begin{array}{c}X\left(z_k\right)=\underset{n=0}{\overset{N-1}{\mathrm{\&Sigma;}}}x\left(n\right)A^{-n}{W}^{[n^2+k^2-{(k-n)}^2]/2}\\ =W^{k^2/2}\underset{n=0}{\overset{N-1}{\mathrm{\&Sigma;}}}x\left(n\right)A^{-n}{W}^{n^2/2}{W}^{-{(k-n)}^2/2}\end{array}$

令

$y\left(n\right)=x\left(n\right)A^{-n}{W}^{n^2/2}$

$h\left(n\right)=W^{-n^2/2}$

则

$X\left(z_k\right)=W^{k^2/2}\underset{n=0}{\overset{N-1}{\mathrm{\&Sigma;}}}y\left(n\right)h(k-n),0\&le;k\&le;M-1$

上式说明，长度为N的序列x(n)的M点谱分析可以通过预乘得到y(n)，再计算y(n)与h(n)的线性卷积，最后乘以三个步骤得到。这种算法的原理框图如图3所示。图3中，看成一个数字滤波器的单位脉冲响应，其输出V(n)＝y(n)*h(n)。

由于本方案中使用的是A₀＝1，并且W₀＝1的情况，所以对此我们重点讨论分析。θ₀为起始取样点的相角，通过提取FFT谱分析得到的信息得到；为取样间隔，M为取样点数，细化的区域由参数θ₀，和M一起决定,所选弧段的弧度范围为当θ₀＝0，M＝N时，Chirp-Z变换就是DFT。θ₀的取值可正可负，但由于实序列谱的对称性,在实际频谱细化中一般参数θ₀，M选择得要使弧段处在[0,π]范围内。X(z_k)(0≤k≤M-1)各谱线对应细化的区域中的数字角频率如下式所示：

与DFT类似，如果序列X(n)的采样率为f_s，则其对应的模拟频率为

由上式可见，Chirp-Z变换的频率分辨率为越小即M越大时，频率分辨率越高。

(6)对加速度补偿后的雷达数字回波信号进行FFT谱分析，得到雷达数字回波多普勒区域频谱；

设加速度补偿后的雷达数字回波时域序列x(n)的长度为N，傅里叶变换(FFT)写成如下形式：

$X\left(k^{\text{'}}\right)=\underset{n=0}{\overset{N-1}{\mathrm{\&Sigma;}}}x\left(n\right)e^{-j2\mathrm{\&pi;n}{k}^{\text{'}}/N},k^{\text{'}}\&Element;[0,N-1]$

(7)对雷达数字回波多普勒频谱进行平滑处理，对平滑处理后的雷达数字回波多普勒频谱进行包络截取，在此基础上计算得到波束中心速度值；

以波束R1为例，图4给出了R1波束照射地面目标及其回波频率分布的示意图。其中，φ₁为速度方向与水平方向的夹角，φ₂为速度方向与R1波束中心的夹角。波束方位向角和俯仰角均为θ。f_min、f_max分别对应回波信号的最小多普勒频率和最大多普勒频率，f_mid为多普勒频率的中值，f_beam对应波束中心点处的多普勒频率值。

回波信号的最小多普勒频率f_min为：

$f_{\min }=\frac{2V}{\mathrm{\&lambda;}}\cos ({\mathrm{\&phi;}}_2+\mathrm{\&theta;})$

回波信号的最大多普勒频率f_max为：

$f_{\max }=\frac{2V}{\mathrm{\&lambda;}}\cos ({\mathrm{\&phi;}}_2-\mathrm{\&theta;})$

根据上式可得多普勒频率的中值f_mid为：

$f_{\mathrm{mid}}=\frac{V}{\mathrm{\&lambda;}}\cos ({\mathrm{\&phi;}}_2+\mathrm{\&theta;})+\frac{V}{\mathrm{\&lambda;}}\cos ({\mathrm{\&phi;}}_2-\mathrm{\&theta;})$

波束中心点处的多普勒频率值f_beam为：

$f_{\mathrm{beam}}=\frac{2V}{\mathrm{\&lambda;}}\cos {\mathrm{\&phi;}}_2$

由此可得：

$\frac{f_{\mathrm{mid}}}{f_{\mathrm{beam}}}=\frac{\cos ({\mathrm{\&phi;}}_2+\mathrm{\&theta;})+\cos ({\mathrm{\&phi;}}_2-\mathrm{\&theta;})}{2\cos {\mathrm{\&phi;}}_2}=\cos \mathrm{\&theta;}$

整理可得：

$f_{\mathrm{beam}}=\frac{f_{\mathrm{mid}}}{\cos \mathrm{\&theta;}}$

因此，我们可以采取平滑处理和包络截取的方法先得到f_mid的值，然后根据上式计算出波束中心对应的频率值。

平滑处理在多普勒功率谱上进行，一方面为了提高多普勒中心估计精度，另一方面避免了地形起伏引起的多普勒频谱不连续的情况。对雷达数字回波多普勒频谱进行平滑处理的方式为：设频谱序列f(n')共N₁个点，平滑点数为2L+1，n'为频谱序列号，那么平滑处理后序列F(n')表示为：

$F\left(n^{\text{'}}\right)=\left\{\begin{array}{cc}f\left(n^{\text{'}}\right)& n^{\text{'}=\mathrm{0,1},...,L-1}\\ \frac{1}{2L+1}\underset{m=n^{\text{'}}-L}{\overset{m=n^{=\text{'}}+L}{\mathrm{\&Sigma;}}}f\left(m\right)& n^{\text{'}}=L,...,N_1-L-1\\ f\left(n^{\text{'}}\right)& n^{\text{'}}=N_1-L,...,N_1-1\end{array}\right\}$

包络截取是从多普勒功率谱的峰值，自顶向下截取天线波束宽度对应门限，这种方法，减小了天线旁瓣引入强散射目标的情况。包络截取的方式为：在平滑处理后的雷达数字回波多普勒区域频谱上选择该段频谱幅度的峰值以及该幅度峰值对应的频率坐标，并根据信道比和天线波束的主干宽度确定幅度门限，在频谱上从幅度峰值对应的频率坐标处向两边进行搜索，当搜索到的频谱幅度小于设定的门限时，得到两个频率坐标值，分别为最小多普勒频率f_min和最大多普勒频率f_max；

计算得到多普勒中心频率估计值的方式为：

计算得到波束中心速度值为其中λ为雷达发射连续波波长。

(8)采用卡尔曼滤波方法对前一个周期波束中心的速度进行跟踪，并对本周期波束中心的加速度信息进行估计，然后重复执行步骤(1)—(7)，实时地得到面目标回波波束中心的速度值。

实施例：

定义AD采集的雷达数字回波信号为s(n)＝s_i(n)+j*s_q(n)，采样率对应两种处理模式分别为600KHz和300KHz，两种采样率下信号积累时间分别为25ms和50ms，数字信号采样点数为15000点；

本发明的处理步骤如下：

(1)天线接收到的射频回波信号经过微波通道变频解调后得到基带模拟信号，经过AD采集模块转换为雷达数字回波信号s(n)＝s_i(n)+j*s_q(n)；

(2)利用卡尔曼滤波对前一个测量周期的速度结果进行跟踪滤波，并对加速度进行估计；

观测矩阵 $H=\left[\begin{array}{c}1\\ 0\\ 0\end{array}\right]$

状态转移矩阵 $\mathrm{\&phi;}=\left[\begin{array}{ccc}1& T& \frac{T^2}{2}\\ 0& 1& T\\ 0& 0& 1\end{array}\right]$

动态噪声协方差矩阵 $Q=\left[\begin{array}{ccc}\frac{T^4}{4}& \frac{T^3}{2}& \frac{T^2}{2}\\ \frac{T^3}{2}& T^2& T\\ \frac{T^2}{2}& T& 1\end{array}\right]{\mathrm{\&sigma;}}_w^2$

误差协方差矩阵 $P=\left[\begin{array}{ccc}1& \frac{1}{T}& \frac{1}{T^2}\\ \frac{1}{T}& \frac{2}{T^2}& \frac{3}{T^3}\\ \frac{2}{T^2}& \frac{3}{T^3}& \frac{6}{T^4}\end{array}\right]q_v^2$

观测噪声协方差矩阵初始值

其中T为周期间隔取128ms，为目标速度动态噪声方差取144，为观测速度方差取9。设速度滤波值三阶矩阵x_filter[]，速度估计值三阶矩阵x_estimate[]，速度测量值为v_m；

第一个周期：x_filter[0]＝v_m；

第二个周期：x_filter[1]＝(v_m-x_filter[0])/T；

x_filter[0]＝v_m；

第三个周期：x_filter[2]＝((v_m-x_filter[0])/T-x_filter[1])/T

x_filter[1]＝(v_m-x_filter[0])/T

x_filter[0]＝v_m；

x_estimate[]＝φ·x_filter[]；

第四个周期以后：

x_filter＝x_estimate+K*z

其中

z＝v_m-H*x_estimate

$K=\left(\frac{1}{H*P*H^T+R}\right)*P*H^T$ 其中P＝φ*P_上一周期*φ^T+Q

P_上一周期为上一个周期更新后的误差协方差矩阵

x_estimate＝φ*x_filter

每个周期最后更新误差协方差矩阵

P＝(I-K*H)*P其中I为单位矩阵；

v＝x_filter[0]为速度滤波值

a＝x_estimate[1]为估计得到的加速度信息。

(3)根据步骤(2)得到的加速度结果，对步骤(1)得到的雷达基带数字回波信号进行加速度补偿；

加速度补偿是给回波信号乘以因子步骤(1)中得到采集的复数数字信号为s(n)＝s_i(n)+j*s_q(n)，即与复数相乘得到补偿后的复数时域序列s'(n)；

(4)根据步骤(2)得到的前一个周期波束中心的速度滤波值v判断该速度滤波值是处于近区测量模式还是远区测量模式，如果是处于近区测量模式，执行步骤(5)；否则执行步骤(6)；

(5)对加速度补偿后的雷达数字回波信号进行FFT谱分析，得到雷达数字回波多普勒频谱X(k')，根据雷达数字回波多普勒频谱对加速度补偿后的雷达数字回波信号进行Chirp-Z变换，得到频谱细化后的雷达数字回波多普勒频谱，然后执行步骤(7)；

对加速度补偿后的时域回波信号，进行FFT谱分析；

补偿后的复数时域序列s'(n)的长度为15000，进行N＝16384点傅里叶变换(FFT)写成如下形式：

$X\left(k^{\text{'}}\right)=\underset{n=0}{\overset{N-1}{\mathrm{\&Sigma;}}}x\left(n\right)e^{-j2\mathrm{\&pi;n}{k}^{\text{'}}/N},k^{\text{'}}\&Element;[0,N-1];$

近区测量模式，测量精度要求高，对回波多普勒区域频谱采用Chirp_Z变换进行频谱细化；利用Chirp-Z变换将频谱峰值坐标为中心的512点频谱作为分析频谱，细化倍数为8倍，得到4096点频谱细化结果作为波束中心分析频谱X'(k')的实现方式为：

补偿后的复数时域序列s'(n)的长度为N＝15000，利用Chirp-Z变换算法得到Z平面上M＝4096点频谱采样值X(z_k)的方式为： $X\left(z_k\right)=W^{k^2/2}\underset{n=0}{\overset{N-1}{\mathrm{\&Sigma;}}}y\left(n\right)h(k-n),0\&le;k\&le;M-1,$ 0≤k≤M-1；

其中

$y\left(n\right)=s^{\text{'}}\left(n\right)A^{-n}{W}^{n^2/2}$

$h(k-n)=W^{-{(k-n)}^2/2}$

A₀和W₀为1；θ₀为频谱序列s'(n)起始取样点的相角，为取样间隔，M为取样点数；X'(k')＝X(z_k)。

(6)对加速度补偿后的雷达数字回波信号进行FFT谱分析，得到雷达数字回波多普勒频谱X(k')，将X(k')频谱中峰值坐标为中心的1024点频谱取出为波束中心分析频谱X'(k')；

(7)对回波频谱X'(k')进行平滑处理和包络截取，从而进行多普勒中心估计值；

对X'(k')的波束中心分析频谱进行平滑处理：

$X^{\text{'}\text{'}}\left(k^{\text{'}}\right)=\left\{\begin{array}{cc}{X}^{\text{'}}\left(k^{\text{'}}\right)& k^{\text{'}=\mathrm{0,1},...,L-1}\\ \frac{1}{2L+1}\underset{m=k^{\text{'}}-L}{\overset{m=k^{\text{'}}+L}{\mathrm{\&Sigma;}}}{X}^{\text{'}}\left(m\right)& k^{\text{'}}=L,...,N_1-L-1\\ X^{\text{'}}\left(k^{\text{'}}\right)& k^{\text{'}}=N_1-L,...,N_1-1\end{array}\right\}$

L取45；

先得到该段频谱的最大幅度A和坐标index，设幅度门限为A/36，从index位置向两边进行搜索，当搜索到的频谱幅值小于A/36的门限时，得到最小多普勒频率f_min和最大多普勒频率f_max，设θ＝6°，根据公式得到多普勒波束中心频率f_beam，再根据公式得到波束中心速度值。

(8)采用Kalman滤波方法对速度进行跟踪，并对下一周期的加速度信息进行估计，然后重复执行步骤(1)—(7)，实时地得到面目标回波波束中心的速度值。

Claims (5)

1.一种面目标回波波束中心速度测量方法，其特征在于步骤如下：

(1)天线接收到的射频回波信号经过微波通道变频解调后得到基带模拟信号，该基带模拟信号经过AD采集模块转换为雷达数字回波信号s(n)；

(2)利用卡尔曼滤波对前一个测量周期面目标回波波束中心的速度进行跟踪滤波，并对本周期面目标回波波束中心的加速度进行估计，得到前一个周期面目标回波波束中心的速度滤波值v和本周期面目标回波波束中心的加速度估计值a；

设速度滤波值的三阶矩阵为x_filter[]，速度估计值的三阶矩阵为x_estimate[]，速度测量值为v_m，则x_filter＝x_estimate+K*z；

其中

z＝v_m-H*x_estimate

P＝φ*P_上一周期*φ^T+Q

x_estimate＝φ*x_filter

P_上一周期为上一个周期更新后的误差协方差矩阵；

每个周期计算完成后，按照下式更新误差协方差矩阵：

P＝(I-K*H)*P，其中I为单位矩阵；

在上述基础上，得到速度滤波值v＝x_filter[0]；加速度估计值a＝x_estimate[1]；

其中，H为观测矩阵，φ为状态转移矩阵，Q为 动态噪声协方差矩阵，P为误差协方差矩阵，R为观测噪声协方差矩阵初始值，T为周期间隔，为目标速度动态噪声方差，为观测速度方差；

(3)利用步骤(2)得到的本周期面目标回波波束中心的加速度估计值a，对步骤(1)得到的雷达数字回波信号s(n)进行加速度补偿；

(4)判断步骤(2)得到的前一个周期面目标回波波束中心的速度滤波值v是处于近区测量模式还是远区测量模式，如果是处于近区测量模式，执行步骤(5)；否则执行步骤(6)；

(5)对加速度补偿后的雷达数字回波信号进行FFT谱分析得到雷达数字回波多普勒频谱，根据雷达数字回波多普勒频谱对加速度补偿后的雷达数字回波信号进行Chirp-Z变换，得到频谱细化后的雷达数字回波多普勒频谱，执行步骤(7)；

(6)对加速度补偿后的雷达数字回波信号进行FFT谱分析，得到雷达数字回波多普勒频谱，执行步骤(7)；

(7)对雷达数字回波多普勒频谱进行平滑处理，对平滑处理后的雷达数字回波多普勒频谱进行包络截取，在此基础上计算得到面目标回波波束中心速度值；

(8)重复执行步骤(1)—(7)，得到面目标回波波束中心的实时速度值。

2.根据权利要求1所述的一种面目标回波波束中心速度测量方法，其特征在于：所述步骤(3)的实现方式为：

加速度补偿后的雷达数字回波信号为

其中A_S为数字回波信号幅度，λ为雷达发射连续波波长，n为离散数字采样序列号，

3.根据权利要求1所述的一种面目标回波波束中心速度测量方法，其特征在于：所述步骤(4)的实现方式为：

当前一个周期波束中心的速度滤波值v∈[80m/s～500m/s]时，为远区测量模式，当v∈[-36m/s～80m/s]时，为近区测量模式。

4.根据权利要求1所述的一种面目标回波波束中心速度测量方法，其特征在于：所述步骤(5)中根据雷达数字回波多普勒频谱对加速度补偿后的雷达数字回波信号进行Chirp-Z变换，得到频谱细化后的雷达数字回波多普勒频谱的实现方式为：

设加速度补偿后的雷达数字回波信号中的时域序列x(n)长度为N，则利用Chirp-Z变换算法得到Z平面上M点频谱采样值X(z_k)的方式为：

其中

A和W为极坐标下的参数，A₀和W₀为1；θ₀为频谱序列x(n)起始取样点的相角，通过从FFT谱分析得到的信息提取得到；为取样间隔，M为取样点数。

5.根据权利要求1所述的一种面目标回波波束中心速度测量方法，其特征在于：所述步骤(7)的实现方式为：

对雷达数字回波多普勒频谱进行平滑处理的方式为：设频谱序列f(n')共N₁个点，n'为频谱序列号，平滑点数为2L+1，那么平滑处理后序列F(n')表示为：

对平滑处理后的雷达数字回波多普勒频谱进行包络截取的方式为：在平滑处理后的雷达数字回波多普勒频谱上选择该段频谱幅度的峰值以及该幅度峰值对应的频率坐标，并根据信噪比和天线波束的主瓣宽度确定幅度门限，在频谱上从幅度峰值对应的频率坐标处向两边进行搜索，当搜索到的频谱幅度小于设定的门限时，得到两个频率坐标值，分别为最小多普勒频率f_min和最大多普勒频率f_max；

计算得到波束中心频率估计值的方式为：其中f_beam为多普勒中心频率，θ表示天线波束主瓣宽度对应角度；从而得到波束中心速度值为其中λ为雷达发射连续波波长。