CN104020466A - 基于变结构多模型的机动目标跟踪方法 - Google Patents

Abstract

本发明提出了一种基于变结构多模型的机动目标跟踪方法，主要解决现有技术中目标跟踪精度低和目标跟踪不及时的问题。其实现步骤是：(1)对雷达探测到的目标运动状态进行采样，得到雷达量测值；(2)根据目标运动特性建立机动目标跟踪数学模型；(3)根据目标跟踪数学模型创建机动目标运动模型，并对其初始化；(4)通过模型集合自适应策略，得到每一时刻的目标运动模型集合；(5)根据雷达量测值和目标运动模型集合，运行变结构交互多模型算法，得到机动目标的状态估计，实现目标的一次跟踪；(6)重复执行步骤(4)和步骤(5)，直到跟踪结束。本发明在保证高跟踪精度的前提下，减小了运算量，提高了机动目标跟踪的及时性。

Description

基于变结构多模型的机动目标跟踪方法

技术领域

本发明属于信号处理技术领域，特别涉及机动目标的跟踪方法。本发明可用于在保证跟踪精度的同时能够提高跟踪的及时性。

背景技术

机动目标跟踪一直是跟踪领域研究的热点。根据所用的目标运动模型数，机动目标跟踪方法可分为单模型跟踪和多模型跟踪。多模型方法主要经历了三代。

第一代多模型方法称为自治多模型AMM，由Magill和Lainiotis首先提出，后来由Maybeck等人广泛应用并推广。这种方法的特点是模型个数固定，各个基础滤波器单独运行并独立于其它滤波器，最后对输出进行融合。由于自治多模型AMM方法没有考虑模型之间的跳变，而且各个基于模型的滤波器之间没有交互，因此称为静态多模型方法。这种方法只有在真实模式未知并且在每一时刻只用一个模型做估计的情况下才是有效的。

第二代多模型方法仍然使用固定的模型个数，但各个模型之间存在交互，具有代表性的是Ackerson和Fu于1970年提出的广义伪贝叶斯方法GPB和Blom于1988年提出的交互多模型方法IMM，由于其在目标跟踪上的优秀表现为多模型方法赢得了无数的荣誉，Bar-Shalom促进了多模型的普及和进一步的发展。

前两代多模型方法都使用固定结构的多模型，而其存在很多固有的局限。因此Li X R于1992年提出了变结构多模型的思想和模型集的自适应方法，并在文献中提出了模型群切换方法MGS，它是第一种可普遍应用于一大类混合估计问题的变结构多模型方法，随后又在相应的文献中提出了可能模型集方法LMS和期望模型集扩展方法EMA。Lan jian于2011年提出了最好模型扩展方法BMA，将期望模型集扩展方法EMA进行了推广，使变结构多模型方法可用于模型结构和参数都存在变化的机动过程。

到目前为止，已经提出了很多变结构多模型方法，但是大部分方法都存在一定的缺陷性。期望模式扩展方法EMA是一种操作相对简单、计算复杂度相对较小的方法，但是该方法跟踪精度较低，同时此方法的跟踪精度对目标的机动方式和模型集合的拓扑结构设计依赖程度较大。实验表明，如果期望模型太接近基本模型时，则会产生模型间竞争，从而会导致跟踪性能下降。最好模型集扩展方法BMA可根据KL准则实时产生与模式匹配的最好模型来提高估计的精度。虽然该方法估计精度很高，但是其运算量较大，实现起来较复杂。

发明内容

本发明的目的在于提出一种基于变结构多模型的机动目标跟踪方法，以在跟踪精度和计算复杂度上取的一个平衡，即在保持高跟踪精度的前提下，减小运算量。此变结构多模型方法对于提高机动目标跟踪精度和降低计算复杂度具有一定的实际意义。

本发明的技术方案是：利用现有的最好模型扩展方法BMA的候选模型集和现有的期望模式扩展EMA中求期望模型的方法，构成一种新的模型集合自适应策略。利用该模型集合自适应策略实时的、自适应的调节模型集的大小，从而获得当前时刻的模型集合，并利用通用的变结构交互多模型算法来估计目标的运动状态，完成对机动目标的实时跟踪。具体步骤包括如下：

