面对接近式威胁的在轨飞行器安全性分析方法

于大腾

(北京跟踪与通信技术研究所，北京 100094)

0 引 言

空间安全在国家安全中具有极其重要的战略地位[1-2]，随着在轨飞行器在未来起到的作用越来越重要，其可能受到的威胁也与日俱增，分析其在轨安全性并据此及时进行下一步规避决策，增强自身生存力显得十分重要。目前，在轨飞行器面临的威胁主要为非合作接近式威胁，不仅包括传统的空间碎片等，更有部分潜在的主动接近式威胁。对于这类接近式威胁，通过规避机动进行防护是可行的。而进行规避防护的第一步便是安全性分析，所得结果可对后续自主规避决策提供支撑。本文研究飞行器面对接近式威胁时的安全性分析方法展开讨论。

在安全性分析方面，诸多学者对空间碎片等传统威胁进行了非常广泛的研究[3-4]。对于空间目标接近分析问题，1984年Hoots等[5]提出了一种通过相位、高度和轨道几何的解析方法，以此来确定接近目标的最小相对距离。文献[6-10]利用轨道根数，对共焦点的开普勒轨道最小相对距离以及距离函数的临界点解析计算方法进行了研究。Alfano[11]和Negron等[12]引入椭球函数和相对距离函数，提出了一种空间目标接近数值分析算法(即A-N算法)。Klinkrad等[13]和Alarcon-Rodriguez等[14]提出了一种多步筛选运动学分析算法，首先进行高度筛选，随后进行粗略距离筛选和细致距离筛选，最后对距离一阶导数求根来确定最小相对距离对应的时刻。Luo等[15]根据椭球距离、碰撞概率、预警门限等安全评价要素对GEO轨道轨迹安全性进行了分析与规划。以上均是对接近问题进行分析，利用两目标的轨道数据，设定接近门限值(给定的距离门限或者安全椭球)和分析起止时间，计算两目标间相对距离到达距离门限或进入安全椭球范围的时间。

在空间碎片碰撞风险评估方面，目前多采用碰撞概率作为评价指标，通过地面对空间目标的观测和轨道预报计算得到评价指标，并据此评估碰撞风险。Chan[16]给出了碰撞概率的计算方法，证明如果飞行器和空间碎片误差的协方差矩阵不相关，则它们的和即是相对位置误差的协方差矩阵，故碰撞概率计算转化为相遇平面二维概率问题。白显宗等[17]推导了碰撞概率显式表达式，对空间碎片的碰撞预警、碰撞危险等问题进行了研究。Patera[18]将二维积分问题转化为一维曲线积分问题，简化了不规则形状飞行器之间的碰撞概率计算方法，并且提高了计算速度。Alfano[19]推导得到了一种级数表达式，该表达式由误差函数和指数函数进行表示，同时还得到了最小级数项数。周威萍等[20]对比分析了圆轨道和椭圆轨道的碰撞通量。殷建丰等[21]基于相对轨道要素的方法，建立了计算碰撞概率的数学模型。王晓伟[22]研究了碰撞概率算法Cube模型参数对空间碎片演化模型的影响。以上方法中，Chan和白显宗的方法是解析法，而Patera和Alfano的方法是数值法。

目前，航天领域的安全性分析多以空间碎片和失效飞行器为威胁目标进行研究，安全性指标主要为碰撞概率和相对距离。未来非合作主动接近飞行器等在轨目标的威胁日趋严重，现有方法在面对这些具有自主机动能力的飞行器时将难以进行有效判断。如何对这些潜在威胁的机动轨迹进行分析，判断其对在轨目标的威胁程度是十分必要的。本文将就这些问题开展研究，尝试给出一种安全性分析解决方案。为便于统一，在后续研究中将我方飞行器称为目标飞行器(简称目标器，Target)，具有自主交会能力的非合作接近飞行器称为追踪飞行器(简称追踪器，Pursuer)。

