0 引　言

无人水面艇的目标跟踪可被视为一个随机过程，而目标跟踪问题的建模是通过数学表达式描述机动目标当前运动模式，并准确估计机动目标位置、速度等运动状态^[1]。常见模型有匀速运动模型（Constant Velocity，CV） 、匀速转弯运动模型（Constant Turn，CT）、匀加速度运动模型（Constant Acceleration，CA） 、Singer模型^[2]以及当前统计（Current Statistical，CS）模型。但单一的运动模型在处理复杂的机动目标跟踪问题上，容易出现较大误差，甚至出现目标丢失情况^[3]。针对单一模型的不足，Bar-Shalom 提出了交互式多模型 （Interactive Multiple Model，IMM）算法^[4]。该算法将目标运动映射到多个已知的模型集中，对每个模型并行滤波。其原理是利用 Markov 概率转移矩阵进行不同子模型之间的交互和切换，并且使得这些模型都相互作用，从各个模型的滤波器得出的状态估计，最后经过加权得到最终状态估计^[5]。

IMM通过组合不同模型进行跟踪，在一定程度上能适应目标多种机动情况，相对于单模型来说，提高了对目标跟踪的自适应性。传统的IMM算法采用的是固定模型集和固定参数，这使得处理复杂场景的机动目标跟踪问题时，只能通过增加模型个数来提高测量的精度。一方面这会降低整体算法的运算速度，在另一方面模型之间的竞争也会降低跟踪的精度问题^[6]。在实际跟踪过程中，以船舶为主的水面障碍物由于强随机性和高动态性，其运动模型会发生频繁切换，直接将具有固定模型集和固定参数的IMM算法应用于无人水面艇目标跟踪场景时，实验表明会造成算法稳定性和跟踪精度不高的问题。

有学者就此提出了一些改进方案。例如变结构多模型（Variable Structure Multiple Model，VSMM）算法^[7]，异构模型集协同滤波（Novel Interactive Multiple Model，Novel-IMM）算法^[8]以及FAIMM（Fuzzy-logic Adaptive Interactive Multiple Model）算法^[9]。尽管VSMM 算法中利用时变模型集来替代IMM中的固定模型集，但是频繁的切换、激活以及初始化模型集会导致滤波估计稳定性下降，影响跟踪效果。Novel-IMM 算法在面对机动目标运动模式时，不需要频繁地切换场景，具有较好的跟踪效果，但是计算量却是成倍增加^[10]。而FAIMM算法经过实验表明，对机动目标运动模式的变化比较敏感，容易输出错误模型，影响跟踪效果。因此，针对无人水面艇目标跟踪场景，除了需要计算结果的实时性，还对跟踪准确性提出了较高要求，以上算法均难以同时满足。

1 交互式多模型算法

IMM算法引入了马尔可夫状态转移矩阵^[11]，内部包含多个模型和与之对应的滤波器，算法将每个模型都参与到状态估计的运算中，各模型之间按照马尔可夫转移概率进行切换^[12]。假设 IMM 采用 N 个运动模型，则模型集$M=\{{M}_{1},\dots {M}_{n}\}$。对于任意k时刻模型${M}_{j}$都有先验概率 ${\mu }_{k}^{j}=P\left\{{M}_{k}^{j}\right\}$。对于下一时刻模型发生切换时，将模型i到模型j的转移概率记作${\pi }^{ij}=P\left\{{M}_{k}^{j}\right|{M}_{k-1}^{i}\Big\}$，而所有模型之间的状态转移概率就组成了模型集的马尔可夫状态转移矩阵。

1.1 输入交互

将上一时刻各运动模型在滤波器中的输出结果进行交互，作为当前迭代过程中的初始状态。${X}_{k-1\mid k-1}^{oj}$为k-1时刻模型j的状态交互值，${X}_{k-1\mid k-1}^{i}$为k-1刻模型i的滤波器状态估计值，${P}_{k-1\mid k-1}^{oj}$为k−1时刻模型j的协方差矩阵交互结果，${P}_{k-1\mid k-1}^{i}$为k−1时刻滤波器i的协方差矩阵输出值。

