CN105182291B - 自适应目标新生强度的phd平滑器的多目标跟踪方法 - Google Patents
自适应目标新生强度的phd平滑器的多目标跟踪方法 Download PDFInfo
- Publication number
- CN105182291B CN105182291B CN201510531102.0A CN201510531102A CN105182291B CN 105182291 B CN105182291 B CN 105182291B CN 201510531102 A CN201510531102 A CN 201510531102A CN 105182291 B CN105182291 B CN 105182291B
- Authority
- CN
- China
- Prior art keywords
- mrow
- target
- msub
- phd
- newborn
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Expired - Fee Related
Links
- 238000000034 method Methods 0.000 title claims abstract description 25
- 230000003044 adaptive effect Effects 0.000 title claims abstract description 16
- 238000005259 measurement Methods 0.000 claims abstract description 25
- 238000001914 filtration Methods 0.000 claims abstract description 17
- 230000003111 delayed effect Effects 0.000 claims abstract description 11
- 238000009966 trimming Methods 0.000 claims abstract description 4
- 230000004083 survival effect Effects 0.000 claims description 26
- 238000001514 detection method Methods 0.000 claims description 19
- 239000000203 mixture Substances 0.000 claims description 9
- 230000000007 visual effect Effects 0.000 claims description 5
- 230000004069 differentiation Effects 0.000 claims description 2
- 244000061408 Eugenia caryophyllata Species 0.000 claims 1
- 235000016639 Syzygium aromaticum Nutrition 0.000 claims 1
- 238000012795 verification Methods 0.000 abstract description 3
- 230000009467 reduction Effects 0.000 abstract description 2
- 230000006870 function Effects 0.000 description 6
- 239000011159 matrix material Substances 0.000 description 6
- 101001086191 Borrelia burgdorferi Outer surface protein A Proteins 0.000 description 5
- 238000000342 Monte Carlo simulation Methods 0.000 description 4
- 230000000694 effects Effects 0.000 description 4
- 238000005516 engineering process Methods 0.000 description 3
- 238000007689 inspection Methods 0.000 description 3
- 238000012545 processing Methods 0.000 description 3
- 238000004458 analytical method Methods 0.000 description 2
- 238000004364 calculation method Methods 0.000 description 2
- 238000012790 confirmation Methods 0.000 description 2
- 230000000977 initiatory effect Effects 0.000 description 2
- 238000004088 simulation Methods 0.000 description 2
- 230000007704 transition Effects 0.000 description 2
- OVIFHFRAOKPVPC-CJMGQXRASA-N (2s)-2-[[(2s)-1-[(2s)-1-[(2s)-2-[[(2s)-2-[[(2s)-2-[[(2s)-2-[[(2s)-6-amino-2-[[(2s)-6-amino-2-[[2-[[(2s)-2-aminopropanoyl]amino]acetyl]amino]hexanoyl]amino]hexanoyl]amino]-3-carboxypropanoyl]amino]-3-carboxypropanoyl]amino]-3-carboxypropanoyl]amino]-3-carb Chemical compound C[C@H](N)C(=O)NCC(=O)N[C@@H](CCCCN)C(=O)N[C@@H](CCCCN)C(=O)N[C@@H](CC(O)=O)C(=O)N[C@@H](CC(O)=O)C(=O)N[C@@H](CC(O)=O)C(=O)N[C@@H](CC(O)=O)C(=O)N1CCC[C@H]1C(=O)N1[C@H](C(=O)N[C@@H](CCC(O)=O)C(O)=O)CCC1 OVIFHFRAOKPVPC-CJMGQXRASA-N 0.000 description 1
- 230000007423 decrease Effects 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 238000009499 grossing Methods 0.000 description 1
- 238000012544 monitoring process Methods 0.000 description 1
- 239000002245 particle Substances 0.000 description 1
- 230000008569 process Effects 0.000 description 1
- 230000000644 propagated effect Effects 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 239000007787 solid Substances 0.000 description 1
- 239000011800 void material Substances 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
Landscapes
- Engineering & Computer Science (AREA)
- Computer Networks & Wireless Communication (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Radar, Positioning & Navigation (AREA)
- Remote Sensing (AREA)
- Radar Systems Or Details Thereof (AREA)
Abstract
本发明公开了一种自适应目标新生强度的PHD平滑器的多目标跟踪方法,主要解决在杂波环境下,自适应目标新生强度的PHD滤波器在目标新生时刻存在目标确认的滞后现象,并给出其线性高斯条件下的实现形式。通过前向滤波和后向平滑,能够准确估计目标数目和状态,减小目标新生时确认滞后对航迹形成带来的影响。其步骤包括目标新生率估计、前向滤波和后向平滑,即首先根据先验杂波数均值来估计k时刻的新生率;其次,采用k时刻的量测对目标进行预测和更新来完成前向滤波;然后,用滞后的L时刻量测对滤波结果进行后向平滑;最后,通过修剪合并和状态提取完成跟踪结果的输出。
Description
技术领域
本发明属于雷达信号处理技术领域,涉及多目标跟踪。具体地说是一种自适应目标新生强度的概率假设密度平滑器(λ-ATBI-PHD Smoother)的多目标跟踪方法,可用于杂波环境下的火控、航管等探测系统中。
背景技术
无论现代的防御探测还是空中交通管制ATC(Air Traffic Control)系统中,多目标跟踪MTT(Mutiple Target Tracking)都是其中的关键技术,也是历来受到最多关注的方向之一。多目标环境下的跟踪问题,存在以下难点:(1)每一时刻都可能存在目标的出现、衍生以及消失,使得目标数目处在一个不断变化的过程中;(2)量测信息不确定,若对漏检、虚警等问题处理不当,将会极大地影响跟踪效果;(3)进行跟踪滤波的同时,还需要完成数据关联。
因此,在复杂环境下,存在的虚警、漏警等量测信息的不确定性和目标生灭带来的数目变化给多目标的稳定跟踪造成很大困难,一直是该领域的研究热点和难点。
传统的多目标跟踪需要进行数据关联,在目标和量测数目较多时计算量会急剧增大。近年来,Mahler等学者提出基于随机有限集(Random Finite Set,RFS)理论的多目标跟踪方法,并推导了PHD(Probability Hypothesis Density,概率假设密度)、CPHD(Cardinalised Probability Hypothesis Density,势概率假设密度)、多伯努利等实现形式。该理论框架下,多目标的量测和状态分别被视为一个随机有限集,避免了复杂的数据关联,因此受到广泛的关注。
PHD滤波是通过递推目标状态的后验概率密度的一阶矩,可从中来估计多目标的状态和个数。PHD的实现方式有两种:高斯混合(Gaussian Mixture,GM)和粒子。在线性高斯条件下,高斯混合PHD(GM-PHD)利用加权的高斯混合分布来拟合目标的后验概率密度函数,状态提取不再需要进行聚类,可以用较小的计算量完成滤波估计。
在传统的PHD和CPHD滤波中,认为新生目标多出现在感兴趣的特定位置(例如:机场,视场的边界处等),新生目标的初始强度是作为先验已知量。而实际应用中,新生目标的初始强度是很难获得的,这给工程应用造成了困难。
B.Ristic等人于2012年提出了基于量测驱动的自适应目标新生强度概率假设密度滤波器(Adaptive Target Birth Intensity Probability Hypothesis DensityFilter,ATBI-PHD Filter),即ATBI-PHD滤波器,能够从量测中估计目标新生强度并同时对新生目标和存活目标进行递推滤波,摆脱对新生强度先验的依赖,并将ATBI方法和多个固定点检测新生的方法进行了对比。Liang Ma等人于2014年针对杂波环境下新生目标的快速生成,提出对PHD滤波器的目标新生率的估计方法,得到λ-ATBI-PHD Filter,即λ-ATBI-PHD滤波器。
但该方法在杂波环境下,自适应目标新生强度的PHD滤波器在目标新生时刻存在目标确认的滞后现象,即对新生目标出现时刻的估计存在延迟,对后续完整的目标航迹形成和关联造成影响。
发明内容
本发明为解决上述现有方法的问题,提出一种自适应目标新生强度的PHD平滑器的多目标跟踪方法,即λ-ATBI-PHD平滑器的多目标跟踪方法。本发明的主要方法:首先,根据先验杂波数均值来估计目标的新生率,能够在得到特定杂波数均值的情况下,尽快检测到新生目标;其次,通过后向平滑能够更加准确的估计目标出现时刻,减小目标新生时刻的确认滞后情况对后续航迹关联的影响,使得目标数目的估计更加准确。
本发明实现上述目的的技术方法包括如下步骤:
1)对k-1时刻跟踪滤波器中存活目标PHD的高斯混合形式进行初始化,设定第i个高斯项均值ms,k-1 (i),高斯项滤波协方差Ps,k-1/k-1 (i)和高斯项权值ws,k-1 (i)的初始值,其中i为高斯项的标号,Dk-1/k-1(y,β)表示全体目标PHD,Dk-1/k-1(y,0)表示存活目标PHD,Dk-1/k-1(y,1)表示新生目标PHD,y为目标状态,N(ms,k-1 (i);Ps,k-1/k-1 (i))表示均值为ms,k-1 (i),方差为Ps,k-1/k-1 (i)的高斯项;
2)根据先验杂波信息对新生目标强度估计:
2a)根据k时刻传感器的测量方差Σk/k-1、杂波检测概率pD (c)、目标检测概率pD (t)、目标存活概率pS (t)、先验杂波数均值N(c)、目标新生权值门限Te和监控区域体积VS信息来估计k时刻的目标新生率λb,k;
2b)根据估计的目标新生率λb,k和量测数据zk,把新生目标的PHD中的高斯项加入到k-1时刻跟踪结果的全体目标的高斯项中,完成对监测区域内新生目标强度的检测,式中wb,k (i)=λb,k,wb,k (i),mb,k (i),Pb,k (i)分别对应新生目标高斯项的权值,均值,方差;
3)前向滤波,包括PHD预测和PHD更新:
3a)PHD预测:对k-1时刻的全体目标PHD向k时刻进行预测,其中全体目标包括存活目标和新生目标;
3b)PHD更新:利用k时刻的量测数据zk,对预测PHD中的新生目标Dk/k-1(y,1)和存活目标Dk/k-1(y,0)分别进行更新,得到前向滤波结果其中wf,k,mf,k,Pf,k分别对应前向滤波后的高斯项权值、均值、方差;
4)后向平滑:利用滞后的L时刻的量测信息zL来平滑前向滤波后的高斯项,得到后向平滑结果其中wk/L,mk/L,Pk/L分别对应后向平滑后的高斯项权值,均值,方差,当L=k+1时为一步后向平滑;
5)删剪合并高斯项:对平滑后权值小于经验门限Tprun的高斯项进行删剪,对均值mk (i)之间的距离小于门限Umerg的高斯项进行合并,得到删剪合并后的状态估计结果其中wk (i),mk (i),Pk (i)分别对应删剪合并后的高斯项权值,均值,方差;
6)估计全体目标数目:对修剪合并后的PHD进行权值求和得到全体目标数目估计
7)输出最终的状态估计和目标数估计结果得到多目标跟踪的结果。
本发明与现有技术相比具有以下优点:
1.针对目标新生强度未知条件下,ATBI-PHD滤波器和λ-ATBI-PHD滤波器在目标新生时刻,对目标数的估计和目标确认存在滞后的这一问题提出了一种自适应目标新生强度的PHD平滑器,即λ-ATBI-PHD平滑器,分别对新生目标和存活目标进行前向滤波和后向平滑,给出了高斯混合条件下的实现形式。对于杂波环境下的多目标跟踪问题,由于PHD平滑器能够利用更多的滞后时间测量,因此相比于PHD滤波器,它在杂波密度较大以及传感器检测概率较低时可以有效提高多目标的个数和状态估计精度。
2.引入目标新生率估计方法,根据先验杂波强度来估计目标新生率,能够在得到特定杂波数均值的情况下,尽快检测到新生目标。
仿真结果表明,相较于ATBI-PHD滤波器、λ-ATBI-PHD滤波器,本发明的λ-ATBI-PHD平滑器能够更好的估计目标状态,更加准确的估计目标出现时刻,可以减小目标新生时确认的滞后对航迹生成的影响并准确估计目标数目。火控和航管等探测系统中常常需要获得准确的航迹信息,本发明的λ-ATBI-PHD平滑器在杂波环境下对于保证多目标跟踪精度和形成正确航迹方面具有积极的意义。
附图说明
图1是本发明方法的流程图
图2是传感器获得的量测数据图
图3是多目标跟踪的目标真实航迹
图4是ATBI-PHD滤波器经过50次蒙特卡洛仿真的跟踪结果
图5是λ-ATBI-PHD滤波器经过50次蒙特卡洛仿真的跟踪结果
图6是λ-ATBI-PHD平滑器经过50次蒙特卡洛仿真的跟踪结果
图7是多目标跟踪的目标数估计图
图8是多目标跟踪的OSPA距离误差图
具体实施方式
依据附图,对本发明的技术方案作具体说明。
本发明使用的量测和状态模型如下:
PHD滤波是将所有量测和目标状态看作两个随机有限集,分别为
Zk={zk,1,........,zk,m}∈F(Ζ) (1)
Ξk=Sk(Xk-1)∪Γk(Xk-1)∈F(χ) (2)
其中Zk为k时刻量测集合,Ξk为全体目标,Sk为存活目标,Γk为新生目标,F(Z)为量测随机有限集,F(χ)为目标状态的随机有限集,Xk-1为k-1时刻的目标状态。
在传感器每一时刻获得的量测是来自于目标或杂波,而全体目标的量测来源又可分为新生目标和存活目标。自适应目标新生强度方法是把目标的状态空间分为存活目标和新生目标,根据每一时刻的量测信息反演出可能的新生目标位置,同时对全体目标进行递推滤波,使得跟踪方法不再依赖已知目标新生强度的先验信息。
本发明中使用的平滑器分为前向滤波、后向平滑两个步骤。在前向滤波时,后验密度是贝叶斯递推的前向传播;在后向平滑时,利用滞后的L时刻的信息来更新k时刻的状态,完成多目标的后向平滑递推,其中k<L。
如附图1的流程图所示,本发明实现步骤如下:
1)对k-1时刻跟踪滤波器中存活目标PHD的高斯混合形式 进行初始化,设定第i个高斯项均值ms,k-1 (i),高斯项滤波协方差Ps,k-1/k-1 (i)和高斯项权值ws,k-1 (i)的初始值,其中i为高斯项的标号,Dk-1/k-1(y,β)表示全体目标PHD,Dk-1/k-1(y,0)表示存活目标PHD,Dk-1/k-1(y,1)表示新生目标PHD,y为目标状态,N(ms,k-1 (i);Ps,k-1/k-1 (i))表示均值为ms,k-1 (i),方差为Ps,k-1/k-1 (i)的高斯项;
用β=0,1来分别表示目标状态x中区分新生目标和存活目标的标志,即
用y表示目标状态x中可以观测到的部分,如目标位置、速度等;
2)根据先验杂波信息对新生目标强度估计:
2a)根据k时刻传感器的测量方差Σk/k-1、杂波检测概率pD (c)、目标检测概率pD (t)、目标存活概率pS (t)、先验杂波数均值N(c)、目标新生权值门限Te和监控区域体积VS信息来估计k时刻的目标新生率λb,k;
在杂波条件下,第j个量测的来源是杂波或目标,故PHD中第i个高斯项的预测值为:
其中pD,k (c),pD,k (t)分别为杂波和目标的检测概率,λb,k为目标新生率;
可以看出,预测的PHD中分为杂波的PHD、新生目标的PHD、存活目标的PHD三部分,pD,k (c)N(c)/VS为杂波部分,λb,k/VS为新生目标部分,pD,k (t)wk/k-1 iqk i,j为存活目标部分,qk (i,j)=N(zk j;Hmk/k-1 (i),R+HPk/k-1 (i)HT)为似然函数;
由于B.Ristic提出的自适应目标新生强度模型是假设新生目标强度在监测视场内满足均匀分布,并没有考虑杂波对目标新生的影响。由(4)式可以看出,当杂波数目均值N(c)增加时,新生部分所占的比例就会相应的下降,即新生目标初始强度降低,会造成目标新生时刻估计的滞后;当新生目标强度部分增加时,又会低估杂波数均值;
Liang Ma提出了目标新生率的估计方法,其目标新生率估计式为:
将引入目标新生率估计的自适应目标新生强度的PHD滤波称为λ-ATBI-PHD滤波,通过对目标新生率的估计,可以根据杂波环境来设置新生目标的初始权值,能够在得到特定杂波数均值的情况下,使得新生目标在连续两帧数据出现后就可以被检测到;
2b)将估计的新生率赋值给对应的新生高斯项权值wb (i),并根据量测数据把新生高斯项的参数wb,k (i),mb,k (i),Pb,k (i)赋值给并把Dk-1/k-1(y,1)加入到原有的中,完成对监测区域内新生目标的检测;
3)前向滤波:
3a)PHD预测:对全体目标的PHD进行预测,其中全体目标分为新生和存活两部分:
Dk/k-1(y,β)为预测的全体目标PHD强度,y和y’分别表示k和k-1时刻目标的可观测状态,ps(y,β)为目标的存活概率,fk/k-1为目标从k-1到k时刻的状态转移函数;
新生目标:Dk/k-1(y,1)=γk(y,1)=Σwb,kN(mb,k;Pb,k) (7)
式中γk(y,1)表示根据k时刻的量测反演的新生目标的强度;
存活目标:
由于γk(y,0)=0,故有
Dk/k-1(y,0)=<Dk-1/k-1(y,1)+Dk-1/k-1(y,0),ps(y',β')fk/k-1(y,β|y',β')>=<ps(y',1)Σwb,k-1N(mb,k-1;Pb,k-1)+ps(y',0)Σws,k-1N(ms,k-1;Ps,k-1),fk/k-1(y,β|y',β')>=Σwb,k/k-1N(mb,k/k-1;Pb,k/k-1)+Σws,k/k-1N(ms,k/k-1;Ps,k/k-1)
(8)
(8)式中
wb,k/k-1=wb,k-1ps(y',1) (8-1)
mb,k/k-1=Fkmb,k-1 (8-2)
Pb,k/k-1=FkPb,k-1Fk T (8-3)
ws,k/k-1=ws,k-1ps(y',0) (8-4)
ms,k/k-1=Fkms,k-1 (8-5)
Ps,k/k-1=FkPs,k-1Fk T+Q (8-6)
wb,k-1,mb,k-1,Pb,k-1和ws,k-1,ms,k-1,Ps,k-1分别为k-1时刻滤波结果中PHD的新生和存活部分的高斯项权值,均值,方差;
wb,k/k-1,mb,k/k-1,Pb,k/k-1和ws,k/k-1,ms,k/k-1,Ps,k/k-1分别为k-1时刻预测PHD的新生和存活部分的高斯项权值,均值,方差,Q为状态噪声协方差,Dk-1/k-1(y,0)=Σws,k-1N(ms,k-1;Ps,k-1)为k-1时刻的存活目标PHD;
3b)PHD更新:根据k时刻的量测数据zk,对预测PHD中的新生目标Dk/k-1(y,1)和存活目标Dk/k-1(y,0)分别进行更新:
其中Dk/k(y,β)为更新后的全体目标的PHD强度,gk(z|y,β)为目标的似然函数,κk(z)=N(i)/VS为监视区域内的杂波强度,为得到的全体目标的PHD前向滤波结果,pD,k(y,β)为目标的检测概率,有
(9)式中
wm,k (i)=(1-PD,k)ws,k/k-1 (i) (9-1)
均值mm,k (i),ms,k (i),mb,k (i)的更新式都形如其中mk (i)=mm,k (i),ms,k (i),mb,k (i);
协方差Pm,k (i),Ps,k (i),Pb,k (i)的更新式都形如Pk (i)=[I-Kk (i)H]Pk/k-1 (i),其中Pk (i)=Pm,k (i),Ps,k (i),Pb,k (i);
Kk (i)=Pk/k-1 (i)HT[Sk/k-1 (i)]-1 (11)
Sk/k-1 (i)=HPk/k-1 (i)HT+R (12)
其中wm,k为更新后漏检部分权值,ws,k为更新后存活部分权值,wb,k为更新后新生部分权值,H为观测矩阵,R为观测噪声协方差,Kk (i)为卡尔曼滤波增益;
4)后向平滑:利用滞后的L时刻的量测信息zL来更新k时刻的状态,完成多目标的后向平滑递推,其中当L=k+1时为一步后向平滑;
平滑器的后向更新方程为:
其中yk为k时刻目标的可观测状态,Z1:k为1到k时刻的量测集合,pS,k+1/k为存活概率,Dk/k(yk|Z1:k)为k时刻存活目标和新生目标的PHD滤波结果,Fk+1/k(yk+1|yk)为k-1到k时刻的状态转移函数,Dk+1/k(yk+1|Z1:k)为k到k+1时刻存活目标和新生目标的PHD预测,Dk+1/L(yk+1|Z1:k)为L向k+1时刻存活目标和新生目标的后向递推结果;
上式可分解为:
Dk/L(yk|Z1:L)=Dk/k(yk|Z1:k)[1-pS,k+1/k+Bk/L(yk|Z1:L)] (14)
其中Bk/L(yk|Z1:L)为L时刻对k时刻的后向平滑算子:
Bk/L(yk|Z1:L)=ps,k+1/k(yk)<Bk+1/kLk+1(Zk+1|yk'),fk+1/k(yk'|yk)>+qs,k+1/k(yk) (15)
自适应目标新生强度的PHD平滑器的高斯混合实现如下:
(16)式中
qk+1/I(zI)=N(zI|mk+1/k+1 (i)(zI),RI/k+HI/kPk/kHI/k (T)) (16-1)
wk+1/k (i)=wk/k (i)pS,k+1/k (16-3)
mk+1/k (i)=Fk+1/kmk/k (i) (16-4)
Pk+1/k (i)=Fk+1/kPk/k (i)Fk+1/k T (16-5)
其中qk+1/I(zI)为量测zI的似然函数,Dk+1/k(yk+1|Z1:k)=wk+1/k (i)N(mk+1/k (i),Pk+1/k (i)),I表示平滑使用的滞后量测的所在时刻,zI为I时刻的量测,Fk+1/k为k到k+1时刻的状态转移矩阵,wk+1/k (i),mk+1/k (i),Pk+1/k (i)分别对应k到k+1时刻的预测高斯项权值,预测均值,预测协方差矩阵,wk/L (i),mk/L (i),Pk/L (i)分别对应L时刻对k时刻平滑结果中的高斯项权值,预测均值,预测协方差矩阵;
5)删剪合并高斯项:
根据经验门限Tprun,对平滑后PHD中权值小于门限的高斯项进行删剪,即若wk/L (i)<Tprun,则将wk/L (i),mk/L (i),Pk/L (i)对应的高斯项从中剔除,得到删剪合并后的状态估计结果其中wk (i),mk (i),Pk (i)分别对应删剪合并后的高斯项权值,均值,方差;
将高斯项均值mk/L (i)之间的距离小于门限Umerg的高斯项进行合并,合并后的权值wk (i)为合并前各高斯项权值wk/L (i)之和,即若有|mk (i)-mk (i+1)|<Umerg,则
wk (i)=wk/L (i)+wk/L (i+1) (17)
mk (i)=wk/L (i)mk/L (i)+wk/L (i+1)mk/L (i+1) (18)
6)估计目标数目:对修剪合并后的PHD进行权值求和得到目标数目估计
7)输出最终的状态估计和目标数估计结果得到多目标跟踪的结果。
下面结合附图对本发明的仿真效果做进一步的描述。
1.仿真条件:
本发明的仿真是在主频3.0GHZ的Intel(R)Pentium(R)CPU G2030、内存4.00GB的硬件环境和MATLAB R2009b的软件环境下进行的。
本实验为对比ATBI-PHD滤波器、λ-ATBI-PHD滤波器和λ-ATBI-PHD平滑器三种方法的跟踪效果,取L=k+1,即一步后向平滑,仿真场景设置如下:
整个观测区域为[-100,100]×[-100,100]㎡,采样周期为1s,观测时刻k=1~40,观测过程持续40帧,相继出现4个目标,不考虑衍生目标情况,杂波数均值λc,k=2;
假设在线性高斯条件下,目标在k时刻的状态向量为状态分别为目标的x轴坐标,x方向速度,y轴坐标,y方向速度,其运动方程为:
Xk=FkXk-1+Gkwk
其中,符号Fk代表目标匀速(CV)运动模型,Gk代表扰动矩阵,wk为状态噪声,wk服从高斯分布,σw为状态噪声根方差,设置σ1,w=0.01,σ2,w=0.01;
观测方程为:Zk=HkXk+nk
其中,观测矩阵nk为观测噪声,nk服从高斯分布,设置σn=0;
假设杂波强度在场景内服从均匀分布,杂波数目服从参数λc,k的泊松分布。目标存活概率Ps=0.9,新生目标的检测概率Pd_birth=1,存活目标的检测概率Pd=0.99;
4个目标在视场内相继出现后运动至仿真结束,新生时刻分别为:目标1在时刻k=1,目标2在时刻k=9,目标3在时刻k=13,目标4在时刻k=20;
4个目标初始状态设定:目标1为[0,0.5,10,1],目标2为[10,0.3,0,-1.1],目标3为[20,0.6,0,0.9],目标4为[30,0.7,0,-1.1];
2.仿真结果及分析
传感器量测如图2所示,目标的真实航迹如图3所示,图中“□”表示目标的起始点。
采用目标状态估计、目标个数均值和最优次模式分配(optimal subpatternassignment,OSPA)距离对各方法的跟踪性能进行评价。
其中,设定参数c=70,p=2,p为距离敏感性参数,c是一个水平调节数,用于调节集合势误差的影响,X和Z为任意两个有限集合,m,n分别对应X,Z两个集合的元素个数,且m≤n,d(xi,zπ(i))表示xi和zπ(i)之间的一阶距离,。
图4,图5,图6分别为ATBI-PHD滤波器,λ-ATBI-PHD滤波器和λ-ATBI-PHD平滑器经过50次蒙特卡洛仿真的跟踪结果。由于λ-ATBI-PHD平滑器采用了一步平滑,所以有效数据处理的时刻为k=2~39。可以看出,在有效时刻上,三者在未知新生目标强度的情况下,都能完成对多目标的跟踪。
图7为三种方法的目标数的估计与真实目标数目。通过对比可以看出,估计目标新生率的λ-ATBI-PHD滤波器比ATBI-PHD滤波器的目标数估计效果更加平稳,这是由于估计目标的新生率,可以得到在得到特定杂波数均值的情况下合适的新生目标初始权值,但是对新生目标的确认滞后也更大。ATBI-PHD滤波器和λ-ATBI-PHD滤波器在k=9,13,20的目标新生的时刻,对目标数目都会产生低估,其原因如(4)式的分析。相比之下,λ-ATBI-PHD平滑器在杂波环境下对目标数的估计更加准确,在目标新生时刻能够及时检测到目标,有利于后续的航迹起始处理。
图8为三种方法的OSPA距离对比,同样,在k=9,13,20的目标新生的时刻,ATBI-PHD滤波器和λ-ATBI-PHD滤波器的OSPA距离明显增大,但本发明的λ-ATBI-PHD平滑器的OSPA距离要明显小于ATBI-PHD滤波器和λ-ATBI-PHD滤波器,表明本发明的λ-ATBI-PHD平滑器有更好的跟踪精度。
综上所述,在杂波环境下,本发明的λ-ATBI-PHD平滑器通过后向平滑能够更准确地估计目标新生时刻,减小目标新生时刻的确认滞后情况对后续航迹形成的影响,获得更好的目标跟踪精度。火控和航管等探测系统中常常需要获得准确的航迹起始信息,本发明的λ-ATBI-PHD平滑器在杂波环境下对于保证多目标跟踪精度和形成正确航迹方面具有积极的意义。
Claims (2)
1.一种自适应目标新生强度的PHD平滑器的多目标跟踪方法,包括如下步骤:
(1)对k-1时刻跟踪滤波器中存活目标的PHD的高斯混合形式 进行初始化,设定第i个高斯项均值ms,k-1 (i),高斯项滤波协方差Ps,k-1/k-1 (i)和高斯项权值ws,k-1 (i)的初始值,其中i为高斯项的标号,Dk-1/k-1(y,β)表示全体目标的PHD,β为区分新生目标和存活目标的标志,Dk-1/k-1(y,0)表示存活目标的PHD,Dk-1/k-1(y,1)表示新生目标的PHD,y为目标的可观测状态,N(ms,k-1 (i);Ps,k-1/k-1 (i))表示均值为ms,k-1 (i),方差为Ps,k-1/k-1 (i)的高斯项;
(2)根据先验杂波信息对新生目标强度估计:
2a)根据k时刻传感器的测量方差杂波检测概率pD (c)、目标检测概率pD (t)、目标存活概率pS (t)、先验杂波数均值N(c)、目标新生权值门限Te和监控区域体积VS信息来估计k时刻的目标新生率λb,k;
2b)根据估计的目标新生率λb,k和量测数据zk,把新生目标的PHD中的高斯项加入到k-1时刻跟踪结果的全体目标的高斯项中,完成对监测区域内新生目标强度的检测,式中wb,k (i)=λb,k,wb,k,mb,k,Pb,k分别对应新生目标高斯项的权值,均值,方差;
(3)前向滤波,包括PHD预测和PHD更新:
3a)PHD预测:对k-1时刻的全体目标PHD向k时刻进行预测,其中全体目标包括存活目标和新生目标;
3b)PHD更新:利用k时刻的量测数据zk,对预测PHD中的新生目标Dk/k-1(y,1)和存活目标Dk/k-1(y,0)分别进行更新,得到前向滤波结果其中wf,k,mf,k,Pf,k分别对应前向滤波后的高斯项权值,均值,方差;
(4)后向平滑:利用滞后的L时刻的量测信息zL来平滑前向滤波后的高斯项,得到后向平滑结果其中wk/L,mk/L,Pk/L分别对应后向平滑后的高斯项的权值,均值,方差;当L=k+1时为一步后向平滑;
其中PHD平滑的计算式为:
<mrow>
<msub>
<mi>D</mi>
<mrow>
<mi>k</mi>
<mo>/</mo>
<mi>L</mi>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mrow>
<msub>
<mi>y</mi>
<mi>k</mi>
</msub>
<mo>|</mo>
<msub>
<mi>Z</mi>
<mrow>
<mn>1</mn>
<mo>:</mo>
<mi>L</mi>
</mrow>
</msub>
</mrow>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msub>
<mi>D</mi>
<mrow>
<mi>k</mi>
<mo>/</mo>
<mi>k</mi>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mrow>
<msub>
<mi>y</mi>
<mi>k</mi>
</msub>
<mo>|</mo>
<msub>
<mi>Z</mi>
<mrow>
<mn>1</mn>
<mo>:</mo>
<mi>k</mi>
</mrow>
</msub>
</mrow>
<mo>)</mo>
</mrow>
<mrow>
<mo>&lsqb;</mo>
<mrow>
<mn>1</mn>
<mo>-</mo>
<msub>
<mi>p</mi>
<mrow>
<mi>S</mi>
<mo>,</mo>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
<mo>/</mo>
<mi>k</mi>
</mrow>
</msub>
<mo>+</mo>
<mo>&Integral;</mo>
<msub>
<mi>D</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
<mo>/</mo>
<mi>L</mi>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mrow>
<msub>
<mi>y</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msub>
<mo>|</mo>
<msub>
<mi>Z</mi>
<mrow>
<mn>1</mn>
<mo>:</mo>
<mi>L</mi>
</mrow>
</msub>
</mrow>
<mo>)</mo>
</mrow>
<mfrac>
<mrow>
<msub>
<mi>F</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
<mo>/</mo>
<mi>k</mi>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mrow>
<msub>
<mi>y</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msub>
<mo>|</mo>
<msub>
<mi>y</mi>
<mi>k</mi>
</msub>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mrow>
<msub>
<mi>D</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
<mo>/</mo>
<mi>k</mi>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mrow>
<msub>
<mi>y</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msub>
<mo>|</mo>
<msub>
<mi>Z</mi>
<mrow>
<mn>1</mn>
<mo>:</mo>
<mi>k</mi>
</mrow>
</msub>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
</mfrac>
<msub>
<mi>&delta;y</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msub>
</mrow>
<mo>&rsqb;</mo>
</mrow>
</mrow>
其中yk为k时刻目标的可观测状态,Z1:k为1到k时刻的量测集合,pS,k+1/k为存活概率,Dk/k(yk|Z1:k)为k时刻存活目标和新生目标的PHD滤波结果,Fk+1/k(yk+1|yk)为k-1到k时刻的状态转移函数,Dk+1/k(yk+1|Z1:k)为k到k+1时刻存活目标和新生目标的PHD预测,Dk+1/L(yk+1|Z1:k)为L向k+1时刻存活目标和新生目标的后向递推结果;
(5)删剪合并高斯项:对平滑后权值小于经验门限Tprun的高斯项进行删剪,对均值mk (i)之间的距离小于门限Umerg的高斯项进行合并,得到删剪合并后的状态估计结果其中wk (i),mk (i),Pk (i)分别对应删剪合并后的高斯项权值,均值,方差;
(6)估计全体目标数目:对修剪合并后的PHD进行权值求和得到全体目标数目估计
(7)输出最终的状态估计和目标数估计结果得到多目标跟踪的结果。
2.根据权利要求1所述的自适应目标新生强度的PHD平滑器的多目标跟踪方法,其特征在于:步骤(2)所述的目标新生率的估计计算式为:
<mrow>
<mover>
<msub>
<mi>&lambda;</mi>
<mrow>
<mi>b</mi>
<mo>,</mo>
<mi>k</mi>
</mrow>
</msub>
<mo>&OverBar;</mo>
</mover>
<mo>=</mo>
<mfrac>
<mrow>
<mn>2</mn>
<mi>&pi;</mi>
<msqrt>
<mrow>
<mo>|</mo>
<msub>
<mi>&Sigma;</mi>
<mrow>
<mi>k</mi>
<mo>/</mo>
<mi>k</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mo>|</mo>
</mrow>
</msqrt>
<msup>
<mrow>
<mo>(</mo>
<msup>
<msub>
<mi>p</mi>
<mrow>
<mi>D</mi>
<mo>,</mo>
<mi>k</mi>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>c</mi>
<mo>)</mo>
</mrow>
</msup>
<msup>
<msub>
<mi>N</mi>
<mrow>
<mi>k</mi>
<mo>/</mo>
<mi>k</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>c</mi>
<mo>)</mo>
</mrow>
</msup>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
</mrow>
<mrow>
<msup>
<msub>
<mi>p</mi>
<mrow>
<mi>D</mi>
<mo>,</mo>
<mi>k</mi>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
</msup>
<msup>
<msub>
<mi>p</mi>
<mrow>
<mi>S</mi>
<mo>,</mo>
<mi>k</mi>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
</msup>
<msub>
<mi>V</mi>
<mi>S</mi>
</msub>
<mrow>
<mo>(</mo>
<msubsup>
<mi>T</mi>
<mi>e</mi>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msubsup>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
</mrow>
</mfrac>
</mrow>
式中∑k/k-1为传感器的测量方差,pD,k (c),pD,k (t)分别为杂波和目标的检测概率,pS,k (t)为目标的存活概率,N(c)为先验杂波数均值,VS为监控区域体积,Te为新生一个目标的新生权值门限,表示目标新生率λb,k的估计值。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201510531102.0A CN105182291B (zh) | 2015-08-26 | 2015-08-26 | 自适应目标新生强度的phd平滑器的多目标跟踪方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201510531102.0A CN105182291B (zh) | 2015-08-26 | 2015-08-26 | 自适应目标新生强度的phd平滑器的多目标跟踪方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN105182291A CN105182291A (zh) | 2015-12-23 |
CN105182291B true CN105182291B (zh) | 2017-08-25 |
Family
ID=54904492
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201510531102.0A Expired - Fee Related CN105182291B (zh) | 2015-08-26 | 2015-08-26 | 自适应目标新生强度的phd平滑器的多目标跟踪方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN105182291B (zh) |
Families Citing this family (16)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105844217A (zh) * | 2016-03-11 | 2016-08-10 | 南京航空航天大学 | 一种基于量测驱动新生目标强度估计的phd多目标跟踪方法 |
CN106054167B (zh) * | 2016-05-27 | 2018-06-26 | 西安电子科技大学 | 基于强度滤波器的多扩展目标跟踪方法 |
CN108121846B (zh) * | 2016-11-29 | 2022-01-11 | 南京航空航天大学 | 一种基于熵惩罚的em未知杂波估计的phd多目标跟踪方法 |
CN107797106A (zh) * | 2017-05-08 | 2018-03-13 | 南京航空航天大学 | 一种加速em未知杂波估计的phd多目标跟踪平滑滤波方法 |
CN108333569B (zh) * | 2018-01-19 | 2021-01-12 | 杭州电子科技大学 | 一种基于phd滤波的异步多传感器融合多目标跟踪方法 |
CN108344981B (zh) * | 2018-01-19 | 2020-08-25 | 杭州电子科技大学 | 面向杂波的多传感器异步检测tsbf多目标跟踪方法 |
CN109297478B (zh) * | 2018-09-19 | 2022-05-17 | 西安汇智信息科技有限公司 | 一种基于GM-CBMeMBer的光纤陀螺导航自适应滤波方法 |
CN109886305B (zh) * | 2019-01-23 | 2021-05-04 | 浙江大学 | 一种基于gm-phd滤波的多传感器非顺序量测异步融合方法 |
CN110376582B (zh) * | 2019-01-24 | 2022-10-04 | 西安电子科技大学 | 自适应gm-phd的机动目标跟踪方法 |
CN109946694A (zh) * | 2019-03-22 | 2019-06-28 | 哈尔滨工业大学 | 基于随机有限集的圆周sar多目标跟踪方法 |
CN110084831B (zh) * | 2019-04-23 | 2021-08-24 | 江南大学 | 基于YOLOv3多伯努利视频多目标检测跟踪方法 |
CN110297221B (zh) * | 2019-06-19 | 2021-04-06 | 西安电子科技大学 | 一种基于高斯混合模型的数据关联方法 |
CN111007487B (zh) * | 2019-12-11 | 2022-11-04 | 西安电子科技大学 | 一种基于时间反演的多基地雷达目标检测方法 |
CN111504327B (zh) * | 2020-04-30 | 2023-10-27 | 江苏理工学院 | 一种基于航迹平滑技术的广义标签多伯努利多目标跟踪方法 |
CN111736145B (zh) * | 2020-06-28 | 2022-04-19 | 电子科技大学 | 一种基于高斯混合概率假设密度滤波的多机动目标多普勒雷达跟踪方法 |
CN111856442A (zh) * | 2020-07-03 | 2020-10-30 | 哈尔滨工程大学 | 一种基于测量值驱动自适应估计新生目标强度的多目标跟踪方法 |
Family Cites Families (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US8332185B2 (en) * | 2009-08-17 | 2012-12-11 | Lockheed Martin Corporation | Method and system for calculating elementary symmetric functions of subsets of a set |
CN101894381A (zh) * | 2010-08-05 | 2010-11-24 | 上海交通大学 | 动态视频序列中多目标跟踪系统 |
US8909589B2 (en) * | 2012-09-12 | 2014-12-09 | Numerica Corp. | Methods and systems for updating a predicted location of an object in a multi-dimensional space |
CN103298156B (zh) * | 2013-06-13 | 2016-01-20 | 北京空间飞行器总体设计部 | 基于无线传感器网络的无源多目标检测跟踪方法 |
CN103310115B (zh) * | 2013-06-27 | 2015-12-23 | 西安电子科技大学 | 一种多目标跟踪的杂波估计方法 |
CN103345577B (zh) * | 2013-06-27 | 2016-05-18 | 江南大学 | 变分贝叶斯概率假设密度多目标跟踪方法 |
-
2015
- 2015-08-26 CN CN201510531102.0A patent/CN105182291B/zh not_active Expired - Fee Related
Also Published As
Publication number | Publication date |
---|---|
CN105182291A (zh) | 2015-12-23 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN105182291B (zh) | 自适应目标新生强度的phd平滑器的多目标跟踪方法 | |
CN103472445B (zh) | 一种针对多目标场景的检测跟踪一体化方法 | |
CN109886305B (zh) | 一种基于gm-phd滤波的多传感器非顺序量测异步融合方法 | |
CN110542885B (zh) | 一种复杂交通环境下的毫米波雷达目标跟踪方法 | |
CN101944234B (zh) | 特征迹驱动的多目标跟踪方法及装置 | |
CN110320512A (zh) | 一种基于带标签的gm-phd平滑滤波多目标跟踪方法 | |
CN111127523B (zh) | 基于量测迭代更新的多传感器gmphd自适应融合方法 | |
CN110596693B (zh) | 一种迭代更新的多传感器gmphd自适应融合方法 | |
CN103678949B (zh) | 基于密度分析和谱聚类的多扩展目标跟踪量测集划分方法 | |
CN105354860B (zh) | 基于箱粒子滤波的扩展目标CBMeMBer跟踪方法 | |
CN104021289B (zh) | 一种非高斯非稳态噪声建模方法 | |
CN109031229B (zh) | 一种杂波环境下目标跟踪的概率假设密度方法 | |
CN104156984A (zh) | 一种不均匀杂波环境下多目标跟踪的概率假设密度方法 | |
CN103729637A (zh) | 基于容积卡尔曼滤波的扩展目标概率假设密度滤波方法 | |
WO2019047455A1 (zh) | 一种适用于闪烁噪声的多机动目标跟踪方法及系统 | |
CN105975772B (zh) | 基于概率假设密度滤波的多目标检测前跟踪方法 | |
CN108344981A (zh) | 面向杂波的多传感器异步检测tsbf多目标跟踪方法 | |
CN107797106A (zh) | 一种加速em未知杂波估计的phd多目标跟踪平滑滤波方法 | |
CN104501812A (zh) | 基于自适应新生目标强度的滤波算法 | |
CN108717702B (zh) | 基于分段rts的概率假设密度滤波平滑方法 | |
CN107102293B (zh) | 基于滑窗累积密度估计的未知杂波无源协同定位方法 | |
CN113673565A (zh) | 多传感器gm-phd自适应序贯融合多目标跟踪方法 | |
CN104777465B (zh) | 基于b样条函数任意扩展目标形状及状态估计方法 | |
Liu et al. | Multi-sensor multi-target tracking using probability hypothesis density filter | |
CN111340853B (zh) | 基于ospa迭代的多传感器gmphd自适应融合方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant | ||
CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20170825 |
|
CF01 | Termination of patent right due to non-payment of annual fee |