CN104330784B - 基于旋翼物理参数估计实现飞机目标分类的方法 - Google Patents

Abstract

本发明公开了一种基于旋翼物理参数估计实现飞机目标分类的方法，主要解决现有基于微动调制特征的目标分类方法需要训练数据库作支撑条件下的问题。其实现过程是：从飞机雷达回波信号中分离旋翼回波信号，并从旋翼回波信号时频域中估计旋翼转速；在旋翼回波信号的时频域中对旋翼回波信号进行旋翼二维成像；对旋翼二维成像结果预处理；从预处理后的旋翼二维成像结果中估计旋翼长度和螺旋桨叶片数；将估计出的旋翼转速ω、旋翼长度L和螺旋桨叶片数N这些旋翼物理参数与飞机类型标准库中飞机旋翼物理参数比对，判断出飞机目标的类型。本发明无需训练数据库作为支撑条件下即能实现对飞机目标的分类，可提高分类速度。

Description

基于旋翼物理参数估计实现飞机目标分类的方法

技术领域

本发明属于雷达技术领域，涉及飞机目标的旋翼物理参数估计，可用于雷达目标的识别。

背景技术

目前，对于飞机目标分类识别，国内外文献提出了多种窄带雷达体制下基于微动调制特征的目标分类识别方法，例如陈凤等发表的《基于特征谱散布特征的低分辨雷达目标识别方法》，根据螺旋桨飞机、喷气式飞机、直升机三类飞机微动调制引起特征谱的差别来提取分类特征，实现飞机目标的分类。这种基于微动调制特征的飞机目标分类方法的突出缺点是需要大量的仿真或实测训练数据离线学习分类器，在没有训练数据库支撑的条件下这种方法无法用于对飞机目标进行实时分类判决。另外飞机目标微动对雷达回波的调制，不仅取决于微动参数，如旋翼转速，还与飞机目标旋翼的尺寸有关。

发明内容

本发明的目的在于针对上述已有方法的不足，提出一种基于旋翼物理参数估计实现飞机目标分类的方法，以在无需大量仿真或实测训练数据作为支撑的条件下，实现对空中飞机目标的分类。

实现本发明目的的技术思路是:通过提取飞机目标螺旋叶片数目、旋翼转速和旋翼长度这些对分类有用的物理参数信息，然后根据不同种类飞机目标已知旋翼参数的差异，对飞机目标进行分类。具体步骤包括如下：

(1)从飞机雷达回波信号中分离出旋翼回波信号；

(2)利用霍夫变换从旋翼回波时频域S(t_m,f)中估计出旋翼转速ω_r；

(3)利用步骤(2)得到的参数在旋翼回波时频域对旋翼进行二维成像；

3a)对旋翼回波信号进行短时傅里叶变换，获得旋翼回波信号的时频域S(t_m,f)，其中，f为t_m时刻对应的瞬时频率；

3b)遍历搜索旋翼上散射点的位置参数[r,θ]^T，根据步骤(2)估计的旋翼转速、位置参数[r,θ]^T及飞机旋翼回波理论模型获得一组随时间t_m变化的瞬时多普勒频率f_d(r,θ)和相位项exp[jφ(r,θ)]；

3c)从旋翼回波信号的时频域S(t_m,f)中提取一组与[t_m，f_d(r,θ)]对应的复值S(t_m,f_d(r,θ))，利用相位项exp[jφ(r,θ)]对复值S(t_m,f_d(r,θ))进行相位补偿，且对相位补偿后的复值S(t_m,f_d(r,θ))进行积分，得到积分函数G(r,θ)：

G(r,θ)＝|∫S(t_m,f_d(r,θ))·exp[jφ(r,θ)]·exp(j2πft_m)dt_m|²；