(1)通过雷达探测机动目标的状态信息即位置信息，对其进行N次采样，得到一个长度为N的量测值序列{z_k},k＝1,2,3...N；

(2)建立机动目标跟踪数学模型：

(2a)用下式对机动目标建立运动状态方程：

x_k＝F_kx_k-1+G_ku_k-1+Γ_kw_k-1

其中x_k表示k时刻机动目标的状态向量，k为采样时刻；F_k表示k时刻的状态转移矩阵；x_k-1表示k-1时刻机动目标的状态向量；G_k表示k时刻的状态输入增益矩阵；u_k-1代表k-1时刻的状态输入；Γ_k表示k时刻的过程噪声增益矩阵；w_k-1表示k-1时刻的过程噪声，其均值和协方差分别为0和Q_k-1的高斯白噪声序列；

(2b)用下式建立机动目标的量测方程：

z_k＝H_kx_k+v_k

其中z_k代表k时刻的雷达量测值，k为采样时刻；H_k表示k时刻的量测矩阵；x_k表示k时刻机动目标的状态向量；v_k表示k时刻的量测噪声，其均值和协方差分别为0和R_k的高斯白噪声序列；

(3)利用机动目标跟踪数学模型进行机动目标的状态估计：

(3a)目标运动模型和模型参数的初始化：

根据机动目标运动特性选择m个目标运动模型作为基本模型集j≤m，选择n个目标运动模型作为候选模型集j≤n，基本模型集合和候选模型集合之间相互独立，每个模型均能作为步骤(2a)状态方程中的一种状态输入；

设在每一个采样周期里均有m+1个运动模型起作用，即每一采样时刻的模型集合均有m个基本模型和1个候选模型。设置m+1个目标运动模型的转移概率矩阵π_ji,i＝1,2,...m+1,j＝1,2,...m+1；

根据目标运动特性设置目标运动模型的初始状态、初始协方差以及其他模型参数：k＝1,2,3,...N,i＝1,2,3,...m+1；

根据机动目标的运动模型个数，初始化目标运动模型概率μ₀＝[1/m,1/m,...1/m,0]_1×(m+1)，其中[·]_1×(m+1)表示一个1行m+1列的矩阵，矩阵中模型概率为0表示此刻这个运动模型不起作用；

(3b)根据模型集合自适应策略，确定每一时刻目标运动模型集合M_k,k＝1,2,3,...N：

(3b1)根据期望模型扩展方法EMA，得出k时刻的期望模型：

${\hat{s}}_k\underset{m^j\∈M_{k-1}}{\mathrm{\Σ}}{m}^j{\mathrm{\μ}}_{k-1}^j,k=\mathrm{1,2,3},...N,$

其中为k-1时刻m^j模型的概率；M_k-1为k-1时刻的模型集合；

(3b2)根据期望模型与候选模型中模型的欧式距离，从候选模型集中选出k时刻与真实运动最匹配的运动模型

${\hat{m}}_k=\underset{m_k^j\∈M_k^c}{\arg \min }\left|\right|{\hat{s}}_k-m_k^j||;$

(3b3)根据最匹配运动模型得到k时刻的目标运动模型集合：

$M_k=M_k^b+{\hat{m}}_k,k=\mathrm{1,2,3}...N;$

(3c)根据步骤(1)中的雷达量测值z_k和步骤(3b)中的目标运动模型集合M_k，运行变结构交互多模型算法VSIMM[M_k,M_k-1]，得到k时刻机动目标的状态估计，完成一个采样点的跟踪；

(3d)判断跟踪是否完成，若k≤N，则跟踪未完成，k递增，进入步骤(3b)继续跟踪，否则目标跟踪过程结束。

本发明与现有技术相比具有以下优点：

1)本发明与现有的期望模式扩展方法EMA相比，由于候选模型与基本模型间不存在模型竞争，所以本发明能够在保证运算量小的前提下提高了目标跟踪精度。

2)本发明与现有的最好模型扩展方法BMA相比，由于在选择最匹配模型时，采用简单的模型距离最小化策略，所以本发明能够在保证跟踪精度的前提下减小了运算量，提高机动目标跟踪的及时性。

附图说明

图1是本发明的实现流程图；

图2是本发明实施例中13个基本模型和12个候选模型的分布图；

