考虑姿轨耦合的航天器高精度实时导航定轨方法

杨盛庆，陈 桦，徐炜莉，刘美师，王 禹，钟 超

（1.上海航天控制技术研究所，上海 201109；2.上海市空间智能控制技术重点实验室，上海 201109；3.上海市宇航学会，上海 200235）

0 引言

随着体系的日趋完善，全球导航卫星系统（GNSS）测量系统已成为近地轨道上航天器轨道确定的主要手段［1-2］。基于GNSS 测量信息能够实现高精度的实时导航与定轨，是实现航天器自主智能轨道控制的基础［3-4］。传统的轨道动力学，通常将航天器视为质点进行分析和轨道积分仿真。实际上，航天器在进行姿轨耦合的刚体运动。作为自主轨控输入参数的实时导航与定轨领域，缺少姿轨耦合方面的研究工作。GNSS 接收机的原始测量数据基于瞬时的运动状态，测量点是GNSS 天线的相位中心［5］。早期，由于GNSS 测量精度在10 m 到100 m量级，GNSS 天线安装位置与航天器质心偏差造成的定轨误差基本可以忽略［6］。与之相反，相对运动的研究内容为编队坐标系下的相对轨道特征参数。相对运动过程中姿态变化对相对测量信息造成的影响较为明显，姿轨耦合运动造成的影响不容忽视［7-8］。此外，相较绝对轨道的参数确定，相对运动可以忽略轨道摄动环境中的长期性和周期性影响。因此，航天器姿轨耦合问题的研究工作，起初主要针对航天器的交汇对接、在轨服务等任务［9-10］。当GNSS 测量达到米级精度，导航与定轨的精度要求也随之提高，单点几何定轨的算法中GNSS 天线的安装关系不能再被忽略。

目前，高精度的绝对轨道定轨工作主要集中在事后精密定轨领域［11-12］。经精密定轨处理的轨道数据在地心轨道系径切法方向已达到厘米级精度，可以作为轨道控制效果评估的手段［13-15］。但是，由于事后精密定轨历史数据充分，数据批处理过程中位速参数的确定可以采用数据拟合插值，不仅仅依赖于测量信息［16］。但是，姿轨控系统的星载导航与定轨算法对实时性要求较高。算法一般采用基于测量信息的贯序数据处理方法，在构造滤波算法时需要针对GNSS 天线的安装关系进行精准的位速修正［17］。研究表明，GNSS 天线相位中心不同于天线的几何参考点（Antenna Reference Point，ARP），需要结合相位中心偏差（Phase Center Offset，PCO）和天线相位中心变化（Phase Center Variations，PCV）进行优化，能够使事后精密定轨得到毫米级的精度提升［18］。因为GNSS 天线相位中心的变化（毫米级）相对于天线的安装关系（米级）是一个小量［13］，本文主要针对天线的几何参考点与航天器质心的安装关系开展研究。

航天器作为刚体，不同质点所受的引力偏差造成的影响表现为重力梯度力矩［19-20］。分析刚体的自由运动情况较为复杂，考虑到航天器在轨的姿态指向与轨道系的关系相对明确，可以围绕其姿态特征提炼出针对安装关系的位速修正项，进而分析由安装关系造成的定轨误差的解析解及其周期变化规律。航天器星载软件的实时导航与定轨，通常采用基于测量信息转换为惯性系状态变量后的扩展卡尔曼滤波算法。滤波算法中的状态估计方程以地球引力中心项为主要因素，滤波过程能够有效修正状态量估计信息和GNSS 测量信息的白噪声。对于在轨的稳态运动，姿态抖动和GNSS 天线安装关系导致的位速偏差可以视为状态量信息与测量信息之间的噪声。但是，针对SAR 卫星载荷工作的多普勒频移原理，往往需要进行偏航为主的姿态导引补偿［21-22］。此刻，GNSS 天线与航天器质心之间的差异对导航与定轨的影响明显。本文结合轨道动力学原理，提出了考虑安装关系的高精度GNSS 测量信息模拟方法，围绕具有姿态导引要求的航天器实时导航与定轨算法，数学仿真证明了GNSS 天线安装关系补偿的必要性。

1 轨道动力学原理

1.1 二体问题的积分常数