3d)当搜索的位置与旋翼上的一个真实散射点位置重合时，相位exp[jφ(r,θ)]与复值S(t_m,f_d(r,θ))的相位相匹配，积分函数G(r,θ)实现相干积累，在(r,θ)位置平面内形成一个峰值点，通过遍历搜索旋翼上散射点位置，在(r,θ)位置平面内形成一组峰值点，该组峰值点的位置分布即为旋翼上散射点的位置，即在(r,θ)位置平面内实现对旋翼的二维成像G；

(4)对步骤(3)获得的旋翼的二维成像G进行图像预处理，获得最终的旋翼图像P；

(5)利用主成分分析PCA投影方法从最终的旋翼图像P中估计出旋翼长度L和螺旋桨叶片数N；

(6)根据步骤(2)和步骤(5)估计出的参数旋翼转速ω_r、旋翼长度L和螺旋桨叶片数N与飞机类型标准库中飞机旋翼转速、旋翼长度及螺旋桨叶片数参数比对，判断出飞机目标的类型。

本发明与现有技术相比的优点

本发明与现有基于微动调制特征的飞机目标分类技术相比的突出优点是不需要大量的仿真或实测训练数据离线学习分类器，即在没有先验数据库支撑的条件下能有效地提取飞机目标螺旋叶片数目、旋翼转速和旋翼长度这些反映飞机物理参数信息，实现对空中飞机目标的分类。

附图说明

图1是本发明的实现流程图；

图2是本发明中对旋翼转速的估计示意图；

图3是本发明中对旋翼的成像结果图；

图4是本发明中对旋翼长度的估计结果图。

具体实施方式

以下结合附图对本发明的实施例做进一步详细描述。

参照附图1，本发明的具体实现步骤如下：

步骤1，从飞机雷达回波信号中分离出旋翼回波信号。

1a)对采集到的飞机回波信号进行离散傅里叶变换，获得信号的多普勒谱；

1b)在多普勒谱中，搜索最大的多普勒谱幅度，并记录该最大的多普勒谱幅度对应的相位φ_body、幅度U_body和多普勒频率f_body；

1c)利用步骤(1b)获得参数，重构该最大的多普勒谱幅度对应的时域信号：S_body＝(U_body/K)exp(j2πf_bodyt_m+jφ_body)，S_body即为机身时域分量，K为离散傅里叶变换的点数；

1d)从飞机雷达回波信号中去除机身分量S_body；

1e)在多普勒谱中，搜索0频附近杂波区域最大的杂波多普勒谱幅度，并记录该最大的杂波多普勒谱幅度对应的相位φ_clutter、幅度U_clutter和多普勒频率f_clutter；

1f)利用步骤(1e)获得的参数，重构该最大的杂波多普勒谱幅度对应的时域信号S_clutter＝(U_clutter/K)exp(j2πf_cluttert_m+jφ_clutter)，S_clutter即为杂波时域分量；

1g)从飞机雷达回波信号中减去机身分量S_clutter；

1h)重复步骤(1e)--(1g)直到迭代次数达到预设迭代门限M，以从飞机雷达回波中去除机身分量和杂波分量，获得旋翼回波信号，其中，M根据杂波谱宽设定。

步骤2，利用霍夫变换从旋翼回波时频矩阵S(t_m,f)中估计出旋翼转速ω_r。

2a)对旋翼回波信号进行短时傅里叶变换，获得旋翼回波信号的时频域S(t_m,f)，其中，f为t_m时刻对应的瞬时频率；

2b)遍历搜索旋翼回波信号时频域S(t_m,f)的频率幅度F和时域周期ω，根据旋翼的瞬时多普勒频率理论模型，获得一条正弦曲线y(F,ω)；

2c)从旋翼回波信号时频域S(t_m,f)中提取一组与正弦曲线y(F,ω)对应的复值S(t_m,y(F,ω))，并对复值S(t_m,y(F,ω))进行求和，得到求和函数Q(F,ω)：

Q(F,ω)＝sum{S[t_m,y(F,ω)]}；