1 安全性分析流程

追踪器通常具有较强的机动变轨能力，一般的单次预测分析往往难以达到目的。本节将对追踪器常见的交会过程进行分析，从而为安全性模型设计相应的算法分析流程。

本节安全性分析的基本思想为化繁为简，在尽可能多的考虑重要安全要素的情况下兼顾计算效率。按照该思想对安全性分析流程进行设计，其基本架构与流程如图1所示。

图1 安全性分析流程图

首先通过情报等渠道给出追踪器当前脉冲机动能力，根据其剩余脉冲冲量计算轨道高度与轨道倾角的变化范围，即威胁范围。将我方目标器轨道参数与威胁范围进行比较，处于威胁范围之外的即可认为威胁度为零并终止后续计算，仅对处于威胁范围内的航天器进行后续安全性分析。这里需要说明的是，为了提高计算效率，本文仅取轨道高度与倾角两项作为威胁范围参考指标。

当我方目标器处于追踪器威胁范围时，进入安全性分析阶段。为了安全性分析的通用性以及合理性，本文选择脉冲能力威胁、交会时间以及相对距离三个参数作为安全性特征指标进行分析。从追踪器与目标器的相对轨道位置来分，可分为共面轨道和异面轨道两种，分别建立相应的威胁评估模型，对更新的轨道数据进行计算，得到各指标的评估威胁度。

最后，在以上威胁度计算结果的基础上，引入多指标加权综合决策模型，将各个安全性特征指标下的威胁度进行归一化，并利用专家权重系数将各项指标下的威胁度相加得到最终的多指标加权综合威胁度，即可利用设置的相应门限值进行最终决策。

下文对已知剩余脉冲冲量情况下的威胁范围快速估算方法进行介绍。将脉冲速度增量视为加速度与时间间隔的乘积，可由轨道摄动相关原理推导得到由脉冲冲量产生的轨道根数瞬时变化为[23]

(1)

式中:Δa和Δi为追踪器与目标器的轨道根数差值，n为轨道角速度，eP,fP和EP分别为偏心率、真近点角以及偏近点角。ΔvPr和ΔvPt分别为所需的径向和迹向脉冲冲量。r为当前位置地心距，u为轨道坐标系相对于赤道直角坐标系绕z轴的旋转角，ΔvPh为所需法向的脉冲冲量。

(2)

则按照式(2)计算结果即可对目标器进行快速筛选，在仅考虑半长轴与倾角改变量时，对处于追踪器威胁范围内的目标器作后续的安全评估分析。

2 安全性分析指标选择与计算模型

对于追踪器，其变轨通常包括升降轨、轨道圆化、轨道倾角调整、以及升交点赤经调整等，其特点为变轨种类多、涉及参数广且不可预测。若按照相应类别分别进行安全性分析，不仅对实际导航的处理要求过高，而且也很难提供一种可以统一处理的分析算法。

本节从能量约束以及交会安全性的角度出发，在涉及的众多计算参数中选取交会脉冲冲量、交会时间以及最小相对距离三种参数作为安全性分析指标。下文对以上三种安全性分析指标构建相应的安全评估计算模型。

2.1 脉冲冲量计算模型

对于追踪器来说，实施变轨机动对目标器进行接近，并最终完成接近是其固定的攻击模式。因此，交会所需脉冲的大小是影响目标器安全性的重要因素。这里以共面和异面两种不同情况来对交会所需脉冲进行快速估计。

当追踪器与目标器共面时，若两飞行器均为圆轨道，采用霍曼变轨进行计算。设追踪器轨道半径为rP，速度为vP，目标器轨道半径为rE，速度为vE，椭圆转移轨道近拱点和远拱点的速度分别为vTP和vTE，则有[23]

(3)

故两次转移所需脉冲冲量为

(4)

式(4)为rPrE时，式(4)右侧取负即可。总脉冲冲量为两者之和，即

