基于火卫二视线测量的火星接近段自主导航算法

唐青原，王晓磊

北京控制工程研究所，北京 100190

目前以火星为主的一系列深空大天体探测中的导航任务大多依然依赖地面测控技术。随着深空探测中任务复杂度的提高，地面测控由于存在大时延、难以全天候观测、测量精度有限等问题，越来越难以满足深空探测的导航要求[1]。未来的深空探测导航关键技术之一，就是在导航需求超过地面测控能力的情况下能够实现星上自主完成较高精度的轨道确定。

深空自主导航常用的测量信息包括测距、测速及测角，其中深空测距的观测信息是脉冲星脉冲间隔时间，深空测速的观测信息是天体光谱。这两种观测方式精度较高[2]，但测量技术尚不成熟，还处于理论研究阶段。测角导航，即深空光学自主导航技术，是一项已在数个深空探测项目中得到成功应用的成熟技术。其中最早的应用出现在阿波罗系列任务中[3]，利用观测月面地标获得位置估计，这一自主导航技术适用于接近目标天体表面的飞行阶段。利用导航相机观测天体视线信息的自主导航方法，其首次利用见于1986年的旅行者2号天王星飞掠任务[4]，此后的伽利略探测器在飞掠一些木星卫星时也利用了该技术[5]，但依然是将自主导航作为地面观测的补充。近些年，小天体探测中光学自主导航已成为热门技术，2004年欧空局发射的“罗塞塔”探测器[6]，2005年美国发射的“深度撞击”探测器[7]，先后在对彗星小行星等小天体的探测活动中，验证或实现了自主导航与控制。2018年抵达小行星101955Benu的OSIRIS-Rex采样返回式探测器也采用了小天体接近段自主导航技术[8]。大行星接近段双曲线轨道由于动力学特性不同，自主导航技术目前还未在轨实现。文献[9]中提出了一种利用大视角导航相机和扩展卡尔曼滤波的双曲线轨道光学自主导航方法，对金星飞越任务的轨道进行了轨道预报研究。文献[10]中研究了利用小型测量仪器进行火星接近段自主导航的情形。然而目前大多数大天体双曲线轨道预报都需要利用中心天体视半径信息，但火星大气稀薄并可能爆发全球性风暴[11]，此时会出现不可忽视的视半径错误测量，严重影响光学导航精度。去除视半径信息的天文自主导航算法此前也多有提出，具有代表性的是文献[12]中提出的一种观测火星本身或其卫星中心视线矢量相对于参考恒星的视线矢量夹角的方法，以此为基础水手9号[13]和火星勘测号对这种方式进行了在轨验证。这种方法需要预先确定双曲线轨道渐近线方向，并且对速度信息的估计能力弱，多作为地面多普勒导航的补充。

因此，本文提出一种基于纯天体视线方向测量的火星接近段自主轨道确定算法，观测信息为一组连续测量的火卫二中心相对探测器的方向矢量，利用最小二乘法，每测量一组火卫二视线矢量得到其中某一时刻探测器完整的轨道速度、位置状态估计量。这种方法不依赖于地面先验信息，可独立完成轨道确定。又考虑到观测过程中火卫二和火星时而同时出现在同一视场下，此时可以结合火星中心视线矢量方向以及火卫二的星历，得到精度较高的飞行器位置估计量，和第一种方法联合可以期望提高导航精度。

1 系统模型

1.1 火星接近段轨道动力学模型

火星接近段轨道为从火星引力球边沿开始，经历双曲线轨道运动，直到接近到火星表面200 km左右的轨道段。这一阶段飞行器动力学模型建立在火星赤道惯性坐标系O-xiyizi下，系统轨道导航的状态向量为：

(1)

(2)

式中：μm为火星引力常数；m为飞行器质量；FS为太阳第三体引力摄动加速度；Fε为所有其他摄动加速度，主要包括火星非球形摄动和太阳光压摄动。

利用EKF对导航方法进行滤波时，该系统的敏感矩阵为：

(3)

式中：aS和aε为太阳第三体引力摄动加速度和其他摄动加速度。

1.2 视线矢量的测量

(1)观测对象的选择