2d)将求和函数Q(F,ω)出现峰值时，其正弦曲线y(F,ω)与旋翼回波时频域S(t_m,f)上的瞬时多普勒频率曲线重合中的正弦曲线y(F,ω)所对应的时域周期ω，作为旋翼转速ω_r。

步骤3，利用步骤(2)得到的参数在旋翼回波时频域对旋翼进行二维成像。

3a)对旋翼回波信号进行短时傅里叶变换，获得旋翼回波信号的时频域S(t_m,f)，其中，f为t_m时刻对应的瞬时频率；

3b)遍历搜索旋翼上散射点的位置参数[r,θ]^T，根据步骤(2)估计的旋翼转速、位置参数[r,θ]^T及飞机旋翼回波理论模型获得一组随时间t_m变化的瞬时多普勒频率f_d(r,θ)和相位项exp[jφ(r,θ)]；

3c)从旋翼回波信号的时频域S(t_m,f)中提取一组与[t_m，f_d(r,θ)]对应的复值S(t_m,f_d(r,θ))，利用相位项exp[jφ(r,θ)]对复值S(t_m,f_d(r,θ))进行相位补偿，且对相位补偿后的复值S(t_m,f_d(r,θ))进行积分，得到积分函数G(r,θ)：

G(r,θ)＝|∫S(t_m,f_d(r,θ))·exp[jφ(r,θ)]·exp(j2πft_m)dt_m|²；

3d)当搜索的位置与旋翼上的一个真实散射点位置重合时，相位exp[jφ(r,θ)]与复值S(t_m,f_d(r,θ))的相位相匹配，积分函数G(r,θ)实现相干积累，在(r,θ)位置平面内形成一个峰值点，通过遍历搜索旋翼上散射点位置，在(r,θ)位置平面内形成一组峰值点，该组峰值点的位置分布即为旋翼上散射点的位置，即在(r,θ)位置平面内实现对旋翼的二维成像G。

步骤4，对旋翼的二维成像G进行图像预处理，获得最终的旋翼图像P。

4a)对旋翼的二维成像G的所有像素值进行对数变换，其表示为：

W(i,j)＝10×In[G(i,j)+0.001]+30；

其中，G(i,j)是旋翼的二维成像G的第(i,j)个像素值，W(i,j)为旋翼的二维成像G所有像素值经对数变换，获得的图像W的第(i,j)个像素值；

4b)利用恒虚警CFAR对上述经过对数变换后的图像W中像素点的坐标进行目标坐标和背景坐标的判断：

其中，μ和σ分别为上述经过对数变换后的图像W所有像素值的均值和方差的估计值，c为虚警率，T代表目标区域像素点的坐标的集合，B代表背景区域像素点的坐标的集合；

4c)将上述经过对数变换后图像W中属于目标区域像素点坐标集合T的像素值置为1，属于背景区域像素点坐标集合B的像素值置为0，得到一个二值图像J；

4d)对二值图像J进行形态学滤波，去除孤立点，获得最终的旋翼图像P。

本发明中，对旋翼的二维成像G的所有像素值进行对数变换，这一处理原因压缩了图像像素值的动态范围，减小图像中大像素值的影响；因为对数变换在零点没有意义，所以在旋翼的二维成像G的所有像素值上任意加一个小的常数0.001，同时为保证对数变换后的图像W的所有像素值大于或者等于零，在最后图像W的所有像素值加一个相应的常数30。

步骤5，利用主成分分析PCA投影方法从最终的旋翼图像P中估计出旋翼长度L和螺旋桨叶片数N；

5a)从最终的旋翼图像P中，提取所有像素值为1的像素点的坐标(i,j)；

5b)对所有像素值为1的像素点的坐标(i,j)进行升序排列，得到坐标矩阵X；

5c)求坐标矩阵X的协方差矩阵：C＝(X^T*X)/M，其中M为像素值为1的像素点数，T代表转置；

