考虑微多普勒特征的鸟类目标多假设跟踪算法∗

（海军工程大学电子工程学院 武汉 430033）

1 引言

鸟击是指航空器飞行和起降过程中和鸟类、蝙蝠等飞行物相撞的事件。鸟击给军队和民航部门带来巨大的财产损失，同时也给机上工作人员和乘客的生命安全造成了巨大威胁［1］。防范鸟击事件的前提是对鸟类目标进行实时探测跟踪，雷达探鸟利用自身辐射的电磁波来探测鸟类，具有自动化程度高、连续工作时间长、数据便于分析存储、受天气和光线条件影响小等诸多优点，故使用雷达来探测跟踪鸟类目标具有天然的优势［2～4］。

鸟类目标具有飞行高度低、目标小、速度慢等特点，其雷达回波信号弱，容易被杂波和噪声信号掩盖，目标出现漏检的概率较大；同时，鸟类经常聚群飞行，数量多，机动性强，空间运动轨迹复杂，航迹交叉与分离频繁，雷达航迹管理难度大［3］。因此，在使用雷达探测跟踪鸟类目标时，对跟踪算法要求非常高，而在大量目标跟踪算法中，多假设跟踪（MHT）算法是能够满足此类目标探测跟踪需求的算法之一［4～5］。

鸟类翅膀的扑动会产生微多普勒效应，且不同的鸟，在不同的姿态和雷达视角下产生的微多普勒特征也不相同［6～7］，因而可利用鸟类产生的微多普勒特征进行鸟类目标探测跟踪。本文将在MHT算法中引入微多普勒特征参量，以提高MHT算法对鸟类目标探测跟踪的准确性。

2 利用鸟类目标微多普勒特征进行参数估计

建立鸟类目标运动模型如图1所示［6］。忽略鸟身的宽度，以翅膀与鸟身的连接处为坐标原点O，鸟飞行方向为X轴正方向，建立如图所示目标坐标系（X，Y，Z）。为简化分析，假设翅膀仅在YOZ平面内作扑翼运动，且扑翼过程中翅膀不弯折，扑翼频率为F0，半翼展为L，扑翼角度ψ（即翼面与XOY平面的夹角）为随时间变化的函数：

式中，A（0＜A＜π）为扑翼幅度，ω0=2πF0表示扑翼角频率，ψ0表示初始时刻扑翼角度的相位，意义不大，为简化分析，令ψ0=0。雷达坐标系（U，V，W）与目标坐标系平行，雷达位于原点S，到O点的距离为R0。雷达视线方向在USV平面上的投影与U轴的夹角为雷达方位角α，与USV平面的夹角为雷达俯仰角β。

图1 鸟类目标运动模型与雷达的位置

通过推导，鸟整体的回波信号s（t）为鸟身和左右翅膀回波信号之和，s（t）的微多普勒特征为鸟身回波的频移fd和左右翅膀回波微多普勒特征fmD1（t）和fmD2（t）的总和，各分量表达式如式（2）～（4）：

利用微多普勒特征进行参数估计的步骤如下：

1）计算回波信号的自相关函数，估计目标扑翼频率ω0。

可知散射点P′和Q′到雷达的距离由初始距离、平动和微动三部分构成，其中平动速度根据回波频谱容易进行估计，微动部分是以T0=2π/ω0为周期的周期函数［6，8］。由复合函数的性质可知，平动补偿后的目标回波s（t）也是周期为T0的时间序列。回波信号的自相关函数和平均幅度差函数均可用于估计信号周期，假设平动补偿后的雷达回波为长度为N的周期序列x（n）（n=1，2，…，N），x（n）的周期延拓为，其自相关函数Ф1（k）和平均幅度差函数Ф2（k）为

当信号长度超过一个周期时，其自相关函数Ф1（k）将在0，±T0，±2T0，…处出现极大值，平均幅度差函数Ф2（k）将0，T0，2T0，…在处出现极小值。目标的扑翼频率和回波信号的频率相等，因此可以通过估计回波信号的周期得到扑翼频率的估计值。

2）对目标回波信号进行时频分析，提取目标微多普勒频率。

3）对目标微多普勒频率进行快速傅里叶变换，得到各谐波分量的幅度绝对值c1～c4。

4）根据雷达视角，选择使用式（6）或者式（7）求得到扑翼幅度估计值。

5）根据式（3），令

则式（3）可以重写为

将扑翼幅度的估计值A和扑翼频率的估计值ω0代入式（9），可以得到g1（t）的估计值。定义微多普勒谱宽BmD1和函数g1（t）的谱宽Bg1为

由式（10）可知：

式（11）中，BmD1由目标回波获取，扑翼幅度的估计值A和扑翼频率的估计值ω0已知，再代入其他参数，即可计算得到半翼展的估计值。

同理，可由式（11）计算得到半翼展的估计值，因此，鸟类目标半翼展的估计值确为和的均值：

3 面向航迹的MHT算法

面向航迹的MHT算法主要分为：1）对现有的航迹进行预测，将新的量测和预测波门进行关联，生成假设得到新的航迹；2）计算新航迹得分，并基于得分进行航迹删除和确认；3）对航迹删除和确认后幸存航迹进行聚类；4）在每个聚类中独立生成全局假设，并寻找最优全局假设；5）基于最优全局假设进行N扫描回溯剪枝，并输出最优航迹［9］。

3.1 航迹生成

