GEO与IGSO卫星轨道数值积分可行性分析

宋传峰，秘金钟，党亚民，薛树强

（中国测绘科学研究院，北京 100830）

1 引言

随着GEO （Geostationary Earth orbit）与IGSO （Inclined Geo Synchronous orbit）卫星在导航定位领域中的应用，对其轨道的精度要求越来越高［1］。利用数值算法积分卫星轨道是人造卫星轨道预报与精密定轨中的基础工作［2］。由于二体问题可以直接得出卫星位置、速度和加速度的解析解，因此本文避开复杂的摄动力计算，在二体意义下分别对GEO和IGSO卫星进行了轨道数值积分，并与其理论值进行了分析比较，得出了有益结论。若考虑摄动力，则不便于结果比较。

目前，在MEO定轨中使用的积分方法有RKF与Admas预测－校正系统结合的方法 （如Gamit软件），Collocation方法 （如Bernese软件）等。本文验证了RKF与Admas预测－校正系统结合的方法在GEO与IGSO卫星轨道积分时的适用性。

2 数值积分方法

2.1 变步长RKF7（8）阶单步积分法

Runge－Kutta－fehlberg方法 （RKF）是一种嵌套的RK方法，即同时给出n阶和n＋1阶两组RK计算公式，用两组公式计算出来的xi＋1之差来估计截断误差，根据截断误差的大小来控制步长［3］。以下给出7阶和8阶的RKF公式

截断误差

式中，h为积分步长，ck，βkλ为积分算法的系数。

2.2 定步长Adams线性多步积分法

在卫星运动方程的数值积分中，单步法只是用于多步法的起步算法，当采用单步法推出足够的步点后，就可采用高精度的多步法往前推算［3］。本文使用变步长的RKF7（8）阶积分法积分10步后，再使用定步长的10阶Adams预测－校正系统进行积分。以下给出显式的Adams－Bashforth公式和隐式的Adams－Moulton公式。

显式

隐式

3 GEO与IGSO卫星二体问题

3.1 GEO卫星二体问题

在二体意义下，GEO地球同步静止轨道卫星在地心地固坐标系中的坐标为一固定值，但轨道积分是在惯性坐标系下进行的，因此需要计算二体问题下地固坐标系与惯性坐标系间的相互转换矩阵。

假设积分初始时刻t0（s），两坐标系坐标轴完全重合，则任意t（s）时刻，地固系向惯性系转换的矩阵为

在惯性系下，GEO卫星绕Z轴旋转的角速度与地球自转角速度相等，因此易得下式

惯性系下GEO卫星在X－Y平面内运动，即Z≡0，由此可得：

式中，X、Y、Z表示卫星在惯性系下的坐标。

在二体意义下，GEO卫星的轨迹为圆，且其在地固坐标系下的速度为0，因此给定初始坐标可得任意时刻GEO卫星在地固和惯性坐标系下的位置、速度和加速度。

3.2 IGSO卫星二体问题

IGSO卫星轨道是特殊的地球同步圆轨道，轨道高度跟GEO一样，均属于地球同步轨道［4］。由于IGSO卫星轨道在惯性系下也是圆轨道，因此可看做由GEO卫星轨道绕X或Y轴或X－Y平面内一条直线旋转i角，i为IGSO卫星轨道倾角，本文中按绕X轴旋转i角计算。

若t0时刻，IGSO卫星在其轨道直角坐标系（惯性系绕X轴旋转i角）下坐标为（，，0）T，速度为（－ωe，ωe，0）T，则t时刻其在轨道直角坐标系中坐标为（X＊，Y＊，0）T＝Rz（α），，0）T，速度为（－ωeY＊，ωeX＊，0）T。ωe，α与3.1节的意义相同。

t时刻IGSO卫星在地心惯性坐标系下位置和速度：

式中，X、Y、Z、VX、VY、VZ为卫星在惯性坐标系下位置和速度，i表示轨道倾角，X＊、Y＊为卫星在轨道直角坐标系下平面坐标。

