韩 伟, 王国师, 晏 凯, 宋亚伟
(空军预警学院雷达士官学校, 湖北 武汉 430345)
由于技术体制的限制,机载预警雷达存在固有的多普勒盲区问题[1-4]。在目标跟踪过程中,多普勒盲区造成的点迹连续丢失现象极易引起目标航迹的中断。因此,国内外很多学者致力于研究多普勒盲区条件下的目标跟踪问题。该问题的实质是将雷达多普勒盲区作为先验信息应用到各类跟踪算法中。Gordon等[5]最早将该先验信息运用到基于粒子滤波的目标跟踪算法中。在此基础上,很多学者针对目标不同的机动场景提出了机动模型下的改进盲区粒子滤波算法[6-9]。Clark等[10-11]则将多普勒盲区先验信息并入高斯和滤波算法中,获得了较好的地面动目标指示雷达对地面目标的跟踪性能。还有一些学者针对地面目标的运动特性,在多模型中增加“停止”状态模型,解决其盲区跟踪问题[12-14]。文献[15-17]则根据目标多普勒盲区先验信息提出了多假设运动模型,解决了空中目标的盲区跟踪问题。但以上研究均是在有杂波或无杂波的单目标环境下开展的,在多目标杂波环境下并不适用。
针对多普勒盲区下的多目标跟踪问题已进行了一些研究。第1类方法是在随机有限集框架下将多普勒盲区信息融入基于高斯混合的跟踪算法中[17-18],获得较好的多目标跟踪性能;第2类方法是基于航迹-航迹关联的思想,在目标航迹已经发生断裂的情况下,研究断裂前航迹与断裂后航迹的关联方法[19-21],但该类方法需要在航迹形成较好的片段下才能有较好的性能;第3类方法是基于边跟踪边关联的思想,研究多普勒盲区下的多目标点迹-航迹数据关联方法[22-24]。文献[22-23]将分配算法应用到点迹-航迹关联中,且综合考虑了多普勒盲区和目标检测概率小于1的两类测量值丢失的原因,但上述研究针对的是GMTI 雷达对地面目标跟踪的背景,不适用于空中目标跟踪。文献[24]根据地面目标和空中目标进入多普勒盲区的运动方式不同,研究了针对空中目标的点迹-航迹关联方法,借用文献[22-23]算法中第2类空量测的思想,进行点迹-航迹的最优二维分配。但二维分配仅利用了当前帧的量测信息完成与前一帧目标航迹的配对,并未充分利用基于多帧量测的航迹演变信息。本文在二维分配算法的基础上,提出了多普勒盲区条件下的基于多维分配算法的目标点迹-航迹关联方法,其思想为通过构建并入两类空量测的多维分配的代价函数,完成雷达当前帧测量数据与多个旧时刻的目标点迹、航迹的配对,代价函数最小的配对即为最终的关联结果。
为减小点迹-航迹二维分配计算量,充分运用目标的速度和位置限制进行初步关联,筛选出候选测量值集合。
(1) 利用速度限制的关联
根据目标的最大速度,可以在k-1时刻状态估计值的条件下预测k-1时刻的测量值范围,经过速度关联的k时刻测量值集合可表示为
(1)
(2) 利用距离误差的关联
在极坐标系下,假设位置测量误差服从高斯分布,则经过无偏转换后,直角坐标下的测量误差近似服从高斯分布,故归一化新息平方服从χ2分布,自由度为测量值的维数。航迹n在k时刻的新息表示如下:
(2)
经过位置关联的测量值集合可表示为
(3)
1.2.1 二维分配算法的代价函数
测量值经过筛选后,剔除了不可能参与关联的点迹,使得点迹-航迹配对的数量大为减小。对于二维分配算法,如图1所示,配对的对象为k-1时刻的目标航迹和k时刻的目标点迹(测量值)。左边航迹集合中的标号0为空航迹,用n=0表示,与该航迹关联的测量值为虚警,右边的测量值集合中的标号0表示空测量值,表示雷达发生漏检,无法获取目标的测量值。其中,由目标检测概率Pd<1造成的空测量值称为第1类空量测,由多普勒盲区引起的空测量值称为第2类空量测,分别用mk=01和mk=02表示。
图1 点迹-航迹二维分配示意图Fig.1 Schematic deagram of two-dimension assignment of plot-to-track
点迹-航迹二维分配目标为全局代价函数最小,因此可等效为求解全局代价函数的最小值问题[25],表示如下:
(4)
(5)
式中:a(k,n,mk)为二进制分配变量;c(k,n,mk)为a(k,n,mk)的代价;Nk-1为k-1时刻航迹的数量;Mk为k时刻测量值的数量。
1.2.2 多维分配算法的代价函数
二维分配算法用于完成当前帧测量值与前一帧目标航迹的关联,这种单测量值与航迹的配对没有充分利用目标航迹演变信息[26]。因此,可以利用雷达从k-S+2时刻到k时刻(当前时刻)的测量序列实现与k-S+1时刻航迹的关联,从而形成点迹-航迹的S维分配,这里S≥3。如图2所示,与k-S+1时刻的航迹关联的测量值除了当前k时刻的测量值,还包括k时刻以前的测量值,即时间窗内的S-1维测量值序列,可表示如下:
图2 点迹-航迹多维分配示意图Fig.2 Schematic diagram multiple-dimension assignment of plot-to-track
Zmk-S+2,…,mk={zmk-S+2(k-S+2),…,zmk(k)}
(6)
k-S+2时刻的航迹则可表示为Tn(k-S+1)。
定义多维分配条件下的二进制变量如下:
(7)
借鉴二维分配算法代价函数的表示方法,多维分配的代价函数也可表示为负对数似然比的形式:
(8)
(9)
(10)
在点迹-航迹关联过程中,将二维分配算法中的两类虚测量值以及与虚警对应的空航迹引入到多维分配的似然函数中,故多维分配的似然函数可以表示为
(11)
式中:f(ms)和g(ms)为二进制函数,分别可表示为
(12)
(13)
(14)
(15)
求解多维分配算法(S≥3)的复杂度会随着点迹-航迹关联数的增加呈现指数增长的规律,计算量成为分配算法可行性的重要制约因素。因此,本文采用文献[28]中的方法,将点迹-航迹的多维分配问题转化为松弛的二维子问题进行处理,再采用广义拍卖算法[29]求得最优解。
跟踪滤波算法:当关联的测量值为第1类空量测值(ms=01)时,采用常规的IMM与扩展卡尔曼滤波(extended Kalman filter, EKF)算法。当关联的量测值为第2类空量测值(ms=02)时,采用文献[30]中的IMM与多普勒盲区粒子滤波(Doppler blind zone particle filter, BDPF)算法。
航迹撤销准则:采用文献[31]中多普勒盲区条件下的改进航迹撤销准则(改进的基于航迹似然比和序列概率比检验的撤销准则)。
进行性能比较的跟踪关联方法:① 基于单类空量测+IMM-EKF算法的二维分配方法(简称为单类空量测二维分配方法);② 基于两类空量测+IMM-BDPF算法的二维分配方法(简称为两类空量测二维分配方法);③ 基于两类空量测+IMM-BDPF的三维分配方法(简称为两类空量测三维分配方法,即本文所提方法)。
运动模型设置:将IMM模型集设为一个恒速度(constant velocity, CV)模型和一个恒转速(constant turn, CT)模型。其中,CT模型的状态方程为
(16)
式中:w(k)为过程噪声;ω为目标角速度。ω为正值时,目标做逆时针圆周运动,ω为负值时,目标做顺时针圆周运动;T为雷达采样间隔。假设目标圆周运动的法向加速度a和线速度V已知,则角速度ω可表示为
(17)
V可以用目标做CT运动的前一时刻的x方向和y方向速度表示:
(18)
参数设置如下。
(2) 目标运动状态参数:目标初始状态x(0)=(-53 km,0 m/s,180 km,-220 m/s)T,雷达状态为(0 km,0 m/s,0 km,0 m/s)T。目标运动过程分为3次匀速运动和两次恒速转弯运动,其中3次匀速运动的持续时间分别为250 s、250 s和350 s,两次恒速转弯运动穿插在3次匀速运动之间,法向加速度分别为0.3g和-0.12g,其中g为重力加速度。正加速度表示目标逆时针转弯,负加速度表示目标顺时针转弯,持续时间分别为120 s和160 s。
这里比较3种方法的跟踪性能,采用位置估计均方根误差(x方向)来衡量,可表示如下:
(19)
图3 目标真实运动航迹和雷达多普勒盲区分布Fig.3 Real motion trajectory of target and distribution of radar Doppler blind zone
由图3可以看到,雷达有3段多普勒盲区,雷达在盲区内丢失的点迹数目分别为14、7和 13。采用3种方法对目标进行跟踪,在运用文献[31]的改进航迹撤销准则后,目标航迹没有发生中断。因此,在航迹不中断的条件下,比较3种算法的滤波估计性能。在不同虚警率条件下经过1 000次蒙特卡罗仿真得到的3种算法的位置估计均方根误差如图4所示。
图4 3种算法下x方向的位置滤波误差Fig.4 x-position filtering errors in three algorithms
由图4可以看到,盲区内的状态估计误差由大到小依次为单类空量测二维分配法,两类空量测二维分配法和两类空量测三维分配法,这是由于两类空量测法充分考虑了多普勒盲区对目标状态的约束信息,减小了多普勒盲区内的状态估计误差,三维分配算法较二维分配算法的位置跟踪误差有一定程度的减小。其中,在高虚警率情况下,误差减小更为明显,这是因为三维分配利用了多帧量测的信息,杂波环境下量测与航迹的误关联概率更小。表1为3种方法在不同虚警率下的平均位置估计均方根误差(多普勒盲区内)。由图3可以看到,在高虚警率和低虚警率条件下,两类空量测三维分配法的盲区内目标位置估计误差明显小于单类空量测法和两类空量测二维分配法,且在虚警率越高的环境下跟踪精度改善越大。
表1 两种算法的多普勒盲区平均位置滤波误差Table 1 Mean position filtering errors of two algorithms in Doppler blind zone m
在多目标杂波环境下,有4批飞机目标编队航行,在航线过程中进行了数次交叉飞行,并数次进入雷达的多普勒盲区。目标1、目标2、目标3和目标4的初始位置分别为(-50 km, 140 km),(-53 km, 140 km),(-56 km, 140 km)和(-59 km, 140 km),4个目标匀速运动持续70 s,速度均为(0 m/s, -220 m/s)。随后,4个目标完成90°的匀速转弯机动,再进行匀速运动。其中,目标2和目标4做加速度为2g的匀速转弯机动,目标1和目标3做加速度为0.8g的匀速转弯机动。图5为多目标的运动轨迹和多普勒盲区分布情况。
图5 多目标真实运动航迹和雷达多普勒盲区分布Fig.5 Real motion trajectory of multi-target and distribution of radar Doppler blind zone
从图5可以看到,目标2、目标4两次进入雷达多普勒盲区,目标1、目标3进入雷达多普勒盲区3次。这里,目标检测概率Pd=0.9,雷达虚警概率Pfa=0.001。跟踪性能的评价指标选用航迹中断率RTBR,可表示为
RTBR=Nbreak/(Nblind·M)
(20)
式中:Nblind为4批目标经历的多普勒盲区的总次数(Nblind=10),M为蒙特卡罗仿真次数;Nbreak为M次蒙特卡罗仿真条件下统计的目标航迹中断总数。
下面比较3种分配方法的目标跟踪关联性能,采用单次蒙特卡罗仿真条件下的计算机运行时间表征算法的复杂度。假设在10 000次蒙特卡罗仿真条件下,统计每次仿真中的航迹中断数量Nb。3种方法的Nb次中断出现频数、航迹中断率和计算时间如表2所示。从表2中的仿真结果可以看到,两类空量测二维分配法由于利用了多普勒盲区引起丢点信息,其航迹中断率RTBR较单类空量测二维分配法有了明显降低,而两类空量测三维分配算法由于利用了连续多帧量测信息,即充分运用了航迹的演变信息,其航迹中断率RTBR较二维分配法有了进一步降低,但其计算复杂度明显高于二维分配算法。
表2 3种算法下的航迹中断率和计算时间的比较(10 000次蒙特卡罗仿真条件下的平均仿真时间)Table 2 Track breakage rate and computing time of three algorithms (mean simulation time of 10 000 Monte Carlo simulations)
本文针对多普勒盲区条件下的机载预警雷达多目标跟踪问题展开了研究。基于边跟踪边关联的思想,在两类空量测二维分配算法的基础上,提出了基于多维分配算法的点迹-航迹关联方法。由于多维分配算法考虑了连续多帧量测信息,充分利用了目标航迹演变信息,进一步提升了杂波环境下目标的跟踪精度,同时在编队目标交叉飞行的复杂多目标场景中,进一步降低了航迹中断率。