5d)对协方差矩阵C进行特征值分解：C＝UΛU^H，其中H代表共轭转置，Λ＝[λ₁,λ₂]为对角矩阵，λ₁,λ₂分别为两个对角线上的元素，表示协方差矩阵C特征值分解后得到的两个特征值，U＝[U₁,U₂]为正交矩阵，U₁,U₂分别为元素λ₁,λ₂对应的特征向量；

5e)对坐标矩阵X进行旋转变换，得到新的坐标矩阵：Y＝(UX^T)^T；

5f)将新的坐标矩阵Y中的所有元素分别向所述的两个特征向量U₁，U₂方向上投影，分别得到新坐标矩阵Y中所有元素在所述特征向量U₁上的长度投影系数向量和在所述特征向量U₂方向上的初相角投影系数向量

$Y_{U_1}=Y^T*U_1,$ Y_U2＝Y^T*U₂；

5h)估计出翼长度L:

(5h1)根据长度投影系数向量得到最小长度投影系数和最大长度投影系数

(5h2)根据初相角投影系数向量得到最小初相角投影系数和最大初相角投影系数

(5h3)根据步骤(5h1)和(5h2)得到的参数估计出翼长度：L＝sqrt[(L₂-L₁)²+(α₂-α₁)²]，其中sqrt为开根号。

步骤6，将步骤(2)和步骤(5)估计出的参数旋翼转速ω_r、旋翼长度L和螺旋桨叶片数N与飞机类型标准库中飞机旋翼转速、旋翼长度及螺旋桨叶片数参数比对，判断出飞机目标的类型。

本发明的效果可以通过以下仿真实验验证：

1.仿真条件：

实验数据为在旋翼电磁仿真回波信号中添加20dB的高斯白噪声的旋翼回波信号，其中，旋翼回波信号中旋翼物理参数为旋翼末端长度L₂＝6.2m，旋翼始端长度L₁＝0.3m，旋翼长度L＝5.9m，旋翼转速ω_r＝5r/s，螺旋桨叶片数N＝4；

2.仿真内容及结果分析：

仿真1，利用本发明对上述物理参数的旋翼回波信号进行转速估计，结果如图2，从图2可以，计算出的旋翼转速ω＝1/[(0.2384-0.03836)]＝4.999r/s。

仿真2，利用本发明对上述物理参数的旋翼回波信号进行二维成像，结果如图3。

仿真3，利用本发明对上述物理参数的旋翼回波信号的螺旋桨叶片数及旋翼长度估计，结果如图4，从图4中可看出估计出的螺旋桨叶片数N＝4，旋翼长度为：L＝{sqrt([(5.916-0.1013)²+(0.3425-0.3678)²])+sqrt([(5.865-0.2554)²+(1.931-1.895)²])}/2＝5.725m

将图2、图3和图4估计的旋翼物理参数旋翼转速、旋翼长度及螺旋桨叶片数与实验中设定的旋翼转速、旋翼长度及螺旋桨叶片数这些物理参数进行对比，可以得出：利用本发明能较准确地估计旋翼转速、旋翼长度和螺旋桨叶片数。

同时将图2、图3和图4估计的旋翼物理参数旋翼转速、旋翼长度及螺旋桨叶片数与表1中列举的常见15架三类飞机的旋翼物理参数进行比对，可以得出图2、图3和图4估计的旋翼参数旋翼转速、旋翼长度及螺旋桨叶片数与直升机的旋翼物理参数相似，判断出实验中具有该旋翼物理数的飞机为直升机。

表1：常见五架直升机、五架螺旋桨飞机及五架喷气式飞机的旋翼参数

Claims (5)

1.一种基于旋翼物理参数估计实现飞机目标分类的方法，包括如下步骤：

(1)从飞机雷达回波信号中分离出旋翼回波信号；

(2)利用霍夫变换从旋翼回波时频域S(t_m,f)中估计出旋翼转速ω_r；

(3)利用步骤(2)得到的参数在旋翼回波时频域对旋翼进行二维成像；

3a)对旋翼回波信号进行短时傅里叶变换，获得旋翼回波信号的时频域S(t_m,f)，其中，f为t_m时刻对应的瞬时频率；