图3是本发明与现有的目标跟踪方法对机动目标跟踪的位置误差和速度误差对比图。

具体实施方式

参照图1，本发明具体的步骤如下：

步骤1，获得雷达量测值。

通过雷达探测机动目标的状态信息即位置信息，以T为采样间隔，对其进行N次采样，获得一个长度为N的量测值序列{z_k},k＝1,2,3...N。

步骤2，建立机动目标跟踪数学模型。

基于机动目标状态变化过程为一个随机过程，可将目标跟踪问题建模为一个离散时间随机混合系统，该混合系统的状态方程和量测方程如下：

$x_k=F_k{x}_{k-1}+G_k{u}_{k-1}+{\mathrm{\Γ}}_k{w}_{k-1},$

$z_k=H_k{x}_k+v_k,$

其中x_k表示k时刻机动目标的状态向量，k为采样时刻；F_k表示k时刻的状态转移矩阵；x_k-1表示k-1时刻机动目标的状态向量；G_k表示k时刻的状态输入增益矩阵；u_k-1代表k-1时刻的状态输入；Γ_k表示k时刻的过程噪声增益矩阵；w_k-1表示k-1时刻的过程噪声，其均值和协方差分别为0和Q_k-1的高斯白噪声序列；z_k代表k时刻的雷达量测值，k为采样时刻；H_k表示k时刻的量测矩阵；x_k表示k时刻机动目标的状态向量；v_k表示k时刻的量测噪声，其均值和协方差分别为0和R_k的高斯白噪声序列。

步骤3，对目标运动模型和模型参数进行初始化。

根据机动目标运动特性选择13个目标运动模型作为基本模型集j≤13，模型的数学表述：m^j＝[a_x,a_y]'，选择12个目标运动模型作为候选模型集j≤12，模型的数学表述：m^cj＝[a_x,a_y]'，基本模型和候选模型的分布如图2所示。其中基本模型集合的数学表述如下：

m¹:  a＝[0,0]',    m²:  a＝[20,0]',

m³:  a＝[0,20]',   m⁴:  a＝[-20,0]',

m⁵:  a＝[0,-20]',  m⁶:  a＝[20,20]',

m⁷:  a＝[-20,20]', m⁸:  a＝[-20,-20]',

m⁹:  a＝[20,-20]', m¹⁰:  a＝[40,0]',

m¹¹:  a＝[0,40]',   m¹²:  a＝[-40,0]',

m¹³:  a＝[0,-40]'；

候选模型集合的数学表述如下：

m^c1:a＝[10,10]',    m^c2:a＝[-10,10]'

m^c3:a＝[-10,-10]',    m^c4:a＝[10,-10]'

m^c5:a＝[30,10]',      m^c6:a＝[10,30]'

m^c7:a＝[-10,30]',     m^c8:a＝[-30,10]'

m^c9:a＝[-30,-10]',    m^c10:a＝[-10,-30]'

m^c11:a＝[10,-30]',    m^c12:a＝[30,-10]'；

设在每一个采样周期里均有14个运动模型起作用，即在每一采样时刻的有效模型集合M_k＝{m¹,m²,...m¹³,m^cj},k＝1,2,3...N,j≤12，在形式上，把m^cj模型作为模型集合M_k的第14个模型m¹⁴，该14个目标运动模型的转移概率矩阵如下：

