0 引言

群目标(如机群编队和弹道导弹群目标), 是由一些具有相似运动特性且在很长一段时间内都处于相互邻近状态的目标组成.同时, 这些目标又会随着时间不断地分裂融合, 构成复杂的运动模式, 给群目标的稳定跟踪带来许多新的挑战[.由于群目标的量测数目是时变的, 传统的多目标跟踪算法往往难以实现对群目标的准确跟踪[.文献[

近年来, 随机有限集(Random finite set, RFS)理论的提出, 为群目标跟踪算法提供了重要的研究思路.文献[[的群目标跟踪算法, 可以获得群结构的动态演化特性.文献[[对群目标状态进行估计, 获得了群内单个目标的状态与目标航迹.但需要说明的是, 文献[

文献[

基于上述问题, 本文提出一种基于均值漂移和双层群结构(Bilayer group structure, BGS)模型的群目标高斯混合概率假设密度滤波(MS-BGS-GMPHD)算法.通过引入均值漂移聚类对群目标量测集进行划分, 以解决群目标合并前和分裂后群数目的漏估问题, 同时降低运算量, 提高算法的实时性.对于群目标交叉时刻的漏估问题, 均值漂移聚类未能解决, 为此, 本文提出一种双层群结构模型, 通过搭建第2层群结构解决该问题.仿真实验表明, 本文所提出算法能够较好地解决群合并、交叉及分裂前后群数目的漏估问题和算法实时性差的问题.

1 双层群结构模型

本文算法主要针对大群目标进行跟踪, 因此, 采用对群整体进行跟踪的办法, 使用椭圆RHM对群目标量测源进行建模, 以实现对群目标质心和扩散形态的跟踪, 并称用质心和扩散形态描述的子群为第1层群结构模型.

1.1 第1层群结构模型

为实现第1层群结构模型的搭建, 本文采用椭圆RHM对群目标量测源进行建模, 同时选用质心坐标(mx, my)、椭圆长半轴a、短半轴b以及长轴与x轴的夹角ϕ(逆时针为正)这5个参数来描述椭圆.群目标状态向量由质心状态和形状参数两部分构成, k时刻群目标状态向量表示为xk=[mxk, myk, vxk, vyk, ak, bk, ϕk], 其中vxk和vyk分别为质心速度在x轴和y轴的分量.

假设每个由目标产生的量测均是由散布于群目标或扩展目标表面的量测源产生的, 对于给定的量测源ykl, 其对应目标产生的量测为

(1)

其中: l=1, 2, ..., nT，nT服从泊松分布, 为目标产生量测的数目; vk为零均值高斯白噪声, 其协方差矩阵为Rk.椭圆RHM如

图 1(Fig. 1)

图 1 椭圆RHM示意图

椭圆RHM是一种描述量测源的特殊模型, 总体而言, 量测源分布于缩小的群目标真实扩散形态上[.假设S为群目标扩散形态的边界, m为其质心, 则量测源y可描述为

(2)

其中s为缩放因子, 且s∈[0, 1].为便于描述和计算, 本文采用极坐标来描述椭圆RHM, 极坐标系下的椭圆边界表示为

(3)

k时刻量测源可表示为

(4)

(5)

其中:质心mk=[mxkmyk]T, [θkl, R(θkl; ak, bk, ϕk)]为k时刻量测源在以极点为中心的极坐标, ek为极坐标到直角坐标系的转换向量.将式(4)代入(1)中, 可得椭圆RHM下的量测方程

(6)

显然, 这里的θkl是未知的, 本文用量测到质心的向量与x轴的夹角来近似[.为了减小对状态向量xk的影响, 此处对式(6)做一简单的数学处理，得到对应的伪量测方程[如下:

(7)

1.2 第2层群结构模型

在群目标交叉时刻, 因不同群目标产生的量测集存在部分重叠, 传统的距离划分算法无法得到准确的量测划分, 导致在该时刻经常发生群目标漏估.本文在此构建第2层群结构以解决群目标交叉时刻量测划分不准确的问题.

构建k时刻的第2层群结构Gk, 用gki表示k时刻各子群目标, mki表示该时刻各子群的质心位置向量, i=1, 2, ..., N, N为子群个数.利用得到的所有子群质心构成一个质心向量集mk1, ..., mkN, 并搭建代表各子群之间关系的邻接矩阵A以表示第2层群结构信息, 同时将之初始化为0矩阵.计算各子群质心间的马氏距离di, j, 并与预设门限ε比较, 若di, j < ε, 则认为子群gki与gkj在k时刻相距较近, 群目标可能发生了交叉运动, 此时将邻接矩阵A的对应位置, 即第i行第j列与第j行第i列, 置为1.最终可以获得代表第2层群结构信息的对称邻接矩阵A.若最终的A为0矩阵, 则表示各子群相距较远; 否则, 认为子群间存在交叉运动.

以g1、g2和g3, 3个点和粗线描绘的椭圆分别表示这3个子群的质心和扩散形态, 3个箭头表示3个子群的运动方向, 2个细线描绘的椭圆代表所划分的第2层群结构.显然, g1与g2和g3相距较远, 在第2层群结构中单独成群, 而g2与g3发生了交叉运动, 相距较近, 在第2层群结构中被划分至同一个群中, 由该图搭建的代表第2层群结构信息的邻接矩阵为

图 2(Fig. 2)

图 2 双层群结构模型

由此可见, 邻接矩阵A很好地反映了各子群间的关系, 可以借此判断群目标是否存在交叉运动的情况, 这将为本文算法的量测划分步中是否进行二次划分提供了直接的依据.

2 均值漂移算法

由于群目标运动模式较为复杂, 存在群合并与分裂等典型的运动特征, 这给群目标的稳定跟踪带来了新的挑战:群合并与分裂前后群数目出现漏估.这是由于群合并与分裂前后群目标量测相距较近, 传统的距离划分算法难以得到准确的量测划分.另外, 因传统的距离划分所得量测划分数较多, 计算量较大, 严重影响了算法的实时性.本文采用均值漂移算法解决这两个问题.

均值漂移算法是从密度函数梯度的非参数估计中获得的, 非参数估计则是从量测集出发对密度函数进行估计.其中, 最常用的是核密度估计, 它根据核函数κ(x)对量测集进行计算得到密度函数.

假定k时刻群目标量测集为Zk, 核函数为κ(·), 量测zk的核密度估计为

(8)

本文采用高斯核函数, 即

(9)

其中: d为量测维数, 本文取2;h为核函数的带宽.依据上述公式可得最终的均值漂移向量为

(10)

令且, 即可求得密度函数的局部极大值点.均值漂移迭代公式为

(11)

均值漂移向量a(zk)总是指向密度大的方向, 即密度梯度增加的方向, 其收敛点为局部密度极大值点.

经均值漂移迭代后, 各量测被划分到不同的类别中, 同时得到多个收敛点, 即不同类别的中心.为获得最终的聚类结果, 需对各收敛点进行相似性度量, 并将相似性大的收敛点合并, 对应的量测类别归为一类.本文选用马氏距离[作为相似性度量的标准.

3 基于均值漂移和双层群结构模型的群目标GMPHD滤波

3.1 算法实现框图与群目标运动模型

3.1.1 算法实现框图

基于均值漂移和双层群结构模型的群目标GMPHD滤波算法的实现框图如

图 3(Fig. 3)

图 3 本文所提出群目标跟踪算法框图

3.1.2 群目标运动模型

假设群目标线性运动且相互独立, 即

(12)

其中:过程噪声ωk(i)服从高斯分布, 其均值为0、协方差矩阵为Qk; Fk为群目标状态转移矩阵.

群目标的量测方程如式(6)所示, 式(7)为其对应的伪量测方程.显然, 这里的群目标量测模型是非线性的.

另外, 需要特别说明的是, 大群目标模型与扩展目标模型基本一致.因此, 与文献[[来解决大群目标跟踪问题, 并在此基础上, 使用均值漂移算法进行量测划分, 给出结合双层群结构模型的PHD更新步骤.在本文中, 将大群目标与扩展目标统称为群目标.

3.2 算法具体流程

3.2.1 剔除杂波量测

引入均值漂移算法中的核密度估计, 使用式(9)所示高斯核函数对每个量测点进行核密度估计, 并求出最大核密度fmax和最小核密度fmin, 将两者N等分, 核密度小于fmin+β(fmax-fmin)/N的量测点作为杂波点剔除, 其中参数β一般取1即可.若杂波强度很大, 则可适当增大β.均值漂移聚类对孤立的杂波点非常敏感, 因此这一步是必须的.当然, 杂波量测的剔除也利于算法整体性能的提升.

3.2.2 均值漂移聚类

对剔除杂波量测后的量测集进行均值漂移聚类, 并依据第2层群结构更新后反馈回的标志位信息判断是否需要进行二次划分.若不需要, 则直接进行下一步; 若需要, 则对经均值漂移聚类所得划分中的每个子集Wj, 计算其中群数目最大似然估计值, 如果, 则使用k-means算法对之进行二次划分, 并在划分完毕后将从第2层群结构更新后反馈回的标志位信息清零.

子集Wj中群数目的最大似然估计计算方法如下:假设每个群目标产生的量测数服从均值为γ的泊松分布, 且各群目标产生的量测相互独立, 则子集Wj中群目标数的似然函数[为

(13)

其中pois(·)表示泊松分布.假设子集Wj中的杂波足够少, 可忽略不计, 则极大似然估计可表示为

(14)

3.2.3 预测

对于群目标的预测与传统的单量测目标的处理方法完全相同, 假设k-1时刻的后验强度具有高斯混合形式, 即

(15)

k时刻的预测强度也是高斯混合形式, 可表示为

(16)

其中: Jk|k-1为预测的高斯项的个数, wk|k-1(j)、mk|k-1(j)和Pk|k-1(j)分别为第j个高斯项预测的权值、均值和协方差.

3.2.4 更新

在单量测目标的量测更新中, 每一个量测都被用来更新每一个高斯分量, 但在群目标跟踪中, 是用划分单元作为一个整体对高斯分量进行更新的.

假设预测PHD具有高斯混合形式, 如式(16)所示.在群目标产生量测个数服从泊松分布假设的前提下, k时刻更新的PHD可表示为对漏检部分的更新和对检测到部分的更新的和的形式, 且同样具有高斯混合形式, 即

(17)

式(17)最终可化简为伪似然函数LZk(x)与预测PHD的乘积, 即

(18)

类似于预测PHD, Jk|k、wk|k(j)、mk|k(j)和Pk|k(j)分别表示更新的高斯项个数、第j个高斯项更新的权值、均值和协方差.

式(18)中的伪似然函数LZk(x)计算方法如下:

(19)

其中: γ(·)为目标产生量测的泊松率; pD(·)为目标检测概率, {1-exp[-γ(x)]}pD(x)为目标至少产生一个量测且被检测到的概率, 1-{1-exp[-γ(x)]}pD(x)为等效的漏检概率; 式(19)加号右侧表达式中p∈Zk为对量测集Zk的p种可能的划分, wp为某种划分在所有可能划分中所占的权重; W为某种划分中各非空集合, |W|为集合中量测个数, dW为该集合的权重; λk和ck(z)分别为杂波率和观测区域内杂波空间分布.

由于这里的群目标量测模型是非线性的, 在此更新步中需对其进行无迹变换(UT)以求解更新后各高斯项的均值mk|k(j)和协方差Pk|k(j)以及单量测伪似然函数ϕz(x), 具体算法可参见文献[

3.2.5 高斯项的修剪与合并

类似于多目标GMPHD滤波器, 群目标GMPHD滤波器后验强度中的高斯项也会随着时间无限增加; 因此, 需要对其进行适当的修剪与合并[.文献[T)的高斯项; 合并相近的(小于合并门限U)高斯项.若最终所得高斯项个数仍然大于最大高斯项的个数Jmax, 则选择权重最大的Jmax个高斯项作为最终修剪合并后的项.

3.2.6 群目标状态提取

由于每个高斯项的均值对应于后验强度的一个局部极值点, 权重表示该高斯项对目标个数期望的贡献, 群目标状态的估计可以根据高斯混合项的权重获得[.首先, 筛选出权重大于0.5的高斯项, 并对权值四舍五入, 值为1;然后, 将该高斯项的均值作为一个群目标状态提取出来, 若值为2, 则说明该高斯项的均值对应两个群的状态, 记录两次; 以此类推, 直至提取出所有群目标的状态.

3.2.7 第2层群结构更新

从式(6)中提取出的群目标状态中再次提取各群的质心信息, 依据状态方程对之进行预测, 得到下一时刻群目标质心的大致位置, 并使用第2层群结构模型对预测的群目标质心状态进行建模和更新, 可以得到反映下一时刻各群目标之间关系的邻接矩阵.根据第2层群结构模型, 若邻接矩阵为0矩阵, 则说明下一时刻各群目标相距较远, 不存在群目标交叉的情况; 若邻接矩阵为非0矩阵, 则认为下一时刻存在群目标交叉运动的情况, 需考虑进行二次划分.

为方便描述与计算, 这里设置一标志位flag, 并将之初始化为0.当最终所得的邻接矩阵为0矩阵时, 标志位flag的值不变, 仍为0;若邻接矩阵为非0矩阵, 则将标志位flag置为1.当标志位为1时, 需进行二次划分.最后, 将此标志位flag的值反馈回量测划分步.

4 仿真实验与分析

为了验证本文所提出基于均值漂移和双层群结构模型的群目标GMPHD滤波算法的性能, 设计了包含群目标融合、交叉及分裂等复杂运动模式的仿真场景进行实验, 并与文献[[, 作为滤波性能的评价指标, 所用PC机平台为Inter(R) Pentium(R) CPU G2030 @ 3.00 GHz, RAM 4.00 GB, matlab(R2010b).

在二维平面上一大小为[-500, 500] m×[-500, 500] m的观测区域内考察多群目标.区域中共有4个不同大小的群目标, 群3在仿真开始9 s后融合进群2, 使得群2扩散形态变大; 第20 s时, 群1和群2交叉; 第32 s后, 群4从群1中分裂出来, 同时, 群1扩散形态变小.各群不同时刻的形状参数如

表 1(Table 1)

表 1 各群目标不同时刻参数表[a/m, b/m, ϕ/(°)]

time/s

group1

group2

group3

group4

1~9

10, 8, -10

8, 6, 45

8, 6, -20

-

9~32

10, 8, -10

9, 7, 10

-

-

32~40

8, 6, -15

9, 7, 10

-

8, 6, -25

表 1 各群目标不同时刻参数表[a/m, b/m, ϕ/(°)]

图 4(Fig. 4)

图 4 群目标仿真场景与单次仿真结果

假设群目标在观测区域内线性运动, 将之建模为CV模型, 各群产生量测的泊松率均为15, 即γ(x)在目标扩散区域内均匀分布; 设缩放因子s服从N(0.7, 0.08);杂波率为5, 在观测区域内均匀分布; 群目标存活概率ps=0.99;检测概率pd=0.99, 采样时间间隔Ts=1 s, OSPA距离参数设置为p=2, c=100;高斯项的修剪门限T=1e-5, 合并门限U=4, 最大高斯项数设置为Jmax=100.

状态方程(12)和量测方程(6)的参数设置为

其中: Id表示d×d维的单位阵, 符号"⊗"表示矩阵的直积.

设观测区域中新生群目标个数Jb=3, 权重为wb, k(j), 方差为Pb, k(j)=diag(5, 5, 10, 10, 2, 2, 0.1), 3个新生群目标状态分别初始化为

群目标分裂强度参数设置为

群目标产生的量测信息与真实质心如

单次蒙特卡洛(MC)仿真结果如

本文采用100次MC仿真的结果来验证算法的跟踪性能.群目标势估计结果如

图 5(Fig. 5)

图 5 100次蒙特卡洛仿真结果

由

两种算法平均运行时间如

表 2(Table 2)

表 2 两种算法平均运行时间对比

算法

RHM-GMPHD

MS-BGS-GMPHD

运行时间/s

247.387 0

143.574 8

表 2 两种算法平均运行时间对比

由

5 结论

本文提出了一种基于均值漂移和双层群结构模型的群目标GMPHD滤波算法.通过引入均值漂移算法, 解决了群目标合并及分裂前后由于量测划分不准确所导致的群目标数漏估问题; 同时, 也大幅减少了量测划分数, 降低了运算量, 提升了算法的实时性.但均值漂移算法未能解决群目标交叉时刻群数目的漏估问题, 因此, 提出一种双层群结构模型, 通过搭建第2层群结构并将其更新后的信息反馈回量测划分步以判断是否需要进行二次划分.双层群结构模型的提出很好地解决了群目标交叉时刻群数目的漏估问题.蒙特卡洛仿真表明, 与传统的基于距离划分的群目标RHM-GMPHD算法相比, 本文所提出算法获得了更高的实时性, 较好地解决了群目标合并、交叉及分裂前后群数目的漏估问题, 实现了对多群目标质心和扩散形态的稳定跟踪.

参考文献

[1] Mihaylova L, Carmi A Y, Septier F. Overview of Bayesian sequential Mote Carlo methods for group and extended object tracking[J]. Digital Signal Processing, 2014, 25(C): 1-16.

[2] Koch J W. Bayesian approach to extended object and cluster tracking using random matrices[J]. IEEE Trans on Aerospace and Electronic Systems, 2008, 44(3): 1042-1059. DOI:10.1109/TAES.2008.4655362

[3] Clark D, Godsill S. Group target tracking with the Gaussian mixture probability hypothesis density filter[C]. The 3rd Int Conf on Intelligent Sensors, Sensor Networks and Information. Melbourne: IEEE, 2007: 149-154.

[4] 李振兴, 刘进忙, 李松, 等. 基于箱式粒子滤波的群目标跟踪算法[J]. 自动化学报, 2015, 41(4): 785-798.

(Li Z X, Liu J M, Li S, et al. Group targets tracking algorithm based on box particle filter[J]. Acta Automatica Sinica, 2015, 41(4): 785-798.)

[5] Abdallah F, Gning A, Bonnifait P. Box particle filtering for nonlinear state estimation using interval analysis[J]. Automatica, 2008, 44(3): 807-815. DOI:10.1016/j.automatica.2007.07.024

[6] Liang M, Song L P, Ji H B, et al. Box-particle CPHD filter for multi-target tracking[C]. 2015 Int Conf on Control, Automation and Information Sciences(ICCAIS). Changshu, 2015: 80-84.

[7] Zhu S J, Liu W F, Wen C L. Multiple group targets tracking using the generalized labeled multi-bernoulli filter[C]. Proc of the 35th Chinese Control Conf. Chengdu, 2016: 4871-4876.

[8] Lin S F, Ba T V, Sven E N. Measurement driven birth model for the generalized labeled multi-bernoulli filter[C]. 2016 Int Conf on Control, Automation and Information Sciences(ICCAIS). Ansan, 2016: 94-99.

[9] Ba-Ngu V, Ba-Tuong V, Hung G H. An efficient implementation of the generalized labeled multi-bernoulli filter[J]. IEEE Trans on Signal Processing, 2017, 65(8): 1975-1987. DOI:10.1109/TSP.2016.2641392

[10] 连峰, 韩崇昭, 刘伟峰, 等. 基于SMC-PHDF的部分可分辨的群目标跟踪算法[J]. 自动化学报, 2010, 36(5): 732-741.

(Lian F, Han C Z, Liu W F, et al. Tracking partly resolvable group targets using SMC-PHDF[J]. Acta Automatica Sinica, 2010, 36(5): 732-741.)

[11] 张慧, 徐晖, 王雪莹, 等. 一种基于椭圆随机超曲面模型的高斯混合PHD滤波器[J]. 光学学报, 2013, 33(9): 14-23.

(Zhang H, Xu H, Wang X Y, et al. A Gaussian mixture PHD filter for group targets tracking based on ellipse random hypersurface models[J]. Acta Optica Sinica, 2013, 33(9): 14-23.)

[12] Baum M, Noack B, Hanebeck U D. Extended object and group tracking with elliptic random hypersurface models[C]. The 13th Conf on Information Fusion. Edinburgh, 2010: 1-8.

[13] Antonio Z, Florian F, Marcus B, et al. Level-set random hypersurface models for tracking nonconvex extended objects[J]. IEEE Trans on Aerospace and Electronic Systems, 2016, 52(6): 2990-3007. DOI:10.1109/TAES.2016.130704

[14] 刘风梅, 葛洪伟, 杨金龙, 等. 基于均值漂移聚类的扩展目标量测集划分算法[J]. 计算机工程, 2014, 40(12): 182-187.

(Liu F M, Ge H W, Yang J L, et al. Extended target measurement set partition algorithm based on mean shift clustering[J]. Computer Engineering, 2014, 40(12): 182-187. DOI:10.3969/j.issn.1000-3428.2014.12.034)

[15] 李翠芸, 桂阳, 刘靳. 基于mean shift迭代的新生未知多扩展目标跟踪[J]. 控制与决策, 2017, 32(3): 521-525.

(Li C Y, Gui Y, Liu J. Unknown newly born multiple extended targets tracking based on mean shift iteration[J]. Control and Decision, 2017, 32(3): 521-525.)

[16] Zhang H, Xu H, Wang X Y, et al. A PHD filter for tracking closely spaced objects with elliptic random hypersurface models[C]. The 16th Int Conf on Information Fusion. Istanbul, 2013: 1558-1565.

[17] Granstrom K, Lundquist C, Orguner O. Extended target tracking using a Gaussian-mixture PHD filter[J]. IEEE Trans on Aerospace and Electronic Systems, 2012, 48(4): 3268-3286. DOI:10.1109/TAES.2012.6324703

[18] Vo B N, Ma W K. The Gaussian mixture probability hypothesis density filter[J]. IEEE Trans on Signal Processing, 2006, 54(11): 4091-4104. DOI:10.1109/TSP.2006.881190

[19] 刘妹琴, 兰剑. 目标跟踪前沿理论与应用[M]. 北京: 科学出版社, 2015: 159-161.

(Liu M Q, Lan J. Leading-edge theory and application of target tracking[M]. Beijing: Science Press, 2015: 159-161.)

[20] Schuhmacher D, Vo B T, Vo B N. A consistent metric for performance evaluation of multi-object filters[J]. IEEE Trans on Signal Processing, 2008, 56(8): 3447-3457. DOI:10.1109/TSP.2008.920469