低轨卫星广播星历参数设计

方善传，杜 兰，高 燕

(1.太原卫星发射中心技术部，山西 太原 030027；2.信息工程大学 地理空间信息学院，河南 郑州 450001)

0 引言

传统的导航增强系统都是采用地球同步轨道(Geosynchronous Earth Orbit，GEO)卫星作为通信媒介，提供差分校正值和完好性数据[1]。虽然GEO仅需要少量的卫星就可以向全球用户发播导航信息，但是其轨道资源有限。按照国际电信联盟规定，GEO卫星每1.5°可部署一颗，总共只有240个轨位[2-3]；其独特的高轨、静地和频繁的轨道定点机动，导致其轨道精度较差[4-5]；较高的轨道导致信号延迟较多，且用户接收到的卫星信号较弱，容易受到干扰。

LEO卫星作为导航增强系统的新兴力量，具有独特的优势。其轨道高度低，轨道资源较为丰富。随着卫星发射技术的进步，LEO卫星发射成本迅速降低，远低于发射中/高轨卫星。随着星载GNSS测定轨技术的成熟和辅助加速度计在摄动力测定中的应用，LEO卫星定轨精度已经优于分米级[6-7]。同时，由于LEO卫星轨道高度低，离地面站较近，卫星信号更强，信号延迟更少[8]。基于LEO的导航增强系统成为目前LEO卫星应用推广发展的趋势之一[9]。

导航电文是卫星导航系统顶层设计的重要组成部分，其结构设计直接决定了首次定位时间、通信效率和可扩展性等关键性能指标。而广播星历是导航电文的核心，其表示方式直接决定了卫星轨道精度和通信效率等[10]。目前，广播星历参数主要有两大类：一类是基于轨道根数型的广播星历，GPS，BDS，GALILEO，QZSS都是使用基于轨道根数型的广播星历[11-12]；另一类是基于轨道状态型的广播星历，GLONASS正在使用这种类型的广播星历[13-14]。对于低轨导航增强卫星，目前还没有成熟的广播星历模型可供使用。

在LEO卫星的星历设计方面，文献[15]参照GPS广播星历16参数拟合法特点，根据约1 000 km高度LEO卫星的轨道特征，设计了一套针对LEO卫星的25参数星历模型，拟合2 h卫星位置误差的均方差小于10 m，部分弧段拟合误差小于25 m，但拟合时间过长，拟合精度较差，参数较多。文献[16]将GPS 16参数星历模型直接用于LEO卫星，星历拟合精度与GPS卫星相当，拟合时长缩短为10～20 min，为确定LEO卫星的星历拟合时间提供了参考。由于LEO导航增强卫星大多数采用近圆轨道，而轨道根数型星历用于近圆轨道拟合时，开普勒根数的强负相关性和拟合算法的矩阵严重病态，会严重影响拟合算法迭代收敛的速度和拟合精度。

本文根据轨道列表型星历模型的无奇点特性，设计了2种基于轨道状态型的16参数星历模型。首先，针对300～1 000 km高度的低轨卫星，根据其摄动力量级和短期动力学特性，在ECEF坐标系下构建经验力模型，利用二次多项式/周期项+常数项来吸收经验力摄动加速度的影响。该方法既保留了轨道列表型星历的无奇点特性，又能够利用较少的参数进行高精度的星历拟合。将拟合时间设计为用户的卫星最大可见时间，降低用户星历更新频率。采用300～1 500 km高度高精度轨道预报数据，通过星历拟合验证了设计的2种星历模型。拟合结果表明，方案1的URE(User Range Error)精度优于11 cm，方案2的URE精度优于9 cm。

1 LEO广播星历参数设计

广播星历的设计原则是利用尽量少的参数实验高精度的卫星位置拟合，减少通信资源损耗，缩短用户首次定位时间。通常广播星历拟合精度要求URE小于10 cm。因此广播星历设计时需要考虑主要摄动力，对于LEO导航增强卫星，其受到的主要摄动力为地球扁率J2摄动、大气阻力摄动、日月引力摄动和其他非球形引力摄动主项。其中，地球扁率摄动力模型简单，可以直接采用。对于力模型计算较为复杂的大气阻力和日月引力摄动，通常简化处理方式是针对其短期内的摄动加速度变化规律，采用简单的数学函数进行模型逼近，即经验力建模。对于短期轨道积分，更小量级的摄动力则可以忽略。

在主要摄动力影响下，存在趋势性和周期性变化的摄动力在短弧段内具有如下规律[17]：

① 由于大气阻力的影响，卫星位置分量呈现逐渐减小的趋势性变化，可考虑用多项式吸收趋势项经验力的影响；

② 由于日月引力和其他摄动主项的影响，摄动力随着卫星运行呈现周期性变化，变化周期与轨道周期相同，可考虑用二次多项式或三角函数吸收周期项经验力的影响。

