CN109444841B - 基于修正切换函数的平滑变结构滤波方法及系统 - Google Patents

基于修正切换函数的平滑变结构滤波方法及系统 Download PDF

Info

Publication number
CN109444841B
CN109444841B CN201811600312.0A CN201811600312A CN109444841B CN 109444841 B CN109444841 B CN 109444841B CN 201811600312 A CN201811600312 A CN 201811600312A CN 109444841 B CN109444841 B CN 109444841B
Authority
CN
China
Prior art keywords
target
radar
state
calculating
measurement
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
CN201811600312.0A
Other languages
English (en)
Other versions
CN109444841A (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.)
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 CN201811600312.0A priority Critical patent/CN109444841B/zh
Publication of CN109444841A publication Critical patent/CN109444841A/zh
Application granted granted Critical
Publication of CN109444841B publication Critical patent/CN109444841B/zh
Active 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
    • G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
    • G01S7/41—Details 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

Landscapes

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

Abstract

本发明公开了一种基于修正切换函数的平滑变结构滤波方法及系统，其中，该方法包括：根据雷达目标的运动过程和量测过程分别建立目标运动模型和雷达量测模型；利用目标运动模型和雷达跟踪的历史数据计算出当前时刻目标状态的预测值；根据雷达目标的运动状态和噪声的先验信息给定出目标运动模型和雷达量测模型的不确定度，根据不确定度的先验估计计算新息增益项；根据目标状态的预测值和新息增益项计算出目标状态的后验更新，并序贯迭代上述步骤对目标状态持续跟踪，直至雷达判断航迹终止。该方法采用双曲正切函数作为切换函数计算新息增益项，抑制滤波过程中的抖振现象，提高滤波器对未建模的系统不确定度的鲁棒性，获得更好的目标状态估计精度。

Description

