CN113128032A - 一种基于轨道解析摄动解的交点时刻及位置求解算法 - Google Patents

Abstract

本发明公开了一种基于轨道解析摄动解的交点时刻及位置求解算法，计算交点时刻采用仅考虑地球扁率J2‑J4项的长期项影响的摄动模型，忽略周期项影响以降低计算量，利用摄动方程计算假定时刻双方航天器轨道面方位，并通过几何关系计算得到交点时刻，将该交点时刻与假定时刻对比，若误差较大，则将该交点时刻作为新的假定时刻以重复计算交点时刻，直至迭代收敛，最后利用同时带有长期项、长周期和短周期振荡的密切轨道要素摄动方程，求得交点时刻对方航天器的位置。上述方法不仅可以精确计算交点时刻及位置，而且算法简单易于实施，节约计算时间，节省计算资源和内存，特别适用于计算资源紧张的星载计算机，具有良好的推广前景。

Description

一种基于轨道解析摄动解的交点时刻及位置求解算法

技术领域

本发明涉及航天器轨道动力学技术领域，尤其涉及一种基于轨道解析摄动解的交点时刻及位置求解算法。

背景技术

为航天器交会任务做准备，需要提前获得对方航天器穿越己方航天器轨道平面的时刻和对应时刻对方航天器的位置。近地航天器(对方航天器和己方航天器都是近地航天器)轨道主要受地球扁率J2到J4项的摄动影响，除了会使平均轨道要素受长期项的摄动影响产生匀速变化外，还会导致航天器密切轨道要素产生长周期和短周期振荡，密切轨道要素包括轨道半长轴a、偏心率e、轨道倾角i、升交点赤经Ω、近地点幅角ω和平近点角M。振荡幅值与轨道半长轴a、轨道倾角i有关，振荡相位与近地点幅角ω、纬度幅角u有关。在计算交点时刻时，如果使用考虑长周期和短周期振荡的密切轨道要素，则密切轨道要素会包含复杂的时间函数，使得交点时刻的计算极其复杂，势必会带来计算量的激增。因此，在平均轨道要素的精度无法满足需求时，才会考虑带有长周期和短周期振荡的密切轨道要素。

目前已有的交点时刻计算方法可大致分为两类。一类是使用考虑周期项的密切轨道要素进行计算，其精度较高，但计算量过于庞大，不利于工程应用；另一类是使用简单的无摄动或仅考虑J2项摄动效果的模型，其计算量较小，但精度较低。现存的算法或精度较低，或需要占用大量计算资源而无法应用于星载计算机，在实际航天任务的应用中都存在不足。

发明内容

有鉴于此，本发明提供了一种基于轨道解析摄动解的交点时刻及位置求解算法，用以在已知双方航天器初始轨道要素的基础上快速求解对方航天器与己方航天器轨道平面的交点时刻以及此时对方航天器的精确位置。

本发明提供的一种基于轨道解析摄动解的交点时刻及位置求解算法，包括如下步骤：

S1：根据双方航天器的平均轨道要素初值，使用仅考虑地球扁率J2-J4项的长期项影响的摄动模型，计算双方航天器的升交点赤经的变化速率、对方航天器的近地点幅角的变化速率和对方航天器的平近点角的变化速率；

S2：根据双方航天器的升交点赤经的变化速率，计算假定时刻双方航天器的轨道面方位；

S3：根据对方航天器的平近点角的变化速率和近地点幅角的变化速率，以及假定时刻双方航天器的轨道面方位，计算交点时刻；

S4：判断交点时刻与假定时刻的误差是否大于阈值；若是，则将交点时刻作为新的假定时刻，返回步骤S2，重复执行步骤S2～步骤S4；若否，则执行步骤S5；

S5：利用带有长周期振荡和短周期振荡的密切轨道要素摄动方程，计算交点时刻对方航天器的位置。

在一种可能的实现方式中，在本发明提供的上述基于轨道解析摄动解的交点时刻及位置求解算法中，步骤S1，根据双方航天器的平均轨道要素初值，使用仅考虑地球扁率J2-J4项的长期项影响的摄动模型，计算双方航天器的升交点赤经的变化速率、对方航天器的近地点幅角的变化速率和对方航天器的平近点角的变化速率，具体包括：

仅考虑地球扁率J2-J4项的长期项影响，双方航天器的升交点赤经Ω_k的变化速率

