CN113641949A - 一种地球同步转移段轨道根数高精度拟合方法 - Google Patents

Abstract

本发明公开了一种地球同步转移段轨道根数高精度拟合方法，具体按照如下步骤实施：步骤1，在轨道拟合时间区间内定步长并计算k组轨道数据，步骤2，计算轨道平均半长轴

和平均轨道角速度n；步骤3，根据k组轨道数据以及轨道平均半长轴

和平均轨道角速度n列写预报方程；步骤4，列写系统状态方程和状态转移矩阵Φ；步骤5，对状态转移矩阵Φ进行SVD奇异值分解；步骤6，计算轨道系数拟合矩阵P，步骤7，将计算的轨道系数拟合矩阵P代入步骤3的预报方程中实现轨道预报。本发明一种地球同步转移段轨道根数高精度拟合方法，采用基函数拟和方法，将复杂的计算工作放在地面完成，实现了高精度的轨道预报。

Description

一种地球同步转移段轨道根数高精度拟合方法

技术领域

本发明属于航天器测量与控制方法技术领域，涉及一种地球同步转移段轨道根数高精度拟合方法。

背景技术

根据地球同步卫星配置要求，在转移轨道段，以星敏感器作为主份方案，实现卫星对地指向、点火姿态建立等姿态控制任务。星敏感器输出卫星相对地心惯性系(通常为J2000)的姿态信息，需要结合轨道信息才能确定当前卫星相对轨道系的三轴姿态。

地球同步卫星转移轨道段，在轨道高度大于20000km以上时，卫星进行对地指向、点火姿态建立等姿态控制任务需要轨道信息。卫星远地点变轨490N发动机工作过程中使用陀螺积分作为主份方案，陀螺积分必须提供精确的初始姿态信息作为陀螺积分初值，才能保证此后变轨过程中的姿态精度。因此，转移段需要高精度的轨道递推方法进行轨道计算。

由于转移段轨道所受的摄动力主要是地球非球形摄动、日月引力摄动、光压摄动。上述三者力的表现形式相对较为固定，经分析三种摄动力最终使得卫星轨道变化的形式主要表现为与时间相关的一次项和周期项。因此，星上采用基于该种表现形式基函数的数据拟合轨道递推方式。充分利用地面高精度轨道预报模型，获取精确的基函数系数，实现星上高精度轨道预报。

在转移段过程中，由于大偏心率的影响，卫星位置速度展开的基函数表达方式极为复杂，各参数之间存在较强的耦合关系，因此不宜使用位置速度作为转移段星上轨道预报模型的参数。

以往星上使用的轨道根数计算方法基本都是采用动力学积分的方法外推实现的，由于受制于星上计算机运算能力，星上采用的动力学模型相对简单，这就使得外推轨道根数不能满足高精度的需求。

发明内容

本发明的目的是提供一种地球同步转移段轨道根数高精度拟合方法，采用基函数拟和方法，将复杂的计算工作放在地面完成，实现了高精度的轨道预报。

本发明所采用的技术方案是，一种地球同步转移段轨道根数高精度拟合方法，采用轨道根数作为星上轨道预报参数的基函数模型，利用SVD奇异值分解对基函数的系数矩阵P进行拟合计算，具体按照如下步骤实施：

步骤1，在轨道拟合时间区间[t₀，t_f]内定步长Δt并计算k组轨道数据，轨道数据包括卫星轨道参数时间t_j＝t₀+(j-1)·Δt，j＝1，2，3，...，k，半长轴a_j，偏心率e_j，倾角i_j，升交点赤经Ω_j，近地点幅角ω_j，平近点角M_j；

步骤2，计算轨道平均半长轴

和平均轨道角速度n；

步骤3，根据k组轨道数据以及轨道平均半长轴

和平均轨道角速度n列写预报方程，在转移段采用轨道六根数作为星上轨道预报参数的基函数模型进行轨道预报；

步骤4，列写系统状态方程，系数矩阵P满足的线性方程组为：

ΦP＝F，其中，

f(t_j)＝(a_j，e_j，i_j，Ω_j，ω_j，M_j)^T，j＝1,2,3,…,k；Φ为状态转移矩阵；

