CN106092099B - 航天器相对位置增量定轨方法 - Google Patents

Abstract

航天器相对位置增量定轨方法，属于航天器轨道确定技术领域。本发明是为了解决现有基于测量历元差分的定轨方法中，需要采用其他定轨方法获得初始位置信息的问题。本发明所述的航天器相对位置增量定轨方法，首先，确定航天器在三维空间中运动的轨道平面；然后，将三维相对位置矢量转换为二维平面内航天器相对位置坐标，在该二维平面内采用最小二乘原理进行航天器运动轨道形状的拟合，采用拉格朗日乘子法求解拟合矩阵的特征值；最后，根据拟合矩阵的特征值判断航天器的轨道类型，并根据中心天体的位置，获得航天器在惯性参考系下的绝对位置信息，实现航天器的初轨确定。本发明适用于对航天器实现定轨。

Description

航天器相对位置增量定轨方法

技术领域

本发明属于航天器轨道确定技术领域。

背景技术

航天器轨道确定的任务就是应用少量的观测数据，快速确定航天器两体问题下的轨道要素。观测数据由于测站数目的不同，测量量(如：位置、速度、角度)属性的差别，相应地定轨方法也有所不同。其中位置观测量可通过全球卫星导航系统或地面雷达获得，测量得到的位置信息均可通过简单的转换获得航天器在地球惯性坐标系下的位置矢量或方向矢量，采用三个地心位置矢量可由吉布斯方法定轨，采用两个位置矢量及他们之间的时间间隔定轨的兰伯特问题。在脉冲星导航系统及全球卫星导航系统中，利用观测数据的历元差分方法可有效消除测量部分的公共误差提供更精确的测量信息，然而在基于测量历元差分的导航方法中，初始位置信息的估计往往先采用地心位置矢量通过其他定轨方法给出。为此使得基于历元差分相对位置信息导航算法无法独立应用，需依靠其他定轨算法提供可靠的初始位置估计后方可进行。

发明内容

本发明是为了解决现有基于测量历元差分的定轨方法中，需要采用其他定轨方法获得初始位置信息的问题，现提供航天器相对位置增量定轨方法。

航天器相对位置增量定轨方法，该方法为：

首先，确定航天器在三维空间中运动的轨道平面；

然后，将三维相对位置矢量转换为二维平面内航天器相对位置坐标，在该二维平面内采用最小二乘原理进行航天器运动轨道形状的拟合，并采用拉格朗日乘子法求解拟合矩阵的特征值；

最后，根据拟合矩阵的特征值判断航天器的轨道类型，若航天器轨道为椭圆轨道，则中心天体位于椭圆的一个焦点处，结合开普勒时间方程与测量相对应的时间序列确定中心天体的位置；若航天器轨道为圆轨道，则中心天体在轨道的中心；根据中心天体的位置，获得航天器在惯性参考系下的绝对位置信息，实现航天器的初轨确定。

本发明所述的航天器相对位置增量定轨方法，是一种直接使用历元差分测量的相对距离增量信息进行定轨的方法，仅依靠航天器相对于空间任意未知位置的相对位置增量及相应的观测时间。该方法使得历元差分测量信息在导航应用时不依赖与其他测量信息，提高了基于历元差分测量信息的导航系统的自主性，本发明定轨方法简单，易于实现。

附图说明

图1为基于历元差分的航天器相对位置示意图；

图2为轨道平面的各参考系示意图；

图3为本发明所述的航天器相对位置增量定轨方法的流程图。

具体实施方式

具体实施方式一：参照图3具体说明本实施方式，本实施方式所述的航天器相对位置增量定轨方法，该方法为：

首先，确定航天器在三维空间中运动的轨道平面；

然后，将三维相对位置矢量转换为二维平面内航天器相对位置坐标，在该二维平面内采用最小二乘原理进行航天器运动轨道形状的拟合，并采用拉格朗日乘子法求解拟合矩阵的特征值；

