弹道导弹自问世以来,以其射程远、威力大、精度高、机动性好和生存能力强等优越性,成为战争中极具价值的武器系统,尤其海湾战争以后,受到世界各国高度重视[1]。应运而生的导弹防御成为各国军事研究不可规避的任务。
关机后,弹道导弹所受空气阻力很小,主要依靠惯性飞行,因此其飞行轨道在理论上是完全确定的。对一组雷达量测数据作平滑优化,选取合适的初值,从而外推出弹道导弹的轨道,是预警、制导等雷达软件系统的关键技术之一。
当前使用较多的弹道外推算法,一类是基于弹道导弹的动力学方程,使用数值法一步一步进行外推[2];另一类是采用解析法[3],求出一组轨道参数,从而求解整个弹道。从理论上,在选取合适步长的情况下,数值法可以达到很高的外推精度,而轨道解析法的精度相对较低,但计算速度较数值法要快。无论是数值法还是解析法,它们的外推精度都极大地依赖于初始点的位置和速度信息的准确性。初始点的选取方法应用比较广的有两种:一种是直接滤波[4-5],取滤波的最后一点为初始点;另一种是用最小二乘法分别在X,Y,Z三个维度上对量测数据进行平滑,取平滑后的中点为初始点[6]。此外,还有一些方法是在这两种方法的基础上作一些调整和组合[7-9]。滤波方法使用了弹道导弹动力学方程,最小二乘法是对弹道作数学近似。然而,弹道导弹作为一种特殊的雷达观测目标,其运行轨道理论上是可完全确定的。但是,上述方法都没有利用到弹道的这种独特特性。
针对上述情况,本文在充分研究弹道导弹轨道特性的基础上,设计并实现了一种基于三维空间几何距离准则的弹道外推方法,对初始点的位置和速度信息进行优化。该方法基于在地球惯性坐标系中,标准弹道为位于过地球质心的一个三维平面内的椭圆这一事实,先确定椭圆平面,使所有量测点到该平面的几何距离平方和最小,然后用各量测点在平面内的投影点代替原始点迹进行最小二乘法平滑,取所有平滑点中点为初始点,进行弹道外推。应用STK软件进行弹道仿真试验,分别使用最小二乘法和基于几何距离准则的方法对量测数据作平滑,选取合适的初始点进行弹道外推,结果表明,该方法的外推精度有显著提高。
弹道导弹由于飞行时间短,可以忽略其轨道摄动,将其弹道看作是绕地球质心的二维椭圆轨道运动。在以地球质心为原点的惯性坐标系(Earth-Centered Inertial,ECI)下,弹道导弹的运动规律满足开普勒定律[10]:1)弹道目标的运行轨道是一个椭圆,地球质心与椭圆的一个焦点重合;2)在相同时间内,弹道目标与地球质心的连线扫过的面积相等;3)弹道目标运行周期的平方与其运行轨道的椭圆长半轴的三次方成正比,比值为地心引力常数的倒数。
根据开普勒定律,可以推导出弹道目标的运行轨道为位于一个过地球质心的平面内的椭圆[7-8]。根据几何知识,弹道目标的椭圆轨道可由以下两步来确定:1)计算轨道所在的平面,使得地球质心位于平面内;2)确定弹道轨道所在的三维椭球体,椭球体的两个焦点都位于第1)步确定的平面内。计算所得的平面和椭球体的交线即为弹道目标的轨道。轨道确定后,弹道目标在任一时刻的位置和速度都可以计算出来。反过来,知道弹道目标在某一时刻的位置和速度,就可以外推出整条轨道。因此,用来计算的这一时刻的位置和速度信息的准确性十分重要。
基于几何距离准则的弹道外推方法分为两步:
1) 通过雷达的一组量测值估算出弹道目标的轨道;
2) 取轨道中点为初始点,用数值法进行弹道外推。
轨道计算的方法是先确定弹道平面,再在平面内计算椭圆的方程。计算中,采用J2000地心惯性坐标系。
1.1.1 计算弹道平面
考虑t1,t2,…,tn时刻的量测数据Xn={Xi},Xi=[xi yi zi]′,i=1,2,…,n。由于弹道目标轨道所在的平面经过坐标系原点(地球质心),可设该平面的方程为a1x+a2y+a3z=0,其中a1,a2,a3为待求参数。
由于随机误差的存在,雷达量测点不会严格在一个平面内,但偏离理论弹道平面不会太远,因此需要拟合一个平面,使所有量测点到拟合平面的距离平方和最小。拟合平面即为待求的弹道平面[11]。
量测点Xi到拟合平面的距离为
(1)
对上式作变换,令
(2)
式中,有以下关系:
(3)
则式(1)化为
(4)
拟合平面方程化为
(5)
t1,t2,…,tn时刻的量测点到拟合平面的距离平方和为
(6)
求取一组合适的值,使得J取值最小,由此确定的三维平面即为弹道平面[12]。引入拉格朗日因子λ,则拉格朗日函数为
(7)
极值点处,对和λ的偏微分均为0,即
(8)
将式(7)代入,化简得
(9)
定义系数矩阵A:
(10)
式中,
式(9)中方程组的前3个方程可表示为
(11)
可以看出,式(9)的解和λ分别对应矩阵A的归一化的特征向量和特征值。矩阵A的各个元素的值由Xi(i=1,2,…,n)确定,A不仅为对角矩阵,且各元素之间有耦合,具有一定的物理意义,因此一定能求出特征值和特征向量。一般情况下,对于三维矩阵A,可以求出3组归一化的特征值和特征向量,对这3组值分别计算其中使F最小的那组值即为待求系数。
1.1.2 野值判定准则
一组量测数据中可能存在偏离弹道平面较远的野值点,判定方法是:对一组量测点,用上述方法求出拟合平面,计算所有量测点到平面的距离的标准差选取合适的常数δ,利用判别式:
|di|<δs
(12)
对每个量测点进行判别,如果上式成立,则判定该量测点为正确值,反之,则判定为野值。δ的取值根据实际情况来选取,经验表明3或者4是不错的选择。
1.1.3 理论弹道计算
量测点在轨道平面的投影为
(13)
对投影点在X,Y,Z三个方向上,分别对时间t作最小二乘法拟合。因为椭球体方程为二次方程,因此采用二项式拟合。综合3个方向的拟合结果,即为理论弹道。
综上所述,弹道目标理论轨道的算法可表示如下:
1) 由一组量测点X=Xn={Xi},其中i=1,2,…,n,拟合平面使所有量测点到平面的距离平方和最小;
2) 计算量测点到平面的距离的标准差s;
3) 选取合适的δ值,计算每个量测点Xi到拟合平面的距离di,若判别式(12)成立,则判该点为正确值,反之,则判定为野值;
4) 由所有判定为正确值的量测点组成新的量测点集Xm(m≤n),令X=Xm,如果m<n,则返回步骤1),否则进入下一步;
5) 计算X中每个点到平面的投影,记投影点集为X′;
6) 对X′,在X,Y,Z三个方向上,分别对时间t作二项式最小二乘法拟合。
研究表明,通常情况下,数值法弹道外推的精度要高于解析法。得到地心惯性坐标系下的理论弹道后,转换到地心地固坐标系,取弹道数据的中点作为弹道外推的初始点,建立弹道导弹动力学方程,采用数值法进行弹道外推。
本文采用STK软件[13-14],仿真了3条标准弹道,分别为近、中、远程导弹,发射时刻均为2018年1月1日0时0分0秒,弹道生成函数为Ballistic,发射地点为同一地点,其他参数如表1所示。
表1 STK模拟导弹参数
参数参数值弹道1弹道2弹道3发射点经度/(°)127.433127.433127.433发射点纬度/(°)39.166739.166739.1667发射点高度/m000落点经度/(°)136.452156.435172.987落点纬度/(°)33.84313.777-7.7311落点高度/m000射程/km100040007000最大高度/km36014002500
雷达站址纬度为33.843°,经度为136.452°,高度为500 m,数据率为1 Hz。雷达测量误差满足高斯分布,测量精度为:距离标准差为100 m,方位角标准差为0.12°,俯仰角标准差为0.12°。试验中,对每组理论弹道进行了100次加噪模拟计算,取100次计算的平均值为最终结果。计算中,分别采用最小二乘法(LS)和基于几何距离准则(CGD)的弹道外推方法进行外推,取中点为初始点。
在量测数据点数分别为50, 100和150的情况下,外推弹道和理论弹道在斜距、方位、俯仰三个纬度上的100次模拟的平均标准差和落点误差分别如表2~表4所示。
表2 50个测量点情况下的试验结果
参数参数值弹道1弹道2弹道3LS斜距标准差/km10.3913.2819.32CGD斜距标准差/km9.7013.1319.44LS方位标准差/(°)22.496.808.16CGD方位标准差/(°) 7.472.654.64LS俯仰标准差/(°) 5.300.820.84CGD俯仰标准差/(°)5.010.800.83LS落点误差/km49.0394.16122.67CGD落点误差/km30.2770.1294.11
表3 100个测量点情况下的试验结果
参数参数值弹道1弹道2弹道3LS斜距标准差/km2.797.4612.84CGD斜距标准差/km2.737.0911.80LS方位标准差/(°)11.272.203.78CGD方位标准差/(°) 7.990.981.44LS俯仰标准差/(°) 3.190.390.33CGD俯仰标准差/(°)3.160.390.33LS落点误差/km12.6836.7940.86CGD落点误差/km9.3230.8933.62
表4 150个测量点情况下的试验结果
参数参数值弹道1弹道2弹道3LS斜距标准差/km2.378.029.24CGD斜距标准差/km2.147.328.95LS方位标准差/(°)7.791.092.16CGD方位标准差/(°) 7.431.021.88LS俯仰标准差/(°) 2.590.340.27CGD俯仰标准差/(°)2.550.340.26LS落点误差/km6.6319.1130.28CGD落点误差/km4.5316.9529.26
从表2~表4中的数据可以看出,对于3种弹道,在量测数据点数分别为50,100和150的情况下,基于几何距离准则的弹道外推方法的外推精度相比传统的最小二乘拟合方法均有所提高,其中,用50组量测数据进行外推时,3种弹道在斜距、方位、俯仰三个纬度上的标准差平均减小2.38%,56.98%和3.03%,落点误差平均减小29.03%;用100组量测数据进行外推时,3个维度的标准差平均减小5.07%,48.82%和0.31%,落点误差平均减小20.08%;用150组量测数据进行外推时,3个维度的标准差平均减小7.19%,8.00%和1.75%,落点误差平均减小15.45%。对比发现,基于几何距离准则的弹道外推方法可以极大提高落点和方位角的预测精度,斜距和俯仰角预测精度提高不明显。这一结果与上文中的算法原理符合。因为在3次仿真试验中,导弹均朝着同一个方向发射,弹道平面基本重合,而雷达正好位于弹道平面内,量测误差表现为对弹道平面的偏离。基于几何距离准则的弹道外推方法首先拟合弹道平面,将偏离弹道平面的量测映射到平面内,因此可以极大提高方位角的外推精度,而斜距和俯仰角的误差表现为平面内的偏移,该方法对二者外推精度的提高有限。可以预计,雷达量测初始值偏离弹道平面越大,该方法对弹道外推精度的提高越明显。
另一方面,仿真试验中,无论是初始量测数据点数为50,100或150,射程为1 000 km和4 000 km的导弹(弹道1和2)外推精度提高的百分比均大于射程为7 000 km的导弹。导致这一结果的原因是,随着导弹射程的增长,导弹飞行时间随之增长,弹道也会随之逐渐偏离同一个平面;同时,现有的弹道外推的数值法都是基于理想的导弹动力学方程,当外推时间增长时误差也不可避免地增大。因此,对于远程导弹,基于几何距离准则的弹道外推方法精度提高不大,但是,对于近、中程导弹,该方法外推精度的提高还是十分明显的。
总的来说,对于近、中程导弹,基于几何距离准则的弹道外推方法可以极大提高弹道外推的精度,对于远程导弹,外推精度也有所提高,但提高幅度较小。特别是在偏离弹道平面的方向,该方法精度提高尤为显著。
本文在充分研究弹道导弹轨道特性的基础上,设计并实现了一种基于三维空间几何距离准则的弹道外推方法。该方法首先在地球惯性坐标系中,确定弹道导弹轨道椭圆所在的平面,使所有量测点到该平面的几何距离平方和最小,然后用各量测点在平面内的投影点代替原始点迹进行最小二乘法平滑,取其中点为初始点进行弹道外推。应用STK软件进行弹道模拟仿真,分别使用传统的最小二乘法和基于几何距离准则的方法进行弹道外推,100次仿真试验结果表明,无论是从雷达量测斜距、方位和俯仰标准差还是预测落点误差进行比较,该方法的外推精度都有所提高,特别是对于近、中程导弹,外推精度提高很大。特别是在偏离弹道平面的方向,该方法精度提高尤为显著。
该方法的核心在于充分利用弹道导弹的特性,从整体上确定弹道平面。如果将该方法与滤波结合起来,也许可以进一步提高弹道外推的精度。此外,该方法也可以应用到空间目标的轨道外推上。
[1] 朱锋. 弹道导弹防御计划与国际安全[M]. 上海: 上海人民出版社, 2001.
[2] 周进. 基于三步数值积分的弹道外推方法[J]. 国外电子测量技术, 2017, 36(9):13-16.
[3] 胡杰,王晟,崔健. 基于定轨方法的无动力飞行段弹道外推过滤算法[J]. 计算机工程与应用, 2017, 53(S2):62-65.
[4] 单奇,钮俊清,李川. 炮位侦校雷达的数据处理研究[J]. 雷达科学与技术, 2010, 8(2):171-176.
[5] 江宝安. 基于UKF的再入大气层弹道导弹跟踪研究[J]. 雷达科学与技术, 2008, 6(2):147-150.
JIANG Baoan. Research on Re-Entry Ballistic Missile Tracking Based on UKF[J]. Radar Science and Technology, 2008, 6(2):147-150. (in Chinese)
[6] 毛艺帆,张多林,王路. 基于轨道根数交接的弹道导弹落点快速预报[J]. 空军工程大学学报(自然科学版), 2017, 18(1):33-38.
[7] 张远,陈勇,吴昊. 弹道导弹落点预报方法研究[J]. 导弹与航天运载技术, 2014(3):5-10.
[8] 陈浩,刘秋生,杨洋. 修正弹落点预测方法[J]. 火力与指挥控制, 2015, 40(12):177-179.
[9] 陈健伟,王良明,李子杰. 基于雷达测量的炮位和落点快速预测方法研究[J]. 弹箭与制导学报, 2017, 37(2):139-142.
[10] 程昊文. 航天器轨道理论在空间目标编目管理中的应用[D]. 南京: 南京大学, 2012.
[11] 林飞,曾五一. 对自变量为随机变量的回归模型估计方法的探讨[J]. 统计研究, 1999, 16(12):22-27.
[12] MALTHOUSE E C, TAMBANE A C, MAH R S H. Nonlinear Partial Least Squares[J]. Computers & Chemical Engineering, 1997, 21(8):875-890.
[13] 沈小仑,刘先博. 基于STK的地面雷达数据处理软件测试平台的设计[J]. 电脑知识与技术, 2015, 11(1):15-18.
[14] 闵艳玲,熊智,邢丽,等. 基于STK的空天飞行器全程航迹模拟仿真研究[J]. 航空计算技术, 2016, 46(2):63-67.