CN109188422B - 一种基于lu分解的卡尔曼滤波目标跟踪方法 - Google Patents

一种基于lu分解的卡尔曼滤波目标跟踪方法 Download PDF

Info

Publication number
CN109188422B
CN109188422B CN201810898258.6A CN201810898258A CN109188422B CN 109188422 B CN109188422 B CN 109188422B CN 201810898258 A CN201810898258 A CN 201810898258A CN 109188422 B CN109188422 B CN 109188422B
Authority
CN
China
Prior art keywords
target
matrix
state
moment
decomposition
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
CN201810898258.6A
Other languages
English (en)
Other versions
CN109188422A (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.)
Leihua Electronic Technology Research Institute Aviation Industry Corp of China
Original Assignee
Leihua Electronic Technology Research Institute Aviation Industry Corp 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 Leihua Electronic Technology Research Institute Aviation Industry Corp of China filed Critical Leihua Electronic Technology Research Institute Aviation Industry Corp of China
Priority to CN201810898258.6A priority Critical patent/CN109188422B/zh
Publication of CN109188422A publication Critical patent/CN109188422A/zh
Application granted granted Critical
Publication of CN109188422B publication Critical patent/CN109188422B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO 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/00Systems 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/66Radar-tracking systems; Analogous systems
    • G01S13/72Radar-tracking systems; Analogous systems for two-dimensional tracking, e.g. combination of angle and range tracking, track-while-scan radar
    • G01S13/723Radar-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

本发明涉及一种基于LU分解的卡尔曼滤波目标跟踪方法,包括:获得当前时刻目标的状态测量信息;设置目标的初始参数,所述初始化参数信息包括目标运动模型和测量噪声;根据上一时刻的目标状态及初始参数预测目标当前时刻的状态信息;对上一时刻的协方差矩阵进行LU分解,并根据分解后的协方差矩阵及初始参数预测当前时刻的协方差矩阵;根据当前时刻目标的状态测量信息与预测的目标状态信息的预测差值计算最优增益矩阵;根据最优增益矩阵更新目标最优状态;对预测的协方差矩阵进行LU分解,并根据最优增益矩阵更新协方差矩阵。本发明解决了目标跟踪过程中卡尔曼滤波的数值稳定性问题,具有速度快,精度高。

Description

