中国科技核心期刊      中国指挥与控制学会会刊     军事装备类重点期刊
指挥控制与仿真

指挥控制与仿真 >

2023 , Vol. 45 >Issue 5: 73 - 83

DOI: https://doi.org/10.3969/j.issn.1673-3819.2023.05.010

多模态信息融合

基于动力特征模板的无源探测主动段弹道估计

  • 张雪松 1, 2 ,
  • 吴楠 1 ,
  • 王锋 1 ,
  • 李东泽 3
展开
  • 1 战略支援部队信息工程大学, 河南 郑州 450001
  • 2 中国人民解放军32035部队, 陕西 西安 710000
  • 3 军事科学院 国防科技创新研究院, 北京 100071

张雪松(1990—),男,硕士研究生,工程师,研究方向为空间目标信息获取与处理。

吴楠(1984—),男,博士研究生,讲师。

Copy editor: 胡前进

收稿日期: 2022-10-12

  修回日期: 2022-12-22

  网络出版日期: 2023-10-13

收起

Trajectory estimation in boost phase of passive detection based on dynamic characteristic template

  • ZHANG Xuesong 1, 2 ,
  • WU Nan 1 ,
  • WANG Feng 1 ,
  • LI Dongze 3
Expand
  • 1 PLA Strategic Support Force Information Engineering University, Zhengzhou 450001
  • 2 Unit 32035 of PLA, Xi’an 710000
  • 3 National Defense Science and Technology Innovation Research Institute, Academy of Military Science, Beijing 100071, China

Received date: 2022-10-12

  Revised date: 2022-12-22

  Online published: 2023-10-13

Fold

摘要

针对单飞行器采用无源探测手段对拦截弹主动段进行弹道估计的问题,提出一种基于动力特征模板的多模型自适应无迹卡尔曼滤波(MMAE-UKF)算法。在观测信息不完备的条件下,通过引入拦截弹动力特征模板,将飞行程序角模型设计参数与发射方位角增广为目标状态分量,基于程序角模型和动力特征模板构建非线性MMAE模型集,采用MMAE-UKF算法对拦截弹的程序角模型和设计参数进行辨识,同时对拦截弹主动段运动状态进行估计。仿真结果表明,算法能够快速准确辨识出拦截弹主动段采用的程序角模型及设计参数;与传统的增广状态滤波算法相比,运动状态估计误差显著减小,滤波稳定速度大幅度提高。

关键词： 主动段; 动力特征模板; 程序角模型; MMAE-UKF算法

本文引用格式

张雪松 , 吴楠 , 王锋 , 李东泽 . 基于动力特征模板的无源探测主动段弹道估计[J]. 指挥控制与仿真, 2023 , 45(5) : 73 -83 . DOI: 10.3969/j.issn.1673-3819.2023.05.010

Abstract

An Unscented Kalman Filter with multi model adaptive estimation algorithm (MMAE-UKF) based on dynamic characteristic template is proposed to solve the problem of trajectory estimation of interceptor in boost phase by using passive detection method for single aircraft. Under the condition of incomplete observation information, by introducing the interceptor dynamic characteristic template, the design parameters of flight program model and launch azimuth angle are expanded to target state components. Based on the flight program model and dynamic characteristic template, a nonlinear MMAE model set is constructed. The flight program model and design parameters of the interceptor are identified using the MMAE-UKF algorithm, and the motion state of the interceptor in boost phase is estimated. The simulation results show that the algorithm can quickly and accurately identify the flight program model and design parameters used in the boost phase of the interceptor. Compared with the traditional expanded state filtering algorithm, the motion state estimation error is significantly reduced, and the filtering stability speed is greatly improved.

Key words： boost phase; dynamic characteristic template; flight program model; MMAE-UKF algorithm

各类新型导弹、助推滑翔飞行器等兵器具有机动能力强、飞行空域大、打击精度高等优势,被广泛运用于实际作战,并得到国内外学者的广泛关注和研究[1-2]。各类飞行器通过搭载的传感器获取拦截器的测量信息进行灵活控制和机动调整,不仅丰富了探测手段、拓展了探测范围,还增强了突防能力,是机动突防和反拦截的发展方向[3⇓-5]。研究人员采用单飞行器搭载红外传感器对拦截弹进行探测,并在弱可观条件下完成对拦截弹主动段的弹道估计,是飞行器信息处理技术研究的关键环节之一。
单飞行器搭载红外探测器探测拦截弹主动段时,同一时刻只能获得两个角度测量值,对拦截弹定位来说属于不完备观测[6]。主动段运动规律复杂,在弱可观测条件的限制下,运动模型复杂会导致模型不可辨识,运动模型简单会带来较大辨识误差[7-8]。解决此类问题有两种思路,一是融合其他先验信息以弥补观测信息的不完备,二是构建主动段的模板信息[9]。
目前,基于模板信息是解决问题的重要方法,主要包括标准弹道模板、推力加速度模板和动力特征模板。基于标准弹道模板或推力加速度模板等导弹的先验信息来提高弹道估计精度,主要受限于模板库的准确性和完备程度,而且运算时会出现不收敛的问题,在实际应用中较为困难[10-11]。动力特征模板基于同型号导弹主动段推力加速度曲线相同的事实,是同型号导弹的本质特征,体现了导弹的运动特性。在动力特征模板已知的条件下,导弹推力随时间的函数是确定的,未知的是所采取的飞行程序角模型[12]。因此,对拦截弹的弹道估计,主要解决的问题是对不同飞行程序角模型的辨识,以及拦截弹位置、速度的准确估计。
因此,在已有文献的研究基础上,本文针对单飞行器对拦截弹主动段弹道估计问题,通过引入拦截弹动力特征模板,将飞行程序角模型设计参数与发射方位角增广为目标状态分量,基于程序角模型和动力特征模板构建非线性MMAE模型集,研究采用MMAE-UKF算法对拦截弹的程序角模型和设计参数进行辨识,同时对主动段运动状态进行估计。最后,通过仿真实验验证了算法的有效性。

1 主动段弹道估计模型

1.1 动力特征模板

对于特定型号的拦截弹,部署区域已知,推进剂质量、发动机性能、燃烧效率和助推时长相对固定,因此推力加速度有着固定的变化规律[13]。也就是说当某型号拦截弹总体参数(包括质量、推力等特性参数)已知时,推力随时间的函数是确定的,未知的是不同弹道采用的程序角φ(t)和发射方位角A0。同时,在实际中为了减少控制设计的复杂度,程序角时间节点也是固定的。图1为某型单级助推拦截弹从不同发射点飞向预测拦截点的弹道和推力加速度曲线,可见对于同种类型的导弹,即使弹道不同,其推力加速度曲线基本相同,动力特征模板正是基于此提出的,这是同种类型导弹的本质特点,具有更强的实用性。
图1 同种类型拦截弹弹道及其推力加速度比较

Fig.1 Trajectory comparison and thrust acceleration comparison of same type interceptor

1.2 主动段动力学模型

在拦截弹的助推段,总加速度a由以下几项组成,包括推力加速度aT、气动力加速度aA、引力加速度ag、离心力加速度acf和科氏力加速度aco,即
a=aT+aA+ag+acf+aco
在地固坐标系下,设拦截弹的状态矢量为 X Y Z X ˙ Y ˙ Z ˙ T,拦截弹动力学方程可表示如下:
a= a X a Y a Z= T L E C F T M B L P / m 0 0+ T L E C F T V L - C X S q / m C Y S q / m C Z S q / m+ - μ X R 3 1 + J a e R 2 1 - 5 Z 2 R 2 - μ Y R 3 1 + J a e R 2 1 - 5 Z 2 R 2 - μ Z R 3 1 + J a e R 2 1 - 5 Z 2 R 2 + 2 μ R 2 J a e R 2 Z R+ ω e 2 X ω e 2 Y 0+ 2 ω e Y ˙ - 2 ω e X ˙ 0
式中,P为推力,m为实时变化的质量;CX为阻力系数,CY为升力系数,CZ为侧力系数,S为拦截弹横截面积,q为动压;μ为地球引力常数,ae为地球的长半轴,J=1.5J2,J2为地球非球形摄动,R= X 2 + Y 2 + Z 2;ωe为地球自转角速度; T M B L为弹体坐标系到发射坐标系的转换矩阵, T M B L= T L M B - 1; T L E C F为发射坐标系到地固坐标系的转换矩阵, T L E C F= T E C F L - 1; T V L为速度坐标系到发射坐标系的转换矩阵, T V L= T L V - 1。
T L M B=R1 γR2 ψR3 φ
式中,R[ ]为初等变换矩阵,φ、ψ、γ分别为俯仰角、偏航角和滚转角,通常偏航角和滚转角为0。
T E C F L=R2 - π 2 + AR1 φR3 - π 2 - λ
式中,A为发射方位角,φ、λ为发射点纬度、经度,可根据作战区域先验信息获得。
T L V=R1 νR2 σR3 θ
式中,ν、σ、θ分别为倾侧角、航迹偏航角、速度倾角,通常倾侧角、航迹偏航角为0。
在拦截弹的总体参数已知的情况下,式(2)中各个作用力的大小是确定的,可以改变的是 T M B L、 T L E C F、 T V L3个转换矩阵,而控制这3个转换矩阵的基本变量是发射方位角A0以及程序角φ(t)。因此,将程序角φ(t)的设计参数φx和A0增广为状态矢量的分量,建立状态矢量微分方程:
d d t X Y X X ˙ Y ˙ Z ˙ φ x A 0 T= X ˙ Y ˙ Z ˙ a X a Y a Z 0 0 T
对式(6)进行离散化处理,可以得到主动段运动模型:
Xk+1=ΦkXk+Ukak
式中,Xk= X k Y k Z k X ˙ k Y ˙ k Z ˙ k φ x k A 0 k T,Φk= I 3 × 3 T I 3 × 3 O 3 × 2 O 3 × 3 I 3 × 3 O 3 × 2 O 2 × 3 O 2 × 3 I 2 × 2,Uk= T 2 I 3 × 3 / 2 T I 3 × 3 O 3 × 2,ak= a X k a Y k a Z k,其中I3×3、I2×2为单位矩阵,O3×2、O2×3和O3×3为零矩阵。