3b)遍历搜索旋翼上散射点的位置参数[r,θ]^T，根据步骤(2)估计的旋翼转速、位置参数[r,θ]^T及飞机旋翼回波理论模型获得一组随时间t_m变化的瞬时多普勒频率f_d(r,θ)和相位项exp[jφ(r,θ)]；

3c)从旋翼回波信号的时频域S(t_m,f)中提取一组与[t_m,f_d(r,θ)]对应的复值S(t_m,f_d(r,θ))，利用相位项exp[jφ(r,θ)]对复值S(t_m,f_d(r,θ))进行相位补偿，且对相位补偿后的复值S(t_m,f_d(r,θ))进行积分，得到积分函数G(r,θ)：

G(r,θ)＝|∫S(t_m,f_d(r,θ))·exp[jφ(r,θ)]·exp(j2πft_m)dt_m|²；

3d)当搜索的位置与旋翼上的一个真实散射点位置重合时，相位exp[jφ(r,θ)]与复值S(t_m,f_d(r,θ))的相位相匹配，积分函数G(r,θ)实现相干积累，在(r,θ)位置平面内形成一个峰值点，通过遍历搜索旋翼上所有散射点位置，在(r,θ)位置平面内形成一组峰值点，该组峰值点的位置分布即为旋翼上散射点的位置，即在(r,θ)位置平面内实现对旋翼的二维成像G；

(4)对步骤(3)获得的旋翼的二维成像G进行图像预处理，获得最终的旋翼图像P；

(5)利用主成分分析PCA投影方法从最终的旋翼图像P中估计出旋翼长度L和螺旋桨叶片数N；

(6)根据步骤(2)和步骤(5)估计出的参数旋翼转速ω_r、旋翼长度L和螺旋桨叶片数N与飞机类型标准库中飞机旋翼转速、旋翼长度及螺旋桨叶片数参数比对，判断出飞机目标的类型。

2.根据权利要求1所述的基于旋翼物理参数估计实现飞机目标分类的方法，其中步骤(1)所述的从飞机雷达回波信号中分离出旋翼回波信号，按如下步骤进行：

(1a)对采集到的飞机回波信号进行离散傅里叶变换，获得信号的多普勒谱；

(1b)在多普勒谱中，搜索最大的多普勒谱幅度，并记录该最大的多普勒谱幅度对应的相位φ_body、幅度U_body和多普勒频率f_body；

(1c)利用步骤(1b)获得参数重构该最大的多普勒谱幅度对应的时域信号S_body＝(U_body/K)exp(j2πf_bodyt_m+jφ_body)，S_body即为机身时域分量，K为离散傅里叶变换的点数；

(1d)从飞机雷达回波信号中去除机身分量S_body；

(1e)在多普勒谱中，搜索0频附近杂波区域最大的杂波多普勒谱幅度，并记录该最大的杂波多普勒谱幅度对应的相位φ_clutter、幅度U_clutter和多普勒频率f_clutter；

(1f)利用步骤(1e)获得的参数重构该最大的杂波多普勒谱幅度对应的时域信号S_clutter＝(U_clutter/K)exp(j2πf_cluttert_m+jφ_clutter)，S_clutter即为杂波时域分量，

(1g)从飞机雷达回波信号中减去机身分量S_clutter；

(1h)重复步骤(1e)--(1g)直到迭代次数达到预设迭代门限M，以从飞机雷达回波中去除机身分量和杂波分量，获得旋翼回波信号，其中，M根据杂波谱宽设定。

3.根据权利要求1所述的基于旋翼物理参数估计实现飞机目标分类的方法，其中步骤(2)所述的利用霍夫变换从旋翼回波时频域S(t_m,f)中估计出旋翼转速ω_r，按如下步骤进行：

(2a)对旋翼回波信号进行短时傅里叶变换，获得旋翼回波信号的时频域S(t_m,f)，其中，f为t_m时刻对应的瞬时频率；