对方航天器的近地点幅角ω_B的变化速率

和对方航天器的平近点角M_B的变化速率

为：

其中，μ表示地球引力常数，R_E表示地球半径，J₂、和J₄表示地球扁率项；i_k表示双方航天器的轨道倾角，a_k表示双方航天器的轨道半长轴，e_k表示双方航天器的离心率，i_k、a_k和e_k在长期项作用下保持不变；k＝A,B，A表示己方航天器，B表示对方航天器。

在一种可能的实现方式中，在本发明提供的上述基于轨道解析摄动解的交点时刻及位置求解算法中，步骤S2，根据双方航天器的升交点赤经的变化速率，计算假定时刻双方航天器的轨道面方位，具体包括：

由双方航天器的升交点赤经的变化速率

求得假定时刻双方航天器的升交点赤经：

其中，Ω_k0表示双方航天器初始的升交点赤经，t表示任意假定时刻；双方航天器的轨道倾角保持不变：

i_k＝i_k0 (5)

其中，i_k0表示双方航天器初始的轨道倾角；由假定时刻双方航天器的升交点赤经Ω_k和轨道倾角i_k表示假定时刻双方航天器的轨道面方位。

在一种可能的实现方式中，在本发明提供的上述基于轨道解析摄动解的交点时刻及位置求解算法中，步骤S3，根据对方航天器的平近点角的变化速率和近地点幅角的变化速率，以及假定时刻双方航天器的轨道面方位，计算交点时刻，具体包括如下步骤：

S31：假定交点时刻为t′_n，n＝1，取t′₁＝t，由对方航天器的近地点幅角的变化速率

计算假定交点时刻对方航天器的近地点幅角：

其中，ω_B0表对方航天器初始的近地点幅角；

S32：计算假定交点时刻对方航天器的真近点角：

θ′_B＝u_B-ω′_B (7)

其中，u_B表示对方航天器穿过己方航天器轨道平面时的纬度幅角；u_B的计算公式如下：

其中，

其中，

表示对方航天器轨道节线的单位矢量；

表示双方航天器轨道面的交线单位矢量；

表示己方航天器的轨道法线方向；

表示对方航天器的轨道法线方向；

S33：计算假定交点时刻对方航天器的偏近点角：

S34：计算假定交点时刻对方航天器的平近点角：

M′_B＝E′_B-e_BsinE′_B (14)

S35：更新后的交点时刻为：

其中，M_B0表示对方航天器初始的平近点角；

S36：判断更新后的交点时刻t′_i+1与假定交点时刻t′_n的误差是否大于阈值；若否，则输出更新后的交点时刻t′_n+1；若是，则将更新后的交点时刻t′_n+1作为新的假定时刻，返回步骤S31，n＝n+1，重复步骤S31～步骤S36。

在一种可能的实现方式中，在本发明提供的上述基于轨道解析摄动解的交点时刻及位置求解算法中，步骤S3中判断更新后的交点时刻t′_n+1与假定交点时刻t′_n的误差，以及步骤S4中判断交点时刻与假定时刻的误差，所依据的阈值的范围为10^-6s～10^-9s。

在一种可能的实现方式中，在本发明提供的上述基于轨道解析摄动解的交点时刻及位置求解算法中，步骤S5，利用带有长周期振荡和短周期振荡的密切轨道要素摄动方程，计算交点时刻对方航天器的位置，具体包括：

利用对方航天器的完整的J2-J4项的解析摄动解，包括长期项、长周期项和短周期项的结果，以及对方航天器的平均轨道要素初值，计算交点时刻t′对方航天器的密切轨道要素：

密切轨道要素(t′)＝长期项(t′)+长周期项(t′)+短周期项(t′)+平均轨道要素初值。

本发明提供的上述基于轨道解析摄动解的交点时刻及位置求解算法，计算交点时刻采用仅考虑地球扁率J2-J4项的长期项影响的摄动模型，忽略周期项影响以降低计算量，利用摄动方程计算假定时刻双方航天器轨道面方位，并通过几何关系计算得到交点时刻，将计算得到的交点时刻与假定时刻对比，若误差较大，则将计算得到的交点时刻作为新的假定时刻以重复计算交点时刻，直至迭代收敛，最后利用同时带有长期项、长周期和短周期振荡的密切轨道要素摄动方程，求得交点时刻对方航天器的位置。上述方法不仅可以精确计算交点时刻及位置，而且算法简单易于实施，节约计算时间，节省计算资源和内存，特别适用于计算资源紧张的星载计算机，具有良好的推广前景。