基于修正切换函数的平滑变结构滤波方法及系统
技术领域
本发明涉及雷达数据处理技术领域，特别涉及一种基于修正切换函数的平滑变结构滤波方法及系统。
背景技术
雷达目标状态估计是指，利用预先假设的目标运动模型和雷达扫描回波中的目标量测点迹，对目标运动的位置、速度、加速度等状态信息进行后验估计，从而抑制过程噪声和量测噪声，提高雷达对目标的跟踪性能。但在现实场景中，用预设的目标运动模型难以准确描述目标真实的状态转移过程，同时模型参数难以实时估计，会导致模型不确定问题。传统的雷达目标状态估计方法，例如卡尔曼滤波(KF)、扩展卡尔曼滤波(EKF)、不敏卡尔曼滤波(UKF)等，基于固定的假设模型和滤波器结构，在模型不确定条件下状态估计误差会急剧增大，甚至滤波发散导致丢失跟踪目标。
平滑变结构滤波(SVSF)是一种近年来兴起的状态估计方法，基于滑膜控制理论和变结构控制理论，能够有效应对模型不确定问题。该方法与传统滤波器的差异主要在于新息增益项的计算：通过引入切换函数，建立了分段的切换控制准则，根据当前时刻和前一时刻的量测误差灵活地控制增益项矢量的方向和大小；在切换函数中引入平滑边界层矢量参数，对模型不确定度的上界进行先验估计，根据该参数和当前时刻量测误差的相对大小确定采用何种切换控制准则，从而实现变结构的基于不确定模型的目标状态估计。符号函数和饱和函数是该方法常用的两种切换函数。
然而，用于描述模型不确定度的上界的平滑层参数是难以准确给定的。如果该参数与实际不确定度相比过大，会造成滤波器收敛慢，状态估计误差增大；如果过小，则会引起系统抖振问题。抖振现象是指，由于平滑层参数低估了真实的模型不确定度，导致切换控制准则的频繁调整，引起增益矢量的计算规则的时序不一致，使得后验状态估计在表征目标真实状态的滑膜面两侧反复“穿越”而不收敛。抖振现象与切换函数的不连续或者不光滑有直接关系，会引起额外的状态估计误差甚至滤波发散，从而导致雷达目标跟踪性能急剧恶化。
因此，如何设计更加鲁棒的切换函数，以应对被平滑层参数低估的模型不确定度导致的抖振现象，从而提高雷达目标状态估计的性能，是SVSF在工程应用中亟待解决的问题。
发明内容
本发明旨在至少在一定程度上解决相关技术中的技术问题之一。
为此，本发明的一个目的在于提出一种基于修正切换函数的平滑变结构滤波方法，该方法采用双曲正切函数作为切换函数计算新息增益项，抑制滤波过程中的抖振现象，提高滤波器对未建模的系统不确定度的鲁棒性，获得更好的目标状态估计精度。
本发明的另一个目的在于提出一种基于修正切换函数的平滑变结构滤波系统。
为达到上述目的，本发明一方面提出了基于修正切换函数的平滑变结构滤波方法，包括：根据雷达目标的运动过程和量测过程分别建立目标运动模型和雷达量测模型；利用所述目标运动模型和雷达跟踪的历史数据计算出当前时刻目标状态的预测值；根据雷达目标的运动状态和噪声的先验信息给定出所述目标运动模型和所述雷达量测模型的不确定度，并根据所述不确定度的先验估计计算新息增益项；根据所述目标状态的预测值和所述新息增益项计算出目标状态的后验更新，并序贯迭代上述步骤对目标状态持续跟踪，直至雷达判断航迹终止。
本发明实施例的基于修正切换函数的平滑变结构滤波方法，通过采用双曲正切函数作为切换函数以计算新息增益项，代替传统的符号函数或者饱和函数，基于“预测—更新”的序贯迭代步骤对雷达目标状态进行后验估计，能够有效抑制滤波过程中的抖振现象，提高滤波器对未建模的系统不确定度的鲁棒性，从而获得更好的目标状态估计精度，具有较高的工程应用价值。
另外，根据本发明上述实施例的基于修正切换函数的平滑变结构滤波方法还可以具有以下附加的技术特征：
进一步地，在本发明的一个实施例中，所述目标运动模型为：
Xk+1＝FkXk+Gkuk+vk
其中，
Figure BDA0001922315880000021
表示第k次雷达扫描回波的目标状态变量，描述目标的位置、速度、加速度等状态信息，
Figure BDA0001922315880000022
表示所述目标运动模型中状态转移矩阵，
Figure BDA0001922315880000023
和
Figure BDA0001922315880000024
分别表示输入项和对应的输入矩阵，
Figure BDA0001922315880000025
是均值为0；
所述雷达量测模型为：
Zk+1＝Hk+1Xk+1+wk+1
其中，
Figure BDA0001922315880000026
和
Figure BDA0001922315880000027
分别是雷达量测变量和量测矩阵，
Figure BDA0001922315880000028
是均值为0，协方差阵为
Figure BDA0001922315880000029
的高斯量测噪声。
进一步地，在本发明的一个实施例中，所述利用所述目标运动模型和雷达跟踪的历史数据计算出当前时刻目标状态的预测值，进一步包括：
根据k时刻的状态变量及其协方差阵的后验估计Xk|k和Pk|k，按照预设的所述目标运动模型和所述雷达量测模型，计算k+1时刻的目标状态估计
Figure BDA0001922315880000031
雷达量测估计
Figure BDA0001922315880000032
和状态估计协方差阵P的先验预测为：
Figure BDA0001922315880000033
Figure BDA0001922315880000034
Figure BDA0001922315880000035
其中，协方差阵为
Figure BDA0001922315880000036
的高斯过程噪声。
进一步地，在本发明的一个实施例中，所述根据雷达目标的运动状态和噪声的先验信息给定出所述目标运动模型和所述雷达量测模型的不确定度，并根据所述不确定度的先验估计计算新息增益项，进一步包括：
根据k+1时刻雷达扫描的目标量测Zk+1，计算先验所述新息增益项：
Figure BDA0001922315880000037
计算混合误差项，式中ez,k|k表示k时刻的后验量测误差，|·|ABS表示对矢量逐元素取绝对值：
EZ＝|ez,k+1|k|ABS+γ|ez,k|k|ABS
计算所述新息增益项，式中上角标+表示伪逆运算，用
Figure BDA00019223158800000310
表示舒尔积算符，γ是衰减因子且满足0<γ<1，diag表示将矢量映射为对角阵：
Figure BDA0001922315880000038
其中，ψ是预设的平滑层矢量参数，tanh(·)表示双曲正切函数，矢量函数的第i维定义为：
tanhi(x,ψ)＝tanh(xi/ψi)，i＝1,2,…,n。
进一步地，在本发明的一个实施例中，所述根据所述目标状态的预测值和所述新息增益项计算出目标状态的后验更新，进一步包括：
利用各项所述目标状态的预测值和所述新息增益项，计算目标状态变量
Figure BDA0001922315880000039
和状态估计协方差阵P的后验更新为：
Figure BDA0001922315880000041
Figure BDA0001922315880000042
其中，
Figure BDA0001922315880000043
是单位阵，Rk+1表示k+1时刻的量测噪声协方差阵。
为达到上述目的，本发明另一方面提出了一种基于修正切换函数的平滑变结构滤波系统，包括：建模模块，用于根据雷达目标的运动过程和量测过程分别建立目标运动模型和雷达量测模型；第一计算模块，用于利用所述目标运动模型和雷达跟踪的历史数据计算出当前时刻目标状态的预测值；第二计算模块，用于根据雷达目标的运动状态和噪声的先验信息给定出所述目标运动模型和所述雷达量测模型的不确定度，并根据所述不确定度的先验估计计算新息增益项；更新迭代模块，用于根据所述目标状态的预测值和所述新息增益项计算出目标状态的后验更新，并序贯迭代上述步骤对目标状态持续跟踪，直至雷达判断航迹终止。
本发明实施例的基于修正切换函数的平滑变结构滤波系统，通过采用双曲正切函数作为切换函数以计算新息增益项，代替传统的符号函数或者饱和函数，基于“预测—更新”的序贯迭代步骤对雷达目标状态进行后验估计，能够有效抑制滤波过程中的抖振现象，提高滤波器对未建模的系统不确定度的鲁棒性，从而获得更好的目标状态估计精度，具有较高的工程应用价值。
另外，根据本发明上述实施例的基于修正切换函数的平滑变结构滤波系统还可以具有以下附加的技术特征：
进一步地，在本发明的一个实施例中，所述目标运动模型为：
Xk+1＝FkXk+Gkuk+vk
其中，
Figure BDA0001922315880000044
表示第k次雷达扫描回波的目标状态变量，描述目标的位置、速度、加速度等状态信息，
Figure BDA0001922315880000045
表示所述目标运动模型中状态转移矩阵，
Figure BDA0001922315880000046
和
Figure BDA0001922315880000047
分别表示输入项和对应的输入矩阵，
Figure BDA0001922315880000048
是均值为0；
所述雷达量测模型为：
Zk+1＝Hk+1Xk+1+wk+1
其中，
Figure BDA0001922315880000049
和
Figure BDA00019223158800000410
分别是雷达量测变量和量测矩阵，
Figure BDA00019223158800000411
是均值为0，协方差阵为
Figure BDA00019223158800000412
的高斯量测噪声。
进一步地，在本发明的一个实施例中，所述第一计算模块进一步包括：
根据k时刻的状态变量及其协方差阵的后验估计Xk|k和Pk|k，按照预设的所述目标运动模型和所述雷达量测模型，计算k+1时刻的目标状态估计
Figure BDA0001922315880000051
雷达量测估计
Figure BDA0001922315880000052
和状态估计协方差阵P的先验预测为：
Figure BDA0001922315880000053
Figure BDA0001922315880000054
Figure BDA0001922315880000055
其中，协方差阵为
Figure BDA0001922315880000056
的高斯过程噪声。
进一步地，在本发明的一个实施例中，所述第二计算模块进一步包括：
根据k+1时刻雷达扫描的目标量测Zk+1，计算先验所述新息增益项：
Figure BDA0001922315880000057
计算混合误差项，式中ez,k|k表示k时刻的后验量测误差，|·|ABS表示对矢量逐元素取绝对值：
EZ＝|ez,k+1|k|ABS+γ|ez,k|k|ABS
计算所述新息增益项，式中上角标+表示伪逆运算，用
Figure BDA00019223158800000513
表示舒尔积算符，γ是衰减因子且满足0<γ<1，diag表示将矢量映射为对角阵：
Figure BDA0001922315880000058
其中，ψ是预设的平滑层矢量参数，tanh(·)表示双曲正切函数，矢量函数的第i维定义为：
tanhi(x,ψ)＝tanh(xi/ψi)，i＝1,2,…,n。
进一步地，在本发明的一个实施例中，所述更新迭代模块进一步包括：用各项所述目标状态的预测值和所述新息增益项，计算目标状态变量
Figure BDA0001922315880000059
和状态估计协方差阵P的后验更新为：
Figure BDA00019223158800000510
Figure BDA00019223158800000511
其中，
Figure BDA00019223158800000512
是单位阵，Rk+1表示k+1时刻的量测噪声协方差阵。
本发明附加的方面和优点将在下面的描述中部分给出，部分将从下面的描述中变得明显，或通过本发明的实践了解到。
附图说明
本发明上述的和/或附加的方面和优点从下面结合附图对实施例的描述中将变得明显和容易理解，其中：
图1为根据本发明实施例的基于修正切换函数的平滑变结构滤波方法流程图；
图2为根据本发明一个实施例的实施例的基于修正切换函数的平滑变结构滤波方法操作流程框图；
图3为2-D近空搜索雷达对单机动目标的滤波跟踪的仿真场景实例图；
图4为平滑变结构滤波的三种切换函数示意图；
图5为仿真场景下基于修正切换函数的平滑变结构滤波方法的目标状态估计误差随时间变化示意图；
图6为根据本发明实施例的基于修正切换函数的平滑变结构滤波系统结构示意图。
具体实施方式
下面详细描述本发明的实施例，所述实施例的示例在附图中示出，其中自始至终相同或类似的标号表示相同或类似的元件或具有相同或类似功能的元件。下面通过参考附图描述的实施例是示例性的，旨在用于解释本发明，而不能理解为对本发明的限制。
下面参照附图描述根据本发明实施例提出的基于修正切换函数的平滑变结构滤波方法及系统，首先将参照附图描述根据本发明实施例提出的基于修正切换函数的平滑变结构滤波方法。
图1为根据本发明实施例的基于修正切换函数的平滑变结构滤波方法流程图。
如图1所示，该基于修正切换函数的平滑变结构滤波方法包括以下步骤：
在步骤S101中，根据雷达目标的运动过程和量测过程分别建立目标运动模型和雷达量测模型。
进一步地，在本发明的一个实施例中，如图2所示，目标运动模型为：
Xk+1＝FkXk+Gkuk+vk
其中，
Figure BDA0001922315880000061
表示第k次雷达扫描回波的目标状态变量，描述目标的位置、速度、加速度等状态信息，
Figure BDA0001922315880000062
表示目标运动模型中状态转移矩阵，
Figure BDA0001922315880000063
和
Figure BDA0001922315880000064
分别表示输入项和对应的输入矩阵，
Figure BDA0001922315880000065
是均值为0，协方差阵为
Figure BDA0001922315880000066
的高斯过程噪声。
雷达量测模型为：
Zk+1＝Hk+1Xk+1+wk+1
其中，
Figure BDA0001922315880000071
和
Figure BDA0001922315880000072
分别是雷达量测变量和量测矩阵，
Figure BDA0001922315880000073
是均值为0，协方差阵为
Figure BDA0001922315880000074
的高斯量测噪声。
进一步地，考虑到模型不确定性问题，用
Figure BDA0001922315880000075
和
Figure BDA0001922315880000076
作为估计的模型矩阵用于滤波。
在步骤S102中，利用目标运动模型和雷达跟踪的历史数据计算出当前时刻目标状态的预测值。
根据目标运动状态和噪声的先验信息，给定平滑层矢量参数
Figure BDA0001922315880000077
用于表示模型不确定度的上界，以k＝0时刻目标状态变量X0和状态估计协方差P0为初始化数据进行迭代，假设已知k时刻的状态变量及其协方差阵的后验估计Xk|k和Pk|k，按照以下步骤操作：
根据k时刻的状态变量及其协方差阵的后验估计Xk|k和Pk|k，按照预设的目标运动模型和雷达量测模型，计算k+1时刻的目标状态估计
Figure BDA0001922315880000078
雷达量测估计
Figure BDA0001922315880000079
和状态估计协方差阵P的先验预测(用下标k+1|k表示)为：
Figure BDA00019223158800000710
Figure BDA00019223158800000711
Figure BDA00019223158800000712
其中，协方差阵为
Figure BDA00019223158800000713
的高斯过程噪声。
在步骤S103中，根据雷达目标的运动状态和噪声的先验信息给定出目标运动模型和雷达量测模型的不确定度，并根据不确定度的先验估计计算新息增益项。
进一步地，根据k+1时刻雷达扫描的目标量测Zk+1，计算先验新息增益项：
Figure BDA00019223158800000714
计算混合误差项，式中ez,k|k表示k时刻的后验量测误差，|·|ABS表示对矢量逐元素取绝对值：
EZ＝|ez,k+1|k|ABS+γ|ez,k|k|ABS
计算新息增益项，式中上角标+表示伪逆运算，用
Figure BDA00019223158800000715
表示舒尔积算符，γ是衰减因子且满足0<γ<1，diag表示将矢量映射为对角阵：
Figure BDA0001922315880000081
其中，ψ是预设的平滑层矢量参数，tanh(·)表示双曲正切函数，矢量函数的第i维定义为：
tanhi(x,ψ)＝tanh(xi/ψi)，i＝1,2,…,n。
在步骤S104中，根据目标状态的预测值和新息增益项计算出目标状态的后验更新，并序贯迭代上述步骤对目标状态持续跟踪，直至雷达判断航迹终止
进一步地，在本发明的一个实施例中，利用步骤S102中各项的先验预测值和S103中的新息增益项，计算目标状态变量
Figure BDA0001922315880000082
和状态估计协方差阵P的后验更新(用下标k+1|k+1表示)为：
Figure BDA0001922315880000083
Figure BDA0001922315880000084
其中，
Figure BDA0001922315880000085
是单位阵，Rk+1表示k+1时刻的量测噪声协方差阵。
最后，对第k＝1,2,…次雷达扫描回波，重复步骤S102、S103和S104的序贯迭代过程，直至雷达系统判定航迹结束。从而实现模型不确定条件下雷达目标状态的实时精确估计。
进一步地，在本发明实施的基于修正切换函数的平滑变结构滤波方法能够有效抑制抖振现象，与传统的基于符号函数或饱和函数的SVSF方法相比，本发明实施例采用的双曲正切函数是一种连续且光滑的切换函数，在保证滤波收敛性的前提下，提供了“软切换”的控制切换准则，避免了控制准则的频繁调换导致的增益矢量计算偏差，从而有效抑制了抖振现象的发生。
进一步地，在本发明实施的基于修正切换函数的平滑变结构滤波方法对目标状态估计的精度高，采用双曲正切函数作为切换函数，建立了“软切换”的控制切换准则，保证了新息增益项计算规则的时序一致性，使得后验状态估计更快地收敛到滑膜面附近，同时降低在平滑层内高频抖振的幅度，从而减小后验估计与真实状态之间的偏差，后验状态估计误差协方差更逼近后验克拉美罗下界(PCRLB)。因此，与传统的两种切换函数相比，本实施例的方法能够取得更高精度的目标状态估计结果。
进一步地，在本发明实施的基于修正切换函数的平滑变结构滤波方法的跟踪鲁棒性好，在模型不确定度难以估计的情况下能够保持更好的目标状态估计精度，从而实现鲁棒的目标跟踪。如果实际的模型不确定度小于预估的平滑层参数，本实施例方法相比于传统的基于饱和函数的SVSF具有相似的滤波输出精度，均好于符号函数；如果实际不确定度与估计参数相比相近或者更大，本实施例方法明显优于传统的两种切换函数。因此，整体而言，本发明实施例的方法能够保持更加鲁棒的状态估计精度。
进一步地，在本发明实施的基于修正切换函数的平滑变结构滤波方法算法复杂度低，相比于传统的卡尔曼滤波(KF)等目标状态估计方法，避免了矩阵求逆等高复杂度的计算过程；相比于传统的SVSF方法，切换函数的替换几乎没有引入额外的计算量。因此，本算法的复杂度较低，在现有硬件技术水平下完全可以实现实时的雷达目标跟踪，符合工程应用的要求。
具体地，下面根据具体实施例对本发明的基于修正切换函数的平滑变结构滤波方法进行详细说明。
如图3所示，机动目标跟踪是模型不确定问题的典型场景。考虑一部2-D近空搜索雷达对单机动目标的跟踪滤波的仿真场景。假定雷达监视区域为50km×50km的平面矩形，雷达建立在坐标原点处，目标机动轨迹如图3所示。本实施例拟比较本发明提出的基于双曲正切函数的SVSF方法(记为Tanh-SVSF)，和传统的基于符号函数和基于饱和函数的SVSF方法(分别记为Sign-SVSF和Sat-SVSF)，三种函数的对比示意图如图4所示，其中，(a)传统的符号函数；(b)传统的饱和函数；(c)本发明方法涉及的双曲正切函数。
其中，采用500次蒙特卡洛实验的平均结果作为对比数据。表1为雷达目标跟踪场景参数，利用表1的仿真参数生成仿真的雷达目标数据。目标状态变量定义为X-Y坐标轴内的位置、速度、加速度，即X＝[x y vx vy ax ay]T；雷达量测变量定义为目标位置，即Z＝[xy]T，采用常加速度模型(CA)对目标状态进行估计，即：
Figure BDA0001922315880000091
量测模型矩阵为：
Figure BDA0001922315880000092
下面给出具体的雷达参数、目标状态参数及滤波算法参数：
表1雷达目标跟踪场景参数
Figure BDA0001922315880000101
以k＝0时刻目标状态变量X0和状态估计协方差P0为初始化数据进行迭代，假设已知k时刻的状态变量及其协方差阵的后验估计Xk|k和Pk|k，按照以下步骤操作：
1)根据k时刻的状态变量及其协方差阵的后验估计Xk|k和Pk|k，按照预设的目标运动模型和雷达量测模型，计算k+1时刻的目标状态估计
Figure BDA0001922315880000102
雷达量测估计
Figure BDA0001922315880000103
和状态估计协方差阵P的先验预测：
Figure BDA0001922315880000104
Figure BDA0001922315880000105
Figure BDA0001922315880000106
2)根据k+1时刻雷达扫描的目标量测Zk+1，计算先验量测误差(即新息)：
Figure BDA0001922315880000107
计算混合误差项，这里用ez,k|k表示k时刻的后验量测误差，|·|ABS表示对矢量逐元素取绝对值：
EZ＝|ez,k+1|k|ABS+γ|ez,k|k|ABS
计算新息增益项，这里用上角标+表示伪逆运算，用
Figure BDA0001922315880000108
表示舒尔积算符，γ是衰减因子且满足0<γ<1，diag表示将矢量映射为对角阵：
Figure BDA0001922315880000109
其中ψ是预设的平滑层矢量参数，tanh(·)表示双曲正切函数，也就是平滑变结构滤波方法的切换函数，该矢量函数的第i维定义为：
tanhi(x,ψ)＝tanh(xi/ψi)
3)利用步骤1)中各项的先验预测值和步骤2)中的新息增益项，计算目标状态变量
Figure BDA0001922315880000114
和状态估计协方差阵P的后验更新：
Figure BDA0001922315880000111
Figure BDA0001922315880000112
其中，
Figure BDA0001922315880000113
是单位阵，Rk+1表示k+1时刻的量测噪声协方差阵。
对第k＝1,2,…,100次雷达扫描回波，重复上述的序贯迭代过程，直至目标航迹终止，从而实现模型不确定条件下雷达目标状态的实时精确估计。
进一步地，理论分析和实验结果均表明，本发明方法在模型不确定度难以准确估计的情况下，取得了很好的抖振抑制效果，提高了雷达目标状态估计的精度和鲁棒性，具有很好的应用价值。
如图5所示，其中，(a)X方向坐标估计误差；(b)X方向速度估计误差；(c)X方向加速度估计误差，给出了X方向目标状态估计误差随扫描时间的收敛情况(Y方向相似)，包括位置、速度和加速度；表2为模型不确定目标的状态估计方均根误差，表2给出了500次蒙特卡洛实验的方均根误差(RMSE)。由图5可知，本发明提出的Tanh-SVSF方法的估计误差随时间收敛更快，速度和加速度误差明显优于传统方法，在模型不确定度难以准确估计的条件下，能够实现更加鲁棒的雷达目标状态估计。
表2模型不确定目标的状态估计方均根误差
Figure BDA0001922315880000121
根据本发明实施例提出的基于修正切换函数的平滑变结构滤波方法，通过采用双曲正切函数作为切换函数以计算新息增益项，代替传统的符号函数或者饱和函数，基于“预测—更新”的序贯迭代步骤对雷达目标状态进行后验估计，能够有效抑制滤波过程中的抖振现象，提高滤波器对未建模的系统不确定度的鲁棒性，从而获得更好的目标状态估计精度，具有较高的工程应用价值。
其次参照附图描述根据本发明实施例提出的基于修正切换函数的平滑变结构滤波系统。
图6为本发明实施例的基于修正切换函数的平滑变结构滤波系统结构示意图。
如图6所示，该基于修正切换函数的平滑变结构滤波系统10包括：建模模块100、第一计算模块200、第二计算模块300和更新迭代模块400。
其中，建模模块100用于根据雷达目标的运动过程和量测过程分别建立目标运动模型和雷达量测模型。第一计算模块200用于利用目标运动模型和雷达跟踪的历史数据计算出当前时刻目标状态的预测值。第二计算模块300根据雷达目标的运动状态和噪声的先验信息给定出目标运动模型和雷达量测模型的不确定度，并根据不确定度的先验估计计算新息增益项。更新迭代模块400用于根据目标状态的预测值和新息增益项计算出目标状态的后验更新，并序贯迭代上述步骤对目标状态持续跟踪，直至雷达判断航迹终止。该平滑变结构滤波系统10采用双曲正切函数作为切换函数计算新息增益项，有效的抑制滤波过程中的抖振现象，提高滤波器对未建模的系统不确定度的鲁棒性，获得更好的目标状态估计精度。
进一步地，在本发明的一个实施例中，目标运动模型为：
Xk+1＝FkXk+Gkuk+vk
其中，
Figure BDA0001922315880000131
表示第k次雷达扫描回波的目标状态变量，描述目标的位置、速度、加速度等状态信息，
Figure BDA0001922315880000132
表示目标运动模型中状态转移矩阵，
Figure BDA0001922315880000133
和
Figure BDA0001922315880000134
分别表示输入项和对应的输入矩阵，
Figure BDA0001922315880000135
是均值为0；
雷达量测模型为：
Zk+1＝Hk+1Xk+1+wk+1
其中，
Figure BDA0001922315880000136
和
Figure BDA0001922315880000137
分别是雷达量测变量和量测矩阵，
Figure BDA0001922315880000138
是均值为0，协方差阵为
Figure BDA0001922315880000139
的高斯量测噪声。
进一步地，在本发明的一个实施例中，第一计算模块进一步包括：根据k时刻的状态变量及其协方差阵的后验估计Xk|k和Pk|k，按照预设的目标运动模型和雷达量测模型，计算k+1时刻的目标状态估计
Figure BDA00019223158800001310
雷达量测估计
Figure BDA00019223158800001311
和状态估计协方差阵P的先验预测为：
Figure BDA00019223158800001312
Figure BDA00019223158800001313
Figure BDA00019223158800001314
其中，协方差阵为
Figure BDA00019223158800001315
的高斯过程噪声。
进一步地，在本发明的一个实施例中，第二计算模块进一步包括：根据k+1时刻雷达扫描的目标量测Zk+1，计算先验新息增益项：
Figure BDA00019223158800001316
计算混合误差项，式中ez,k|k表示k时刻的后验量测误差，|·|ABS表示对矢量逐元素取绝对值：
EZ＝|ez,k+1|k|ABS+γ|ez,k|k|ABS
计算新息增益项，式中上角标+表示伪逆运算，用
Figure BDA00019223158800001318
表示舒尔积算符，γ是衰减因子且满足0<γ<1，diag表示将矢量映射为对角阵：
Figure BDA00019223158800001317
其中，ψ是预设的平滑层矢量参数，tanh(·)表示双曲正切函数，矢量函数的第i维定义为：
tanhi(x,ψ)＝tanh(xi/ψi)，i＝1,2,…,n。
进一步地，在本发明的一个实施例中，更新迭代模块进一步包括：用各项目标状态的预测值和新息增益项，计算目标状态变量
Figure BDA0001922315880000141
和状态估计协方差阵P的后验更新为：
Figure BDA0001922315880000143
Figure BDA0001922315880000144
其中，
Figure BDA0001922315880000145
是单位阵，Rk+1表示k+1时刻的量测噪声协方差阵。
需要说明的是，前述对基于修正切换函数的平滑变结构滤波方法实施例的解释说明也适用于该基于修正切换函数的平滑变结构滤波系统，此处不再赘述。
根据本发明实施例提出的基于修正切换函数的平滑变结构滤波系统，通过采用双曲正切函数作为切换函数以计算新息增益项，代替传统的符号函数或者饱和函数，基于“预测—更新”的序贯迭代步骤对雷达目标状态进行后验估计，能够有效抑制滤波过程中的抖振现象，提高滤波器对未建模的系统不确定度的鲁棒性，从而获得更好的目标状态估计精度，具有较高的工程应用价值。
此外，术语“第一”、“第二”仅用于描述目的，而不能理解为指示或暗示相对重要性或者隐含指明所指示的技术特征的数量。由此，限定有“第一”、“第二”的特征可以明示或者隐含地包括至少一个该特征。在本发明的描述中，“多个”的含义是至少两个，例如两个，三个等，除非另有明确具体的限定。
在本发明中，除非另有明确的规定和限定，术语“安装”、“相连”、“连接”、“固定”等术语应做广义理解，例如，可以是固定连接，也可以是可拆卸连接，或成一体；可以是机械连接，也可以是电连接；可以是直接相连，也可以通过中间媒介间接相连，可以是两个元件内部的连通或两个元件的相互作用关系，除非另有明确的限定。对于本领域的普通技术人员而言，可以根据具体情况理解上述术语在本发明中的具体含义。
在本发明中，除非另有明确的规定和限定，第一特征在第二特征“上”或“下”可以是第一和第二特征直接接触，或第一和第二特征通过中间媒介间接接触。而且，第一特征在第二特征“之上”、“上方”和“上面”可是第一特征在第二特征正上方或斜上方，或仅仅表示第一特征水平高度高于第二特征。第一特征在第二特征“之下”、“下方”和“下面”可以是第一特征在第二特征正下方或斜下方，或仅仅表示第一特征水平高度小于第二特征。
在本说明书的描述中，参考术语“一个实施例”、“一些实施例”、“示例”、“具体示例”、或“一些示例”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施例或示例中。在本说明书中，对上述术语的示意性表述不必须针对的是相同的实施例或示例。而且，描述的具体特征、结构、材料或者特点可以在任一个或多个实施例或示例中以合适的方式结合。此外，在不相互矛盾的情况下，本领域的技术人员可以将本说明书中描述的不同实施例或示例以及不同实施例或示例的特征进行结合和组合。
尽管上面已经示出和描述了本发明的实施例，可以理解的是，上述实施例是示例性的，不能理解为对本发明的限制，本领域的普通技术人员在本发明的范围内可以对上述实施例进行变化、修改、替换和变型。

