CN104793201A - 一种跟踪临近空间高超声速目标的修正变结构网格交互多模型滤波方法 - Google Patents

Abstract

一种跟踪临近空间高超声速目标的修正变结构网格交互多模型滤波方法，本发明涉及修正变结构网格交互多模型滤波方法。本发明的目的是为了解决现有单模型滤波算法、固定结构交互式多模型算法以及传统变结构交互多模型算法无法实现高精度、快速跟踪临近空间高超声速机动目标的问题。通过以下技术方案实现：步骤一、建立惯性参考坐标系，并在惯性参考坐标系中建立目标机动运动的状态方程；步骤二、中心模型采用机动目标当前统计模型，左转弯模型和右转弯模型采用匀速转弯模型；步骤三、基于惯性参考坐标系确定目标跟踪系统测量模型；步骤四、进行状态估计和误差协方差矩阵融合。本发明应用于航行器领域。

Description

一种跟踪临近空间高超声速目标的修正变结构网格交互多模型滤波方法

技术领域

本发明涉及修正变结构网格交互多模型滤波方法。

背景技术

临近空间高超声速滑翔弹头是一种能够在临近空间无动力滑翔数千至一万余公里，并具有较强的机动能力和突防性能，它能够携带核弹头或常规弹头实施远距离快速打击，具有较高的升阻比，可在大气层内长时间飞行，其运动轨迹往往呈现出“跳跃”特征。临近空间高超声速巡航导弹在临近空间巡航飞行时，其机动模式为侧向摆式机动，高度和速度基本保持不变。

从国内外目前机动目标跟踪问题的研究情况来看，多采用基于单模型，例如Singer模型和Jerk模型的跟踪滤波算法，以及固定结构交互式多模型算法或传统变结构交互多模型算法。单模型滤波算法难以覆盖目标可能采取的机动形式，跟踪误差大，固定结构交互式多模型算法存在模型需求量大、计算量大、模型之间转换过程计算效率低等问题，传统变结构交互多模型算法存在模型后验概率精确度较低且模型之间的切换需要一定时延的问题，导致单模型滤波算法、固定结构交互式多模型算法以及传统变结构交互多模型算法无法实现高精度、快速跟踪临近空间高超声速机动目标，因此基于单模型滤波算法、固定结构交互式多模型算法以及传统变结构交互多模型算法无法实现高精度、快速跟踪临近空间高超声速机动目标，需要提出一种新的跟踪滤波方法。

发明内容

本发明的目的是为了解决现有单模型滤波算法、固定结构交互式多模型算法以及传统变结构交互多模型算法无法实现高精度、快速跟踪临近空间高超声速机动目标的问题，而提出了一种跟踪临近空间高超声速目标的修正变结构网格交互多模型滤波方法。

上述的发明目的是通过以下技术方案实现的：

步骤一、建立惯性参考坐标系，并在惯性参考坐标系中建立目标机动运动的状态方程，即匀速转弯模型和机动目标当前统计模型，构成目标运动模型库；

步骤二、中心模型采用步骤一中的机动目标当前统计模型，左转弯模型和右转弯模型采用步骤一中的匀速转弯模型；

步骤三、基于步骤一中的惯性参考坐标系确定目标跟踪系统测量模型；

步骤四、将步骤二中的左转弯模型、中心模型和右转弯模型分别作为状态方程与步骤三中的目标跟踪系统测量模型相结合，采用离散型推广卡尔曼滤波进行计算，并根据计算结果对左转弯模型和右转弯模型的网格中心和网格距离重新调整，然后进行状态估计和误差协方差矩阵融合。

发明效果

采用本发明的一种跟踪临近空间高超声速目标的修正变结构网格交互多模型滤波方法，基于修正变结构网格算法，提供一种跟踪临近空间高超声速机动目标的交互多模型滤波方法，该方法采用机动目标当前统计模型，左右两边采用自适应转弯模型，从而得到基于修正变结构网格算法的临近空间机动目标交互多模型跟踪滤波方法。该方法采用推广卡尔曼滤波算法得到各状态变量估计，随后进行网格中心和网格距离重新调整，根据没有跳跃、向左跳变、向右跳变的转弯速率公式，得出状态估计和误差协方差矩阵融合计算结果。相比于传统变结构多模型方法，该方法对临近空间机动目标的跟踪精度更高。解决以往单模型滤波算法存在的难以覆盖目标可能机动形式，跟踪误差大等问题，固定结构交互式多模型算法存在的模型需求量大、计算量大、模型之间转换过程计算效率低的问题，以及传统变结构交互多模型算法存在模型后验概率精确度较低且模型之间的切换需要一定时延，导致单模型滤波算法、固定结构交互式多模型算法以及传统变结构交互多模型算法无法实现高精度、快速跟踪临近空间高超声速机动目标的问题。

在X-51A等高飞行时横向机动试验中，从表1可以看出，三个方向MVSIMM算法加速度估计误差均比其余两种滤波算法的误差要小，说明在该种机动模式下使用MVSIMM算法能够得到较高的跟踪滤波精度，可以实现较高精度跟踪临近空间目标飞行器。MVSIMM算法在X轴方向的加速度均方根误差比IMM-UKF算法少0.41，MVSIMM算法在X轴方向的加速度均方根误差比IMM-PF算法少1.21，MVSIMM算法在Y轴方向的加速度均方根误差比IMM-UKF算法少2.82，MVSIMM算法在Y轴方向的加速度均方根误差比IMM-PF算法少4.81，MVSIMM算法在Z轴方向的加速度均方根误差比IMM-UKF算法少1.84，，MVSIMM算法在Z轴方向的加速度均方根误差比IMM-PF算法少4.07，MVSIMM算法在Z轴方向的加速度均方根误差率比IMM-PF算法减少了29％。

附图说明

图1为本发明流程图；

图2为三维双站测量系统示意图，X为惯性参考坐标系x轴坐标，Y为惯性参考坐标系y轴坐标，Z为惯性参考坐标系z轴坐标，α为高低角，β为方位角，r为观测站到目标的距离，T为采样周期，V_T为目标飞行速度，O为惯性参考坐标系中心；

图3为X-51A纵向等高飞行轨迹图；

图4为X-51A横侧向机动飞行轨迹图，横坐标为惯性参考坐标系x轴坐标，纵坐标为惯性参考坐标系z轴坐标。

图5为X-51A等高飞行且横向机动时x方向加速度估计图；

图6为X-51A等高飞行且横向机动时y方向加速度估计图；

图7为X-51A等高飞行且横向机动时z方向加速度估计图；

图8为HTV-2常值攻角和常值倾斜角飞行轨迹图；

图9为HTV-2常值攻角和常值倾斜角飞行时x方向加速度估计图；

图10为HTV-2常值攻角和常值倾斜角飞行时y方向加速度估计图；

图11为HTV-2常值攻角和常值倾斜角飞行时z方向加速度估计图；

图12为HTV-2最大升阻比和常值倾斜角飞行轨迹图；

图13为HTV-2最大升阻比和常值倾斜角飞行时x方向加速度估计图；

图14为HTV-2最大升阻比和常值倾斜角飞行时y方向加速度估计图；

图15为HTV-2最大升阻比和常值倾斜角飞行时z方向加速度估计图。

具体实施方式

具体实施方式一：结合图1说明本实施方式，一种跟踪临近空间高超声速目标的修正变结构网格交互多模型滤波方法，具体是按照以下步骤进行的：

步骤一、建立惯性参考坐标系，并在惯性参考坐标系中建立目标机动运动的状态方程，即匀速转弯模型和机动目标当前统计模型，构成目标运动模型库；

步骤二、中心模型采用步骤一中的机动目标当前统计模型，左转弯模型和右转弯模型采用步骤一中的匀速转弯模型，通过中心模型、左转弯模型和右转弯模型构成修正变结构滤波器；

步骤三、基于步骤一中的惯性参考坐标系确定目标跟踪系统测量模型；

步骤四、将步骤二中的左转弯模型、中心模型和右转弯模型分别作为状态方程与步骤三中的目标跟踪系统测量模型相结合，采用离散型推广卡尔曼滤波进行计算，并根据计算结果对左转弯模型和右转弯模型的网格中心和网格距离重新调整，然后进行状态估计和误差协方差矩阵融合。

具体实施方式二：本实施方式与具体实施方式一不同的是：所述步骤一中建立惯性参考坐标系，并在惯性参考坐标系中建立目标机动运动的状态方程，即匀速转弯模型和机动目标当前统计模型，构成目标运动模型库；具体过程为：

建立惯性参考坐标系中临近空间目标运动状态向量取为：

