一种测量数据丢失情况下的多目标跟踪方法与流程

文档序号:11134407阅读:1868来源:国知局
一种测量数据丢失情况下的多目标跟踪方法与制造工艺

本发明涉及目标跟踪技术领域,特别是一种测量数据丢失情况下的多目标跟踪方法。



背景技术:

随着非传统防务及安全挑战的不断涌现,多目标跟踪算法的研究成为一个热点。而在多目标跟踪算法中,有两类主要方法,数据关联(如PDA、JPDA)和绕过关联直接处理(如PHD、CPHD)。在目标个数较多且含杂波时,数据关联(PDA等)的计算量会非常大,不利于工程应用。而近年来,多目标跟踪研究专家Ronald P.S.Mahler教授提出了基于有限集统计学(FISST)的随机有限集(random finite set,RFS)理论,在此基础上推出概率假设密度(Probability Hypothesis Density,PHD)滤波器。PHD滤波方法将集函数积分方法变换为单个目标的积分,它首先跟踪整个目标群,随后再去检测每个变量,然而PHD滤波也存在一些问题,如漏检敏感,无数目分布等。针对这一问题,Ronald P.S.Mahler教授提出了势概率假设密度(CPHD)滤波器,相比较PHD滤波器,CPHD滤波器能在传递概率密度假设函数的基础上对目标的数目分布也进行更新,在目标的估计方面做的比PHD更好。

在实际雷达量测信息中,由于大多数采用机载脉冲多普勒雷达,而多普勒盲区不可避免的会存在并导致部分目标测量信息的丢失,从而影响滤波精度。因此,需要一个鲁棒的滤波算法,能够在目标量测数据丢失(雷达盲区)情况下保持其滤波精度。在滤波过程中,改进算法更新中的增益矩阵是一种常见的方法,通过引入比例因子调节滤波算法的增益矩阵,降低算法在目标数据丢失时给滤波精度带来的影响,以提高滤波器的鲁棒性。虽然CPHD滤波器的性能优于PHD滤波器,但是CPHD滤波器的计算量也比PHD滤波器多得多,在CPHD滤波器中,计算复杂度为O(NM3),而相比之下在PHD滤波器中的计算复杂度只有O(NM),其中N为跟踪的目标数目,M为当前观测集中的观测数目。从计算复杂度可以看出,CPHD滤波器比PHD滤波器大,而当前观测集中的观测数目M作为计算复杂度的关键部分也比跟踪的目标数目N要大。目前,高斯混合概率假设密度滤波器(GM-CPHD)进行目标跟踪时,在雷达多普勒存在盲区,且计算复杂度过大、计算效率低。



技术实现要素:

本发明的目的在于提供一种测量数据丢失情况下的多目标跟踪方法,通过引入比例因子提高滤波器的鲁棒性,通过自适应门限减小计算量。

实现本发明目的的技术解决方案为:一种测量数据丢失情况下的多目标跟踪方法,包括以下步骤:

步骤1,对于多目标跟踪,目标状态集Xk={xk,1,…,xk,m(k)},m(k)是目标状态向量个数,下标k表示k时刻;目标状态随机有限集Ξk=Sk(Xk-1)∪Nk(Xk-1),其中Sk(Xk-1)、Nk(Xk-1)分别为原保存和新产生的目标随机有限集;k时刻新生目标强度函数其中分别代表第i个新生目标的权值、均值和协方差矩阵,Jγ,k为新生目标的总数;真实目标和杂波源的新生概率假设密度为势分布为

步骤2,初始化初始目标的概率假设密度D0(x)及势分布p0(n);

步骤3,预测:对目标状态集Xk在k+1时刻的概率假设密度及势分布进行预测,k≥1,得到k+1时刻的概率假设密度Dk+1|k(x)及势分布pk+1|k(n);

步骤4,更新:对目标状态集Xk在k+1时刻的概率假设密度及势分布进行更新,k≥1,得到此时刻的概率假设密度及势分布Dk+1(x)、pk+1(n);

步骤5,修剪合并:对目标集强度函数υk+1(x)的高斯项进行修剪合并,提取目标状态估计进行性能评估;

步骤6,重复步骤3~5,对目标进行跟踪直至目标消失。

