CN111711432B - 一种基于ukf和pf混合滤波的目标跟踪算法 - Google Patents

一种基于ukf和pf混合滤波的目标跟踪算法 Download PDF

Info

Publication number
CN111711432B
CN111711432B CN202010547096.9A CN202010547096A CN111711432B CN 111711432 B CN111711432 B CN 111711432B CN 202010547096 A CN202010547096 A CN 202010547096A CN 111711432 B CN111711432 B CN 111711432B
Authority
CN
China
Prior art keywords
target
observation
algorithm
ukf
time
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.)
Active
Application number
CN202010547096.9A
Other languages
English (en)
Other versions
CN111711432A (zh
Inventor
张烈平
彭慧豪
唐柳
于滟琳
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Guilin University of Technology
Original Assignee
Guilin University of Technology
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Guilin University of Technology filed Critical Guilin University of Technology
Priority to CN202010547096.9A priority Critical patent/CN111711432B/zh
Publication of CN111711432A publication Critical patent/CN111711432A/zh
Application granted granted Critical
Publication of CN111711432B publication Critical patent/CN111711432B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • H—ELECTRICITY
    • H03—ELECTRONIC CIRCUITRY
    • H03H—IMPEDANCE NETWORKS, e.g. RESONANT CIRCUITS; RESONATORS
    • H03H21/00—Adaptive networks
    • H03H21/0012—Digital adaptive filters
    • H03H21/0043—Adaptive algorithms

Landscapes

  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明提供一种基于UKF和PF混合滤波的目标跟踪算法，属于目标跟踪算法技术领域。本发明是一种基于UKF和PF混合滤波的多观测站适应度加权目标跟踪算法，其中UKF和PF混合滤波算法，是先用粒子滤波对状态向量进行初始估计，这一滤波过程中，系统不受非线性系统的制约；为了消除粒子退化、奇异值对估计结果造成的影响，将上一步所求的估计结果再进行一次无迹卡尔曼滤波，从而提高估计的精度。本发明提出的UKF和PF混合滤波算法不仅能克服非线性系统的制约，也可以粒子退化、奇异值对估计结果造成的影响，滤波精度高，具有较高应用前景。相比于但观测站目标跟踪算法，有效改善了跟踪效果，跟踪精度更高。

Description

一种基于UKF和PF混合滤波的目标跟踪算法
【技术领域】
本发明涉及目标跟踪算法技术领域，具体涉及一种基于UKF和PF混合滤波的目标跟踪算法。
【背景技术】
无论是在军事领域还是民用领域，目标跟踪技术均发挥着至关重要的作用，其在敌情监控、海陆空防御、智能交通、货物仓库的物流跟踪等领域具有重要应用价值，因此一直以来都是国内外学者的研究热点。目标跟踪实际上就是应用滤波方法处理目标状态估计问题，不同的滤波方法直接影响着目标跟踪的精度，也很大程度上决定了跟踪的实时性，由于实际应用中环境和条件的复杂性，对不同滤波方法的研究与改进一直是目标跟踪研究的难点和关键问题。解决目标跟踪问题的传统方法有：卡尔曼滤波(Kalman Filter，KF)、扩展卡尔曼滤波(Extend Kalman Filter，EKF)、无迹卡尔曼滤波(Unscented Kalman Filter,UKF)等，但是对于强非线性系统，或者过程噪声、观测噪声非高斯白噪声时，上述滤波算法则无法保证稳定的跟踪性能，会造成较大跟踪误差，严重时甚至会发散。Gordon等在1993年提出一种新的理论，即重采样技术，进而粒子滤波(Particle Filter，PF)方法应运而生，并且被应用于单站纯角度目标跟踪，取得了优于卡尔曼滤波器的跟踪效果。由于PF对非线性非高斯系统具有较强的适用性，与传统方法相比其跟踪性能更加稳健，近年来，其在目标跟踪、机器人定位等领域均取得了成功的应用。
尽管基于PF的跟踪算法可以应用于非线性非高斯的系统，但是其也面临粒子退化、易受奇异值影响等问题。因此，研究卡尔曼与粒子混合滤波的目标跟踪算法，利用粒子滤波克服非线性非高斯环境问题，然后再用卡尔曼滤波方法消除奇异值的影响，对于进一步提高跟踪精度具有十分重要的意义。
随着目标机动性增强、电子对抗能力的不断提高，单个观测站在测量精度、范围、稳定性和可靠性方面存在明显不足，因此仅依靠单观测站进行目标跟踪已经不能满足跟踪精度的要求。近年来，传感器技术不断进步，形式多样的侦测手段也相继出现，人们开始探索利用多个传感器对目标进行协同跟踪，并进一步运用信息融合技术进行融合滤波处理，最终形成统一的目标状态估计。时至今日，多传感器目标跟踪问题仍然涉及多学科研究，是一个复杂的研究课题，当前还没有合适的可以应用于各种跟踪情况的通用解决方案。所以研究适合的多观测站信息融合方法，以获得更好的跟踪效果，已成为目标跟踪领域中研究和关注的重点，并具有重要的战略意义和应用前景。
【发明内容】
本发明的发明目的在于：针对上述存在的问题，提供一种基于UKF和PF混合滤波的目标跟踪算法，本发明提出的适应度加权融合方法有效改善了多观测站目标跟踪融合效果，并且基于UKF和PF混合滤波的多观测站适应度加权目标跟踪算法取得了较好的跟踪效果，对于提高目标跟踪精度的研究具有较高的参考价值。
为了实现上述目的，本发明采用的技术方案如下：
一种基于UKF和PF混合滤波的目标跟踪算法，其具体计算步骤如下：
设定目标运动的状态方程和目标的观测方程为：
xk＝Fxk-1+ωk-1 (1)
zk＝h(xk)+vk (2)
式中，
Figure GDA0002577317530000021
代表k时刻的目标状态变量；x′k、y′k表示k时刻目标在X方向、Y方向的位置；/>
Figure GDA0002577317530000022
表示k时刻目标在X方向、Y方向的速度，/>
Figure GDA0002577317530000023
表示k时刻目标在X方向、Y方向的加速度，T是采样时间；F代表状态转移矩阵；xk-1表示k-1时刻目标的状态变量；ωk-1表示k-1时刻系统的过程噪声，其是均值0，方差分别为Qk-1的高斯白噪声；代表zk为k时刻的目标观测量；vk是k时刻均值为0，方差为Rk的过观测噪声序列；对xk的滤波估计值来自于所有1～k时刻的目标观测量/>
Figure GDA0002577317530000024
h()表示线性或非线性的观测模型函数；按照以下步骤进行运算：
1.1：应用PF算法，根据已知的先验分布p(x0)采样得到N个粒子
Figure GDA0002577317530000025
1.2：根据目标运动的状态方程得先验密度
Figure GDA0002577317530000026
选用先验密度作为重要性密度函数，即/>
Figure GDA0002577317530000027
对该函数采样得到粒子/>
Figure GDA0002577317530000028
其中，
Figure GDA0002577317530000029
表示在k-1时刻第i个粒子的状态变量；/>
Figure GDA00025773175300000210
表示在k时刻第i个粒子的状态变量；
1.3：根据目标的观测方程得
Figure GDA00025773175300000211
进一步计算权值/>
Figure GDA00025773175300000212
并进行权值归一化
Figure GDA0002577317530000031
式中，
Figure GDA0002577317530000032
表示在k时刻第i个粒子的权值；/>
Figure GDA0002577317530000033
表示在k时刻第i个粒子归一化后的权值；
1.4：当有效粒子个数估计值
Figure GDA0002577317530000034
时，Nth为阙值，进行重采样，生成新的粒子并重新赋予权值/>
Figure GDA0002577317530000035
Figure GDA0002577317530000036
1.5：进而可以得到k时刻应用PF算法得到的目标状态变量估计值
Figure GDA0002577317530000037
1.6：进一步应用UKF算法进行混合滤波消除奇异值影响；根据k-1时刻的
Figure GDA0002577317530000038
利用Sigma采样点的选取规则获得2n+1个sigma点/>
Figure GDA0002577317530000039
和对应权值w0,w1,…,w2n，并根据/>
Figure GDA00025773175300000310
接着计算sigma点集的一步预测/>
Figure GDA00025773175300000311
其中，n代表随机变量xk的状态的维数；/>
1.7：根据式(5)、(6)计算目标的状态变量的预测均值
Figure GDA00025773175300000312
及协方差矩阵Pk|k-1；
Figure GDA00025773175300000314
Figure GDA00025773175300000315
式中，
Figure GDA00025773175300000316
和/>
Figure GDA00025773175300000317
为Sigma采样点的权重系数，下标m表示均值，下标c表示协方差；Qk-1为过程噪声方差；
1.8：根据计算所得目标的状态变量的预测均值
Figure GDA00025773175300000318
和协方差矩阵Pk|k-1，再次使用UT变换，产生新的sigma点/>
Figure GDA00025773175300000319
和对应权值w0,w1,…,w2n；
1.9：根据目标的观测方程进行sigma点非线性变换，根据式(7)计算观测量预测采样点
Figure GDA00025773175300000320
并根据式(8)计算系统观测量的预测均值/>
Figure GDA00025773175300000321
Figure GDA00025773175300000322
Figure GDA0002577317530000041
1.10：根据式(9)计算新息协方差矩阵Pzz，k|k-1，根据式(10)计算状态与观测间的互协方差矩阵Pxz，k|k-1，并根据式(11)计算滤波增益矩阵Kk；
Figure GDA0002577317530000042
Figure GDA0002577317530000043
Figure GDA0002577317530000044
式中，Rk为该时刻的观测噪声方差；
1.11：最后进行状态更新；根据式(12)
Figure GDA0002577317530000045
计算k时刻后验状态估计均值/>
Figure GDA0002577317530000046
根据式(13)/>
Figure GDA0002577317530000047
计算协方差矩阵Pk；
Figure GDA0002577317530000048
Figure GDA0002577317530000049
1.12：循环步骤1.2至步骤1.11，继续进行k+1时刻的混合滤波计算。
本发明中，进一步地，所述F的取值如下：
Figure GDA00025773175300000410
本发明中，进一步地，设观测站的位置为[x0,y0]，采用的传感器为距离传感器，则
Figure GDA00025773175300000411
本发明还提供一种基于上述目标跟踪算法的多观测站适应度加权目标跟踪算法，即采用a个观测站进行目标跟踪，利用不同观测站在时间、空间及功能上的互补，通过融合这些观测站的观测信息，提高目标跟踪的精度，a为正整数，所述算法包括以下步骤：
步骤1，各观测站采集观测信息，利用UKF和PF混合滤波算法进行各观测站的目标状态估计，并分别记为
Figure GDA0002577317530000051
步骤2：根据步骤1求出的混合滤波算法的目标观测量的预测均值
Figure GDA0002577317530000052
根据公式(15)计算当前k时刻第i个观测站的适应度函数fk (i)；
Figure GDA0002577317530000054
式中：exp表示以e为底的指数函数，
Figure GDA0002577317530000055
表示k时刻各观测站的实际目标观测量，
Figure GDA0002577317530000056
表示各观测站的观测预测均值，R(i)表示各观测站的观测噪声方差；
步骤3：根据fk (i)和公式(16)计算k时刻各观测站进行融合时所占的权重Wk (i)；
Figure GDA0002577317530000059
步骤4：根据Wk (i)对k时刻各观测站生成的目标状态变量估计值进行加权融合，获取相应的目标状态变量融合估计值
Figure GDA00025773175300000511
Figure GDA00025773175300000512
步骤5：置k＝k+1，循环进行步骤1到步骤4，进行下一时刻跟踪计算。
综上所述，由于采用了上述技术方案，本发明的有益效果是：
1、本发明提出的UKF和PF混合滤波算法，先用粒子滤波对状态向量进行初始估计，这一滤波过程中，系统不受非线性系统的制约；然后，为了消除粒子退化、奇异值对估计结果造成的影响，将上一步所求的估计结果再进行一次无迹卡尔曼滤波，从而提高估计的精度。相较于单独采用UKF算法或PF算法，本发明提出的UKF和PF混合滤波算法不仅能克服非线性系统的制约，也可以粒子退化、奇异值对估计结果造成的影响，滤波精度高，具有较高应用前景。
2、本发明提出的基于UKF和PF混合滤波算法的多观测站目标跟踪算法在整个运动过程中均能较好地跟踪真实轨迹，而单观测站目标跟踪算法在运动中期、后期与目标真实轨迹有一定偏差。因此，基于本发明提出的的多观测站目标跟踪算法有效改善了跟踪效果，相较于单观测站跟踪时，跟踪精度更高。
【附图说明】
图1是UKF和PF混合滤波算法流程图；
图2是目标真实轨迹与UKF算法、PF算法、UKF和PF混合滤波算法的跟踪轨迹曲线；
图3是UKF算法、PF算法、UKF和PF混合滤波算法仿真计算的X轴方向位置RMSE曲线；
图4是UKF算法、PF算法、UKF和PF混合滤波算法仿真计算的Y轴方向位置RMSE曲线；
图5是目标真实轨迹与UKF的单观测站目标跟踪算法、UKF的多观测站目标跟踪算法的跟踪轨迹曲线；
图6是UKF的单观测站目标跟踪算法、UKF的多观测站目标跟踪算法仿真计算的X轴方向位置RMSE曲线；
图7是UKF的单观测站目标跟踪算法、UKF的多观测站目标跟踪算法仿真计算的Y轴方向位置RMSE曲线；
图8是目标真实轨迹与UKF的多观测站适应度加权目标跟踪算法、UKF和PF混合滤波的多观测站目标跟踪算法(采用等权重加权)、UKF和PF混合滤波的多观测站适应度加权目标跟踪算法的跟踪轨迹；
图9是UKF的多观测站适应度加权目标跟踪算法、UKF和PF混合滤波的多观测站目标跟踪算法(采用等权重加权)、UKF和PF混合滤波的多观测站适应度加权目标跟踪算法仿真计算的X轴方向位置RMSE曲线；
图10是UKF的多观测站适应度加权目标跟踪算法、UKF和PF混合滤波的多观测站目标跟踪算法(采用等权重加权)、UKF和PF混合滤波的多观测站适应度加权目标跟踪算法仿真计算的Y轴方向位置RMSE曲线。
【具体实施方式】
为了更清楚地表达本发明，以下通过具体实施例对本发明作进一步说明。
一、技术方案介绍
本发明所提供的一种基于UKF和PF混合滤波的目标跟踪算法，是采用a个部署在不同方位的的观测站进行多观测站目标跟踪，各个观测站的加权因子值分别为W(i)，(i＝1,2,…,a)，a为正整数，各个观测站进行观测时彼此独立，即他们对同一目标进行跟踪时，观测误差是独立的，且假设各观测站局部状态估计的误差也是不相关的，本发明利用不同观测站在时间、空间及功能上的互补，通过融合这些观测站的观测信息，可提高目标跟踪的精度，各观测站均采用UKF和PF混合滤波算法，则基于UKF和PF混合滤波的多观测站适应度加权目标跟踪算法步骤为：
步骤1，各观测站采集观测信息，利用UKF和PF混合滤波算法计算各观测站的目标状态变量估计值，其流程图参见图1，首先采用PF算法得到k时刻的状态估计
Figure GDA0002577317530000071
然后应用UKF算法进行混合滤波消除奇异值影响，计算目标状态变量的预测均值/>
Figure GDA0002577317530000072
目标观测量的预测均值/>
Figure GDA0002577317530000073
最后得到k时刻a个观测站的目标状态变量估计值分别为/>
Figure GDA0002577317530000074
具体步骤如下：
设定目标运动的状态方程和目标的观测方程为：
xk＝Fxk-1+ωk-1 (1)
zk＝h(xk)+vk (2)
式中，
Figure GDA0002577317530000075
代表k时刻的目标状态变量；x′k、y′k表示k时刻目标在X方向、Y方向的位置；/>
Figure GDA0002577317530000076
表示k时刻目标在X方向、Y方向的速度，/>
Figure GDA0002577317530000077
表示k时刻目标在X方向、Y方向的加速度，T是采样时间；F代表状态转移矩阵，其取值如式(17)；xk-1表示k-1时刻目标的状态变量；ωk-1表示k-1时刻系统的过程噪声，其是均值0，方差分别为Qk-1的高斯白噪声；代表zk为k时刻的目标观测量；vk是k时刻均值为0，方差为Rk的过观测噪声序列；对xk的滤波估计值来自于所有1～k时刻的目标观测量/>
Figure GDA0002577317530000078
h()表示非线性的观测模型函数；设观测站的位置为[x0,y0]，采用的传感器为距离传感器，则
Figure GDA0002577317530000079
按照以下步骤进行运算：
1.1：应用PF算法，根据已知的先验分布p(x0)采样得到N个粒子
Figure GDA00025773175300000710
1.2：根据目标运动的状态方程得先验密度
Figure GDA00025773175300000711
选用先验密度作为重要性密度函数，即/>
Figure GDA00025773175300000712
对该函数采样得到粒子/>
Figure GDA00025773175300000713
其中，
Figure GDA00025773175300000714
表示在k-1时刻第i个粒子的状态变量；/>
Figure GDA00025773175300000715
表示在k时刻第i个粒子的状态变量；
1.3：根据目标的观测方程得
Figure GDA0002577317530000081
进一步计算权值/>
Figure GDA0002577317530000082
并进行权值归一化
Figure GDA0002577317530000083
式中，
Figure GDA0002577317530000084
表示在k时刻第i个粒子的权值；/>
Figure GDA0002577317530000085
表示在k时刻第i个粒子归一化后的权值；
1.4：当有效粒子个数估计值
Figure GDA0002577317530000086
时，Nth为阙值，进行重采样，生成新的粒子并重新赋予权值/>
Figure GDA0002577317530000087
Figure GDA0002577317530000088
1.5：进而可以得到k时刻应用PF算法得到的目标状态变量估计值
Figure GDA0002577317530000089
1.6：进一步应用UKF算法进行混合滤波消除奇异值影响；根据k-1时刻的
Figure GDA00025773175300000810
利用式(18)和(19)获得2n+1个sigma点/>
Figure GDA00025773175300000811
和对应权值w0,w1,…,w2n；
Figure GDA00025773175300000812
式中，
Figure GDA00025773175300000813
n代表随机变量x的状态维数，/>
Figure GDA00025773175300000814
表示矩阵平方根的第i列；/>
Figure GDA00025773175300000815
和Px为x的均值和方差；
这些Sigma采样点对应的权重系数为：
Figure GDA00025773175300000816
式中，下标m表示均值，下标c表示协方差，其中：λ＝α2(n+κ)-n表示缩放参数，决定了sigma采样点与均值
Figure GDA00025773175300000818
的距离，可以降低总的预测误差。式中参数α确定了/>
Figure GDA00025773175300000817
周围sigma采样点的分布状态，通常设置为较小的正值，且10-4≤α＜1；κ表示辅助缩放参数，通常设置为0或3-n；参数β是调整参数(对于高斯分布，β＝2是最优的，而如果状态变量是单变量，则最优的选择是β＝0)。
这里是把sigma点集看做一个向量
Figure GDA0002577317530000091
接着根据式(20)计算sigma点集的一步预测
Figure GDA0002577317530000092
Figure GDA0002577317530000093
1.7：根据式(5)、(6)计算目标的状态变量的预测均值
Figure GDA0002577317530000094
及协方差矩阵Pk|k-1；
Figure GDA0002577317530000096
Figure GDA0002577317530000097
式中，
Figure GDA0002577317530000098
和/>
Figure GDA0002577317530000099
为Sigma采样点的权重系数，下标m表示均值，下标c表示协方差；Qk-1为过程噪声方差；
1.8：根据计算所得目标的状态变量的预测均值
Figure GDA00025773175300000910
和协方差矩阵Pk|k-1，再次使用UT变换，产生新的sigma点/>
Figure GDA00025773175300000911
和对应权值w0,w1,…,w2n；
1.9：根据目标的观测方程进行sigma点非线性变换，根据式(7)计算观测量预测采样点
Figure GDA00025773175300000912
并根据式(8)计算系统观测量的预测均值/>
Figure GDA00025773175300000913
Figure GDA00025773175300000914
Figure GDA00025773175300000915
1.10：根据式(9)计算新息协方差矩阵Pzz，k|k-1，根据式(10)计算状态与观测间的互协方差矩阵Pxz，k|k-1，并根据式(11)计算滤波增益矩阵Kk；
Figure GDA00025773175300000916
Figure GDA0002577317530000101
Figure GDA0002577317530000102
式中，Rk为该时刻的观测噪声方差；
1.11：最后进行状态更新；根据式(12)
Figure GDA0002577317530000103
计算k时刻后验状态估计均值/>
Figure GDA0002577317530000104
根据式(13)/>
Figure GDA0002577317530000105
计算协方差矩阵Pk；
Figure GDA0002577317530000106
Figure GDA0002577317530000107
1.12：循环步骤1.2至步骤1.11，继续进行k+1时刻的混合滤波计算。
步骤2：根据步骤1求出的混合滤波算法的目标观测量的预测均值
Figure GDA0002577317530000108
根据公式(15)计算当前k时刻第i个观测站的适应度函数fk (i)；
Figure GDA00025773175300001010
式中：exp表示以e为底的指数函数，
Figure GDA00025773175300001011
表示k时刻各观测站的实际目标观测量，
Figure GDA00025773175300001012
表示各观测站的观测预测均值，R(i)表示各观测站的观测噪声方差；
步骤3：根据fk (i)和公式(16)计算k时刻各观测站进行融合时所占的权重Wk (i)；
Figure GDA00025773175300001015
步骤4：根据Wk (i)对k时刻各观测站生成的目标状态变量估计值进行加权融合，获取相应的目标状态变量融合估计值
Figure GDA00025773175300001017
Figure GDA00025773175300001018
步骤5：置k＝k+1，循环进行步骤1到步骤4，进行下一时刻跟踪计算。
二、实验仿真
为了讨论本发明在多观测站目标跟踪中的应用效果，进行以下仿真实验。
本次实验的硬件环境为英特尔Intel Corei-7200U处理器、4GB内存，在Windows10环境下，仿真软件采用MATLAB(R2018a)版本。假设目标运动状态方程和观测方程与上文相同，假设跟踪区域内部署有15个观测站，观测站1～观测站15的位置分别为[1200m,5100m]、[1100m,5000m]、[1400m,5150m]、[1700m,5050m]、[1100m,5150m]、[1800m,5150m]、[1300m,5375m]、[1600m,5300m]、[1000m,5250m]、[1350m,5250m]、[1725m,5250m]、[1400m,5350m]、[1250m,5200m]、[1150m,5350m]、[1600m,5100m]，各观测站通过距离传感器获得目标与观测站之间的距离信息，实际观测过程中具有加性观测噪声vk，且vkN(0,R)，并假设各观测站噪声方差R的取值为：观测站1～观测站5：R＝52；观测站6～观测站10：R＝82；观测站11～观测站15：R＝102。各参数设置如下：目标运动的初始值x0＝[1000m 5000m 10m/s 50m/s 2m/s2 -4m/s2]T，采样时间T＝0.5s，时间步长为50，UT变换的相关系数设置为：α＝0.01，κ＝0，β＝2，以及采样粒子数目N＝500。
三、结果分析
1、为了验证算法的性能，采用UKF算法、PF算法、UKF和PF混合滤波算法做单观测站混合滤波目标跟踪实验仿真。图2给出了三种算法的目标跟踪效果图，基于UKF和PF混合滤波的目标跟踪算法效果最好，其跟踪轨迹最接近目标真实轨迹；PF算法跟踪效果次之，后期其跟踪轨迹稍微偏离目标真实轨迹；UKF算法跟踪效果最差，其在坐标[1300m 5300m]附近偏离目标真实估计较大，且在跟踪后期与目标运动真实轨迹也偏离较大。
为了直观展现不同算法的跟踪效果，这里采用目标位置均方根误差以及目标速度均方根误差进行比较，图3到图4给出了经过100次蒙特卡洛仿真后计算所得的RMSE。从图3到图4可以看出：虽然三种算法均能大致估计出目标的位置和速度变化趋势，但本发明提出的UKF和PF混合滤波的跟踪性能最优，位置和速度RMSE均最小；PF滤波跟踪效果次之，虽然不是每一次的估计结果都比UKF准确，但整体上位置和速度RMSE比UKF算法小；UKF算法跟踪效果在三种算法中排最后。
2、为验证单观测站和多观测站的跟踪精度，将基于UKF的单观测站目标跟踪算法、基于UKF的多观测站目标跟踪算法做分析对比，且假设基于UKF的单观测站目标跟踪算法只采用观测站1进行跟踪。
图5是两种算法的跟踪轨迹比较，可以看出基于UKF的多观测站目标跟踪算法在整个运动过程中均能较好地跟踪真实轨迹，而基于UKF的单观测站目标跟踪算法在运动中期、后期与目标真实轨迹有一定偏差。因此，基于UKF的多观测站目标跟踪算法有效改善了跟踪效果，相较于单观测站跟踪时，跟踪精度更高。
为了更直观分析两种算法的跟踪效果，还需要进一步比较两种算法的跟踪误差。采用均方根误差作为评价指标，对比X轴方向和Y轴方向的RMSE，图6和图7分别是两种算法经过100次蒙特卡洛实验后，X轴方向和Y轴方向的RMSE。可以看出，基于UKF的多观测站目标跟踪算法性能更优，其在X轴方向和Y轴方向的RMSE最大值不超过6米，而采用单观测站跟踪时，X轴方向和Y轴方向的RMSE最大值达到16米。这是因为当只用单观测站进行跟踪时，目标在运动过程中逐渐远离观测站，随着观测距离的增加，获得的观测数据变差，因此跟踪精度下降。但是，多观测站目标跟踪算法由于融合了15个观测站的数据，对观测站的观测位置不敏感，能够在目标运动过程获得较优的跟踪精度。
3、为了验证不同融合方法的性能，将基于UKF的多观测站适应度加权目标跟踪算法、基于UKF和PF混合滤波的多观测站目标跟踪算法(采用等权重加权)、以及本发明提出的基于UKF和PF混合滤波的多观测站适应度加权目标跟踪算法做分析对比。图8给出了三种算法的目标跟踪效果图，可以看出：在上述三种算法中，基于UKF和PF混合滤波的多观测站适应度加权目标跟踪算法效果最好、基于UKF和PF混合滤波的多观测站目标跟踪算法(采用等权重加权)效果次之，基于UKF的多观测站适应度加权目标跟踪算法效果最差。
为了直观展现不同算法的跟踪效果，这里同样采用目标位置均方根误差以及目标速度均方根误差进行比较，图9到图10给出了经过100次蒙特卡洛仿真后计算所得的RMSE。
分析图9到图10可以得出进一步结论：(1)从X轴方向、Y轴方向位置RMSE曲线图中可以看出：当同样采用混合滤波算法时，本发明的多观测站适应度加权与对比组的等权重加权相比，虽然有一定的改进效果，但改进效果不大，这是因为本文提出的混合滤波算法对跟踪效果有显著改进，再结合多观测站进行跟踪，仅是这两种方案已经很大程度上提高了跟踪精度，所以即使进一步应用适用度加权进行改进，对跟踪效果的提高也贡献不大了。(2)基于UKF的多观测站适应度加权目标跟踪算法在三种算法中效果最差，说明本文提出的混合滤波算法、适应度加权这两种方案相比较而言，混合滤波算法对提高跟踪精度贡献更大，效果更明显。
上述说明是针对本发明较佳可行实施例的详细说明，但实施例并非用以限定本发明的专利申请范围，凡本发明所提示的技术精神下所完成的同等变化或修饰变更，均应属于本发明所涵盖专利范围。

Claims (4)

1.一种基于UKF和PF混合滤波的目标跟踪算法，其特征在于，具体计算步骤如下：
设定目标运动的状态方程和目标的观测方程为
xk＝Fxk-1+ωk-1
zk＝h(xk)+vk
式中，
Figure FDA0002541092890000011
代表k时刻的目标状态变量；x′k、y′k表示k时刻目标在X方向、Y方向的位置；/>
Figure FDA0002541092890000012
表示k时刻目标在X方向、Y方向的速度，/>
Figure FDA0002541092890000013
表示k时刻目标在X方向、Y方向的加速度，T是采样时间；F代表状态转移矩阵；xk-1表示k-1时刻目标的状态变量；ωk-1表示k-1时刻系统的过程噪声，其是均值0，方差分别为Qk-1的高斯白噪声；zk为k时刻的目标观测量；vk是k时刻均值为0，方差为Rk的过程观测噪声序列；对xk的滤波估计值来自于所有1～k时刻的目标观测量/>
Figure FDA0002541092890000014
h()表示线性或非线性的观测模型函数；
1.1：应用PF算法，根据已知的先验分布p(x0)采样得到N个粒子
Figure FDA0002541092890000015
1.2：根据目标运动的状态方程得先验密度
Figure FDA0002541092890000016
选用先验密度作为重要性密度函数，即/>
Figure FDA0002541092890000017
对该函数采样得到粒子/>
Figure FDA0002541092890000018
其中，/>
Figure FDA0002541092890000019
表示在k-1时刻第i个粒子的状态变量；/>
Figure FDA00025410928900000110
表示在k时刻第i个粒子的状态变量；
1.3：根据目标的观测方程得
Figure FDA00025410928900000111
进一步计算权值/>
Figure FDA00025410928900000112
并进行权值归一化/>
Figure FDA00025410928900000113
Figure FDA00025410928900000114
表示在k时刻第i个粒子的权值；/>
Figure FDA00025410928900000115
表示在k时刻第i个粒子归一化后的权值；
1.4：当有效粒子个数估计值
Figure FDA00025410928900000116
时，Nth为阙值，进行重采样，生成新的粒子并重新赋予权值/>
Figure FDA00025410928900000117
1.5：进而可以得到k时刻应用PF算法得到的目标状态变量估计值
Figure FDA00025410928900000118
1.6：进一步应用UKF算法进行混合滤波消除奇异值影响；根据k-1时刻的
Figure FDA00025410928900000119
利用Sigma采样点的选取规则获得2n+1个sigma点/>
Figure FDA00025410928900000120
和对应权值w0,w1,…,w2n，并根据/>
Figure FDA00025410928900000121
接着计算sigma点集的一步预测/>
Figure FDA00025410928900000122
其中，n代表随机变量xk的状态的维数；
1.7：计算目标的状态变量的预测均值
Figure FDA0002541092890000021
及协方差矩阵Pk|k-1；
Figure FDA0002541092890000022
Figure FDA0002541092890000023
/>
式中，
Figure FDA0002541092890000024
和/>
Figure FDA0002541092890000025
为Sigma采样点的权重系数，下标m表示均值，下标c表示协方差；Qk-1为过程噪声方差；
1.8：根据计算所得目标的状态变量的预测均值
Figure FDA0002541092890000026
和协方差矩阵Pk|k-1，再次使用UT变换，产生新的sigma点/>
Figure FDA0002541092890000027
和对应权值w0,w1,…,w2n；
1.9：根据目标的观测方程进行sigma点非线性变换，根据式
Figure FDA0002541092890000028
计算观测量预测采样点/>
Figure FDA0002541092890000029
并根据式/>
Figure FDA00025410928900000210
计算系统观测量的预测均值
Figure FDA00025410928900000211
1.10：根据式
Figure FDA00025410928900000212
计算新息协方差矩阵Pzz，k|k-1、根据式/>
Figure FDA00025410928900000213
计算状态与观测间的互协方差矩阵Pxz，k|k-1，并根据式/>
Figure FDA00025410928900000214
计算滤波增益矩阵Kk；
1.11：最后进行状态更新；根据式
Figure FDA00025410928900000215
计算k时刻后验状态估计均值/>
Figure FDA00025410928900000216
根据式/>
Figure FDA00025410928900000217
计算协方差矩阵Pk；
1.12：循环步骤1.2至步骤1.11，继续进行k+1时刻的混合滤波计算。
2.根据权利要求1所述的目标跟踪算法，其特征在于：所述F的取值如下：
Figure FDA00025410928900000218
3.根据权利要求1所述的目标跟踪算法，其特征在于：设观测站的位置为[x0,y0]，采用的传感器为距离传感器，则
Figure FDA0002541092890000031
4.根据权利要求1-3任一所述的目标跟踪算法的多观测站适应度加权目标跟踪算法，其特征在于：采用a个观测站进行目标跟踪，利用不同观测站在时间、空间及功能上的互补，通过融合这些观测站的观测信息，提高目标跟踪的精度，a为正整数，所述算法包括以下步骤：
步骤1：各观测站采集观测信息，利用UKF和PF混合滤波算法进行各观测站的目标状态估计，并分别记为
Figure FDA0002541092890000032
步骤2：根据上述求出的混合滤波算法的目标观测量的预测均值
Figure FDA0002541092890000033
计算当前k时刻第i个观测站的适应度函数/>
Figure FDA0002541092890000034
上式中：exp表示以e为底的指数函数，/>
Figure FDA0002541092890000035
表示k时刻各观测站的实际目标观测量，/>
Figure FDA0002541092890000036
表示各观测站的观测预测均值，R(i)表示各观测站的观测噪声方差；
步骤3：根据
Figure FDA0002541092890000037
计算k时刻各观测站进行融合时所占的权重/>
Figure FDA0002541092890000038
步骤4：根据
Figure FDA0002541092890000039
对k时刻各观测站生成的目标状态变量估计值进行加权融合，获取相应的目标状态变量融合估计值/>
Figure FDA00025410928900000310
的值；
步骤5：置k＝k+1，循环进行步骤1到步骤4，进行下一时刻跟踪计算。
CN202010547096.9A 2020-06-16 2020-06-16 一种基于ukf和pf混合滤波的目标跟踪算法 Active CN111711432B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010547096.9A CN111711432B (zh) 2020-06-16 2020-06-16 一种基于ukf和pf混合滤波的目标跟踪算法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010547096.9A CN111711432B (zh) 2020-06-16 2020-06-16 一种基于ukf和pf混合滤波的目标跟踪算法

Publications (2)

Publication Number Publication Date
CN111711432A CN111711432A (zh) 2020-09-25
CN111711432B true CN111711432B (zh) 2023-03-28

Family

ID=72540626

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010547096.9A Active CN111711432B (zh) 2020-06-16 2020-06-16 一种基于ukf和pf混合滤波的目标跟踪算法

Country Status (1)

Country Link
CN (1) CN111711432B (zh)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113504012B (zh) * 2021-07-14 2023-12-22 国网福建省电力有限公司检修分公司 基于卡尔曼滤波算法的sf6压力状态判别方法及装置
CN113630106A (zh) * 2021-08-02 2021-11-09 杭州电子科技大学 基于强跟踪滤波的高阶扩展卡尔曼滤波器设计方法
CN115248415B (zh) * 2022-06-27 2024-05-14 中国矿业大学 一种用于矿井非视距场景的目标定位跟踪方法

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102830402A (zh) * 2012-09-10 2012-12-19 江苏科技大学 水下传感器网络目标跟踪系统及方法
WO2015026902A1 (en) * 2013-08-22 2015-02-26 Amazon Technologies, Inc. Multi-tracker object tracking
CN108802692A (zh) * 2018-05-25 2018-11-13 哈尔滨工程大学 一种基于最大互相关熵容积粒子滤波的目标跟踪方法
CN108983215A (zh) * 2018-05-25 2018-12-11 哈尔滨工程大学 一种基于最大互相关熵自适应无迹粒子滤波器的目标跟踪方法
CN111047627A (zh) * 2019-11-14 2020-04-21 中山大学 一种平滑约束无迹卡尔曼滤波方法及目标跟踪方法

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102830402A (zh) * 2012-09-10 2012-12-19 江苏科技大学 水下传感器网络目标跟踪系统及方法
WO2015026902A1 (en) * 2013-08-22 2015-02-26 Amazon Technologies, Inc. Multi-tracker object tracking
CN108802692A (zh) * 2018-05-25 2018-11-13 哈尔滨工程大学 一种基于最大互相关熵容积粒子滤波的目标跟踪方法
CN108983215A (zh) * 2018-05-25 2018-12-11 哈尔滨工程大学 一种基于最大互相关熵自适应无迹粒子滤波器的目标跟踪方法
CN111047627A (zh) * 2019-11-14 2020-04-21 中山大学 一种平滑约束无迹卡尔曼滤波方法及目标跟踪方法

Also Published As

Publication number Publication date
CN111711432A (zh) 2020-09-25

Similar Documents

Publication Publication Date Title
CN111985093B (zh) 一种带噪声估计器的自适应无迹卡尔曼滤波状态估计方法
CN111711432B (zh) 一种基于ukf和pf混合滤波的目标跟踪算法
CN110503071B (zh) 基于变分贝叶斯标签多伯努利叠加模型的多目标跟踪方法
CN109002835B (zh) 一种基于最大熵模糊聚类的粒子滤波数据关联方法
CN111722214B (zh) 雷达多目标跟踪phd实现方法
CN101980044B (zh) 未知测量噪声分布下的多目标跟踪方法
CN109990786A (zh) 机动目标跟踪方法及装置
CN108871365B (zh) 一种航向约束下的状态估计方法及系统
CN106021697A (zh) 一种快速相控阵雷达时间-能量资源联合管理方法
CN107797106A (zh) 一种加速em未知杂波估计的phd多目标跟踪平滑滤波方法
CN111259332B (zh) 一种杂波环境下的模糊数据关联方法及多目标跟踪方法
CN109509207A (zh) 一种对点目标和扩展目标进行无缝跟踪的方法
CN109115228A (zh) 一种基于加权最小二乘容积卡尔曼滤波的目标定位方法
CN116047495B (zh) 一种用于三坐标雷达的状态变换融合滤波跟踪方法
CN111262556B (zh) 一种同时估计未知高斯测量噪声统计量的多目标跟踪方法
Wang et al. Design of mean shift tracking algorithm based on target position prediction
CN104750998B (zh) 基于强度滤波器的被动多传感器目标跟踪方法
CN115905986A (zh) 一种基于联合策略的稳健卡尔曼滤波方法
CN115859626A (zh) 针对周期运动目标的自适应无迹卡尔曼滤波器设计方法
CN115114985A (zh) 一种基于集合理论的传感器系统分布式融合方法
Luo et al. An improved two-phase robust distributed Kalman filter
CN111722213B (zh) 一种机动目标运动参数的纯距离提取方法
CN115169136A (zh) 三维空间下快速uk-gmphd多目标跟踪方法
Wei-Feng et al. Multi-measurement target tracking by using random sampling approach
US20140095565A1 (en) Constrained Grid-Based Filter

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant