近年来，临近空间高超声速飞行器的快速发展对中国的空天防御安全构成了极大的威胁，开展针对此类高速目标的拦截策略研究迫在眉睫^[1]。与传统的飞机类目标以及弹道类目标不同，高超声速再入滑翔式飞行器(Hypersonic Reentry-Glide Vehicle, HRGV)具有高机动、高速、打击范围广等特点，现阶段几乎无防御系统能够完成对其的拦截任务，故研究HRGV对中国的国防事业有着重大意义。

HRGV的机动范围在距地面20~60 km、空气稀薄的临近空间内，飞行速度极快，与大气相互作用形成等离子鞘套并产生湍流尾迹等特殊的电磁/光学现象，对拦截方的探测预警系统的要求极高^[1-2]，高超声速目标的跟踪难度大。但同时热流、动压、过载等硬性条件将高超目标的机动限制在一定范围内^[3-5]，且高超声速滑翔式目标的机动能力有限。这些限制都为高超声速目标的跟踪提供了基础。HRGV的典型代表是美国军队研制的高超声速技术飞行器(Hypersonic Technology Vehicle, HTV-2)和先进高超声速武器(Advanced Hypersonic Weapon, AHW)，其中AHW因其外形结构更适合有效负载，实用性更高。本文以AHW为研究对象展开讨论。

高超声速目标的运动模式多变，运动模式与跟踪模型难以匹配，且目标运动速度快、范围大，因此对目标的跟踪精度难以满足拦截方的要求。为了提高跟踪精度，国内外相关文献主要围绕跟踪方法与目标运动模型展开研究。

在目标运动模型方面，文献[6]针对作周期性的滑跃机动的临近空间高超声速目标，设计了一种运动学模型，该模型将加速度建模为正弦波自相关的零均值随机过程。文献[7]基于“当前”统计模型的思想，提出了一种可实时修正转弯角速度的交互式多模型(IMM)算法，用于高超声速目标跟踪。但以上研究采用了较为单一的模型进行目标跟踪，无法完全覆盖高超声速目标的机动模式，故多采用多模型对高超目标进行跟踪。文献[8-9]将白噪声模型(如常速(CV)模型、常加速度(CA)模型)与Morkov过程模型(如Singer模型)相结合，提出了基于CV-CA-Singer的IMM算法，分别对X-51飞行器的平飞加速巡航和俯冲攻击段以及滑跃式机动飞行的高超声速目标进行跟踪。文献[10]则在文献[6]的基础上针对滑跃式机动的高超声速目标，将正弦波模型与CV^[11]、CA^[12-13]模型组成IMM算法。

然而若采用固定结构的交互式多模型(FSIMM)算法，模型集中很多模型在特定的时间内与系统的有效模式差别较大，且这些“多余”模型的不必要的竞争降低了跟踪精度^[14]，对跟踪模型数量和类型的综合选取限制了FSIMM算法的应用。为此，国外学者提出了基于变结构的多模型估计(VSMM)理论的跟踪方法，以突破FSIMM算法的局限性并提高了跟踪算法的精度。其中，文献[15]提出了3种模型集自适应方法：激活有向图(AD)方法，自适应网格(AG)方法，有向图切换(DS)方法。文献[16]基于DS方法利用加速度与角速度为模型的参数，同时依据模型的后验概率对子模型实时地进行调整，从而实现模型集自适应的变结构算法。但模型可跳转的模式是有限的，即有向图结构是固定，算法的自适应程度受限于有向图结构的复杂程度。与文献[16]相比，文献[17-19]均增强了算法的自适应性，并选取角速度为模型参数。其中，文献[17]提出了利用自适应网格交互式多模型(AGIMM)算法跟踪高机动模型。文献[18]则借鉴了文献[17]的思想，将AGIMM算法运用到对临近空间高超声速目标的跟踪，通过转弯角速度自适应，结合了多模型粒子滤波算法，调整模型集中3个模型的参数，从而完成滤波跟踪算法的自适应。文献[19]基于有向图的切换，根据卡尔曼滤波公式，提出了一种角速度滤波方法，并根据角速度的估计值修正有向图。但文献[17-19]的算法中均存在运动学滤波模型较为单一，无法完全覆盖所有机动模式的不足。AGIMM算法中每个模型的运动学跟踪模型始终采用的是匀速转弯(CT^[20])模型，然而当转弯角速率很小，接近于0时，可将目标运动近似为直线运动，此时，若仍然只采用CT模型为跟踪系统的运动学模型进行滤波，则会对滤波精度造成影响。针对上述AGIMM算法运动学模型单一，无法适用于所有运动模式的问题，本文提出了一种自适应变结构交互多式模型(AVSIMM)算法跟踪高超声速目标，同时选择更优的跟踪运动模型并进行比较，通过对滤波算法与跟踪模型的双重自适应调整，最终改善了目标跟踪算法的跟踪性能。