Δv=Δv1+Δv2

(5)

当追踪器与目标器为共面非圆轨道时，两者的轨道形状有很多种，若根据每种特例均进行分析计算将十分复杂且计算时间成本过高。为了简化分析，提高计算速度，本节根据追踪器的交会特点，即其最终目的是要与目标器交会，故两者轨道的半长轴及偏心率会逐渐接近并最终相等。因此，这里取轨道根数中半长轴和偏心率的变化作为所需脉冲计算依据。易知，冲量与半长轴和偏心率的变化关系为[23-24]

(6)

式中:Δa和Δe为追踪器与目标器的轨道根数差值，n为轨道角速度，eP,fP和EP分别为偏心率、真近点角以及偏近点角。ΔvPr和ΔvPt分别为所需的径向和迹向脉冲冲量。

式(6)中，Δa和Δe为已知量，联立即可求得ΔvPr和ΔvPt，总脉冲冲量Δv为两者之和。需要指出的是，为寻找最优脉冲，应对fP进行简单遍历，从而提高计算效率。

当追踪器与目标器异面时，本节选用双脉冲Lambert算法对变轨接近所需脉冲进行简单估计。如图2所示，rP和vP0为追踪器某一时刻的初始位置和初始速度，vPf为转移轨道的初始速度，rE和vEf为目标器在期望交会时刻的位置和速度，vE0为转移轨道的终端速度。双脉冲Lambert问题可以描述为给定追踪器某一时刻的位置rP和速度vPf，目标器在期望交会时刻的位置rE和速度vEf，以及追踪器沿转移轨道从rP到rE的飞行时间ΔtTr，进而确定追踪器的两次脉冲冲量Δv1和Δv2。

图2 异面交会示意

显然，求解Lambert问题实质上就是求解高斯问题。而高斯问题的求解，研究者们早已推导给出一组经典方程，该方程对rP，rE，rPf和vE0之间的关系进行了诠释[25]

(7)

为求解轨道转移的最优脉冲，本节将构建双重遍历模型。首先，外层遍历模型为对目标器的轨道参数进行遍历，以一定步长计算其一个轨道周期内各时刻的位置速度。其次，内部遍历模型为对转移轨道的转移时长进行遍历，转移时长遍历范围为零到目标器一个轨道周期。轨道转移时长可代入Lambert问题进行求解得到所需脉冲总冲量，在遍历的数据中取最小值即是Lambert交会最优脉冲冲量解。计算流程如图3所示。需要指出的是，由于上节中威胁范围的筛选为了计算的便捷性设计得较为简单，因此这里计算得出的最优脉冲冲量解集中的解仍有可能超出追踪器的脉冲冲量上限。这种情况将在下一节的多指标加权综合评价中加以考虑，以使综合评价更合理。

图3 Lambert最优转移脉冲计算流程

至此，便将共面和异面情况下的脉冲冲量计算模型构建完成，实际应用中可根据不同的情况选择相应模型进行脉冲冲量计算。所得脉冲冲量可代入下节中的多指标加权综合评价模型中，为最终的安全性评估提供支撑。

2.2 交会时间计算模型

对于追踪器来说，执行接近任务时，双方很可能已经处于交战状态，时间的重要性大大提升，因此交会时间是除了脉冲冲量外另一个重要的参考要素。交会时间计算模型与脉冲计算模型是相匹配的，两者的计算结果将在下节中进行综合评价。

需要指出的是，本节所作工作是在追踪器与目标器相位角理想情况下进行的，即最有利于追踪器的情况。实际情况中若存在相位差，且追踪器不进行相位等待直接机动，则需要结合当时相位差进行相应计算，所需要脉冲消耗将大于本章中情况。

当追踪器与目标器为共面轨道时，假设采用霍曼转移，最短交会时间即转移时长为转移轨道的半个周期，则有

(8)

