CN108896986B - 一种基于预测值的量测转换序贯滤波机动目标跟踪方法 - Google Patents

一种基于预测值的量测转换序贯滤波机动目标跟踪方法 Download PDF

Info

Publication number
CN108896986B
CN108896986B CN201810364436.7A CN201810364436A CN108896986B CN 108896986 B CN108896986 B CN 108896986B CN 201810364436 A CN201810364436 A CN 201810364436A CN 108896986 B CN108896986 B CN 108896986B
Authority
CN
China
Prior art keywords
measurement
model
calculating
pseudo
representing
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
Application number
CN201810364436.7A
Other languages
English (en)
Other versions
CN108896986A (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.)
University of Electronic Science and Technology of China
Original Assignee
University of Electronic Science and Technology of China
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 University of Electronic Science and Technology of China filed Critical University of Electronic Science and Technology of China
Priority to CN201810364436.7A priority Critical patent/CN108896986B/zh
Publication of CN108896986A publication Critical patent/CN108896986A/zh
Application granted granted Critical
Publication of CN108896986B publication Critical patent/CN108896986B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

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
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/66—Radar-tracking systems; Analogous systems
    • G01S13/72—Radar-tracking systems; Analogous systems for two-dimensional tracking, e.g. combination of angle and range tracking, track-while-scan radar
    • G01S13/723—Radar-tracking systems; Analogous systems for two-dimensional tracking, e.g. combination of angle and range tracking, track-while-scan radar by using numerical data