已知惯性系位置r、速度v，二体相对运动的微分方程［23-24］为

矢量叉积满足

式中：r×v为积分常量，定义为轨道的动量矩为［1］

使用点积式（1），可得

作为积分常数，可以推导得到活力公式为

式中：a为轨道半长轴。

利用椭圆方程半通径的定义p=a(1−e2)及，轨道半长轴满足

1.2 惯性系下的刚体运动

考虑空间中的刚体运动，刚体中任意质点可以描述为惯性空间oIxIyIzI的平移运动和旋转运动，如图1 所示，惯性空间中质点mi的位置为

图1 空间刚体的运动学Fig.1 Kinematics of rigid body in space

式中：A为本体系相对于惯性空间的旋转矩阵；A⋅Δrb为质点mi关于质心的修正项。安装关系Δrb不变的情况下，惯性空间中质点mi的速度为

因此，轨道半长轴的瞬态偏差为

1.3 轨道运动特征

航天器作为刚体，在轨道运动过程中姿态特征的一种为保证星体的稳态对地指向。理想情况是航天器的本体系与轨道系重合，即零姿态运动［25］。在理想情况的假设下，分析无摄轨道运动在惯性空间中的闭合曲线运动。如图2 所示，本体系在一个轨道周期内的指向变化相对于绕其yb旋转1 周。

图2 瞬时轨道运动的坐标系Fig.2 Coordinate systems describing instantaneous orbital motion

不妨假设纬度幅角u=0 处的本体系到惯性系转换矩阵为A0，不同纬度幅角处A(u)=Ry(u) ⋅A0。导数形式，其中，n为轨道角速度。

在零姿态运动假设的基础之上，考虑本体系三轴独立偏离质心±2 m 的质点，进行速度项补偿与不补偿的数值仿真。仿真分析如图3 和图4 所示，可以发现在进行式（8）所示速度项修正的基础之上，本体系X轴（对应轨道系航向）安装与质心的位置偏值对半长轴确定的影响轨道周期明显，修正后半长轴确定精度改良明显。在本体系Y轴（对应轨道系法向）安装与质心的位置偏值对半长轴确定的影响较小。在本体系Z轴（对应轨道系径向指向地心）安装与质心的位置偏值对半长轴确定存在相同极性的偏值现象，可以理解为在零姿态运动情况下质点扫过的行程不同。此外，不同轴偏离质心量的极性会影响波动周期的相位特性。

图3 不进行速度项修正的半长轴确定偏差Fig.3 Semi⁃axis determing errors without velocity correction

图4 进行速度项修正的半长轴确定偏差Fig.4 Semi⁃axis determining errors with velocity correction

2 GNSS 天线安装与航天器质心

2.1 GNSS 天线的安装关系

在半长轴确定偏差的力学解析解基础之上，进一步分析姿态运动对半长轴确定的影响。已知GNSS 天线几何中心相对航天器质心的安装关系为Δrb=[Δx，Δy，Δz]T，使用刻画轨道运动的坐标转换惯例，将航天器本体系到惯性系的转换分为轨道系到惯性系、本体系到轨道系两个步骤，物理含义清晰。如图5 所示，假设GNSS 天线几何中心相对航天器质心的安装关系在惯性系的表述为

图5 GNSS 天线与航天器的安装关系Fig.5 Installation relationship of GNSS antenna and spacecraft

已知惯性系航天器的质心矢量rsat，GNSS 天线几何中心的位置参数在惯性系下的描述为

记航天器质心的速度为vsat，vsat=ṙsat，GNSS 天线几何中心的速度参数在惯性系下的描述为［26］

2.2 轨道系与惯性系的转换

区别于理想情况的封闭轨道假设，实际情况轨道受摄运动情况复杂，需要根据实时确定的轨道参数确定轨道系到惯性系的坐标转换矩阵Aoi，

2.3 本体系与轨道系的转换

区别于理想情况的零姿态假设，实际情况需要考虑航天器在主动姿态运动控制情况下，由实时测量的姿态角确定本体系到轨道系的坐标转换矩阵：

其导数项满足：

3 导航算法的设计与验证

3.1 GNSS 接收机测量信息的模拟

