一种电推进GEO卫星自学习定轨方法

一种电推进geo卫星自学习定轨方法
技术领域
1.本发明涉及卫星定轨领域，具体涉及一种电推进geo卫星自学习定轨方法。


背景技术：

2.传统的卫星定轨主要依赖地面测控，然而，电推进geo卫星在入轨过程中需要进行长期持续轨道机动，这使得基于地面测控的定轨方式具有诸多局限性。其一，地面测控站资源有限，对于电推进geo卫星的轨道转移全过程，难以做到实时全程监控；其二，地面测控对于长期持续机动卫星的定轨精度不够高，已经严重制约了电推进geo卫星的最终入轨定点捕获精度；其三，一旦地面测控站受到干扰或破坏，卫星将难以正常的飞行和执行预定的任务。因此，卫星自主定轨技术是未来电推进geo卫星定轨的必然发展方向。自主定轨技术是航天器自主运行的核心技术，是航天器实现轨道姿态自主控制、精确定点捕获和执行在轨服务等空间任务的前提。拥有自主定轨能力的电推进geo卫星不需要通过全球布站来提供全过程定轨支持，从而可以大幅减轻地面测控的负担，提高定轨的实时性、准确性和系统的安全性。
3.卫星自主定轨是指卫星利用自身携带的测量设备来实时确定自己的位置和速度的一项技术。惯性导航系统依靠携带的imu(惯性导航单元)，通过测量卫星的加速度，并根据轨道动力学模型进行积分运算，获得卫星瞬时速度和瞬时位置数据。惯性导航具有完全自主的优势，但是随着运行时间增加，轨道动力学模型积分导致的系统误差将不断增大，因此需要借助其他导航信息修正。利用全球卫星导航系统，可通过geo卫星携带的gnss接收机接受gnss信号(如北斗、gps的信号等)，从而解算获得卫星的位置和速度信息。利用gnss的定轨能够提供高精度的导航信息，但是相比惯性导航，其信息输出频率较低，且受到导航卫星可见性限制，gnss信号对高轨卫星难以全覆盖，在高动态条件下不能保证接收到足够数量的gnss信号。
4.然而，当前的基于组合导航系统的自主定轨方法大多适用于不进行机动的卫星。若卫星进行长期持续轨道机动，则轨道动力学模型积分所得数据的误差累积，严重影响定轨精度。此外，当前导航滤波算法大多基于卡尔曼滤波理论设计，例如ekf(扩展卡尔曼滤波)、ukf(无迹卡尔曼滤波)等。基于卡尔曼滤波理论设计的导航滤波器依赖于先验已知的系统噪声方差和测量噪声方差。但在实际的卫星运行环境中，受地球非球形引力摄动、大气阻力摄动、太阳光压摄动等多种摄动、空间热环境变化、卫星机动导致的观测几何变化等影响，系统噪声和测量噪声的统计特性存在不确定性，这导致传统导航滤波算法的精度上限难以进一步提升。


技术实现要素：

