1、第 45 卷第 2 期 2024 年 2 月宇 航 学 报Journal of AstronauticsNo.22024FebruaryVol.45面向高超声速滑翔目标的多模型多意图融合轨迹预测李佳丽1,2,郭杰1,2,唐胜景1,2(1.北京理工大学宇航学院,北京 100081;2.飞行器动力学与控制教育部重点实验室(北京理工大学),北京 100081)摘要:针对高超声速滑翔目标速度快、机动强导致轨迹预测精度低且窗口有限的问题,提出了一种多模型多意图融合轨迹预测方法。首先从目标的控制规律出发,构建了纵向、侧向解耦的时变预测模型集,以升力参数、阻力参数与倾侧角幅值为特征参数,通过最小二乘法与高斯
2、假设设计纵向模型;以倾侧角符号为特征参数,针对攻击意图、规避禁飞区以及C/S型机动等典型动作设计侧向模型。然后考虑意图距离、航向角与规避禁飞区等因素设计意图代价函数,基于贝叶斯理论推导预测模型与意图的后验概率递推公式,通过多机动模式、多意图的双层次融合提高轨迹预测准确性与鲁棒性。仿真结果表明,该方法能够应对高超声速滑翔目标的复杂机动,有效提高了长时轨迹预测精度。关键词:高超声速滑翔飞行器;轨迹预测;机动模式推断;意图推断;贝叶斯理论中图分类号:V249.1 文献标识码:A 文章编号:1000-1328(2024)02-0167-14 DOI:10.3873/j.issn.1000-1328.2
3、024.02.002Trajectory Prediction Based on Multi-model and Multi-intent Fusion for Hypersonic Gliding TargetsLI Jiali1,2,GUO Jie1,2,TANG Shengjing1,2(1.School of Aerospace Engineering,Beijing Institute of Technology,Beijing 100081,China;2.Key Laboratory of Dynamics and Control of Flight Vehicle(Beijin
4、g Institute of Technology),Ministry of Education,Beijing 100081,China)Abstract:To address the problem of low trajectory prediction accuracy and limited prediction window caused by the high speed and strong maneuverability of hypersonic gliding vehicles,a trajectory prediction method based on multi-m
5、odel and multi-intent fusion is proposed.By analyzing the control mechanism of targets,a time-varying prediction model set with longitudinal and lateral decoupling is constructed.Taking the lift variable,the drag variable and the amplitude of bank angle as the characteristic parameters,longitudinal
6、maneuver models are designed through the least square method and Gaussian assumption.Taking the sign of bank angle as the characteristic parameter,lateral maneuver models are designed for typical actions such as attacking intents,avoiding no-fly zones,C-type and S-type maneuvers.Furthermore,the inte
7、nt cost function is designed by combining the intents distance,heading angles and avoiding no-fly zones.The posterior probability recursive formulas of maneuver patterns and intents are derived based on the Bayesian theory.Through double-level fusion of multiple maneuver patterns and multiple intent
8、s,the accuracy and robustness of trajectory prediction are improved.Simulation results show that the proposed method can cope with complex maneuvers of hypersonic gliding targets and effectively improve the performance of long-term trajectory prediction.Key words:Hypersonic gliding vehicles;Trajecto
9、ry prediction;Maneuver model inference;Intent inference;Bayesian theory收稿日期:2023-06-18;修回日期:2023-09-09基金项目:国家自然科学基金(11572036);上海航天科技创新基金(SAST201711)宇航学报第 45 卷0引言高超声速滑翔飞行器指飞行马赫数超过5的飞行器,主要在距离海平面20100 km的临近空间飞行。它具有飞行速度快、作战范围广、机动能力强等特点,因而与传统目标相比,其轨迹可预测性大大降低,轨迹预测难度显著增加1-2。根据目标潜在信息挖掘利用方式的不同,可将轨迹预测方法分为基于统计
10、分析的轨迹预测、基于运动机理的轨迹预测与基于意图推断的轨迹预测3类。基于统计分析的轨迹预测方法通过目标历史信息分析,获取飞行轨迹的统计学特征,并据此预测未来轨迹。韩春耀等1将轨迹序列分解为趋势项、周期项与随机项,通过加法或乘法模型集成获取预测轨迹。Li等2通过Hough变换将目标横向、纵向轨迹分别描述为二次函数与复合函数(线性项+周期项)进行预测。张博伦等3将侧向加速度拟合为三次函数,假定其他方向加速度为常值,通过数值积分解算预测轨迹。随着机器学习技术的发展,部分学者将神经网络应用于轨迹预测领域。Xie等4将轨迹预测考虑为多元时间序列预测问题,构造包含非线性模式与线性模式的双通道双向神经网络挖
11、掘轨迹规律。类似地,杨春伟等5利用长短期记忆网络解决序列预测问题,基于注意力机制构建序列到序列预测模型。Sun等6在无跟踪误差的假设下推导了预测模型,并利用长短期记忆网络估计参数范围,实现飞行管道预测。基于运动机理的轨迹预测方法通过挖掘和利用目标的运动机理信息解算预测轨迹。李世杰等7将控制变量拟合为多项式和三角函数的复合函数,将其代入动力学模型积分得到预测轨迹。Moon等8采用多模型交互算法估计目标状态,通过正则化线性回归、支持向量回归和高斯过程回归等方法实现弹道系数预测。Li等9利用经验模态分解方法对气动加速度的多个子项进行分解与聚合,通过先验基函数对聚类数据进行拟合实现轨迹预测。程云鹏等1
12、0对高超声速滑翔飞行器的纵向、侧向典型机动模式进行标定,利用支持向量机进行模式识别从而实现轨迹预测。张君彪等11与李明杰等12-13利用长短期记忆神经网络对再入滑翔目标控制参数进行预测与重构,代入动力学模型解算预测轨迹。Zhang等14通过经验小波变换与灰色关联度对控制参数进行去噪处理,利用注意力卷积长短期记忆网络进行轨迹预测。基于统计分析与基于运动机理的预测方法均依赖目标历史信息,一旦目标在预测过程中发生机动,预测精度便会大大降低,且预测窗口非常有限,无法满足高超声速滑翔目标防御作战需求。基于意图推断的轨迹预测方法通过挖掘攻防博弈过程中的目标意图信息,能够提高目标运动趋势预测的可靠度,提高长
13、时轨迹预测精度,被视为高超声速滑翔目标轨迹预测的有效手段。张凯等15对传统轨迹预测方法与高超声速滑翔飞行器特点进行了分析,从预测模型、预测误差、制导规律和意图信息4个角度探讨了长时轨迹预测关键技术,并阐述了目标状态与其意图推断相结合的轨迹预测思路。基于此思路,王开园等16在确定攻击某意图的前提下,通过离散化设计H-V剖面二次项系数、升阻比与终端需求速度构建机动模型集,利用跟踪阶段获取的模型概率生成预测轨迹簇。张凯等17通过意图代价函数与贝叶斯理论推断机动模式,利用蒙特卡洛采样近似飞行器状态后验概率分布实现轨迹预测。Hu等18假设高度与弹道倾角为常值简化动力学模型,快速解算可达区并离散化,基于贝
14、叶斯理论与蒙特卡洛采样实现轨迹预测,但是该文献将参数概率转化为区域概率估计问题,未考虑状态转移的动力学约束。目前基于意图推断的轨迹预测方法研究成果较少,且已有方法中,多采用蒙特卡洛方法,导致预测模型误差大,计算复杂。为进一步提高强机动目标轨迹预测精度,满足长时轨迹预测需求,本文提出了一种基于贝叶斯理论的多模型多意图融合轨迹预测方法。首先,根据典型控制规律分析,以升力参数、阻力参数、倾侧角幅值与倾侧角符号为特征参数设计纵向、侧向预测模型,综合构成时变预测模型集,模型集根据目标历史运动规律、典型动作与任务环境实时更新,提高了未知机动适应能力。然后综合攻击意图与行为意图设计意图代价函数,在此基础上,
15、通过贝叶斯理论推导机动模式与攻击意图的后验概率,并据此对多模型、多意图进行融合实现轨迹预测,提高了轨迹预测准确性与鲁棒性。1高超声速滑翔飞行器模型假设高超声速滑翔目标侧滑角为零,忽略地球168第 2 期李佳丽等:面向高超声速滑翔目标的多模型多意图融合轨迹预测自转与扁率,其在半速度坐标系下的动力学模型为:r=Vsin=Vcossinrcos=VcoscosrV=-D-gsin=1V()Lcos+()V2r-g cos=1V()Lsincos+V2rcossintan(1)式中:r为地心距;和分别为经度、纬度;V为飞行速度大小;和分别为弹道倾角与航向角;g为重力加速度大小;为倾侧角;L和D分别为升
16、力加速度与阻力加速度大小,其计算公式为:L=12mV2CL(),Ma SD=12mV2CD(),Ma S(2)式中:m为飞行器质量;为大气密度;S为飞行器特征面积;CL和CD分别为升力系数与阻力系数,是关于攻角和马赫数Ma的函数。高超声速滑翔飞行器控制量为攻角与倾侧角。攻角的主要作用是提供气动力,维持飞行器的稳定飞行;而倾侧角的主要作用是通过改变升力方向,实现飞行器机动飞行。由于非合作目标的气动系数、质量与参考面积是未知的,在滤波模型中,气动参数按下式计算:uLuD=CL(),Ma SmCD(),Ma Sm(3)式中:uL为升力参数;uD为阻力参数。二者与倾侧角共同构成高超声速滑翔目标的控制变
17、量。跟踪过程中,将3个控制变量建模为一阶马尔可夫过程:uLuD=-LuL+L-DuD+D-u+(4)式中:L,D与分别表示升力参数、阻力参数与倾侧角的机动频率;L,D与表示对应的过程噪声。综合式(1)与式(4)构成滤波模型,通过扩展卡尔曼滤波(EKF)对控制规律进行辨识,实现对目标的跟踪估计。2预测模型集高超声速滑翔飞行器一般采用纵向、侧向解耦控制。纵向需满足动压、过载、热流以及航程、终端高度等约束,通过设计攻角方案与倾侧角剖面来实现。侧向需满足规避禁飞区、指向意图以及机动动作等约束,一般通过设计倾侧角符号翻转逻辑来实现。侧向机动不消耗能量,也无法增加航程,它通过倾侧角符号翻转来满足地理约束,
18、提升突防能力。综上,以升力参数、阻力参数与倾侧角幅值为特征参数,设计纵向机动模型;以倾侧角符号为特征参数设计侧向机动模型,综合构成高超声速滑翔飞行器预测模型集,模型集根据目标飞行状态与任务环境在线更新,从而提高机动目标预测能力。2.1纵向预测模型高超声速滑翔飞行器的纵向机动特性由攻角方案与倾侧角剖面共同确定。攻角一般设计为关于最大攻角和最大升阻比攻角的分段方案,变化幅度较小。根据式(3),升力参数与阻力参数也呈现出分段特点,且变化幅度较小。倾侧角剖面中倾侧角幅值变化连续且变化幅度较小。因此,目标历史控制规律辨识是其控制规律预测的重要参考。基于目标控制规律的辨识结果,通过最小二乘法进行滑窗拟合得
19、到:ufitL,k=fL(uL,kest-nw+1,uL,kest)ufitD,k=fD(uD,kest-nw+1,uD,kest)|fitk=f(|kest-nw+1,|kest)(5)式中:ufitL,k,ufitD,k与|fitk表示纵向特征参数的拟合预测结果;fL(),fD()与f()表示对应的拟合函数;nw为滑窗容量;kest为跟踪结束时刻。由于辨识过程中存在滤波误差,在预测过程中假设升力参数、阻力参数与倾侧角幅值均满足以拟合预测结果为均值的高斯分布,即:uL,kN(ufitL,k,2L)uD,kN(ufitD,k,2D)|kN(|fitk,2|)(6)式中:2L、2D与2|分别表示
20、升力参数、阻力参数与倾侧角幅值的方差,它们根据跟踪结束时刻的滤波协方差矩阵确定,即:169宇航学报第 45 卷2L=Pkest772D=Pkest882|=Pkest99(7)式中:Pkest为kest时刻的滤波协方差矩阵,表征了各状态量的滤波估计精度;Pkestij表示Pkest的第i行、第j列元素;为放缩因子。综上,根据式(6)可生成nv组纵向特征参数,从而构成高超声速滑翔飞行器的纵向预测模型集 Mv1,Mvnv。2.2侧向预测模型高超声速滑翔目标再入过程中主要有规避禁飞区(地缘政治限制或威胁规避)等地理约束以及攻击意图等终端约束。在满足以上约束的前提下,高超声速滑翔目标通常采用C型与S型
21、2种侧向机动模式,以增强突防能力。因此针对以上4种典型动作设计侧向预测模型。(1)攻击意图模型当高超声速滑翔飞行器处于攻击意图状态时,倾侧角符号由高超声速滑翔飞行器与其意图之间的视线角T确定,其计算公式如下:T=arctan(sin(T-k)cosktanT-sinkcos(T-k)(8)式中:(k,k)表示k时刻高超声速滑翔目标的经度、纬度;(T,T)表示意图T的经度、纬度。根据视线角设计航向角走廊为:maxT=T+minT=T-(9)式中:表示航向角误差门限。根据意图航向角走廊确定倾侧角符号:signk+1=-1,1,sgn()k,k maxTk minTminT k maxT(10)式中
22、:signk+1表示k+1时刻的倾侧角符号;k表示k时刻的倾侧角;k为高超声速滑翔目标的航向角。(2)规避禁飞区模型当高超声速滑翔目标规避禁飞区时,航向角走廊根据目标当前位置、意图与禁飞区的相对关系确定。如果目标与意图在禁飞区的同一侧,即满足:k M-rMR0且T M+rMR0且T M+rMR0k M-rMR0且T M+rMR0且T M+rMR0(11)式中:(M,M)表示禁飞区中心的经度、纬度;rM表示禁飞区半径;R0表示地球半径。此时目标不受禁飞区影响,由禁飞区确定的航向角走廊为:maxM=minM=-(12)式(11)中第1种情况表示目标与意图均在禁飞区的西侧,此时目标攻击意图不需经过禁
23、飞区,因此禁飞区对目标状态无影响。同理,当目标与意图均在禁飞区的东侧、南侧与北侧时,分别对应式(11)中第2、3、4种情况,禁飞区对目标状态同样无影响。当不满足式(11)中4种情况,高超声速滑翔目标受到禁飞区影响。考虑禁飞区的航向角走廊如图1所示,图中+M与-M分别表示目标当前位置与禁飞区M的2个切线角,其计算公式如下:+M=M+arcsin()sin(rM/R0)sin(LM/R0)-M=M-arcsin()sin(rM/R0)sin(LM/R0)(13)式中:M表示目标与禁飞区M的视线角;LM表示目标与禁飞区M中心的距离,计算公式为:LM=R0arccos(sinMsink+cosMcos
24、kcos(M-k)(14)当|-M-T|+M-T|时,高超声速滑翔飞行器从禁飞区南侧绕过,此时航向角走廊为:maxM=+MminM=(16)考虑所有禁飞区M,航向角走廊为:170第 2 期李佳丽等:面向高超声速滑翔目标的多模型多意图融合轨迹预测 max=max(maxM)min=min(minM),M(17)则根据禁飞区航向角走廊设计倾侧角符号翻转逻辑:signk+1=-1,1,sgn()k,k maxk LM时,目标瞬时转弯圆与禁飞区相交,目标从禁飞区内部穿过的风险增大。此时禁飞区意图代价IM,k 0-。相反,当相对距离比小于 1,即LminM 0(LM)(28)综上,禁飞区意图代价为:I,
25、k=I*M,k(29)综合攻击意图与行为意图,定义意图代价函数:g(T,xk)Ip exp(-1Idis,k-2Iangle,k-3I,k)(30)式中:1、2与3为代价系数。考虑到目标在整个飞行过程中攻击意图与行为意图具有连贯性,将意图代价函数构建为衰减自相关递推序列,即:g(T,xk)=g()T,xk-1+g()T,xk+1(31)式中:为衰减系数。3.2多模型多意图融合轨迹预测当攻击意图T时,根据意图代价函数,其在机动模式uj下状态转移的似然概率定义为:Pr(xk|xk-1,T,uj)K-1g(T,xk)(32)式中:K为归一化系数:K=uj g()T,xk(33)根据贝叶斯定理,机动模
26、式uj的后验概率为:Pr()uj|x1:k,T=Pr()uj,x1:k,TPr()x1:k,T=Pr()xk|x1:k-1,T,ujPr()x1:k-1,T,ujPr()xk|x1:k-1,TPr()x1:k-1,T=Pr()xk|x1:k-1,T,ujPr()uj|x1:k-1,TPr()xk|x1:k-1,T(34)由于状态转移过程具有马尔可夫性,上式可进一步简化为:Pr()uj|x1:k,T=Pr()xk|xk-1,T,ujPr()uj|x1:k-1,TPr()xk|x1:k-1,T(35)式中:Pr(xk|xk-1,T,uj)为机动模式uj的似然概率,根据其定义式(32)计算得到;Pr
27、(uj|x1:k-1,T)为上一时刻机动模式uj的后验概率,通过迭代计算得到;Pr(xk|x1:k-1,T)表示意图T的似然概率,其计算公式如下:Pr()xk|x1:k-1,T=uj Pr()xk,uj|x1:k-1,T=uj Pr()xk,uj,x1:k-1,TPr()x1:k-1,T=uj Pr()xk|x1:k-1,T,ujPr()uj|x1:k-1,T(36)图2禁飞区约束示意图Fig.2Schematic diagram of no-fly zone constraints172第 2 期李佳丽等:面向高超声速滑翔目标的多模型多意图融合轨迹预测进一步,根据贝叶斯定理,意图T的后验概率
28、计算公式如下:Pr()T|x1:k=Pr()xk,x1:k-1,TPr()x1:k=Pr()xk|xk-1,TPr()T|x1:k-1Pr()xk|x1:k-1(37)根据全概率公式,利用意图T的似然概率与后验概率对高超声速滑翔目标状态进行一步预测:Pr()xk+1|x1:k=T Pr()xk+1,T|x1:k=T Pr()xk+1,T,x1:kPr()x1:k=T Pr()xk+1|xk,TPr()T|x1:k(38)同理,递推得到多步预测的计算公式:Pr()xk+j|x1:k=xk+j-1 Pr()xk+j,xk+j-1|x1:k=xk+j-1 Pr()xk+j|xk+j-1Pr()xk+
29、j-1|x1:k(39)其中,状态转移似然概率计算公式为:Pr()xk+j|xk+j-1=T Pr()xk+j,T|xk+j-1=T Pr()xk+j|xk+j-1,TPr()T|xk+j-1(40)式中:第 1 项为似然概率Pr(xk+j|xk+j-1,T),根据式(36)计算得到;第2项为意图后验概率。在预测阶段,观测信息缺失,依靠误差较大的预测信息进行意图推断难以保证准确性与可靠性,同时考虑到目标意图在飞行过程中不会频繁改变,因而高超声速滑翔目标的意图后验概率分布在预测阶段不再迭代更新,即:Pr(T|xk+j-1)=Pr(T|x1:kest)(41)注3.多模型多意图融合轨迹预测方法综合
30、了基于运动机理与基于意图推断 2类轨迹预测方法的优势,以意图代价函数作为评价规则,基于贝叶斯理论推导多模型、多意图的融合规则,提高了长时轨迹预测的准确性与鲁棒性。3.3轨迹预测流程由于目标轨迹连续且光滑,轨迹序列拟合预测短时精度较高。因此,为进一步提高数据利用率,在预测初期将高超声速滑翔目标的轨迹拟合预测结果引入预测方法进行修正,修正方案如下:rfk+j=wrk+j+(1-w)rfitk+jfk+j=wk+j+(1-w)fitk+jfk+j=wk+j+(1-w)fitk+j(42)式中:rfk+j、fk+j与fk+j表示最终轨迹预测结果;rk+j、k+j与k+j表示多模型多意图融合轨迹预测结果
31、;rfitk+j、fitk+j与fitk+j为基于滤波结果,通过最小二乘法拟合预测得到的目标位置信息;w为修正系数,设计为指数函数形式:w=1-exp(-j)(43)式中:1为一常数。轨迹拟合预测结果短时精度较高,但预测误差会急剧增加,因而w在初期为较大值,之后快速减小,以剥离拟合预测结果,避免引入拟合误差。多模型多意图融合轨迹预测框图如图3所示。轨迹预测步骤如下:1)初始化a)读取传感器参数、观测数据、禁飞区以及意图分布等数据;b)初始化滤波状态及协方差矩阵、机动模式后验概率与意图后验概率等;2)预测模型集构建a)通过EKF辨识高超声速滑翔目标控制规律;b)根据式(6)可生成nv组纵向特征参
32、数,构建纵向预测模型;c)针对2.2节中4种典型动作设计侧向特征参数,构建侧向预测模型;3)意图推断a)解算目标可达区,确定可达意图;b)基于预测模型集进行预测,并根据式(31)计算局部预测结果与滤波结果的意图代价函数;高超声速滑翔目标滤波结果预测模型1预测模型2预测模型n局部预测1局部预测2局部预测n任务环境意图代价函数多模型融合多意图融合图3多模型多意图融合轨迹预测框图Fig.3Block diagram of trajectory prediction173宇航学报第 45 卷c)根据式(32)(37)迭代更新意图后验概率;4)轨迹预测a)将局部预测的意图代价函数代入式(32)计算状态转
33、移似然概率;b)根据式(35)计算预测模型的后验概率,并据此对多模型预测轨迹进行融合;c)根据意图后验概率对多意图预测轨迹进行融合;d)根据式(42)引入轨迹序列拟合预测结果进行修正,获取目标轨迹预测结果;e)保存数据并返回2)进行下一次迭代更新。4仿真校验为验证所提预测方法的有效性与优越性,以美国CAV-H19为典型目标开展仿真实验。该飞行器参考面积为 0.483 9 m2,质量为 907 kg。目标初始条件如表1所示。存在3个意图与2个禁飞区,其参数设置如表2所示,3个意图价值相同,其中意图1为飞行器真实意图。高超声速滑翔飞行器采用预测校正与倾侧角符号翻转逻辑进行制导20-21。传感器位于
34、东经60、北纬20,测量量为距离、高低角与方位角。距离测量精度为10 m,高低角与方位角测量精度为0.01 rad,数据采样率为2 Hz。仿真中滑窗容量nw=400,放缩因子=0.8,纵向机动模型数nv=5,航向角误差门限=5,代价系数1=0.000 02,2=1,3=100,衰减系数=0.1。设定倾侧角取值范围为-90,90,利用常值倾侧角法22解算高超声速滑翔飞行器可达区,结果如图4所示,3个意图均在可达区内。为验证本文所提预测方法的有效性与优越性,对文献 7 与文献 18 的轨迹预测方法进行仿真以作对比。为便于标记,以上2种轨迹预测方法分别记为方法1与方法2,本文所提轨迹预测方法记为方法
35、3。方法2中代价系数与方法3设置相同。4.1算例1高超声速滑翔目标纵向采用准平衡滑翔模式,侧向做S型机动,同时规避禁飞区并飞向意图,其真实轨迹如图5所示。分别以1 000 s、1 300 s与1 700 s为起始点进行预测,以验证轨迹预测方法在接近禁飞区、经过禁飞区、攻击意图3种典型情况下的预测性能,图5中三角形标记了预测起始点。意图后验概率如表3所示,表中tpre表示轨迹预测起始时间。在tpre=1 000 s与tpre=1 300 s时,由于受到禁飞区影响,飞行器在距离与角度方面都更倾向于意图2,因此意图2的后验概率相对较高。根据图5(a)中所示目标的轨迹变化趋势,此时将意图2作为主要意图
36、进行轨迹预测是合理的。当飞行器绕过禁飞区继续朝其真实表1目标初始状态Table 1Initial states of the target项目取值h0/km800/()100/()-20V0/(ms-1)7 1000/()-10/()45表2意图与禁飞区参数Table 2Parameters of intents and no-fly zones项目意图 1意图 2意图 3禁飞区 1禁飞区 2经度/()9095706070纬度/()3025351225半径/km40050010080604002040608010012014016020020纬度/()经度/()意图1意图2意图3图4高超声速滑
37、翔飞行器可达区Fig.4Reachable region of the hypersonic gliding vehicle4030201010203040506070809010001020201001020304080706050403020050100纬度/()高度/km纬度/()经度/()经度/()(a)水平面轨迹(b)三维轨迹禁飞区1禁飞区2意图1意图2意图3图5高超声速滑翔飞行器真实轨迹Fig.5Real trajectory of the hypersonic gliding vehicle174第 2 期李佳丽等:面向高超声速滑翔目标的多模型多意图融合轨迹预测意图(意图1)运动
38、时,意图1的后验概率逐渐增大,在tpre=1 700 s超过80%,明显高于其他意图。从tpre=1 000 s开始预测,预测时长200 s,水平面与竖直面内轨迹预测结果如图6所示,对应的预测误差如图7所示。根据预测结果可知,方法1对控制变量的历史规律进行外推实现轨迹预测,在20 s内预测误差较小,之后随着时间推移,目标实际运动与历史运动偏差增大,预测误差显著增加。由于意图推断的修正作用,方法2与方法3在水平面内预测精度较高,相对方法1分别提高了61.39%、81.05%。在高度预测方面,方法3结合历史规律与高斯过程构建纵向机动模型进行预测,同时叠加轨迹拟合预测结果的修正作用,获得了较高的预测
39、精度。方法1 完全依赖历史运动规律,预测效果稳定性较差。由于方法2假设飞行器高度为常值,而此时飞行器处于平衡滑翔状态,高度预测误差逐渐增大。当tpre=1 300 s时,预测结果如图8图9所示。预测初期飞行器主要任务为绕开禁飞区,之后朝其真实意图飞行。根据图9(a)可得,3种预测方法在前80 s的预测误差相对较小,这是由于此阶段飞行器经过了禁飞区最近点,倾侧角符号发生偏转,因此预测误差呈现出先增大后减小的特点,但是飞行器侧向机动范围较小,因而整体预测误差尚可。在80 s后,飞行器开始朝向真实意图飞行,观察图8发现,方法1不具备意图推断能力,依赖历史信息进行预测,显然其预测模型偏离了飞行器真实运
40、动模型,导致在水平与高度上均出现了较大偏差。方法2在飞行器可达集范围内通过蒙特卡洛采样进行预测,根据图8(a)可以看出,在规避禁飞区后未能及时调整倾侧角符号,导致侧向预测误差逐渐增大。基于高度不变假设,正值此时飞行器下降缓慢,方法2在高度上预测精度较高。相比之下,方法3预测精度优势明显,尤其是侧向预测方面,体现出了较强的机动自适应能力。以tpre=1 700 s为起始时刻的预测结果如图10图11所示,仿真结果表明,本文提出的多模型多意图融合轨迹预测方法的预测精度高于方法1与方法2。以上3种情况的综合预测误差与计算时间如表4所示。方法1类似于开环控制,计算速度快,但精度无法保证。方法2与方法3通
41、过意图推断提高了预测精度,但增加了计算时间。该类方法依赖于意图推断的准确性。由于方法2基于可达集区域离散进行预测,计算量较大,耗时相对较长。504030201001 0001 0501 1001 1501 200预测误差/km时间/s(a)水平误差意图365432101 0001 0501 1001 1501 200预测误差/km时间/s(b)高度误差意图3方法1方法2方法3方法1方法2方法3图7预测误差(tpre=1 000 s)Fig.7Prediction errors(tpre=1 000 s)3530252020191815105505055606560708090100纬度/()经
42、度/()(a)水平面意图3意图1意图2意图3真实轨迹方法1方法2方法348464442403850556065高度/km经度/()(b)竖直面意图3真实轨迹方法1方法2方法3禁飞区2禁飞区1图6预测轨迹(tpre=1 000 s)Fig.6Prediction trajectory(tpre=1 000 s)表3意图后验概率Table 3Posterior probability of intents项目意图 1意图 2意图 3tpre=1 000 s11.34%88.66%0%tpre=1 300 s3.44%96.56%0%tpre=1 700 s80.98%19.02%0%175宇航学报
43、第 45 卷4.2算例2目标纵向采用跳跃滑翔模式,侧向做C型与S型机动,同时规避禁飞区并指向意图,其真实轨迹如图12所示。分别在1 200 s、1 400 s与1 800 s开始预测,对应的意图后验概率如表5所示。在1 500 s之前意图1与意图2尚未分明,因此概率几乎持平。之后飞行器朝向其真实意图飞行,意图1相对距离、航向角度优势逐渐明显,对应概率逐渐提高。从tpre=1 200 s开始向后预测 200 s,轨迹预测结果如图13所示,对应预测误差如图14所示。由353025202019181510550707560708090100纬度/()经度/()(a)水平面意图3意图1意图265意图3
44、真实轨迹方法1方法2方法35040302010010206466687072747678高度/km经度/()(b)竖直面意图3真实轨迹方法1方法2方法3禁飞区2禁飞区1图8预测轨迹(tpre=1 300 s)Fig.8Prediction trajectory(tpre=1 300 s)1015501 3201 360 1 3801 3401 3003502503002001505010001 3001 3501 4001 4501 500预测误差/km时间/s(a)水平误差意图360504020301001 3001 3501 4001 4501 500预测误差/km时间/s(b)高度误差意
45、图3方法1方法2方法3方法1方法2方法3图9预测误差(tpre=1 300 s)Fig.9Prediction errors(tpre=1 300 s)3530252026241510550848860708090100纬度/()经度/()(a)水平面意图3意图1意图28682意图3真实轨迹方法1方法2方法3403530252015818283848587868889高度/km经度/()(b)竖直面意图3真实轨迹方法1方法2方法3禁飞区2禁飞区1图10预测轨迹(tpre=1 700 s)Fig.10Prediction trajectory(tpre=1 700 s)4003003502502
46、001505010001 7001 7501 8001 8501 900预测误差/km时间/s(a)水平误差意图3141210846201 7001 7501 8001 8501 900预测误差/km时间/s(b)高度误差意图3方法1方法2方法3方法1方法2方法3图11预测误差(tpre=1 700 s)Fig.11Prediction errors(tpre=1 700 s)176第 2 期李佳丽等:面向高超声速滑翔目标的多模型多意图融合轨迹预测图可得,方法1在前20 s预测误差较小,之后随着预测模型的偏离,预测误差急剧增大。方法2与方法3预测精度远高于方法1。在高度方面,方法3有效预测了飞
47、行器的跳跃弹道特性,但是预测模型误差依然存在。而方法2高度维持不变,与高超声速滑翔目标的真实跳跃弹道先远离后接近。以 1 400 s、1 800 s 为起始时刻的预测结果如图15图16所示。根据图16(a),在攻击意图过程中,目标侧向做S型机动,方法1与方法2未能及时调整倾侧角符号,导致逐渐偏离真实轨迹,水平预测误差显著增大。方法3通过侧向机动模式集,实现了倾侧角符号的在线调整,预测误差因延迟呈现出波动特征,但误差始终较小。表4预测方法对比Table 4Comparison of prediction methodstpre/s1 0001 3001 700指标预测误差/km计算时间/s预测误
48、差/km计算时间/s预测误差/km计算时间/s方法 123.452 30.386 0123.940 50.339 160.947 20.355 3方法 29.046 322.106 392.542 021.982 7164.370 522.289 8方法 34.418 913.359 818.512 813.328 39.246 513.202 8403020101020304050607080901000102010080604020807060504030202020040纬度/()高度/km经度/()经度/()纬度/()(a)水平面轨迹(b)三维轨迹禁飞区1禁飞区2意图1意图2意图3图1
49、2高超声速滑翔飞行器真实轨迹Fig.12Real trajectory of the hypersonic gliding vehicle表5意图后验概率Table 5Posterior probability of intents项目意图 1意图 2意图 3tpre=1 200 s49.07%50.76%0.17%tpre=1 400 s50.98%48.99%0.03%tpre=1 800 s80.35%19.64%0.01%40303525202830261510560706050708090100纬度/()经度/()(a)水平面意图3意图1意图26555意图3真实轨迹方法1方法2方法3
50、504030201001020565860626468667072高度/km经度/()(b)竖直面意图3真实轨迹方法1方法2方法3禁飞区2禁飞区1真实轨迹方法1方法2方法3图13预测轨迹(tpre=1 200 s)Fig.13Prediction trajectory(tpre=1 200 s)70050060040030010020001 2001 2501 3001 3501 400预测误差/km时间/s(a)水平误差意图360504020301001 2001 2501 3001 3501 400预测误差/km时间/s(b)高度误差意图3方法1方法2方法3方法1方法2方法3图14预测误差