火星拥有两颗形状不规则的卫星，其中火卫二Deimos是距离火星较远、质量较小的一颗，表1列出了火星两颗卫星的基本参数[4-6]。

表1 火星卫星基本参数Table 1 Physical parameters of Mars moons

可以看出火卫二的轨道半径较大，适合作为测量的基线；火卫二和火卫一都是不规则天体，但火卫二有平均直径，而火卫一的不规则程度过高，不能用直径来衡量；火卫二比火卫一体积小，是太阳系内最小的卫星，对其形心的判断更精确；且火卫二、火卫一的轨道偏心率和轨道倾角都非常小，可以近似为没有倾角的圆轨道。这些特性都使得火卫二适合作为辅助观测对象。

(2)视线矢量测量原理

导航相机对被观测天体拍照，经过图像处理可提取目标形心，结合星敏感器数据测量出该星光在星敏感器测量坐标系下的方向矢量s，再通过星敏感器安装矩阵Cs进行坐标转化，得到该星光在飞行器本体坐标系下的方向矢量sb，即本体坐标系下星光相对飞行器的视线矢量。进一步通过姿态矩阵Cib可得到惯性系下星光相对飞行器的视线矢量。

sb=Cs·s

(4)

si=Cib·sb

(5)

通过测量一个或多个天体的视线矢量，并结合探测器和天体间的几何关系以及最优估计方法，可估计出探测器在惯性系下的位置速度信息，从而实现轨道确定，具体方法将在下节详细介绍。

2 观测方法

2.1 连续观测火卫二

这里介绍的连续观测火卫二自主导航方法是在单纯角度观测的高斯初轨确定方法[17]基础上，将观测视角的位置由地面转移到火星探测器上，火星是引力中心，火卫二围绕火星转动，因只观测火卫二视线信息，故后文简称为“单目标导航方法”。探测器进入火星影响球后沿双曲线接近轨道飞行，此时探测器、火星、火卫二间的位置几何关系可抽象成图1。

图1 火卫二观测几何原理Fig.1 Geometric relation of Deimos observation

(6)

这3个矢量方程也是9个标量方程，包含12个未知数，又由于3个矢量r1、r2和r3属于同一轨道，根据角动量守恒原理有：

r2=c1r1+c3r3

(7)

根据二体运动方程，利用拉格朗日系数，可知探测器轨道某一时刻的状态矢量可由任意时刻(t2)的状态矢量(r2,v2)来表示，故有：

(8)

式中：fk和gk为拉格朗日系数在tk时刻的值，若3次观测的时间间隔足够小，fk和gk仅与r2有关。这样得到18个方程解18个未知数，方程组可解。传统测定轨所用的高斯解析解法对测量误差容忍度极低，因此本文利用一次完成最小二乘法，放弃求取以上方程组式(6)～(8)的唯一解，转而增加观测次数。当解不唯一时，求取最小二乘意义下的最优解。

最小二乘格式的系统输入/输出关系可以用下式表示：

z(k)=hT(k)θ+n(k)

(9)

式中：z(k)是待辨识系统(过程)的输出；hT(k)为观测向量；n(k)为随机噪声，均值为零；θ为待辨识的参数。则有以下准则函数：

(10)

zL=HLθ+nL

(11)

其中：

(12)

(13)

(14)

式中：na为输出向量的阶次。用nb表示输入向量的阶次，则上述线性方程组式(14)有解的条件是L≥na+nb。而L=na+nb。只有nL=[0,0,…,0]T时，θ有唯一解。但只有当L>na+nb，nL≠[0,0,…,0]T时，问题才有意义。方程组有不止一个解，可求取最优解，所以为了保证辨识精度，数据长度L可以取得充分大。

最小二乘算法分为一次完成算法和递推算法，其中一次完成算法的解为：

(15)

具体到本节利用火卫二中心视线矢量作为测量量的轨道参数预报问题，将式(8)代入式(6)中，有：

(16)

待辨识的参数为t2时刻位置矢量r2和速度矢量v2以及3个测量时刻火卫二到探测器的距离ρ1、ρ2、ρ3。将拉格朗日系数在t2作泰勒展开，截至Δt的4次方：