5.针对现有技术中的上述不足，本发明提供的一种电推进geo卫星自学习定轨方法解决了现有电推进geo卫星自主定轨精度低的问题。
6.为了达到上述发明目的，本发明采用的技术方案为：
7.提供一种电推进geo卫星自学习定轨方法，其包括以下步骤：
8.s1、对电推进geo卫星建立带动力改进春分点轨道根数模型；获取电推进geo卫星的加速度信息；
9.s2、根据带动力改进春分点轨道根数模型获取电推进geo卫星以位置和速度矢量为状态变量的系统状态方程；
10.s3、根据加速度信息和系统状态方程对位置和速度矢量进行时间更新，获取电推进geo卫星的理论轨道信息；
11.s4、判断电推进geo卫星是否接收到足够的gnss信号，若是则进入步骤s5；否则进入步骤s8；
12.s5、基于gnss信号，采用ukf算法进行第一次多源信息融合，获取电推进geo卫星的第一状态变量估计值和第一状态变量协方差估计矩阵；
13.s6、将上一时刻的状态变量值和加速度信息输入lstm网络，获取对应的误差预测值，并结合第一状态变量估计值和第一状态变量协方差估计矩阵进行第二次多源信息融合，获取电推进geo卫星的第二状态变量估计值，并将其作为电推进geo卫星的最终状态估计值；
14.s7、判断是否结束自学习定轨，若是则结束；否则返回步骤s3；
15.s8、根据星历和电推进geo卫星上的太阳帆板输出功率信息，计算当前卫星飞行周期内的地影信息；
16.s9、通过地影信息分段补偿法修正带动力改进春分点轨道根数模型，得到修正后的系统状态方程；
17.s10、通过修正后的系统状态方程进行轨道递推，对位置和速度矢量进行时间更新，获取电推进geo卫星的修正后理论轨道信息，并将其作为电推进geo卫星的最终状态估计值；返回步骤s7。
18.进一步地，步骤s1中带动力改进春分点轨道根数模型为：
[0019][0020][0021][0022][0023]
[0024][0025][0026]
其中[p,f,g,h,k,l]
t
为带动力改进春分点轨道根数模型的状态变量，p为轨道半通径，f和g为偏心率矢量分量，h和k为轨道倾角矢量分量，l为真经度，[
·
]
t
表示矩阵的转置；μ为地球引力常数；w为常数，w＝1+fcosl+gsinl；f
t
为电推进geo卫星的切向加速度；fr为电推进geo卫星的径向加速度；fn为电推进geo卫星的法向加速度；s为常数，s2＝1+h2+k2；m为电推进geo卫星质量；t为电推进geo卫星的推力大小；g0为海平面处重力加速度；i
sp
为电推进系统的比冲；字符上一点表示求导。
[0027]
进一步地，步骤s2的具体方法为：
[0028]
根据变量转换关系将改进春分点轨道根数模型的状态变量[p,f,g,h,k,l]
t
转换为位置、速度矢量[r,v]
t
，得到以位置和速度矢量为状态变量的系统状态方程。
[0029]
进一步地，步骤s5的具体方法包括以下子步骤：
[0030]
s5-1、根据公式：
[0031][0032][0033]
获取状态初始值x0的期望值和状态变量协方差矩阵初始值p0，并将期望值作为初始状态变量的均值；其中(
·
)
t
表示矩阵的转置；e(
·
)表示求期望；
[0034]
s5-2、根据公式：
[0035][0036]
构造2n+1个sigma点；其中χk为k时刻的sigma矢量；为k时刻状态变量的均值；p
kxx
为k时刻状态变量的协方差矩阵；n为状态向量长度；λ为比例参数；表示p
kxx
的cholesky分解；
[0037]
s5-3、根据公式：
[0038][0039]
计算每个sigma点对应的权重；其中w
im
为计算均值的权重系数；w
ic
为计算协方差的权重系数；α和κ均为尺度参数；
[0040]
s5-4、根据公式：
[0041][0042][0043][0044]
建立预测方程；其中为k+1时刻状态变量的sigma矢量；φ(χ
i,k
,u
i,k
)表示非线性系统状态转移方程；χ
i,k
为k时刻状态变量的sigma矢量；u
i,k
表示k时刻的控制参数；为k+1时刻通过sigma矢量得到的状态变量的先验均值；为k+1时刻通过sigma矢量得到的状态变量的先验协方差矩阵；q为系统噪声协方差矩阵；
[0045]
s5-5、基于gnss信号获取k+1时刻的状态变量的先验均值
[0046]
s5-6、根据公式：
[0047][0048][0049]
[0050][0051][0052][0053]
建立测量修正方程；其中h(
·
)为观测方程，y
i,k+1
为k+1时刻对测量得到的状态变量进行转移得到的sigma矢量；为k+1时刻测量得到的状态变量的先验协方差矩阵；r为量测噪声协方差矩阵；为k+1时刻的交叉协方差阵；k
k+1
为k+1时刻的卡尔曼增益系数；为k+1时刻更新后的状态变量值，即第一状态变量估计值；y
k+1
为k+1时刻通过gnss测量得到的状态变量；为k+1时刻更新后的状态变量协方差矩阵；为k
k+1
的转置。
[0054]
进一步地，步骤s6的具体方法包括以下子步骤：
[0055]
s6-1、将上一时刻的卫星位置、速度矢量和加速度信息输入至训练好的lstm模型中，获取lstm模型预测的状态变量值和lstm模型预测的定轨误差估计向量；
[0056]
s6-2、根据公式：
[0057][0058]
进行第二次多源信息融合，获取电推进geo卫星的第二状态变量估计值并将其作为电推进geo卫星的最终状态估计值；其中为电推进geo卫星的第一状态变量估计值；p
ukf
(tk)为tk时刻下ukf算法中的状态变量的协方差矩阵；p
lstm
(tk)为tk时刻下lstm模型预测得到的定轨误差估计向量；为tk时刻下lstm模型预测得到的状态变量值。
[0059]
进一步地，步骤s8的具体方法包括以下子步骤：
[0060]
s8-1、以地心为起点，以太阳光线照射方向为方向作射线，记为地影中线；
[0061]
s8-2、根据电推进geo卫星上的太阳帆板输出功率信息，获取真实的地影进入时刻t
en
和离开时刻t
ex
；
[0062]
s8-3、根据公式：
[0063][0064]
获取电推进geo卫星到达地影中线的时刻t
mid
；
[0065]
s8-4、根据公式：
[0066][0067]
获取地影中线真经度其中l0为初始状态的真经度；w
*
为地影中线转动的角速度。
[0068]
进一步地，步骤s9的具体方法包括以下子步骤：
[0069]
s9-1、根据公式：
[0070][0071][0072]
获取电推进geo卫星绕地球旋转的当前圈的真实真经度平均角速度和电推进geo卫星绕地球旋转的当前圈的估计真经度角速度；其中为电推进geo卫星绕地球旋转的上一圈数据对应的地影中线真经度；电推进geo卫星绕地球旋转的上一圈数据对应的地影中线时刻；l
es
为通过带动力改进春分点轨道根数模型递推获得的对应时刻的卫星真经度估计值；
[0073]
s9-2、根据公式：
[0074][0075]
获取电推进geo卫星绕地球旋转的当前圈的真经度角速度补偿值wc；其中为电推进geo卫星绕地球旋转的上一圈的真经度角速度补偿值；
[0076]
s9-3、根据公式：
[0077][0078][0079][0080][0081][0082][0083][0084]
修正带动力改进春分点轨道根数模型，得到修正后的系统状态方程；其中p为轨道半通径，f和g为偏心率矢量分量，h和k为轨道倾角矢量分量，l为真经度，[
·
]
t
表示矩阵的转置；μ为地球引力常数；w为常数，w＝1+fcosl+gsinl；f
t
为电推进geo卫星的切向加速度；fr为电推进geo卫星的径向加速度；fn为电推进geo卫星的法向加速度；s为常数，s2＝1+h2+k2；m为电推进geo卫星的质量；t为电推进geo卫星的推力大小；g0为海平面处重力加速度；i
sp
为电推进系统的比冲；字符上一点表示求导。
[0085]
本发明的有益效果为：
[0086]
1、本发明对电推进geo卫星通过加速度计对卫星的加速度进行测量，将卫星的控制加速度信息引入系统状态方程，减小轨道动力学模型误差，使轨道递推结果更精确。
[0087]
2、当有gnss定轨信息时，本发明利用lstm算法对gnss导航的定轨误差进行估计，然后将ukf融合结果进行二次融合，进一步提升卫星自主定轨精度。
[0088]
3、当gnss不可用时，本发明在不额外增加星载敏感器的条件下，利用卫星上的太阳帆板输出功率信息，获取真实的地影进入时刻和离开时刻，从而构建地影信息。根据卫星飞行每圈的进出地影时间，构建卫星真经度角速度补偿量，从而减小轨道动力学模型误差。
附图说明
[0089]
图1为本方法的流程示意图；
[0090]
图2为地球阴影区域示意图。
具体实施方式
[0091]
下面对本发明的具体实施方式进行描述，以便于本技术领域的技术人员理解本发明，但应该清楚，本发明不限于具体实施方式的范围，对本技术领域的普通技术人员来讲，只要各种变化在所附的权利要求限定和确定的本发明的精神和范围内，这些变化是显而易见的，一切利用本发明构思的发明创造均在保护之列。
[0092]
如图1所示，该电推进geo卫星自学习定轨方法包括以下步骤：
[0093]
s1、对电推进geo卫星建立带动力改进春分点轨道根数模型；获取电推进geo卫星的加速度信息；
[0094]
s2、根据带动力改进春分点轨道根数模型获取电推进geo卫星以位置和速度矢量为状态变量的系统状态方程；
[0095]
s3、根据加速度信息和系统状态方程对位置和速度矢量进行时间更新，获取电推进geo卫星的理论轨道信息；
[0096]
s4、判断电推进geo卫星是否接收到足够的gnss信号，若是则进入步骤s5；否则进入步骤s8；
[0097]
s5、基于gnss信号，采用ukf算法进行第一次多源信息融合，获取电推进geo卫星的第一状态变量估计值和第一状态变量协方差估计矩阵；
[0098]
s6、将上一时刻的状态变量值和加速度信息输入lstm网络，获取对应的误差预测值，并结合第一状态变量估计值和第一状态变量协方差估计矩阵进行第二次多源信息融合，获取电推进geo卫星的第二状态变量估计值，并将其作为电推进geo卫星的最终状态估计值；
[0099]
s7、判断是否结束自学习定轨，若是则结束；否则返回步骤s3；
[0100]
s8、根据星历和电推进geo卫星上的太阳帆板输出功率信息，计算当前卫星飞行周期内的地影信息；
[0101]
s9、通过地影信息分段补偿法修正带动力改进春分点轨道根数模型，得到修正后的系统状态方程；
[0102]
s10、通过修正后的系统状态方程进行轨道递推，对位置和速度矢量进行时间更新，获取电推进geo卫星的修正后理论轨道信息，并将其作为电推进geo卫星的最终状态估计值；返回步骤s7。
[0103]
步骤s1中带动力改进春分点轨道根数模型为：
[0104][0105][0106][0107][0108][0109][0110][0111]
其中[p,f,g,h,k,l]
t
为带动力改进春分点轨道根数模型的状态变量，p为轨道半通径，f和g为偏心率矢量分量，h和k为轨道倾角矢量分量，l为真经度，[
·
]
t
表示矩阵的转置；μ为地球引力常数；w为常数，w＝1+fcosl+gsinl；f
t
为电推进geo卫星的切向加速度；fr为电推进geo卫星的径向加速度；fn为电推进geo卫星的法向加速度；s为常数，s2＝1+h2+k2；m为电推进geo卫星质量；t为电推进geo卫星的推力大小；g0为海平面处重力加速度；i
sp
为电推进系统的比冲；字符上一点表示求导。
[0112]
步骤s2的具体方法为：根据变量转换关系将改进春分点轨道根数模型的状态变量[p,f,g,h,k,l]
t
转换为位置、速度矢量[r,v]
t
，得到以位置和速度矢量为状态变量的系统状态方程。
[0113]
步骤s5的具体方法包括以下子步骤：
[0114]
s5-1、根据公式：
[0115][0116]
[0117]
获取状态初始值x0的期望值和状态变量协方差矩阵初始值p0，并将期望值作为初始状态变量的均值；其中(
·
)
t
表示矩阵的转置；e(
·
)表示求期望；
[0118]
s5-2、根据公式：
[0119][0120]
构造2n+1个sigma点；其中χk为k时刻的sigma矢量；为k时刻状态变量的均值； p
kxx
为k时刻状态变量的协方差矩阵；n为状态向量长度； λ为比例参数；表示p
kxx
的cholesky分解；
[0121]
s5-3、根据公式：
[0122][0123]
计算每个sigma点对应的权重；其中w
im
为计算均值的权重系数；w
ic
为计算协方差的权重系数；α和κ均为尺度参数；β为常数；i表示sigma点的序号；α用于确定周围sigma点的分布情况，一般取一小正数(0.0001≤α≤1)；κ通常为0或3-n；β是一个与状态的先验分布信息有关的参数，对于高斯分布β＝2最优。
[0124]
s5-4、根据公式：
[0125][0126][0127][0128]
建立预测方程；其中为k+1时刻状态变量的sigma矢量；φ(χ
i,k
,u
i,k
)表示非线性系统状态转移方程；χ
i,k
为k时刻状态变量的sigma矢量；u
i,k
表示k时刻的控制参数；为k+1时刻通过sigma矢量得到的状态变量的先验均值；为k+1时刻通过sigma矢量得到的
状态变量的先验协方差矩阵；q为系统噪声协方差矩阵；
[0129]
s5-5、基于gnss信号获取k+1时刻的状态变量的先验均值
[0130]
s5-6、根据公式：
[0131][0132][0133][0134][0135][0136][0137]
建立测量修正方程；其中h(
·
)为观测方程，y
i,k+1
为k+1时刻对测量得到的状态变量进行转移得到的sigma矢量；为k+1时刻测量得到的状态变量的先验协方差矩阵；r为量测噪声协方差矩阵；为k+1时刻的交叉协方差阵；k
k+1
为k+1时刻的卡尔曼增益系数；为k+1时刻更新后的状态变量值，即第一状态变量估计值；y
k+1
为k+1时刻通过gnss测量得到的状态变量；为k+1时刻更新后的状态变量协方差矩阵；为k
k+1
的转置。
[0138]
步骤s6的具体方法包括以下子步骤：
[0139]
s6-1、将上一时刻的卫星位置、速度矢量和加速度信息输入至训练好的lstm模型中，获取lstm模型预测的状态变量值和lstm模型预测的定轨误差估计向量；
[0140]
s6-2、根据公式：
[0141][0142]
进行第二次多源信息融合，获取电推进geo卫星的第二状态变量估计值并将其作为电推进geo卫星的最终状态估计值；其中为电推进geo卫星的第一状态变量估计值；p
ukf
(tk)为tk时刻下ukf算法中的状态变量的协方差矩阵；p
lstm
(tk)为tk时刻下lstm模型预测得到的定轨误差估计向量；为tk时刻下lstm模型预测得到的状态变量值。
[0143]
步骤s8的具体方法包括以下子步骤：
[0144]
s8-1、以地心为起点，以太阳光线照射方向为方向作射线，记为地影中线；
[0145]
s8-2、根据电推进geo卫星上的太阳帆板输出功率信息，获取真实的地影进入时刻t
en
和离开时刻t
ex
；
[0146]
s8-3、根据公式：
[0147][0148]
获取电推进geo卫星到达地影中线的时刻t
mid
；
[0149]
s8-4、根据公式：
[0150][0151]
获取地影中线真经度其中l0为初始状态的真经度；w
*
为地影中线转动的角速度。
[0152]
步骤s9的具体方法包括以下子步骤：
[0153]
s9-1、根据公式：
[0154][0155][0156]
获取电推进geo卫星绕地球旋转的当前圈的真实真经度平均角速度和电推进geo卫星绕地球旋转的当前圈的估计真经度角速度；其中为电推进geo卫星绕地球旋转的上一圈数据对应的地影中线真经度；电推进geo卫星绕地球旋转的上一圈数据对应的地影中线时刻；l
es
为通过带动力改进春分点轨道根数模型递推获得的对应时刻的卫星真经度估计值；
[0157]
s9-2、根据公式：
[0158][0159]
获取电推进geo卫星绕地球旋转的当前圈的真经度角速度补偿值wc；其中为电推进geo卫星绕地球旋转的上一圈的真经度角速度补偿值；
[0160]
s9-3、根据公式：
[0161][0162][0163][0164]
[0165][0166][0167][0168]
修正带动力改进春分点轨道根数模型，得到修正后的系统状态方程；其中p为轨道半通径，f和g为偏心率矢量分量，h和k为轨道倾角矢量分量，l为真经度，[
·
]
t
表示矩阵的转置；μ为地球引力常数；w为常数，w＝1+fcosl+gsinl；f
t
为电推进geo卫星的切向加速度；fr为电推进geo卫星的径向加速度；fn为电推进geo卫星的法向加速度；s为常数，s2＝1+h2+k2；m为电推进geo卫星的质量；t为电推进geo卫星的推力大小；g0为海平面处重力加速度；i
sp
为电推进系统的比冲；字符上一点表示求导。
[0169]
在具体实施过程中，电推进geo卫星三个方向上的加速度可表示如下：
[0170]fr
＝f
re
+f
rj2
[0171]ft
＝f
te
+f
tj2
[0172]fn
＝f
ne
+f
nj2
[0173]
其中fe＝[f
re
,f
te
,f
ne
]
t
表示电推力产生的加速度矢量，f
j2
＝[f
rj2
,f
tj2
,f
nj2
]
t
表示j2摄动产生的加速度。
[0174]
步骤s3中加速度信息的量测方程为：
[0175]ab
＝a
t
+ac+as[0176]
其中ab为卫星本体系下测量的加速度矢量；a
t
为真实的加速度矢量；ac为加速度计的常值漂移，可以通过在轨标定的方式获得较为准确的常漂估计值；as为测量噪声。
[0177]
步骤s3中时间更新的方法可以采用龙格库塔45数值积分方法，也可以采用其他数值积分方法。
[0178]
步骤s4中是否接收到足够的gnss信号的判断标准为是否接收到以解算gnss定位数据的gnss信号。
[0179]
lstm网络在使用之前可以先进行训练，训练的方式为：构建样本数据集，包括卫星的位置矢量(3维)、速度矢量(3维)和加速度计测量得到的卫星测量加速度矢量(3维)，采用样本数据集对lstm网络进行训练，并使lstm网络的输出状态变量值定轨误差估计向量。
[0180]
如图2所示，地影信息是地球围绕太阳公转而产生的周期性的阴影区域，l
en
和l
ex
分别为卫星进入地影区域和离开地影区域的真赤经。本方法定义一条地影中线，是以地心为起点，太阳光线照射方向为方向的一条射线，l
mid
为卫星到达地影中线的真赤经。由于电推进geo卫星在阴影区域的时间短，轨道变化不明显，因此视为真经度l均匀变化。由于地球绕太阳公转，周期为1年，则可视为均速绕行转动。地影中线的真赤经可简化为绕太阳公转，周期为1年，则可视为均速绕行转动。地影中线的真赤经可简化为
[0181]
综上所述，本发明提出了一种结合lstm算法的基于gnss和惯性导航的组合导航方法，该方法利用ukf算法和lstm进行两次多源信息融合，以减小电推进geo卫星自主定轨误
差。并针对gnss信号弱或不可用的场景，提出了利用地影信息修正轨道积分模型的方法，实现在gnss信号弱或gnss卫星可见性差条件下的电推进geo卫星自主定轨。