当追踪器与目标器为共面非圆轨道时，为了简化计算模型种类，这里将其归到异面轨道情况下，统一采用Lambert算法进行计算。在采用Lambert算法进行计算时，遍历的转移轨道时长即是交会时间，其与所得到的脉冲冲量解集是一一对应的。在第3节，将利用多指标加权综合评价模型对脉冲冲量及其对应的交会时间进行安全性评估，并选择综合性威胁较大者作为追踪器的交会策略进行安全性评估。

2.3 最小相对距离计算模型

相对距离在传统空间碎片安全性分析中常作为一项核心指标，由其来对目标器的安全性进行评估具有简单客观和表征性强等特点。无论追踪器进行何种复杂机动，其最终目的都是接近目标器，因此两者相对距离在其路径规划的终端必然呈逐渐减少的趋势，其对目标器的威胁也将越来越大。

本文利用这一特点将相对距离引入安全性分析模型，以当前输入的追踪器与目标器两者的轨道状态为时间零点，以目标器轨道周期为时长分别对追踪器和目标器进行轨道外推，计算两者的相对距离，并取一个周期内的最小值作为最小相对距离计算模型的输出。

3 多指标加权综合评价与决策模型

在第2节中对安全性指标进行选择与计算后，后续安全性综合评价需要解决各项指标权重的确定以及如何综合的问题。对多指标进行综合的方法目前有很多种，但应用最广、意义最直观的是加权和法，该方法将各个指标评价值与相应的指标权重相乘后累加求和，即可得到整个方案综合评价值，从而根据评价值的大小对方案进行排序或者决策。当然也可以采用其它非线性形式，但经验判断、理论推导和仿真计算都表明，虽然加权和法较非线性形式更为简单，但两者计算的结果非常接近，而前者具有简单以及容易理解的优势[27]。

设空间飞行器安全性分析的各评价指标最终权重向量为

σ=(σ1,σ2,…,σm)σj∈[0,1],

(9)

则采用加权和法的安全性分析综合评价与决策模型可以表示为

(10)

式中:Zi即当前轨迹xi的综合评价值，根据Zi的大小即可给出追踪器对目标器i的威胁评价。qij为轨迹i的指标j标准化值。

本节接下来重点对指标结果标准化以及指标权重的确定进行介绍。

3.1 指标标准化

在应用加权和法之前，首先需要将第2节指标计算结果进行标准化。本文采用常见的向量规范化方法对所选择的三个指标进行标准化。

对于脉冲冲量指标，其当前位置对目标器进行交会所需脉冲越小，交会越容易实现，则相应的威胁就越大。考虑极限情况，当Δv→0时，即两者已趋近于交会状态；当Δv>VPmax时，即目标器超出追踪器可达域，暂无威胁，记指标值为0。因此以追踪器当前剩余脉冲估计值为参考量进行标准化，则第i条轨迹的脉冲冲量指标为

(11)

由于前述威胁范围计算较为简单，故此处计算仍存在超出当前剩余燃料的情况，当Δv>VPmax时，计qi1=0。

对于交会时间指标，易知其当前位置对目标器进行交会所需时间越短则相应的威胁就越大，当交会时间为0时，即已经交会，则指标值为1；当交会所需时间超过当前轨道轨道周期时间时，认为暂时没有威胁，指标值为0。因此以追踪器当前轨道周期为参考量进行标准化，第i条轨迹的交会时间指标为

(12)

对于相对距离指标，考虑未来一个轨道周期内追踪器与目标器外推值的相对距离越大则相应的威胁就越小，在相对距离过大时继续进行威胁考虑将没有意义，因此引入一个相对距离参考值进行标准化，小于该参考值时进行威胁计算，相对距离为0时指标值最大为1，相对距离大于该参考值时，则认为暂无威胁，指标值为0，则第i条轨迹的相对距离指标为

(13)

通过以上计算便完成了指标结果的标准化，随后需要对各个指标的权重进行计算。