步骤5，对状态转移矩阵Φ进行SVD奇异值分解，分解后为：Φ_m×l＝U_m×mD_m×lV_l×l，计算矩阵U_m×m、D_m×l、V_l×l，m为系统状态方程维数，l为系数矩阵P维数，其中，U_m×m为m维正交矩阵，V_l×l为l维数正交矩阵，

为分块矩阵，S_m×m为m维对角矩阵，对角线线元素为系数矩阵P的非零特征值；

步骤6，计算轨道系数拟合矩阵P，计算方法为：

步骤7，将计算的轨道系数拟合矩阵P代入步骤3的预报方程中实现轨道预报。

本发明的特征还在于，

步骤2的轨道平均半长轴

和平均轨道角速度n按照如下公式计算：

其中，μ为地球引力常数。

步骤3中的预报方程如下：

a_j＝P_a0+P_a1t_j+P_a1s sin nt_j+P_a1c cos nt_j+P_a2s sin 2nt_j+P_a2c cos 2nt_j+P_a3s sin3nt_j+P_a3c cos 3nt_j

e_j＝P_e0+P_e1t_j+P_e1s sin nt_j+P_e1c cos nt_j+P_e2s sin 2nt_j+P_e2c cos 2nt_j+P_e3s sin3nt_j+P_e3c cos 3nt_j

i_j＝P_i0+P_i1t_j+P_i1s sin nt_j+P_i1c cos nt_j+P_i2s sin 2nt_j+P_i2c cos 2nt_j+P_i3s sin3nt_j+P_i3c cos 3nt_j

Ω_j＝P_Ω0+P_Ω1t_j+P_Ω1s sin nt_j+P_Ω1c cos nt_j+P_Ω2s sin 2nt_j+P_Ω2c cos 2nt_j+P_Ω3ssin 3nt_j+P_Ω3c cos 3nt_j

ω_j＝P_w0+P_w1t_j+P_w1s sin nt_j+P_w1c cos nt_j+P_w2s sin 2nt_j+P_w2c cos 2nt_j+P_w3s sin3nt_j+P_w3c cos 3nt_j

M_j＝P_M0+P_M1t_j+P_M1s sin nt_j+P_M1c cos nt_j+P_M2s sin 2nt_j+P_M2c cos 2nt_j+P_M3s sin3nt_j+P_M3c cos 3nt_j

其中，P＝[P_h0,P_h1,P_h1s,P_h1c,P_h2s,P_h2c,P_h3s,P_h3c]^T，h代表a,e,i,Ω,ω,M，P为待求解的轨道拟合系数矩阵。

状态转移矩阵Φ为：

其中，Φ^T(t_j)＝(1，t_j，sinnt_j，cosnt_j，sin2nt_j，cos2nt_j，sin3nt_j，cos3nt_j)，j＝1,2,3,…，k。

步骤1中半长轴a_j单位为米，角度单位为弧度，当升交点赤经、近地点幅角、平近点角出现2π跳转时，需要将跳转后的数据加2π或减2π，使得原始数据连续。

本发明的有益效果是：

本发明采用的基函数拟和方法，将复杂的计算工作放在地面完成，充分利用地面高精度轨道预报模型，获取精确的基函数系数，避开了星上计算机运算能力的限制，实现了高精度的轨道预报。本发明的方法可靠性好、可操作性强、易推广和使用。

附图说明

图1是本发明一种地球同步转移段轨道根数高精度拟合方法的流程图。

具体实施方式

下面结合附图和具体实施方式对本发明进行详细说明。

本发明一种地球同步转移段轨道根数高精度拟合方法，采用轨道根数作为星上轨道预报参数的基函数模型，利用SVD奇异值分解对基函数的系数矩阵P进行拟合计算，如图1所示，具体按照如下步骤实施：

步骤1，在轨道拟合时间区间[t₀,t_f]内定步长Δt并计算k组轨道数据，轨道数据包括卫星轨道参数时间t_j＝t₀+(j-1)·Δt，j＝1,2,3,…,k，半长轴a_j，偏心率e_j，倾角i_j，升交点赤经Ω_j，近地点幅角ω_j，平近点角M_j，半长轴a_j单位为米，角度单位为弧度，当升交点赤经、近地点幅角、平近点角出现2π跳转时，需要将跳转后的数据加2π或减2π，使得原始数据连续；