附图说明

图1为本发明实施例1中基于轨道解析摄动解的交点时刻及位置求解算法的流程示意图；

图2为定轨方式示意图；

图3为距离误差变化图；

图4为位置误差变化图；

图5为方法二和方法三分别与方法一的交点时刻的差值变化图；

图6为方法一和方法二的交点时刻真实距离变化图；

图7为方法一和方法三的交点时刻真实距离变化图；

图8为方法一和方法二的交点时刻真实夹角变化图；

图9为方法一和方法三的交点时刻真实夹角变化图；

图10为J2-J4模型与对比模型在一天内位置误差比较图；

图11为J2-J4模型在一天内位置误差变化图。

具体实施方式

下面将结合本发明实施方式中的附图，对本发明实施方式中的技术方案进行清楚、完整的描述，显然，所描述的实施方式仅仅是作为例示，并非用于限制本发明。

本发明提供的一种基于轨道解析摄动解的交点时刻及位置求解算法，包括如下步骤：

S1：根据双方航天器的平均轨道要素初值，使用仅考虑地球扁率J2-J4项的长期项影响的摄动模型，计算双方航天器的升交点赤经的变化速率、对方航天器的近地点幅角的变化速率和对方航天器的平近点角的变化速率；

S2：根据双方航天器的升交点赤经的变化速率，计算假定时刻双方航天器的轨道面方位；

S3：根据对方航天器的平近点角的变化速率和近地点幅角的变化速率，以及假定时刻双方航天器的轨道面方位，计算交点时刻；

S4：判断交点时刻与假定时刻的误差是否大于阈值；若是，则将交点时刻作为新的假定时刻，返回步骤S2，重复执行步骤S2～步骤S4；若否，则执行步骤S5；

S5：利用带有长周期振荡和短周期振荡的密切轨道要素摄动方程，计算交点时刻对方航天器的位置。

本发明提供的上述基于轨道解析摄动解的交点时刻及位置求解算法，所依托的设备包括星敏感器、视觉测量敏感器和星载计算机，对象是交会任务的目标航天器。通过视觉测量敏感器了解对方航天器的初始位置，使用本发明提供的上述基于轨道解析摄动解的交点时刻及位置求解算法求解对方航天器穿越己方航天器轨道平面的时刻(即交点时刻)，以及交点时刻对方航天器的位置，从而为接下来的交会操作做准备。计算交点时刻采用基于轨道长期摄动项的迭代算法，计算位置采用基于完整解析摄动解的精确算法。

下面通过一个具体的实施例并结合图1对本发明提供的上述基于轨道解析摄动解的交点时刻及位置求解算法的具体实施进行详细说明。

实施例1：

第一步，根据双方航天器的平均轨道要素初值，使用仅考虑地球扁率J2-J4项的长期项影响的摄动模型，计算双方航天器的升交点赤经的变化速率、对方航天器的近地点幅角的变化速率和对方航天器的平近点角的变化速率。

考虑到轨道为大气密度稀薄的近地轨道，且计算交点时刻考虑的时间段仅为24小时，因此，轨道摄动力仅考虑地球扁率J2到J4项，忽略太阳和月球的引力摄动、大气阻力以及地球非球形引力的高阶项等其他摄动的影响。地球扁率J2-J4项的长期项影响，会使双方航天器的升交点赤经Ω_k、对方航天器的近地点幅角ω_B以及对方航天器的平近点角M_B按照如下的变化速率进行变化：

其中，μ表示地球引力常数，R_E表示地球半径，J₂、和J₄表示地球扁率项；i_k表示双方航天器的轨道倾角，a_k表示双方航天器的轨道半长轴，e_k表示双方航天器的离心率，i_k、a_k和e_k在长期项作用下保持不变；k＝A,B，A表示己方航天器，B表示对方航天器。