一种基于LU分解的卡尔曼滤波目标跟踪方法
技术领域
本发明属于雷达信号处理技术领域,尤其涉及一种基于LU分解的卡尔曼滤波目标跟踪方法。
背景技术
卡尔曼滤波(Kalman filtering)一种利用线性系统状态方程,通过系统输入/输出观测数据,对系统状态进行最优估计的算法。经过多年的研究实践,卡尔曼滤波已经广泛地应用到了工程中,但是卡尔曼滤波在线性的前提假设下是一个线性无偏、最小方差估计器,从而可以为线性滤波问题提供精确的解析解。自卡尔曼滤波提出以来,它已成为控制、信号处理与通信工程等领域最基本、最重要的计算方法和工具之一,并成功地应用到航空、航天、电力系统及社会经济等不同领域。
然而随着微型计算机的普及应用,对卡尔曼滤波的数值稳定性、实用性和有效性的要求越来越高,为此人们在如何改善卡尔曼滤波的计算复杂性和数值稳定性方面做了大量的探索工作。卡尔曼滤波算法是目标跟踪领域和数据处理领域研究比较多的问题,也相应形成了不少方法。虽然卡尔曼滤波算法能够取得比较准确的结果,但仍然可能出现滤波发散的问题。这是由于计算机的字长有限而引起的舍入误差使得矩阵P(k|k)失去正定性或对称性,从而导致矩阵K(k)计算失真,造成发散。因此,需要不同的平方根滤波算法来解决上述问题。
发明内容
本发明的目的是提供一种基于LU分解的卡尔曼滤波目标跟踪方法,提高多功能雷达目标跟踪的稳健性,降低目标跟踪的运算量。
为达到上述目的,本发明采用的技术方案是:一种基于LU分解的卡尔曼滤波目标跟踪方法,所述目标跟踪方法包括:
S1:获得当前时刻目标的状态测量信息,所述状态测量信息包括目标距离、目标方位角、目标俯仰角和目标速度;
S2:设置目标的初始参数,所述初始参数信息包括目标运动模型和测量噪声;
S3:根据上一时刻的目标状态及初始参数预测目标当前时刻的状态信息;
S4:对上一时刻的协方差矩阵进行LU分解,并根据分解后的协方差矩阵及初始参数预测当前时刻的协方差矩阵;
S5:根据当前时刻目标的状态测量信息与预测的目标状态信息的预测差值计算最优增益矩阵;
S6:根据最优增益矩阵更新目标最优状态;
S7:对预测的协方差矩阵进行LU分解,并根据最优增益矩阵更新协方差矩阵。
其中,步骤S3中,目标当前时刻的状态信息预测方法为:
X(k|k-1)=F(k|k-1)×X(k-1|k-1)
式中:X(k|k-1)为k时刻预测状态向量;
F(k|k-1)为k时刻转移状态矩阵;
X(k-1|k-1)为k-1时刻的最优状态向量。
其中,步骤S4中,协方差矩阵在时间上的预测和LU分解方法为:
P(k|k-1)=F(k|k-1)×P(k-1|k-1)×F(k|k-1)T+Q(k|k)=F(k|k-1)×{L(k|k-1)×U(k|k-1)×L(k|k-1)T}×F(k|k-1)T+Q(k|k)
式中:P(k|k-1)为k时刻的预测误差协方差;
F(k|k-1)为k时刻转移状态矩阵;
P(k-1|k-1)为k-1时刻的预测误差协方差;
Q(k|k)为k时刻过程噪声;
L(k|k-1)为P(k|k-1)分解后的上三角矩阵;
L(k|k-1)T为P(k|k-1)分解后的下三角矩阵;
U(k|k-1)为P(k|k-1)分解后的对称矩阵。
其中,步骤S5中,最优增益矩阵为:
K(k|k)=P(k|k-1)×H(k|k)T×{H(k|k)×P(k|k-1)×H(k|k)T×R(k|k)}1=G(k|k)×S(k|k)-1
式中:H(k|k)为观测矩阵;R(k|k)为k时刻测量噪声。
其中,矩阵F(k|k)=L(k|k-1)×L(k|k-1)T×H(k|k)T
矩阵G(k|k)=L(k|k-1)×F(k|k);
矩阵S(k|k)=H(k|k)×G(k|k)+R(k|k)。
其中,步骤S6中,目标状态更新方法为:
X(k|k)=X(k|k-1)+K(k|k)×{Z(k|k)-H(k|k)×X(k|k-1)}=X(k|k-1)+G(k|k)×S(k|k)-1×{Z(k|k)-H(k|k)×X(k|k-1)}
式中:Z(k|k)为k时刻测量信息。
其中,步骤S7中,协方差矩阵更新方法为:
P(k|k)={I-K(k|k)×H(k|k)}×P(k|k-1)=L(k|k-1)×{L(k|k-1)-F(k|k)×S(k|k)-1×F(k|k)T}×L(k|k-1)T
本发明的基于LU分解的卡尔曼滤波目标跟踪方法在已定目标运动模型的标准卡尔曼滤波算法的基础上,引入LU矩阵分解,构造出一种协方差平方根自适应滤波算法,保证协方差矩阵在递推计算过程中始终保持正定性,解决了目标跟踪过程中卡尔曼滤波的数值稳定性问题,具有速度快,精度高,计算量少的特点。
附图说明
此处的附图被并入说明书中并构成本说明书的一部分,示出了符合本发明的实施例,并与说明书一起用于解释本发明的原理。
图1为卡尔曼滤波示意图;
图2为本发明的卡尔曼滤波目标跟踪方法流程图;
图3为计算机舍入误差未采取措施的目标跟踪效果图;
图4为计算机舍入误差采取LU分解的目标跟踪效果图。
具体实施方式
为使本发明实施的目的、技术方案和优点更加清楚,下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行更加详细的描述。
如图1所示为卡尔曼滤波示意图,在雷达跟踪过程中,根据上一个时刻目标状态的最优状态向量X(k-1|k-1)和预测的误差协方差P(k-1|k-1)来预测下一时刻目标状态的预测状态向量X(k|k-1)和预测误差协方差P(k|k-1),但由于计算机的字长有限而引起的舍入误差使得矩阵P(k|k)失去正定性或对称性,从而导致矩阵K(k)计算失真,产生预测误差。
为了解决上述问题,如图2所示,本发明提出了一种基于LU分解的卡尔曼滤波目标跟踪方法,包括如下步骤:
1)雷达发射天线辐射电磁波,接收天线获取当前时刻目标的测量回波信息,上述回波信息包括目标距离、目标方位角、目标俯仰角和目标速度等信息;
2)初始化目标运动模型和测量噪声(系数),以及目标状态信息和状态协方差矩阵;
其中,所述运动模型是描述目标速度、目标距离和目标加速度等变化信息的一组函数,不同的运动模型,其速度与距离的变化具有差别。上述目标运动模型一般包括匀速模型、匀加速模型和转圈模型等。测量噪声是一种描述测量信息的不确定性的干扰。
3)状态量在时间上的预测:
X(k|k-1)=F(k|k-1)×X(k-1|k-1)
4)协方差矩阵在时间上的预测和它的LU分解
P(k|k-1)=F(k|k-1)×P(k-1|k-1)×F(k|k-1)T+Q(k|k)=F(k|k-1)×{L(k|k-1)×U(k|k-1)×L(k|k-1)T}×F(k|k-1)T+Q(k|k)
5)卡尔曼滤波的最优增益矩阵
K(k|k)=P(k|k-1)×H(k|k)T×{H(k|k)×P(k|k-1)×H(k|k)T×R(k|k)}-1=G(k|k)×S(k|k)-1
其中,矩阵F(k|k)=L(k|k-1)×L(k|k-1)T×H(k|k)T
矩阵G(k|k)=L(k|k-1)×F(k|k);
矩阵S(k|k)=H(k|k)×G(k|k)+R(k|k)。
6)利用最优增益矩阵计算最优状态量
X(k|k)=X(k|k-1)+K(k|k)×{Z(k|k)-H(k|k)×X(k|k-1)}=X(k|k-1)+G(k|k)×S(k|k)-1×{Z(k|k)-H(k|k)×X(k|k-1)}
7)最后,更新协方差矩阵
P(k|k)={I-K(k|k)×H(k|k)}×P(k|k-1)=L(k|k-1)×{L(k|k-1)-F(k|k)×S(k|k)-1×F(k|k)T}×L(k|k-1)T
上述公式中的含义:
X(k|k-1)为k时刻预测状态向量;
F(k|k-1)为k时刻转移状态矩阵;
X(k|k)为k时刻的最优状态向量;
X(k-1|k-1)为k-1时刻的最优状态向量;
P(k|k)为k时刻的最优预测误差协方差;
P(k|k-1)为k时刻的预测误差协方差;
P(k-1|k-1)为k-1时刻的预测误差协方差;
K(k|k)为最优增益矩阵;
Q(k|k)为k时刻过程噪声;
L(k|k-1)为P(k|k-1)分解后的上三角矩阵;
L(k|k-1)T为P(k|k-1)分解后的下三角矩阵;
U(k|k-1)为P(k|k-1)分解后的对称矩阵;
H(k|k)为观测矩阵;
R(k|k)为k时刻测量噪声;
Z(k|k)为k时刻测量信息。
根据上述方法,并结合具体案例对本发明的优点作进一步阐述:
首先雷达向地面辐射电磁波,获取当前时刻目标的测量回波集合S,包括目标距离、速度、方位角和俯仰角,假设获得了三个测量回波集合S1[80000,200,1.2,0.1]、S2[80000,205,1.21,0.101]和S3[80000,202,1.23,0.11],集合内元素单位分别对应为m、m/s、rad和rad。
之后,对初始化目标运动模型、测量噪声水平、杂波个数、状态提取门限、探测概率、目标幸存概率、初始化目标状态向量和状态协方差;初始化参数设置参见表1。
表1初始参数设置
参数名称 初始值设置
目标运动模型 匀加速
距离测量噪声水平 X,Y,Z都是50m的标准差值
速度测量噪声水平 X,Y,Z都是2m/s的标准差值
方位角噪声水平 0.0005rad
俯仰角噪声水平 0.00012rad
杂波个数 5个
状态提取门限 Thrnd=0.6
探测概率 DetRatio=0.9
目标幸存概率 TarSurv=0.95
其中,目标距离、目标速度、目标方位角、目标俯仰角和测量噪声水平对应公式中的R矩阵。
通过计算机程序处理之后,对比图3和图4可见,在经历较长睡觉时间后,未采用本发明中的基于LU分解的卡尔曼滤波目标跟踪方法中,后续阶段具有较为明显的发散现象,导致计算结果失真、计算精度降低。而采用了本发明中的基于LU分解的卡尔曼滤波目标跟踪方法中,后续阶段的卡尔曼滤波数值稳定性较好,未发生具有较为明显的发散现象,计算结果未失真、计算精度提高。
综上,本发明的基于LU分解的卡尔曼滤波目标跟踪方法在已定目标运动模型的标准卡尔曼滤波算法的基础上,引入LU矩阵分解,构造出一种协方差平方根自适应滤波算法,保证协方差矩阵在递推计算过程中始终保持正定性,解决了目标跟踪过程中卡尔曼滤波的数值稳定性问题,具有速度快,精度高,计算量少的特点。
本发明最适用于机载多功能雷达目标跟踪的场景,也可以应用于其它平台雷达对目标跟踪的数据处理中。
以上所述,仅为本发明的最优具体实施方式,但本发明的保护范围并不局限于此,任何熟悉本技术领域的技术人员在本发明揭露的技术范围内,可轻易想到的变化或替换,都应涵盖在本发明的保护范围之内。因此,本发明的保护范围应以所述权利要求的保护范围为准。