步骤2，计算轨道平均半长轴

和平均轨道角速度n，按照如下公式计算：

其中，μ为地球引力常数；

步骤3，根据k组轨道数据以及轨道平均半长轴

和平均轨道角速度n列写预报方程，在转移段采用轨道六根数作为星上轨道预报参数的基函数模型进行轨道预报；预报方程如下：

a_j＝P_a0+P_a1t_j+P_a1s sin nt_j+P_a1c cos nt_j+P_a2s sin 2nt_j+P_a2c cos 2nt_j+P_a3s sin3nt_j+P_a3c cos 3nt_j

e_j＝P_e0+P_e1t_j+P_e1s sin nt_j+P_e1c cos nt_j+P_e2s sin 2nt_j+P_e2c cos 2nt_j+P_e3s sin3nt_j+P_e3c cos 3nt_j

i_j＝P_i0+P_i1t_j+P_i1s sin nt_j+P_i1c cos nt_j+P_i2s sin 2nt_j+P_i2c cos 2nt_j+P_i3s sin3nt_j+P_i3c cos 3nt_j

Ω_j＝P_Ω0+P_Ω1t_j+P_Ω1s sin nt_j+P_Ω1c cos nt_j+P_Ω2s sin 2nt_j+P_Ω2c cos 2nt_j+P_Ω3ssin 3nt_j+P_Ω3c cos 3nt_j

ω_j＝P_w0+P_w1t_j+P_w1s sin nt_j+P_w1c cos nt_j+P_w2s sin 2nt_j+P_w2c cos 2nt_j+P_w3s sin3nt_j+P_w3c cos 3nt_j

M_j＝P_M0+P_M1t_j+P_M1s sin nt_j+P_M1c cos nt_j+P_M2s sin 2nt_j+P_M2c cos 2nt_j+P_M3s sin3nt_j+P_M3c cos 3nt_j

其中，P＝[P_h0,P_h1,P_h1s,P_h1c,P_h2s,P_h2c,P_h3s,P_h3c]^T，h代表a,e,i,Ω,ω,M，P为待求解的轨道拟合系数矩阵；

步骤4，列写系统状态方程，系数矩阵P满足的线性方程组为：

ΦP＝F，其中，

f(t_j)＝(a_j，e_j，i_j，Ω_j，ω_j，M_j)^T，j＝1,2,3,…,k；Φ为状态转移矩阵；状态转移矩阵Φ为：

其中，Φ^T(t_j)＝(1，t_j，sinnt_j，cosnt_j，sin2nt_j，cos2nt_j，sin3nt_j，cos3nt_j)，j＝1,2,3,…，k；

步骤5，对状态转移矩阵Φ进行SVD奇异值分解，Φ_m×l＝U_m×mD_m×lV_l×l，计算矩阵U_m×m、D_m×l、V_l×l，m为系统状态方程维数，l为系数矩阵P维数，其中，U_m×m为m维正交矩阵，V_l×l为l维数正交矩阵，

为分块矩阵，S_m×m为m维对角矩阵，对角线线元素为系数矩阵P的非零特征值；

步骤6，计算轨道系数拟合矩阵P，计算方法为：

步骤7，将计算的轨道系数拟合矩阵P代入步骤3的预报方程中实现轨道预报。

实施例

本实施例采用的航天器为同步卫星，该同步卫星运行在36000km高度的地球同步轨道上，主要任务是进行通信技术试验。

如图1所示，具体步骤如下：

步骤1：在轨道拟合时间区间[t₀,t_f]内定步长Δt计算k组轨道数据，数据包括卫星轨道参数时间t_j＝t₀+(j-1)·Δt,j＝1,2,L k，半长轴a_j，偏心率e_j，倾角i_j，升交点赤经Ω_j，近地点幅角ω_j，平近点角M_j；

步骤2：计算轨道平半长轴

和平均轨道角速度n，其求法如下：