${\mathrm{\π}}_{\mathrm{ji}}=\begin{array}{}\end{array}\left[\begin{array}{cccccccccccccc}287/300& 1/120& 1/120& 1/120& 1/120& 0& 0& 0& 0& 0& 0& 0& 0& 0.01\\ 0.02& 0.9& 0& 0& 0& 0.01& 0& 0& 0.01& 0.01& 0& 0& 0& 0.05\\ 0.02& 0& 0.9& 0& 0& 0.01& 0.01& 0& 0& 0& 0.01& 0& 0& 0.05\\ 0.02& 0& 0& 0.9& 0& 0& 0.01& 0.01& 0& 0& 0& 0.01& 0& 0.05\\ 0.02& 0& 0& 0& 0.9& 0& 0& 0.01& 0.01& 0& 0& 0& 0.01& 0.05\\ 0& 1/30& 1/30& 0& 0& 53/60& 0& 0& 0& 0& 0& 0& 0& 0.05\\ 0& 0& 1/30& 1/30& 0& 0& 53/60& 0& 0& 0& 0& 0& 0& 0.05\\ 0& 0& 0& 1/30& 1/30& 0& 0& 53/60& 0& 0& 0& 0& 0& 0.05\\ 0& 1/30& 0& 0& 1/30& 0& 0& 0& 53/60& 0& 0& 0& 0& 0.05\\ 0& 0.1& 0& 0& 0& 0& 0& 0& 0& 0.85& 0& 0& 0& 0.05\\ 0& 0& 0.1& 0& 0& 0& 0& 0& 0& 0& 0.85& 0& 0& 0.05\\ 0& 0& 0& 0.1& 0& 0& 0& 0& 0& 0& 0& 0.85& 0& 0.05\\ 0& 0& 0& 0& 0.1& 0& 0& 0& 0& 0& 0& 0& 0.85& 0.05\\ 0.01& 0.01& 0.01& 0.01& 0.01& 0.01& 0.01& 0.01& 0.01& 0.01& 0.01& 0.01& 0.01& 0.87\end{array}\right];$

根据目标运动模型数，初始化目标运动模型概率μ₀＝[1/13,1/13,...1/13,0]_1×14，其中[·]_1×14表示一个1行14列的矩阵，矩阵中模型概率为0表示此刻这个运动模型不起作用；

设置目标运动模型参数如下：

$F_k^j=\mathrm{diag}[F,F],G_k^j={\mathrm{\Γ}}_k^j=\mathrm{diag}[G,G]$

$F=\left[\begin{array}{cc}1& T\\ 0& 1\end{array}\right],G=\left[\begin{array}{c}{T}^2/2\\ T\end{array}\right],H_k^j=\left[\begin{array}{cccc}1& 0& 0& 0\\ 0& 0& 1& 0\end{array}\right]$

$w_{k-1}^j~N\left[0{,Q}_{k-1}^j\right],v_k^j~N\[0,R_k^i\]$

其中T为采样间隔；为量测噪声协方差；为过程噪声协方差；分别设置为：T＝1.0s, $R_k^j=1250I{m}^2,Q_{k-1}^j=0.001m^2/s^4.$

步骤4，根据模型集合自适应策略，确定每一时刻目标运动模型集合M_k,k＝1,2,3,...N。

4a)根据期望模型扩展方法EMA，得出k时刻的期望模型：

${\hat{s}}_k\underset{m^j\∈M_{k-1}}{\mathrm{\Σ}}{m}^j{\mathrm{\μ}}_{k-1}^j,k=\mathrm{1,2,3},...N,$

其中为k-1时刻m^j模型的概率；M_k-1为k-1时刻的模型集合；

4b)根据期望模型与候选模型中模型的欧式距离，从候选模型集中选出k时刻与真实运动最匹配的运动模型

${\hat{m}}_k=\underset{m_k^j\∈M_k^c}{\arg \min }\left|\right|{\hat{s}}_k-m_k^j||;$

4c)根据最匹配运动模型得到k时刻的目标运动模型集合：

$M_k=M_k^b+{\hat{m}}_k,k=\mathrm{1,2,3},...N;$

步骤5，根据步骤1中的雷达量测值z_k和步骤4中的目标运动模型集合M_k，运行变结构交互多模型算法VSIMM[M_k,M_k-1]，得到k时刻机动目标的状态估计，完成一个采样点的跟踪。

5a)模型集合的条件初始化：

5a1)根据k-1时刻m^j模型的概率和模型转移概率矩阵π_ji，得出k时刻模型的预测概率和混合权重u^j|i：

${\hat{u}}_{k|k-1}^i=\underset{m^j\∈M_{k-1}}{\mathrm{\Σ}}{\mathrm{\π}}_{\mathrm{ji}}{u}_{k-1}^j,k=\mathrm{1,2,3},...N,$

$u^{j|i}={\mathrm{\π}}_{\mathrm{ji}}{u}_{k-1}^j/{\hat{u}}_{k|k-1}^i;$

5a2)根据步骤5a1)中的混合权重u^j|i，得出混合估计和混合协方差