1.3 程序角模型

对于一级助推的拦截弹来说,助推段主要包括垂直上升段(t0<t<t1)、程序转弯段(t1<t<t2)、重力转弯段(t2<t<t3,姿态保持),随后进入瞄准段,在接近目标达到一定距离时进入制导段。在助推段共有4种程序角模型[14]。
1)直线型程序角
导弹在程序转弯段是直接给出直线斜率形式变化的程序角。程序角变化规律为
φ(t)= π 2 　 　 　 　 　 　 t 0 ≤ t ≤ t 1 π 2 - φ ˙ t - t 1 t 1 < t ≤ t 2 φ ( t 2 ) t 2 < t ≤ t 3
式中, φ ˙为俯仰角变化率,是程序角的设计参数。
2)攻角控制型程序角
拦截弹在程序转弯段采用攻角控制程序角变化,可分为指数型攻角和三角函数型攻角。
φ(t)= π 2 　 　 　 　 　 　 t 0 ≤ t < t 1 θ ( t ) - α ( t ) t 1 ≤ t < t 2 θ ( t ) t 2 ≤ t < t 3
式中,θ(t)为速度倾角;α(t)为攻角。
当α(t)为指数型攻角时
α(t)=-4αmz 1 - z
式中,z= e - 0.28 t - t 1,αm为攻角绝对值的最大值。
当α(t)为三角函数型攻角
α(t)=-αmsin2 f ( t )
式中,
f(t)= π t - t 1 0.34 t 2 - t + t - t 1
采用该种程序角的设计参数为αm。
3)抛物线型程序角
工程上设计出一种实用的抛物线形式的程序角,变化规律为
φ(t)= π 2 　 　 　 　 　 　 　 　 t 0 ≤ t < t 1 π 2 + π 2 - C f ( t ) t 1 ≤ t < t 2 θ ( t ) t 2 ≤ t < t 3
其中,
f(t)= t - t 1 t 2 - t 1 2-2 t - t 1 t 2 - t 1
式中,C为选择的某一常数,是程序角的设计参数。

1.4 飞行器探测模型

利用红外探测器进行探测时,观测方程为角度观测量 A , E与目标状态Xk间的关系。
Zk(A,E)=h(Xk)+wk
式中,A为观测方位角,E为观测俯仰角,Xk为每一个时刻目标在地固坐标系下的状态矢量,wk为等效测量噪声,协方差矩阵R=E w k w T k。为了简化,以下公式省略下标k。假设拦截弹在地固坐标系中的位置矢量为 X Y Z T M,飞行器在地固坐标系中的位置矢量为 X Y Z T T,则拦截弹在飞行器视场坐标系中的位置矢量为
x y z= T L q T E C F L X Y Z M - X Y Z T
式中, T L q为发射系到视场坐标系的转换矩阵,具体公式参见式(3),取飞行器的俯仰角、偏航角和滚转角进行计算; T E C F L为地固坐标系到发射系的转换矩阵,具体公式参见式(4),取飞行器的发射方位角,发射点纬度、经度进行计算。
由拦截弹位置矢量可以计算红外探测器的两个角度测量值:定义方位角A为位置矢量在xoy平面内的投影与oy轴的夹角,仰角E为位置矢量与xoy平面的夹角。则位置矢量与红外探测器角度测量值的关系为
A = a r c t a n ( y x ) E = a r c t a n ( z x 2 + y 2 )
式(16)~(17)构成了观测方程(12)的解析表达式。

2 主动段弹道估计多模型滤波器设计

MMAE算法可以对参数不确定模型的状态进行估计,具有模型匹配程度高、参数辨识准确、易于工程实现的优点。在本文中,算法的核心步骤有以下3步:一是根据拦截弹的动力特征模型和程序角设计参数构建模型集;二是基于各个模型分别对拦截弹运动状态、程序角设计参数和发射方位角进行滤波初始化;三是基于各个模型分别设计一组单元滤波器,各个滤波器并行计算,获得各自状态估计值、程序角设计参数和发射方位角。通过各个单元滤波器的测量残差计算得到相应滤波器的权值(预测概率),能反映各个模型与实际系统的一致程度,从而实现拦截弹程序角模型确定及弹道参数的辨识,同时运动状态估计值是各个并行滤波器状态估计值的加权和。

2.1 MMAE模型集构建

在本文中,假设拦截弹的程序角模型为1.3节中的4种模型之一,并根据每种程序角的设计参数,与拦截弹的动力特征模型组合,构成MMAE算法的所需模型集β= β j | j = 1,2 , 3,4,βj表示模型j。

