CN114217283A - 多普勒雷达静态融合平滑变结构滤波方法及装置 - Google Patents

多普勒雷达静态融合平滑变结构滤波方法及装置 Download PDF

Info

Publication number
CN114217283A
CN114217283A CN202111300430.1A CN202111300430A CN114217283A CN 114217283 A CN114217283 A CN 114217283A CN 202111300430 A CN202111300430 A CN 202111300430A CN 114217283 A CN114217283 A CN 114217283A
Authority
CN
China
Prior art keywords
target
state
estimation
pseudo
posterior
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.)
Pending
Application number
CN202111300430.1A
Other languages
English (en)
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.)
Tsinghua University
Naval Aeronautical University
Original Assignee
Tsinghua University
Naval Aeronautical University
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 Tsinghua University, Naval Aeronautical University filed Critical Tsinghua University
Priority to CN202111300430.1A priority Critical patent/CN114217283A/zh
Publication of CN114217283A publication Critical patent/CN114217283A/zh
Pending legal-status Critical Current

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
    • G01S7/00Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/02Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
    • G01S7/41Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
    • 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
    • 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
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/16Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization

Landscapes

  • Engineering & Computer Science (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Mathematical Physics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Computational Mathematics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Theoretical Computer Science (AREA)
  • Computing Systems (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

本申请公开了一种多普勒雷达静态融合平滑变结构滤波方法及装置,其中,方法包括以下步骤:利用第一广义平滑变结构滤波器对目标状态进行更新估计,获取目标状态后验估计;利用第二广义平滑变结构滤波器对目标伪状态进行更新估计,获取目标伪状态后验估计;利用静态融合器对目标状态后验估计和目标伪状态后验估计进行融合处理,获得目标最终的后验状态估计。本申请在保持模型不确定条件下的鲁棒估计性能的前提下,解决了非线性观测模型和非满秩观测矩阵导致的欠定问题,从而提高多普勒雷达系统的目标跟踪性能。可用于多普勒雷达的非线性观测模型,并在目标运动模型不确定条件下实现鲁棒的雷达目标状态估计。

Description

多普勒雷达静态融合平滑变结构滤波方法及装置
技术领域
本申请涉及雷达数据处理技术领域,特别涉及一种多普勒雷达静态融合平滑变结构滤波方法及装置。
背景技术
雷达目标跟踪系统,利用雷达回波中包含的目标距离、方位角、多普勒频偏等信息,对目标数量和位置、速度、加速度等运动状态进行持续有效的估计。雷达跟踪系统广泛应用于自动驾驶、空中交通管制、气象监测、安防、国防军事等领域。其中,多普勒雷达能够观测目标运动的径向速度信息,在雷达目标状态估计环节增加有效信息维度,从而大幅提高目标跟踪性能。
多普勒雷达目标跟踪算法和系统首先需要建立雷达观测模型和目标状态转移模型,它们分别面临两个重要挑战:1)高度非线性的观测模型,包含距离、方位角、俯仰角和径向速度;高度非线性模型会造成卡尔曼滤波器(KF)等常用的雷达跟踪方法产生较大的非线性系统误差,甚至导致滤波发散。2)目标状态转移模型的不确定性;传统贝叶斯滤波方法,包括KF等,依赖于对目标运动模型的精确描述;如果系统状态方程和真实目标的运动模型偏差较大,例如目标发生强机动运动,会导致跟踪误差急剧增大、甚至滤波发散丢失目标。
用于解决问题1)非线性观测模型的方法包括一系列非线性滤波方法,例如扩展卡尔曼滤波(EKF)、不敏卡尔曼滤波(UKF)、粒子滤波(PF)等。特别地,针对多普勒雷达的数学模型,静态融合卡尔曼滤波(SF-CMKF)方法,对非线性的位置量测和径向速度量测解耦处理,分别建立目标状态方程和伪状态方程,利用静态融合方法避免了非线性截断误差在迭代过程中的传递和累积,从而显著提高了估计精度,受到广泛欢迎。但是这些传统的贝叶斯方法,包括静态融合卡尔曼滤波器及其大量的改进形式,都不能解决目标状态转移模型的不确定性问题
用于解决问题2)目标状态转移模型的不确定性问题的方法包括鲁棒卡尔曼滤波,H-infinite滤波,以及近年来备受关注的平滑变结构滤波(SVSF)等方法。特别地,平滑变结构滤波(SVSF)是一种模型不确定条件下的鲁棒滤波方法,基于滑模控制理论和变结构控制方法,能够保证状态估计误差的有界性,同时计算复杂度较低。但值得注意的是,标准的平滑变结构滤波方法要求线性的观测模型和满秩的观测矩阵,从而实现目标状态空间到观测空间的双射映射。为了解决非线性观测模型问题,相关技术试图用雅克比矩阵代替非线性函数代入标准平滑变结构滤波器结构,但这会引入新的欠定矩阵求逆的难题,从而导致目标速度的估计误差很大。因此现有的平滑变结构滤波方法仍然无法应用于多普勒雷达系统。
发明内容
本申请提供一种多普勒雷达静态融合平滑变结构滤波方法及装置,解决了多普勒雷达的非线性欠定观测模型问题和雷达目标的运动模型不确定条件下的鲁棒跟踪问题,可用于多普勒雷达的非线性观测模型,并在目标运动模型不确定条件下实现鲁棒的雷达目标状态估计。
第一方面,本申请实施例提供一种多普勒雷达静态融合平滑变结构滤波方法,包括以下步骤:利用第一广义平滑变结构滤波器对目标状态进行更新估计,获取目标状态后验估计;利用第二广义平滑变结构滤波器对目标伪状态进行更新估计,获取目标伪状态后验估计;利用静态融合器对所述目标状态后验估计和所述目标伪状态后验估计进行融合处理,获得目标最终的后验状态估计。
在本申请的一个实施例中,所述利用第一广义平滑变结构滤波器对目标状态进行更新估计,获取目标状态后验估计,包括:采用无偏量测转换方法将极坐标下或球坐标系下多普勒雷达的非线性位置量测转换到笛卡尔坐标系下,获得线性位置量测及其转换偏差和转换量测误差协方差;利用所述第一广义平滑变结构滤波器和所述线性位置量测及其转换偏差和转换量测误差协方差得到所述目标状态后验估计。
在本申请的一个实施例中,所述利用所述第一广义平滑变结构滤波器和所述线性位置量测及其转换偏差和转换量测误差协方差得到所述目标状态后验估计,包括:计算k+1时刻的目标状态先验估计、雷达量测先验估计和先验状态估计误差协方差阵;计算位置量测新息,并根据所述位置量测新息计算第一混合误差项和第一新息增益项;计算基于目标位置量测的目标状态后验估计及其协方差阵。
在本申请的一个实施例中,所述利用第二广义平滑变结构滤波器对目标伪状态进行更新估计,获取目标伪状态后验估计,包括:构造距离-径向速度积伪状态矢量和伪量测变量;采用无偏量测转换方法计算所述伪量测变量的转换偏差、量测误差方差和协方差;利用所述第二广义平滑变结构滤波器和所述伪量测变量的转换偏差、量测误差方差和协方差得到所述目标伪状态后验估计。
在本申请的一个实施例中,所述利用所述第二广义平滑变结构滤波器和所述伪量测变量的转换偏差、量测误差方差和协方差得到所述目标伪状态后验估计,包括:计算k+1时刻的目标伪状态先验估计、雷达伪量测先验估计和伪状态先验估计协方差阵:计算伪量测新息,并根据所述伪量测新息计算第二混合误差项和第二新息增益项;计算目标伪状态的后验估计及其协方差阵。
在本申请的一个实施例中,利用所述第一广义平滑变结构滤波器输出的所述目标状态后验估计及其协方差阵和所述第二广义平滑变结构滤波器输出的所述目标伪状态后验估计及其协方差阵,采用基于线性最小均方误差准则的静态融合器计算目标最终的后验状态估计及其协方差。
第二方面,本申请实施例提供一种多普勒雷达静态融合平滑变结构滤波装置,包括:第一估计模块,用于利用第一广义平滑变结构滤波器对目标状态进行更新估计,获取目标状态后验估计;第二估计模块,用于利用第二广义平滑变结构滤波器对目标伪状态进行更新估计,获取目标伪状态后验估计;融合模块,用于利用静态融合器对所述目标状态后验估计和所述目标伪状态后验估计进行融合处理,获得目标最终的后验状态估计。
第三方面,本申请实施例提供一种电子设备,包括:存储器、处理器及存储在所述存储器上并可在所述处理器上运行的计算机程序,所述处理器执行所述程序,以实现如上述实施例所述的多普勒雷达静态融合平滑变结构滤波方法。
第四方面,本申请实施例提供一种计算机可读存储介质,所述计算机可读存储介质存储计算机指令,所述计算机指令用于使所述计算机执行如上述实施例所述的多普勒雷达静态融合平滑变结构滤波方法。
本申请实施例的多普勒雷达静态融合平滑变结构滤波方法及装置,能够显著提高多普勒雷达目标状态估计的精度,保证模型不确定条件下的滤波鲁棒性,具有很好的应用价值,具体为:
1)目标跟踪的鲁棒性好
本申请的静态融合平滑变结构滤波方法,基于滑模变结构控制理论,能够保证目标运动模型未知条件下的跟踪误差有界性,解决了传统的非线性贝叶斯滤波器(例如卡尔曼滤波器和静态融合卡尔曼滤波器等)面对目标状态转移模型存在建模误差时跟踪滤波误差急剧增大、甚至滤波发散的问题。在第一广义平滑变结构滤波器(CPM-SVSF)中,本申请的方法采用无偏转换方法将多普勒雷达的位置量测转换到笛卡尔坐标系下,从而可以构造线性形式的经典平滑变结构滤波器,实现对目标状态的鲁棒估计;在第二广义平滑变结构滤波器(CDM-SVSF)中,本申请利用径向速度及其导数构造线性的伪状态方程和伪量测方程,从而利用线性形式的经典平滑变结构滤波器,实现对目标伪状态的鲁棒估计;在静态融合处理时,采用最小均方误差准则策略将CPM-SVSF输出的目标状态的鲁棒估计和CDM-SVSF输出的目标伪状态的鲁棒估计进行融合,获得目标最终的状态估计,从而保证了模型不确定条件下的鲁棒状态估计性能。
2)目标跟踪的精度高
本申请所提出的静态融合平滑变结构滤波方法能够显著提升跟踪精度,这来自于两方面的优势:
第一,本申请所提出的静态融合平滑变结构滤波方法,解决了现有的平滑变结构滤波方法面临的多普勒雷达非线性模型、非满秩观测矩阵的欠定问题,从而充分利用多普勒雷达的径向速度量测信息有效提高了平滑变结构滤波器的目标状态估计精度。具体而言,静态融合平滑变结构滤波器解决了标准的平滑变结构滤波器无法利用非线性的多普勒量测的问题,也解决了现有文献中改进的平滑变结构滤波器采用雅克比矩阵局部线性化操作时面临的非满秩观测矩阵求逆的欠定性难题。因此,本申请的方法相比于现有的平滑变结构滤波方法更有效地利用了多普勒速度信息,从而显著提高了目标状态跟踪的精度。
第二,本申请所提出的静态融合平滑变结构滤波方法,采用位置量测和多普勒量测解耦处理、并行鲁棒滤波后静态融合的方法,避免了非线性量测导致的截断误差在迭代估计过程中的传递和累积。具体来说,为了降低非线性误差的影响,第一广义平滑变结构滤波器(CPM-SVSF)和第二广义平滑变结构滤波器(CDM-SVSF)均建模为线性的状态/伪状态方程和量测/伪量测方程,迭代估计过程仅在线性的子滤波器内部完成,不存在非线性系统误差的传递和累积过程。因此,本申请的方法相比于现有的非线性平滑变结构滤波方法有效降低了非线性系统误差的影响,从而显著提高了目标状态的跟踪精度。
本申请附加的方面和优点将在下面的描述中部分给出,部分将从下面的描述中变得明显,或通过本申请的实践了解到。
附图说明
本申请上述的和/或附加的方面和优点从下面结合附图对实施例的描述中将变得明显和容易理解,其中:
图1为根据本申请实施例提供的一种多普勒雷达静态融合平滑变结构滤波方法的流程示意图;
图2为根据本申请实施例提供的一种多普勒雷达静态融合平滑变结构滤波方法流程逻辑框图;
图3为根据本申请实施例提供的一种自动驾驶场景下2-D多普勒毫米波雷达对车辆目标跟踪的仿真场景实例图;
图4为根据本申请实施例提供的一种仿真场景下本申请方法的目标状态估计误差随时间变化示意图;图4的(a)为坐标估计误差,图4的(b)为速度估计误差;
图5为根据本申请实施例的多普勒雷达静态融合平滑变结构滤波装置的示例图;
图6为根据本申请实施例的一种电子设备示意图。
具体实施方式
下面详细描述本申请的实施例,所述实施例的示例在附图中示出,其中自始至终相同或类似的标号表示相同或类似的元件或具有相同或类似功能的元件。下面通过参考附图描述的实施例是示例性的,旨在用于解释本申请,而不能理解为对本申请的限制。
图1为根据本申请实施例提供的一种多普勒雷达静态融合平滑变结构滤波方法的流程示意图。
在本申请的实施例中,静态融合平滑变结构滤波器(SF-SVSF)包含两个并行的子滤波器(第一广义平滑变结构滤波器和第二广义平滑变结构滤波器)和一个静态融合器。第一广义平滑变结构滤波器(称为CPM-SVSF)中,先利用多普勒雷达的非线性位置量测(包含目标距离、方位角、俯仰角)无偏转换到笛卡尔坐标系下,再利用笛卡尔坐标系下的线性位置量测构造线性的第一广义平滑变结构滤波器,迭代更新基于位置量测的目标状态估计。第二广义平滑变结构滤波器(称为CDM-SVSF)中,先利用多普勒雷达的径向速度量测构造距离-径向速度积伪状态方程和伪观测方程,在此基础上利用线性的第二广义平滑变结构滤波器迭代更新伪状态估计。静态融合器对目标状态后验估计和目标伪状态后验估计进行融合处理,获得目标最终的状态估计。
如图1所示,该多普勒雷达静态融合平滑变结构滤波方法包括以下步骤:
在步骤S101中,利用第一广义平滑变结构滤波器对目标状态进行更新估计,获取目标状态后验估计。
在本申请的一个实施例中,利用第一广义平滑变结构滤波器对目标状态进行更新估计,获取目标状态后验估计,包括:采用无偏量测转换方法将极坐标下或球坐标系下多普勒雷达的非线性位置量测转换到笛卡尔坐标系下,获得线性位置量测及其转换偏差和转换量测误差协方差;利用第一广义平滑变结构滤波器和线性位置量测及其转换偏差和转换量测误差协方差得到目标状态后验估计。
具体地,利用第一广义平滑变结构滤波器(CPM-SVSF)对目标状态进行更新估计,首先采用无偏量测转换的方法将极坐标下或球坐标系下多普勒雷达的目标位置量测转换到笛卡尔坐标系下,然后利用广义平滑变结构滤波器和无偏转换的位置量测获得目标状态估计。具体步骤为:
步骤1-1:采用无偏量测转换(UCM)方法,将多普勒雷达球坐标系或极坐标下的非线性位置量测
Figure BDA0003338160690000061
(分别表示距离、方位角和俯仰角)转换为笛卡尔坐标系下的线性位置量测
Figure BDA0003338160690000062
及其转换偏差μp=[μxyz]T和转换量测误差协方差Rp(k)。其中,UCM方法可参考相关技术文献中的方法,不作具体限定。
在本申请的一个实施例中,利用第一广义平滑变结构滤波器和线性位置量测及其转换偏差和转换量测误差协方差得到目标状态后验估计,包括:计算k+1时刻的目标状态先验估计、雷达量测先验估计和先验状态估计误差协方差阵;计算位置量测新息,并根据位置量测新息计算第一混合误差项和第一新息增益项;计算基于目标位置量测的目标状态后验估计及其协方差阵。
具体地,将上述的无偏转换后的目标位置量测
Figure BDA0003338160690000063
作为输入,采用如下的广义平滑变结构滤波器获得目标运动状态估计,即目标运动状态估计矢量
Figure BDA0003338160690000064
和估计误差协方差矩阵Pp(k|k)。具体步骤如下:
步骤1-2:计算k+1时刻的目标状态先验估计
Figure BDA0003338160690000065
雷达量测先验估计
Figure BDA0003338160690000066
和先验状态估计误差协方差阵Pp(k|k-1):
Figure BDA0003338160690000067
Figure BDA0003338160690000068
Figure BDA0003338160690000069
其中
Figure BDA00033381606900000610
是包含三维位置、速度的状态估计变量,
Figure BDA00033381606900000611
是预设的状态转移矩阵,
Figure BDA00033381606900000612
是观测矩阵,Q和Γp是过程噪声协方差及其系数矩阵。
Figure BDA00033381606900000613
和Γp分别定义为:
Figure BDA00033381606900000614
需要说明的是,不同于传统贝叶斯滤波器(例如SF-CMKF)要求准确已知的目标状态转移模型Fp,这里的
Figure BDA00033381606900000615
相比于真实值允许存在范数有界的模型误差,即本申请方法具有模型不确定条件下的鲁棒性优势;本申请引入了无偏量测转换的步骤,因此
Figure BDA00033381606900000616
是线性时不变的,优于直接使用非线性位置量测的传统的平滑变结构滤波器。
步骤1-3:计算位置量测新息:
Figure BDA0003338160690000071
计算第一混合误差项:
Figure BDA0003338160690000072
Figure BDA0003338160690000073
其中ep(k-1|k-1)表示子滤波器1上一时刻的后验量测误差,·ABS表示对矢量逐元素取绝对值,γp,z∈i3×1和γp,y∈i3×1是取值在(0,1)区间的衰减因子;
Figure BDA0003338160690000074
表示主对角元素为矢量a的对角阵;
Figure BDA0003338160690000075
其中Tp是变换矩阵,在本方法中可选择为单位阵。
计算基于位置量测的第一新息增益项Kp(k):
Figure BDA0003338160690000076
Figure BDA0003338160690000077
Figure BDA0003338160690000078
其中sat(.)表示饱和函数;Ψp,z∈i3×1和Ψp,y∈i3×1分别是预设的平滑层参数矢量;Hp,1=I3×3是步骤1-2中
Figure BDA0003338160690000079
的满秩分块子矩阵,是线性时不变的:
步骤:1-4:计算目标状态的后验估计
Figure BDA00033381606900000710
和状态协方差阵Pp(k|k):
Figure BDA00033381606900000711
Figure BDA00033381606900000712
其中,
Figure BDA00033381606900000714
是单位阵,Rp(k)是步骤1-1中获得的位置量测误差协方差。
第一广义平滑变结构滤波器(CPM-SVSF)输出目标状态估计
Figure BDA00033381606900000713
及其协方差阵Pp(k|k)。
步骤S102,利用第二广义平滑变结构滤波器对目标伪状态进行更新估计,获取目标伪状态后验估计。
在本申请的实施例中,利用第二广义平滑变结构滤波器对目标伪状态进行更新估计,获取目标伪状态后验估计,包括:构造距离-径向速度积伪状态矢量和伪量测变量;采用无偏量测转换方法计算伪量测变量的转换偏差、量测误差方差和协方差;利用第二广义平滑变结构滤波器和伪量测变量的转换偏差、量测误差方差和协方差得到目标伪状态后验估计。
对第二广义平滑变结构滤波器进行状态更新,首先,定义伪状态矢量:
Figure BDA0003338160690000081
其中,标量η(k)是距离-径向速度积,定义为
Figure BDA0003338160690000082
相应的伪量测定义为
Figure BDA0003338160690000083
采用无偏量测转换(UCM)方法计算伪量测变量的转换偏差μη,量测误差方差Rη和协方差R;计算方法可以使用相关技术中的方法,不作具体限定。
在本申请的实施例中,利用第二广义平滑变结构滤波器和伪量测变量的转换偏差、量测误差方差和协方差得到目标伪状态后验估计,包括:计算k+1时刻的目标伪状态先验估计、雷达伪量测先验估计和伪状态先验估计协方差阵:计算伪量测新息,并根据伪量测新息计算第二混合误差项和第二新息增益项;计算目标伪状态的后验估计及其协方差阵。
具体地,将无偏转换后的伪量测zη作为输入,采用如下的第二广义平滑变结构滤波器获得目标伪状态估计
Figure BDA0003338160690000084
和估计误差协方差矩阵Pη(k|k)。具体步骤如下:
步骤2-1:计算k+1时刻的目标伪状态先验估计
Figure BDA0003338160690000085
雷达伪量测先验估计
Figure BDA0003338160690000086
和伪状态先验估计协方差阵Pη(k|k-1):
Figure BDA0003338160690000087
Figure BDA0003338160690000088
Figure BDA0003338160690000089
其中
Figure BDA00033381606900000810
是伪状态估计变量,
Figure BDA00033381606900000811
是预设的伪状态转移矩阵,Gη和uη是过程噪声高阶量引入的输入项系数矩阵和输入项,Hη是伪量测矩阵,Qx和Γx是与目标状态相关而引入的过程噪声协方差及其系数矩阵,Qs和Γs是与过程噪声高阶量有关的过程噪声协方差及其系数矩阵。分别定义如下为:
Figure BDA0003338160690000091
Figure BDA0003338160690000092
其中q是笛卡尔坐标系下的过程噪声(高斯白噪声)的方差,T是雷达扫描周期。
需要说明的是,不同于传统贝叶斯滤波器(例如SF-CMKF)要求准确已知的目标伪状态转移模型Fη,这里的
Figure BDA0003338160690000093
相比于真实值允许存在范数有界的模型误差,即本申请的方法具有模型不确定条件下的鲁棒性优势;本申请方法引入了无偏量测转换的步骤,因此Hη(k)=[1,0]是线性时不变的,优于直接使用非线性多普勒量测的传统的平滑变结构滤波器。
步骤2-2:计算伪量测新息:
Figure BDA0003338160690000094
计算第二混合误差项:
Eη,z=|eη(k|k-1)|+γη,z|eη(k-1|k-1)|
Figure BDA0003338160690000095
其中eη(k-1|k-1)表示第二广义平滑变结构滤波器上一时刻的后验量测误差,γη,z和γη,y是取值在(0,1)区间的衰减因子;
Figure BDA0003338160690000096
其中Tη是变换矩阵,在本申请方法中可选择为单位阵。
计算基于伪量测的第二新息增益项Kη(k):
Figure BDA0003338160690000097
Kη,u(k)=Hη,1 -1·sat(eη(k|k-1),Ψη,z)·Eη,z·ep(k|k-1)-1
Figure BDA0003338160690000098
其中sat(.)表示饱和函数;Ψη,z和Ψη,y分别是预设的平滑层参数矢量;Hη,1=1是步骤2-1中Hη的满秩分块子矩阵,是线性时不变的。
步骤2-3:计算目标状态的后验估计
Figure BDA0003338160690000101
和状态协方差阵Pη(k|k):
Figure BDA0003338160690000102
Pη(k|k)=[I-Kη(k)Hη(k)]Pη(k|k-1)[I-Kη(k)Hη(k)]T+Kη(k)Rη(k)Kη(k)T
其中,
Figure BDA0003338160690000109
是单位阵,Rη(k)是上述步骤获得的伪量测方差。
第二广义平滑变结构滤波器(CDM-SVSF)输出伪状态估计
Figure BDA0003338160690000103
及其协方差Pη(k|k)。
在步骤S103中,利用静态融合器对目标状态后验估计和目标伪状态后验估计进行融合处理,获得目标最终的后验状态估计。
在本申请的实施例中,利用第一广义平滑变结构滤波器输出的目标状态后验估计及其协方差阵和第二广义平滑变结构滤波器输出的目标伪状态后验估计及其协方差阵,采用基于线性最小均方误差准则的静态融合器计算目标最终的后验状态估计及其协方差。
静态融合估计器是一个局部近似的线性最小均方误差估计器,利用上述两个广义平滑变结构滤波器的输出融合处理,得到当前帧目标状态的最终估计,利用上述过程迭代更新以实现对目标状态的持续跟踪。
具体地,将第一广义平滑变结构滤波器(CPM-SVSF)和第二广义平滑变结构滤波器(CDM-SVSF)的输出作为输入,利用静态融合器获得当前时刻目标运动状态的最终估计
Figure BDA0003338160690000104
和P(k|k),计算方法如下:
Figure BDA0003338160690000105
Figure BDA0003338160690000106
其中,
Figure BDA0003338160690000107
和Pp(k|k)是第一广义平滑变结构滤波器(CPM-SVSF)输出的状态估计及其协方差,
Figure BDA0003338160690000108
是第一广义平滑变结构滤波器(CDM-SVSF)输出的伪状态估计。伪状态均值E(η(k)),静态融合器的互相关矩阵Pxz和自相关矩阵Pzz的计算方法可以利用相关技术中的计算方法,不作具体限定。
对第k=1,2,…次雷达扫描回波,重复步骤S101至S103的迭代估计过程,直至多普勒雷达系统判定航迹结束,从而实现模型不确定条件下多普勒雷达目标状态的非线性鲁棒估计。
下面结合附图和具体实施例对本申请实施例的多普勒雷达静态融合平滑变结构滤波方法进行说明。
如图2所示,为多普勒雷达静态融合平滑变结构滤波方法的执行步骤,其中,子滤波器1为第一广义平滑变结构滤波器,子滤波器2为第二广义平滑变结构滤波器。机动目标跟踪是模型不确定问题的典型场景,考虑自动驾驶场景下2-D多普勒毫米波雷达对单机动目标的跟踪滤波的仿真场景。雷达建立在坐标原点处,目标机动轨迹如图3所示。本申请实施例拟比较本申请提出的静态融合平滑变结构滤波方法(记为SF-SVSF),和扩展卡尔曼滤波方法(EKF),仅利用位置量测的平滑变结构滤波方法(PO-SVSF),现有文献中采用雅克比矩阵局部线性化的平滑变结构滤波方法(记为L-SVSF)。采用2000次蒙特卡洛实验的平均结果作为对比数据。
利用表1的仿真参数生成仿真的雷达目标数据。目标状态变量定义为X-Y坐标轴内的位置、速度,即x=[x y vx vy]T;目标伪状态变量定义为
Figure BDA0003338160690000111
雷达量测变量定义
Figure BDA0003338160690000112
在子滤波器1(CPM-SVSF)中采用二维匀速运动模型(CV)对目标状态进行估计,即:
Figure BDA0003338160690000113
子滤波器2(CDM-SVSF)的相关设置通过上述步骤记载的方法进行设置。
表1雷达目标跟踪场景参数
Figure BDA0003338160690000114
按照上述实施例记载的步骤进行目标状态的迭代估计。
图4给出了目标状态估计的坐标误差和速度误差随扫描时间的变化情况,其中,图4的(b)中L-SVSF方法面临欠定矩阵求逆的数值稳定性问题,导致速度估计误差超出展示范围。表2给出了2000次蒙特卡洛实验的方均根误差(RMSE)。由图4可知,本申请提出的静态融合平滑变结构滤波(SF-SVSF)方法取得了最好的跟踪精度。相比于扩展卡尔曼滤波(EKF),SF-SVSF方法在目标机动运动时跟踪性能更加鲁棒(10-15s和23-28s,峰值坐标误差降低了超过一半)。相比于只利用位置量测的PO-SVSF,所提出的SF-SVSF能够有效利用多普勒速度信息,大幅提高了目标状态估计精度。现有文献的局部线性化L-SVSF则面临非满秩观测矩阵引起的欠定问题,目标状态估计性能(尤其是速度估计性能)最差。综上,本申请的方法解决于多普勒雷达的非线性欠定观测模型问题,在目标运动模型不确定条件下,能够充分利用多普勒速度量测信息,实现更加鲁棒和精确的雷达目标状态估计。
表2运动模型不确定目标的状态估计方均根误差
Figure BDA0003338160690000121
根据本申请实施例提出的多普勒雷达静态融合平滑变结构滤波方法,将状态空间静态融合策略与平滑变结构滤波方法相结合,解决了标准的平滑变结构滤波器无法应用于非线性、非满秩观测矩阵的多普勒雷达系统的问题,同时解决了传统的非线性贝叶斯滤波器在目标运动模型不确定条件下容易滤波发散的问题;此外,目标位置量测和多普勒量测解耦处理、并行滤波和静态融合的结构能够避免非线性截断误差在迭代过程中的传递和累积,从而进一步降低估计误差。因此,本申请能够有效利用多普勒观测信息提高目标状态的跟踪精度,并且在目标运动模型不确定情况下保持鲁棒估计性能,具有较高的工程应用价值。
其次参照附图描述根据本申请实施例提出的多普勒雷达静态融合平滑变结构滤波装置。
图5是本申请实施例的多普勒雷达静态融合平滑变结构滤波装置的方框示意图。
如图5所示,该多普勒雷达静态融合平滑变结构滤波装置10包括:第一估计模块100、第二估计模块200和融合模块300。
其中,第一估计模块100,用于利用第一广义平滑变结构滤波器对目标状态进行更新估计,获取目标状态后验估计。第二估计模块200,用于利用第二广义平滑变结构滤波器对目标伪状态进行更新估计,获取目标伪状态后验估计。融合模块300,用于利用静态融合器对目标状态后验估计和目标伪状态后验估计进行融合处理,获得目标最终的后验状态估计。
在本申请的一个实施例中,第一估计模块100,进一步用于,采用无偏量测转换方法将极坐标下或球坐标系下多普勒雷达的非线性位置量测转换到笛卡尔坐标系下,获得线性位置量测及其转换偏差和转换量测误差协方差;利用第一广义平滑变结构滤波器和线性位置量测及其转换偏差和转换量测误差协方差得到目标状态后验估计。
在本申请的一个实施例中,第一估计模块100,具体用于,计算k+1时刻的目标状态先验估计、雷达量测先验估计和先验状态估计误差协方差阵;计算位置量测新息,并根据位置量测新息计算第一混合误差项和第一新息增益项;计算基于目标位置量测的目标状态后验估计及其协方差阵。
在本申请的一个实施例中,第二估计模块200,进一步用于,构造距离-径向速度积伪状态矢量和伪量测变量;采用无偏量测转换方法计算伪量测变量的转换偏差、量测误差方差和协方差;利用第二广义平滑变结构滤波器和伪量测变量的转换偏差、量测误差方差和协方差得到目标伪状态后验估计。
在本申请的一个实施例中,第二估计模块200,具体用于,计算k+1时刻的目标伪状态先验估计、雷达伪量测先验估计和伪状态先验估计协方差阵;计算伪量测新息,并根据伪量测新息计算第二混合误差项和第二新息增益项;计算目标伪状态后验估计及其协方差阵。
在本申请的一个实施例中,融合模块300,进一步用于,利用第一广义平滑变结构滤波器输出的目标状态后验估计及其协方差阵和第二广义平滑变结构滤波器输出的目标伪状态后验估计及其协方差阵,采用基于线性最小均方误差准则的静态融合器计算目标最终的后验状态估计及其协方差。
需要说明的是,前述对多普勒雷达静态融合平滑变结构滤波方法实施例的解释说明也适用于该实施例的多普勒雷达静态融合平滑变结构滤波装置,此处不再赘述。
根据本申请实施例提出的多普勒雷达静态融合平滑变结构滤波装置,将状态空间静态融合策略与平滑变结构滤波方法相结合,解决了标准的平滑变结构滤波器无法应用于非线性、非满秩观测矩阵的多普勒雷达系统的问题,同时解决了传统的非线性贝叶斯滤波器在目标运动模型不确定条件下容易滤波发散的问题;此外,目标位置量测和多普勒量测解耦处理、并行滤波和静态融合的结构能够避免非线性截断误差在迭代过程中的传递和累积,从而进一步降低估计误差。因此,本申请能够有效利用多普勒观测信息提高目标状态的跟踪精度,并且在目标运动模型不确定情况下保持鲁棒估计性能,具有较高的工程应用价值。
图6为本申请实施例提供的电子设备的结构示意图。该电子设备可以包括:
存储器601、处理器602及存储在存储器601上并可在处理器602上运行的计算机程序。
处理器602执行程序时实现上述实施例中提供的多普勒雷达静态融合平滑变结构滤波方法。
进一步地,电子设备还包括:
通信接口603,用于存储器601和处理器602之间的通信。
存储器601,用于存放可在处理器602上运行的计算机程序。
存储器601可能包含高速RAM存储器,也可能还包括非易失性存储器(non-volatile memory),例如至少一个磁盘存储器。
如果存储器601、处理器602和通信接口603独立实现,则通信接口603、存储器601和处理器602可以通过总线相互连接并完成相互间的通信。总线可以是工业标准体系结构(Industry Standard Architecture,简称为ISA)总线、外部设备互连(PeripheralComponent,简称为PCI)总线或扩展工业标准体系结构(Extended Industry StandardArchitecture,简称为EISA)总线等。总线可以分为地址总线、数据总线、控制总线等。为便于表示,图6中仅用一条粗线表示,但并不表示仅有一根总线或一种类型的总线。
在具体实现上,如果存储器601、处理器602及通信接口603,集成在一块芯片上实现,则存储器601、处理器602及通信接口603可以通过内部接口完成相互间的通信。
处理器602可能是一个中央处理器(Central Processing Unit,简称为CPU),或者是特定集成电路(Application Specific Integrated Circuit,简称为ASIC),或者是被配置成实施本申请实施例的一个或多个集成电路。
本实施例还提供一种计算机可读存储介质,其上存储有计算机程序,该程序被处理器执行时实现如上的多普勒雷达静态融合平滑变结构滤波方法。
在本说明书的描述中,参考术语“一个实施例”、“一些实施例”、“示例”、“具体示例”、或“一些示例”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本申请的至少一个实施例或示例中。在本说明书中,对上述术语的示意性表述不必须针对的是相同的实施例或示例。而且,描述的具体特征、结构、材料或者特点可以在任一个或N个实施例或示例中以合适的方式结合。此外,在不相互矛盾的情况下,本领域的技术人员可以将本说明书中描述的不同实施例或示例以及不同实施例或示例的特征进行结合和组合。
此外,术语“第一”、“第二”仅用于描述目的,而不能理解为指示或暗示相对重要性或者隐含指明所指示的技术特征的数量。由此,限定有“第一”、“第二”的特征可以明示或者隐含地包括至少一个该特征。在本申请的描述中,“N个”的含义是至少两个,例如两个,三个等,除非另有明确具体的限定。
流程图中或在此以其他方式描述的任何过程或方法描述可以被理解为,表示包括一个或更N个用于实现定制逻辑功能或过程的步骤的可执行指令的代码的模块、片段或部分,并且本申请的优选实施方式的范围包括另外的实现,其中可以不按所示出或讨论的顺序,包括根据所涉及的功能按基本同时的方式或按相反的顺序,来执行功能,这应被本申请的实施例所属技术领域的技术人员所理解。
在流程图中表示或在此以其他方式描述的逻辑和/或步骤,例如,可以被认为是用于实现逻辑功能的可执行指令的定序列表,可以具体实现在任何计算机可读介质中,以供指令执行系统、装置或设备(如基于计算机的系统、包括处理器的系统或其他可以从指令执行系统、装置或设备取指令并执行指令的系统)使用,或结合这些指令执行系统、装置或设备而使用。就本说明书而言,"计算机可读介质"可以是任何可以包含、存储、通信、传播或传输程序以供指令执行系统、装置或设备或结合这些指令执行系统、装置或设备而使用的装置。计算机可读介质的更具体的示例(非穷尽性列表)包括以下:具有一个或N个布线的电连接部(电子装置),便携式计算机盘盒(磁装置),随机存取存储器(RAM),只读存储器(ROM),可擦除可编辑只读存储器(EPROM或闪速存储器),光纤装置,以及便携式光盘只读存储器(CDROM)。另外,计算机可读介质甚至可以是可在其上打印所述程序的纸或其他合适的介质,因为可以例如通过对纸或其他介质进行光学扫描,接着进行编辑、解译或必要时以其他合适方式进行处理来以电子方式获得所述程序,然后将其存储在计算机存储器中。
应当理解,本申请的各部分可以用硬件、软件、固件或它们的组合来实现。在上述实施方式中,N个步骤或方法可以用存储在存储器中且由合适的指令执行系统执行的软件或固件来实现。如,如果用硬件来实现和在另一实施方式中一样,可用本领域公知的下列技术中的任一项或他们的组合来实现:具有用于对数据信号实现逻辑功能的逻辑门电路的离散逻辑电路,具有合适的组合逻辑门电路的专用集成电路,可编程门阵列(PGA),现场可编程门阵列(FPGA)等。
本技术领域的普通技术人员可以理解实现上述实施例方法携带的全部或部分步骤是可以通过程序来指令相关的硬件完成,所述的程序可以存储于一种计算机可读存储介质中,该程序在执行时,包括方法实施例的步骤之一或其组合。
此外,在本申请各个实施例中的各功能单元可以集成在一个处理模块中,也可以是各个单元单独物理存在,也可以两个或两个以上单元集成在一个模块中。上述集成的模块既可以采用硬件的形式实现,也可以采用软件功能模块的形式实现。所述集成的模块如果以软件功能模块的形式实现并作为独立的产品销售或使用时,也可以存储在一个计算机可读取存储介质中。
上述提到的存储介质可以是只读存储器,磁盘或光盘等。尽管上面已经示出和描述了本申请的实施例,可以理解的是,上述实施例是示例性的,不能理解为对本申请的限制,本领域的普通技术人员在本申请的范围内可以对上述实施例进行变化、修改、替换和变型。