$X_k={\left[\begin{array}{ccccccccc}{x}_k& {\stackrel{\·}{x}}_k& {\stackrel{\·\·}{x}}_k& y_k& {\stackrel{\·}{y}}_k& {\stackrel{\·\·}{y}}_k& z_k& {\stackrel{\·}{z}}_k& {\stackrel{\·\·}{z}}_k\end{array}\right]}^T---\left(1\right)$

式中，

X_k为k时刻的临近空间目标运动状态向量,k取正整数；

x_k为k时刻的目标位置在惯性参考坐标系x轴中的分量；

y_k为k时刻的目标位置在惯性参考坐标系y轴中的分量；

z_k为k时刻的目标位置在惯性参考坐标系z轴中的分量；

为k时刻的目标速度在惯性参考坐标系x轴中的分量；

为k时刻的目标速度在惯性参考坐标系y轴中的分量；

为k时刻的目标速度在惯性参考坐标系z轴中的分量；

为k时刻的目标加速度在惯性参考坐标系x轴中的分量；

为k时刻的目标加速度在惯性参考坐标系y轴中的分量；

为k时刻的目标加速度在惯性参考坐标系z轴中的分量；

匀速转弯模型(CT)的离散化动态方程为：

X_k＝F_CTX_k-1+Gw_k-1   (2)

式中，

$F_{\mathrm{CT}}=\left[\begin{array}{ccc}{A}_{\mathrm{CT}}& 0_{3\×3}& 0_{3\×3}\\ 0_{3\×3}& A_{\mathrm{CT}}& 0_{3\×3}\\ 0_{3\×3}& 0_{3\×3}& A_{\mathrm{CT}}\end{array}\right],A_{\mathrm{CT}}=\left[\begin{array}{ccc}1& \frac{\sin \mathrm{\ω\Δt}}{\mathrm{\ω}}& \frac{1-\cos \mathrm{\ω\Δt}}{{\mathrm{\ω}}^2}\\ 0& \cos \mathrm{\ω\Δt}& \frac{\sin \mathrm{\ω\Δt}}{\mathrm{\ω}}\\ 0& -\mathrm{\ω}\sin \mathrm{\ω\Δt}& \cos \mathrm{\ω\Δt}\end{array}\right],G=\left[\begin{array}{ccc}{\mathrm{\Δt}}^2/2& 0& 0\\ \mathrm{\Δt}& 0& 0\\ 1& 0& 0\\ 0& {\mathrm{\Δt}}^2/2& 0\\ 0& \mathrm{\Δt}& 0\\ 0& 1& 0\\ 0& 0& {\mathrm{\Δt}}^2/2\\ 0& 0& \mathrm{\Δt}\\ 0& 0& 1\end{array}\right]$

X_k为k时刻的临近空间目标运动状态向量，k取正整数；

ω为转弯速率，当向左转弯时取ω>0，向右转弯时取ω<0，ω＝0时模型以常值速度运动；

Δt为采样周期，一般取0.1s到0.5s；

w_k-1为k-1时刻的状态噪声向量；

X_k-1为k-1时刻的临近空间目标运动状态向量，在前一步滤波算法中得到的；

G为状态噪声向量转移矩阵；

A_CT为单通道匀速转弯模型状态转移矩阵；

F_CT为匀速转弯模型状态向量转移矩阵；

0_3×3为3×3零矩阵；

其中，所述，k-1时刻的状态噪声向量为根据目标实际运动状态来设定，

w_k-1＝diag(300m,20m/s,10m/s²,300m,20m/s,10m/s²,300m,20m/s,10m/s²)；

机动目标当前统计模型(CS)的离散化动态方程为：

$X_k=F_{\mathrm{CS}}{X}_{k-1}+B{\stackrel{\&OverBar;}{a}}_k+w_{k-1}---\left(3\right)$

式中，

$F_{\mathrm{CS}}=\left[\begin{array}{ccc}{A}_{\mathrm{CS}}& 0_{3\&times;3}& 0_{3\&times;3}\\ 0_{3\&times;3}& A_{\mathrm{CS}}& 0_{3\&times;3}\\ 0_{3\&times;3}& 0_{3\&times;3}& A_{\mathrm{CS}}\end{array}\right],A_{\mathrm{CS}}=\left[\begin{array}{ccc}1& \mathrm{\&Delta;t}& (\mathrm{\&lambda;\&Delta;t}-1+e^{-\mathrm{\&lambda;\&Delta;t}})/{\mathrm{\&lambda;}}^2\\ 0& 1& (1-e^{-\mathrm{\&lambda;\&Delta;t}})/\mathrm{\&lambda;}\\ 0& 0& e^{-\mathrm{\&lambda;\&Delta;t}}\end{array}\right],$

$B=(\left[\begin{array}{c}{\mathrm{\&Delta;t}}^2/2\\ \mathrm{\&Delta;t}\\ 1\end{array}\right]-\left[\begin{array}{c}(\mathrm{\&lambda;\&Delta;t}-1+e^{-\mathrm{\&lambda;\&Delta;t}})/{\mathrm{\&lambda;}}^2\\ (1-e^{-\mathrm{\&lambda;\&Delta;t}})/\mathrm{\&lambda;}\\ e^{-\mathrm{\&lambda;\&Delta;t}}\end{array}\right])$

X_k为k时刻的临近空间目标运动状态向量，k取正整数；

λ为机动时间常数的倒数，即机动频率，取值范围为：转弯机动λ＝1/60，逃避机动λ＝1/20，大气扰动λ＝1；

e为指数函数；

A_CS为单通道当前统计模型状态转移矩阵；

Δt为采样周期，一般取0.1s到0.5s；

F_CS为当前统计模型状态转移矩阵；

B为机动加速度当前均值转移矩阵；

w_k-1为k-1时刻的状态噪声向量；

X_k-1为k-1时刻的临近空间目标运动状态向量，在前一步滤波算法中得到的；

为机动加速度“当前”均值，在每一个采样周期内设定为常数。

其它步骤及参数与具体实施方式一相同。

具体实施方式三：本实施方式与具体实施方式一或二不同的是：所述步骤二中的中心模型采用机动目标当前统计模型，左转弯模型和右转弯模型采用匀速转弯模型，通过中心模型、左转弯模型和右转弯模型构成修正变结构滤波器；具体过程为：

(1)中心模型采用机动目标当前统计模型(CS)，

解决目标直线运动中匀速运动和匀加速运动的跟踪，

一种方法是同时增加CV和CA模型，但这种增加模型个数的方法会导致模型集中和模型个数增多，造成模型之间竞争而使起决策作用的各个模型后验概率的值分散，从而造成模型切换滞后；CV为匀速运动；CA为匀加速运动；

另一种方法是采用机动目标当前统计模型：机动目标当前统计模型同时具有跟踪匀速运动和匀加速运动的能力，因此可以代替CV和CA模型来跟踪这两种运动，同时并不增加模型个数；

机动目标当前统计模型(CS)的跟踪精度与最大加速度参数±a_max设置有关，为了在匀速运动和弱机动情况下有较高的跟踪精度，加速度误差小于2m/s²，最大加速度参数±a_max设置为±40m/s²；在滤波过程中，中心模型保持不变；

(2)左转弯模型和右转弯模型采用匀速转弯模型(CT)

目标在机动结束时需要延时，即延时大于0.05s才能切回到由机动目标当前统计模型和左右匀速转弯模型构成的模型集中的正确模型上，在跟踪过程中，中心模型始终保持为机动目标当前统计模型；而目标机动较大的运动则通过左右模型与中心模型的间隔来控制；左右模型为左转弯模型和右转弯模型；

左转弯模型和右转弯模型采用匀速转弯模型(CT)，分别设为左转弯模型(CTL)和右转弯模型(CTR)，左转弯模型转弯速率为右转弯模型转弯速率为其中，且 的取值为-0.098rad/s，，的取值为0rad/s；且 的取值为0rad/s，的取值为0.098rad/s；

为了缩小模型间隔，提高跟踪精度，借鉴有向图转换多模型算法(SGIMM算法)，以右转弯模型为例，将区间离散为N段而滤波器中的转弯速率则在此区间进行自适应切换。为了保证模型间的独立性和良好的模型转换速度，离散化间隔(即模型间隔)应不小于0.025rad/s，且离散化个数不易过多。

其它步骤及参数与具体实施方式一或二相同。

具体实施方式四：本实施方式与具体实施方式一、二或三不同的是：所述步骤三中基于步骤一中的惯性参考坐标系确定目标跟踪系统测量模型；具体过程为：