1.1 广播星历参数模型

在经验力的趋势性和周期性变化分析基础上，综合考虑参数个数、拟合精度以及用户卫星位置速度的计算复杂度等，提出了2种星历模型方案。

1.1.1 方案1

采用二次多项式来拟合经验力的摄动加速度，即

(1)

设计的16参数 LEO广播星历在toe时刻表示为：

x16=(toe，x，y，z，vx，vy，vz，a0X，a0Y，a0Z，a1X，a1Y，a1Z，a2X，a2Y，a2Z)T。

1.1.2 方案2

采用常数项加周期项来拟合经验力的摄动加速度，即

(2)

设计的16参数 LEO广播星历在toe时刻表示为：

x16=(toe，x，y，z，vx，vy，vz，a0X，a0Y，a0Z，AX，AY，AZ，BX，BY，BZ)T。

1.2 拟合算法

广播星历的参数x需要通过一段时间内精密外推的卫星位置序列R(默认为真值，m×3阶，含m组卫星位置序列)拟合得来，其计算方法有很多，本文使用最小二乘法进行拟合[18-20]。

设参考历元星历参数的状态向量为x(t0)，其概略初值可以取为：

x(t0)=(t0，x，y，z，vx，vy，vz，0，0，0，0，0，0，0，0，0)T，

式中，(x，y，z，vx，vy，vz)为t0时刻卫星位置分量，其他摄动参数初值取为0。

卫星位置向量序列R与参考历元星历参数x(t0)之间的关系为：

R=Y(x(t0)，t)，

(3)

式中，t为m组卫星位置序列对应的时间向量序列；函数Y(x(t0)，t)的算法即是前面介绍的用户算法。使用最小二乘法对式(3)进行线性化之后得到：

(4)

令y为卫星位置真值与计算值Y(xi(t0)，t)之差，即

y=R-Y(xi(t0)，t)，

(5)

令

(6)

Δx=x-xi，

(7)

略去式(4)中的二阶及其以上的高阶项，可得：

y=AΔx，

(8)

根据最小二乘法原理，可得Δx的最优估值为

Δx=(ATA)-1ATy，

(9)

迭代公式为：

xi=x(t0)+Δx，

(10)

迭代结束的条件为：

|δi-δi-1|/δi-1<ε，

(11)

式中，ε是预先给定的收敛条件，其值越小，迭代次数越多，迭代精度越高；δi为第i次迭代时y中的误差统计值，其计算方法为

(12)

1.3 用户算法

用户接收到卫星播发的轨道列表型星历后，需要通过星历提供的广播星历参数，计算出卫星在可见时间内某一时刻的位置速度[22-23]。

(13)

(14)

1.4 LEO的星历拟合精度指标

URE性能指标反映了导航星历和预报钟差的精度，是计算星历和钟差误差对定位影响的一个重要评价指标[25]。在卫星导航系统中，广播星历参数拟合的限制性要求是拟合结果等效URE损失小于0.1 m[26]。URE的来源主要是轨道外推误差和钟差误差，这里仅分析星历参数拟合误差的单项影响，计算公式为[27]：

(15)

式中，ΔR，ΔT，ΔN分别为卫星径向、沿迹和外法向位置分量的拟合误差；AR，AT，AN分别对应3个方向的贡献因子。表1列出了利用上述公式计算出的不同高度近圆轨道的URE贡献因子值大小(用户最小观测仰角设为0°)。随着轨道升高，T方向和N方向的影响逐渐减小，R方向影响逐渐增大。

表1 500～1 000 km高度URE贡献因子

轨道高度/kmARATAN3000.374 30.655 70.655 75000.455 70.629 40.629 47000.515 30.606 00.606 09000.562 20.584 80.584 81 1000.600 70.565 30.565 31 3000.633 10.547 30.547 31 5000.661 00.530 60.530 620 2000.979 40.142 70.142 7

2 拟合实验

为了验证星历参数模型的适用性，将设计的2种星历方案用于300～1 500 km高度的低轨卫星轨道拟合。首先对卫星在ECI坐标系下进行复杂力模型下的轨道外推(数据点间隔30 s)，然后采用IERS发布的地球定向参数将轨道数据转换至ECEF下，外推过程中的摄动力模型设置如表2所示。按照拟合时段长度对实验数据进行分组，历元时刻取每组数据的中间时刻。

表2 摄动力模型参数设置

摄动力模型模型参数地球非球形引力WGS84_EGM96① 阶数:21∗21② 地球半径r=6 378.1 363 km③ 地球自转角速度ω=7.292 158 553×10-5 rad/sec大气阻力大气密度模型Harris-priester① CD=2.2② A/M=0.02 m2/kg③ F10.7(2016)=150三体摄动太阳、月球① 太阳引力常数: 1.327 122 000 000×1011 km3/sec2② 月球引力常数: 4.902 801 076 000×103 km3/sec2太阳光压阴影模型:Dual Cone① Cr=1② A/M=0.02 m2/kg