传统的轨道仿真中，普遍采用轨道动力学积分得到的位置、速度叠加误差噪声作为GNSS 接收机测量数据的模拟量。本文提出了考虑姿轨耦合的GNSS接收机测量信息模拟方法，要求GNSS 接收机测量数据的模拟考虑航天器的刚体运动特征、GNSS 天线的几何中心rGNSS、vGNSS与航天器质心rsat、vsat差别。

3.2 测量信息与导航定轨的滤波算法

航天器的姿轨控系统一般会对GNSS 接收机的测量数据进行导航滤波与轨道确定。扩展卡尔曼滤波（EKF）导航滤波算法流程可分为状态估计、测量修正两大步骤。滤波算法的状态变量是惯性系的位置、速度X=[x，y，z，vx，vy，vz]T，滤波过程中使用轨道的二阶微分方程进行状态估计，其单拍积分的有限差分形式为

式中：加速度项

式中：μ为地球引力常数；aJ2、aJ3、aJ4分别为地球引力场J2、J3和J4项摄动加速度矢量；axJ2、ayJ2、azJ2为aJ2的三轴分量；axJ3、ayJ3、azJ3为aJ3的三轴分量；axJ4、ayJ4、azJ4为aJ4的三轴分量。

3.3 导航算法的天线安装修正

在模拟刚体运动特征的GNSS 接收机测量信息的基础之上，星载软件的导航算法同样需要结合姿态测量信息，将GNSS 天线测量的惯性系位置、速度转换到航天器质心。

以航天器质心的惯性系位置、速度作为状态变量和滤波观测量，进行EKF 滤波。滤波输出的航天器质心惯性系位置、速度与轨道动力学积分理论值的偏差，即为导航滤波的精度。

在轨导航系统已知rGNSS、vGNSS，可得

3.4 姿轨耦合导航滤波的算法验证

在假设航天器质心满足轨道运动特征的基础之上，天线几何中心相对航天器质心的惯性系位置、速度需要补偿修正项。在上述两类质点的惯性系位置、速度的基础之上，分析轨道半长轴的确定偏差。

同时，结合导航滤波算法，分析不同质点惯性系信息接入滤波时的定轨差别。为了验证高精度姿轨耦合的导航算法，设计如图6 所示的导航精度验证系统。

图6 导航精度分析的逻辑说明Fig.6 Logical illumination of the navigation precision analysis

仿真结果如图7～图9 所示。仿真结果显示，当航天器不进行姿态导引时，使用GNSS 天线测量的惯性系位速信息、补偿安装关系后的航天器质心惯性系位速信息作为EKF 的滤波输入，定轨偏差较小。当航天器进行姿态导引时，定轨偏差出现与导引规律相应的周期特性。

图7 非导引状态下两种测量输入的轨道平根Fig.7 Mean orbital elements without steering

图8 导引状态下两种测量输入的轨道平根Fig.8 Mean orbital elements with steering

图9 仿真过程中导引状态下半长轴平根变化情况Fig.9 Variations of the mean semi-major axis under the steering state during the simulation process

根据前文分析，本体系X轴安装关系导致的周期特性明显，设置Δrb=[±0.5，0.2，−1.5]T，单位m。如图8（c）和图8（d）所示，两者安装关系会导致周期性规律的极性偏差。

4 结束语

本文针对航天器自主轨控使用的轨道参数，围绕具有偏航导引姿态控制特征的合成孔径雷达（SAR）卫星存在的导引与定轨数据的特征关联问题，考虑采用不同质点测量信息进行导航滤波对定轨精度的影响，并进行了仿真验证。本文的研究基于轨道动力学原理，将航天器视为刚体运动，分析了GNSS 天线几何中心与航天器质心在惯性系中的相对运动关系，明确了实时定轨误差与惯性系位置、速度的关系。姿轨控制系统在使用GNSS 接收机测量数据时，需求辨明测量数据属于几何法还是经滤波处理的动力学版本。后续研究将围绕在轨大角度机动状态下的实时定轨精度与导航滤波参数优化设计。结果表明，适当地增大误差协方差矩阵主对角线元素的速度相关项（认为速度测量置信度低），能有效减小机动对速度项造成的影响，从而保证实时导航轨道平半长轴的稳定性。