设测量系统的站址就在惯性参考坐标系的原点，定义方位角为β，高低角为α，方位角会受到量测方位角噪声v_β的污染，高低角会受到量测高低角噪声v_α的污染，v_β和v_α均是零均值高斯白噪声，标准差分别为σ_β和σ_α；

根据观察站以及目标S位置，得出相对位置关系如图2所示，由图2可见，方位角和高低角对于相对位置矢量的非线性函数为：

$\left\{\begin{array}{c}\mathrm{\&alpha;}={\tan }^{-1}\frac{x}{\sqrt{x^2+z^2}}+v_{\mathrm{\&alpha;}}\\ \mathrm{\&beta;}=-{\tan }^{-1}\frac{z}{x}+v_{\mathrm{\&beta;}}\end{array}\right.---\left(4\right)$

式中，

x为惯性参考坐标系x轴坐标；

z为惯性参考坐标系z轴坐标；

v_α为量测高低角噪声；

v_β为量测方位角噪声；

r为观察站到目标S位置的距离，距离测量方程为：

$r=\sqrt{x^2+y^2+z^2}+v_r---\left(5\right)$

式中，y为惯性参考坐标系y轴坐标；；

v_r为距离测量噪声，是零均值高斯白噪声，标准差为σ_r；

则目标跟踪系统测量模型方程为：

$Z_k=\left[\begin{array}{c}{\mathrm{\&alpha;}}_k\\ {\mathrm{\&beta;}}_k\\ r_k\end{array}\right]=\left[\begin{array}{c}{\tan }^{-1}\frac{x_k}{{\sqrt{x_k^2+z_k^2}}^{}}\\ -{\tan }^{-1}\frac{z_k}{x_k}\\ \sqrt{x_k^2+y_k^2+z_k^2}\end{array}\right]+\left[\begin{array}{c}{v}_{\mathrm{\&alpha;},k}\\ v_{\mathrm{\&beta;},k}\\ v_{r,k}\end{array}\right]=h\left(X_k\right)+\left[\begin{array}{c}{v}_{\mathrm{\&alpha;},k}\\ v_{\mathrm{\&beta;},k}\\ v_{r,k}\end{array}\right]---\left(6\right)$

式中，Z_k为k时刻的测量向量；

α_k为α在k时刻的值；

β_k为β在k时刻的值；

r_k为r在k时刻的值；

v_α,k为v_α在k时刻的值；

v_β,k为v_β在k时刻的值；

v_r,k为v_r在k时刻的值；

h(X_k)为测量函数在k时刻的值，

$h\left(X_k\right)=\left[\begin{array}{c}{\tan }^{-1}\frac{x_k}{\sqrt{x_k^2+z_k^2}}\\ -{\tan }^{-1}\frac{z_k}{x_k}\\ \sqrt{x_k^2+y_k^2+z_k^2}\end{array}\right]---\left(7\right).$

其它步骤及参数与具体实施方式一、二或三相同。

具体实施方式五：本实施方式与具体实施方式一、二、三或四不同的是：所述步骤四中将步骤二中的左转弯模型、中心模型和右转弯模型分别作为状态方程与步骤三中的目标跟踪系统测量模型相结合，采用离散型推广卡尔曼滤波进行计算，并根据计算结果对左转弯模型和右转弯模型的网格中心和网格距离重新调整，然后进行状态估计和误差协方差矩阵融合；具体过程为：

以作为k时刻的系统模型集，经过输入交互及预测估计，求得各模型的似然函数j＝L,C,R，有

${\mathrm{\&Lambda;}}_k^j=N({\tilde{z}}_k^j,0,S_k^j)=\left|2\mathrm{\&pi;}{S}_k^j{|}^{-\frac{1}{2}}\exp \right[-\frac{1}{2}{\tilde{z}}_k^{\mathrm{jT}}{\left(S_k^j\right)}^{-1}{\tilde{z}}_k^j]---\left(8\right)$

其中，j＝L,C,R，为左转弯模型、中心模型和右转弯模型的新息；

为左转弯模型、中心模型和右转弯模型的状态一步预报估计；

为左转弯模型、中心模型和右转弯模型的测量函数一步预报值；

为左转弯模型、中心模型和右转弯模型的新息协方差矩阵；

Z_k为k时刻的测量向量；

上标T为转置运算；

L为左转弯模型；

C为中心模型；

R为右转弯模型；

代表以为变量，均值为0，方差为的高斯分布函数；

k时刻的左转弯模型、中心模型和右转弯模型的交互式多模型后验概率μ_CTL,k、μ_CS,k和μ_CTR,k的计算公式为

${\mathrm{\&mu;}}_{\mathrm{CTL},k}=\frac{1}{C_{\mathrm{CTL}}}{\mathrm{\&Lambda;}}_k^j{C}_{\mathrm{CTL},j},C_{\mathrm{CTL},j}=\underset{i=L}{\overset{R}{\mathrm{\&Sigma;}}}{p}_{\mathrm{ij}}{\mathrm{\&mu;}}_{\mathrm{CTL},k-1}---\left(9\right)$

${\mathrm{\&mu;}}_{\mathrm{CS},k}=\frac{1}{C_{\mathrm{CS}}}{\mathrm{\&Lambda;}}_k^j{C}_{\mathrm{CS},j},C_{\mathrm{CS},j}=\underset{i=L}{\overset{R}{\mathrm{\&Sigma;}}}{p}_{\mathrm{ij}}{\mathrm{\&mu;}}_{\mathrm{CS},k-1}---\left(10\right)$

${\mathrm{\&mu;}}_{\mathrm{CTR},k}=\frac{1}{C_{\mathrm{CTR}}}{\mathrm{\&Lambda;}}_k^j{C}_{\mathrm{CTR},j},C_{\mathrm{CTR},j}=\underset{i=L}{\overset{R}{\mathrm{\&Sigma;}}}{p}_{\mathrm{ij}}{\mathrm{\&mu;}}_{\mathrm{CTR},k-1}---\left(11\right)$

式中，μ_CTL,k为k时刻的左转弯模型的交互式多模型后验概率，k取正整数；

μ_CS,k为k时刻的中心模型的交互式多模型后验概率；

μ_CTR,k为k时刻的右转弯模型的交互式多模型后验概率；

p_ij为从模型i向模型j转换的概率，i＝L,C,R，j＝L,C,R；i为模型i，j为模型j；

μ_CTL,k-1为k-1时刻的左转弯模型的交互式多模型后验概率；

μ_CS,k-1为k-1时刻的中心模型的交互式多模型后验概率；

μ_CTR,k-1为k-1时刻的右转弯模型的交互式多模型后验概率；

c_CTL为左转弯模型的归一化常数；

c_CS为中心模型的归一化常数；

c_CTR为右转弯模型的归一化常数；

L为左转弯模型；

R为右转弯模型；

c_CTL、c_CS和c_CTR分别定义为：

$C_{\mathrm{CTL}}=\underset{j=L}{\overset{R}{\mathrm{\&Sigma;}}}{\mathrm{\&Lambda;}}_k^j{C}_{\mathrm{CTL},j}---\left(12\right)$

$C_{\mathrm{CS}}=\underset{j=L}{\overset{R}{\mathrm{\&Sigma;}}}{\mathrm{\&Lambda;}}_k^j{C}_{\mathrm{CS},j}---\left(13\right)$

$C_{\mathrm{CTR}}=\underset{j=L}{\overset{R}{\mathrm{\&Sigma;}}}{\mathrm{\&Lambda;}}_k^j{C}_{\mathrm{CTR},j}---\left(14\right)$

L为左转弯模型；

R为右转弯模型；

采用离散型推广卡尔曼滤波，离散型推广卡尔曼滤波算法中模型转换方法如下：

(a)若右转弯模型的后验概率μ_CTR>t₂，其中的t₂＝0.9为探测重要模式的阈值,则目标发生了向右的机动，将当前右转弯模型转弯速率向右扩大0.025rad/s，即其中，和分别代表k时刻和k-1时刻的右转弯模型的转弯速率；

若左转弯模型的后验概率μ_CTL>t₂，则目标发生了向左的机动，将当前左转弯模型转弯速率向左扩大0.025rad/s,即其中，和分别代表k时刻和k-1时刻的左转弯模型的转弯速率；

(b)当左转弯模型和右转弯模型的后验概率都小于探测模式的阈值t₁时，t₁＝0.1，即μ_CTL<t₁,μ_CTR<t₁，则目标转弯速率下降，将当前左转弯模型转弯速率和右转弯模型转弯速率一共往中间靠拢一个单位，即 ${\mathrm{\&omega;}}_k^L={\mathrm{\&omega;}}_{k-1}^L+0.0125\mathrm{rad}/s,{\mathrm{\&omega;}}_k^R={\mathrm{\&omega;}}_{k-1}^R-0.0125\mathrm{rad}/s;$

定义有向图D_k是由点与点之间的连线组成的一种数据结构，即

D_k＝(V(D_k),E(D_k))

式中，V(D_k)是非空的顶点集合，E(D_k)是不与V(D_k)相交的边集合，并具有方向性，在有向图D_k中，若从顶点U到顶点V之间有路径，则称U和V是连通的；

设目标机动转弯速率ω在连续范围[-ω_max,ω_max]内变化，针对时变支撑的有向图D_k的交互式多模型算法，网格值k＝0,1,…N,N为正整数，为机动目标当前统计中心模型的转弯速率，-ω_max的取值为-0.098rad/s，ω_max的取值为0.098rad/s，

模型集合即包含左转弯模型、机动目标当前统计中心模型、右转弯模型转弯速率集合；

状态集合 $S_{X_k}={\{X}_{L,k},X_{C,k},X_{R,k}\},$

式中，为k时刻左转弯模型、中心模型、右转弯模型构成的的运动状态向量集合；

X_L,k为左转弯模型估计出来的k时刻状态；

X_C,k为中心模型估计出来的k时刻状态；

X_R,k为右转弯模型估计出来的k时刻状态；

设定误差协方差矩阵集合为 $S_{P_k}=\{P_{L,k},P_{C,k},P_{R,k}\},$

式中，为k时刻左转弯模型、中心模型、右转弯模型的状态估计误差协方差矩阵构成的集合；

P_L,k为左转弯模型k时刻状态误差协方差矩阵；

P_R,k为右转弯模型k时刻状态误差协方差矩阵；

P_C,k为中心模型k时刻状态误差协方差矩阵；

模式转换概率矩阵P_LCR为：

$P_{\mathrm{LCR}}=\left[\begin{array}{ccc}{p}_{\mathrm{LL}}& p_{\mathrm{LC}}& p_{\text{LR}}\\ p_{\mathrm{CL}}& p_{\mathrm{CC}}& p_{\mathrm{CR}}\\ p_{\mathrm{RL}}& p_{\mathrm{RC}}& p_{\mathrm{RR}}\end{array}\right]---\left(15\right)$

式中，p_LL为左转弯模型继续保持的概率；

p_LC为由左转弯模型转化为中心模型的概率；

p_LR为由左转弯模型转换为右转弯模型的概率；

p_CL为由中心模型转化为左转弯模型的概率；

p_CC为中心模型继续保持的概率；

p_CR为由中心模型转化为右转弯模型的概率；

p_RL为由右转弯模型转化为左转弯模型的概率；

p_RC为由右转弯模型转化为中心模型的概率；

p_RR为右转弯模型继续保持的概率；

修正变结构交互多模型滤波算法以粗网格：

$D_0=\{{\mathrm{\&omega;}}_0^L=-{\mathrm{\&omega;}}_{\max },{\mathrm{\&omega;}}_0^C=0,{\mathrm{\&omega;}}_0^R={\mathrm{\&omega;}}_{\max }\}---\left(16\right)$

式中，D₀为初始时刻修正变结构交互多模型滤波算法左转弯模型、机动目标当前统计中心模型、右转弯模型的转弯速率集合；

为初始时刻修正变结构交互多模型滤波算法左转弯模型的转弯速率；

ω_max为初始时刻模型转弯速率最大值；

为初始时刻修正变结构交互多模型滤波算法机动目标当前统计中心模型的转弯速率；

为初始时刻修正变结构交互多模型滤波算法右转弯模型的转弯速率；

在每一个循环时间步(k→k+1)网格按照如下四步进行调整：

(1)左转弯模型、中心模型，以及右转弯模型，均采用离散型推广卡尔曼滤波：

$\left\{\begin{array}{c}{\stackrel{\&OverBar;}{X}}_{k|k-1}=F{\hat{X}}_{k-1}+B{\stackrel{\&OverBar;}{a}}_{k-1}\\ P_{k|k-1}=F{R}_{k-1}{F}^T+Q_{k-1}\\ {\hat{X}}_k={\stackrel{\&OverBar;}{X}}_{k|k-1}+K_k[Z_k-h({\hat{X}}_{k|k-1},k)]\\ P_k=(I-K_k{H}_k)P_{k|k-1}{(I-K_k{H}_k)}^T+K_k{R}_k{K}_k^T\\ K_k=P_{k|k-1}{H}_k^T[H_k{P}_{k|k-1}{H}_k^T+R_k{]}^{-1}\end{array}\right.---\left(17\right)$

式中，为从k-1时刻第k时刻的状态一步预报值；

F为子模型状态转移矩阵；

B为机动加速度均值转移矩阵；

为k-1时刻机动加速度均值，在每一个采样周期内设定为常数；

F^T为子模型状态转移矩阵的转置矩阵；

h为系统量测函数；

上标_T为转置矩阵符号；

H_k为根据测量方程计算出来的Jacobi矩阵；

为第k-1步的状态估计值；

P_k|k-1为的状态预报误差协方差矩阵；

P_k-1为的状态估计误差协方差矩阵；

Q_k-1为状态噪声协方差矩阵；

为第k步的状态估计值；

K_k为卡尔曼滤波增益矩阵；

Z_k为第k步的测量向量；

为根据测量方程计算出来的测量向量的一步预报值；

P_k为的状态估计误差协方差矩阵；

I为单位矩阵；

R_k为测量噪声协方差矩阵；

(2)对左转弯模型和右转弯模型的网格中心和网格距离重新调整：

网格中心调整按照下式计算：

${\mathrm{\&omega;}}_k^C={\mathrm{\&mu;}}_{\mathrm{CTL},k}{\mathrm{\&omega;}}_{k-1}^L+{\mathrm{\&mu;}}_{\mathrm{CS},k}{\mathrm{\&omega;}}_{k-1}^C+{\mathrm{\&mu;}}_{\mathrm{CTR},k}{\mathrm{\&omega;}}_{k-1}^R---\left(18\right)$

式中，和分别为k时刻和k-1时刻的中心模型转弯速率；

μ_CTL,k为k时刻左转弯模型的后验概率；

μ_CS,k为k时刻中心模型的后验概率；

μ_CTR,k为k时刻右转弯模型的后验概率；

为k-1时刻左转弯模型转弯速率；

为k-1时刻右转弯模型转弯速率；

左转弯模型和右转弯模型网格距离重新调整分为三种情况：

(a)没有跳跃

k时刻中心模型的后验概率最大，即μ_CS,k＝max{μ_CTL,k,μ_CS,k,μ_CTR,k}，

式中， ${\mathrm{\&lambda;}}_{L,k}=\max \{{\mathrm{\&omega;}}_k^C-{\mathrm{\&omega;}}_{k-1}^L,{\mathrm{\&delta;}}_{\mathrm{\&omega;}}\},$

${\mathrm{\&lambda;}}_{R,k}=\max \{{\mathrm{\&omega;}}_{k-1}^R-{\mathrm{\&omega;}}_k^C,{\mathrm{\&delta;}}_{\mathrm{\&omega;}}\};$

t₁＝0.1为探测模式的阈值；

δ_ω为模型分离距离0.025rad/s；

为k时刻左转弯模型的转弯速率；

为k时刻右转弯模型的转弯速率；

为k时刻中心模型的转弯速率；

λ_L,k为k时刻左转弯模型分离距离最大值；

λ_R,k为k时刻右转弯模型分离距离最大值；

(b)向左跳变

k时刻左转弯模型的后验概率最大，即μ_CTL,k＝max{μ_CTL,k,μ_CS,k,μ_CTR,k}，

${\mathrm{\&omega;}}_k^R={\mathrm{\&omega;}}_k^C+{\mathrm{\&lambda;}}_{R,k}---\left(22\right)$

式中，t₂＝0.9为探测重要模式的阈值；

(c)向右跳变

k时刻右转弯模型的后验概率最大，即μ_CTR,k＝max{μ_CTL,k,μ_CS,k,μ_CTR,k}，

${\mathrm{\&omega;}}_k^L={\mathrm{\&omega;}}_k^C+{\mathrm{\&lambda;}}_{L,k}---\left(23\right)$

(3)状态估计和误差协方差矩阵融合

状态向量X_k取作左转弯模型、中心模型，以及右转弯模型的运动状态的融合；

$X_k=X_{L,k}{\mathrm{\&omega;}}_k^L+X_{C,k}{\mathrm{\&omega;}}_k^C+X_{R,k}{\mathrm{\&omega;}}_k^R---\left(25\right)$

状态向量误差协方差矩阵P_k取作左转弯模型、中心模型，以及右转弯模型的状态的误差协方差矩阵P_Lk、P_Ck和P_Rk的融合；

$P_k=P_{L,k}{\mathrm{\&omega;}}_k^L+P_{C,k}{\mathrm{\&omega;}}_k^C+P_{R,k}{\mathrm{\&omega;}}_k^R---(26)$

式中，为第k时刻左转弯模型的转弯速率；

为第k时刻中心模型的转弯速率；

为第k时刻右转弯模型的转弯速率；

每一次循环对的子模型转移矩阵F需重新计算。

其它步骤及参数与具体实施方式一、二、三或四相同。

采用以下实施例验证本发明的有益效果：

实施例1

下面，通过数值仿真验证跟踪临近空间高超声速目标的修正变结构网格交互多模型滤波方法的效果。

假定地面雷达站采样周期为T＝0.1s，距离量测标准差为σ_r＝300m，角度测量标准差为σ_β＝0.001rad和σ_α＝0.001rad。

设仿真初始时刻t₀＝82s，仿真结束时刻为t_f＝400s。为了进行对比分析，我们还利用Singer模型分别设计了交互式多模型无迹卡尔曼滤波算法和交互式多模型粒子滤波算法，用来与本专利提出的修正变结构交互式多模型算法进行对比仿真分析。

目标初始位置在发射点惯性坐标系中的分量为x₀＝1176km,y₀＝-76km,z₀＝3578m，目标初始速度在发射点惯性坐标系中的分量为

目标跟踪滤波器的初始估计值取为

$X_0={\left[\begin{array}{ccccccccc}{x}_0+300m& {\stackrel{\&CenterDot;}{x}}_0+20m/s& 0& y_0+300m& {\stackrel{\&CenterDot;}{y}}_0+20m/s& 0& z_0+300m& {\stackrel{\&CenterDot;}{z}}_0+20m/s& 0\end{array}\right]}^T$

基于Singer模型的IMM-UKF算法和IMM-PF算法的参数设置如下：

Singer模型1：机动频率α₁＝1，表示大气扰动，最大加速度a_1,max＝0.1m/s²，最大加速度发生概率P_1,max＝0.1，非机动概率P_1,0＝0.9；

Singer模型2：机动频率α₂＝1/60，表示慢速转弯，最大加速度a_2,max＝20m/s²，最大加速度发生概率P_2,max＝0.9，非机动概率P_2,0＝0.1；

Singer模型3：机动频率α₃＝1/20，表示逃避机动，最大加速度a_3,max＝60m/s²，最大加速度发生概率P_3,max＝0.9，非机动概率P_3,0＝0.1；

IMM-UKF算法和IMM-PF算法：采用三种Singer模型，模型先验概率u_i(0)＝1/3(i＝1,2,3),Markov模型转移概率矩阵：

$\mathrm{\&Pi;}=\left[\begin{array}{ccc}0.90& 0.05& 0.05\\ 0.20& 0.75& 0.05\\ 0.20& 0.05& 0.75\end{array}\right]$

初始状态误差协方差矩阵取为

$P_0=\left[\begin{array}{ccc}{P}_{x0}& 0_{3\&times;3}& 0_{3\&times;3}\\ 0_{3\&times;3}& P_{y0}& 0_{3\&times;3}\\ 0_{3\&times;3}& 0_{3\&times;3}& P_{z0}\end{array}\right],P_{x0}=P_{y0}=P_{z0}\left[\begin{array}{ccc}{300}^2& 0& 0\\ 0& {20}^2& 0\\ 0& 0& {10}^2\end{array}\right]$

下面重点讨论X-51A和HTV-2加速度的估计结果，使用三种滤波算法分别对两类飞行器可能存在的四种机动模式进行跟踪滤波，得到三轴方向加速度估计：

a)X-51A等高飞行时横向机动

该机动模式下的目标等高飞行、横侧向机动运动轨迹如图3和图4所示，三个方向加速度估计如图5-7所示。

计算目标跟踪加速度均方根误差公式如下：

$a_{\mathrm{RMSE}}=\left\{\frac{1}{M}\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}\right[{({\stackrel{\&CenterDot;\&CenterDot;}{x}}_k-{\widehat{\stackrel{\&CenterDot;\&CenterDot;}{x}}}_k)}^2+{({\stackrel{\&CenterDot;\&CenterDot;}{y}}_k-{\widehat{\stackrel{\&CenterDot;\&CenterDot;}{y}}}_k)}^2+{({\stackrel{\&CenterDot;\&CenterDot;}{z}}_k-{\widehat{\stackrel{\&CenterDot;\&CenterDot;}{z}}}_k)}^2]{\}}^{\frac{1}{2}}---\left(27\right)$

根据上式计算得到加速度均方根误差如表1所示。

表1 三轴方向加速度均方根误差比较

从表1可以看出，三个方向MVSIMM算法加速度估计误差均比其余两种滤波算法的误差要小，说明在该种机动模式下使用MVSIMM算法能够得到较高的跟踪滤波精度，可以实现较高精度跟踪临近空间目标飞行器。

b)HTV-2常值攻角和常值倾斜角飞行

该机动模式下的目标飞行轨迹如图8所示。三个方向加速度估计如图9-11所示。

根据式(11)，得到加速度均方根误差如表2所示。从表2可以看出，三个方向MVSIMM算法加速度估计误差均比其余两种滤波算法的误差要小，说明在该种机动模式下MVSIMM算法跟踪滤波精度更高，跟踪效果更好。

表2 三轴方向加速度均方根误差比较

c)HTV-2最大升阻比和常值倾斜角飞行

该机动模式下的目标飞行轨迹如图12所示。三个方向加速度估计如图13-15所示。根据式(11)，计算得到加速度均方根误差如表4所示。从表4可以看出，三个方向MVSIMM算法加速度均方根误差均比其余两种滤波算法的误差要小，说明在该种机动模式下使用MVSIMM算法能够得到较高的跟踪滤波精度，具有一定的参考作用。

表4 三轴方向加速度均方根误差比较

综上所述，从四种机动模式的仿真结果来看，由于“等高飞行时进行横向机动”机动模式要求弹道倾角为0，目标的机动性相对较高，而交互式多模型无迹卡尔曼滤波算法主要是在非线性高斯情况下利用最小均方误差准则获得目标的动态估计，但是在用于强非线性、非高斯系统时会产生极大误差，因此该算法滤波得到的加速度均方根误差最大，跟踪滤波效果最差。交互式多模型粒子滤波算法计算量庞大，而且迭代后会产生退化现象。MVSIMM算法相比IMM-UKF算法和IMM-PF算法来说，该算法的适应环境最为广泛，更适合于复杂的非高斯环境。所以经过MVSIMM算法跟踪滤波后的加速度估计精度最好，跟踪滤波误差最小；当以“常值攻角和常值倾斜角”机动模式飞行时，由于该模式运动较为平缓，所以使用三种滤波算法均方根误差相差不多，没有出现加速度反向的情况，但是相比来说，MVSIMM算法的均方根误差最小；当以“等动压飞行，常值倾斜角”的机动模式飞行时，由于该模式相比其他模式机动性最强，会产生强非线性问题，因此使用交互式多模型无迹卡尔曼滤波算法时加速度均方根误差最大，相比而言，MVSIMM算法的加速度均方根误差最小，说明该算法在强机动性情况下仍然可以保持相对较好的估计精度；当以“最大升阻比飞行，常值倾斜角”的机动模式飞行时，该机动模式仍然保持较强的机动性，从y向加速度均方根误差较大可以看出，MVSIMM算法的加速度估计精度最好，说明该算法具有一定的工程实用价值。

使用本专利提出的修正变结构交互式多模型算法对X-51A和HTV-2等高超声速飞行器跟踪滤波时，跟踪误差较小，精度较高，误差在允许的范围之内，可见修正变结构交互式多模型算法相比固定结构交互式多模型算法在跟踪临近空间目标方面更具优势，具有工程实用价值。

Claims (5)

1.一种跟踪临近空间高超声速目标的修正变结构网格交互多模型滤波方法，其特征在于，一种跟踪临近空间高超声速目标的修正变结构网格交互多模型滤波方法，具体是按照以下步骤进行的：

步骤一、建立惯性参考坐标系，并在惯性参考坐标系中建立目标机动运动的状态方程，即匀速转弯模型和机动目标当前统计模型，构成目标运动模型库；