除长期项带来的三种轨道要素(升交点赤经Ω、近地点幅角ω和平近点角M)的平均变化外，J2到J4项还会使得航天器密切轨道要素产生长周期和短周期振荡。如果考虑长周期和短周期振荡，轨道要素会变得极为复杂，交点时刻的计算时间也会大大增加。故放宽对方航天器严格位于己方航天器轨道面这一要求，使用仅考虑长期项的平均轨道要素摄动方程以降低计算量。在计算双方航天器的轨道面参数时，半长轴a、偏心率e和轨道倾角i不受长期项影响，因此，仅考虑另外三种轨道要素(升交点赤经Ω、近地点幅角ω和平近点角M)的长期项影响以计算变化速率，即可基于初始轨道参数得到假定时刻的轨道参数。

第二步，根据双方航天器的升交点赤经的变化速率，计算假定时刻双方航天器的轨道面方位。

由双方航天器的升交点赤经的变化速率

可以求得假定时刻双方航天器的升交点赤经：

其中，Ω_k0表示双方航天器初始的升交点赤经，t表示任意假定时刻；双方航天器的轨道倾角保持不变：

i_k＝i_k0 (5)

其中，i_k0表示双方航天器初始的轨道倾角；由假定时刻双方航天器的升交点赤经Ω_k和轨道倾角i_k表示假定时刻双方航天器的轨道面方位。

第三步，根据对方航天器的平近点角的变化速率和近地点幅角的变化速率，以及假定时刻双方航天器的轨道面方位，计算交点时刻。

(1)假定交点时刻为t′_n，n＝1，取t′₁＝t，由对方航天器的近地点幅角的变化速率

可以计算假定交点时刻对方航天器的近地点幅角：

其中，ω_B0表对方航天器初始的近地点幅角；

(2)计算假定交点时刻对方航天器的真近点角：

θ′_B＝u_B-ω′_B (7)

其中，u_B表示对方航天器穿过己方航天器轨道平面时的纬度幅角；u_B的计算公式如下：

其中，

其中，

表示对方航天器轨道节线的单位矢量；

表示双方航天器轨道面的交线单位矢量；

表示己方航天器的轨道法线方向；

表示对方航天器的轨道法线方向；

(3)根据假定交点时刻对方航天器的真近点角θ′_B，可以计算假定交点时刻对方航天器的偏近点角E′_B：

(4)根据假定交点时刻对方航天器的偏近点角E′_B，可以计算假定交点时刻对方航天器的平近点角：

M′_B＝E′_B-e_BsinE′_B (14)

(5)根据假定交点时刻对方航天器的平近点角M′_B，可以得到更新后的交点时刻为：

其中，M_B0表示对方航天器初始的平近点角；

(6)判断更新后的交点时刻t′_i+1与假定交点时刻t′_n的误差是否大于阈值，阈值可以在10^-6s～10^-9s范围内选择；当更新后的交点时刻t′_i+1与假定交点时刻t′_n的误差小于或等于阈值时，将更新后的交点时刻t′_n+1作为交点时刻；当更新后的交点时刻t′_i+1与假定交点时刻t′_n的误差大于阈值时，将更新后的交点时刻t′_n+1作为新的假定时刻，返回第(1)步，n＝n+1，重复第(1)步～第(6)步，直至迭代收敛(即更新后的交点时刻t′_i+1与假定交点时刻t′_n的误差小于或等于阈值)。

第四步，判断交点时刻与假定时刻的误差是否大于阈值，阈值可以在10^-6s～10^-9s范围内选择；当交点时刻与假定时刻的误差大于阈值时，将交点时刻作为新的假定时刻，返回第二步，重复执行第二步～第四步；当交点时刻与假定时刻的误差小于或等于阈值时，执行第五步。

第五步，利用带有长周期振荡和短周期振荡的密切轨道要素摄动方程，计算交点时刻对方航天器的位置。

在已知交点时刻的基础上，利用对方航天器的完整的J2-J4项的解析摄动解，包括长期项、长周期项和短周期项的结果，以及对方航天器的平均轨道要素初值，计算交点时刻t′对方航天器的密切轨道要素：

密切轨道要素(t′)＝长期项(t′)+长周期项(t′)+短周期项(t′)+平均轨道要素初值。

使用完整的考虑J2-J4项影响的解析摄动解，即同时包括长期项、长周期项和短周期项的结果，根据对方航天器的初始平均轨道要素得到交点时刻的密切轨道要素，因其解析性可以较快地求得交点时刻对方航天器的位置，并且得到的位置非常精确。