(17)

若测量间隔时间Δt足够短可略去2次以上的项，这样拉格朗日系数中只出现未知量r2，而更高次项会含有v2，故这一步计算存在近似处理。

将(16)式写成最小二乘格式：

(18)

式中：I3为3阶单位阵，tri为第i列为全1列，其余列为全0列的3×3阶矩阵，此时对应最小二乘各个量有：

(19)

从式(19)可以看出至少3个时刻的测量量才能解出一个θ，数据长度为3。为了求出“最优”辨识参数，需要取更大的数据长度，即更多的测量量。而增加测量量意味着HL矩阵维数增加，式(15)中的求逆运算计算量急剧增加，所以数据长度亦不能取过大。考虑到计算的时间成本，数据长度取6～7在可接受范围。又如式(17)所示，若按截至Δt的2次方求拉格朗日系数，要求Δt足够小，若Δt过大会导致系统误差过大，因此待辨识的速度位置矢量应避免选择处于一组测量值两头时刻的值。本文选择t2时刻的轨道半径作为待估计量，则有：

(20)

为了弥补数据长度的不足，可以进一步求取辨识参数的均值，具体作法可将以下两种方法结合：

1)测量参数不变，随机误差变化，反复求取同一组辨识参数，并取均值。

2)除待辨识参数所在时刻的所有参数不变，其他测量时刻每次加1，再次测量，即更新测量量，再辨识同一组参数。如此进行n步，再对辨识出的n组参数求均值。

图2 准则函数变化规律曲线Fig.2 Criterion function curve

这一变化规律与轨道预报误差的变化规律一致，故选取在准则函数“拐点”处停止观测，即可表征估计误差最小处。此后只进行轨道递推，此时距近点有足够的时间制定后续制动计划。

2.2 同时观测火星和火卫二

从火星影响球边界到轨道确定完成的过程中，火卫二和火星有时会出现在同一视场下，尤其是当探测器距离火星较远时，火卫二和火星长时间保持在同一视场内。因此这种情况下可以考虑加入火星中心相对探测器的视线矢量作为补充观测量，以期提高轨道预报的精度。

同时对火卫二和火星拍照可以同时得到方向矢量ρ和r的测量值，方便直接计算出探测器到火星的距离r，因同时观测火卫二和火星视线信息，故后文简称为“双目标导航方法”。根据火星、火卫二以及探测器间的几何关系有：

rr=Rd+ρρ

(21)

(22)

根据以上原理，观测方程可写为：

(23)

观测方程的敏感矩阵为：

(24)

3 联合导航策略

以上两种导航方法在适用范围、预报精度等方面存在差异。只观测火卫二的单目标导航方法可以得到完整的观测信息，其导航精度随着探测器沿双曲线轨道接近火星呈先提高再降低的趋势，需在轨计算停止估计的时刻；同时观测火卫二和火星的方法位置估计精度很高，没有速度的直接观测信息，速度估计精度较低。因此，将两种导航方法结合起来，根据马尔可夫最优估计理论可期得到完整的观测信息和更高精度的轨道估计结果。

具体联合导航的流程如图3所示。

图3中初始轨道参数指进入火星影响球以后双曲线轨道的初始值，这个值可由上一阶段巡航段轨道的导航数据给出，也可由地面测控给出一个误差较大的初始测定轨道。T1代表火卫二观测时刻，T2代表同时观测火卫二及火星时刻。当未达到T2或者仅能在视场内观测到火卫二时,只执行连续观测火卫二的自主导航方法;当同时观测火卫二及火星时,结合该方法预报的轨道半径模值和只观测火卫二方法预报的速度矢量及位置单位矢量，得到更高精度的三自由度轨道速度位置预报值。

图3 联合导航方法计算流程Fig.3 Flow chart of combined navigation algorithm

4 仿真校验

以常见火星接近段双曲线轨道为例，从进入火星影响球不久，距离火星约5.8×105km处开始，分别对单目标和双目标视线导航以及联合导航方法进行仿真分析。

4.1 仿真条件

仿真条件如表2所示。

表2 仿真条件Table 2 Simulation condition