最后，根据拟合矩阵的特征值判断航天器的轨道类型，若航天器轨道为椭圆轨道，则中心天体位于椭圆的一个焦点处，结合开普勒时间方程与测量相对应的时间序列确定中心天体的位置；若航天器轨道为圆轨道，则中心天体在轨道的中心；根据中心天体的位置，获得航天器在惯性参考系下的绝对位置信息，实现航天器的初轨确定。

具体实施方式二：参照图1具体说明本实施方式，本实施方式是对具体实施方式一所述的航天器相对位置增量定轨方法作进一步说明，本实施方式中，确定航天器在三维空间中运动的轨道平面的具体方法为：

测量航天器在地心惯性坐标系下历元时刻t_j与t_j-1之间航天器移动的相对位置矢量Δr(t_j-1,t_j)，并根据该相对位置矢量Δr(t_j-1,t_j)获得航天器的初始位置r(t₀)，其中j＝1,2,3,...；

根据相对位置矢量Δr(t_j-1,t_j)，采用最小二乘原理确定单位角动量矢量h；

其中，h_x、h_y、h_z分别为单位角动量矢量h在空间三维参考系中的x轴、y轴、z轴的分量，△r_x(t_j-1,t_j)、△r_y(t_j-1,t_j)、△r_z(t_j-1,t_j)分别为相对位置矢量Δr(t_j-1,t_j)在三维参考系中x轴、y轴、z轴的分量，

根据下式利用单位角动量矢量h获得轨道倾角i和升交点赤经Ω，

i＝arccos(h_z) (2)

利用轨道倾角i和升交点赤经Ω确定三维轨道平面坐标。

图1描述了相对位置测量信息，r为航天器在地心惯性坐标系下相对于地心的位置矢量，Δr(t_j-1,t_j)和Δr(t_j,t_j+1)分别为历元时刻(t_j-1,t_j)和(t_j,t_j+1)间航天器移动的相对位置矢量，定轨就是根据一系列的测量信息Δr(t_j-1,t_j)及t_j(j＝1,2,3,...)获得航天器的初始位置r(t₀)及速度或轨道六根数。

设相对位置Δr为地心惯性J2000参考系下的矢量。根据向量间关系：任意两个矢量差乘即可获得轨道平面的法向量，即单位角动量矢量，这里考虑测量噪声的存在，利用所有观测量，采用最小二乘原理确定单位角动量矢量h。上述式(1)中可通过求解目标函数：

来获得约束条件G(h_x,h_y,h_z)下的最小值。并根据非线性优化理论，选取适当的优化算法(如牛顿法、最速下降法等)获得单位角动量矢量h。

具体实施方式三：本实施方式是对具体实施方式二所述的航天器相对位置增量定轨方法作进一步说明，本实施方式中，将三维相对位置矢量转换为二维平面内航天器相对位置坐标的具体方法为：

设参考系P_xyz以航天器的初始位置为原点，z轴垂直于轨道平面，x,y轴由地球J2000参考系根据坐标转换矩阵Q_Xx获得；

根据轨道倾角i和升交点赤经Ω获得航天器相对位置在参考系P_xyz中的矢量x_k：

x_k＝Q_Xxm_k (4)

式中，为k时刻航天器位置相对于初始位置的矢量，Q_Xx为坐标转换矩阵，具体公式如下：

利用矢量x_k在x,y轴上的投影获得航天器各点在二维平面内的位置坐标(x_k,y_k)。

图2给出了在轨道平面内各参考系之间的关系，其中参考系P_xyz以航天器未知的初始位置为原点，z轴垂直于轨道平面，x,y轴由地球J2000参考系根据坐标转换矩阵Q_Xx获得。参考系O_x'y'z'以椭圆轨道中心为原点，z轴垂直于轨道平面，x,y轴分别沿椭圆长轴和短轴方向，满足右手定则。参考系O'_xyz以地球中心为原点，为轨道平面内的近焦点坐标系。

具体实施方式四：本实施方式是对具体实施方式一所述的航天器相对位置增量定轨方法作进一步说明，本实施方式中，采用最小二乘原理进行航天器运动轨道形状拟合的具体方法为：

设平面内任意轨道的表达式为：

其中，U＝[A,B,C,D,E,M]，X＝[x²,xy,y²,x,y,1]，A、B、C、D、E、M均为多项式各项的系数，x和y分别为二维平面参考系P_xyz中轨道上各点坐标参数，