1 高超声速目标运动学模型

高超声速飞行器在纵向平面内发生机动时，描述其运动的状态方程的形式为

当目标的状态变量为时, 状态量分别为二维平面内x、y方向的位置、速度；F为系统的状态转移矩阵，而Γ则为噪声矩阵；取w_CV~N(0, Q_CV)，w_CT2~N(0, Q_CT2)表示过程演化噪声，Q为协方差噪声，下标CV和CT₂分别表示CV模型和二维CT模型。

此时当系统采用常速(CV)模型和匀速转弯(CT)模型进行目标跟踪时，假设系统的采样时间间隔为T，目标的转弯角速度为ω，则系统的离散化状态转移矩阵为

当目标的状态变量为时，状态量分别为二维平面内x、y方向的位置、速度、加速度，取w_Singer~N(0, Q_Singer)，w_CT3~N(0, Q_CT3)，下标Singer和CT₃分别表示Singer模型和三维CT模型。此时当系统采用CT模型和Singer^[21]模型时，假设系统的采样时间间隔为T，目标的转弯角速度为ω，此时系统的离散化状态转移矩阵为

式中：α为机动时间常数的倒数，而机动时间常数又与目标的机动时间持续的长短有关。

2 AVSIMM算法

为解决AGIMM算法中运动学滤波模型单一的问题，考虑添加CV、CA、Singer等运动学模型，用于描述直线运动的模型作为系统的。但在直接增添运动学模型的同时，会增加模型集中的模型数量，进而加强模型之间的竞争关系，降低跟踪精度，故设计了一种AVSIMM算法：设置一个合适的门限值，当滤波更新得到的角速度值小于门限值时，自适应选择跟踪系统的运动学模型为CV、CA或Singer模型。至此，在模型集参数自适应的基础上，实现了算法结构的自适应。

初始时刻，取能够覆盖所有机动空间的粗网格模型集的3个模型参数为D₀={ω_L, ω_C, ω_R}=。假设k(k=0, 1, …)时刻，目标以角速度ω(未知)在[ω_min, ω_max]范围内机动。模型集中状态向量的集合为X_k={x_k^L, x_k^C, x_k^R}，误差协方差矩阵集合为P_k={p_k^L, p_k^C, p_k^R}，模型转换概率为

在每个循环时间步长(k→k+1)内，模型集按照以下步骤进行自适应调整。

步骤1  更新模型集网格中心的ω_k^C参数:

式中：{ω_k^L, ω_k^C, ω_k^R}为模型集，μ_k^L、μ_k^C和μ_k^R为模型的交互后验概率。

步骤2  根据更新的网格中心ω_k+1^C，调整模型间距。

1) 若μ_k^C=max{μ_k^L, μ_k^C, μ_k^R}，即模型集无跳跃

式中：λ_k^L=max{ω_k^C－ω_k^L, δ_ω}，λ_k^R=max{ω_k^R－ω_k^C, δ_ω}，δ_ω为模型最小间距；t₁为检测阈值。

2) 若μ_k^L=max{μ_k^L, μ_k^C, μ_k^R}，即模型集向左跳跃

式中：t₂为检测阈值。

3) 若μ_k^R=max{μ_k^L, μ_k^C, μ_k^R}，即模型集向右跳跃

至此，获得更新的模型集D_k+1={ω_k+1^L, ω_k+1^C, ω_k+1^R}。

步骤3  k时刻误差协方差矩阵与状态估计值融合。

状态与误差协方差矩阵采用的是交互式多模型算法，具体过程如下：

1) 获得k时刻的混合概率μ_k|k^i|j

式中：为归一化常数; p_ij为模型i转移到模型j的状态转移概率。此处因模型集中有3个模型，故取r=3，下文同理。

2) k时刻每个模型的状态估计和误差协方差矩阵初始化

步骤4  更新状态转移矩阵F_k+1。

根据更新的D_k+1={ω_k+1^L, ω_k+1^C, ω_k+1^R}，将转弯角速度与门限值作比较，再选择合适的状态转移矩阵。取门限G=0.002 rad/s。

1) 若ω_k+1^L < G，则更新状态转移矩阵F_k+1为可描述直线运动的运动学模型，如CV、CA、Singer模型等，否则取状态转移矩阵F_k+1为CT模型。