本发明提供的上述基于轨道解析摄动解的交点时刻及位置求解算法，在不同阶段使用不同的模型：在求解交点时刻过程中，轨道摄动只考虑地球扁率J2-J4项中的长期项，忽略周期项影响，通过放宽对方航天器严格位于己方航天器轨道面这一要求来大大降低交点时刻和位置的计算量，具有非常高的计算速度，可以极大地节约计算资源，虽然所得结果为近似交点时刻，但此时对方航天器仍在己方航天器轨道平面附近，在后续交会问题求解时，己方航天器只需很小的轨道面外速度增量就可到达对方位置，可以满足任务要求，预期对己方航天器的可行弧段影响很小；最后根据所得较低时刻求解航天器位置时再使用考虑周期项的完整模型，以此得到精确的对方航天器位置。

下面从四个方面说明本发明实施例1提供的上述基于轨道解析摄动解的交点时刻及位置求解算法的性能。

(1)为了判断本发明计算得到的交点时刻和位置的准确性，下面将本发明所得结果与stk中HPOP模块的轨道仿真结果进行对比。使用HPOP模块进行仿真时，其重力模型选用EGM2008，且考虑海洋潮汐的影响。航天器空气阻力系数为2.2，面质比为0.01m^2/kg，太阳光压系数为1。同时考虑太阳与月球三体引力，其余使用默认条件。考虑所有轨道摄动项，以此作为真实轨道，以分析本发明计算结果的误差。

本发明所得结果与stk仿真结果的对比，需要明确轨道要素初值的确定方式，即定轨方式，图2为定轨方式的说明图。由于HPOP模块中轨道要素为密切轨道要素即瞬根，因此，首先，使用瞬根初值在HPOP模块仿真一天，利用全天的瞬根数据，以弧段定轨的方式生成该天最后时刻的TLE(两行轨道要素，为平均轨道要素即平根)数据，以此作为交点时刻计算和位置预报的初值；然后，以该天的HPOP瞬根末值，以相同的条件在HPOP中仿真一天，作为真实轨道，和本发明得到的结果进行对比，即弧段平根初值和HPOP瞬根末值相对应。双方航天器的轨道要素初值如表1所示。

表1

(2)下面通过求解交点时刻的距离误差和位置误差，分析本发明实施例1提供的上述基于轨道解析摄动解的交点时刻及位置求解算法的精度。

距离误差即在所得交点时刻对方航天器到己方航天器轨道平面的真实距离，可以反映本发明得到的交点时刻的准确性。位置误差即在所得交点时刻使用同时考虑长期项与周期项的J2-J4模型得到的对方航天器位置与stk中仿真得到的真实位置之间的误差，可以反映对方航天器位置预报的精度。图3和图4分别为距离误差变化图和位置误差变化图，由图3可以看出，距离误差在数千米内，说明得到的交点时刻是比较准确的，由图4可以看出，位置误差最高约1km，说明J2-J4模型的精度很高，得到的对方航天器的位置比较准确。

(3)下面对比本发明实施例1提供的上述基于轨道解析摄动解的交点时刻及位置求解算法与已有的另外两种计算交点时刻方法，以说明本发明实施例1提供的上述基于轨道长期项和解析摄动解计算交点时刻及位置的方法的准确性。

分别使用三种交点时刻计算方法得到不同的交点时刻，并计算对方航天器与己方航天器轨道面的真实距离、对方航天器位置矢量与己方航天器轨道面的真实夹角，以分析每种交点时刻计算方法的计算精度。下面对三种方法原理作简单介绍：方法一即本发明实施例1提供的上述基于轨道解析摄动解的交点时刻及位置求解算法；方法二，己方航天器仍使用长期项递推得到轨道要素，对方航天器根据HPOP的离散结果，插值得到某个时间的密切轨道要素，根据密切轨道要素初值计算得到大致的交点时刻，在大致的交点时刻附近使用真近角的数据进行插值迭代，求得交点时刻；方法三，在不考虑摄动的情况下，计算相交时纬度幅角，以得到交点时刻。

下面对比三种方法的结果。由于三种方法得到的交点时刻比较接近，因此，以方法一为对照，计算方法二与方法三分别同方法一的交点时刻的差值变化图如图5所示。由图5可以看出，方法二与方法一的结果比较接近，但会有10s以内的波动，方法三和方法一的结果差别比较大，最高有60s以上的差别。