轨道平半长轴

平均轨道角速度

其中μ为地球引力常数；

步骤3：列写预报方程，在转移段采用轨道六根数作为星上轨道预报参数的基函数模型进行轨道预报，预报方程如下：

a_j＝P_a0+P_a1t_j+P_a1s sin nt_j+P_a1c cos nt_j+P_a2s sin 2nt_j+P_a2ccos 2nt_j+P_a3s sin3nt_j+P_a3ccos 3nt_j

e_j＝P_e0+P_e1t_j+P_e1s sin nt_j+P_e1ccos nt_j+P_e2s sin 2nt_j+P_e2ccos 2nt_j+P_e3s sin3nt_j+P_e3ccos 3nt_j

i_j＝P_i0+P_i1t_j+P_i1s sin nt_j+P_i1ccos nt_j+P_i2s sin 2nt_j+P_i2c cos 2nt_j+P_i3s sin3nt_j+P_i3ccos 3nt_j

Ω_j＝P_Ω0+P_Ω1t_j+P_Ω1s sin nt_j+P_Ω1ccos nt_j+P_Ω2s sin 2nt_j+P_Ω2ccos 2nt_j+P_Ω3s sin3nt_j+P_Ω3ccos 3nt_j

ω_j＝P_w0+P_w1t_j+P_w1s sin nt_j+P_w1ccos nt_j+P_w2s sin 2nt_j+P_w2ccos 2nt_j+P_w3s sin3nt_j+P_w3ccos 3nt_j

M_j＝P_M0+P_M1t_j+P_M1s sin nt_j+P_M1ccos nt_j+P_M2s sin 2nt_j+P_M2ccos 2nt_j+P_M3s sin3nt_j+P_M3ccos 3nt_j

其中，P＝[P_h0,P_h1,P_h1s,P_h1c,P_h2s,P_h2c,P_h3s,P_h3c]^T，h代表轨道六根数a,e,i,Ω,ω,M，P为待求解的轨道拟合系数矩阵。半长轴单位为米，角度单位为弧度，当升交点赤经、近地点幅角、平近点角出现2π跳转时，需要将跳转后的数据加2π或减2π，使得原始数据连续；

步骤4：计算状态转移矩阵Φ，

其中，Φ^T(t_j)＝(1，t_j，sinnt_j，cosnt_j，sin2nt_j，cos2nt_j，sin3nt_j，cos3nt_j)，j＝1,2,3,…，k；

步骤5：列写系统状态方程，系数矩阵P满足的线性方程组为：

ΦP＝F，其中，

f(t_j)＝(a_j，e_j，i_j，Ω_j，ω_j，M_j)^T，j＝1,2,3,…,k；

步骤6：对状态转移矩阵Φ进行SVD奇异值分解，计算矩阵U_m×m、D_m×l、V_l×l,m为系统状态方程维数，l为系数矩阵P维数，Φ_m×l＝U_m×mD_m×lV_l×l，其中U_m×m为m维正交矩阵，V_l×l为l维数正交矩阵，

为分块矩阵，S_m×m为m维对角矩阵，对角线线元素为系数矩阵P的非零特征值；

步骤7：计算轨道系数拟合矩阵P，计算方法为：

将计算的轨道系数拟合矩阵P代入步骤3的预报方程中实现轨道预报。

上述同步卫星地球同步转移段轨道根数高精度拟合方法计算结果如下：

从拟和结果中的回代原始数据检核信息可以看出，星上采用这种拟和方法进行轨道预报产生的轨道根数与地面预报值相比，对应的误差均方根(RMS)与平均误差(MEAN)都很小，能够显著提高星上轨道预报精度。

Claims (5)

1.一种地球同步转移段轨道根数高精度拟合方法，其特征在于，采用轨道根数作为星上轨道预报参数的基函数模型，利用SVD奇异值分解对基函数的系数矩阵P进行拟合计算，具体按照如下步骤实施：

步骤1，在轨道拟合时间区间[t₀,t_f]内定步长Δt并计算k组轨道数据，所述轨道数据包括卫星轨道参数时间t_j＝t₀+(j-1)·Δt，j＝1,2,3,…,k，半长轴a_j，偏心率e_j，倾角i_j，升交点赤经Ω_j，近地点幅角ω_j，平近点角M_j；