步骤二、中心模型采用步骤一中的机动目标当前统计模型，左转弯模型和右转弯模型采用步骤一中的匀速转弯模型；

步骤三、基于步骤一中的惯性参考坐标系确定目标跟踪系统测量模型；

步骤四、将步骤二中的左转弯模型、中心模型和右转弯模型分别作为状态方程与步骤三中的目标跟踪系统测量模型相结合，采用离散型推广卡尔曼滤波进行计算，并根据计算结果对左转弯模型和右转弯模型的网格中心和网格距离重新调整，然后进行状态估计和误差协方差矩阵融合。

2.根据权利要求1所述一种跟踪临近空间高超声速目标的修正变结构网格交互多模型滤波方法，其特征在于，所述步骤一中建立惯性参考坐标系，并在惯性参考坐标系中建立目标机动运动的状态方程，即匀速转弯模型和机动目标当前统计模型，构成目标运动模型库；具体过程为：

建立惯性参考坐标系中临近空间目标运动状态向量取为：

$X_k={\left[\begin{array}{ccccccccc}{x}_k& {\stackrel{\&CenterDot;}{x}}_k& {\stackrel{\&CenterDot;\&CenterDot;}{x}}_k& y_k& {\stackrel{\&CenterDot;}{y}}_k& {\stackrel{\&CenterDot;\&CenterDot;}{y}}_k& z_k& {\stackrel{\&CenterDot;}{z}}_k& {\stackrel{\&CenterDot;\&CenterDot;}{z}}_k\end{array}\right]}^T---\left(1\right)$

式中，

X_k为k时刻的临近空间目标运动状态向量,k为正整数；

x_k为k时刻的目标位置在惯性参考坐标系x轴中的分量；

y_k为k时刻的目标位置在惯性参考坐标系y轴中的分量；

z_k为k时刻的目标位置在惯性参考坐标系z轴中的分量；

为k时刻的目标速度在惯性参考坐标系x轴中的分量；

为k时刻的目标速度在惯性参考坐标系y轴中的分量；

为k时刻的目标速度在惯性参考坐标系z轴中的分量；

为k时刻的目标加速度在惯性参考坐标系x轴中的分量；

为k时刻的目标加速度在惯性参考坐标系y轴中的分量；

为k时刻的目标加速度在惯性参考坐标系z轴中的分量；

匀速转弯模型的离散化动态方程为：

X_k＝F_CTX_k-1+Gw_k-1            (2)

式中，

$F_{\mathrm{GT}}=\left[\begin{array}{ccc}{A}_{\mathrm{CT}}& 0_{3\&times;3}& 0_{3\&times;3}\\ 0_{3\&times;3}& A_{\mathrm{CT}}& 0_{3\&times;3}\\ 0_{3\&times;3}& 0_{3\&times;3}& A_{\mathrm{CT}}\end{array}\right],A_{\mathrm{CT}}=\left[\begin{array}{ccc}1& \frac{\sin \mathrm{\&omega;\&Delta;t}}{\mathrm{\&omega;}}& \frac{1-\cos \mathrm{\&omega;\&Delta;t}}{{\mathrm{\&omega;}}^2}\\ 0& \cos \mathrm{\&omega;\&Delta;t}& \frac{\sin \mathrm{\&omega;\&Delta;t}}{\mathrm{\&omega;}}\\ 0& -\mathrm{\&omega;}\sin \mathrm{\&omega;\&Delta;t}& \cos \mathrm{\&omega;\&Delta;t}\end{array}\right],G=\left[\begin{array}{ccc}\mathrm{\&Delta;}{t}^2/2& 0& 0\\ \mathrm{\&Delta;t}& 0& 0\\ 1& 0& 0\\ 0& \mathrm{\&Delta;}{t}^2/2& 0\\ 0& \mathrm{\&Delta;t}& 0\\ 0& 1& 0\\ 0& 0& \mathrm{\&Delta;}{t}^2/2\\ 0& 0& \mathrm{\&Delta;t}\\ 0& 0& 1\end{array}\right]$

X_k为k时刻的临近空间目标运动状态向量,k为正整数；

ω为转弯速率，当向左转弯时取ω>0，向右转弯时取ω<0，ω＝0时模型以常值速度运动；

CT为匀速转弯模型；

Δt为采样周期，取值范围为0.1s～0.5s；

w_k-1为k-1时刻的状态噪声向量；

X_k-1为k-1时刻的临近空间目标运动状态向量；

G为状态噪声向量转移矩阵；

A_CT为单通道匀速转弯模型状态转移矩阵；

F_CT为匀速转弯模型状态向量转移矩阵；

0_3×3为3×3零矩阵；

机动目标当前统计模型的离散化动态方程为：

$X_k=F_{\mathrm{CS}}{X}_{k-1}+B{\stackrel{\&OverBar;}{a}}_k+w_{k-1}---\left(3\right)$

式中，

$F_{\mathrm{CS}}=\left[\begin{array}{ccc}{A}_{\mathrm{CS}}& 0_{3\&times;3}& 0_{3\&times;3}\\ 0_{3\&times;3}& A_{\mathrm{CS}}& 0_{3\&times;3}\\ 0_{3\&times;3}& 0_{3\&times;3}& A_{\mathrm{CS}}\end{array}\right],A_{\mathrm{CS}}=\left[\begin{array}{ccc}1& \mathrm{\&Delta;t}& (\mathrm{\&lambda;\&Delta;t}-1+e^{-\mathrm{\&lambda;\&Delta;t}})/{\mathrm{\&lambda;}}^2\\ 0& 1& (1-e^{-\mathrm{\&lambda;\&Delta;t}})/\mathrm{\&lambda;}\\ 0& 0& e^{-\mathrm{\&lambda;\&Delta;t}}\end{array}\right],$

$B=(\left[\begin{array}{c}\mathrm{\&Delta;}{t}^2/2\\ \mathrm{\&Delta;t}\\ 1\end{array}\right]-\left[\begin{array}{c}(\mathrm{\&lambda;\&Delta;t}-1+e^{-\mathrm{\&lambda;\&Delta;t}})/{\mathrm{\&lambda;}}^2\\ (1-e^{-\mathrm{\&lambda;\&Delta;t}})/\mathrm{\&lambda;}\\ e^{-\mathrm{\&lambda;\&Delta;t}}\end{array}\right])$

X_k为k时刻的临近空间目标运动状态向量,k为正整数；

λ为机动时间常数的倒数，即机动频率，取值范围为：转弯机动λ＝1/60，逃避机动λ＝1/20，大气扰动λ＝1；

CS为机动目标当前统计模型；

e为指数函数；

A_CS为单通道当前统计模型状态转移矩阵；

△t为采样周期，取值范围为0.1s～0.5s；

F_CS为当前统计模型状态转移矩阵；

B为机动加速度当前均值转移矩阵；

w_k-1为k-1时刻的状态噪声向量；

X_k-1为k-1时刻的临近空间目标运动状态向量；

为机动加速度“当前”均值，在每一个采样周期内设定为常数。

3.根据权利要求2所述一种跟踪临近空间高超声速目标的修正变结构网格交互多模型滤波方法，其特征在于，所述步骤二中的中心模型采用步骤一中的机动目标当前统计模型，左转弯模型和右转弯模型采用步骤一中的匀速转弯模型；具体过程为：

(1)机动目标当前统计模型速度误差小于2m/s²，最大加速度参数±a_max设置为±40m/s²；

(2)左转弯模型和右转弯模型采用匀速转弯模型

左转弯模型和右转弯模型采用匀速转弯模型，左转弯模型转弯速率为右转弯模型转弯速率为其中，且 的取值为-0.098rad/s，的取值为0rad/s；且 的取值为0rad/s，的取值为0.098rad/s；为左转弯模型转弯速率的最小值；为左转弯模型转弯速率的最大值；为右转弯模型转弯速率的最小值；为右转弯模型转弯速率的最大值。

4.根据权利要求1所述一种跟踪临近空间高超声速目标的修正变结构网格交互多模型滤波方法，其特征在于，所述步骤三中基于步骤一中的惯性参考坐标系确定目标跟踪系统测量模型；具体过程为：

定义方位角为β，高低角为α，方位角会受到量测方位角噪声v_β的污染，高低角会受到量测高低角噪声v_α的污染，v_β和v_α均是零均值高斯白噪声，标准差分别为σ_β和σ_α；

方位角和高低角对于相对位置矢量的非线性函数为：