2) 若ω_k+1^C < G，则更新F_k+1为可描述直线运动的运动学模型，否则取状态转移矩阵F_k+1为CT模型。

3) 若ω_k+1^R < G，则更新F_k+1为可描述直线运动的运动学模型，否则取状态转移矩阵F_k+1为CT模型。

步骤5  扩展卡尔曼滤波。

分别对模型集中的每个模型进行扩展卡尔曼滤波(EKF)，得到每个模型的状态估计值与误差协方差矩阵估计值，即。

系统的状态与误差协方差估计为

系统的量测预测为

量测矩阵的雅可比矩阵为

新息协方差矩阵为

量测一步预测为

状态与误差协方差估计为

式中：i=1, 2, …, r，代表模型集中的第i个模型。

步骤6  k+1时刻误差协方差矩阵与状态估计值融合。

1) 模型似然函数更新：

2) 模型概率更新：

3) 误差协方差矩阵与状态估计值融合：

AVSIMM算法采用了直线运动模型与协调转弯模型自适应，与AGIMM算法相比，转弯协调模型仍采用CT模型，而在描述直线运动时，可选择CV、CA、Singer模型。

Singer模型是一个Markov过程模型，通过选择Singer模型中的参数0 < αT < ∞，能够使Singer模型介于CV模型和CA模型之间，较之二者，Singer模型具有更宽的覆盖面^[14]。同时，在关于目标机动的在线信息不足或不准确的情况下，相比较CV、CA模型，Singer模型是最理想的模型。

下文将把基于CT、Singer模型的AVSIMMⅠ算法与基于CT、CV模型的AVSIMMⅡ算法以及文献[18]中只基于CT模型的AGIMM算法的仿真结果进行比较。同时，考虑到临近空间高超声速目标与普通高机动目标存在一定的差异性，而本文设计的方法的理论基础不是针对的临近空间目标的，故为了体现设计方法的性能，还将与文献[10]中针对临近空间高超声速目标提出的IMM算法作比较。

3 算法对比与仿真

本文以跟踪临近空间高超声速再入滑翔弹AHW为研究背景。因现阶段无AHW的真实实验数据，故利用通用飞行器(Common Aerial Vehicle, CAV)模型^[22]作为参考，结合基本的动力学方程，利用GPOPS^[23]软件近似模拟为AHW的机动轨迹。考虑到AHW主要是在纵向平面内大范围的跳跃机动，本文针对该机动方式研究平面内的目标跟踪问题^{[16, 24-25]}。同时考虑到假设在纵向二维平面内对高超声速飞行器进行跟踪，设目标在(1 500, 60) km处开始机动，仿真得到目标的真实轨迹如图 1所示。

在二维平面内，设机动目标与基站的相对位置关系如图 1所示。图中d为目标与基站之间的相对距离，φ为目标俯仰角，规定逆时针为正。故可得高超声速模型的量测方程为

离散化得

式中：

目标在(1 500, 60) km处开始机动，其真实轨迹如图 2所示。

设定观测基站坐标为(0, 0)，跟踪时间为277 s，采样周期为T=0.5 s，量测距离标准差为50 m，俯仰角的标准差为0.05°，即量测协方差矩阵为diag(0.05 km, 0.05π/180 rad)。同时，考虑到EKF算法的局限性，所有的位置、速度、加速度单位均为km、km/s、km/s²。

AVSIMMⅠ算法中，目标的状态变量为X=，w_CT~N(0, Q_CT)，w_Singer~N(0, Q_Singer)。模型的初始状态估计值与状态更新值均为[1 500+0.1, －3+0.01, 2×10^－4+10^－3, 60+0.1, －0.1+0.01, －9×10^－3+10^－3]^T, 初始模型更新协方差矩阵为diag(0.1, 10^－2, 10^－3, 0.1, 10^－2, 10^－3)。

AVSIMMⅡ算法中，目标的状态变量为X=，w_CT~N(0, Q_CT)，w_CV~N(0, Q_CV)。模型的初始状态估计值与状态更新值均为[1 500+0.1, －3+0.01, 60+0.1, －0.1+0.01]^T，初始模型更新协方差矩阵为diag(0.1, 10^－2, 0.1, 10^－2)。