Landscapes

  • Engineering & Computer Science (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明公开了一种基于预测值的量测转换序贯滤波机动目标跟踪方法，涉及雷达目标跟踪领域，主要针对机动目标跟踪情况下的多普勒雷达跟踪领域。本发明在交互多模型算法的框架下结合基于预测值的序贯滤波方法实现目标跟踪。其中，模型概率由位置模型概率和多普勒模型概率共同决定，提高了模型概率估计的准确性；同时，基于预测值的量测转换，消除了量测转换误差统计特性与量测值之间的相关性。该算法是一种有效的可同时处理多普勒量测的机动目标跟踪算法。

Description

一种基于预测值的量测转换序贯滤波机动目标跟踪方法
技术领域
本发明涉及雷达目标跟踪领域，主要针对机动目标跟踪情况下的多普勒雷达跟踪领域。
背景技术
将卡尔曼滤波(Kalman Filter,KF)算法运用到目标跟踪时，其状态方程通常用来描述目标的运动特性，往往是在直角坐标系中建立的，其状态向量主要由目标位置、速度及加速度等参量构成；在观测坐标系下，目标的量测信息通常表现为斜距、方位角及俯仰角。在绝大多数跟踪场景中，目标与观测者处于不同的坐标系时，导致了两者的非线性关系。为解决这类非线性滤波估计问题，最简单的方法是使用传统转换量测方式(ConventionalConversion)，即将极(球)坐标系下的量测值变换到直角坐标系下，使得目标状态与量测呈线性关系。然而在互距离(Cross-range)较大时，即目标位于远距离且量测角度误差较大的情况下。传统转换量测的方法会产生较大的偏差。因此，解决转换偏差所带来的跟踪性能下降的问题是非线性滤波中的首要工作。
为了消除非线性变换带来的偏差，许多学者在此提出了一些改进方法。Lerro D等(见文献Lerro D,Bar-Shalom Y.Tracking with debiased consistent convertedmeasurements versus EKF[J].Aerospace&Electronic Systems IEEE Transactions on,1993,29(3):1015-1022.)提出了一种加性去偏的量测转换方式(DCMKF)，即利用传统转换后的量测值去除转换误差，该方法在方位角度量测误差较小时有较好的去偏效果和量测转换一致性，结果表明其跟踪性能优于EKF算法。但由于量测转换误差的均值和方差是采用嵌套方式得到，Mo L.等提出一种无偏转换量测卡尔曼滤波(UCMKF)算法(见文献：Song X,Zhou Y,Bar-Shalom Y.Unbiased converted measurements for tracking[J].Aerospace&Electronic Systems IEEE Transactions on,1998,34(3):1023-1027.)，该方法采用乘性补偿方式去除偏差，弥补了传统转换方式的缺陷，然而进一步研究发现，该方法在量测转换误差均值和方差的推导上存在兼容性问题。修正的无偏转换量测卡尔曼滤波(MUCMKF)算法直接基于量测值对量测转换方差进行估计，解决了UCMKF的兼容性问题。但在推导量测误差协方差矩阵时，修正无偏量测转换卡尔曼滤波算法(MUCMKF)是在量测值已知的条件下进行的，这会导致量测误差协方差和量测值具有一定的相关性，因此其状态估计也是有偏的。对此进一步研究的相关学者提出了基于预测位置的去相关无偏量测转换卡尔曼滤波算法(DUCMKF)解决了以上算法存在的问题。(见文献：Bordonaro S,Willett P,Barshalom Y.Decorrelated unbiased converted measurement Kalman filter[J].IEEETransactions on Aerospace&Electronic Systems,2014,50(2):1431-1444.)
在目标跟踪领域内，多普勒雷达除了可以获得位置量测信息，还能够提供多普勒量测信息。理论计算与实践已经证明，处理位置量测信息的同时加入多普勒量测，能够更加有效地提升目标的跟踪精度(见文献：Farina A,Studer F A.Radar data processing I-Introduction and Tracking[J].Memorie Della Societa Astronomica Italiana,1985,56.)。以上量测转换算法都仅仅考虑了雷达的位置量测。先前有关多普勒量测的算法中，通常假设多普勒与距离之间的量测误差是统计独立的。然而Bar-Shalom在研究中发现，对于线性调频情况下的某些波形而言，距离和多普勒的量测误差之间存在一定的负相关性(见文献：Bar-Shalom Y.Negative correlation and optimal tracking with Dopplermeasurements[J].Aerospace&Electronic Systems IEEE Transactions on,2001,37(3):1117-1120.)。为了提升雷达的跟踪精度，又需要减小引入的多普勒量测时导致与位置之间的较强非线性，可将斜距和多普勒量测的乘积构造伪量测，来减小这种相关性。(见文献：ZhanshengDuan,Chongzhao Han and X.Rong Li,Sequential Nonlinear TrackingFilter with Range-rateMeasurements in Spherical Coordinates,7th InternationalConference on Information Fusion,Stockholm,2004,131-138.)段战胜等将多普勒量测信息引入到量测转换卡尔曼滤波算法中，将其推广为序贯滤波跟踪算法，但在求解转换量测误差统计特性时采用了嵌套的条件期望，且误差与量测具有相关性也会引入偏差。
当目标具有机动特性时，以上算法采用固定的系统状态方程无法描述目标的运动状态。交互多模型(IMM，interacting multiple model)算法能够有效的解决机动目标跟踪时存在的问题(见文献：Blom H AP,Barshalom Y.The interacting multiple modelalgorithm for systems with Markovian switching coefficients[J].IEEETransactions on Automatic Control,1988,33(8):780-783.)。该算法假设有限个目标模型存在，每个模型对应于不同的机动输入水平。在计算出各模型为正确的后验概率之后，就可以通过对各模型正确时的状态估计加权求和来给出最终的目标状态估计，其中的加权因子为模型正确的后验概率。针对机动目标情况，将基于预测值的序贯滤波算法与交互多模型算法相结合。本发明首先对量测值进行无偏转换处理，基于预测值推导转换误差统计特性，将序贯滤波算法与交互多模型算法框架相结合，利用基于线性KF算法的位置量测滤波器和基于非线性二阶EKF的伪量测滤波器所获得的模型概率平均值进行状态加权，从而获得一种基于预测值量测转换的序贯滤波机动目标跟踪算法。
发明内容
本发明针对背景技术的不足之处，考虑在机动目标跟踪情况下，利用多普勒雷达所提取的目标位置信息和多普勒信息进行序贯滤波，与交互多模型算法相结合，提出了一种针对机动目标的非线性量测转换序贯滤波跟踪算法。
本发明由混合去偏(位置量测采用相乘去偏方式；伪量测采用相减去偏方式)量测转换模块获得各模型的位置状态信息和多普勒信息的转换量测值，输入到滤波器输入端的交互模块，计算出交互后的各模型输入估计值，输出给各目标运动模型的滤波模块进行滤波：由位置量测滤波器产生的估计结果作为下一时刻的输入交互计算进行迭代，同时输出到伪量测滤波器中进行序贯处理，由伪量测滤波器产生的估计结果作为该模型的滤波估计。所述各目标运动模型滤波模块内的位置状态滤波器和伪量测滤波器，两者的新息计算结果共同输出给各模型概率计算模块，分别得到位置量测概率计算结果和伪量测概率计算结果后，再计算平均概率，作为当前时刻更新的模型概率。输出端的估计组合模块将各模型的滤波估计与平均概率进行加权求和计算得到目标位置状态，实现目标跟踪。
现对下文将用到的部分符号进行说明。目标在直角坐标系下的状态可描述为
Figure BDA0001636797590000031
系统在x,y方向上的过程噪声假设为零均值高斯噪声，其协方差Qk＝diag(qx,qy)，q表示各方向上的噪声方差；设多普勒雷达采样周期为T，获取的量测数据：
Figure BDA0001636797590000032
和
Figure BDA0001636797590000033
分别为距离量测、方位角量测和多普勒量测(下标k表示时刻点，上标m量测数据)；其对应的量测噪声
Figure BDA0001636797590000034
和
Figure BDA0001636797590000035
是零均值高斯白噪声序列，方差分别为
Figure BDA0001636797590000036
和
Figure BDA0001636797590000037
且
Figure BDA0001636797590000038
和
Figure BDA0001636797590000039
不相关，
Figure BDA00016367975900000310
和
Figure BDA00016367975900000311
的相关系数为ρ。
本发明的书方案为一种基于预测值的量测转换序贯滤波机动目标跟踪方法，该方法包括：
步骤1：计算量测转换值；
采用无偏量测转换方式计算
Figure BDA00016367975900000312
计算方法如下：
Figure BDA00016367975900000313
其中，
Figure BDA00016367975900000314
表示k时刻转换后的量测值，
Figure BDA00016367975900000315
表示k时刻转换后的量测值在x方向上的分量，
Figure BDA00016367975900000316
表示k时刻转换后的量测值在y方向上的分量，
Figure BDA00016367975900000317
表示k时刻转换后伪量测值，
Figure BDA00016367975900000318
为乘性偏差补偿因子；
Figure BDA00016367975900000319
和
Figure BDA00016367975900000320
分别为距离量测、方位角量测和多普勒量测，下标k表示时刻点，上标m量测数据；ρ相关系数，σr、
Figure BDA0001636797590000041
距离和方位角的量测误差标准差；
步骤2：计算各模型的输入交互；
首先假设目标在k时刻运动模式表示为Mk，若多模型系统中存在N个模型，这N个模型能包含任何可能的目标运动模式，那么：
Figure BDA0001636797590000042
模型之间由模型i转移到模型j的转移概率为πij；
步骤2.1：计算混合概率
Figure BDA0001636797590000043
k时刻运动模型为j条件下，k-1时刻运动模型为i的混合概率
Figure BDA0001636797590000044
为：
Figure BDA0001636797590000045
其中，
Figure BDA0001636797590000046
为k-1时刻运动模型i的更新概率，归一化常数：
Figure BDA0001636797590000047
步骤2.2：计算各模型目标状态估计的交互输入；
通过对上一时刻估计值
Figure BDA0001636797590000048
和协方差
Figure BDA0001636797590000049
的交互计算得到每个模型滤波器输入的状态估计和协方差：
Figure BDA00016367975900000410
Figure BDA00016367975900000411
其中：
Figure BDA00016367975900000412
和
Figure BDA00016367975900000413
表示每个模型滤波器输入的状态估计和协方差；
步骤3：模型j的位置量测滤波估计；
基于预测值计算量测转换误差统计特性，并进行序贯滤波，包括以下几个子步骤：现为了公式表达简洁，将原处于每个符号表示当前目标状态为模型j的上标(j)省略，后续中每个符号均隐含有状态模型为j的含义；
步骤3.1:由卡尔曼一步预测方程计算预测状态
Figure BDA00016367975900000414
Figure BDA0001636797590000051
其中，
Figure BDA0001636797590000052
分别表示目标在x和y方向上的预测位置，预测速度和预测加速度值；
步骤3.2：计算预测误差协方差Pk,k-1；
步骤3.3：计算极坐标系下的距离rt、方位角θt和多普勒
Figure BDA0001636797590000053
的预测值：
Figure BDA0001636797590000054
Figure BDA0001636797590000055
Figure BDA0001636797590000056
步骤3.4：计算预测值对应的预测误差方差
Figure BDA0001636797590000057
Figure BDA0001636797590000058
Figure BDA0001636797590000059
Figure BDA00016367975900000510
其中：
Figure BDA00016367975900000511
表示距离的预测方差，02×1表示2行1列的零矩阵，
Figure BDA00016367975900000512
表示方位角的预测方差，
Figure BDA0001636797590000061
表示径向速度的预测方差；
步骤3.5：计算量测转换误差协方差Rk；
在上一步计算出预测距离、方位角、多普勒信息rt、θt、
Figure BDA0001636797590000062
及其预测误差方差条件下,计算去偏量测转换误差协方差Rk，下标k表示第k时刻；
Figure BDA0001636797590000063
其中，协方差矩阵各元素计算方式如下：
Figure BDA0001636797590000064
Figure BDA0001636797590000065
Figure BDA0001636797590000066
Figure BDA0001636797590000067
Figure BDA0001636797590000068
Figure BDA0001636797590000069
步骤3.6：分块量测协方差矩阵，构造位置量测协方差矩阵和伪量测互协方差矩阵；
Figure BDA00016367975900000610
步骤3.7：计算位置量测新息及新息误差协方差矩阵；
Figure BDA0001636797590000071
Figure BDA0001636797590000072
其中，
Figure BDA0001636797590000073
为位置量测观测矩阵，
Figure BDA0001636797590000074
表示k时刻转换后的位置量测值，
Figure BDA0001636797590000075
表示位置的一步预测协方差矩阵；
步骤3.8：计算位置量测滤波器的增益；
Figure BDA0001636797590000076
步骤3.9：计算位置状态估计
Figure BDA0001636797590000077
和状态估计误差协方差矩阵
Figure BDA0001636797590000078
步骤4：模型j的伪量测滤波估计；
为了公式表达的简洁性，将原处于每个符号表示当前目标状态为模型j的上标(j)省略，后续中每个符号均隐含有状态模型为j的含义；
步骤4.1：消除相关性；
去除位置量测与伪量测的相关性；将基于预测值求取的量测误差协方差矩阵
Figure BDA0001636797590000079
进行分块：
Figure BDA00016367975900000710
其中，令
Figure BDA00016367975900000711
则有：
Figure BDA00016367975900000712
其中：Lk表示由量测误差协方差
Figure BDA00016367975900000713
构造的分解向量，
Figure BDA00016367975900000714
表示Lk中的第一列元素，
Figure BDA00016367975900000715
表示Lk中的第二列元素，Bk表示由单位阵和Lk构成的分解矩阵，I2表示2阶单位向量；
由矩阵的Cholesky分解得到：
Figure BDA00016367975900000716
Figure BDA0001636797590000081
其中，
Figure BDA0001636797590000082
表示位置量测噪声，
Figure BDA0001636797590000083
表示位置量测观测矩阵；
Figure BDA0001636797590000084
Figure BDA0001636797590000085
Figure BDA0001636797590000086
其中：
Figure BDA0001636797590000087
分别表示k时刻的x,y方向的转换量测值和转换后的伪量测值，
Figure BDA0001636797590000088
Figure BDA0001636797590000089
分别表示k时刻的x,y方向的转换量测误差均值和伪量测转换误差均值；
步骤4.2：计算与伪量测信息相关的量测转换误差协方差
Figure BDA00016367975900000810
Figure BDA00016367975900000811
步骤4.3：这里使用二阶EKF对伪量测进行序贯滤波估计，计算伪量测新息及其协方差方法如下；
Figure BDA00016367975900000812
Figure BDA00016367975900000813
其中，
Figure BDA00016367975900000814
表示与伪量测转换误差的均值，
Figure BDA00016367975900000815
表示由位置估计结果
Figure BDA00016367975900000816
变换得到的伪量测真实值，
Figure BDA00016367975900000817
表示与伪量测相关的量测转换误差协方差矩阵，
Figure BDA00016367975900000818
表示位置滤波的估计协方差，
Figure BDA00016367975900000819
为
Figure BDA00016367975900000820
在
Figure BDA00016367975900000821
处的Jacobian矩阵；
Figure BDA00016367975900000822
用
Figure BDA00016367975900000823
来表示位置估计误差协方差矩阵
Figure BDA00016367975900000824
的第i行第j列的元素，则
Figure BDA00016367975900000825
步骤4.4：计算伪量测滤波增益
Figure BDA00016367975900000826
Figure BDA00016367975900000827
步骤4.5：计算伪量测估计结果及估计误差协方差:
Figure BDA0001636797590000091
Figure BDA0001636797590000092
步骤4.6：最终模型j的滤波估计及误差自相关矩阵为：
Figure BDA0001636797590000093
步骤5：模型概率计算模块；
步骤5.1:计算位置量测滤波器似然函数；
Figure BDA0001636797590000094
其中：
Figure BDA0001636797590000095
表示位置新息自相关矩阵，
Figure BDA0001636797590000096
表示位置新息，上标j表示运动模型，p表示与位置量测相关的信息，下标k表示当前时刻；
步骤5.2:计算位置模型概率；
Figure BDA0001636797590000097
其中，
Figure BDA0001636797590000098
C(j),p表示模型j的位置归一化常数；
步骤5.3:计算伪量测滤波器似然函数；
Figure BDA0001636797590000099
其中，
Figure BDA00016367975900000910
表示伪量测新息自相关矩阵，
Figure BDA00016367975900000911
表示伪量测新息，上标j表示运动模型，ε表示与伪量测相关的信息，下标k表示当前时刻；
步骤5.4:计算伪量测模型概率；
Figure BDA00016367975900000912
其中，
Figure BDA00016367975900000913
C(j),ε表示模型j的伪量测归一化常数；
步骤5.5:计算运动模型j平均概率
Figure BDA00016367975900000914
Figure BDA00016367975900000915
步骤6：计算加权的状态估计
Figure BDA0001636797590000101
和估计误差自相关矩阵Pk,k；
Figure BDA0001636797590000102
Figure BDA0001636797590000103
步骤7：重复步骤1-6进行递推滤波算法，实现目标跟踪。
从上述步骤可见，本发明在交互多模型算法的框架下结合基于预测值的序贯滤波方法实现目标跟踪。其中，模型概率由位置模型概率和多普勒模型概率共同决定，提高了模型概率估计的准确性；同时，基于预测值的量测转换，消除了量测转换误差统计特性与量测值之间的相关性。该算法是一种有效的可同时处理多普勒量测的机动目标跟踪算法。
在本发明中每个模型可看作是一个序贯滤波器，都包含有一个位置量测滤波器和伪量测滤波器。每个序贯滤波都有两个计算模型概率的部分，分别是位置滤波估计新息及其协方差矩阵计算位置信息的模型概率和多普勒量测滤波估计新息及其协方差矩阵计算多普勒信息的模型概率估计。为了提高模型概率的准确性，本发明分别计算位置信息和多普勒信息的模型概率估计后计算它们的均值。
附图说明
图1为本发明一种基于预测值的序贯滤波机动目标跟踪流程图；
图2为本发明交互多模型算法示意图；
图3为本发明具体实施方式中场景1 RMSE曲线；
图4为本发明具体实施方式中场景2 RMSE曲线；
图5为本发明具体实施方式中场景3 RMSE曲线。
其中，图3a，图4a，图5a是SQU-PRE和SQU-PRE-IMM的RMSE对比结果；图3b，图4b，图5b是SQ-COVN-IMM和SQU-PRE-IMM的RMSE对比结果。
具体实施方式
在此次实施方案中，考虑对二维平面内机动场景下的某一运动目标进行跟踪，目标初始值位置为(x＝5000m,y＝1000m)，初始速度为(vx＝25m/s,vy＝2m/s)。0～60s目标做匀速运动，60s～85s做匀加速运动，加速度为(-2m/s2,2m/s2)，85s～135s做匀速运动，135s～150s做匀加速运动，加速度为(3m/s2,-0.5m/s2)，150s～200s做匀速运动。径向速度和距离的相关系数为0.5。针对距离、方位角、多普勒量测噪声标准差取不同的值，做如下场景设置：
场景一：传感器的距离量测误差为
Figure BDA0001636797590000111
方位角的量测误差为
Figure BDA0001636797590000112
度，多普勒径向速度的量测误差为0.1m/s；
场景二：传感器的距离量测误差为
Figure BDA0001636797590000113
方位角的量测误差为
Figure BDA0001636797590000114
度，多普勒径向速度的量测误差为0.3m/s；
场景三：传感器的距离量测误差为
Figure BDA0001636797590000115
方位角的量测误差为
Figure BDA0001636797590000116
度，多普勒径向速度的量测误差为0.5m/s；
在上述3个场景下，采用基于传统量测转换序贯滤波的交互多模型方法(SQ-COVN-IMM)和基于预测值的序贯滤波方法(SQU-PRE)与本发明的基于预测值的量测转换序贯滤波机动目标跟踪方法(SQU-PRE-IMM)进行500次蒙特卡罗仿真。
图3、图4和图5分别采取上述三种机动目标跟踪方法的仿真结果的位置估计均方误差(RMSE)曲线。其中，图(3a)，(4a)，(5a)是SQU-PRE和SQU-PRE-IMM的RMSE对比结果，图(3b)，(4b)，(5b)是SQ-COVN-IMM和SQU-PRE-IMM的RMSE对比结果。
实施结果分析：一种基于预测值的序贯滤波机动目标跟踪方法可以实现带有多普勒量测的机动场景下的目标跟踪，并且跟踪精度比基于传统量测转换序贯滤波的交互多模型方法高；由于考虑了目标的机动特性，其跟踪精度也比没使用交互多模型的基于预测值的序贯滤波算法高。
三个场景下，根据SQU-PRE和SQU-PRE-IMM的对比结果可以看出，由于目标进行机动运动，因此SQU-PRE-IMM在交互多模型框架下进行目标跟踪，能获得较好的跟踪效果，而SQU-PRE在匀速运动阶段能够保持较好的跟踪效果，但从目标开始做匀加速运动时，滤波器中的模型与目标真实运动状态不相匹配，滤波误差大幅度增大，直至85秒时目标恢复匀速运动，滤波误差才开始下降，在135秒进行第二次匀加速运动时，曲线又出现发散情况。另一方面，对比SQ-CONV-IMM可以看出，SQU-PRE-IMM具有较好的跟踪性能和更高的跟踪精度，三个场景下随着量测误差的增加更加明显。
综上所述，基于预测值的序贯滤波机动目标跟踪算法SQU-PRE-IMM可以有效的实现机动目标跟踪，且具有较高的跟踪精度。

Claims (1)

1.一种基于预测值的量测转换序贯滤波机动目标跟踪方法，该方法包括：
步骤1：计算量测转换值；
采用无偏量测转换方式计算
Figure FDA0003558649940000011
计算方法如下：
Figure FDA0003558649940000012
其中，
Figure FDA0003558649940000013
表示k时刻转换后的量测值，
Figure FDA0003558649940000014
表示k时刻转换后的量测值在x方向上的分量，
Figure FDA0003558649940000015
表示k时刻转换后的量测值在y方向上的分量，
Figure FDA0003558649940000016
表示k时刻转换后伪量测值，
Figure FDA0003558649940000017
为乘性偏差补偿因子；
Figure FDA0003558649940000018
和
Figure FDA0003558649940000019
分别为距离量测、方位角量测和多普勒量测，下标k表示时刻点，上标m表示量测数据；ρ表示相关系数，σr和
Figure FDA00035586499400000118
分别表示距离和方位角的量测误差标准差；
步骤2：计算各模型的输入交互；
首先假设目标在k时刻运动模式表示为Mk，若多模型系统中存在N个模型，这N个模型能包含任何可能的目标运动模式，那么：
Figure FDA00035586499400000110
模型之间由模型i转移到模型j的转移概率为πij；
步骤2.1：计算混合概率
Figure FDA00035586499400000111
k时刻运动模型为j条件下，k-1时刻运动模型为i的混合概率
Figure FDA00035586499400000112
为：
Figure FDA00035586499400000113
其中，
Figure FDA00035586499400000114
为k-1时刻运动模型i的更新概率，归一化常数：
Figure FDA00035586499400000115
步骤2.2：计算各模型目标状态估计的交互输入；
通过对上一时刻估计值
Figure FDA00035586499400000116
和协方差
Figure FDA00035586499400000117
的交互计算得到每个模型滤波器输入的状态估计和协方差：
Figure FDA0003558649940000021
Figure FDA0003558649940000022
其中：
Figure FDA0003558649940000023
和
Figure FDA0003558649940000024
表示每个模型滤波器输入的状态估计和协方差；
步骤3：模型j的位置量测滤波估计；
基于预测值计算量测转换误差统计特性，并进行序贯滤波，包括以下几个子步骤：现为了公式表达简洁，将原处于每个符号表示当前目标状态为模型j的上标(j)省略，后续中每个符号均隐含有状态模型为j的含义；
步骤3.1:由卡尔曼一步预测方程计算预测状态
Figure FDA0003558649940000025
Figure FDA0003558649940000026
其中，
Figure FDA0003558649940000027
分别表示目标在x和y方向上的预测位置，预测速度和预测加速度值；
步骤3.2：计算预测误差协方差Pk,k-1；
步骤3.3：计算极坐标系下的距离rt、方位角θt和多普勒信息
Figure FDA0003558649940000028
的预测值：
Figure FDA0003558649940000029
Figure FDA00035586499400000210
Figure FDA00035586499400000211
步骤3.4：计算预测值对应的预测误差方差
Figure FDA00035586499400000212
Figure FDA00035586499400000213
Figure FDA0003558649940000031
Figure FDA0003558649940000032
其中：
Figure FDA0003558649940000033
表示距离的预测方差，02×1表示2行1列的零矩阵，
Figure FDA0003558649940000034
表示方位角的预测方差，
Figure FDA0003558649940000035
表示径向速度的预测方差；
步骤3.5：计算量测转换误差协方差Rk；
在上一步计算出预测距离rt、方位角θt、多普勒信息
Figure FDA0003558649940000036
及其预测误差方差条件下,计算去偏量测转换误差协方差Rk，下标k表示第k时刻；
Figure FDA0003558649940000037
其中，协方差矩阵各元素计算方式如下：
Figure FDA0003558649940000038
Figure FDA0003558649940000039
Figure FDA0003558649940000041
Figure FDA0003558649940000042
Figure FDA0003558649940000043
Figure FDA0003558649940000044
步骤3.6：分块量测协方差矩阵，构造位置量测协方差矩阵和伪量测互协方差矩阵；
Figure FDA0003558649940000045
步骤3.7：计算位置量测新息及新息误差协方差矩阵；
Figure FDA0003558649940000046
Figure FDA0003558649940000047
其中，
Figure FDA0003558649940000048
为位置量测观测矩阵，
Figure FDA0003558649940000049
表示k时刻转换后的位置量测值，
Figure FDA00035586499400000410
表示位置的一步预测协方差矩阵；
步骤3.8：计算位置量测滤波器的增益；
Figure FDA00035586499400000411
步骤3.9：计算位置状态估计
Figure FDA00035586499400000412
和状态估计误差协方差矩阵
Figure FDA00035586499400000413
步骤4：模型j的伪量测滤波估计；
为了公式表达的简洁性，将原处于每个符号表示当前目标状态为模型j的上标(j)省略，后续中每个符号均隐含有状态模型为j的含义；
步骤4.1：消除相关性；
去除位置量测与伪量测的相关性；将基于预测值求取的量测误差协方差矩阵
Figure FDA00035586499400000414
进行分块：
Figure FDA0003558649940000051
其中，令
Figure FDA0003558649940000052
则有：
Figure FDA0003558649940000053
其中：Lk表示由量测误差协方差
Figure FDA0003558649940000054
构造的分解向量，
Figure FDA0003558649940000055
表示Lk中的第一列元素，
Figure FDA0003558649940000056
表示Lk中的第二列元素，Bk表示由单位阵和Lk构成的分解矩阵，I2表示2阶单位向量；
由矩阵的Cholesky分解得到：
Figure FDA0003558649940000057
Figure FDA0003558649940000058
其中，
Figure FDA0003558649940000059
表示位置量测噪声，
Figure FDA00035586499400000510
表示位置量测观测矩阵；
Figure FDA00035586499400000511
Figure FDA00035586499400000512
Figure FDA00035586499400000513
其中：
Figure FDA00035586499400000514
分别表示k时刻的x,y方向的转换量测值和转换后的伪量测值，
Figure FDA00035586499400000515
Figure FDA00035586499400000516
分别表示k时刻的x,y方向的转换量测误差均值和伪量测转换误差均值；
步骤4.2：计算与伪量测信息相关的量测转换误差协方差
Figure FDA00035586499400000517
Figure FDA00035586499400000518
步骤4.3：这里使用二阶EKF对伪量测进行序贯滤波估计，计算伪量测新息及其协方差方法如下；
Figure FDA00035586499400000519
Figure FDA0003558649940000061
其中，
Figure FDA0003558649940000062
表示与伪量测转换误差的均值，
Figure FDA0003558649940000063
表示由位置估计结果
Figure FDA0003558649940000064
变换得到的伪量测真实值，
Figure FDA0003558649940000065
表示与伪量测相关的量测转换误差协方差矩阵，
Figure FDA0003558649940000066
表示位置滤波的估计协方差，
Figure FDA0003558649940000067
为
Figure FDA0003558649940000068
在
Figure FDA0003558649940000069
处的Jacobian矩阵；
Figure FDA00035586499400000610
用
Figure FDA00035586499400000611
来表示位置估计误差协方差矩阵
Figure FDA00035586499400000612
的第i行第j列的元素，则
Figure FDA00035586499400000613
步骤4.4：计算伪量测滤波增益
Figure FDA00035586499400000614
Figure FDA00035586499400000615
步骤4.5：计算伪量测估计结果及估计误差协方差:
Figure FDA00035586499400000616
Figure FDA00035586499400000617
步骤4.6：最终模型j的滤波估计及误差自相关矩阵为：
Figure FDA00035586499400000618
步骤5：模型概率计算模块；
步骤5.1:计算位置量测滤波器似然函数；
Figure FDA00035586499400000619
其中：
Figure FDA00035586499400000620
表示位置新息自相关矩阵，
Figure FDA00035586499400000621
表示位置新息，上标j表示运动模型，p表示与位置量测相关的信息，下标k表示当前时刻；
步骤5.2:计算位置模型概率；
Figure FDA00035586499400000622
其中，
Figure FDA00035586499400000623
C(j),p表示模型j的位置归一化常数；
步骤5.3:计算伪量测滤波器似然函数；
Figure FDA0003558649940000071
其中，
Figure FDA0003558649940000072
表示伪量测新息自相关矩阵，
Figure FDA0003558649940000073
表示伪量测新息，上标j表示运动模型，ε表示与伪量测相关的信息，下标k表示当前时刻；
步骤5.4:计算伪量测模型概率；
Figure FDA0003558649940000074
其中，
Figure FDA0003558649940000075
C(j),ε表示模型j的伪量测归一化常数；
步骤5.5:计算运动模型j平均概率
Figure FDA0003558649940000076
Figure FDA0003558649940000077
步骤6：计算加权的状态估计
Figure FDA0003558649940000078
和估计误差自相关矩阵Pk,k；
Figure FDA0003558649940000079
Figure FDA00035586499400000710
步骤7：重复步骤1-6进行递推滤波算法，实现目标跟踪。
CN201810364436.7A 2018-04-23 2018-04-23 一种基于预测值的量测转换序贯滤波机动目标跟踪方法 Expired - Fee Related CN108896986B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810364436.7A CN108896986B (zh) 2018-04-23 2018-04-23 一种基于预测值的量测转换序贯滤波机动目标跟踪方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810364436.7A CN108896986B (zh) 2018-04-23 2018-04-23 一种基于预测值的量测转换序贯滤波机动目标跟踪方法

Publications (2)

Publication Number Publication Date
CN108896986A CN108896986A (zh) 2018-11-27
CN108896986B true CN108896986B (zh) 2022-06-03

Family

ID=64342364

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810364436.7A Expired - Fee Related CN108896986B (zh) 2018-04-23 2018-04-23 一种基于预测值的量测转换序贯滤波机动目标跟踪方法

Country Status (1)

Country Link
CN (1) CN108896986B (zh)

Families Citing this family (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109655826B (zh) * 2018-12-16 2021-06-08 成都汇蓉国科微系统技术有限公司 一种低慢小目标轨迹滤波方法及装置
CN110501696B (zh) * 2019-06-28 2022-05-31 电子科技大学 一种基于多普勒量测自适应处理的雷达目标跟踪方法
CN110516193B (zh) * 2019-07-08 2023-03-21 杭州电子科技大学 笛卡尔坐标系下基于转换瑞利滤波器的机动目标跟踪方法
CN110426702B (zh) * 2019-07-24 2021-06-08 中国人民解放军海军航空大学 基于多雷达的系统偏差与目标状态的序贯滤波方法及系统
CN110716196B (zh) * 2019-11-04 2023-04-25 广东中科四创科技有限公司 一种多点低慢小空中目标跟踪临视方法
CN111008364B (zh) * 2019-12-09 2023-05-09 北京壹氢科技有限公司 一种双观测者协同被动定位方法及系统
CN111077518B (zh) * 2019-12-20 2020-11-06 哈尔滨工业大学 一种基于距离-多普勒量测的跟踪滤波方法及装置
CN111190173B (zh) * 2020-01-14 2023-02-14 电子科技大学 一种基于预测值量测转换的相控阵雷达目标跟踪方法
CN111650577B (zh) * 2020-06-12 2022-05-24 电子科技大学 极坐标系下含多普勒量测的机动目标跟踪方法
CN111965618B (zh) * 2020-08-18 2022-09-23 中国电子科技集团公司第三十八研究所 一种融合多普勒量测的转换量测跟踪方法及系统

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105785359A (zh) * 2016-05-19 2016-07-20 哈尔滨工业大学 一种多约束机动目标跟踪方法
CN106646453A (zh) * 2016-11-17 2017-05-10 电子科技大学 一种基于预测值量测转换的多普勒雷达目标跟踪方法
CN106950562A (zh) * 2017-03-30 2017-07-14 电子科技大学 一种基于预测值量测转换的状态融合目标跟踪方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
SG101444A1 (en) * 2001-04-20 2004-01-30 Singapore Tech Aerospace Ltd A maneuvering target tracking method via modifying the interacting multiple model (imm) and the interacting acceleration compensation (iac) algorithms

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105785359A (zh) * 2016-05-19 2016-07-20 哈尔滨工业大学 一种多约束机动目标跟踪方法
CN106646453A (zh) * 2016-11-17 2017-05-10 电子科技大学 一种基于预测值量测转换的多普勒雷达目标跟踪方法
CN106950562A (zh) * 2017-03-30 2017-07-14 电子科技大学 一种基于预测值量测转换的状态融合目标跟踪方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
A nonlinear tracking algorithm with range-rate measurements based on unbiased measurement conversion;Lianmeng Jiao et al.;《2012 15th International Conference on Information Fusion》;20120831;第1400-1405页 *
一种改进的自适应网格交互多模型跟踪算法;唐婷 等;《信号处理》;20090531;第25卷(第5期);第816-819页 *
基于量测转换与输入估计的机动目标跟踪算法;盛琥 等;《系统工程与电子技术》;20150131;第37卷(第1期);第31-36页 *

Also Published As

Publication number Publication date
CN108896986A (zh) 2018-11-27

Similar Documents

Publication Publication Date Title
CN108896986B (zh) 一种基于预测值的量测转换序贯滤波机动目标跟踪方法
CN108226920B (zh) 一种基于预测值处理多普勒量测的机动目标跟踪系统及方法
CN111985093B (zh) 一种带噪声估计器的自适应无迹卡尔曼滤波状态估计方法
CN106950562B (zh) 一种基于预测值量测转换的状态融合目标跟踪方法
CN107045125B (zh) 一种基于预测值量测转换的交互多模型雷达目标跟踪方法
CN110208792B (zh) 同时估计目标状态和轨迹参数的任意直线约束跟踪方法
CN110501696B (zh) 一种基于多普勒量测自适应处理的雷达目标跟踪方法
CN106372646B (zh) 基于srck-gmcphd滤波的多目标跟踪方法
Gao et al. Adaptive unscented Kalman filter based on maximum posterior and random weighting
KR100816269B1 (ko) 언센티드 필터를 적용한 강인한 동시 위치 추정 및 지도작성 방법
Aftab et al. A learning Gaussian process approach for maneuvering target tracking and smoothing
CN104035083B (zh) 一种基于量测转换的雷达目标跟踪方法
CN104199022B (zh) 一种基于目标模态估计的临近空间高超声速目标跟踪方法
CN111965618B (zh) 一种融合多普勒量测的转换量测跟踪方法及系统
CN110231620B (zh) 一种噪声相关系统跟踪滤波方法
CN107688179A (zh) 基于多普勒信息辅助的综合概率数据互联方法
Kim et al. Exactly Rao-Blackwellized unscented particle filters for SLAM
CN116047498A (zh) 基于最大相关熵扩展卡尔曼滤波的机动目标跟踪方法
CN111693984A (zh) 一种改进的ekf-ukf动目标跟踪方法
CN113030945B (zh) 一种基于线性序贯滤波的相控阵雷达目标跟踪方法
CN107886058B (zh) 噪声相关的两阶段容积Kalman滤波估计方法及系统
CN111190173B (zh) 一种基于预测值量测转换的相控阵雷达目标跟踪方法
Ma et al. Tracking of maneuvering star-convex extended target using modified adaptive extended Kalman Filter
CN114898202A (zh) 基于多模型融合的水下视频目标尺度空间判别式跟踪系统及方法
CN114445459A (zh) 基于变分贝叶斯理论的连续-离散最大相关熵目标跟踪方法

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
CF01 Termination of patent right due to non-payment of annual fee
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20220603