本发明与现有技术相比,其显著优点为:(1)在传递PHD函数的同时传递目标数分布,并对概率假设密度及势分布进行预测和更新,可以在杂波环境下对目标的状态和数目准确的估计;(2)比例因子的引入,可以通过调节算法的增益矩阵提高滤波器的鲁棒性,解决雷达多普勒盲区的数据丢失问题;(3)自适应门限的设定,对减小滤波器的计算量起到良好的作用,使GM-CPHD滤波器的工程应用成为可能。

附图说明

图1是本发明测量数据丢失情况下的多目标跟踪方法的流程图。

图2是机载雷达与目标相对位置示意图。

图3是目标量测信息图。

图4是本发明滤波结果图。

图5是本发明与传统方法估计的目标个数图。

图6是OSPA距离图。

图7是三个算法OSPA距离比较差值图。

图8是各目标多普勒频率变化图。

图9是算法运行耗费时间对比图。

具体实施方式

下面结合附图及具体实施例对本发明做进一步详细说明。

结合图1,本发明测量数据丢失情况下的多目标跟踪方法,具体步骤如下:

步骤1,对于多目标跟踪,它的目标状态集Xk={xk,1,…,xk,m(k)},m(k)是目标状态向量个数,下标k表示k时刻;目标状态随机有限集Ξk=Sk(Xk-1)∪Nk(Xk-1),其中Sk(Xk-1)、Nk(Xk-1)分别为原保存和新产生的目标随机有限集;k时刻新生目标强度函数其中分别代表第i个新生目标的权值、均值和协方差矩阵,Jγ,k为新生目标的总数;真实目标和杂波源的新生概率假设密度为势分布为

步骤2,初始化:主要包括初始目标的概率假设密度D0(x)及势分布p0(n),具体如下:

初始目标的概率假设密度D0(x)及势分布p0(n)的关系为:

D0(x)=n0·s0(x)

其中,s0(x)为概率密度,s0(x)峰值对应先验的目标位置;初始目标的势分布p0(n)是目标数n的概率分布,p0(n)期望值为n0,即:

在高斯势概率假设密度方法中,初始概率假设密度D0(x)符合高斯分布,D0(x)由每个目标的正态分布概率和表示;而初始势分布选择为二项分布,则:

其中,[0,L]为目标满足均匀分布的区间,n0为初始目标数的推测值,n0=Lq0,q0为二项式分布发生概率,为伯努利分布、CL,n为分布系数。

步骤3,预测:对目标状态集Xk在k+1时刻的概率假设密度及势分布进行预测,k≥1,得到k+1时刻的概率假设密度Dk+1|k(x)及势分布pk+1|k(n),具体如下:

在势概率假设密度滤波器中:目标运动是独立的、不相关的,即目标x在k时刻雷达中出现的概率为bk(x)是确定的,与目标个数、状态等无关,同理目标消失的概率也一样。

在k时刻,已知的参数有:概率假设密度Dk(x)、目标数的期望nk、势分布Pk(x),k时刻存活下来的目标状态集Xk的概率假设密度如下:

Dk+1|k(ξ)=∫ps(x')·fk+1|k(x|x')·Dk|k(x')dx'

其中:ps(x')表示目标存活概率,取0.9;fk+1|k(x|x')表示单目标马儿可夫转移密度;Dk|k(x')表示前一时刻目标状态集Xk的概率假设密度,则k时刻目标状态集Xk的概率假设密度Dk+1|k(x)=b(x)+∫ps(x')·fk+1|k(x|x')·Dk|k(x')dx',b(x)为衍生目标概率;

势分布:

其中:pk+1|k(n)为k+1时刻的预测值、j为样本数、ps为目标存活概率、Jk为k个高斯成分、为k时刻期望目标数,代表第j个目标的权值,Nmax代表势分布的最大可能数,pk(l)代表前一时刻即k时刻的目标存活概率,代表二项式系数;

目标数的期望值预测:

其中:分别代表期望的新生目标数和存活目标数。

步骤4,更新:对目标状态集Xk在k+1时刻的概率假设密度及势分布进行更新,k≥1,得到此时刻的概率假设密度及势分布Dk+1(x)、pk+1(n)。包括对真实协方差矩阵和真实偏差进行无偏转换并解耦合;分析多普勒频移背景下的传感器测量数据丢失现象,引入比例因子调节系统增益矩阵;设置自适应门限对量测集合进行简化,减小当前观测集中的观测数目M,具体如下:

在预测目标状态的强度函数和预测目标状态集的势分布pk+1|k已知,且高斯混合的情况下,可以得到CPHD滤波器的更新方程如下:

势分布更新:

目标状态的强度函数更新:

目标数更新:

其中

其中,δj(·)为均衡函数,κk(·)为杂波强度函数。

4.1)盲区内的跟踪算法

在实际雷达量测信息中,由于大多数采用机载脉冲多普勒雷达,而多普勒盲区不可避免的会存在并导致部分目标测量信息的丢失,从而影响滤波精度。因此,需要一个鲁棒的滤波算法,能够在目标量测数据丢失(雷达盲区)情况下保持其滤波精度。在滤波过程中,改进算法更新中的增益矩阵是一种常见的方法,通过引入比例因子调节滤波算法的增益矩阵,降低算法在目标数据丢失时给滤波精度带来的影响,以提高滤波器的鲁棒性。

在三维坐标系中,目标运动时的多普勒频移可以表示为:

其中Vt和Va分别是目标和载机的运动速度,φ是载机航线与雷达间的偏角,β是目标与载机间的航向偏角。

机载雷达上静止目标的多普勒频移为:

多普勒盲区为|fdt|≤Δf,其中fdt=fd-fdc。此时的多普勒盲区就为[-Δf,Δf],目标多普勒频移为:

其中,vr是目标相对于传感器的径向速度,f0是目标辐射源的发射频率,c是目标辐射源信号的传播速度,分别是目标速度矢量在x,y,z方向上的分量,分别是载机速度矢量在x,y,z方向上的分量。

基于此,本发明提出了一种鲁棒的UCM-CPHD滤波算法,具体如下:

首先,定义量测误差向量为e(k+1):

其中,y(k+1)是第k+1时刻的目标笛卡尔坐标系下的量测值,是第k+1时刻的预测值,它们是3×1的向量。

然后,在增益矩阵求取公式中添加一个n维方阵S(k)用来调节GM-CPHD滤波器的增益矩阵,可得:

E(k+1)=H(k+1)P(k+1|k)HT(k+1)+(S(k)-I)R0(k+1)+R(k+1)

其中,I是单位阵,H为观测矩阵,R0(k+1)是对角阵,其对角元素的取值为平均真实偏差R(k+1)的对角元素,即在MATLAB中计算公式如下:

R0(k+1)=diag(diag(R(k+1)))

且E(k+1)计算公式如下:

ξ是滑动窗口数目,也就是求取E(k+1)所需数据的数量。

由此可得S(k)的计算公式:

设增益调节矩阵为S*(k),它是用来调节增益矩阵,求取公式如下:

Sii(k)是S(k)的i行i列元素,即S(k)对角线上的第i个元素。

整理后,改进的R-UCM-CPHD滤波算法迭代步骤如下:

步骤一:求取去偏转换量测卡尔曼滤波状态初值和估计协方差矩阵初值P(0)。设需要滤波的目标状态变量包含目标在笛卡尔坐标系上的位置、速度和加速度,因为是三维的,故包含三个坐标轴方向的位置和速度共九个变量,组成9×1向量。将接收到的前两组目标测量值代入公式,求得前三组数据的量测噪声协方差矩阵:R(1)、R(2)和R(3),并根据下式求取卡尔曼滤波初值和P(0|0):

其中

设当接收到第二组量测数据时,目标状态变量的真实值为X,则状态量测误差为:

则估计协方差矩阵的初值P(0|0)为:

步骤二:计算目标状态预测值

步骤三:计算预测估计值协方差矩阵

P(k+1|k)=F(k)P(k|k)FT(k)+Γ(k)Q(k)ΓT(k)

步骤四:求出解耦后的量测转换的均值偏差u(k+1),协方差R(k+1)和无偏量测值Z(k+1),特别的,当目标数据丢失时,将使用预测值来代替量测值。

真实协方差矩阵和真实偏差中添加偏差补偿因子,进行无偏估计并解耦合,得到改进的和代入更新方程,求取真实协方差矩阵和真实偏差公式如下:

μ=E[υm|rmm]=[μxyz]T

分别是和的协方差,μx、μy、μz分别为偏差μ在x、y、z轴上的投影,为协方差矩阵中的量。

噪声协方差矩阵是一个3×3的对称矩阵,其中非对角线上的元素Rxy,Ryz,Rxz代表着x,y,z轴的噪声耦合相,在协方差矩阵R中,相对于主对角线上的元素Rxx,Ryy,Rzz,非对角线上的元素的影响可以忽略,此时的协方差矩阵R可以简化成对角矩阵,即Rnew=diag(Rxx,Ryy,Rzz),则解耦后的协方差矩阵Rnew的计算量是未解耦前的33/(3×13)=9倍。具体解耦细节如下:

在笛卡尔坐标系上,噪声协方差矩阵R可以描述为;

Rnew=MRMT

其中,M是转移矩阵。

此时,转换后的噪声协方差矩阵Rnew可以写成:

Rnew=diag(Rxx,Ryy,Rzz)

步骤五:求出增益调节矩阵S*(k)

步骤六:计算增益矩阵

K(k+1)=P(k+1|k)HT(k+1)[H(k+1)P(k+1|k)HT(k+1)+(S*(k)-I)R0(k+1)+R(k+1)]-1

步骤七:计算滤波估计值.

步骤八:计算滤波估计值协方差矩阵.

P(k+1|k+1)=P(k+1|k)-K(k+1)H(k+1)P(k+1|k)

步骤九:返回步骤二,并使k+1后进行迭代。

与传统的UCM-CPHD算法不同的是,R-UCM-CPHD算法主要是对计算增益矩阵的公式进行修改,加入了增益调节矩阵S*(k),当预测数据与量测数据误差很大时,此时很可能产生了奇异值甚至数据丢失。此时增益调节矩阵S*(k)通过之前的目标数据信息,对增益矩阵进行调节,减小盲区导致的数据丢失对滤波估计值和协方差的影响。这样,R-UCM-CPHD算法能够改进传统的UCM-CPHD算法在盲区导致的数据丢失时的滤波性能,具有很好的鲁棒性。

4.2)自适应门限

在CPHD滤波器中,算法的计算量主要来自每个更新周期都要计算M+1个元素的均衡函数(即前文提到的δj(·)),其计算复杂度为O(NM3),而相比之下在PHD滤波器中的计算复杂度只有O(NM),由此可知,减小M值能更有效的降低计算复杂度。

设γ为跟踪门限的大小,观测维数M,残差协方差矩阵S(k),则残差向量d(k)的范数为g(k)=dT(k)S-1(k)d(k),g(k)为服从自由度为M的分布,当g(k)≤γ时,跟踪门规则满足。

首先计算出目标量测数据落入跟踪门内的概率PG,结合跟踪门规则可得:其中跟踪门的体积为其中系数由观测空间的维数nz决定(c1=2,c2=π,c3=4π/3)。

设置自适应门限去掉所有和已知被探测目标不相关的量测数据,设γ为自适应跟踪门限的大小,观测维数M,残差协方差矩阵S(k),此时门限γ满足(nz=3):

其中,PD为检测概率,β为新源密度;根据公式,门限γ与残差协方差矩阵S(k)有关,而S(k)是一个与目标状态集有关的变量。从公式中可以看出,门限γ与残差协方差矩阵S有关,而S是一个与目标状态集有关的变量。

门限所包含的区域落入联合门限区域的量测集合就可以表示为:与传统算法相比,采用自适应门限可能会去掉所有和已知被探测目标不相关的量测数据,减小计算量。

最后更新目标状态的强度函数υk+1|k(x)及势分布pk+1|k(n),由于量测集合由Zk变化为引起杂波强度函数变化为则更新公式

杂波强度的变化直接影响到势概率假设密度滤波器中高斯分量对应权值和量测状态集的变化,而后者的大小M恰恰是CPHD滤波器计算复杂度O(NM3)的关键部分,所以通过自适应门限可以有效的减小滤波器的计算量。

步骤5,修剪合并:对目标集强度函数υk+1(x)的高斯项进行修剪合并,提取目标状态估计进行性能评估。首先修剪,该步骤需要将小于权值τ的高斯成分滤掉。

式中是大于阈值的权值,为单个高斯分量的权值,为GM-CPHD函数中单个高斯分量的均值,为目标状态协方差矩阵,x为目标状态集,为与目标状态、均值和协方差矩阵有关的函数,为k+1时刻的概率假设密度函数;

然后合并,当一些高斯成分间的距离小于阈值U时,需要将这些高斯成分合并;最后提取权值大于τ1的高斯成分。

最后状态估计,目标的状态估计是提取权值大于τ1的高斯成分,其提取公式如下:其中,为单个高斯分量的权值,为GM-CPHD函数中单个高斯分量的均值,τ1为设定的门限。