技术特征：
1.一种电推进geo卫星自学习定轨方法，其特征在于，包括以下步骤：s1、对电推进geo卫星建立带动力改进春分点轨道根数模型；获取电推进geo卫星的加速度信息；s2、根据带动力改进春分点轨道根数模型获取电推进geo卫星以位置和速度矢量为状态变量的系统状态方程；s3、根据加速度信息和系统状态方程对位置和速度矢量进行时间更新，获取电推进geo卫星的理论轨道信息；s4、判断电推进geo卫星是否接收到足够的gnss信号，若是则进入步骤s5；否则进入步骤s8；s5、基于gnss信号，采用ukf算法进行第一次多源信息融合，获取电推进geo卫星的第一状态变量估计值和第一状态变量协方差估计矩阵；s6、将上一时刻的状态变量值和加速度信息输入lstm网络，获取对应的误差预测值，并结合第一状态变量估计值和第一状态变量协方差估计矩阵进行第二次多源信息融合，获取电推进geo卫星的第二状态变量估计值，并将其作为电推进geo卫星的最终状态估计值；s7、判断是否结束自学习定轨，若是则结束；否则返回步骤s3；s8、根据星历和电推进geo卫星上的太阳帆板输出功率信息，计算当前卫星飞行周期内的地影信息；s9、通过地影信息分段补偿法修正带动力改进春分点轨道根数模型，得到修正后的系统状态方程；s10、通过修正后的系统状态方程进行轨道递推，对位置和速度矢量进行时间更新，获取电推进geo卫星的修正后理论轨道信息，并将其作为电推进geo卫星的最终状态估计值；返回步骤s7。2.根据权利要求1所述的电推进geo卫星自学习定轨方法，其特征在于，步骤s1中带动力改进春分点轨道根数模型为：力改进春分点轨道根数模型为：力改进春分点轨道根数模型为：力改进春分点轨道根数模型为：力改进春分点轨道根数模型为：
其中[p,f,g,h,k,l]
t
为带动力改进春分点轨道根数模型的状态变量，p为轨道半通径，f和g为偏心率矢量分量，h和k为轨道倾角矢量分量，l为真经度，[
·
]
t
表示矩阵的转置；μ为地球引力常数；w为常数，w＝1+fcosl+gsinl；f
t
为电推进geo卫星的切向加速度；f
r
为电推进geo卫星的径向加速度；f
n
为电推进geo卫星的法向加速度；s为常数，s2＝1+h2+k2；m为电推进geo卫星质量；t为电推进geo卫星的推力大小；g0为海平面处重力加速度；i
sp
为电推进系统的比冲；字符上一点表示求导。3.根据权利要求2所述的电推进geo卫星自学习定轨方法，其特征在于，步骤s2的具体方法为：根据变量转换关系将改进春分点轨道根数模型的状态变量[p,f,g,h,k,l]
t
转换为位置、速度矢量[r,v]
t
，得到以位置和速度矢量为状态变量的系统状态方程。4.根据权利要求1所述的电推进geo卫星自学习定轨方法，其特征在于，步骤s5的具体方法包括以下子步骤：s5-1、根据公式：1、根据公式：获取状态初始值x0的期望值和状态变量协方差矩阵初始值p0，并将期望值作为初始状态变量的均值；其中(
·
)
t
表示矩阵的转置；e(
·
)表示求期望；s5-2、根据公式：构造2n+1个sigma点；其中χ
k
为k时刻的sigma矢量；为k时刻状态变量的均值；p
kxx
为k时刻状态变量的协方差矩阵；n为状态向量长度；λ为比例参数；表示p
kxx
的cholesky分解；s5-3、根据公式：
计算每个sigma点对应的权重；其中w
im
为计算均值的权重系数；w
ic
为计算协方差的权重系数；α和κ均为尺度参数；s5-4、根据公式：4、根据公式：4、根据公式：建立预测方程；其中为k+1时刻状态变量的sigma矢量；φ(χ
i,k
,u
i,k
)表示非线性系统状态转移方程；χ
i,k
为k时刻状态变量的sigma矢量；u
i,k
表示k时刻的控制参数；为k+1时刻通过sigma矢量得到的状态变量的先验均值；为k+1时刻通过sigma矢量得到的状态变量的先验协方差矩阵；q为系统噪声协方差矩阵；s5-5、基于gnss信号获取k+1时刻的状态变量的先验均值s5-6、根据公式：6、根据公式：6、根据公式：6、根据公式：
建立测量修正方程；其中h(
·
)为观测方程，y
i,k+1
为k+1时刻对测量得到的状态变量进行转移得到的sigma矢量；为k+1时刻测量得到的状态变量的先验协方差矩阵；r为量测噪声协方差矩阵；为k+1时刻的交叉协方差阵；k
k+1
为k+1时刻的卡尔曼增益系数；为k+1时刻更新后的状态变量值，即第一状态变量估计值；y
k+1
为k+1时刻通过gnss测量得到的状态变量；为k+1时刻更新后的状态变量协方差矩阵；为k
k+1
的转置。5.根据权利要求1所述的电推进geo卫星自学习定轨方法，其特征在于，步骤s6的具体方法包括以下子步骤：s6-1、将上一时刻的卫星位置、速度矢量和加速度信息输入至训练好的lstm模型中，获取lstm模型预测的状态变量值和lstm模型预测的定轨误差估计向量；s6-2、根据公式：进行第二次多源信息融合，获取电推进geo卫星的第二状态变量估计值并将其作为电推进geo卫星的最终状态估计值；其中为电推进geo卫星的第一状态变量估计值；p
ukf
(t
k
)为t
k
时刻下ukf算法中的状态变量的协方差矩阵；p
lstm
(t
k
)为t
k
时刻下lstm模型预测得到的定轨误差估计向量；为t
k
时刻下lstm模型预测得到的状态变量值。6.根据权利要求1所述的电推进geo卫星自学习定轨方法，其特征在于，步骤s8的具体方法包括以下子步骤：s8-1、以地心为起点，以太阳光线照射方向为方向作射线，记为地影中线；s8-2、根据电推进geo卫星上的太阳帆板输出功率信息，获取真实的地影进入时刻t
en
和离开时刻t
ex
；s8-3、根据公式：获取电推进geo卫星到达地影中线的时刻t
mid
；s8-4、根据公式：获取地影中线真经度其中l0为初始状态的真经度；w
*
为地影中线转动的角速度。7.根据权利要求6所述的电推进geo卫星自学习定轨方法，其特征在于，步骤s9的具体方法包括以下子步骤：
s9-1、根据公式：1、根据公式：获取电推进geo卫星绕地球旋转的当前圈的真实真经度平均角速度和电推进geo卫星绕地球旋转的当前圈的估计真经度角速度；其中为电推进geo卫星绕地球旋转的上一圈数据对应的地影中线真经度；电推进geo卫星绕地球旋转的上一圈数据对应的地影中线时刻；l
es
为通过带动力改进春分点轨道根数模型递推获得的对应时刻的卫星真经度估计值；s9-2、根据公式：获取电推进geo卫星绕地球旋转的当前圈的真经度角速度补偿值w
c
；其中为电推进geo卫星绕地球旋转的上一圈的真经度角速度补偿值；s9-3、根据公式：3、根据公式：3、根据公式：3、根据公式：3、根据公式：3、根据公式：3、根据公式：修正带动力改进春分点轨道根数模型，得到修正后的系统状态方程；其中p为轨道半通径，f和g为偏心率矢量分量，h和k为轨道倾角矢量分量，l为真经度，[
·
]
t
表示矩阵的转置；μ为地球引力常数；w为常数，w＝1+fcosl+gsinl；f
t
为电推进geo卫星的切向加速度；f
r
为电
推进geo卫星的径向加速度；f
n
为电推进geo卫星的法向加速度；s为常数，s2＝1+h2+k2；m为电推进geo卫星的质量；t为电推进geo卫星的推力大小；g0为海平面处重力加速度；i
sp
为电推进系统的比冲；字符上一点表示求导。

技术总结
本发明公开了一种电推进GEO卫星自学习定轨方法，本发明对电推进GEO卫星通过加速度计对卫星的加速度进行测量，将卫星的控制加速度信息引入系统状态方程，减小轨道动力学模型误差，使轨道递推结果更精确。当有GNSS定轨信息时，本发明利用LSTM算法对GNSS导航的定轨误差进行估计，然后与UKF定轨结果进行二次融合，进一步提升卫星自主定轨精度。当GNSS不可用时，本发明在不额外增加星载敏感器的条件下，利用卫星上的太阳帆板输出功率信息，获取真实的地影进入时刻和离开时刻，从而构建地影信息。根据卫星飞行每圈的进出地影时间，构建卫星真经度角速度补偿量，从而减小轨道动力学模型误差。差。差。


技术研发人员：江秀强 吕文杰 罗跃龙 周创 季袁冬 钟苏川 孙国皓
受保护的技术使用者：四川大学
技术研发日：2023.04.12
技术公布日：2023/8/9