3.2 权重计算模型

权重体现了安全性分析评价指标的重要程度，一般来说最终权重是主观权重与客观权重的综合，比较常用的方法是采用加权求和的方式。设σs与σo分别为指标主观权重和客观权重的向量，则组合权重可表示如下

σ=γσs+(1-γ)σo

(14)

式中：σ为组合权重向量，γ和1-γ分别为主观权重和客观权重的权系数，γ可根据实际需要确定。

主观权重常由专家打分得到，下文主要介绍客观权重的计算方法。通常情况下，某个指标计算结果的标准差越大，表明该指标值的差异程度越大，则提供的信息量就越多，在综合评价中其起的作用越大，因此其权重也应该越大。相反，某个指标计算结果的标准差越小，表明该指标值的差异程度越小，则提供的信息量也越小，在综合评价中其所起的作用越小，其权重也应越小。

本节采用指标标准差来计算相对应的客观权重，设第i条轨迹的指标qij的标准差为δj，可由多组不同情况下的追踪器轨迹实际计算得到。则其对应的客观权重计算公式为

(15)

至此，再结合式(14)即可确定组合权重。需要指出的是，由于脉冲冲量与交会时间具有一定的相关性，因此需要先利用前几节的相关内容对标准化指标进行计算，计算脉冲冲量与交会时间各自的权重并进行指标融合，第i条轨迹的计算公式如下

Ji1=σ′1qi1+σ′2qi2

(16)

式中:σ′1为脉冲冲量的组合权重，且σ′1=γσ′s1+(1-γ)σ′o1，σ′2为交会时间的组合权重，且σ′2=γσ′s2+(1-γ)σ′o2。

随后将融合后的指标Ji1作为一个整体再与最小相对距离指标qi3一起按照前述步骤计算各自的组合权重σ1和σ2，其计算方式与前述权重计算方式相同，并由此进行融合，由式(10)可得到最终的安全性综合指标值Zi，具体公式如下

Zi=σ1Ji1+σ2qi3

(17)

式中：σ1为脉冲冲量与交会时间指标融合后的组合权重，且σ1=γσs1+(1-γ)σo1，σ2为相对距离的组合权重，且σ2=γσs2+(1-γ)σo2。

根据实际情况，可对该综合性指标设置一定的阈值，超过该阈值即认定为威胁不可避免，准备采取相应的规避机动方法进行规避。

4 仿真校验

本节用对前面建立的模型进行仿真校验，为保证仿真与实际情况相接近，仿真取美国2015年发射的通信卫星AprizeSat-7/8以及电视卫星AlSat-2A作为假想目标器，相关轨道参数由国外公开网站查询获得。同时假设追踪器初始位于与AprizeSat-7共面的停泊轨道上，其余参数随机设置以保证仿真的普适性。以上飞行器各自的初始轨道根数如表1所示，轨道历元为UTC时间2019年4月17日04时0分0秒。

表1 各飞行器初始轨道根数

以轨道历元为初始时刻，设追踪器采用三脉冲变轨策略对AprizeSat-8进行接近。为检验安全性分析模型是否具有普适性，这里取的三脉冲变轨策略为非最优策略。三次脉冲在J2000坐标系下的矢量分别为

vm1=[248.584,634.62,364.293]Tm/s

vm2=[1084.415,1.221,628.763]Tm/s

vm3=[-541.598,302.413,-348.289]Tm/s

以初始时刻为零点，第二和第三次脉冲施加时间为tm1=1030 s和tm2=2660 s。设数据观测频率为每5 min一次，总时长80 min，共16个观测数据，则该变轨轨迹与上述三个目标器轨道参数可按照上节中建立的模型计算相应的指标参数。设相对距离参考值Rref=800000 m，追踪器剩余脉冲估计矩阵为

VPmax=[3000, 2200, 2200, 2200, 2200, 1800,1800, 1800, 1800, 1800, 1800, 1800,1800, 1800, 1100, 1100]m/s

