具体实施方式

下面结合具体实施例对本发明做进一步详细的描述，但本发明的实施方式不限于此。

实施例一

首先，介绍一下箱粒子(Box Particle，BP)滤波和泊松多伯努利混合(PoissonMulti-Bernoulli Mixture，PMBM)滤波方法。

BP滤波是种广义的粒子滤波算法，其结合了区间分析与序贯蒙特卡洛方法，使用一种具有可控大小且最大误差已知的箱粒子来代替粒子滤波中所用的点粒子并使用误差界限模型来描述系统的随机不确定性。

PMBM滤波方法是一种RFS框架下的多目标跟踪方法，该方法将多目标的状态建模为Poisson随机集和多Bernoulli随机集的混合，即Poisson和MBM两部分，其中Poisson部分用来表示所有未被探测到的目标，MBM部分则用来处理所有的数据关联假设。

本实施例提供的基于BP-PMBM滤波算法的多目标跟踪方法是在PMBM滤波算法的基础上，将其推广到标签随机有限集(Labeled RFS，LRFS)中以实现对目标航迹的跟踪，并提出算法的箱粒子实现方式，即BP-PMBM滤波算法，其滤波结构如下：

算法可分为Poisson与MBM两部分，其中，Poisson部分可用其强度描述，可表示为一组加权箱粒子集表示状态区间(箱)，为其权值，N为箱粒子的数量。MBM部分为多个带标签的Bernoulli分量的混合，每个Bernoulli分量表示一种目标—量测关联假设，可由参数集{l,r,w,p}描述，其中l表示标签，r和w分别为Bernoulli分量的存在概率和权重，p表示Bernoulli分量的概率密度，可以由一组加权箱粒子集表示。

具体地，请参见图1，图1是本发明实施例提供的一种基于BP-PMBM滤波算法的多目标跟踪方法流程示意图，包括：

(1)初始化算法；

具体地，PMBM概率密度包括Poisson与MBM两部分，Poisson部分的初始强度可由一组加权箱粒子表示，为一组从初始状态采样所得的箱粒子，N₀为初始时刻箱粒子的数量，每个箱粒子的初始权值为1/N₀。

由于初始时刻尚未产生Bernoulli分量，因此MBM的初始概率密度为空集。

(2)对目标状态进行预测；

(21)对所述目标状态中的Poisson分量进行预测，得到所述Poisson分量的预测强度并用箱粒子集表示；

进一步地，步骤(21)包括：

(21-1)设k-1时刻Poisson分量的强度由箱粒子集表示为其中，表示第i个箱粒子在k-1时刻的状态，表示第i个箱粒子在k-1时刻的重要性权值，表示在k-1时刻所有箱粒子的数量；

(21-2)相应地，k时刻所述Poisson分量的预测强度箱粒子集为其中，表示预测箱粒子的数量。

进一步地，所述Poisson分量的预测强度箱粒子集包括存活目标的预测箱粒子集和新生目标箱粒子集，其中，

所述存活目标的预测箱粒子表示为其中，表示第i个存活目标的预测箱粒子的状态，表示第i个存活目标的预测箱粒子的重要性权值，且

式中，p_s表示目标的存活概率，[ω]表示过程噪声区间，[f_k|k-1]表示状态转移函数的包含函数；

所述新生目标箱粒子集表示为其中，B_k表示新生箱粒子的数量，表示第i个新生目标箱粒子的状态，表示第i个新生目标箱粒子的重要性权值。

则所述Poisson分量的预测强度箱粒子集表示为：

(22)对所述目标状态中的MBM分量进行预测，得到所述MBM分量的预测参数集；

进一步地，步骤(22)包括：

(22-1)设k-1时刻MBM分量由参数集表示为其中，h表示Bernoulli分量的索引，H_k-1表示所有Bernoulli分量的数量，表示Bernoulli分量的标签，和分别表示Bernoulli分量的存在概率和权值，表示Bernoulli分量的概率密度；

具体地，可由一组加权箱粒子集描述。

(22-2)相应地，k时刻MBM分量预测参数集为其中，其中各项参数计算如下：

H_k|k-1＝H_k-1 (8)

进一步地，预测概率密度可由一组加权箱粒子集表示，其中为预测箱粒子的数量，其值与相同。箱粒子的权值与状态分别为：

(3)更新所述目标状态；

(31)更新所述Poisson分量，根据所述Poisson分量的预测强度得到所述Poisson分量的后验强度；

具体地，由步骤(21)已得到，k时刻Poisson分量的预测强度由箱粒子集表示，则k时刻Poisson分量的后验强度可表示为箱粒子集其中，

(32)更新所述MBM分量，根据所述MBM分量的预测参数集和量测信息得到更新后的MBM分量参数集；

在本实施例中，由步骤(22)已得到，k时刻预测MBM分量可表示为参数集的形式，其中Bernoulli分量的概率密度可由一组加权箱粒子集表示，为箱粒子的数量。