2种AVSIMM算法对模型集中3个模型的初始化均为D₀={ω₀^L, ω₀^C, ω₀^R}={ω_min, ω_min+={－0.3, －0.05, 0.2}。根据式(2)~式(7)，设置变量初始误差协方差矩阵为diag(0.1, 10^－2, 0.1, 10^－2)。最小模型间距为，检测门限t₁=0.2 s，t₂=0.9 s。

参照文献[18]，AGIMM算法中，目标的状态量为，w_CT~N(0, Q_CT)。状态变量初值和状态更新均为[1500+0.1, －3+0.01, 60+0.1, －0.1+0.01]^T。

参照文献[10]，设计的CV-CA-Singer模型的IMM算法中，目标的状态量为，w_sin~N(0, Q_sin)，w_CV~N(0, Q_CV)，w_CA~N(0, Q_CA)，其中Q_sin表达式如文献[10]所述，取目标的机动频率为1/200 Hz，设目标的加速度方差为0.09²。状态变量初值为[1 500+0.1, －3×10^－4+10^－5, 3×10^－6+10^－6, 60+0.1, －0.1+0.01, －10^－2+10^－3, 1.6×10^－5+10^－5]。

4种算法中，每个模型集包含3个模型，设定初始时刻3个模型的初始概率为[0.01, 0.98, 0.01]^T。初始模型转移概率矩阵为

跟踪目标的位置和速度的均方根误差(RMSE)的计算公式为

3.1 跟踪轨迹对比

根据以上设定的仿真条件，考虑到观测矩阵的非线性性，利用EKF算法得到基于AVSIMMⅠ算法、AVSIMMⅡ算法、AGIMM算法、IMM算法的目标跟踪轨迹。进行100次蒙特卡罗仿真，结合跟踪目标的真实运动轨迹，得到对比如图 3所示。

由图 3易知，4种算法都能够实现对临近空间高超声速目标的有效跟踪。为了完成对4种算法的进一步比较，下面将对4种算法的速度误差、位置误差、速度均方根误差以及位置均方根误差的比较进行具体讨论与分析。

3.2 位置与速度的估计误差结果对比

基于AVSIMMⅠ算法、AVSIMMⅡ算法、AGIMM算法以及IMM算法对目标位置的估计误差如图 4所示，速度误差估计如图 5所示。

结合表 1和图 4可知，AVSIMM算法中x方向与y方向的位置误差要小于AGIMM算法与IMM算法中x方向与y方向的位置误差。由表 1与图 5可知，AVSIMM算法x方向与y方向的速度误差要小于AGIMM算法的x方向与y方向的速度误差，而AGIMM算法在x方向上的速度误差大于IMM算法，而其在y方向上的速度误差则小于IMM算法。而结合表 1与图 4、图 5比较AVSIMMⅠ算法与AVSIMMⅡ算法的位置与速度误差可知，前者的速度与位置误差均小于后者。综上，AVSIMMⅠ算法的速度与位置误差最小，跟踪性能最优。

根据式(34)，可以得到基于AVSIMMⅠ、AVSIMMⅡ、AGIMM算法以及IMM算法的目标位置跟踪的均方根误差(RMSE)，如图 6所示。

根据式(35)，可以得到基于AVSIMMⅠ、AVSIMMⅡ、AGIMM算法、IMM算法的目标速度跟踪的RMSE，如图 7所示。

结合表 2与图 6可知，AVSIMM算法的位置均方根误差要小于AGIMM算法，而AGIMM算法的位置均方根误差则小于IMM算法，结合表 2与图 7可知，AVSIMM算法与AGIMM算法、IMM算法相比，AVSIMM算法的速度均方根误差更小。同时由表 2与图 6、图 7可知，AVSIMMⅠ算法的位置均方根误差与速度均方根误差均小于与AVSIMMⅡ算法。同时由表 2易知，AVSIMM算法仿真的时长虽略长于AGIMM算法，但AVSIMM算法的速度与位置均方根误差更小，且IMM算法的仿真时间最长，故AVSIMM算法在单位时间内的跟踪效率更高, 实用性更强。同样的，AVSIMMⅠ算法的跟踪效率高于AVSIMMⅡ算法。

4 结论

本文设计了一种同时具有结构与参数自适应的交互式多模型的高精度滤波算法，研究了对临近空间高超声速再入滑翔目标的跟踪问题。

通过仿真结果可以得到如下结论：

1) AVSIMM算法在算法结构上明显优于AGIMM算法与文献[10]中设计的针对临近空间高超声速目标的IMM算法。

2) AVSIMMⅠ算法在运动学模型上明显优于AVSIMMⅡ算法。故设计的2种AVSIMM算法中，AVSIMMⅠ算法的跟踪效果最佳。

本文算法重点讨论了对高超声速目标的跟踪情况，但考虑到目标的加速度气动参数信息也与目标机动情况紧密相关，故下一步工作考虑将加速度与气动参数扩维到状态向量中进行滤波算法设计。