$X_{k-1\mid k-1}^{oj}=\sum _{i=1}^{N}{X}_{k-1\mid k-1}^{i}\cdot {\mu }_{k-1}^{ij},i,j=1,2,\dots ,n，$

(1)

$\begin{split} {P}_{k-1\mid k-1}^{oj}=& \sum _{i=1}^{N}{\mu }_{k-1}^{ij}\cdot [{P}_{k-1\mid k-1}^{i}+ ({X}_{k-1\mid k-1}^{i}-{X}_{k-1\mid k-1}^{oj})\times\\ & \left.{\left({X}_{k-1\mid k-1}^{i}-{X}_{k-1\mid k-1}^{oj}\right)}^{\mathrm{T}}\right]。\end{split}$

(2)

其中：

${\mu }_{k-1}^{ij}=\frac{1}{\overline{{C}_{j}}}{{\text{π}} }^{ij}{\mu }_{k-1}^{i}，$

(3)

$\overline{{C}_{j}}=\sum _{i=1}^{N}{{\text{π}} }^{ij}\cdot {\mu }_{k-1}^{i}。$

(4)

式中：${\mu }_{k-1}^{ij}$为模型i到模型j的混合转移概率；$\overline{{C}_{j}}$为标准化常数。

1.2 滤波处理

容积卡尔曼滤波器（Cubature Kalman Filter，CKF）^[13]，相比其他非线性滤波器如EKF、UKF等，其数值精度更高，稳定性更好。本文采用文献[14]提出的SRCKF滤波器进行目标状态估计。其过程如下：

步骤1　初始化。给定初始条件${\widehat{\boldsymbol{x}}}_{0}、{\boldsymbol{P}}_{0}、{\boldsymbol{Q}}_{0}、{\boldsymbol{R}}_{0}$。

${\widehat{\boldsymbol{x}}}_{0}=E\left({\boldsymbol{x}}_{0}\right)，$

(5)

${\boldsymbol{P}}_{0}=\mathrm{cov}\left({\boldsymbol{x}}_{0}\right)。$

(6)

步骤2　计算容积点${\xi }_{i}$和响应权值${\omega }_{i}$。三阶容积标准下的 容积点和相应权重为：