进一步地，步骤(32)包括：

(32-1)设k时刻量测随机集为Z＝{[z₁],…,[z_m]}；将所述MBM分量更新后的参数集表示为其中，所述更新后的MBM分量包括三种单目标假设，分别为潜在目标首次被检测到的假设、存活目标的漏检假设以及存活目标与每个量测匹配所形成的假设；

具体地，每个量测都可以被认为是来源于一个首次被检测到的潜在目标，这种假设的数量与量测数m相同，每个存活目标在当前时刻都有可能漏检，漏检假设的数量与预测单目标假设数H_k|k-1相同，每个存活目标在当前时刻可能和任何一个量测匹配，这种假设的数量为预测假设数与量测数的乘积H_k|k-1×m，更新后假设的数量H_k＝m+H_k|k-1+H_k|k-1×m，可由加权箱粒子集表示。

(32-2)根据所述量测随机集得到所述潜在目标首次被检测到所形成的Bernoulli分量的参数集；

具体地，对于任意一量测，其可以被认为是一个首次被检测到的潜在目标，会形成一个Bernoulli分量，设该分量为(l_n,r_n,w_n,p_n)，p_n可由加权箱粒子集描述。

在本实施例中，对于潜在目标首次被检测到所形成的Bernoulli分量，如果Poisson分量的预测强度中所有箱粒子关于某一量测[z]的似然和高于某一预设门限，则为中所有似然不为0的每个箱粒子的权值为：

其中，N为似然不为0的箱粒子数。

Bernoulli分量的存在概率和权值分别为：

标签为l_n＝(k,j)，k为当前时刻，j为量测的索引。

若关于量测[z]的似然和未达到预设门限，则r_n＝0，w_n＝c(z)，l_n＝0。

(32-3)建立漏检假设并计算所述存活目标的漏检假设的参数集；

具体地，对于每个存活目标的单目标假设首先建立一个漏检假设(l_mis,r_mis,w_mis,p_mis)，该假设的概率密度可由一组加权箱粒子集表示，其中每个箱粒子的状态与预测Bernoulli概率密度的加权箱粒子集中的箱粒子状态相同，权值为：

假设的存在概率和权值分别为：

标签为：

(32-4)根据所述量测随机集得到所述存活目标与每个量测匹配所形成的假设的参数集。

具体地，遍历所有的量测，对于每个量测[z]，都将其与进行匹配从而形成一个新的假设(l_det,r_det,w_det,p_det)，描述该假设概率密度的箱粒子集中每个箱粒子的状态与预测Bernoulli概率密度中的箱粒子状态相同，权值为：

假设的存在概率和权值分别为：

r_det＝1 (22)

标签为：

(4)对所述Poisson分量和所述MBM分量进行修剪；

由于随着时间的递推，Poisson分量中的箱粒子数和MBM分量中的假设的数量会越来越多，滤波器的运行速率也会不断降低，因此需要对Poisson分量和MBM分量进行修剪。

具体地，可以通过设定两个门限T_P和T_B分别对应于Poisson分量和MBM分量，对于Poisson分量，裁减掉箱粒子权值低于预设门限值的箱粒子，对于MBM分量，则裁减掉其中Bernoulli权值低于门限值的Bernoulli项。

(5)对所述Poisson分量和MBM分量进行箱粒子重采样；

具体地，采用随机子划分策略，根据每个箱粒子的权值确定采样次数，假设某箱粒子的重采样次数为n，则该箱粒子在重采样后会被划分为n个子区间，最后将这些子区间作为重采样之后的箱粒子。

(6)对所述目标状态进行估计，得到全局目标状态估计值。

(61)将所述单目标假设的权值之积作为全局假设的权重，并找出所有全局假设中权值最大的全局假设；

具体地，将全局假设j中所有的单目标假设的权值之积即作为全局假设j的权重，找出权值最大的全局假设，设其为j^*，即

(62)从所述权值最大的全局假设中筛选出其中存在概率大于某一预设门限的单目标假设，得到目标的状态估计区间；

具体地，遍历全局假设j^*中所有的单目标假设并筛选出其中存在概率大于某一预设门限T的单目标假设，这些单目标假设的概率密度可由粒子集表示，则目标的状态估计区间可由每个假设的加权粒子和来表示，即

(63)根据所述目标的状态估计区间得到目标状态估计值。

具体地，目标的状态估计值可表示为：

其中，mid(·)表示区间的中心点。

本实施例提供的基于BP-PMBM滤波算法的多目标跟踪方法在PMBM滤波算法的基础上，将其推广到标签随机有限集(Labeled RFS，LRFS)中以实现对目标航迹的跟踪，并提出算法的箱粒子实现方式，具有跟踪精度高、运算速度快、可区分航迹等优点。

实施例二

下面结合MATLAB仿真对本发明做出进一步说明。