Claims (8)

1.一种基于修正切换函数的平滑变结构滤波方法，其特征在于，包括以下步骤：
根据雷达目标的运动过程和量测过程分别建立目标运动模型和雷达量测模型；
利用所述目标运动模型和雷达跟踪的历史数据计算出当前时刻目标状态的预测值；
根据雷达目标的运动状态和噪声的先验信息给定出所述目标运动模型和所述雷达量测模型的不确定度，并根据所述不确定度的先验估计计算新息增益项，具体地，根据k+1时刻雷达扫描的目标量测Zk+1，计算先验所述新息增益项：
Figure FDA0002536391210000011
计算混合误差项，式中ez,k|k表示k时刻的后验量测误差，|·|ABS表示对矢量逐元素取绝对值：
EZ＝|ez,k+1|k|ABS+γ|ez,k|k|ABS
计算所述新息增益项，式中上角标+表示伪逆运算，用
Figure FDA0002536391210000012
表示舒尔积算符，γ是衰减因子且满足0<γ<1，diag表示将矢量映射为对角阵：
Figure FDA0002536391210000013
其中，ψ是预设的平滑层矢量参数，tanh(·)表示双曲正切函数，矢量函数的第i维定义为：
tanhi(x,ψ)＝tanh(xi/ψi)，i＝1,2,…,n；以及
根据所述目标状态的预测值和所述新息增益项计算出目标状态的后验更新，并序贯迭代上述步骤对目标状态持续跟踪，直至雷达判断航迹终止。
2.根据权利要求1所述的基于修正切换函数的平滑变结构滤波方法，其特征在于，所述目标运动模型为：
Xk+1＝FkXk+Gkuk+vk
其中，
Figure FDA0002536391210000014
表示第k次雷达扫描回波的目标状态变量，描述目标的位置、速度、加速度等状态信息，
Figure FDA0002536391210000015
表示所述目标运动模型中状态转移矩阵，
Figure FDA0002536391210000016
和
Figure FDA0002536391210000017
分别表示输入项和对应的输入矩阵，
Figure FDA0002536391210000018
是均值为0，协方差阵为
Figure FDA0002536391210000019
的高斯过程噪声；
所述雷达量测模型为：
Zk+1＝Hk+1Xk+1+wk+1
其中，
Figure FDA0002536391210000021
和
Figure FDA0002536391210000022
分别是雷达量测变量和量测矩阵，
Figure FDA0002536391210000023
是均值为0，协方差阵为
Figure FDA0002536391210000024
的高斯量测噪声。
3.根据权利要求1所述的基于修正切换函数的平滑变结构滤波方法，其特征在于，所述利用所述目标运动模型和雷达跟踪的历史数据计算出当前时刻目标状态的预测值，进一步包括：
根据k时刻的状态变量及其协方差阵的后验估计Xk|k和Pk|k，按照预设的所述目标运动模型和所述雷达量测模型，计算k+1时刻的目标状态估计
Figure FDA0002536391210000025
雷达量测估计
Figure FDA0002536391210000026
和状态估计协方差阵P的先验预测为：
Figure FDA0002536391210000027
Figure FDA0002536391210000028
Figure FDA0002536391210000029
其中，协方差阵为
Figure FDA00025363912100000210
的高斯过程噪声。
4.根据权利要求1所述的基于修正切换函数的平滑变结构滤波方法，其特征在于，所述根据所述目标状态的预测值和所述新息增益项计算出目标状态的后验更新，进一步包括：
利用各项所述目标状态的预测值和所述新息增益项，计算目标状态变量
Figure FDA00025363912100000211
和状态估计协方差阵P的后验更新为：
Figure FDA00025363912100000212
Figure FDA00025363912100000213
其中，
Figure FDA00025363912100000214
是单位阵，Rk+1表示k+1时刻的量测噪声协方差阵。
5.一种基于修正切换函数的平滑变结构滤波系统，其特征在于，包括：
建模模块，用于根据雷达目标的运动过程和量测过程分别建立目标运动模型和雷达量测模型；
第一计算模块，用于利用所述目标运动模型和雷达跟踪的历史数据计算出当前时刻目标状态的预测值；
第二计算模块，用于根据雷达目标的运动状态和噪声的先验信息给定出所述目标运动模型和所述雷达量测模型的不确定度，并根据所述不确定度的先验估计计算新息增益项，具体地，根据k+1时刻雷达扫描的目标量测Zk+1，计算先验所述新息增益项：
Figure FDA0002536391210000031
计算混合误差项，式中ez,k|k表示k时刻的后验量测误差，|·|ABS表示对矢量逐元素取绝对值：
EZ＝|ez,k+1|k|ABS+γ|ez,k|k|ABS
计算所述新息增益项，式中上角标+表示伪逆运算，用
Figure FDA0002536391210000032
表示舒尔积算符，γ是衰减因子且满足0<γ<1，diag表示将矢量映射为对角阵：
Figure FDA0002536391210000033
其中，ψ是预设的平滑层矢量参数，tanh(·)表示双曲正切函数，矢量函数的第i维定义为：
tanhi(x,ψ)＝tanh(xi/ψi)，i＝1,2,…,n；以及
更新迭代模块，用于根据所述目标状态的预测值和所述新息增益项计算出目标状态的后验更新，并序贯迭代上述步骤对目标状态持续跟踪，直至雷达判断航迹终止。
6.根据权利要求5所述的基于修正切换函数的平滑变结构滤波系统，其特征在于，所述目标运动模型为：
Xk+1＝FkXk+Gkuk+vk
其中，
Figure FDA0002536391210000034
表示第k次雷达扫描回波的目标状态变量，描述目标的位置、速度、加速度等状态信息，
Figure FDA0002536391210000035
表示所述目标运动模型中状态转移矩阵，
Figure FDA0002536391210000036
和
Figure FDA0002536391210000037
分别表示输入项和对应的输入矩阵，
Figure FDA0002536391210000038
是均值为0；
所述雷达量测模型为：
Zk+1＝Hk+1Xk+1+wk+1
其中，
Figure FDA0002536391210000039
和
Figure FDA00025363912100000310
分别是雷达量测变量和量测矩阵，
Figure FDA00025363912100000311
是均值为0，协方差阵为
Figure FDA00025363912100000312
的高斯量测噪声。
7.根据权利要求5所述的基于修正切换函数的平滑变结构滤波系统，其特征在于，所述第一计算模块进一步包括：
根据k时刻的状态变量及其协方差阵的后验估计Xk|k和Pk|k，按照预设的所述目标运动模型和所述雷达量测模型，计算k+1时刻的目标状态估计
Figure FDA0002536391210000041
雷达量测估计
Figure FDA0002536391210000042
和状态估计协方差阵P的先验预测为：
Figure FDA0002536391210000043
Figure FDA0002536391210000044
Figure FDA0002536391210000045
其中，协方差阵为
Figure FDA0002536391210000046
的高斯过程噪声。
8.根据权利要求5所述的基于修正切换函数的平滑变结构滤波系统，其特征在于，所述更新迭代模块进一步包括：用各项所述目标状态的预测值和所述新息增益项，计算目标状态变量
Figure FDA0002536391210000047
和状态估计协方差阵P的后验更新为：
Figure FDA0002536391210000048
Figure FDA0002536391210000049
其中，
Figure FDA00025363912100000410
是单位阵，Rk+1表示k+1时刻的量测噪声协方差阵。
CN201811600312.0A 2018-12-26 2018-12-26 基于修正切换函数的平滑变结构滤波方法及系统 Active CN109444841B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201811600312.0A CN109444841B (zh) 2018-12-26 2018-12-26 基于修正切换函数的平滑变结构滤波方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201811600312.0A CN109444841B (zh) 2018-12-26 2018-12-26 基于修正切换函数的平滑变结构滤波方法及系统