Claims (1)

1.一种基于LU分解的卡尔曼滤波目标跟踪方法,其特征在于,所述目标跟踪方法包括:
S1:获得当前时刻目标的状态测量信息,所述状态测量信息包括目标距离、目标方位角、目标俯仰角和目标速度;
S2:设置目标的初始参数,所述初始参数信息包括目标运动模型和测量噪声;
S3:根据上一时刻的目标状态及初始参数,预测目标当前时刻的状态信息,其中,目标当前时刻的状态信息的预测方法为:
X(k|k-1)=F(k|k-1)×X(k-1|k-1)
式中:X(k|k-1)为k时刻预测状态向量;
F(k|k-1)为k时刻转移状态矩阵;
X(k-1|k-1)为k-1时刻的最优状态向量;
S4:对上一时刻的协方差矩阵进行LU分解,并根据分解后的协方差矩阵及初始参数,预测当前时刻的协方差矩阵,其中,协方差矩阵在时间上的预测和LU分解的方法为:
P(k|k-1)=F(k|k-1)×P(k-1|k-1)×F(k|k-1)T+Q(k|k)
=F(k|k-1)×{L(k|k-1)×U(k|k-1)×L(k|k-1)T}×F(k|k-1)T+Q(k|k)
式中:P(k|k-1)为k时刻的预测误差协方差;
F(k|k-1)为k时刻转移状态矩阵;
P(k-1|k-1)为k-1时刻的预测误差协方差;
Q(k|k)为k时刻过程噪声;
L(k|k-1)为P(k|k-1)分解后的上三角矩阵;
L(k|k-1)T为P(k|k-1)分解后的下三角矩阵;
U(k|k-1)为P(k|k-1)分解后的对称矩阵;
S5:根据当前时刻目标的状态测量信息与预测的目标状态信息的预测差值,计算最优增益矩阵,最优增益矩阵K(k|k)为:
K(k|k)=P(k|k-1)×H(k|k)T×{H(k|k)×P(k|k-1)×H(k|k)T×R(k|k)}-1
=G(k|k)×S(k|k)-1
式中:H(k|k)为观测矩阵;R(k|k)为k时刻测量噪声;
S6:根据最优增益矩阵更新目标最优状态,其中,目标状态的更新方法为:
X(k|k)=X(k|k-1)+K(k|k)×{Z(k|k)-H(k|k)×X(k|k-1)}
=X(k|k-1)+G(k|k)×S(k|k)-1×{Z(k|k)-H(k|k)×X(k|k-1)}
式中:X(k|k)为k时刻的最优状态向量,Z(k|k)为k时刻测量信息;
S7:对预测的协方差矩阵进行LU分解,并根据最优增益矩阵更新协方差矩阵,其中,协方差矩阵的更新方法为:
P(k|k)={I-K(k|k)×H(k|k)}×P(k|k-1)
=L(k|k-1)×{L(k|k-1)-F(k|k)×S(k|k)-1×F(k|k)T}×L(k|k-1)T
式中,P(k|k)为k时刻的最优预测误差协方差。
CN201810898258.6A 2018-08-08 2018-08-08 一种基于lu分解的卡尔曼滤波目标跟踪方法 Active CN109188422B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810898258.6A CN109188422B (zh) 2018-08-08 2018-08-08 一种基于lu分解的卡尔曼滤波目标跟踪方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810898258.6A CN109188422B (zh) 2018-08-08 2018-08-08 一种基于lu分解的卡尔曼滤波目标跟踪方法