仿真实验一：线性场景下BP-PMBM的滤波性能表现

假设在区域大小为[-250m,250m]×[-250m,250]m的具有随机噪声的二维模拟监控区域内有6个目标在先后50个观测时刻内做匀速转弯运动，目标状态为其中(x,y)表示目标的位置，为目标在x方向和y方向上的速度，表示目标的角速度。表1所示为6个目标的初始状态和存活时间。

表1目标初始状态与存活时间

假设新生目标的概率密度服从高斯分布，即p_b＝N(x；m_b,P_b)，目标的新生概率为r_b＝0.01，新生目标初始状态分别为 方差为P_B＝diag([10,10,10,10,3π/180])²。假设采样周期T＝1s，目标的状态方程和量测方程为：

x_k＝Fx_k-1+Gv_k (27)

z_k＝Hx_k+w_k (28)

其中，

设量测区间长度为Δ＝[15m,15m]^T，其噪声协方差和过程协方差分别为R＝diag([1.5²,1.5²])m²和其中σ_ω＝3m/s²，σ_u＝π/180rad/s。杂波的数量服从均值为r＝10的Poisson分布，其在监控区域内均匀分布。目标存活概率和检测概率分别为p_s＝0.99，p_d＝0.98，OSPA参数为p＝1，c＝300。6个目标的真实运动轨迹和量测和跟踪结果如图2～3所示。

其中，图2是本发明实施例提供的线性场景下目标运动轨迹与量测图，图2中实线所示为6个目标的真实运动轨迹，方形区域则表示区间量测。图3是本发明实施例提供的线性场景下BP-PMBM滤波器一次蒙特卡洛的跟踪效果图，图3中不同目标的航迹使用不同符号标记，从图3中可以看出，所提算法可以准确跟踪多目标的状态与航迹。

请参见图4和图5，图4是本发明实施例提供的线性场景下100次蒙特卡洛仿真的平均目标数估计；图5是本发明实施例提供的线性场景下100次蒙特卡洛仿真的平均OSPA距离。从图4中可以看出，滤波器可以较为准确的估计出目标的数量，且具有较好的跟踪精度。由于滤波器均基于多假设跟踪思想，当目标消失时，消失目标假设分量的权值变化存在一定的延迟，因此在图4中，在24、27、30、35、37时刻，当有目标消失时，真实目标数减少，估计目标数会延迟一个时刻减少，并且在图5中对应时刻的OSPA距离会出现尖峰。

实验中单次跟踪用时约35s，箱粒子数目为40，每个新生目标的箱粒子数目为5。

仿真实验二：非线性场景下SMC-PMBM的滤波性能表现

本实验在非线性雷达观测系统下对BP-PMBM滤波算法进行MATLAB仿真并评估及跟踪性能及表现。假设在50个观测时刻内，区域大小为[-250m,250m]×[-250m,250m]的模拟监控区域内先后新生6个目标做匀速转弯运动，目标状态为表2所示为6个目标的初始状态与起始、终止时刻。

表2目标初始状态与起止时刻

目标的状态方程与非线性量测方程为：

x_k＝Fx_k-1+Gv_k (29)

其中，

传感器采样周期为T＝1s，过程噪声协方差为其中σ_ω＝3m/s²，σ_u＝4π/180rad/s；量测噪声协方差为其中σ_α＝π/(4×180)rad，σ_ρ＝1m，量测区间为Δ＝[Δα,Δρ]^T，其中Δα＝6π/180rad，Δρ＝20m。目标的新生概率为r_b＝0.01，新生目标的初始状态分别为 方差为P_b＝diag([10,10,10,10,3π/180])²。

杂波的数量服从均值为r＝10的Poisson分布且杂波均匀分布在监控区域内。目标存活概率为p_s＝0.99，检测概率为p_d＝0.98，OSPA参数设置为c＝300,p＝1。

请参见图6和图7，图6是本发明实施例提供的非线性场景下观测区域内6个目标的真实运动轨迹以及杂波的量测分布状况，图6中实线为目标的真实运动轨迹，方形区域则表示目标与杂波的区间量测。图7是本发明实施例提供的非线性场景下BP-PMBM滤波算法的一次蒙特卡洛跟踪表现，从图7中可以看出，滤波器基本可以有效跟踪多目标的状态并区分出不同目标的航迹。请参见图8和图9，图8是本发明实施例提供的非线性场景下BP-PMBM滤波算法在平均杂波数r＝10时100次蒙特卡洛仿真的平均目标数估计；图9是本发明实施例提供的非线性场景下BP-PMBM滤波算法在平均杂波数r＝10时100次蒙特卡洛仿真的平均OSPA距离。

从图8和图9中可以看出，BP-PMBM滤波算法的数目估计较为准确且具有较好的OSPA距离表现。

以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明，不能认定本发明的具体实施只局限于这些说明。对于本发明所属技术领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干简单推演或替换，都应当视为属于本发明的保护范围。