基于Viterbi 算法和ROMP 的多弹道目标分离与特征提取*

2020-05-11冯存前许旭光

火力与指挥控制 2020年3期

陈帅，冯存前，许旭光

（1.空军工程大学研究生院，西安 710051；2.空军工程大学防空反导学院，西安 710051）

0 引言

如今，随着弹道目标技术的发展，为保证杀伤效果，提高弹头生存能力，在弹道导弹飞行中段，往往有若干真弹头和伴飞的假弹头等干扰目标。由于这些目标的尺寸相近，在窄带体制雷达信号下，目标回波信号复杂混叠，难以分离并提取各个目标的参数［1-3］。邵长宇等人利用目标跟踪技术，将目标回波的时频曲线看作机动目标的航迹，从而完成了散射点频率信息的提取。高昭昭等人利用时间距离像上微多普勒曲线的正弦形式特点实现对目标的分离与快速参数提取。束长勇等人利用matlab 自带函数提取时频图像中的微多普勒曲线，完成了目标的微动类型分类［4-6］。上述方法使用的微动模型单一，参数提取时利用的先验信息少，不适合用于复杂微动目标的分离及参数提取。

本文以单频窄带脉冲作为雷达的发射信号，同时对旋转和滑动两种不同的微动目标建模，分析二者的运动特征。

1 多弹道目标模型

图1 弹头微动模型图

建立如图1 所示的微动模型，OXYZ 为雷达坐标系，o1x1y1z1为旋转目标参考坐标系，o2x2y2z2为滑动目标参考坐标系，分别以弹头A 和弹头B 的质心o1、o2为坐标原点，各对应坐标轴与雷达坐标轴相平行。弹头A 绕自旋轴以角速度ω1旋转，初始时刻o1在雷达坐标系中的位置为（α1，β1）。θ 为弹头B 对称轴与锥旋轴的夹角即进动角。初始时刻o2在雷达坐标系中位置为（α2，β2）。则LOSi为

设弹头A 上p1点散射中心。根据罗德里戈旋转公式，已知原向量、旋转轴以及它们的叉乘积可以得到旋转后的向量［3］。文献［3］推导出已知旋转轴方位角和俯仰角以及散射点绕轴旋转角速度下的旋转矩阵，则旋转微动参考坐标系下的旋转矩阵可以表示为

由于实际目标的微动距离远小于目标到雷达的径向距离，因此，在t 时刻散射点p1到雷达坐标系原点的微距离为

由于目标的遮挡效应，单部雷达观测滑动目标可等效为锥顶p2和锥底p3两个散射点。通过推导计算可知［7-9］，p2、p3两点到雷达坐标系原点的微距离为

由式（4）、式（5）可知，锥顶等效散射中心p2点的微距离满足正弦形式特点，锥底等效散射点p3点多出了复杂非正弦调制项，因而不再满足正弦变化规律，但仍具有周期变化规律。

2 目标回波处理

2.1 回波模型分析

设窄带雷达发射载频为fc的脉冲，单个脉冲宽度为，脉冲重复间隔为T。则雷达发射信号表达式为

其中，门函数

假设在t 时刻，等效散射中心pi到雷达的微距离为Rpi，则雷达接收到的回波为

式中，σ 为散射点散射系数，经过正交双通道解调后，目标散射点的回波为

通过上面的推导，已经得到目标等效散射中心pi的雷达回波信号和微多普勒模型，可以看出微动目标的回波是典型的非平稳信号，简单的频域分析，难以获得更多有效信息。利用时频分析从时域和频域上同时描述信号的能量密度，可以确定在特定时间信号的能量集中分布的频点，可以用时频分析得到回波信号的微多普勒频率变化曲线。

2.2 重排Gabor 分布

时域信号可以写成

其中，amn为分解系数，与信号自身相关。正交解调后的回波信号Gabor 分布为

在Gabor 变换后，将每个时频点处的幅值移动到该点附近的能量中心，从而得到时频信号sb（t）的重排Gabor 幅值。

2.3 目标分离与分类

根据能量法则，可以通过时频图像处理来获得频率变化曲线，最常用的方法是按照能量最大原则，寻找每个时刻能量最大的频点，没有考虑时频曲线上相邻两个时刻的幅值相关性，在处理多分量信号时效果不佳［11-12］。因此，提取多分量时频曲线时应当同时考虑信号时频幅值的连续性和信号时频分布能量的大小。这样分离曲线的任务就等效为在时频图像中寻找能量之和最大，并且任意两点的变化率要尽可能小的路径。构造目标寻优函数

利用Viterbi 算法求解目标寻优函数，在时频分布图上依次寻找出各个散射点时频曲线对应的最优路径。Viterbi 算法是常用的估计信号最可能的隐状态动态规划算法，在每个时刻寻找当前最优状态，有效减小了运算量并提取出各个散射点的时频曲线。

图2 时频曲线分离示意图

由文献［13］的结论，通过对分析散射点微多普勒变化曲线的频谱可以判断出目标的微动形式，故本文利用这一结论，对分离出的目标散射点曲线做频域分析，按照频域的分布特点判断出目标的运动状态，由于分离多分量信号在交叉点处难以判断，可能会产生一些误差甚至是断点，本文采取有理数拟合的方法，获取平滑的特征曲线，这样能更好地提取出目标的特征。