${\stackrel{\‾}{x}}^i=\underset{m^j\∈M_{k-1}}{\mathrm{\Σ}}{\hat{x}}_{k-1}^j{u}^{j|i},k=\mathrm{1,2,3},...N,$

${\stackrel{\‾}{P}}^i=\underset{m^j\∈M_{k-1}}{\mathrm{\Σ}}[P_{k-1}^j+({\stackrel{\‾}{x}}^i-{\hat{x}}_{k-1}^j){({\stackrel{\‾}{x}}^i-{\hat{x}}_{k-1}^j)}^{\′}]\·u^{j|i},k=\mathrm{1,2,3},...N,$

其中为m^j模型在k-1时刻的状态估计；为m^j模型在k-1时刻的状态估计协方差；M_k-1为k-1时刻的模型集合；(·)′表示矩阵或向量的转置；

5b)模型条件滤波：

5b1)根据步骤5a2)中的混合估计得出机动目标的状态预测 ${\hat{x}}_{k|k-1}^i:$

${\hat{x}}_{k|k-1}^i=F_k^i{\stackrel{\‾}{x}}^i+G_k^i{m}_k^i+{\mathrm{\Γ}}_k^i{w}_{k-1}^i,$

其中为k时刻的状态转移矩阵；为k时刻的状态输入增益矩阵；为k时刻的机动运动模型；为k时刻的过程噪声增益矩阵；为k-1时刻的过程噪声；

5b2)根据步骤5a2)中的混合协方差得出状态协方差预测

$P_{k|k-1}^i=F_k^i{\stackrel{\‾}{P}\left(F_k^i\right)}^{\′}+G_k^i{Q}_{k-1}^i{\left(G_k^i\right)}^{\′},$

其中为k时刻的状态转移矩阵；为k时刻的状态输入增益矩阵；为k-1时刻的过程噪声协方差；

5b3)根据步骤1)中的雷达量测值z_k和步骤5b1)中的机动目标的状态预测得出量测残差

${\tilde{z}}^i=z_k-H_k^i{\hat{x}}_{k|k-1}^i-v_k^i,$

其中为k时刻的量测矩阵；为k时刻的量测噪声；

5b4)根据步骤5b2)中的预测协方差得出量测残差协方差Sⁱ：

$S^i=H_k^i{P}_{k|k-i}^i{\left(H_k^i\right)}^{\′}+R_k^i,$

其中为k时刻的量测矩阵；为k时刻的量测噪声协方差；

5b5)根据步骤5b2)中的预测协方差和步骤5b4)中的量测残差协方差Sⁱ，得出滤波增益Kⁱ：

$K^i=P_{k|k-1}^i{\left(H_k^i\right)}^{\′}{\left(S^i\right)}^{-1};$

5b6)根据步骤5b1)中的状态预测步骤5b3)中的量测残差和步骤5b5)中的滤波增益Kⁱ，得出k时刻的目标运动状态：

${\hat{x}}_k^i={\hat{x}}_{k|k-1}^i+K\stackrel{i~i}{z};$

5b7)根据步骤5b2)中的预测协方差步骤5b4)中的量测残差协方差Sⁱ和步骤5b5)中的滤波增益Kⁱ，得出k时刻的目标运动状态协方差：

$P_k^i=P_{k|k-1}^i-K^i{S}^i{\left(K^i\right)}^{\′};$

5c)模型概率更新：

5c1)根据步骤5b4)中的量测残差协方差Sⁱ和步骤5b5)中的滤波增益Kⁱ，得出模型的似然函数：

$L^i={\left|{2\mathrm{\πS}}^i\right|}^{-1/2}\exp [-(1/2){\left({\tilde{z}}^i\right)}^{\′}{\left(S^i\right)}^{-1}{\tilde{z}}^i];$

5c2)根据步骤5a1)中的mⁱ模型的预测概率和步骤5c1)中的似然函数Lⁱ，得出k时刻的模型概率

$u_k^i=\frac{1}{c}{\hat{u}}_{k|k-1}^i{L}^i,$

$c=\underset{m^j\∈M_k}{\mathrm{\Σ}}{\hat{u}}_{k|k-1}^j{L}^j,k=\mathrm{1,2,3}...N,$

其中L^j为m^j模型的似然函数；为m^j模型的预测概率；M_k为k时刻的目标运动模型集合；