Publications (2)

Publication Number Publication Date
CN109444841A CN109444841A (zh) 2019-03-08
CN109444841B true CN109444841B (zh) 2020-08-04

Family

ID=65535789

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201811600312.0A Active CN109444841B (zh) 2018-12-26 2018-12-26 基于修正切换函数的平滑变结构滤波方法及系统

Country Status (1)

Country Link
CN (1) CN109444841B (zh)

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109974698B (zh) * 2019-04-10 2020-11-17 清华大学深圳研究生院 一种室内物流小车定位方法和终端设备
CN112230195A (zh) * 2020-09-02 2021-01-15 清华大学 基于非线性最优平滑层策略的平滑变结构滤波方法及系统
CN115451952B (zh) * 2022-08-29 2023-11-07 南京航空航天大学 一种故障检测及抗差自适应滤波的多系统组合导航方法和装置
CN117406590B (zh) * 2023-10-08 2024-04-02 哈尔滨工业大学 一种基于平滑变结构滤波的鲁棒机动目标跟踪方法及系统

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101622669A (zh) * 2007-02-26 2010-01-06 高通股份有限公司 用于信号分离的系统、方法及设备
CN104007423A (zh) * 2014-05-27 2014-08-27 电子科技大学 基于混沌序列预测的天波雷达海杂波抑制方法
CN108319666A (zh) * 2018-01-19 2018-07-24 国网浙江省电力有限公司电力科学研究院 一种基于多模态舆情分析的供电服务评估方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101622669A (zh) * 2007-02-26 2010-01-06 高通股份有限公司 用于信号分离的系统、方法及设备
CN104007423A (zh) * 2014-05-27 2014-08-27 电子科技大学 基于混沌序列预测的天波雷达海杂波抑制方法
CN108319666A (zh) * 2018-01-19 2018-07-24 国网浙江省电力有限公司电力科学研究院 一种基于多模态舆情分析的供电服务评估方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
A SMOOTH VARIABLE STRUCTURE FILTER FOR STATE ESTIMATION;S.Wang等;《Control and Intelligent Systems》;20071231;第2页左列第1段-第4页右列倒数第2段 *
平滑变结构滤波算法研究及其在惯性导航系统初始对准中的应用;陈帅;《中国博士学位论文全文数据库 信息科技辑》;20180615;第1页第1段、第11页第5段、第14页最后1段-第26页第3段 *
陈帅.平滑变结构滤波算法研究及其在惯性导航系统初始对准中的应用.《中国博士学位论文全文数据库 信息科技辑》.2018, *