如图2所示，假定在某个区域1～k时刻收到了多批量测。k-1时刻，已经生成了多个航迹，可能存在的航迹如图中实线，并且通过卡尔曼预测，得到了量测②的跟踪波门和量测③的跟踪波门。在k时刻，收到了五个量测，其中量测①、②和④落入波门中，量测①和③落入波门中，量测⑤没有落入任何波门（如图3所示）。因此，将k时刻的量测和k-1时刻的航迹关联，可得到如图2中虚线所示假设：1）k时刻五个量测都是新目标；2）k时刻五个量测都是虚警；3）k时刻量测①、②和④是k-1时刻量测③的继续，k时刻量测①和③是k-1时刻量测②的继续，其中，k-1时刻量测②可能是k-2时刻两个量测的继续，因此从该量测到k时刻量测①和③都分别有两个假设。同时，可以得到1～k时刻的航迹树，如图4所示。

图2 k-1时刻航迹

图3 k时刻量测关联

图4 k时刻航迹树

3.2 航迹得分计算

航迹得分由递归累积产生。每一个航迹的得分等于它的上一次的值加上一个得分增量Δlk，即：

其中，Zk是目标测量值，Xk|k-1是目标预测位置，S是残差协方差矩阵；Pd为目标的探测概率；βf为虚警的空间密度；Pf为虚警概率；M为量测的维数。且航迹初始得分为l1=ln[(1-βf)/βf]。

通过航迹得分，可以进行航迹的删除和确认，即lk≤ TL，则删除航迹；lk≥ TU，则确认航迹；TL≤ lk≤TU，等待更多的数据更新航迹，其中TL定义为航迹删除阈值值，TU为航迹确认阈值。这样很多得分低的航迹在此过程中被删除，大大减少了后面生成假设的数量［10～11］。

3.3 聚类

由于某些目标在时空中的可区分性很高，因此可将所有的假设航迹划分成子集，在每个子集中独立地进行假设生成、全局级航迹删减和航迹更新等操作，即可将大规模跟踪问题化解为若干个小规模跟踪问题来处理，降低了算法复杂度和计算量，假设航迹子集即为航迹类。

图5 k时刻聚类及全局假设生成

在新扫描周期内分别将新接收到的量测与以前的类进行互联，若上一周期处理中两个独立的类与同一个量测相关，则这两个类合并成一个类，不与任何旧聚类相关的量测形成新的类。将图4中的航迹树进行聚类，由于航迹1和航迹2有共同量测，必须合并成一个类，航迹3为一个新类，最终得到两个类，如图5所示。

3.4 全局假设生成

在每个类中，利用图论算法，将假设航迹映射为图论中图的顶点，航迹共享量测则对应顶点是相邻的（即顶点间存在边）。全局假设生成即为顶点集的划分操作，根据聚类原则，相邻顶点必在同一集合中（如图5），而寻找全局最优假设则为寻求极大连通子集的过程，此过程可由深度优先搜索算法实现。图5中，深色部分即为全局最优假设［12］。

3.5 N-scan剪枝

N-scan剪枝法通过限制航迹树深度来控制假设数量，强制在k-N时刻产生的不确定性在k时刻延迟解决。当轨迹树的深度大于N时，N-scan剪枝法将保留全局最优假设中的根节点分枝，删除其余分枝［10］。如图6所示，N=2时在k时刻航迹树执行N-scan剪枝操作，上一步已得到最优假设为132、223、005，保留该分枝，而删除其他分枝，最终得到k时刻的航迹（如图7所示）。

图6 k时刻剪枝

图7 k时刻航迹

4 引入微多普勒特征的MHT算法

考虑在MHT算法的第二步，即计算航迹得分时引入微多普勒频率。利用雷达历史回波得到目标的参数估计，根据MHT算法中与目标相关联量测的坐标得到目标的雷达视角，将目标的参数估计和雷达视角代入式（2）～（4），计算目标的多普勒特征理论值ft。根据当前回波提取目标多普勒特征观测值fa，并与理论值ft做对比，将两者差值的均方根作为目标多普勒特征得分lf。给式（13）中MHT算法航迹得分lk分配权重ak，目标多普勒特征得分lf分配权重af，则假设航迹最终得分为

其中，权重ak、af与雷达性能、目标特征、雷达工作环境有关，需经过多次实验确定。

5 仿真计算

仿真中，模拟产生两只不同种类的鸟目标，其半翼展、扑翼幅度和扑翼频率均不相同，目标1做直线运动，目标2做曲线运动，两个目标在空间交汇两次（如图8所示）。杂波点个数服从泊松分别，杂波密度3×10-5，N-scan剪枝深度为3，量测噪声服从高斯分布，标准差为100m。分别用本文算法和传统MHT算法进行100次Monte Carlo仿真计算，得到两种算法失跟率如表1所示，两种算法对两个目标跟踪的距离均方根误差曲线分别如图9、10所示。

图8 目标轨迹

表1 算法失跟率统计表

图9 目标1跟踪误差RMSE曲线

从表1可以看出，在目标跟踪稳定跟踪率上，本文算法要优于传统MHT算法。而从图9、10可以看出，在位置的跟踪精度上，两种算法差距不大。另外，将图8放大，可以看到在两个目标交汇处，点迹距离很近，航迹关联难度更大，对算法的性能提出了挑战。通过统计，传统MHT算法将交汇处点迹错误关联的概率高达30%，而本文算法错误关联的概率仅1%。因此，在航迹精确关联方面，本文算法性能要优于传统MHT算法。

图10 目标2跟踪误差RMSE曲线

图11 目标轨迹局部图

6 结语

本文介绍了根据微多普勒特征进行鸟类目标参数估计的方法和面向对象的MHT算法基本流程，提出利用参数估计的结果计算鸟类目标微多普勒特征得分，并将该得分引入MHT算法，作为MHT算法中航迹得分的组成部分，从而实现了考虑微多普勒特征的鸟类目标MHT算法。通过仿真计算表明，该算法与传统MHT算法相比，失跟率更低，且在处理交叉航迹时，航迹关联精度更高。