5d)估计融合：

根据步骤5b6)中的目标运动状态和步骤5b7)中的状态协方差和步骤5c2)中的k时刻mⁱ模型的概率得出k时刻目标运动状态的总体估计和总体协方差P_k：

${\hat{x}}_k=\underset{m^i\∈M_k}{\mathrm{\Σ}}{\hat{x}}_k^i{u}_k^i,k=\mathrm{1,2,3},...N,$

$P_k=\underset{m^i\∈M_k}{\mathrm{\Σ}}[P_k^i+({\hat{x}}_k-{\hat{x}}_k^i){({\hat{x}}_k-{\hat{x}}_k^i)}^{\′}]\·u_k^i,k=\mathrm{1,2,3},...N,$

其中为mⁱ模型在k时刻的目标运动状态；为mⁱ模型在k时刻的目标运动状态协方差；M_k为k时刻的目标运动模型集合。

步骤6，判断跟踪是否完成，若k≤N，则跟踪未完成，令k＝k+1，返回步骤4继续跟踪，否则目标跟踪过程结束。

本发明的效果可以通过以下实验进一步说明：

1.仿真条件：

设机动目标真实运动轨迹的初始状态为：x₀＝[8000m,600m/s,8000m,600m/s]，对目标进行连续160s的连续量测。目标机动过程如下：在1～20s做匀速运动，在21～110s做恒等转弯运动，转弯率为1/30rad/s，在111～160s做匀速运动。

在笛卡尔坐标系下，采取200次的Monte Carlo仿真实验。

采用均方根误差RMSE来评估机动目标的跟踪性能。

2.仿真内容：

仿真1，采用本发明方法和现有的交互多模型IMM、期望模型扩展EMA、最好模型扩展BMA目标跟踪方法分别对目标的位置和速度进行跟踪估计，并比较跟踪效果，结果如图3，其中：

图3(a)是本发明和交互多模型IMM、期望模型扩展EMA、最好模型扩展BMA目标跟踪方法在位置上的跟踪效果示意图。从图3(a)中可以看出，本发明的位置跟踪误差明显小于现有的交互多模型IMM和期望模型扩展EMA目标跟踪方法，与现有的最好模型扩展BMA目标跟踪方法的跟踪性能接近。

图3(b)是本发明和交互多模型IMM、期望模型扩展EMA、最好模型扩展BMA目标跟踪方法在速度上的跟踪效果示意图。从图3(b)中可以看出，本发明的速度跟踪误差明显小于现有的交互多模型IMM和期望模型扩展EMA目标跟踪方法，与现有的最好模型扩展BMA目标跟踪方法的跟踪性能接近。

仿真2，将本发明与现有的交互多模型IMM、期望模型扩展EMA、最好模型扩展BMA机动目标跟踪方法的运算量对比，结果如表1：

表1.本发明与现有的机动目标跟踪方法的运算量对比

从表1可以看出，本发明与交互多模型IMM和期望模型扩展EMA方法相比，目标跟踪精度提高，运算量增加不多。本发明与最好模型扩展方法BMA相比，跟踪精度保持相当，运算量比最好模型扩展方法BMA小。由此表明了本发明能在保持高跟踪精度的同时减少运算量，提高目标跟踪的及时性。

Claims (2)

1.一种基于变结构多模型的机动目标跟踪方法，包括如下步骤：

(1)通过雷达探测机动目标的状态信息即位置信息，对其进行N次采样，得到一个长度为N的量测值序列{z_k},k＝1,2,3...N；

(2)建立机动目标跟踪数学模型：

(2a)用下式对机动目标建立运动状态方程：

x_k＝F_kx_k-1+G_ku_k-1+Γ_kw_k-1

其中x_k表示k时刻机动目标的状态向量，k为采样时刻；F_k表示k时刻的状态转移矩阵；x_k-1表示k-1时刻机动目标的状态向量；G_k表示k时刻的状态输入增益矩阵；u_k-1代表k-1时刻的状态输入；Γ_k表示k时刻的过程噪声增益矩阵；w_k-1表示k-1时刻的过程噪声，其均值和协方差分别为0和Q_k-1的高斯白噪声序列；

(2b)用下式建立机动目标的量测方程：