Also Published As

Publication number Publication date
CN109444841A (zh) 2019-03-08

Similar Documents

Publication Publication Date Title
CN109444841B (zh) 基于修正切换函数的平滑变结构滤波方法及系统
CN111985093B (zh) 一种带噪声估计器的自适应无迹卡尔曼滤波状态估计方法
CN109100714B (zh) 一种基于极坐标系的低慢小目标跟踪方法
CN109901153B (zh) 基于信息熵权和最近邻域数据关联的目标航迹优化方法
Janabi-Sharifi et al. A kalman-filter-based method for pose estimation in visual servoing
CN109633590B (zh) 基于gp-vsmm-jpda的扩展目标跟踪方法
CN111178385A (zh) 一种鲁棒在线多传感器融合的目标跟踪方法
CN109782240B (zh) 一种基于递推修正的多传感器系统误差配准方法和系统
CN112198503A (zh) 一种目标航迹预测优化方法、装置及雷达系统
CN110376582B (zh) 自适应gm-phd的机动目标跟踪方法
CN111291471B (zh) 一种基于l1正则无迹变换的约束多模型滤波方法
CN108710125A (zh) 用于目标跟踪的距离方位滤波方法
CN108871365B (zh) 一种航向约束下的状态估计方法及系统
CN110261859B (zh) 一种水下机动静止交替状态目标跟踪方法
CN114217283A (zh) 多普勒雷达静态融合平滑变结构滤波方法及装置
CN109613477B (zh) 一种复杂环境下的tdoa定位追踪方法
CN111679269B (zh) 一种基于变分的多雷达融合航迹状态估计方法
CN111262556B (zh) 一种同时估计未知高斯测量噪声统计量的多目标跟踪方法
CN107391446A (zh) 基于随机矩阵的不规则形状多扩展目标形状和状态估计方法
CN111274529A (zh) 一种鲁棒的高斯逆威沙特phd多扩展目标跟踪算法
CN110007298B (zh) 一种目标超前预测跟踪方法
CN116224320A (zh) 一种极坐标系下处理多普勒量测的雷达目标跟踪方法
CN112230195A (zh) 基于非线性最优平滑层策略的平滑变结构滤波方法及系统
CN115633306A (zh) 一种关于多区域uwb信号的定位修正方法及装置
CN110208791B (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