3 ROMP 算法在特征提取中的应用

贪婪类算法由于其原理简单且易于实现，在稀疏表示的系数求解中得到了广泛应用［14-15］。为了将一组信号分解到非正交基上，并求解信号的稀疏系数，Mallat 在1993 年提出了匹配追踪算法。由于匹配追踪算法信号在选中的原子中进行垂直投影是非正交性的，因此，该算法收敛需要多次迭代。为了解决这个问题，正交匹配追踪算法被提出，它的改进在于分解的每一步都对选择的所有原子进行正交化处理。在上面的基础上，Deann 提出了正则化正交匹配追踪算法，将正则化应用到原子的选择过程中，每次选出与信号或残差内积最大的k 个原子，在这k 个原子中，按照正则化标准选出本次迭代的原子，即所选原子与信号的内积最大值不能超过最小值的两倍且能量最大。选出原子后，再用最小二乘法求解分解系数，最终将信号迭代分解。

由于ROMP 算法分解效率高且简单易行，本文将它引入到雷达信号处理中。对弹道目标的雷达回波进行ROMP 分解，通过选出的原子确定回波信号的参数，从而实现弹道目标旋转微动的特征提取。由ROMP 算法的原理可知，原子库的构造必须考虑到待分解信号的各个特征。由式（14）、式（15）可知，旋转目标散射点和滑动目标锥顶散射点的雷达回波基带信号为正弦调频信号，故构造原子

滑动目标锥底散射点的雷达回波基带信号为比较复杂的调制形式，构造原子为

表1 ROMP 特征提取算法流程

本文提取的整个信号处理流程如图3 所示。将获得雷达回波信号进行预处理后得到降噪后的基带信号，首先对信号进行重排Gabor 时频变换，从变换得到的时频图像中提取各个散射点的时频曲线。对曲线进行拟合和频域分析判断得到散射点的运动类型，针对不同类型构造相应的原子库，利用正则化正交匹配分解对曲线进行特征提取，最终完成目标分离和特征提取。

图3 算法流程示意图

4 仿真分析

本文对两个处于不同微动状态下的锥形弹头进行仿真。初始时刻，弹头A 质心在雷达坐标系中的方位角为60°，俯仰角为60°，弹头锥顶到质心的长度为1.5 m。弹头A 绕自旋轴以角速度ω=3π rad/s旋转，自旋轴在参考坐标系中的方位角为35π/18，俯仰角为4π/9。初始时刻，弹头B 弹体对称轴与参考坐标系z1轴（锥旋轴）的夹角为2π/9，进动角速度为1.5π rad/s，弹头锥顶到质心的长度为1.8 m，弹头总高度为2.4 m。SNR=5 dB 时，对接收信号进行重排Gabor 分析，得到信号的时频分布图像，可以看到重排Gabor 分布得到的时频图的时频聚焦性要明显好于文献［12］提到的Margenau-Hill 变换和短时傅里叶变换。

由下页图4 可知，本文中的重排Gabor 变换在非平稳信号回波的处理上有最好的时频聚焦性和连续性，每个散射点对应的时频曲线都清晰可见，对接下来的分离提取非常有利。对得到的图像使用奇异值分解降噪处理后，利用Viterbi 算法依次寻找3 个散射点对应的时频图中的微动曲线最优路径。可以看出在初始时刻交叉点附近算法的分离效果不够理想，在初期很小的时间段内甚至出现了误判。为了使分离后的曲线也具有连续性，便于频谱分析，本文利用有理数拟合的方法分别拟合分离后的3 个散射点所对应的曲线。得到的结果如图6 所示。

图4 时频曲线图

图5 Viterbi 分离结果

图6 散射点1、2、3 微多普勒曲线拟合图

图7 散射点1、2、3 微多普勒曲线频谱图

对各个散射点的微多普勒曲线进行频谱分析，由文献［13］归纳的各种微动类型的频谱特点，可以判断出散射点3 的频谱单一，为旋转目标的等效散射点。散射点2 的频点分布在n 倍的第一频点上，判断为锥底散射点。散射点1 的频点与散射点2 第一频点相等，判断为滑动目标的锥顶散射点。从而可以判断散射点1 为弹头A 的等效散射点，散射点2、3 为弹头B 的等效散射点。

对不同运动类型的散射点构造不同的匹配原子库，利用ROMP 对多普勒特征曲线进行特征提取，设迭代次数为5。提取出的目标的参数如表2 所示。

将提取的部分参数代入方程组可以求解出弹头A 的自旋速度为9.61 rad/s，弹头B 的弹头高度为2.37 m，进动角为0.767，进动角速度为4.73 rad/s，与理论值均十分接近，充分验证了算法的有效性。

5 结论

本文主要提出了一套从分离、分类到特征提取的算法流程。在分离方面，应用时频分析对信号回波进行时频域复合分析，利用Viterbi 算法优化分离的计算量，算法优点在于使用重排Gabor 分布，有效提高了时频曲线的时频聚焦性。在特征提取方面，本算法的优点在于先判断出各个回波分量的运动形式，根据不同形式构造不同匹配原子库，使用正则化正交匹配算法提取特征。通过构造最合适的匹配原子库，有效提高特征提取的精度和效率。