2.2 滤波初始化方法

1)运动状态初始化
由于红外探测器仅能测得视线矢量的两个角度,无法测量距离,无法对目标形成三维立体定位,需要通过引入其他约束条件进行初始化定位。如图2所示,假设拦截弹在P处被飞行器S搭载的红外探测器首次探测,此时拦截弹距离地面很近,可认为位于地球椭球表面,视线射线与地球椭球表面的交点便是目标位置,可以利用红外探测器首次测量数据的两个角度测量值对目标进行定位。
图2 飞行器与拦截弹的空间几何关系

Fig.2 Space geometric relationship between aircraft and interceptor

在地固坐标系中,O为地心,P为拦截弹初始位置,S为飞行器位置。飞行器地心矢量为 O S → E C F,长度为b,拦截弹地心矢量 O P → E C F,长度为d,探测射线矢量为 S P → E C F,长度f。拦截弹初始位置定位方程推导如下:
在飞行器视场坐标系中,设三角形OSP中∠OSP=S,则
cos S=e S P → q·e S O → q
其中,
e S P → q= c o s E s i n A c o s E c o s A s i n E
e S O → q=- 1 b O S → q
O S → q= T L q T E C F L O S → E C F
式中,A为观测方位角,E为观测俯仰角。
根据余弦定理有
d2=f2+b2-2bfcos S
整理可得
f2+ - 2 b c o s Sf+ b 2 - d 2=0
式中,拦截弹地心矢量 O P →长度d未知,可将地球平均半径作为初值,即d=6 371 110 m。计算出P点大地坐标后,通过判断计算出的大地高度的正负值,进行迭代计算,解出f有两个数值,即视线矢量与地球椭球面有两个交点,较小的值为所求。
因此,矢量 S P →在视场坐标系中的表达式为
S P → q= f c o s E s i n A f c o s E c o s A f s i n E
转换到地固坐标系中:
S P → E C F= T L E C F T q L S P → q
则 O P →矢量在地固坐标系中的表达式为
O P → E C F= O S → E C F+ S P → E C F
由此,式(18)~式(26)便构成了拦截弹初始位置的定位方程,可以求得初始位置作为滤波初始状态估计值。根据拦截弹初始位置的定位方程以及测量数据的均值和协方差,采用无迹变换的方法可以得到滤波初始协方差矩阵。
2)弹道参数初始化
拦截弹弹道参数的初始化主要包括对程序角设计参数和发射方位角的初始化。程序角模型设计参数反映了拦截弹助推段控制策略,主要影响拦截弹的转弯速度,决定了拦截弹飞行的轨迹和射程。对于程序角设计参数,不同程序角模型有不同的选取范围,一般根据经验,直线型程序角设计参数 φ ˙初始化取值范围为[0.01,0.1],单位为rad/sec;指数型程序角、三角函数型程序角设计参数αm的初始化取值范围为[5,20],单位为度;抛物线型程序角设计参数C的初始化取值范围为[0.1,1.3]。根据上面的程序角设计参数取值范围,本文的程序角设计参数初始化估计值为取值范围的中值,并对应设置较大的初始协方差。
发射方位角主要决定拦截弹纵向射面,即拦截弹拦截目标时的射向。对于发射方位角的初始化,本文根据飞行器当前位置 X Y Z T T及飞行速度 V - T、拦截弹初始化位置 X Y Z T M及飞行平均速度 V - M,并假设两者相向飞行,计算出相遇时间t,即
t=R/ V T + V M
式中,R= ( X T - X M ) 2 + ( Y T - Y M ) 2 + ( Z T - Z M ) 2。
根据相遇时间t及飞行器当前的位置速度,可以大致估算经过t时间后的飞行器位置 X t Y t Z t T T,结合拦截弹初始化位置,采用文献[15]方法计算发射方位角作为拦截弹的初始化发射方位角,同时对应设置较大的初始协方差,以保证滤波快速收敛。

2.3 MMAE-UKF算法流程