根据最小二乘原理将公式(6)转换为以下形式：

式中可由矢量x_k在x,y轴上的投影获得航天器各点在二维平面内的位置坐标(x_k,y_k)计算获得。

本实施方式公式(6)中含有的：

4AC-B²＝1 (7)

公式(7)是为了避免全零解及抛物线、双曲线情况，所增加的约束条件。

平面椭圆拟合问题可通过Hough变换、Kalman滤波及最小二乘法等求解，考虑到计算的复杂度与鲁棒性，本实施方式中选择最小二乘法求解椭圆拟合问题。

具体实施方式五：本实施方式是对具体实施方式四所述的航天器相对位置增量定轨方法作进一步说明，本实施方式中，采用拉格朗日乘子法求解拟合矩阵的特征值的具体方法为：

将公式(8)改写成矩阵形式：

其中，

采用拉格朗日乘子求解式(9)，则有：

上式中W＝V^TV记为散布矩阵，采用拉格朗日乘子求解矩阵W特征值(λ_i,U_i)，其中i＝1,2,...,6，λ_i为特征值，U_i为特征向量，最小的正实数解即为轨道拟合问题的解，即获得参数U的估计值。

Fitzgibbon(菲茨吉本)给出了基于最小二乘原理的轨道拟合问题的解，在此基础上，Halir(海立)进一步给出了快速非迭代的椭圆拟合算法，Nikolai(尼古拉)在Matlab的CodeExchange(程序交流论坛)中公布了稳定可靠的拟合程序。上述对式(9)的求解就基于Nikolai的程序。将式(6)与椭圆几何参数模型相比较，从而可根据计算得到的参数A获得图2中参数x₀，y₀，a，b及α。

具体实施方式六：本实施方式是对具体实施方式五所述的航天器相对位置增量定轨方法作进一步说明，本实施方式中，根据拟合矩阵的特征值判断航天器的轨道类型的具体方法为：

在无噪声的情况下参数A满足如下条件：

考虑噪声的影响，设置门限eps用以区分圆轨道情况，

将参数A,B,C代入下式，

(A/C-1)²+B²<eps (22)

若结果满足式(22)，则轨道为圆轨道，

若结果不满足式(22)，则轨道为椭圆轨道。

本实施方式中，由于圆轨道可视为偏心率为零时的一种特殊的椭圆轨道，因此才有式(21)的条件。

具体实施方式七：本实施方式是对具体实施方式六所述的航天器相对位置增量定轨方法作进一步说明，本实施方式中，航天器轨道为椭圆轨道时，航天器初轨确定的具体方法为：

根据以下椭圆的几何参数方程确定椭圆在参考系P_xyz下的几何参数α、x₀、y0、a和b：

利用下式获得偏心率e：

设参考系O_x'y'z'以椭圆轨道中心为原点，z'轴垂直于轨道平面，x'，y'轴分别与椭圆长轴方向和短轴方向同向，将椭圆轨道上各点坐标转换到参考系O_x'y'z'下，满足右手定则，有：

p_k＝Q_xx'(x_k+S_xx') (14)

其中，p_k为航天器相对于参考系O_x'y'z'原点的位置矢量，Q_xx'为坐标转换矩阵，S_xx'为坐标平移向量，解得：

设地球位于椭圆轨道的右焦点处，根据位置矢量p_k在参考系O_x'y'z'下各方向分量获得轨道各点的真近点角θ_k：

将真近点角θ_k代入到开普勒时间方程中，获得平近点角M_e(k)：

根据平近点角与时间的关系获得相应位置相对于近地点时刻的时间t'_k，其中以近地点时刻为时间起点，若t'_k与测量量记录的时间t_k具有相同的变化趋势，则地球位于椭圆轨道的右焦点处，且仅有常值偏差；若t'_k与t_k变化趋势不符，则地球位于椭圆轨道的左焦点处，

根据地球位于椭圆轨道的位置获得近地点幅角ω和真近点角θ₀：

其中，q₀＝x₀sinω-y₀cosω，(p₀,q₀)为航天器初始位置在近焦点坐标系O'_xyz下的坐标；