步骤2，计算轨道平均半长轴

和平均轨道角速度n；

步骤3，根据k组轨道数据以及轨道平均半长轴

和平均轨道角速度n列写预报方程，在转移段采用轨道六根数作为星上轨道预报参数的基函数模型进行轨道预报；

步骤4，列写系统状态方程，系数矩阵P满足的线性方程组为：

ΦP＝F，其中，

f(t_j)＝(a_j，e_j，i_j，Ω_j，ω_j，M_j)^T，j＝1,2,3,…,k；Φ为状态转移矩阵；

步骤5，对状态转移矩阵Φ进行SVD奇异值分解，Φ_m×l＝U_m×mD_m×lV_l×l，计算矩阵U_m×m、D_m×l、V_l×l，m为系统状态方程维数，l为系数矩阵P维数，其中，U_m×m为m维正交矩阵，V_l×l为l维数正交矩阵，

为分块矩阵，S_m×m为m维对角矩阵，对角线线元素为系数矩阵P的非零特征值；

步骤6，计算轨道系数拟合矩阵P，计算方法为：

步骤7，将计算的轨道系数拟合矩阵P代入步骤3的预报方程中实现轨道预报。

2.根据权利要求1所述的一种地球同步转移段轨道根数高精度拟合方法，其特征在于，所述步骤2所述的轨道平均半长轴

和平均轨道角速度n按照如下公式计算：

其中，μ为地球引力常数。

3.根据权利要求2所述的一种地球同步转移段轨道根数高精度拟合方法，其特征在于，所述步骤3中的预报方程如下：

a_j＝P_a0+P_a1t_j+P_a1ssinnt_j+P_a1ccosnt_j+P_a2ssin2nt_j+P_a2ccos2nt_j+P_a3ssin3nt_j+P_a3ccos3nt_j

e_j＝P_e0+P_e1t_j+P_e1ssinnt_j+P_e1ccosnt_j+P_e2ssin2nt_j+P_e2ccos2nt_j+P_e3ssin3nt_j+P_e3ccos3nt_j

i_j＝P_i0+P_i1t_j+P_i1ssinnt_j+P_i1ccosnt_j+P_i2ssin2nt_j+P_i2ccos2nt_j+P_i3ssin3nt_j+P_i3ccos3nt_j

Ω_j＝P_Ω0+P_Ω1t_j+P_Ω1ssinnt_j+P_Ω1ccosnt_j+P_Ω2ssin2nt_j+P_Ω2ccos2nt_j+P_Ω3ssin3nt_j+P_{Ω 3c}cos3nt_j

ω_j＝P_w0+P_w1t_j+P_w1ssinnt_j+P_w1ccosnt_j+P_w2ssin2nt_j+P_w2ccos2nt_j+P_w3ssin3nt_j+P_w3ccos3nt_j

M_j＝P_M0+P_M1t_j+P_M1ssinnt_j+P_M1ccosnt_j+P_M2ssin2nt_j+P_M2ccos2nt_j+P_M3ssin3nt_j+P_M3ccos3nt_j

其中，P＝[P_h0,P_h1,P_h1s,P_h1c,P_h2s,P_h2c,P_h3s,P_h3c]^T，h代表a,e,i,Ω,ω,M，P为待求解的轨道拟合系数矩阵。

4.根据权利要求3所述的一种地球同步转移段轨道根数高精度拟合方法，其特征在于，所述状态转移矩阵Φ为：

其中，Φ^T(t_j)＝(1，t_j，sinnt_j，cosnt_j，sin2nt_j，cos2nt_j，sin3nt_j，cos3nt_j)，j＝1，2，3，…，k。

5.根据权利要求1所述的一种地球同步转移段轨道根数高精度拟合方法，其特征在于，所述步骤1中半长轴a_j单位为米，角度单位为弧度，当升交点赤经、近地点幅角、平近点角出现2π跳转时，需要将跳转后的数据加2π或减2π，使得原始数据连续。

Images