MMAE-UKF算法主要分为3个步骤,即并行UKF滤波、模型概率概新和状态及弹道参数估计输出。
1)并行UKF滤波
由于各模型的滤波算法相同,均为UKF滤波算法,因此这里仅给出第j个模型βj的UKF滤波算法,设k时刻的状态为
X ^ k / k j= X k j Y k j Z k j X ^ k j Y ^ k j Z ^ k j φ k j A 0 k j T
式中, φ k j、 A 0 k j分别为第j个模型在k时刻的程序角设计参数估计值、发射方位角估计值,j=1,2,3,4。
具体算法流程如下:
a.状态预测
X ^ k + 1 / k j= Φ k + 1 / k j X ^ k / k j
b.协方差预测
P k + 1 / k j= Φ k + 1 / k j P k / k j Φ k + 1 / k j T+ Q k j
c.采用修正比例采样规则构造Sigma采样点χk+1/k(i)和相应权重Wi
χ k + 1 / k j ( i ) = X ^ k + 1 / k j W i = κ n + κ , i = 0 χ k + 1 / k j i = X ^ k + 1 / k j + n + κ P k + 1 / k j i W i = 1 2 n + κ , i = 1,2 , … , n χ k + 1 / k j n + i = X ^ k + 1 / k j - n + κ P k + 1 / k j i W n + i = 1 2 n + κ , i = 1,2 , … , n
d.基于观测方程的非线性传播
Z k + 1 / k j i=h χ k + 1 / k j i,i=0,1,…,2n
e.计算观测值的预估值和协方差矩阵
Z ^ k + 1 / k j = ∑ i = 0 2 n W i Z k + 1 / k j i P k + 1 / k Z Z , j = ∑ i = 0 2 n W i Z k + 1 / k j ( i ) - Z ^ k + 1 / k j [ Z k + 1 / k j ( i ) - Z ^ k + 1 / k j ] T + R k + 1 j P k + 1 / k X Z , j = ∑ i = 0 2 n W i χ k + 1 / k j i - X ^ k + 1 / k j Z k + 1 / k j i - Z ^ k + 1 / k j T
f.计算新息矢量和新息协方差矩阵
N k + 1 j = Z ^ k + 1 / k j - Z k + 1 j S k + 1 j = P k + 1 / k Z Z , j
g.状态估计与协方差更新
X ^ k + 1 | k + 1 j = X ^ k + 1 | k j + P k + 1 | k X Z , j S k + 1 j - 1 N k + 1 j P k + 1 | k + 1 j = P k + 1 | k j - P k + 1 | k X Z , j S k + 1 j - 1 P k + 1 | k X Z , j T
2) 模型概率更新
第j个模型对应的单元滤波器在k+1时刻的状态估计输出值为 X ^ k + 1 | k + 1 j,对应的测量残差为 N k + 1 j。
采用Gaussian density函数计算在k+1时刻模型j为匹配模型的似然概率为
f k + 1 j= 2 π n S k + 1 j - 1 / 2×exp - 1 2 N k + 1 j T S k + 1 j - 1 N k + 1 j
根据似然概率进行模型匹配概率的更新,各模型的预测概率为
μ k + 1 | k + 1 j= 1 C μ k + 1 | k j f k + 1 j
式中,C为归一化系数。
C= ∑ i = 1 n μ k + 1 | k j f k + 1 j
由上述公式可知,各个模型的预测概率之和为1,且由于各模型之间的程序角设计参数的差异,导致单元滤波器对模型匹配的程度不同,因而滤波得到的预测概率存在差异,即模型集中能够与实际弹道采用的程序角模型相匹配的,预测概率快速收敛至1,而不能够匹配的,预测概率收敛于0,据此可以判断实际弹道所采用的程序角模型,对应单元滤波器输出的程序角设计参数、发射方位角估计值即为所求。
3)运动状态及弹道参数估计
根据多个模型的估计输出,则估计的总输出为
X ^ k + 1 | k + 1= ∑ i = 1 r μ k + 1 | k + 1 j X ^ k + 1 | k + 1 j
总的协方差矩阵为
P ^ k + 1 | k + 1= ∑ i = 1 r μ k + 1 | k + 1 i{ P ^ k + 1 | k + 1 i+[ X ^ k + 1 | k + 1 i- X ^ k + 1 | k + 1- X ^ k + 1 | k + 1]T}
由模型概率更新结果确定实际拦截弹所采用的程序角设计模型,则对应滤波器输出的程序角设计参数、发射方位角即为弹道参数辨识结果。
算法的整体流程如图3所示。在获得首次探测数据Z0时,可以根据2.2节滤波初始化方法对k时刻各模型的运动状态、程序角设计参数φx和发射方位角A0进行初始化;后续基于各单元滤波器进行滤波计算,得到对应的状态估计值 X k + 1 | k + 1 j、程序角设计参数 φ ^ x j、发射方位角 A 0 j以及测量残差 N k + 1 j,并计算相应滤波器对实际弹道采用程序角模型的预测概率 μ k j;最后根据预测概率判断实际弹道所采用的程序角设计模型βj,输出对应单元滤波器的程序角设计参数 φ ^ x、发射方位角估计值 A ^ 0,运动状态估计值为各状态估计值的加权和。
图3 算法的整体流程

Fig.3 The overall flow of the algorithm

3 仿真实验与结果分析

为验证算法有效性,本文构建三维仿真场景,一方面验证对于不同程序角模型的辨识以及对增广参数估计的能力,另一方面验证算法对于运动状态的滤波性能,并与传统的增广状态的滤波算法进行了对比。

3.1 仿真场景