性能评价:对于多目标跟踪算法的性能评价指标,一般采用均方误差(MSE)、均方根误差(RMSE)、圆丢失概率(CPEP)、Wasserstein距离和OSPA距离等。一般OSPA距离可以通过水平参数c来很好地体现多目标跟踪算法对位置和目标数目跟踪的性能。OSPA距离定义如下:

其中,

步骤6,重复步骤3~5,对目标进行跟踪直至目标消失。

实施例1

本发明的效果可以通过以下仿真实验进一步说明:

1.仿真条件

假设目标状态其中位置单位是m,速度单位是m/s。本仿真有四个目标,各个目标运动模型是CA模型,且四个目标的初始状态如下:

Xt1=[300m,300m,40m,-20m/s,-20m/s,-1m/s,0m/s2,0m/s2,-2m/s2]T

Xt2=[40m,-300m,-300m,1m/s,20m/s,20m/s,2m/s2,0m/s2,0m/s2]T

Xt3=[200m,300m,300m,-1m/s,-20m/s,-20m/s,-2m/s2,0m/s2,0m/s2]T

Xt4=[-200m,40m,-200m,20m/s,-1m/s,20m/s,0m/s2,-2m/s2,0m/s2]T

目标的运动方程:

仿真实验假设目标1的出生时刻为第0s,死亡时刻为第7s,目标2的出生时刻为第7s,死亡时刻为25s,目标3的出生时刻为第11s,死亡时刻为37s,目标4的出生时刻为25s,死亡时刻为第40s。设雷达采样周期T=1s,径向距离量测方差方位角与高低角量测方差取检测概率PD=0.99,目标存活概率PS=0.9,合并阈值U=4,裁剪阈值τ=1e-5,状态估计阈值τ1=0.5,最大高斯数Jmax=100。

新出现的目标符合泊松过程,bk(x)=N(x,mb,Pb),

mb=[300m,300m,40m,-20m/s,-20m/s,-1m/s,0m/s2,0m/s2,-2m/s2]T,

Pb=diag(104×[5,5,5,1,1,1,0.5,0.5,0.5]),Q=diag(10-2×[1,1,1])。

2.仿真内容和结果分析

生成的目标量测信息如图3所示,量测信息中含有目标信息和杂波信息。

各算法目标位置估计与去除杂波后的量测值和真实值对比如图4所示。图中CPHD代表正常情况下的UCM-AG-CPHD,CPHD2代表加盲区情况下的UCM-AG-CPHD,CPHD3代表改进后的加盲区R-UCM-AG-CPHD。从图中可以看出,各CPHD算法可以对多个目标同时跟踪。

算法各个时刻估计的目标个数结果如图5所示,可以看出,各算法的估计精度大于90%,且盲区对目标数目估计的影响较小。

由OSPA距离结果图6可以看出,各CPHD算法不管是对目标数目的估计,还是对目标位置的估计都具有很好地性能。

三个算法OSPA距离比较如图7,从图中可以看出,R-UCM-AG-CPHD算法在多目标跟踪方面有着良好的性能,改进的算法(CPHD3)的跟踪精度优于不改进的算法(CPHD2),说明改进方法具有良好的效果,但显然会比没有盲区的算法(CPHD1)精度低。蓝色部分代表CPHD3-2,即改进的盲区算法与未改进的盲区算法的OSPA距离差,在16~18s目标1发生机动,数据丢失;在26~28s目标2数据丢失;在31~33s目标3数据丢失;在36~38s目标4数据丢失。对应图中可以看到相应的时间内改进的算法OSPA距离小于未改进的(即蓝色线条<0),仿真表明,改进的算法对目标机动产生的盲区有良好的效果。

目标多普勒频率(径向速度)随仿真次数的变化情况如图8所示。根据公式可知,目标多普勒盲区区间是[-100 100]HZ,其中各目标运动轨迹有一段落入多普勒盲区内,分别是在16~18s目标1;在26~28s目标2;在31~33s目标3;在36~38s目标4数据丢失。在这段区域内,雷达的目标数据将会受到严重影响。

通过多次实验,得到本发明改进的算法和传统UCM-PHD算法耗费时间如图9中,从可知在保证算法性能的前提下,本算法比改进前的计算量相比降低8.6%左右,新算法在计算效率上的优势十分明显,具有更好的工程应用前景。

当前第1页1 2 3 
网友询问留言 已有0条留言
  • 还没有人留言评论。精彩留言会获得点赞!
1