Claims (9)

1.一种多普勒雷达静态融合平滑变结构滤波方法,其特征在于,包括以下步骤:
利用第一广义平滑变结构滤波器对目标状态进行更新估计,获取目标状态后验估计;
利用第二广义平滑变结构滤波器对目标伪状态进行更新估计,获取目标伪状态后验估计;
利用静态融合器对所述目标状态后验估计和所述目标伪状态后验估计进行融合处理,获得目标最终的后验状态估计。
2.根据权利要求1所述的方法,其特征在于,所述利用第一广义平滑变结构滤波器对目标状态进行更新估计,获取目标状态后验估计,包括:
采用无偏量测转换方法将极坐标下或球坐标系下多普勒雷达的非线性位置量测转换到笛卡尔坐标系下,获得线性位置量测及其转换偏差和转换量测误差协方差;
利用所述第一广义平滑变结构滤波器和所述线性位置量测及其转换偏差和转换量测误差协方差得到所述目标状态后验估计。
3.根据权利要求2所述的方法,其特征在于,利用所述第一广义平滑变结构滤波器和所述线性位置量测及其转换偏差和转换量测误差协方差得到所述目标状态后验估计,包括:
计算k+1时刻的目标状态先验估计、雷达量测先验估计和先验状态估计误差协方差阵;
计算位置量测新息,并根据所述位置量测新息计算第一混合误差项和第一新息增益项;
计算基于目标位置量测的所述目标状态后验估计及其协方差阵。
4.根据权利要求1所述的方法,其特征在于,利用所述第二广义平滑变结构滤波器对目标伪状态进行更新估计,获取目标伪状态后验估计,包括:
构造距离-径向速度积伪状态矢量和伪量测变量;
采用无偏量测转换方法计算所述伪量测变量的转换偏差、量测误差方差和协方差;
利用所述第二广义平滑变结构滤波器和所述伪量测变量的转换偏差、量测误差方差和协方差得到所述目标伪状态后验估计。
5.根据权利要求4所述的方法,其特征在于,利用所述第二广义平滑变结构滤波器和所述伪量测变量的转换偏差、量测误差方差和协方差得到所述目标伪状态后验估计,包括:
计算k+1时刻的目标伪状态先验估计、雷达伪量测先验估计和伪状态先验估计协方差阵;
计算伪量测新息,并根据所述伪量测新息计算第二混合误差项和第二新息增益项;
计算所述目标伪状态后验估计及其协方差阵。
6.根据权利要求1-5任一项所述的方法,其特征在于,利用所述第一广义平滑变结构滤波器输出的所述目标状态后验估计及其协方差阵和所述第二广义平滑变结构滤波器输出的所述目标伪状态后验估计及其协方差阵,采用基于线性最小均方误差准则的静态融合器计算目标最终的后验状态估计及其协方差。
7.一种多普勒雷达静态融合平滑变结构滤波装置,其特征在于,包括:
第一估计模块,用于利用第一广义平滑变结构滤波器对目标状态进行更新估计,获取目标状态后验估计;
第二估计模块,用于利用第二广义平滑变结构滤波器对目标伪状态进行更新估计,获取目标伪状态后验估计;以及
融合模块,用于利用静态融合器对所述目标状态后验估计和所述目标伪状态后验估计进行融合处理,获得目标最终的后验状态估计。
8.一种电子设备,其特征在于,包括:存储器、处理器及存储在所述存储器上并可在所述处理器上运行的计算机程序,所述处理器执行所述程序,以实现如权利要求1-6任一项所述的多普勒雷达静态融合平滑变结构滤波方法。
9.一种计算机可读存储介质,其上存储有计算机程序,其特征在于,该程序被处理器执行,以用于实现如权利要求1-6任一项所述的多普勒雷达静态融合平滑变结构滤波方法。
CN202111300430.1A 2021-11-04 2021-11-04 多普勒雷达静态融合平滑变结构滤波方法及装置 Pending CN114217283A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202111300430.1A CN114217283A (zh) 2021-11-04 2021-11-04 多普勒雷达静态融合平滑变结构滤波方法及装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202111300430.1A CN114217283A (zh) 2021-11-04 2021-11-04 多普勒雷达静态融合平滑变结构滤波方法及装置