即能够获得航天器在地心J2000参考系下相对地心的绝对位置矢量r_k，

其中，为坐标转换矩阵的转置矩阵，为坐标平移向量，

具体实施方式八：本实施方式是对具体实施方式六所述的航天器相对位置增量定轨方法作进一步说明，本实施方式中，航天器轨道为圆轨道时，航天器初轨确定的具体方法为：

利用下式确定圆轨道的半径R_c及中心位置坐标O_c：

利用圆轨道在平面参考系下的中心位置坐标O_c，确定航天器各点在地心J2000参考系下航天器的绝对位置矢量r_k：

根据角动量矢量h获得位置矢量r′_k的速度信息v_k：

其中，μ为中心天体的引力常数。

本实施方式中所述的引力常数μ，对于地球一般取μ＝398600km³/s²。本实施方式对圆轨道确定了轨道各点的绝对位置及速度信息，即完成轨道的确定。

本发明讨论了根据相对位置矢量确定航天器椭圆、圆轨道的方法。下面对上述定轨算法分别采用椭圆轨道和圆轨道进行检验。设观测的相对位置矢量的信噪比为20dB，每15s获得一次观测信息，总观测数据为200时，分别进行1000次蒙特卡洛仿真，结果如表1和表2所示。

表1 ST5-A卫星椭圆轨道定轨结果

表2 PICARD卫星圆轨道定轨结果

仿真结果验证了本发明所述方法的可行性。

Claims (8)

1.航天器相对位置增量定轨方法，其特征在于，该方法为：

首先，确定航天器在三维空间中运动的轨道平面；

然后，将三维相对位置矢量转换为二维平面内航天器相对位置坐标，在该二维平面内采用最小二乘原理进行航天器运动轨道形状的拟合，并采用拉格朗日乘子法求解拟合矩阵的特征值；

最后，根据拟合矩阵的特征值判断航天器的轨道类型，若航天器轨道为椭圆轨道，则中心天体位于椭圆的一个焦点处，结合开普勒时间方程与测量相对应的时间序列确定中心天体的位置；若航天器轨道为圆轨道，则中心天体在轨道的中心；根据中心天体的位置，获得航天器在惯性参考系下的绝对位置信息，实现航天器的初轨确定。

2.根据权利要求1所述的航天器相对位置增量定轨方法，其特征在于，确定航天器在三维空间中运动的轨道平面的具体方法为：

测量航天器在地心惯性坐标系下历元时刻t_j与t_j-1之间航天器移动的相对位置矢量Vr(t_j-1,t_j)，并根据该相对位置矢量Vr(t_j-1,t_j)获得航天器的初始位置r(t₀)，其中j＝1,2,3,...；

根据相对位置矢量Vr(t_j-1,t_j)，采用最小二乘原理确定单位角动量矢量h；

其中，h_x、h_y、h_z分别为单位角动量矢量h在空间三维参考系中的x轴、y轴、z轴的分量，△r_x(t_j-1,t_j)、△r_y(t_j-1,t_j)、△r_z(t_j-1,t_j)分别为相对位置矢量Vr(t_j-1,t_j)在三维参考系中x轴、y轴、z轴的分量，

根据下式利用单位角动量矢量h获得轨道倾角i和升交点赤经Ω，

i＝arccos(h_z) (2)

利用轨道倾角i和升交点赤经Ω确定三维轨道平面坐标。

3.根据权利要求2所述的航天器相对位置增量定轨方法，其特征在于，将三维相对位置矢量转换为二维平面内航天器相对位置坐标的具体方法为：

设参考系P_xyz以航天器的初始位置为原点，z轴垂直于轨道平面，x,y轴由地球J2000参考系根据参考系P_xyz的坐标转换矩阵Q_Xx获得；

根据轨道倾角i和升交点赤经Ω获得航天器相对位置在参考系P_xyz中的矢量x_k：

x_k＝Q_Xxm_k (4)

式中，为k时刻航天器位置相对于初始位置的矢量，Q_Xx为参考系P_xyz的坐标转换矩阵，具体公式如下：