根据飞行器动力学方程和数值积分算法,优化一条飞行器滑翔弹道。拦截弹动力特征已知,飞行程序角节点为t1=2 s,t2=20 s,即t1时刻前垂直上升,t1~t2进行程序转弯。发射点地理坐标为19.6°E,0.4°N,高程为5 m,发射方位角为-117.65°,分别采用1.3节的程序角模型优化出一条飞向预测拦截点的弹道,4种程序角模型的参数取值如表1所示。
表1 程序角参数

Tab.1 The flight program parameters

程序角模型 设计参数 取值
直线型程序角 φ ^ 0.053 8 rad/sec
指数型程序角 αm 14.14°
三角函数型程序角 αm 14.14°
抛物线型程序角 C 0.66
构建三维场景如图4所示,采用不同程序角模型的拦截弹射程-高度曲线如图5所示。飞行器搭载红外探测器对拦截弹进行探测,并假设拦截弹发射后即被飞行器发现,根据飞行器和拦截弹的相对关系,计算出助推段观测数据,并加入测量噪声,生成测量数据,数据率为T=0.1 s。
图4 三维弹道曲线

Fig.4 The 3D curve of trajectory

图5 射程-高度曲线

Fig.5 Range-altitude curve

3.2 仿真环境及参数设置

实验平台如下:移动工作站,操作系统为Windows 10,CPU为Intel Core(TM) i7-9850(2.60 GHz),内存为32 GB,显卡为NVIDIA Quadro RTX 3000;编程环境为Matlab R2020b。
本文算法的滤波初始化采用2.2节方法,利用首次测量数据得到运动状态的初始估计值为
X0= 6008434.5,2139507.3,44229.417,0 , 0,0 T
初始协方差设为
P0=diag 1 e 4,1 e 4,1 e 4,1 e 2,1 e 2,1 e 2
程序角模型设计参数初值为直线型程序角0.05,指数型程序角、三角函数型程序角12.5,抛物线型程序角0.7。首次探测到拦截弹时,飞行器位置为 6155618.11780988.1 - 0.00046042 T,速度为2 752 m/s,拦截弹初始化位置为[6 008 434.5,2 139 507.3,44 229.417]T,飞行平均速度取750 m/s,经计算发射方位角初值为115.1°。上述状态分量的过程噪声方差设计采用文献[16]的方法。传统的增广状态的滤波算法采用UKF,将时变参数程序角φ(t)与发射方位角A0作为增广状态分量进行滤波,采用同样的滤波初始化方法。观测角度标准差为1×10-4 rad。

3.3 仿真结果及分析

1)程序角模型辨识
由飞行器测量数据和2.2节中的滤波算法对拦截弹采用的不同程序角模型进行辨识。直线型程序角,设计参数 φ ^=0.053 8 rad/sec,辨识结果如图6所示;指数型程序角,设计参数αm=14.14°,辨识结果如图7所示;三角函数型程序角,设计参数αm=14.14°,辨识结果如图8所示;抛物线型程序角,设计参数C=0.66,辨识结果如图9所示。
图6 直线型程序角辨识结果

Fig.6 The identified result for flight program of linear form

图7 指数型程序角辨识结果

Fig.7 The identified result for flight program of exponential form

图8 三角函数型程序角辨识结果

Fig.8 The identified result for flight program of trigonometric form

图9 抛物线型程序角辨识结果

Fig.9 The identified result for flight program of parabolic form

从上面的辨识结果可以看出,对于4种不同程序角模型的拦截弹助推段弹道,MMAE-UKF算法可以有效实现对拦截弹程序角模型的辨识。如图6a)所示,由于直线型程序角模型采用俯仰角变化率直接控制俯仰角的变化,与指数型程序角模型、三角函数型程序角模型及抛物线型程序角模型差异较大,能很快地被排除或者确认,因此在拦截弹577.4 s进入程序转弯后匹配概率便不断上升,在583 s前便很快被辨识匹配。对于指数型程序角模型、三角函数型程序角模型来说,两者属于攻角控制型程序角模型,对于程序俯仰角控制具有一定的相似性,如图7a)、图8a)所示,在辨识过程中两者在程序转弯时初始时会彼此互相影响,且同时受到抛物线型程序角模型的影响,但两者分别能在586.5 s和588 s前辨识出所匹配的程序角模型。对于抛物线型程序角来说,如图9a)所示,主要受两种攻角控制型程序角模型影响,且受三角函数型程序角模型影响较大,因此匹配概率在585 s才开始缓慢上升,约在590 s才完全匹配。
从4种程序角模型对应的拦截弹攻角来分析,可以更为显著地说明产生这种结果的原因,如图10所示。可以看出,在577.4 s进入程序转弯后,直线型程序角模型与其余模型的攻角差异明显,这也导致弹道有明显差异,因而能够较为快速地辨识匹配。同理,两种攻角控制型程序角模型虽然具有相同的最大负攻角参数,但两者到达最大负攻角的时间是不同的,因此虽然两者在辨识的起始阶段互有影响,但仍然能够明显辨识匹配。而抛物线型程序角模型的攻角在起始阶段介于指数型程序角模型与三角函数型程序角模型之间,且到达最大负攻角的时间与三角函数型程序角模型相近,因此在辨识起始时受三角函数型程序角模型影响较大,匹配概率上升缓慢。
图10 四种程序角模型对应的攻角