Publications (1)

Publication Number Publication Date
CN114217283A true CN114217283A (zh) 2022-03-22

Family

ID=80696672

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202111300430.1A Pending CN114217283A (zh) 2021-11-04 2021-11-04 多普勒雷达静态融合平滑变结构滤波方法及装置

Country Status (1)

Country Link
CN (1) CN114217283A (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117406590A (zh) * 2023-10-08 2024-01-16 哈尔滨工业大学 一种基于平滑变结构滤波的鲁棒机动目标跟踪方法及系统
CN117521018A (zh) * 2024-01-08 2024-02-06 鹏城实验室 基于扩展观测的融合估计方法、装置、设备及存储介质

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117406590A (zh) * 2023-10-08 2024-01-16 哈尔滨工业大学 一种基于平滑变结构滤波的鲁棒机动目标跟踪方法及系统
CN117406590B (zh) * 2023-10-08 2024-04-02 哈尔滨工业大学 一种基于平滑变结构滤波的鲁棒机动目标跟踪方法及系统
CN117521018A (zh) * 2024-01-08 2024-02-06 鹏城实验室 基于扩展观测的融合估计方法、装置、设备及存储介质
CN117521018B (zh) * 2024-01-08 2024-03-26 鹏城实验室 基于扩展观测的融合估计方法、装置、设备及存储介质

Similar Documents

Publication Publication Date Title
CN114217283A (zh) 多普勒雷达静态融合平滑变结构滤波方法及装置
CN110501696B (zh) 一种基于多普勒量测自适应处理的雷达目标跟踪方法
CN109444841B (zh) 基于修正切换函数的平滑变结构滤波方法及系统
CN113189561B (zh) 一种海杂波参数估计方法、系统、设备及存储介质
CN110865343A (zh) 基于lmb的粒子滤波检测前跟踪方法及系统
CN108871365B (zh) 一种航向约束下的状态估计方法及系统
He et al. A spline filter for multidimensional nonlinear state estimation
Govaers On independent axes estimation for extended target tracking
Frery et al. Analysis of minute features in speckled imagery with maximum likelihood estimation
CN111965618A (zh) 一种融合多普勒量测的转换量测跟踪方法及系统
Sliwka et al. Using interval methods in the context of robust localization of underwater robots
Hirsenkorn et al. Learning sensor models for virtual test and development
CN112328959B (zh) 一种基于自适应扩展卡尔曼概率假设密度滤波器的多目标跟踪方法
CN111736144B (zh) 一种仅用距离观测的机动转弯目标状态估计方法
Fu et al. Maneuvering target tracking with improved unbiased FIR filter
CN111708013B (zh) 一种距离坐标系目标跟踪滤波方法
Konatowski et al. A comparison of estimation accuracy by the use of KF, EKF & UKF filters
Ebert et al. Deep radar sensor models for accurate and robust object tracking
Kamthe et al. Iterative state estimation in non-linear dynamical systems using approximate expectation propagation
Mohammed et al. Reduced cubature Kalman filtering applied to target tracking
CN111190173B (zh) 一种基于预测值量测转换的相控阵雷达目标跟踪方法
Baerveldt et al. Extended target PMBM tracker with a Gaussian Process target model on LiDAR data
CN115047444A (zh) 基于广义时变平滑层策略的多普勒雷达滤波方法及系统
CN117784114B (zh) 异常噪声下基于混合熵的不规则扩展目标跟踪方法
CN113009468B (zh) 一种视线坐标系内的解耦cmkf跟踪方法及系统

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