对高度为300～1 500 km、倾角为45°、偏心率为0.001的卫星轨道进行了2种星历方案的实验验证，拟合时长设置为卫星最大可见时间附近的整分钟数，如图1所示。

图1 不同高度卫星的最大可见时间

2.1 方案1星历拟合实验

方案1的星历拟合实验结果如图2和图3所示。

图2 方案1星历R/T/N方向拟合误差统计

图3 方案1星历URE统计

实验过程中，LEO的拟合成功率均为100%，拟合迭代次数均稳定为3或4次，拟合消耗总时间376 s(处理器：i5-5200U，内存：8 GB)。对于300～1 500 km高度的卫星，随着轨道高度增加，拟合时间变长，R/T/N三个方向位置分量误差以及URE都会增大，结合文献[17]的结论(相同拟合时段条件下，拟合误差随轨道升高而减小)可知，导致拟合误差增大的原因是拟合时间的延长；R/T/N三个方向平均拟合精度分别保持在12/11/3 cm以内，N方向拟合误差比R，T方向小8～9 cm；URE均方根误差可以保持在11 cm以内。若要保证星历URE精度小于10 cm的要求，轨道高度应控制在1 400 km以下。

星历拟合过程中的900 km高度轨道URE序列如图4所示。除了个别组URE大于20 cm外，大多数结果波动范围在0～20 cm。

图4 方案1星历URE序列(900 km高度轨道)

2.2 方案2星历拟合实验

方案2的拟合实验结果如图5和图6所示。实验过程中，LEO的拟合成功率均为100%，拟合迭代次数均稳定为3或4次，拟合消耗总时间525 s(硬件配置与方案1相同)。可以看出：R/T/N三个方向拟合误差分别保持在10/10/3 cm以内；URE的均方根误差保持在9 cm以内。

图5 方案2星历RTN方向拟合误差统计

图6 方案2星历URE统计

星历拟合过程中的900 km高度轨道URE序列如图7所示。URE波动范围在0～22 cm之间，波动范围相对于方案1更小，数据集中程度比方案1更高，说明了方案2的算法稳定性更好。

图7 方案2星历URE序列(900 km高度轨道)

2种方案的实验结果表明，方案2星历拟合精度优于方案1，周期项+常数项能够更好地拟合低轨卫星摄动力影响，但方案2算法时间复杂度比方案1高。此外，设计的2种轨道列表型16参数星历模型拟合精度稍优于GPS 16参数星历模型(GPS 16参数星历模型精度直接用于低轨卫星拟合10～20 min，精度为分米级[28])。

2.3 IRIDIUM卫星数据拟合实验

为验证真实轨道的数据拟合效果，对采用IRIDIUM增强系统的卫星轨道进行了验证。IRIDIUM系统卫星设计轨道高度为780 km，倾角86.4°，卫星总数66+6(备用)个，均匀分布在6个圆轨道面上。在每个轨道面各取一颗卫星进行实验，选取2010年3月1-10日共10天的卫星轨道预报数据，力模型设置与前面相同。方案1和方案2的星历拟合实验结果如图8和图9所示。

图8 方案1星历拟合误差统计

图9 方案2星历拟合误差统计

拟合时间均为16 min(卫星用户最大可见时间)。实验过程中，LEO的拟合成功率均为100%，拟合迭代次数稳定为3次。可以看出：方案1中R/T/N三个方向拟合误差基本保持在9/8/3 cm以内；URE的均方根误差保持在8 cm以内；方案2中R/T/N三个方向拟合误差基本分别保持在8/8/3 cm以内；URE的均方根误差保持在8 cm以内。2种方案拟合结果精度与前面进行的实验结果精度基本相同，方案2精度稍高于方案1。

3 结束语

针对低轨卫星摄动力的趋势项和周期项特征，利用二次多项式/周期项+常数项吸收经验力摄动加速度影响，设计了2种16参数星历模型。设计的星历模型具有参数个数少、拟合精度高和无奇点的优点。

基于300～1 500 km高度的LEO拟合实验表明：① 方案1中R/T/N三个方向平均拟合精度分别优于12/11/3 cm，URE精度优于11 cm。方案2中R/T/N三个方向拟合误差基本分别优于10/10/3 cm，URE精度优于9 cm；② 2种拟合方案的参数个数相同，方案2拟合精度更高，拟合结果更稳定，但方案2时间复杂度更高。若对时间复杂度有较高要求，建议采用方案1的星历模型；若对星历精度有较高要求，建议采用方案2的星历模型。