(2b)遍历搜索旋翼回波信号时频域S(t_m,f)的频率幅度F和时域周期ω，根据旋翼的瞬时多普勒频率理论模型，获得一条正弦曲线y(F,ω)；

(2c)从旋翼回波信号时频域S(t_m,f)中提取一组与正弦曲线y(F,ω)对应的复值S(t_m,y(F,ω))，并对复值S(t_m,y(F,ω))进行求和，得到求和函数Q(F,ω)：

Q(F,ω)＝sum{S[t_m,y(F,ω)]}；

(2d)当求和函数Q(F,ω)出现峰值时，将正弦曲线y(F,ω)与旋翼回波时频域S(t_m,f)上的瞬时多普勒频率曲线重合中的正弦曲线y(F,ω)所对应的时域周期ω，作为旋翼转速ω_r。

4.根据权利要求1所述的基于旋翼物理参数估计实现飞机目标分类的方法，其中步骤(4)所述的对旋翼的二维成像G进行图像预处理，获得最终的旋翼图像P，按如下步骤进行：

(4a)对旋翼的二维成像结果G的所有像素值进行对数变换，其表示为：

W(i,j)＝10×In[G(i,j)+0.001]+30；

其中，G(i,j)是旋翼的二维成像结果G的第(i,j)个像素值，W(i,j)为旋翼的二维成像结果G所有像素值经对数变换获得的图像W的第(i,j)个像素值；

(4b)利用恒虚警CFAR对上述经过对数变换后的图像W中像素点的坐标进行目标坐标和背景坐标的判断：

其中，μ和σ分别为上述经过对数变换后的图像W所有像素值的均值和方差的估计值，c为虚警率，T代表目标区域像素点的坐标的集合，B代表背景区域像素点的坐标的集合；

(4c)将上述经过对数变换后图像W中属于目标区域像素点坐标集合T的像素值置为1，属于背景区域像素点坐标集合B的像素值置为0，得到一个二值图像J；

(4d)对二值图像J进行形态学滤波，去除孤立点，获得最终的旋翼图像P。

5.根据权利要求1所述的基于旋翼物理参数估计实现飞机目标分类的方法，其中步骤(5)所述的利用主成分分析PCA投影方法估计旋翼长度L，按如下步骤进行：

(5a)从最终的旋翼图像P中，提取所有像素值为1的像素点的坐标(i,j)；

(5b)对所有像素值为1的像素点的坐标(i,j)进行升序排列，得到坐标矩阵X；

(5c)求坐标矩阵X的协方差矩阵：C＝(X^T*X)/M，其中M为像素值为1的像素点数，T代表转置；

(5d)对协方差矩阵C进行特征值分解：C＝UΛU^H，其中H代表共轭转置，Λ＝[λ₁,λ₂]为对角矩阵，λ₁,λ₂分别为两个对角线上的元素，为对协方差矩阵C特征值分解后得到的特征值，正交矩阵U＝[U₁,U₂]，U₁,U₂分别为元素λ₁,λ₂对应的特征向量；

(5e)对坐标矩阵X进行旋转变换，得到新的坐标矩阵：Y＝(UX^T)^T；

(5f)将新的坐标矩阵Y中的所有元素分别向所述的两个特征向量U₁，U₂方向上投影，分别得到新坐标矩阵Y中所有元素在所述特征向量U₁上的长度投影系数向量和在所述特征向量U₂方向上的初相角投影系数向量

Y_U2＝Y^T*U₂；

(5h)估计旋翼长度L:

(5h1)根据长度投影系数向量得到最小长度投影系数和最大长度投影系数

(5h2)根据初相角投影系数向量得到最小初相角投影系数和最大初相角投影系数

(5h3)根据步骤(5h1)和(5h2)得到的参数估计出旋翼长度：L＝sqrt[(L₂-L₁)²+(α₂-α₁)²]，其中sqrt为开根号。