z_k＝H_kx_k+v_k

其中z_k代表k时刻的雷达量测值，k为采样时刻；H_k表示k时刻的量测矩阵；x_k表示k时刻机动目标的状态向量；v_k表示k时刻的量测噪声，其均值和协方差分别为0和R_k的高斯白噪声序列；

(3)利用机动目标跟踪数学模型进行机动目标的状态估计：

(3a)目标运动模型和模型参数的初始化：

根据机动目标运动特性选择m个目标运动模型作为基本模型集j≤m，选择n个目标运动模型作为候选模型集j≤n，基本模型集和候选模型集之间相互独立，每个模型均能作为步骤(2a)状态方程中的一种状态输入；

设在每一个采样周期里均有m+1个运动模型起作用，即m个基本模型和1个候选模型。设置m+1个目标运动模型的转移概率矩阵π_ji,i＝1,2,...m+1,j＝1,2,...m+1；

根据目标运动特性设置目标运动模型的初始状态、初始协方差以及其他模型参数：k＝1,2,3,...N,i＝1,2,3,...m+1；

根据目标运动模型数，初始化目标运动模型概率μ₀＝[1/m,1/m,...1/m,0]_1×(m+1)，其中[·]_1×(m+1)表示一个1行m+1列的矩阵，矩阵中模型概率为0表示此刻这个运动模型不起作用；

(3b)根据模型集合自适应策略，确定每一时刻目标运动模型集合M_k,k＝1,2,3,...N：

(3b1)根据期望模型扩展方法EMA，得出k时刻的期望模型：

${\hat{S}}_k\underset{m^j\∈M_{k-1}}{\mathrm{\Σ}}{m}^j{\mathrm{\μ}}_{k-1}^j,k=\mathrm{1,2,3},...N,$

其中为k-1时刻m^j模型的概率；M_k-1为k-1时刻的模型集合；

(3b2)根据期望模型与候选模型中模型的欧式距离，从候选模型集中选出k时刻与真实运动最匹配的运动模型

${\hat{m}}_k=\underset{m_k^j\∈M_k^c}{\arg \min }\left|\right|{\hat{s}}_k-m_k^j||;$

(3b3)根据最匹配运动模型得到k时刻的目标运动模型集合：

$M_k=M_k^b+{\hat{m}}_k,k=\mathrm{1,2,3},...N;$

(3c)根据步骤(1)中的雷达量测值z_k和步骤(3b)中的目标运动模型集合M_k，运行变结构交互多模型算法VSIMM[M_k,M_k-1]，得到k时刻机动目标的状态估计，完成一个采样点的跟踪；

(3d)判断跟踪是否完成，若k≤N，则跟踪未完成，k递增，进入步骤(3b)继续跟踪，否则目标跟踪过程结束。

2.根据权利要求1所述的基于变结构多模型的机动目标跟踪方法，其中步骤(3c)所述的运行变结构交互多模型算法VSIMM[M_k,M_k-1]，得到k时刻机动目标的状态估计，按如下步骤进行：

(3c1)模型集合的条件初始化

根据k-1时刻m^j模型的概率和模型转移概率矩阵π_ji，得出k时刻模型的预测概率和混合权重u^j|i：

${\hat{u}}_{k|k-1}^i=\underset{m^j\∈M_{k-1}}{\mathrm{\Σ}}{\mathrm{\π}}_{\mathrm{ji}}{u}_{k-1}^j,k=\mathrm{1,2,3},...N,$

$u^{j|i}={\mathrm{\π}}_{\mathrm{ji}}{u}_{k-1}^j/{\hat{u}}_{k|k-1}^i;$

根据混合权重u^j|i，得出混合估计和混合协方差

${\stackrel{\‾}{x}}^i\underset{m^j\∈M_{k-1}}{\mathrm{\Σ}}{\hat{x}}_{k-1}^j{u}^{j|i},k=\mathrm{1,2,3},...N,$

${\stackrel{\‾}{P}}^i\underset{m^j\∈M_{k-1}}{\mathrm{\Σ}}[P_{k-1}^j+({\stackrel{\‾}{x}}^i-{\hat{x}}_{k-1}^j){({\stackrel{\‾}{x}}^i-{\hat{x}}_{k-1}^j)}^{\′}]\·u^{j|i},k=\mathrm{1,2,3},...N,$