下面计算交点时刻对方航天器与己方航天器轨道面的真实距离、对方航天器位置矢量与己方航天器轨道面的真实夹角，这可以代表方法的结果精度。方法一和方法二的交点时刻真实距离变化图如图6所示，方法一和方法三的交点时刻真实距离变化图如图7所示，由图6和图7可以看出，方法一的距离误差最高约6km，方法二的距离误差最高约48km，方法三的距离误差最高约379km。方法一和方法二的交点时刻真实夹角变化图如图8所示，方法一和方法三的交点时刻真实夹角变化图如图9所示，由图8和图9可以看出，方法一的真实夹角最高约0.03°，方法二的真实夹角最高约0.30°，方法三的真实夹角最高约2.40°。综上可以得出，方法一的精度最高，方法二稍次之，方法三的精度要远低于方法一与方法二。

(4)下面分别计算本发明使用的J2-J4模型与较为简单的解析递推模型(后简称对比模型)在24小时内的预报误差演变情况，以分析本发明使用的J2-J4模型的精确程度。定轨方式与前面相同，仍为STK弧段定轨，误差的参照是对应瞬根初值在HPOP中仿真得到的轨道(步长为3s)。以对方航天器为例，其轨道要素初值如表2所示。

表2

两种模型(J2-J4模型与对比模型)在一天内位置误差比较结果如图10所示，由于J2-J4模型的误差过小，因此单独在图11中给出。从以上结果可以看出，对比模型的位置误差在一天内大致呈周期变化，最高可达数百甚至上千千米，而J2-J4模型的位置误差最高约2km左右，因此，J2-J4模型的精度远远高于对比模型。从位置误差的变化情况来看，虽然对比模型的计算速度要高于本发明的J2-J4模型(一天内计算了28800组轨道数据)，但对比模型的计算精度远不如J2-J4模型。J2-J4模型的计算速度与精度都非常适用于实际轨道预报。

综上，通过将本发明所得结果与stk中精确仿真结果进行对比，验证了本发明实施例1提供的上述基于轨道解析摄动解的交点时刻及位置求解算法的正确性和有效性。并且，本发明与不考虑轨道摄动，或是对方航天器使用真实轨道数据插值得轨道要素的计算方法相比，计算速度相差较小，但计算精度明显占优势。

本发明提供的上述基于轨道解析摄动解的交点时刻及位置求解算法，计算交点时刻采用仅考虑地球扁率J2-J4项的长期项影响的摄动模型，忽略周期项影响以降低计算量，利用摄动方程计算假定时刻双方航天器轨道面方位，并通过几何关系计算得到交点时刻，将计算得到的交点时刻与假定时刻对比，若误差较大，则将计算得到的交点时刻作为新的假定时刻以重复计算交点时刻，直至迭代收敛，最后利用同时带有长期项、长周期和短周期振荡的密切轨道要素摄动方程，求得交点时刻对方航天器的位置。上述方法不仅可以精确计算交点时刻及位置，而且算法简单易于实施，节约计算时间，节省计算资源和内存，特别适用于计算资源紧张的星载计算机，具有良好的推广前景。

显然，本领域的技术人员可以对本发明进行各种改动和变型而不脱离本发明的精神和范围。这样，倘若本发明的这些修改和变型属于本发明权利要求及其等同技术的范围之内，则本发明也意图包含这些改动和变型在内。

Claims (6)

1.一种基于轨道解析摄动解的交点时刻及位置求解算法，其特征在于，包括如下步骤：

S1：根据双方航天器的平均轨道要素初值，使用仅考虑地球扁率J2-J4项的长期项影响的摄动模型，计算双方航天器的升交点赤经的变化速率、对方航天器的近地点幅角的变化速率和对方航天器的平近点角的变化速率；

S2：根据双方航天器的升交点赤经的变化速率，计算假定时刻双方航天器的轨道面方位；

S3：根据对方航天器的平近点角的变化速率和近地点幅角的变化速率，以及假定时刻双方航天器的轨道面方位，计算交点时刻；

S4：判断交点时刻与假定时刻的误差是否大于阈值；若是，则将交点时刻作为新的假定时刻，返回步骤S2，重复执行步骤S2～步骤S4；若否，则执行步骤S5；