Publications (2)

Publication Number Publication Date
CN109188422A CN109188422A (zh) 2019-01-11
CN109188422B true CN109188422B (zh) 2023-01-06

Family

ID=64920666

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810898258.6A Active CN109188422B (zh) 2018-08-08 2018-08-08 一种基于lu分解的卡尔曼滤波目标跟踪方法

Country Status (1)

Country Link
CN (1) CN109188422B (zh)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112444266B (zh) * 2019-08-27 2023-04-14 北京魔门塔科技有限公司 一种参数的标定方法及装置
CN116781115B (zh) * 2022-03-08 2024-03-26 南方科技大学 扩展目标动态波束跟踪方法、装置及终端设备

Family Cites Families (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP5229126B2 (ja) * 2009-06-17 2013-07-03 日本電気株式会社 目標追尾処理器及びそれに用いる誤差共分散行列の補正方法
US20160013773A1 (en) * 2012-11-06 2016-01-14 Pavel Dourbal Method and apparatus for fast digital filtering and signal processing
CN105352487B (zh) * 2015-10-13 2018-06-15 上海华测导航技术股份有限公司 一种姿态测量系统的精度校准方法
CN106646356B (zh) * 2016-11-23 2019-07-26 西安电子科技大学 一种基于卡尔曼滤波定位的非线性系统状态估计方法
CN107045125B (zh) * 2017-03-17 2020-01-14 电子科技大学 一种基于预测值量测转换的交互多模型雷达目标跟踪方法
CN107506332B (zh) * 2017-07-25 2021-01-29 四川航天系统工程研究所 Kalman滤波器快速实现方法
CN107561503B (zh) * 2017-08-28 2020-08-14 哈尔滨工业大学 一种基于多重渐消因子的自适应目标跟踪滤波方法

Also Published As

Publication number Publication date
CN109188422A (zh) 2019-01-11

Similar Documents

Publication Publication Date Title
CN110503071B (zh) 基于变分贝叶斯标签多伯努利叠加模型的多目标跟踪方法
Granström et al. New prediction for extended targets with random matrices
CN103902812B (zh) 一种粒子滤波方法、装置及目标跟踪方法、装置
CN107688179A (zh) 基于多普勒信息辅助的综合概率数据互联方法
CN110209180B (zh) 一种基于HuberM-Cubature卡尔曼滤波的无人水下航行器目标跟踪方法
CN109188422B (zh) 一种基于lu分解的卡尔曼滤波目标跟踪方法
CN103605886A (zh) 一种船舶动力定位系统多模型自适应融合滤波方法
CN115981148A (zh) 无人机地面移动目标跟踪方法
CN116500575A (zh) 一种基于变分贝叶斯理论的扩展目标跟踪方法和装置
CN111121778B (zh) 一种导航系统初始化方法
CN112328959A (zh) 一种基于自适应扩展卡尔曼概率假设密度滤波器的多目标跟踪方法
CN114236480A (zh) 一种机载平台传感器系统误差配准算法
CN103888100B (zh) 一种基于负熵的非高斯线性随机系统滤波方法
CN116047495B (zh) 一种用于三坐标雷达的状态变换融合滤波跟踪方法
CN116224320B (zh) 一种极坐标系下处理多普勒量测的雷达目标跟踪方法
Ruslan et al. Parameters effect in Sampling Importance Resampling (SIR) particle filter prediction and tracking of flood water level performance
CN115469314A (zh) 一种均匀圆环阵稳健水下目标方位跟踪方法及系统
CN109375160B (zh) 纯方位无源定位中一种测角误差估计方法
CN102707278B (zh) 奇异值分解的多目标跟踪方法
Kumar et al. Performance comparision of α-β-γ filter and kalman filter for CA, NCA target tracking using bistatic range and range rate measurements
CN113009468B (zh) 一种视线坐标系内的解耦cmkf跟踪方法及系统
Ma et al. Variational Bayesian cubature Kalman filter for bearing-only tracking in glint noise environment
CN112491393B (zh) 一种基于观测噪声协方差矩阵未知的线性滤波方法
CN111504326B (zh) 一种基于t分布的鲁棒glmb多目标跟踪方法
CN110208790B (zh) 一种基于mgekf的多传感器目标跟踪方法

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