Fig.10 Angle of attack corresponding to four flight program

对于设计参数估计来说,如图6b)、图7b)、图8b)和图9b)所示,指数型程序角、三角函数型程序角以及抛物线型程序角辨识结果较好,三者均能在585 s前实现参数估计的收敛,且参数估计值稳定收敛于真值附近。而直线型俯仰角虽然能够辨识出设计参数,但在585 s后参数才收敛,且参数估计值不够平稳,这是因为直线型设计参数值较小,在滤波过程中受噪声影响明显。
传统增广状态滤波算法直接以程序角这一时变参数作为增广状态分量进行滤波,这里以直线型程序角模型对应的弹道为例进行分析,其余几种程序角模型的弹道类似,辨识结果如图11所示。图中黑色实线为利用本文算法辨识出的设计参数重建后得到的程序角曲线,黑色虚线为传统增广状态滤波算法的辨识结果,红色虚线为程序角的真实值。可以看出,与真实程序角曲线相比,传统增广状态滤波算法的辨识结果精度较低,而本文算法采用对设计参数这一固定值进行辨识,能得到较为准确的设计参数估计值,重建后能得到比较准确的程序角变化曲线。
图11 程序角辨识结果对比

Fig.11 The identified results comparison of flight program angle

2)发射方位角辨识
由于拦截弹的发射点和预测拦截点都是固定的,采用4种程序角模型的弹道,其发射方位角是相同的。这里仅对采用直线型程序角模型的弹道发射方位角进行辨识分析,其余几种程序角模型的弹道类似。如图12所示,本文算法和传统增广状态滤波算法对发射方位角的辨识估计结果总体上比较准确,但本文算法收敛速度更快,辨识更为准确。本文算法在583 s时滤波估计结果开始平稳且收敛于真值,但在初始阶段发射方位角估计值产生剧烈震荡,这是因为拦截弹在垂直上升段和程序转弯起始的一个阶段射向不明显,因此辨识估计结果较差。
图12 发射方位角辨识结果

Fig.12 The identified results for launching azimuth

3)运动状态估计
对于运动状态估计,这里以直线型程序角模型对应的弹道进行分析,其余几种程序角模型的弹道类似。采用均方根误差RMSE表示跟踪过程中的误差变化。定义如下:
RMSE= 1 N ∑ i = 1 N x i k - X ^ i k | k 2
式中,xi k和xi k | k分别表示第i次仿真中第k时刻真值和估计值,N为蒙特卡洛仿真次数。
以x轴方向状态估计均方根误差为例,其他方向基本相同,进行500次蒙特卡罗仿真,结果如图13和图14所示。实线描述本文算法的弹道估计结果,虚线描述基于传统增广状态滤波算法的弹道估计结果。
图13 X方向位置均方根误差

Fig.13 The RMSE of position in X direction

图14 X方向速度均方根误差

Fig.14 The RMSE of velocity in X direction

从图13和图14可以看出,与传统增广状态滤波算法相比,本文算法估计精度更高。本文算法滤波稳定时的位置均方根误差约为5 m,而传统增广状态滤波算法的位置均方根误差约为10 m,本文算法的位置均方根误差比传统增广状态滤波算法减小50%;本文算法速度估计的均方根误差约为1.2 m/s,而传统增广状态滤波算法的速度均方根误差约为1.8 m/s,本文算法速度估计的均方根误差约比传统增广状态滤波算法减小30%。此外,本文算法在滤波时能更快到达稳定状态。

4 结束语

本文以单飞行器搭载红外探测器对拦截弹主动段探测跟踪为背景,基于动力特征模板和程序角模型构建非线性MMAE模型集,提出一种可以对拦截弹主动段程序角模型及设计参数进行辨识,同时对运动状态进行估计的MMAE-UKF算法。仿真结果表明,本文算法能够快速准确对拦截弹主动段程序角模型及设计参数进行辨识。同时,与传统的增广状态的滤波算法相比,运动状态估计误差显著减小,滤波稳定速度大幅度提高。后续,可以在此基础上进行弹道外推,对拦截弹的预测拦截点进行初步估计,评估拦截弹中末制导所需修正误差的量级大小,对飞行器规避拦截具有理论研究意义。

参考文献

原文顺序 | 文献年度倒序 | 文中引用次数倒序
[1]
孟二龙, 高桂清, 吴鹏程, 等. 智能导弹战发展趋势研究[J]. 飞航导弹, 2021(9): 26-31.

MENG E L, GAO G Q, WU P C, et al. Research on development trend of intelligent missile warfare[J]. Aerodynamic Missile Journal, 2021(9): 26-31.

[2]
程进, 齐航, 袁健全, 等. 关于导弹武器智能化发展的思考[J]. 航空兵器, 2019, 26(1): 20-24.