将惯性系下位置左乘矩阵RT，R同式 （6），即可得地固系下卫星位置，再转成大地坐标系 （B，L，H）即可得IGSO卫星的星下点轨迹，图1为经度为45°，倾角为10°和30°的IGSO卫星星下点轨迹 （二体意义下）。

4 算例分析与比较

对不同经度的GEO卫星、不同经度和倾角的IGSO卫星，先使用变步长的RKF7（8）阶积分法积分10步，截断误差设定为10－12，再使用定步长75s的10阶Adams预测－校正系统进行积分48h（2圈），每隔15min输出一组结果。不同经度的GEO卫星在地心地固坐标系下的积分初值见表1，其速度均为0；不同经度和倾角的IGSO卫星在惯性坐标系下的积分初值见表2。

图1 IGSO卫星星下点轨迹 （二体意义）

表1 GEO卫星积分初值

表2 IGSO卫星积分初值

由于篇幅原因，仅给出了Y方向差异的统计信息和经度为30°的GEO卫星差异变化图。图2是经度为30°的GEO卫星轨道数值积分结果与二体意义下结果在地心地固坐标系下的差异，历元间隔为15分钟，表3是GEO卫星轨道积分结果与真值 （二体意义）在地固系下Y方向差异的统计信息。

图2 GEO卫星轨道数值积分差异

表3 GEO卫星轨道数值积分Y方向差异统计信息

由于篇幅原因，仅给出了X方向差异的统计信息和经度为45°倾角为30°的IGSO卫星差异变化图。图3是经度为45°倾角为30°的IGSO卫星轨道数值积分结果与二体意义下结果在地心地固坐标系下的差异，历元间隔为15min，表4是IGSO卫星轨道积分结果与真值 （二体意义）在地固系下X方向差异的统计信息。

表4 IGSO卫星轨道数值积分X方向差异统计信息

图3 IGSO卫星轨道数值积分差异

5 结论

由以上图表可以得出，二体意义下GEO、IGSO卫星轨道数值积分结果与真值相差甚小，均在亚毫米级或以内，与现有精密定轨精度相比可以忽略，适当减小积分步长可获得更小的理论值差异，但这将明显增加计算耗时。结果表明对MEO（如GPS）卫星适用的轨道数值积分方法对GEO和IGSO卫星也同样适用。

本文仅在二体意义下验证了数值积分方法的适用性，现实中导航卫星摄动力复杂多变，因此建立合适的力学模型再结合数值积分方法进行轨道积分是统计定轨中的重要步骤。

［1］刘吉华，欧吉坤，钟世民，等.GEO与IGSO卫星星间差分的精密定轨［C］//中国卫星导航学术年会组委会.第一届中国卫星导航学术年会文集.北京：出版社不详，2010：1－7.

［2］李得海，袁运斌，欧吉坤，等.12阶Runge－kutta2次算法的卫星轨道积分研究［J］.武汉大学学报：信息科学版，2010，35（11）：1335－1338.

［3］赵齐乐.GPS导航星座及低轨卫星的精密定轨理论和软件研究［D］.武汉：武汉大学，2004.

［4］杨 洋，范 丽，董绪荣.IGSO星座分析与优化设计［J］.兵工自动化，2008，27（11）：53－55.

［5］ES－HAGH M.Step Variable Numerical Orbit Integration of a Low Earth Orbiting Satellite［J］.Journal of the Earth ＆Space Physics，2005，31（1）：1－12.

［6］XU Guo－chang.GPS Theory Algorithms and Applications［M］.Berlin：Springer－Verlag，2003.

［7］IERS Convention Centre.IERS Conventions：2010［M］.Frankfurt am Main：Verlag des Bundesamts für Kartographie und Geodäsie，2010.

［8］LEICK A.GPS Satellite Surveying［M］.Hoboken：John Wiley ＆Sons Inc.，2004.