其中为m^j模型在k-1时刻的状态估计；为m^j模型在k-1时刻的状态估计协方差；M_k-1为k-1时刻的模型集合；(·)′表示矩阵或向量的转置；

(3c2)模型条件滤波

根据步骤(3c1)中的混合估计，得出机动目标的状态预测

${\hat{x}}_{k|k-1}^i=F_k^i{\stackrel{\‾}{x}}^i+G_k^i{m}_k^i+{\mathrm{\Γ}}_k^i{w}_{k-1}^i,$

其中为k时刻的状态转移矩阵；为k时刻的状态输入增益矩阵；为k时刻的机动运动模型；为k时刻的过程噪声增益矩阵；为k-1时刻的过程噪声；

根据步骤(3c1)中的混合协方差得出状态协方差预测

$P_{k|k-1}^i=F_k^i{\stackrel{\‾}{P}\left(F_k^i\right)}^{\′}+G_k^i{Q}_{k-1}^i{\left(G_k^i\right)}^{\′},$

其中为k时刻的状态转移矩阵；为k时刻的状态输入增益矩阵；为k-1时刻的过程噪声协方差；

根据雷达量测值z_k和机动目标的状态预测得出量测残差

${\tilde{z}}^i=z_k-H_k^i{\hat{x}}_{k|k-1}^i-v_k^i,$

其中为k时刻的量测矩阵；为k时刻的量测噪声；

根据预测协方差得出量测残差协方差Sⁱ：

$S^i=H_k^i{P}_{k|k-i}^i{\left(H_k^i\right)}^{\′}+R_k^i,$

其中为k时刻的量测矩阵；为k时刻的量测噪声协方差；

根据预测协方差和量测残差协方差Sⁱ，得出滤波增益Kⁱ：

$K^i{=P}_{k|k-1}^i{\left(H_k^i\right)}^{\′}{\left(S^i\right)}^{-1};$

根据状态预测滤波增益Kⁱ和量测残差得出k时刻的目标运动状态：

${\hat{x}}_k^i={\hat{x}}_{k|k-1}^i+K\stackrel{i~i}{z};$

根据预测协方差滤波增益Kⁱ和量测残差协方差Sⁱ，得出k时刻的目标运动状态协方差：

$P_k^i=P_{k|k-1}^i-K^i{S}^i{\left(K^i\right)}^{\′};$

(3c3)模型概率更新

根据步骤(3c2)中的滤波增益Kⁱ和量测残差协方差Sⁱ，得出模型的似然函数：

$L^i={\left|{2\mathrm{\πS}}^i\right|}^{-1/2}\exp [-(1/2){\left({\tilde{z}}^i\right)}^{\′}{\left(S^i\right)}^{-1}{\tilde{z}}^i];$

根据似然函数Li和步骤(3c1)中mⁱ模型的预测概率得出k时刻的模型概率：

$u_k^i=\frac{1}{c}{\hat{u}}_{k|k-1}^i{L}^i,$

$c=\underset{m^j\∈M_k}{\mathrm{\Σ}}{\hat{u}}_{k|k-1}^j{L}^j,k=\mathrm{1,2,3}...N,$

其中L^j为m^j模型的似然函数；为m^j模型的预测概率；M_k为k时刻的目标运动模型集合；

(3c4)估计融合

根据步骤(3c2)中的目标运动状态和状态协方差和步骤(3c3)中的k时刻mⁱ模型的概率得出k时刻目标运动状态的总体估计和总体协方差P_k：

${\hat{x}}_k=\underset{m^i\∈M_k}{\mathrm{\Σ}}{\hat{x}}_k^i{u}_k^i,k=\mathrm{1,2,3},...N,$

$P_k=\underset{m^i\∈M_k}{\mathrm{\Σ}}[P_k^i+({\hat{x}}_k-{\hat{x}}_k^i){({\hat{x}}_k-{\hat{x}}_k^i)}^{\′}]\·u_k^i,k=\mathrm{1,2,3},...N,$

其中为mⁱ模型在k时刻的目标运动状态；为mⁱ模型在k时刻的目标运动状态协方差；M_k为k时刻的目标运动模型集合。