$\left\{\begin{split}&{\xi }_{j}=\sqrt{\dfrac{m}{2}}[l{]}_{j}，\\& {\omega }_{j}=\dfrac{1}{m}，\end{split}\right.j=1,2,\cdots m。$

(7)

式中：${m}={n}_{x}$，${n}_{x}$为状态向量的维数，$\left[{l}\right]=[{\boldsymbol{I}}_{n\times n},-{\boldsymbol{I}}_{n\times n}]$，I为n维单位阵；${\left[{l}\right]}_{j}$为第j列向量。

步骤3　时间更新

1）对k-1时刻的协方差矩阵进行Cholesky 分解

${\boldsymbol{P}}_{k-1}={{{\boldsymbol{S}}}}_{k-1}{{{\boldsymbol{S}}}}_{k-1}^{{\rm{T}}}。$

(8)

2）计算状态预测容积点和传播容积点

${\boldsymbol{\zeta }}_{j,k-1}={\boldsymbol{S}}_{k-1}{\boldsymbol{\xi }}_{j}+{\widehat{\boldsymbol{x}}}_{k-1}，$

(9)

${\boldsymbol{X}}_{j,\frac{k}{k}-1}^{\mathrm{*}}=f\left({\zeta }_{j,k-1}\right)。$

(10)

3）状态预测估计和预测估计协方差

$\left\{\begin{aligned} & {\widehat{\boldsymbol{x}}}_{k|k-1}=\sum _{j=1}^{m}{\omega }_{j}{\boldsymbol{X}}_{j,\frac{k}{k}-1}^{\mathrm{*}}，\\ & {\boldsymbol{S}}_{k\mid k-1}=\mathrm{Tria}\left(\left[{\chi }_{k\mid k-1}^{\mathrm{*}},{\boldsymbol{S}}_{Q,k-1}\right]\right)。\end{aligned}\right.$

(11)

式中：$\mathrm{T}\mathrm{r}\mathrm{i}\mathrm{a}(\cdot )$为利用QR分解取下三交矩阵，${\boldsymbol{S}}_{Q,k-1}= \mathrm{c}\mathrm{h}\mathrm{o}\mathrm{l}\left({\boldsymbol{Q}}_{k-1}\right)$为取矩阵${\boldsymbol{Q}}_{k-1}$的平方根。

步骤4　量测更新

1）计算量测容积点和传播容积点

$\left\{\begin{aligned} & {\zeta }_{j,k|k-1}={\boldsymbol{S}}_{k|k-1}{\boldsymbol{\xi }}_{j}+{\widehat{\mathbf{x}}}_{k|k-1}，\\ & {\boldsymbol{y}}_{j,k|k-1}^{\mathrm{*}}=h\left({\zeta }_{j,k|k-1}\right)。\end{aligned}\right.$

(12)

2）量测估计和量测估计协方差

$\left\{\begin{aligned} & {\widehat{\boldsymbol{y}}}_{k|k-1}=\sum _{j=1}^{m}{\omega }_{j}{\boldsymbol{y}}_{\mathrm{j},\mathrm{k}|\mathrm{k}-1}^{\mathrm{*}}，\\ & {\boldsymbol{S}}_{yy}=Tria\left(\right[{\eta }_{k\mid k-1},\sqrt{{\boldsymbol{R}}_{k-1}}\left]\right)。\end{aligned}\right.$

(13)

$\begin{split} {\eta }_{k\mid k-1}=& \frac{1}{\sqrt{2n}}\left[\left({\boldsymbol{y}}_{1,{k}|{k}-1}^{\mathrm{*}}-{\widehat{\boldsymbol{y}}}_{k|k-1}\right),\left({\boldsymbol{y}}_{2,{k}|{k}-1}^{\mathrm{*}}-{\widehat{\boldsymbol{y}}}_{k|k-1}\right),\cdots ,\right.\\ & \left.\left({\boldsymbol{y}}_{2n,k|k-1}^{\mathrm{*}}-{\widehat{\boldsymbol{y}}}_{k|k-1}\right)\right]。\\[-10pt] \end{split}$

(14)

3）获得互协方差矩阵

${\boldsymbol{P}}_{xy}={\boldsymbol{\gamma }}_{k\mid k-1}{\left({\boldsymbol{\eta }}_{k\mid k-1}\right)}^{{\rm{T}}} ，$

(15)

$\begin{split} {\gamma }_{k\mid k-1}=& \frac{1}{\sqrt{2n}}\left[\left({\zeta }_{j,k|k-1}-{\widehat{\mathbf{x}}}_{k|k-1}\right),\left({\zeta }_{j,k|k-1}-{\widehat{\mathbf{x}}}_{k|k-1}\right),\cdots ,\right.\\ & \left.\left({\zeta }_{j,k|k-1}-{\widehat{\mathbf{x}}}_{k|k-1}\right)\right]。\\[-10pt] \end{split}$

(16)

4）增益更新

${\boldsymbol{K}}_{k}=({\boldsymbol{P}}_{xy}∕{\boldsymbol{S}}_{yy}^{{\rm{T}}})∕{\boldsymbol{S}}_{yy} 。$

(17)

5）计算状态估计值和平方根形式的估计协方差

$\left\{\begin{aligned} & {\widehat{\boldsymbol{x}}}_{k}={\widehat{\boldsymbol{x}}}_{k|k-1}+{\boldsymbol{K}}_{k}\left({\boldsymbol{y}}_{k}-{\widehat{\boldsymbol{y}}}_{k|k-1}\right)，\\ & {\boldsymbol{S}}_{k}=\mathrm{Tria}\left(\left[{\gamma }_{k\mid k-1}-{\boldsymbol{K}}_{k}{\eta }_{k\mid k-1},{\boldsymbol{K}}_{k}\sqrt{{\boldsymbol{R}}_{k-1}}\right]\right)。\end{aligned}\right.$

(18)

1.3 概率更新

通常假设模型的似然函数服从正态分布。因此k时刻模型j匹配的似然函数以及更新后的模型概率分别为

$\left\{\begin{aligned} & {\mathrm{\Lambda }}_{k}^{j}=\frac{1}{\sqrt{\left|2{\text{π}} {\mathbf{S}}_{k}^{j}\right|}}\mathrm{exp}\left[-\frac{1}{2}{\left({\boldsymbol{e}}_{k}^{j}\right)}^{{\rm{T}}}{\left({\boldsymbol{S}}_{k}^{j}\right)}^{-1}{\boldsymbol{e}}_{k}^{j}\right]，\\ & {u}_{k}^{j}=\frac{{\mathrm{\Lambda }}_{k}^{j}\overline{{C}_{j}}}{\displaystyle\sum _{i=1}^{n}{\mathrm{\Lambda }}_{k}^{i}\overline{{C}_{i}}}。\end{aligned}\right.$

(19)

其中，滤波残差${e}_{k}^{j}$及其协方差${\boldsymbol{S}}_{k}^{j}$可以表示为：

$\left\{\begin{aligned} & {e}_{k}^{j}={\boldsymbol{y}}_{k}-{\widehat{\boldsymbol{y}}}_{k|k-1}，\\ & {\boldsymbol{S}}_{k}^{j}={\boldsymbol{S}}_{zz,k|k-1}{\boldsymbol{S}}_{zz,k|k-1}^{{\rm{T}}}。\end{aligned}\right.$

(20)

1.4 数据融合

将各个滤波器输出的估计值根据模型概率进行加权融合，得到IMM算法的最终状态估计输出${\widehat{\boldsymbol{x}}}_{k}$和对应的协方差矩阵${\boldsymbol{P}}_{k}$。

$\left\{\begin{aligned} & {\widehat{\boldsymbol{x}}}_{k}=\sum _{j=1}^{r}{u}_{k}^{j}{\widehat{\boldsymbol{x}}}_{k}^{j}，\\ & {\boldsymbol{P}}_{k}=\sum _{j=1}^{r}{u}_{k}^{j}\left[{\boldsymbol{P}}_{k}^{j}+\left({\widehat{\boldsymbol{x}}}_{k}^{j}-{\boldsymbol{x}}_{k}\right){\left({\widehat{\boldsymbol{x}}}_{k}^{j}-{\boldsymbol{x}}_{k}\right)}^{{\rm{T}}}\right]。\end{aligned}\right.$

(21)

假设N=3，则IMM算法结构图如图1所示。

IMM算法综合考虑了各模型之间的切换，相比传统单模型跟踪，提供了更高的跟踪精度和更快的收敛速度。当机动目标包含多种运动模式时，IMM模型需包含尽可能多的运动模型才能更好匹配目标运动状态，而过多的运动模型也意味着非匹配模型的竞争，引入误差，导致跟踪精度下降的同时，增大算法运算量。

2 SWA-IMM

为解决IMM算法固定模型集带来的局限性问题，同时满足无人艇在水面机动目标跟踪场景下对跟踪精度的要求，水面机动目标常见且高频的3种运动模型分别是CV、CT、CA模型，令模型数量n=3，${M}= \{ \mathrm{CV},\mathrm{CT},\mathrm{CA}\}$，其中${M}^{1}=\mathrm{C}\mathrm{V},{M}^{2}=\mathrm{C}\mathrm{T},{M}^{3}=\mathrm{C}\mathrm{A}$。因此本文提出的SWA-IMM算法设立3种模型权重不同的滑动窗口，分别为${\mathrm{C}\mathrm{V}}^{*}-\mathrm{C}\mathrm{T}-\mathrm{C}\mathrm{A}$、${\mathrm{C}\mathrm{T}}^{*}- \mathrm{C}\mathrm{V}-\mathrm{C}\mathrm{A}$、${\mathrm{C}\mathrm{A}}^{*}-\mathrm{C}\mathrm{V}-\mathrm{C}\mathrm{T}$窗口，其中${\boldsymbol{M}}^{\boldsymbol{*}}\in \{{\boldsymbol{M}}^{1},{\boldsymbol{M}}^{2},{\boldsymbol{M}}^{3}\}$代表窗口中权重最大的模型，窗口中其余模型的权重相同。通过维护一定数据长度的窗口，根据滑动窗口中输出的模型集权重增益，动态更新主模型集的概率转移矩阵，达到模型集自适应更新的同时，改善跟踪效果。

根据模型权重，初始化每个独立滑动窗口模型集的模型概率转移矩阵${{PP}}\left({{w}}\right)$，将k-1时刻的融合状态估计及对应的协方差矩阵作为先验输入，并经过IMM算法得到k时刻目标估计值${\widehat{\boldsymbol{x}}}_{k}\left(w\right)$以及协方差矩阵${P}_{k}\left(w\right)$。其中$w=\mathrm{1,2},3$ 代表不同的窗口。3个滑动窗口并行且独立的工作，利用滑窗判定方法，从3个滑动窗口输出值中选择目标机动状态最优匹配估计结果。

滑窗判定方法核心是m/n滑窗逻辑法^[15]，序列 ${\widehat{\boldsymbol{x}}}_{k+1}\left(w\right),{\widehat{\boldsymbol{x}}}_{k+2}\left(w\right)$,…,${\widehat{\boldsymbol{x}}}_{k+i}\left(w\right)$,…,${\widehat{\boldsymbol{x}}}_{k+n}\left(w\right)$ 表示含n次扫描的时间窗输入，即当前模型集经过IMM算法得到的目标估计值序列，如图2所示。

SWA-IMM算法流程图如图3所示。

2.1 初始化

由于滑动窗口中，模型集前几个时刻先验信息及后验信息较少，而SWA-IMM算法的实现依赖于窗口内长度为n的序列。因此前n-1时刻$\{0\;\mathrm{s},1\;\mathrm{s},\cdots ,(\mathrm{n}-2)\;\mathrm{s}, ({n}-1\left)\;\mathrm{s}\right\}$取3个窗口状态估计值的算数平均值作为估计输出,则在前n−1时间序列的任意k时刻，SWA-IMM算法输出的状态估计值${\overline{\boldsymbol{X}}}_{k}$及其误差协方差${\overline{\boldsymbol{P}}}_{k}$可表示为：

$\begin{split} & {\overline{\boldsymbol{X}}}_{k}=\frac{1}{3}{\left[\begin{array}{ccc}{\lambda }_{cv}& {\lambda }_{ct}& {\lambda }_{ca}\end{array}\right]}^{{\rm{T}}}\left[\begin{array}{ccc}{\widehat{\boldsymbol{x}}}_{k}\left(1\right)& {\widehat{\boldsymbol{x}}}_{k}\left(2\right)& {\widehat{\boldsymbol{x}}}_{k}\left(3\right)\end{array}\right]，\\ & {\overline{\boldsymbol{P}}}_{k}=\frac{1}{3}{\left[\begin{array}{ccc}{\lambda }_{cv}& {\lambda }_{ct}& {\lambda }_{ca}\end{array}\right]}^{{\rm{T}}}\left[\begin{array}{ccc}{\boldsymbol{P}}_{k}\left(1\right)& {\boldsymbol{P}}_{k}\left(2\right)& {\boldsymbol{P}}_{k}\left(3\right)\end{array}\right]。\end{split}$

(22)

式中，${\left[\begin{array}{ccc}{\lambda }_{cv}& {\lambda }_{ct}& {\lambda }_{ca}\end{array}\right]}^{{\rm{T}}}={\left[\begin{array}{ccc}1& 1& 1\end{array}\right]}^{{\rm{T}}}$，${\lambda }_{cv}$、${\lambda }_{ct}$、${\lambda }_{ca}$分别代表不同窗口的滑窗判定系数。

2.2 滑窗判定

当窗口内序列长度大于n时，可进行滑窗判定。输入目标状态以及状态协方差矩，并执行IMM算法得到状态估计${\widehat{\boldsymbol{x}}}_{k}\left(w\right)$及其协方差矩阵${\boldsymbol{P}}_{k}\left(w\right)$；计算各模型概率值${\mathrm{\mu }}_{{w}{j}}$，记${\mathrm{\mu }}_{{w}}$为模型集中模型概率的最大值。其中${M}_{wj}$为第w个滑动窗口的第j个模型；${\mu }_{wj}$为第w个滑动窗口的第j个模型的概率值；${\mu }_{w}$为第w个滑动窗口模型集中模型概率的最大值。

滑窗判定主要是将窗口中，概率值最大值对应的模型${M}_{wj}$与当前窗口权重最大的模型${\boldsymbol{M}}_{\boldsymbol{w}}^{\boldsymbol{*}}$进行比较。当模型匹配时，说明当前目标运动状态正处于或者倾向于运动模型${\boldsymbol{M}}^{\boldsymbol{*}}$。此时应该将模型集设定为$\left\{{M}^{*},{M}^{i}\right\}$，其中${M}^{i}\in M$且${M}^{i}\ne {M}^{*}$，相比于传统的固定模型集$M$，能够根据目标运动状态的先验、后验信息，实现动态调整模型集。滑窗判定方法如下：

1）获取最佳匹配模型

检测数m和序列n一起构成滑动窗口的逻辑判定。对于序列中的n个时刻，获取每个时刻${\mu }_{w}$所对应的模型，作为最佳匹配模型，记作G(${\mu }_{w}$)。统计n个时刻的最佳匹配模型数量分布，如果当前窗口中最佳匹配模型为权重模型的同时，其统计数量大于等于检测数m，则当前滑动窗口的权重模型应判定有效，记作$F\left(w,m\right)$，$F\left(w,m\right)$保留了最佳匹配模型结果。

$F\left(w,m\right)=\left\{\begin{array}{cc}\varnothing，& count\left({\boldsymbol{M}}_{\boldsymbol{w}}^{\boldsymbol{*}}\right) < m，\\ {\boldsymbol{M}}^{\boldsymbol{*}}，& count\left({\boldsymbol{M}}_{\boldsymbol{w}}^{\boldsymbol{*}}\right)\geqslant m。\end{array}\right.$

(23)

2）设定滑窗判定系数

${\lambda }_{cv},{\lambda }_{ct},{\lambda }_{ca}=\left\{\begin{array}{c}1,F\left(w,m\right)={\boldsymbol{M}}^{\boldsymbol{*}}，\\ 0 ,F\left(w,m\right)\ne {\boldsymbol{M}}^{\boldsymbol{*}}。\end{array}\right.$

(24)

式中：${\lambda }_{cv}、{\lambda }_{ct}、{\lambda }_{ca}$对应的最大权重模型${\boldsymbol{M}}^{\boldsymbol{*}}$分别为CV、CT、CA。

2.3 融合输出结果

当多个滑动窗口均判定其权重模型有效时，取所有窗口结果的算术平均值作为结果。

$\begin{split} {\overline{\boldsymbol{X}}}_{k}=& \dfrac{1}{{\lambda }_{cv}+{\lambda }_{ct}+{\lambda }_{ca}}{\left[\begin{array}{ccc}{\lambda }_{cv}& {\lambda }_{ct}& {\lambda }_{ca}\end{array}\right]}^{{\rm{T}}}\times\\ & \left[\begin{array}{ccc}{\widehat{\boldsymbol{x}}}_{k}\left(1\right)& {\widehat{\boldsymbol{x}}}_{k}\left(2\right)& {\widehat{\boldsymbol{x}}}_{k}\left(3\right)\end{array}\right],\\ {\overline{\boldsymbol{P}}}_{k}=& \dfrac{1}{{\lambda }_{cv}+{\lambda }_{ct}+{\lambda }_{ca}}{\left[\begin{array}{ccc}{\lambda }_{cv}& {\lambda }_{ct}& {\lambda }_{ca}\end{array}\right]}^{{\rm{T}}}\times\\ & \left[\begin{array}{ccc}{\boldsymbol{P}}_{k}\left(1\right)& {\boldsymbol{P}}_{k}\left(2\right)& {\boldsymbol{P}}_{k}\left(3\right)\end{array}\right]。\end{split}$

(25)

3 仿真与分析 3.1 目标运动轨迹

假设目标在x和y方向的初始位置分别为[100 m, 100 m]，初始加速度为[0.8 m/s², 0.8 m/s²]。目标的轨迹由多段不同运动状态下的路径组成：该目标在0～20 s时做加速度[0.8 m/s², 0.8 m/s²]的匀加速度直线运动；20～150 s 时做角速度为0.05 rad/s的匀速转弯运动；150～180 s时做加速度为[1 m/s², 1 m/s²]的匀加速度运动。产生的机动目标运动仿真轨迹如图4所示。

3.2 算法验证

仿真中目标运动存在匀速直线运动、匀速转弯运动、匀加速直线运动3种状态。假设采样间隔T=1 s，设定量测向量${Z}={\left[{p}_{x},{p}_{y}\right]}^{{\rm{T}}}$, 其中${p}_{x}$、${p}_{y}$ 分别为目标x方向和y方向上的位置；量测矩阵${H}=[\mathrm{1,0};\mathrm{0,1}]$，假定量测为白噪声序列，在真实值中加入均值为 0、方差为${\left(10\;{\rm{m}}\right)}^{2}$的随机噪声得到仿真的测量值，即${R}={\left(10\;{\rm{m}}\right)}^{2}$。假设3个模型的初始时刻概率值为1/3，马尔科夫概率转移矩阵为：

$P=\left[\begin{array}{ccc}0.90& 0.25& 0.25\\ 0.25& 0.90& 0.25\\ 0.25& 0.25& 0.90\end{array}\right]。$

(26)

对于SAC-IMM算法，设定滑动窗口中模型概率初值为：

$\varPi =\left[\begin{array}{ccc}0.6& 0.2& 0.2\\ 0.2& 0.6& 0.2\\ 0.2& 0.2& 0.6\end{array}\right]。$

(27)

设定滑动窗口中概率转移矩阵$PP=[{PP}_{1},{PP}_{2},{PP}_{3}]$，其中${{P}{P}}_{i}$ 为：

${PP}_{i}=\left[\begin{array}{ccc}0.9& 0.05& 0.05\\ 0.05& 0.9& 0.05\\ 0.05& 0.05& 0.9\end{array}\right]，i=1,2,3。$

(28)

图5为笛卡尔坐标系下 IMM、SWA-IMM 算法的跟踪轨迹。可知，SWA-IMM算法相比传统IMM算法在机动目标发生运动状态改变时，能够更加有效地结合量测结果，对目标进行轨迹跟踪。

图6所示为目标真实运动轨迹、量测轨迹以及2种滤波算法轨迹的局部对比图。可以发现，无论目标处于匀速转弯运动阶段还是匀加速运动阶段，本文提出的SWA-IMM滤波跟踪算法相比传统IMM算法，在滤波轨迹上更接近于真实轨迹。

3.3 性能分析

对于2种算法性能，使用目标位置、速度、加速度的均方根误差作为评价指标。

${\boldsymbol{RMSE}}=\sqrt{\frac{1}{M}\sum _{j=1}^{M}{\left({\boldsymbol{X}}_{i}\left(k\right)-{\overline{\boldsymbol{X}}}_{i}\left(k\right)\right)}^{2}} 。$

(29)

式中：${\overline{\boldsymbol{X}}}_{i}\left(k\right)$和${\boldsymbol{X}}_{i}\left(k\right)$分别表示k时刻的状态估计值和真实值，M表示蒙特卡罗次数。

根据表1中的数据分析可知，SWA-IMM算法相较于固定模型集IMM算法，在位置、速度和加速度性能上分别提高了33%、39%和17%。

分析图7～图9可知，目标处于弱机动或强机动状态时，本文提出的算法相比较传统IMM算法无论在位置还是速度亦或是加速度上的跟踪误差都有所提升，在匹配目标运动模型时效果良好。这是因为SWA-IMM算法通过记录历史信息，当滑动窗口判定目标运动状态发生改变时，及时调窗口模型集中失配模型和配准模型的权重比，即提高配准模型的概率，从而比固定模型集得到更加准确的匹配结果。由图10可看出，在目标运动的0～30 s和150～180 s内，权重模型的概率分布基本接近目标实际运动状态；在目标运动的130～150 s内，即目标处于匀加速转弯状态时，由于滤波结果不稳定，多个窗口满足滑窗判定，最终的模型集为多个窗口模型集融合结果，因此最佳模型的匹配度相比于其他运动阶段会有所下降。

4 结　语

本文针对无人水面艇目标跟踪场景，对固定模型集的交互式多模型算法进行改进，提出了SWA-IMM算法。实验表明，相比传统IMM算法，SWA-IMM 算法能更加准确获取运动目标速度、加速度等信息，使得无人艇在跟踪水面机动目标时能更加快速的匹配当前目标的运动状态，并完成目标的跟踪定位以及运动信息获取。