CHENG J, QI H, YUAN J Q, et al. Discussion on the development of intelligent missile technology[J]. Aero Weaponry, 2019, 26(1): 20-24.

[3]
宗凯彬, 张承龙, 卓志敏. 智能导弹武器系统发展综述[J]. 现代防御技术, 2020, 48(3): 49-55, 85.

DOI

ZONG K B, ZHANG C L, ZHUO Z M. Survey of the development of intelligent missile weapon systems[J]. Modern Defence Technology, 2020, 48(3): 49-55, 85.

[4]
雷虎民, 骆长鑫, 周池军, 等. 临近空间防御作战拦截弹制导与控制关键技术综述[J]. 航空兵器, 2021, 28(2): 1-10.

LEI H M, LUO C X, ZHOU C J, et al. Summary of key technologies of interceptor guidance and control in near space defense operations[J]. Aero Weaponry, 2021, 28(2): 1-10.

[5]
REIM G. Hypersonic arms race accelerates[J]. Flight International, 2020, 197(5717): 31.

[6]
于大腾, 王华, 尤岳, 等. 基于单星观测的弹道导弹参数估计方法综述[J]. 现代防御技术, 2014, 42(4): 62-68.

YU D T, WANG H, YOU Y, et al. Review for the method of single satellite early warning ballistic missile parameter estimation[J]. Modern Defence Technology, 2014, 42(4): 62-68.

[7]
WU N, CHEN L, LEI Y J, et al. Adaptive estimation algorithm of boost-phase trajectory using binary asynchronous observation[J]. Proceedings of the Institution of Mechanical Engineers, Part G: Journal of Aerospace Engineering, 2016, 230(14): 2661-2672.

DOI

[8]
杨锐. 单星观测下弹道导弹主动段参数估计技术研究[D]. 南京: 南京航空航天大学, 2020.

YANG R. Research on parameter estimation technology of ballistic missile boost phase under single star observation[D]. Nanjing: Nanjing University of Aeronautics and Astronautics, 2020.

[9]
王栋. 基于先验信息融合的单星预警算法及误差分析[D]. 长沙: 国防科学技术大学, 2015.

WANG D. Single satellite early warning algorithm based on prior information fusion and error analysis[D]. Changsha: National University of Defense Technology, 2015.

[10]
李晓宇, 田康生, 郑玉军, 等. 单星观测下弹道导弹状态估计与预测误差分析[J]. 火力与指挥控制, 2015, 40(8): 5-8, 13.

LI X Y, TIAN K S, ZHENG Y J, et al. Error analysis of single satellite observe ballistic missile state estimation and forecast[J]. Fire Control & Command Control, 2015, 40(8): 5-8, 13.

[11]
储雪峰, 吴楠, 王锋, 等. 基于不完备推力加速度模板的主动段弹道估计[J]. 信息工程大学学报, 2019, 20(4): 507-512.

CHU X F, WU N, WANG F, et al. Trajectory estimation in boost phase based on incomplete thrust acceleration profile[J]. Journal of Information Engineering University, 2019, 20(4): 507-512.

[12]
申镇, 强胜, 易东云. 单星探测弹道特征信息提取方法研究[J]. 系统仿真学报, 2010, 22(2): 295-297, 356.

SHEN Z, QIANG S, YI D Y. Research on estimating missile trajectory characteristic based on single satellite[J]. Journal of System Simulation, 2010, 22(2): 295-297, 356.

[13]
强胜, 申镇, 易东云. 基于导弹动力特征的单星预警算法[J]. 系统工程与电子技术, 2011, 33(10): 2234-2238.

QIANG S, SHEN Z, YI D Y. Method of early warning for a single satellite based on missile dynamic characteristics[J]. Systems Engineering and Electronics, 2011, 33(10): 2234-2238.

[14]
杨希祥, 江振宇, 张为华. 小型运载火箭大气层飞行段飞行程序设计研究[J]. 飞行力学, 2010, 28(4): 68-72.

YANG X X, JIANG Z Y, ZHANG W H. Flight program design for small launch vehicle in atmosphere flight stage[J]. Flight Dynamics, 2010, 28(4): 68-72.

[15]
王宝和. 运载火箭射击方位角计算方法研究[J]. 舰船电子工程, 2019, 39(2): 95-98.

WANG B H. Research on the calculation method of launch vehicle launching azimuth[J]. Ship Electronic Engineering, 2019, 39(2): 95-98.

[16]
丁力全, 吴楠, 孟凡坤, 等. 基于无迹变换的协同探测与速度方向控制制导律一体化设计[J]. 弹道学报, 2022, 34(2): 25-32.

DOI

DING L Q, WU N, MENG F K, et al. Integrated design of collaborative detection and velocity direction control guidance law based on unscented transformation[J]. Journal of Ballistics, 2022, 34(2): 25-32.

DOI

Options
文章导航

/

〈 〉