S5：利用带有长周期振荡和短周期振荡的密切轨道要素摄动方程，计算交点时刻对方航天器的位置。

2.如权利要求1所述的基于轨道解析摄动解的交点时刻及位置求解算法，其特征在于，步骤S1，根据双方航天器的平均轨道要素初值，使用仅考虑地球扁率J2-J4项的长期项影响的摄动模型，计算双方航天器的升交点赤经的变化速率、对方航天器的近地点幅角的变化速率和对方航天器的平近点角的变化速率，具体包括：

仅考虑地球扁率J2-J4项的长期项影响，双方航天器的升交点赤经Ω_k的变化速率

对方航天器的近地点幅角ω_B的变化速率

和对方航天器的平近点角M_B的变化速率

为：

其中，μ表示地球引力常数，R_E表示地球半径，J₂、和J₄表示地球扁率项；i_k表示双方航天器的轨道倾角，a_k表示双方航天器的轨道半长轴，e_k表示双方航天器的离心率，i_k、a_k和e_k在长期项作用下保持不变；k＝A,B，A表示己方航天器，B表示对方航天器。

3.如权利要求2所述的基于轨道解析摄动解的交点时刻及位置求解算法，其特征在于，步骤S2，根据双方航天器的升交点赤经的变化速率，计算假定时刻双方航天器的轨道面方位，具体包括：

由双方航天器的升交点赤经的变化速率

求得假定时刻双方航天器的升交点赤经：

其中，Ω_k0表示双方航天器初始的升交点赤经，t表示任意假定时刻；双方航天器的轨道倾角保持不变：

i_k＝i_k0 (5)

其中，i_k0表示双方航天器初始的轨道倾角；由假定时刻双方航天器的升交点赤经Ω_k和轨道倾角i_k表示假定时刻双方航天器的轨道面方位。

4.如权利要求3所述的基于轨道解析摄动解的交点时刻及位置求解算法，其特征在于，步骤S3，根据对方航天器的平近点角的变化速率和近地点幅角的变化速率，以及假定时刻双方航天器的轨道面方位，计算交点时刻，具体包括如下步骤：

S31：假定交点时刻为t'_n，n＝1，取t'₁＝t，由对方航天器的近地点幅角的变化速率

计算假定交点时刻对方航天器的近地点幅角：

其中，ω_B0表对方航天器初始的近地点幅角；

S32：计算假定交点时刻对方航天器的真近点角：

θ'_B＝u_B-ω'_B (7)

其中，u_B表示对方航天器穿过己方航天器轨道平面时的纬度幅角；u_B的计算公式如下：

其中，

其中，

表示对方航天器轨道节线的单位矢量；

表示双方航天器轨道面的交线单位矢量；

表示己方航天器的轨道法线方向；

表示对方航天器的轨道法线方向；

S33：计算假定交点时刻对方航天器的偏近点角：

S34：计算假定交点时刻对方航天器的平近点角：

M'_B＝E'_B-e_Bsin E'_B (14)

S35：更新后的交点时刻为：

其中，M_B0表示对方航天器初始的平近点角；

S36：判断更新后的交点时刻t'_i+1与假定交点时刻t'_n的误差是否大于阈值；若否，则输出更新后的交点时刻t'_n+1；若是，则将更新后的交点时刻t'_n+1作为新的假定时刻，返回步骤S31，n＝n+1，重复步骤S31～步骤S36。

5.如权利要求4所述的基于轨道解析摄动解的交点时刻及位置求解算法，其特征在于，步骤S3中判断更新后的交点时刻t'_n+1与假定交点时刻t'_n的误差，以及步骤S4中判断交点时刻与假定时刻的误差，所依据的阈值的范围为10^-6s～10^-9s。

6.如权利要求1～5任一项所述的基于轨道解析摄动解的交点时刻及位置求解算法，其特征在于，步骤S5，利用带有长周期振荡和短周期振荡的密切轨道要素摄动方程，计算交点时刻对方航天器的位置，具体包括：

利用对方航天器的完整的J2-J4项的解析摄动解，包括长期项、长周期项和短周期项的结果，以及对方航天器的平均轨道要素初值，计算交点时刻t′对方航天器的密切轨道要素：

密切轨道要素(t′)＝长期项(t′)+长周期项(t′)+短周期项(t′)+平均轨道要素初值。

Images