现假设专家打分得到的脉冲冲量与交会时间的主观权重系数为σs1=[0.7, 0.3]，主观权重的权系数为γ=0.75，客观权重系数可利用轨迹数据结合,上文公式计算得到，经计算可得σo1=[0.549, 0.451]，则脉冲冲量的组合权重σ′1=0.663，交会时间的组合权重σ′2=0.337。根据各自的组合权重可计算得到脉冲冲量与交会时间的融合指标，类似的，该融合指标与最小相对距离的主观权重系数为σs2=[0.45, 0.55]，主观权重的权系数不变，则上述融合指标与最小相对距离的组合权重分别为σ1=0.441和σ2=0.559。为便于比较观察，在进行指标综合评价时将各指标计算结果均乘以100。

根据以上数据可对追踪器安全性指标进行计算，脉冲和交会时间的融合指标如图4所示，最小相对距离安全指标如所图5示。

图4 脉冲冲量与交会时间的安全性融合指标

图5 最小相对距离安全性指标

由图4可知，初始时刻追踪器对AlSat-2A的脉冲冲量与交会时间融合指标最大，在第一次脉冲施加后，其对AlSat-2A的脉冲冲量与交会时间融合指标快速减小，其对AprizeSat-8的指标值逐渐增大，可见安全性模型可以判断出了脉冲意图。由图5可知，第一次脉冲使最小相对距离指标快速增大后又开始逐渐衰减，在第二次脉冲施加后该指标开始稳步增加，直至第三次脉冲后指标达到峰值。

图6为多指标加权后的综合评价结果，由图可见，AprizeSat-7全程不受任何威胁，AlSat-2A初始时刻受到的威胁与AprizeSat-8相近，而在追踪器施加脉冲后，AlSat-2A所受威胁逐渐减小，AprizeSat-8所受威胁呈不规则增大趋势。之所以增大趋势不规则，是由于本文为考虑普适性而选择的交会策略并非最优，但本文算法仍对每次脉冲意图进行了准确判断。可见该模型对交会轨迹具有一定的鲁棒性。

图6 多指标加权安全性综合评价结果

为进一步验证本文模型的合理性，下文假设追踪器在上述过程中不施加任何机动，即处于待命状态，对上述指标进行重新计算，所得结果如图7、图8和图9所示。

图7 脉冲冲量与交会时间的安全性融合指标(待命状态)

图8 最小相对距离安全性指标(待命状态)

图9 多指标加权安全性综合评价结果(待命状态)

从待命状态下的计算结果可知，多数指标均维持在一个较低的范围内上下波动，且综合评价最大值也不超过30。待命状态下虽然评价值均很低，但通过比较仍可以发现追踪器潜伏状态下对我方哪些飞行器具有威胁，并可根据综合评价值的大小对潜在的机动位置进行有效估计。

5 结 论

本文提出的分析方法可推动在轨飞行器面对各种主动接近威胁时的安全性分析研究，对未来卫星防护体系及相关技术的研究具有一定的参考意义。主要结论如下：

1)对于在轨飞行器的安全性分析问题，本文以主动接近式威胁为例，分析给出了一种安全性分析综合评价模型，拓展了目前以空间碎片等非自主接近威胁的研究目标集，使飞行器进行安全分析时具有更多选择。

2)本文所建立的模型可对在轨飞行器的安全性进行有效分析。若对综合评价值设置安全阈值，受威胁的飞行器可据此对规避机动窗口进行自主决策，提升复杂态势下的自身生存能力。

3)在本文分析模型下，追踪器施加多脉冲机动变轨对目标器进行接近，其每次脉冲均会在安全性评价值中得到体现，算法可对脉冲真实意图进行准确判断。甚至在追踪器未机动状态下，仍可发现其对我方哪些飞行器具有威胁，并对潜在机动位置进行有效估计。