需要说明，表中给出的轨道初始位置速度信息因篇幅原因保留小数点后4位，实际仿真中保留14位。

此外，单目标测量间隔时间为2 400 s，双目标测量时间间隔100 s，为节省计算量，联合导航算法中双目标测量时间间隔放宽至1 000 s。

4.2 仿真结果

图4、图5是连续单视线矢量测量导航方法的位置、速度预报误差曲线，仿真中根据准则函数变化规律大概于近点前19 h停止观测，仅进行轨道递推。

图4 单目标观测轨道位置预报误差曲线Fig.4 Altitude error estimation of Deimos-only method

可以看出，随着探测器接近火星，该方法位置预报误差最初阶段经历了大幅度震荡，甚至大于初始轨道误差，而后迅速收敛，在距近点33 h时位置估计误差大小已明显小于初始误差，最终位置预报误差逐渐收敛至20 km以内，在近点处误差稍有增大。该方法的速度预报误差逐渐收敛，除近点跳变区域外，速度预报精度可以达到0.2 m/s以内，近点处预报精度在5 m/s以内。

图6、图7分别是双视线矢量测量导航的轨道位置、速度预报误差曲线。

图6 双目标观测轨道位置预报误差曲线Fig.6 Altitude error estimation of Two-target observation method

图7 双目标观测轨道速度预报误差曲线Fig.7 Velocity error estimation of Two-target observation method

图6中呈现脉冲较大的区域是火卫二和火星在导航相机中成像较接近的位置，二者间角度很小，测量误差影响较大，可设计导航流程忽略这部分测量信息，只进行轨道递推，此处保留该区域是为说明这一特性。其余区域的位置预报精度良好，误差不到10 km，并随接近火星而逐渐收敛。

根据第2.2小节的分析，可以看出由于缺乏速度观测信息，速度估计误差先是缓慢累积，在近点位置误差累积速度急剧增大，速度预报精度不理想。

加入火星视线矢量这一补充观测量后联合导航位置预报误差曲线如图8所示。

图8 联合导航方法轨道位置速度预报误差曲线Fig.8 Altitude & velocity error estimation of combined navigation algorithm

联合后轨道位置预报精度进一步提高，即便在近点脉冲区域内，位置导航精度也达到了3 km左右；速度预报误差也在单目标导航方法的导航结果上进一步提高，近点跳变区域精度达到1 m/s以内，平缓区域误差小于0.2 m/s。

进一步对火星近点轨道要素信息进行预报，主要仿真了和后续轨道机动任务息息相关的近点轨道半径预报误差和近点轨道速度预报误差，如图9、图10所示。

图9 联合导航方法近点半径预报误差曲线Fig.9 Perigee altitude error estimation of combined navigation algorithm

图10 联合导航方法近点速度预报误差曲线Fig.10 Perigee velocity error estimation of combined navigation algorithm

近点轨道半径预报误差等价于近点轨道高度预报误差，可以看出利用初始轨道信息估计的近点半径误差接近1 000 km，不能满足任何轨道捕获任务需求。经过导航预报后，近点高度的预报精度达到了2 km以内，近点速度预报精度收敛至0.8 m/s以内，满足轨道捕获任务需求。

5 结束语

通过理论分析与仿真校验可知，文中提出的基于火卫二视线矢量测量的火星接近段自主导航方法能够通过单纯依赖视线矢量测量，在缺少速度直接测量信息的情况下达到接近段导航后续任务所需的导航精度。其中火卫二单视线测量法可以实现通过一组视线矢量测量得到探测器轨道速度信息，而双视线法可以辅助提高状态量的估计精度。联合导航的最终结果显示，近点轨道高度预报精度达到了2 km，近点轨道速度预报精度达到了0.8 m/s。

大天体接近段轨道前后速度特性变化剧烈，非线性较强，难以通过单一观测手段实现高精度的自主导航目标。本文所讨论的基于单纯视线矢量测量的自主导航方法在接近火星后因轨道动力学非线性程度增加，在近点附近预报精度有所下降。在后续的研究中，可尝试对导航相机拍摄图像进行更多地信息挖掘，利用如火星日照面积等信息进一步提高轨道预报精度。