利用矢量x_k在x,y轴上的投影获得航天器各点在二维平面内的位置坐标(x_k,y_k)。

4.根据权利要求1所述的航天器相对位置增量定轨方法，其特征在于，采用最小二乘原理进行航天器运动轨道形状拟合的具体方法为：

设平面内任意轨道的表达式为：

其中，U＝[A,B,C,D,E,M]，X＝[x²,xy,y²,x,y,1]，A、B、C、D、E、M均为多项式各项的系数，x和y分别为二维平面参考系P_xyz中轨道上各点坐标参数，

根据最小二乘原理将公式(6)转换为以下形式：

式中，

5.根据权利要求4所述的航天器相对位置增量定轨方法，其特征在于，采用拉格朗日乘子法求解拟合矩阵的特征值的具体方法为：

将公式(8)改写成矩阵形式：

其中，

采用拉格朗日乘子求解式(9)，则有：

上式中W＝V^TV记为散布矩阵，采用拉格朗日乘子求解矩阵W特征值(λ₁,U_i)，其中i＝1,2,...,6，λ₁为特征值，U_i为特征向量，最小的正实数解即为轨道拟合问题的解，即获得参数U的估计值。

6.根据权利要求5所述的航天器相对位置增量定轨方法，其特征在于，根据拟合矩阵的特征值判断航天器的轨道类型的具体方法为：

在无噪声的情况下参数A满足如下条件：

考虑噪声的影响，设置门限eps用以区分圆轨道情况，

将参数A,B,C代入下式，

(A/C-1)²+B²<eps (22)

若结果满足式(22)，则轨道为圆轨道，

若结果不满足式(22)，则轨道为椭圆轨道。

7.根据权利要求6所述的航天器相对位置增量定轨方法，其特征在于，航天器轨道为椭圆轨道时，航天器初轨确定的具体方法为：

根据以下椭圆的几何参数方程确定椭圆在参考系P_xyz下的几何参数α、x₀、y₀、a和b：

利用下式获得偏心率e：

设参考系O_x'y'z'以椭圆轨道中心为原点，z'轴垂直于轨道平面，x'，y'轴分别与椭圆长轴方向和短轴方向同向，将椭圆轨道上各点坐标转换到参考系O_x'y'z'下，满足右手定则，有：

p_k＝Q_xx'(x_k+S_xx') (14)

其中，p_k为航天器相对于参考系O_x'y'z'原点的位置矢量，Q_xx'为参考系O_x'y'z'的坐标转换矩阵，S_xx'为坐标平移向量，解得：

设地球位于椭圆轨道的右焦点处，根据位置矢量p_k在参考系O_x'y'z'下各方向分量获得轨道各点的真近点角θ_k：

将真近点角θ_k代入到开普勒时间方程中，获得平近点角M_e(k)：

μ为中心天体的引力常数，

根据平近点角与时间的关系获得相应位置相对于近地点时刻的时间t'_k，其中以近地点时刻为时间起点，若t'_k与测量量记录的时间t_k具有相同的变化趋势，则地球位于椭圆轨道的右焦点处，且仅有常值偏差；若t'_k与t_k变化趋势不符，则地球位于椭圆轨道的左焦点处，

根据地球位于椭圆轨道的位置获得近地点幅角ω和真近点角θ₀：

其中，q₀＝x₀sin ω-y₀cos ω，(p₀,q₀)为航天器初始位置在近焦点坐标系下的坐标；

即能够获得航天器在地心J2000参考系下相对地心的绝对位置矢量r_k，

其中，为坐标系的坐标转换矩阵的转置矩阵，为坐标平移向量，

8.根据权利要求6所述的航天器相对位置增量定轨方法，其特征在于，航天器轨道为圆轨道时，航天器初轨确定的具体方法为：

利用下式确定圆轨道的半径R_c及中心位置坐标O_c：

利用圆轨道在平面参考系下的中心位置坐标O_c，确定圆轨道时航天器各点在地心J2000参考系下航天器的绝对位置矢量r'_k：

根据角动量矢量h获得位置矢量r'_k的速度信息v_k：

其中，μ为中心天体的引力常数。