$\left\{\begin{array}{c}\mathrm{\&alpha;}={\tan }^{-1}\frac{x}{\sqrt{x^2+z^2}}+v_{\mathrm{\&alpha;}}\\ \mathrm{\&beta;}=-{\tan }^{-1}\frac{z}{x}+v_{\mathrm{\&beta;}}\end{array}\right.---\left(4\right)$

式中，

x为惯性参考坐标系x轴坐标；

z为惯性参考坐标系z轴坐标；

v_α为量测高低角噪声；

v_β为量测方位角噪声；

r为观察站到目标S位置的距离，距离测量方程为：

$r=\sqrt{x^2+y^2+z^2}+v_r---\left(5\right)$

式中，y为惯性参考坐标系y轴坐标；；

v_r为距离测量噪声，是零均值高斯白噪声，标准差为σ_r；

则目标跟踪系统测量模型方程为：

$Z_k=\left[\begin{array}{c}{\mathrm{\&alpha;}}_k\\ {\mathrm{\&beta;}}_k\\ r_k\end{array}\right]=\left[\begin{array}{c}{\tan }^{-1}\frac{x_k}{\sqrt{x_k^2+z_k^2}}\\ -{\tan }^{-1}\frac{z_k}{x_k}\\ \sqrt{x_k^2+y_k^2+z_k^2}\end{array}\right]+\left[\begin{array}{c}{v}_{\mathrm{\&alpha;},k}\\ v_{\mathrm{\&beta;},k}\\ v_{r,k}\end{array}\right]=h\left(X_k\right)+\left[\begin{array}{c}{v}_{\mathrm{\&alpha;},k}\\ v_{\mathrm{\&beta;},k}\\ v_{r,k}\end{array}\right]---\left(6\right)$

式中，

Z_k为k时刻的测量向量；

α_k为α在k时刻的值；

β_k为β在k时刻的值；

r_k为r在k时刻的值；

v_α,k为v_α在k时刻的值；

v_β,k为v_β在k时刻的值；

v_r,k为v_r在k时刻的值；

h(X_k)为测量函数在k时刻的值，

$h\left(X_k\right)=\left[\begin{array}{c}{\tan }^{-1}\frac{x_k}{\sqrt{x_k^2+z_k^2}}\\ -{\tan }^{-1}\frac{z_k}{x_k}\\ \sqrt{x_k^2+y_k^2+z_k^2}\end{array}\right]---\left(7\right).$

5.根据权利要求4所述一种跟踪临近空间高超声速目标的修正变结构网格交互多模型滤波方法，其特征在于，所述步骤四中将步骤二中的左转弯模型、中心模型和右转弯模型分别作为状态方程与步骤三中的目标跟踪系统测量模型相结合，采用离散型推广卡尔曼滤波进行计算，并根据计算结果对左转弯模型和右转弯模型的网格中心和网格距离重新调整，然后进行状态估计和误差协方差矩阵融合；具体过程为：

以作为k时刻的系统模型集，经过输入交互及预测估计，求得各模型的似然函数j＝L,C,R，有

${\mathrm{\&Lambda;}}_k^j=N({\tilde{z}}_k^j,0,S_k^j)={\left|2\mathrm{\&pi;}{S}_k^j\right|}^{-\frac{1}{2}}\exp [-\frac{1}{2}{{\tilde{z}}_k^j}^T{\left(S_k^j\right)}^{-1}{\tilde{z}}_k^j]---\left(8\right)$

式中，

j＝L,C,R，为左转弯模型、中心模型和右转弯模型的新息；

为左转弯模型、中心模型和右转弯模型的状态一步预报估计；

为左转弯模型、中心模型和右转弯模型的测量函数一步预报值；

为左转弯模型、中心模型和右转弯模型的新息协方差矩阵；

Z_k为k时刻的测量向量；

上标T为转置运算；

L为左转弯模型；

C为中心模型；

R为右转弯模型；

代表以为变量，均值为0，方差为的高斯分布函数；

k时刻的左转弯模型、中心模型和右转弯模型的交互式多模型后验概率μ_CTL,k、μ_CS,k和μ_CTR,k的计算公式为：

${\mathrm{\&mu;}}_{\mathrm{CTL},k}=\frac{1}{c_{\mathrm{CTL}}}{\mathrm{\&Lambda;}}_k^j{c}_{\mathrm{CTL},j},c_{\mathrm{CTL},j}=\underset{i=L}{\overset{R}{\mathrm{\&Sigma;}}}{p}_{\mathrm{ij}}{\mathrm{\&mu;}}_{\mathrm{CTL},k-1}---\left(9\right)$

${\mathrm{\&mu;}}_{\mathrm{CS},k}=\frac{1}{c_{\mathrm{CS}}}{\mathrm{\&Lambda;}}_k^j{c}_{\mathrm{CS},j},c_{\mathrm{CS},j}=\underset{i=L}{\overset{R}{\mathrm{\&Sigma;}}}{p}_{\mathrm{ij}}{\mathrm{\&mu;}}_{\mathrm{CS},k-1}---\left(10\right)$

${\mathrm{\&mu;}}_{\mathrm{CTR},k}=\frac{1}{c_{\mathrm{CTR}}}{\mathrm{\&Lambda;}}_k^j{c}_{\mathrm{CTR},j},c_{\mathrm{CTR},j}=\underset{i=L}{\overset{R}{\mathrm{\&Sigma;}}}{p}_{\mathrm{ij}}{\mathrm{\&mu;}}_{\mathrm{CTR},k-1}---\left(11\right)$

式中，μ_CTL,k为k时刻的左转弯模型的交互式多模型后验概率，k为正整数；

μ_CS,k为k时刻的中心模型的交互式多模型后验概率；

μ_CTR,k为k时刻的右转弯模型的交互式多模型后验概率；

p_ij为从模型i向模型j转换的概率，i＝L,C,R，j＝L,C,R；i为模型i，j为模型j；

μ_CTL,k-1为k-1时刻的左转弯模型的交互式多模型后验概率；

μ_CS,k-1为k-1时刻的中心模型的交互式多模型后验概率；

μ_CTR,k-1为k-1时刻的右转弯模型的交互式多模型后验概率；

c_CTL为左转弯模型的归一化常数；

c_CS为中心模型的归一化常数；

c_CTR为右转弯模型的归一化常数；

L为左转弯模型；

C为中心模型；

c_CTL、c_CS和c_CTR分别定义为：

$c_{\mathrm{CTL}}=\underset{j=L}{\overset{R}{\mathrm{\&Sigma;}}}{\mathrm{\&Lambda;}}_k^j{c}_{\mathrm{CTL},j}---\left(12\right)$

$c_{\mathrm{CS}}=\underset{j=L}{\overset{R}{\mathrm{\&Sigma;}}}{\mathrm{\&Lambda;}}_k^j{c}_{\mathrm{CS},j}---\left(13\right)$

$c_{\mathrm{CTR}}=\underset{j=L}{\overset{R}{\mathrm{\&Sigma;}}}{\mathrm{\&Lambda;}}_k^j{c}_{\mathrm{CTR},j}---\left(14\right)$

L为左转弯模型；

C为中心模型；

采用离散型推广卡尔曼滤波，离散型推广卡尔曼滤波算法中模型转换方法如下：

(a)若右转弯模型的后验概率μ_CTR>t₂，其中的t₂＝0.9为探测重要模式的阈值,则目标发生了向右的机动，将当前右转弯模型转弯速率向右扩大0.025rad/s，即其中，和分别代表k时刻和k-1时刻的右转弯模型的转弯速率；

若左转弯模型的后验概率μ_CTL>t₂，则目标发生了向左的机动，将当前左转弯模型转弯速率向左扩大0.025rad/s,即其中，和分别代表k时刻和k-1时刻的左转弯模型的转弯速率；

(b)当左转弯模型和右转弯模型的后验概率都小于探测模式的阈值t₁时，t₁＝0.1，即μ_CTL<t₁,μ_CTR<t₁，则目标转弯速率下降，将当前左转弯模型转弯速率和右转弯模型转弯速率一共往中间靠拢一个单位，即 ${\mathrm{\&omega;}}_k^L={\mathrm{\&omega;}}_{k-1}^L+0.0125\mathrm{rad}/s,{\mathrm{\&omega;}}_k^R={\mathrm{\&omega;}}_{k-1}^R-0.0125\mathrm{rad}/s;$

定义有向图D_k是由点与点之间的连线组成的一种数据结构，即

D_k＝(V(D_k),E(D_k))

式中，V(D_k)是非空的顶点集合，E(D_k)为不与V(D_k)相交的边集合，并具有方向性，在有向图D_k中，若从顶点U到顶点V之间有路径，则称U和V是连通的；

设目标机动转弯速率ω在连续范围[-ω_max,ω_max]内变化，针对时变支撑的有向图D_k的交互式多模型算法，网格值k＝0,1,…N,N为正整数，为机动目标当前统计中心模型的转弯速率，-ω_max的取值为-0.098rad/s，ω_max的取值为0.098rad/s，

模型集合即包含左转弯模型、机动目标当前统计中心模型、右转弯模型转弯速率集合；

状态集合 $S_{X_k}=\{X_{L,k},X_{C,k},X_{R,k}\},$

式中，为k时刻左转弯模型、中心模型、右转弯模型构成的的运动状态向量集合；X_L,k为左转弯模型估计出来的k时刻状态，X_C,k为中心模型估计出来的k时刻状态，X_R,k为右转弯模型估计出来的k时刻状态；

设定误差协方差矩阵集合为

式中，为k时刻左转弯模型、中心模型、右转弯模型的状态估计误差协方差矩阵构成的集合，P_L,k为左转弯模型k时刻状态误差协方差矩阵，P_R,k为右转弯模型k时刻状态误差协方差矩阵，P_C,k为中心模型k时刻状态误差协方差矩阵；

模式转换概率矩阵P_LCR为：

$P_{\mathrm{LCR}}=\left[\begin{array}{ccc}{p}_{\mathrm{LL}}& p_{\mathrm{LC}}& p_{\mathrm{LR}}\\ p_{\mathrm{CL}}& p_{\mathrm{CC}}& p_{\mathrm{CR}}\\ p_{\mathrm{RL}}& p_{\mathrm{RC}}& p_{\mathrm{RR}}\end{array}\right]---\left(15\right)$

式中，p_LL为左转弯模型继续保持的概率；

p_LC为由左转弯模型转化为中心模型的概率；

p_LR为由左转弯模型转换为右转弯模型的概率；

p_CL为由中心模型转化为左转弯模型的概率；

p_CC为中心模型继续保持的概率；

p_CR为由中心模型转化为右转弯模型的概率；

p_RL为由右转弯模型转化为左转弯模型的概率；

p_RC为由右转弯模型转化为中心模型的概率；

p_RR为右转弯模型继续保持的概率；

修正变结构交互多模型滤波算法以粗网格：

$D_0=\{{\mathrm{\&omega;}}_0^L=-{\mathrm{\&omega;}}_{\max },{\mathrm{\&omega;}}_0^C=0,{\mathrm{\&omega;}}_0^R={\mathrm{\&omega;}}_{\max }\}---\left(16\right)$

式中，D₀为初始时刻修正变结构交互多模型滤波算法左转弯模型、机动目标当前统计中心模型、右转弯模型的转弯速率集合；

为初始时刻修正变结构交互多模型滤波算法左转弯模型的转弯速率；

ω_max为初始时刻模型转弯速率最大值；

为初始时刻修正变结构交互多模型滤波算法机动目标当前统计中心模型的转弯速率；

为初始时刻修正变结构交互多模型滤波算法右转弯模型的转弯速率；

在每一个循环时间步(k→k+1)网格按照如下四步进行调整：

(1)左转弯模型、中心模型，以及右转弯模型，均采用离散型推广卡尔曼滤波：

$\left\{\begin{array}{c}{\stackrel{\&OverBar;}{X}}_{k|k-1}=F{\hat{X}}_{k-1}+B{\stackrel{\&OverBar;}{a}}_{k-1}\\ P_{k|k-1}=F{P}_{k-1}{F}^T+Q_{k-1}\\ {\hat{X}}_k={\stackrel{\&OverBar;}{X}}_{k|k-1}+K_k[Z_k-h({\hat{X}}_{k|k-1},k)]\\ P_k=(I-K_k{H}_k)P_{k|k-1}{(I-K_k{H}_k)}^T+K_k{R}_k{K_k}^T\\ K_k=P_{k|k-1}{H_k}^T{[H_k{P}_{k|k-1}{H_k}^T+R_k]}^{-1}\end{array}\right.---\left(17\right)$

式中，为从k-1时刻第k时刻的状态一步预报值；

F为子模型状态转移矩阵；

B为机动加速度均值转移矩阵；

为k-1时刻机动加速度均值，在每一个采样周期内设定为常数；

F^T为子模型状态转移矩阵的转置矩阵；

h为系统量测函数；

上标T为转置矩阵符号；

H_k为根据测量方程计算出来的Jacobi矩阵的转置矩阵；

为第k-1步的状态估计值；

P_k|k-1为的状态预报误差协方差矩阵；

P_k-1为的状态估计误差协方差矩阵；

Q_k-1为状态噪声协方差矩阵；

为第k步的状态估计值；

K_k为卡尔曼滤波增益矩阵；

Z_k为第k步的测量向量；

为根据测量方程计算出来的测量向量的一步预报值；

P_k为的状态估计误差协方差矩阵；

I为单位矩阵；

R_k为测量噪声协方差矩阵；

(2)对左转弯模型和右转弯模型的网格中心和网格距离重新调整：

网格中心调整按照下式计算：

${\mathrm{\&omega;}}_k^C={\mathrm{\&mu;}}_{\mathrm{CTL},k}{\mathrm{\&omega;}}_{k-1}^L+{\mathrm{\&mu;}}_{\mathrm{CS},k}{\mathrm{\&omega;}}_{k-1}^C+{\mathrm{\&mu;}}_{\mathrm{CTR},k}{\mathrm{\&omega;}}_{k-1}^R---\left(18\right)$

式中，和分别为k时刻和k-1时刻的中心模型转弯速率；

μ_CTL,k为k时刻左转弯模型的后验概率；

μ_CS,k为k时刻中心模型的后验概率；

μ_CTR,k为k时刻右转弯模型的后验概率；

为k-1时刻左转弯模型转弯速率；

为k-1时刻右转弯模型转弯速率；

左转弯模型和右转弯模型网格距离重新调整分为三种情况：

(a)没有跳跃

k时刻中心模型的后验概率最大，即μ_CS,k＝max{μ_CTL,k,μ_CS,k,μ_CTR,k}，

式中， ${\mathrm{\&lambda;}}_{L,k}=\max \{{\mathrm{\&omega;}}_k^C-{\mathrm{\&omega;}}_{k-1}^L,{\mathrm{\&delta;}}_{\mathrm{\&omega;}}\};$

${\mathrm{\&lambda;}}_{R,k}=\max \{{\mathrm{\&omega;}}_{k-1}^R-{\mathrm{\&omega;}}_k^C,{\mathrm{\&delta;}}_{\mathrm{\&omega;}}\};$

t₁＝0.1为探测模式的阈值；

δ_ω为模型分离距离0.025rad/s；

为k时刻左转弯模型的转弯速率；

为k时刻右转弯模型的转弯速率；

为k时刻中心模型的转弯速率；

λ_L,k为k时刻左转弯模型分离距离最大值；

λ_R,k为k时刻右转弯模型分离距离最大值；

(b)向左跳变

k时刻左转弯模型的后验概率最大，即μ_CTL,k＝max{μ_CTL,k,μ_CS,k,μ_CTR,k}，

${\mathrm{\&omega;}}_k^R={\mathrm{\&omega;}}_k^C+{\mathrm{\&lambda;}}_{R,k}---\left(22\right)$

式中，t₂＝0.9为探测重要模式的阈值；

(c)向右跳变

k时刻右转弯模型的后验概率最大，即μ_CTR,k＝max{μ_CTL,k,μ_CS,k,μ_CTR,k}，

${\mathrm{\&omega;}}_k^L={\mathrm{\&omega;}}_k^C-{\mathrm{\&lambda;}}_{L,k}---\left(23\right)$

(3)状态估计融合X_k和误差协方差矩阵融合P_k

状态估计融合X_k取作左转弯模型、中心模型，以及右转弯模型的运动状态的融合；

$X_k=X_{L,k}{\mathrm{\&omega;}}_k^L+X_{C,k}{\mathrm{\&omega;}}_k^C+X_{R,k}{\mathrm{\&omega;}}_k^R---\left(25\right)$

误差协方差矩阵P_k融合取作左转弯模型、中心模型，以及右转弯模型的状态的误差协方差矩阵P_L,k、P_C,k和P_R,k的融合；

$P_k=P_{L,k}{\mathrm{\&omega;}}_k^L+P_{C,k}{\mathrm{\&omega;}}_k^C+P_{R,k}{\mathrm{\&omega;}}_k^R---\left(26\right)$

式中，为第k时刻左转弯模型的转弯速率；

为第k时刻中心模型的转弯速率；

为第k时刻右转弯模型的转弯速率；

每一次循环对的子模型转移矩阵F需重新计算。