CN103728980A - 航天器相对轨道的控制方法 - Google Patents

Abstract

航天器相对轨道的控制方法,本发明涉及航天器的近距离相对轨道控制方法。以实现航天器的掠飞模式，即追踪航天器在进入与目标航天器相关的指定空间范围后按自身轨道运行，只需要进行姿态控制；从而克服传统的悬停、伴随飞行、绕飞等方法可能出现计算复杂、姿轨控耦合导致指向精度不高、易暴露身份、时间难以保持等问题。本发明的方法通过下述步骤实现：一、追踪航天器进入目标航天器的视线角范围内而且追踪航天器进入二者之间确定的距离范围内；二、计算并确定追踪航天器期望轨道的起点、末点和初始入轨速度，并确定主飘方向；三、追踪航天器在期望轨道的起点，以上述计算并确定的初始入轨速度进入轨道，并在期望轨道的末点脱离轨道。

Description

航天器相对轨道的控制方法

技术领域

本发明涉及航天器的近距离相对轨道控制方法。 

背景技术

当今航天领域的一个重要研究热点就是航天器的近距离相对轨道运动控制。航天器的相对轨道运动是研究一个航天器（追踪航天器）处于另一个航天器（目标航天器）周围的持续运动规律。常被应用到编队飞行、在轨维护、交会对接、跟踪监视等空间任务。目前最常用的相对轨道运动形式有悬停（追踪航天器与目标航天器保持相对位置不变）、伴随飞行（追踪航天器围绕目标航天器附近某点进行封闭轨迹飞行）和绕飞（伴随飞行的一种特殊情况，封闭轨迹的中心是目标航天器）等。 

其中，对空间非合作目标的跟踪监视问题更加体现了相对轨道控制的重要性，随着航天器机动性的增强，跟踪监视的精度、范围等要求也越来越高，当要求追踪航天器在目标航天器附近一个限制的空间范围内保持一定时间的姿态指向时，传统的悬停、伴随飞行、绕飞等方法可能出现计算复杂、姿轨控耦合导致指向精度不高、易暴露身份、时间难以保持等问题。 

发明内容

本发明提供一种航天器相对轨道的控制方法，以实现航天器的掠飞模式，即追踪航天器在进入与目标航天器相关的指定空间范围后按自身轨道运行，只需要进行姿态控制；从而克服传统的悬停、伴随飞行、绕飞等方法可能出现计算复杂、姿轨控耦合导致指向精度不高、易暴露身份、时间难以保持等问题。 

本发明的方法通过下述步骤实现：一、追踪航天器进入目标航天器的视线角范围内而且追踪航天器进入二者之间确定的距离范围内；二、计算并确定追踪航天器期望轨道的起点、末点和初始入轨速度，并确定主飘方向；三、追踪航天器在期望轨道的起点，以上述计算并确定的初始入轨速度进入轨道，并在期望轨道的末点脱离轨道，从而完成对目标航天器的掠飞动作。 

本发明在确定期望轨道的起点、末点和初始入轨速度和主飘方向的时候采用了“主飘方向”分析法，在找到主飘方向的期望起点与终点情况下，这种分析法将一个相对hill坐标系的三维掠飞问题，化成了一个一维初速度值求解问题，进而再通过一些简单的原则化成一个三维初速度值确定问题，分析思路简单。与无轨道控制的伴随飞行和悬停方案相比，具有视场方向可以任意的优点；与无轨道控制的绕飞方案相比，具有掠飞时间可设计的优点；与脉冲控制的“水滴”形轨迹飞行方案相比，具有掠飞区域中完全不进行轨道控制，以实现足够长时间内完成高精度姿态指向控制等任务需求的优点；与连续控制的悬停方案相比，除了避免姿轨耦合控制，利于提高姿态指向控制等任务精度外，还具有不易暴露身份的优点。 

附图说明

图1本发明方法的流程示意图，图2掠飞区域说明图，图3是z₀的确定方式示意图，图4圆锥台侧视图，图5是近似圆柱区域示意图，图6是z=z₀在圆锥体表面截出的近似椭圆区域分析图，图7是c的变化趋势图，图8是c在x轴上的情况示意图，图9是c在y轴上的情况示意图，图10是c与x轴夹角小于45°的情况示意图，图11是c与x轴夹角大于45°的情况示意图，图12是主飘方向分析图，图13是允许掠飞的圆锥台区域的后视图，图14是地心惯性系与轨道坐标系示意图，图15是轨道平面的示意图，图16是相对运动坐标系与地心惯性坐标系的关系图，图17是两星相对位置示意图，图18是轨迹投影在xy面形成水滴形的示意图，图19是追踪星B与目标星A的空间关系示意图，图20是追踪星受力分析图，图21是掠飞轨迹示意图，图22是掠飞轨迹示意图放大图，图23目标星坐标系相对轨迹三轴分量时间变化曲线对比图，图24是掠飞时间触发信号示意图。 

具体实施方式

具体实施方式一：下面结合图1具体说明本实施方式。本实施方式通过下述步骤实现：一、追踪航天器进入目标航天器的视线角范围内而且追踪航天器进入二者之间确定的距离范围内；二、计算并确定追踪航天器期望轨道的起点、末点和初始入轨速度，并确定主飘方向；三、追踪航天器在期望轨道的起点，以上述计算并确定的初始入轨速度进入轨道，并在期望轨道的末点脱离轨道，从而完成对目标航天器的掠飞动作。 

（一）相关概念的解释和自定义： 

地心惯性坐标系（O-XYZ）：坐标原点O在地球质心，Z轴向北指向平赤道面北极，X轴指向平春分点，Y轴与Z轴、X轴组成直角右手坐标系。 

轨道坐标系（s-xyz）：坐标原点s与航天器的质心固连并随其沿轨道运动，x轴与航天器的地心矢量

重合，由地心指向s，y轴在航天器的轨道面内垂直于x轴，并指向运动方向，z轴与x轴，y轴组成直角右手坐标系。 

开普勒方程：两体运动问题中，对于椭圆轨道，偏心率为e，偏近点角为E，平近点角为M，沿逆时针方向为正，开普勒方程可以表示为 

E-esin(E)=M    (1) 

Hill方程：假定两航天器仅受地球的引力作用，以目标航天器的轨道坐标系作为相对运动坐标系，通过一阶线性化，将相对运动动力学方程化为一组常系数线性微分方程 

$\left\{\begin{array}{c}\stackrel{\·\·}{x}-2n\stackrel{\·}{y}-2n^{2}x=a_x\\ \stackrel{\·\·}{y}+2n\stackrel{\·}{x}=a_y\\ \stackrel{\·\·}{z}+n^{2}z=a_z\end{array}\right.---\left(2\right)$

脉冲控制：发动机短时工作情况下的轨道控制，由于发动机工作的时间比轨道飞行周期短得多，因而可以视为是脉冲作用。 

相对轨道掠飞运动：追踪航天器在进入与目标航天器相关的指定空间范围后按自身轨道运行，不需要进行主动轨道控制下的一种相对轨道运动。 

“主飘方向”分析法：如果要求追踪航天器在目标航天器的视场范围内，hill坐标系三轴分量同时满足保持相同掠飞时间的要求，对初始位置信息

的要求将过于严格。由于hill相对坐标系的xy平面和z方向的运动是独立的，因此可以简单进行如下设计，即x或y一个方向的分量的掠飞时间满足相应的时间要求，另外两个方向的分量的掠飞时间只要长于该方向的时间即可。这种分析方法称之为“主飘方向”分析法。 

（二）掠飞原理分析 

非合作目标处在圆形高轨轨道，以地球静止轨道GEO为例。 

追踪星保持在非合作目标的视线角±θ小角度(一般小于1°)范围内且距离目标星有上下边界的要求。 

追踪星在要求的掠飞区域内掠飞几百秒时间，以便于通过姿态控制完成对非合作目标高精度指向保持等任务需求。 

在基于hill方程的相对运动坐标系下考虑，如，图2所示 

其中e=[e₁;e₂;e₃]是非合作目标星体视线矢量l的单位矢量，e_i(i=1,2,3)是在hill方程坐标系下的分量。由此可见视线方向和两星距离决定了追踪星允许掠飞的范围，是hill方程坐标系下的一段圆锥台。设计的目标是基于hill方程坐标系选取合适的进入点A信息[r₀,v₀](r₀=[x₀y₀z₀]^T，v₀=[v_x0v_y0v_z0]^T)和离开点B信息[r_f,v_f](r_f=[x_f y_f z_f]^T，v_f=[v_xf v_yf v_zf]^T)，使得追踪星在不加任何主动轨道控制的基础上，在圆锥台内自由掠飞，并且保持满足任务要求的时间。 

（三）基于“主飘方向”分析法的掠飞轨迹设计 

3.1对于z₀的选取 

如0-3所示，z₀的选取原则是使得平面z=z₀与圆锥台所截取的xy平面尽量大，且z₀容易确定。对于本设计，取圆锥台下底面圆心对应的z轴坐标为z₀。由e的定义，如图3所示， 

z₀=lsinα=l·e₃    (3) 

3.2关于掠飞段进入速度的说明 

对于不加轨道控制的自由掠飞过程，存在两种情况，即顺着目标航天器轨道角速度和逆着轨道角速度两种情况。对于本发明申请，研究顺着目标航天器轨道角速度的情况，即选取相对hill系的初始速度的y方向分量为正进行设计。 

（四）关于xy面的图形分析 

对于本发明，由于视线角较小，即圆锥台的顶角非常小，只有2×θ，因此圆锥的母线围成的图形的投影形成的等腰三角形的底角较大，近似为90°（90°-θ）。因此可以看成是一个矩形，如图4、图5所示。从图5可以看出目标星视线轴范围内的圆锥区域。 

因此，实际上在相对较短的区域内，追踪星掠飞的区域近似是一段圆柱。因此z=z₀在圆柱表面截出的区域是一段椭圆，如图6所示。 

（五）关于椭圆相对xy面位置的分析与“主飘方向”分析法 

可以容易看出，椭圆的一个半轴长度为圆锥台底面半径R，另一半轴长度c由e与xy面的夹角α决定，即 

c=R/sinα    (4) 

可见随着α的变化，c也发生着变化，即 

$\mathrm{\α}:0\→\frac{\mathrm{\π}}{2},c:\∞\→R---\left(5\right)$

也就是，随着α的变大，c逐渐变小，最后将使得截出的区域是一个半圆，如图7所示。 

对于该椭圆，根据R与c的关系，可以分成两大类进行分析，即c＞＞2R和c<2R两种情况。以下对第一种情况进行分析，对于第二种情况可类比分析。 

对于该截出的椭圆区域（c＞＞2R），根据其对称轴在xy平面投影方向的不同，可以分成以下四种主要情况： 

5.1对称轴投影在x轴上情况 

对于该情况，如图8所示。这时很明显，2R是处在y轴方向上，并且由于2R较小，因此，选择y轴为主飘方向。x、z方向的初始速度保证一个较小值即可。 

5.2对称轴投影在y轴上情况 

对于本情况，如图9所示。分析方法与第一种情况类似，且选x轴为主飘方向。y、 z方向的初始速度保证一个较小值即可。 

5.3对称轴投影与x轴夹角小于45°情况 

对于本情况，如图10所示。可以看出，此时各个轨迹y轴方向的分量总是小于x轴方向的分量，因此相当于y轴更敏感，因此选取y轴为主飘方向。x、z方向的初始速度保证一个较小值即可。 

5.4对称轴投影与x轴夹角大于45°情况 

对于本情况，如图11所示。可以看出，此时各个轨迹x轴方向的分量总是小于y轴方向的分量，因此相当于x轴更敏感，因此选取x轴为主飘方向。y、z方向的初始速度保证一个较小值即可。 

（六）关于主飘方向理想主飘轨迹始末点的计算以及主飘方向初始速度计算 

以对称轴与x轴夹角小于45°情况为例进行分析，如图12所示。 

期望的主飘方向是从Y₀点到Y_f的轨迹。各点的hill方程下的坐标表示计算如下： 

首先Y₀和Y_f的x坐标相同，为 

x_f=X₀-Rsinβ    (6) 

其中X₀为圆锥体底面中心x坐标，R为底面圆半径，β的计算满足： 

$\mathrm{\&beta;}=\arctan \left(\frac{e_2}{e_1}\right)---\left(7\right)$

由于Y₀和Y_f两点都在圆锥台的表面上，因此其与原点的连线与目标星视线方向l成θ，通过求解如下方程即可确定两点的y方向坐标值，即 

$\mathrm{ar}\cos ((x_f{e}_1+y_{o,f}{e}_2+z_0{z}_3)/\sqrt{x_f^2+y_{0,f}^2+z_0^2})=\mathrm{\&theta;}---\left(8\right)$

该方程有两个解，其中y_o,f为两点的y方向坐标。 

z方向的坐标已由前述式(3)确定。 

在给出了期望主飘轨线的进入点和离开点后，给出主飘方向的期望速度分量的计算原则和方法。 

设期望掠飞时间是t_lve，则主飘方向的期望速度分量计算为： 

${\stackrel{.}{y}}_0=\frac{|Y_f-Y_0|}{t_{\mathrm{lve}}}---\left(9\right)$

（七）关于非主飘方向的初始速度说明 

理想的主飘轨迹是希望能沿着直线，可由于这时相当于两个非主飘方向的相对速度严格为0，这种情况显然过于理想，无法达到。并且，如果两个非主飘方向选取得当，可以起到减少期望掠飞时间损失或延长期望掠飞时间的效果。 

依然以图12为例进行说明。 

假设y方向的初始位置和速度都已给出，观察图12，如果x方向存在过大的正的初始速度，那么将使得最终轨迹可能沿着Y₀Y_f1飞行，在y的投影小于预计的主飘轨迹长度；如果x方向的速度为负，将沿着轨迹Y₀Y_f2提前飞出允许的区域。因此在设计好的y的初始速度基础上，都将使实际的掠飞时间受到严重损失。因此x方向的速度应为正，且不应过大。 

对于z方向的分析，如图13所示。z方向如果具有一定的负方向的初始速度，那么实际轨迹将沿着Y₀Y_f3飞行，在xy方向的投影，将大于预计的主飘轨迹长度，因此在y方向投影将得到延长，在y方向初始速度设计好前提下，将使掠飞时间得到延长。相反，类似的分析方法，如果z方向的初始速度具有正方向，将使得掠飞时间受到损失。 

（八）掠飞初始点变轨所需速度脉冲坐标系分量转换方法 

假设在任意时刻t，目标星的位置和速度信息在惯性系表示分别为r_t，v_t，则有： 

$e_x=\frac{r_t}{\left|r_t\right|},e_y=\frac{v_t}{\left|v_t\right|},e_z=e_x\&times;e_y$

则从惯性系到hill方程坐标系的余弦矩阵为： 

$\mathrm{DCMTi}=\left[\begin{array}{c}{e}_x^T\\ e_y^T\\ e_z^T\end{array}\right]---\left(10\right)$

追踪星在进入掠飞区域初始点t时刻的变轨所需要的相对惯性系的速度脉冲为： 

$\begin{array}{c}\mathrm{\&Delta;v}=\left[\begin{array}{c}\mathrm{\&Delta;}\frac{\mathrm{dx}}{{\mathrm{dt}}^{+}}\\ \mathrm{\&Delta;}\frac{\mathrm{dy}}{{\mathrm{dt}}^{+}}\\ \mathrm{\&Delta;}\frac{\mathrm{dz}}{{\mathrm{dt}}^{+}}\end{array}\right]-\left[\begin{array}{c}\mathrm{\&Delta;}\frac{\mathrm{dx}}{{\mathrm{dt}}^{-}}\\ \mathrm{\&Delta;}\frac{\mathrm{dy}}{{\mathrm{dt}}^{-}}\\ \mathrm{\&Delta;}\frac{\mathrm{dz}}{{\mathrm{dt}}^{-}}\end{array}\right]={\mathrm{DCMTi}}^T\{{\left[\begin{array}{c}\stackrel{\&CenterDot;}{x}\\ \stackrel{\&CenterDot;}{y}\\ \stackrel{\&CenterDot;}{z}\end{array}\right]}_{t+}-{\left[\begin{array}{c}\stackrel{\&CenterDot;}{x}\\ \stackrel{\&CenterDot;}{y}\\ \stackrel{\&CenterDot;}{z}\end{array}\right]}_{t^{-}}\}\\ ={\mathrm{DCMTi}}^T{\mathrm{\&Delta;v}}_h\end{array}---\left(11\right)$

其中△v_h为hill方程系下的速度增量。 

（九）理论基础 

9.1两体问题轨道模型 

航天器轨道是指航天器在地心惯性系下的位置，描述的是航天器的平动运动。 

9.1.1Kepler轨道及轨道要素 

作为二体问题，两个航天器均作为质点对待。考虑质量分别为m和M的两个天体间运动方程为： 

$\stackrel{\&CenterDot;\&CenterDot;}{r}=-\frac{G(M+m)}{r^3}r=-\frac{\mathrm{\&mu;}}{r^3}r---\left(12\right)$

通过对二体问题的求解，天体（航天器）在惯性空间中的运动可以用六个经典的轨道要素（也称轨道根数，如图14所示）来描述： 

a:轨道半长轴； 

e:偏心率； 

Ω:升交点赤经，在J2000地心惯性系中自Ox轴方向在xy平面内逆时针度量到升交点的角度，0≤Ω<2π； 

i:轨道倾角，轨道正法向h与J2000地心惯性系Oz轴的夹角，0≤i≤π。若0≤i<π/2，为顺行轨道，航天器偏向东飞行；若π/2≤i<π，为逆行轨道，航天器偏向西飞行；若i=π/2，为极轨道； 

ω:近地点幅角，在轨道平面内自轨道升交点方向到偏心率矢量e之间的夹角，顺航天器方向度量，0≤ω<2π； 

f:真近点角，在轨道平面内从e到r之间的夹角。 

9.1.2轨道计算与星历表计算 

给定初始条件，从观测的航天器位置和速度数据出发，可以计算轨道根数，这一过程称为轨道计算。同时，轨道根数给定后，也可以计算出任意时刻航天器的位置矢量和速度矢量，该过程称为星历表计算。 

以下简单介绍一下两种过程的计算方法： 

9.1.2.1轨道计算 

求半长轴a： 

$a=\frac{\mathrm{\&mu;r}}{2\mathrm{\&mu;}-{\mathrm{rv}}^2}---\left(13\right)$

其中v为航天器轨道速度。 

求偏心率e： 

$e=\sqrt{1-p/a}---\left(14\right)$

其中p为椭圆半通径。 

轨道倾角i与升交点赤经计算如下： 

i=arccos(h_z/h)    (15) 

Ω=-arctan(2h_x/h_y)    (16) 

求真近点角f： 

cosf=(p/r-1)/e    (17) 

求近地点幅角ω： 

令i_Ω为地心指向升交点的单位矢量，在地心惯性系中表示为i_Ω=cosΩ·x_i+sinΩ·y_i。故可由下式确定轨道角u 

$\cos u=\frac{1}{r}{r}^T{i}_{\mathrm{\&Omega;}}---\left(18\right)$

进而可得近地点幅角为 

ω=u-f    (19) 

注：当e=0时，会出现奇异。因此一般考虑当e<10^-6时，轨道称为圆轨道，有ω=0，f=u。 

9.1.2.2星历表计算 

在轨道平面内定义如图15所示的坐标系，有 

$\begin{array}{c}r=a(\cos E-e)\&CenterDot;P+a\sqrt{1-e^2}\sin E\&CenterDot;Q\\ v=-\sqrt{\mathrm{\&mu;}/p}\sin f\&CenterDot;P+\sqrt{\mathrm{\&mu;}/p}(e+\cos f)\&CenterDot;Q\end{array}---\left(20\right)$

其中，P和Q分别表示卫星近地点和半通径方向的单位矢量。通过3-1-3的旋转，P和Q的表达式如下所示： 

$P=R_3(-\mathrm{\&Omega;})R_1(-i)R_3(-\mathrm{\&omega;})\left[\begin{array}{c}1\\ 0\\ 0\end{array}\right]\left[\begin{array}{c}\cos \mathrm{\&Omega;}\cos \mathrm{\&omega;}-\sin \mathrm{\&Omega;}\sin \mathrm{\&omega;}\cos i\\ \sin \mathrm{\&Omega;}\cos \mathrm{\&omega;}+\cos \mathrm{\&Omega;}\sin \mathrm{\&omega;}\cos i\\ \sin \mathrm{\&omega;}\sin i\end{array}\right]---\left(21\right)$

$Q=R_3(-\mathrm{\&Omega;})R_1(-i)R_3(-\mathrm{\&omega;})\left[\begin{array}{c}0\\ 1\\ 0\end{array}\right]\left[\begin{array}{c}-\cos \mathrm{\&Omega;}\sin \mathrm{\&omega;}-\sin \mathrm{\&Omega;}\cos \mathrm{\&omega;}\cos i\\ -\sin \mathrm{\&Omega;}\sin \mathrm{\&omega;}+\cos \mathrm{\&Omega;}\cos \mathrm{\&omega;}\cos i\\ \cos \mathrm{\&omega;}\sin i\end{array}\right]---\left(22\right)$

9.2卫星相对轨道动力学方程与分析 

9.2.1相对轨道动力学方程 

记目标星为s，追踪星为c。取目标星的轨道坐标系s-xyz作为相对运动坐标系，其原点与目标星的质心固连并随其沿轨道运动，x轴与目标星的地心矢量r_s重合，由地心指向s，y轴于目标星的轨道面内垂直于x轴，并指向运动方向，z轴由右手规则确定，亦即z轴与目标星轨道动量矩矢量的方向一致。轨道坐标系s-xyz与地心惯性坐标系O_E-XYZ的关系如图16所示。 

在轨道坐标系中有 

$r_s=\left[\begin{array}{c}{r}_s\\ 0\\ 0\end{array}\right]---\left(23\right)$

设追踪星的地心位置为r_c，则其对于目标星的位置矢量ρ为 

$\mathrm{\&rho;}=r_c-r_s=\left[\begin{array}{c}x\\ y\\ z\end{array}\right]---\left(24\right)$

在地心惯性坐标中，目标星和追踪星的动力学方程如下 

$\frac{d^2{r}_s}{{\mathrm{dt}}^2}=-\frac{{\mathrm{\&mu;r}}_s}{{r_s}^3}+a_s---\left(25\right)$

$\frac{d^2{r}_c}{{\mathrm{dt}}^2}=-\frac{{\mathrm{\&mu;r}}_c}{{r_s}^3}+a_c---\left(26\right)$

其中a_s和a_c分别为目标星和追踪星除地球中心重力以外的其他作用力的合力的加速度矢量，即对推力和摄动力(包括地球形状摄动，大气阻力摄动和太阳光压摄动)的加速度矢量。 

由式(24)、式(25)和式(26)可得追踪星与目标星的绝对加速度之差

为 

$\frac{d^2\mathrm{\&rho;}}{{\mathrm{dt}}^2}=\frac{d^2{r}_c}{{\mathrm{dt}}^2}=\frac{d^2{r}_s}{{\mathrm{dt}}^2}=-\frac{{\mathrm{\&mu;r}}_c}{{r_c}^3}+\frac{{\mathrm{\&mu;r}}_s}{{r_s}^3}+a_c-a_s---\left(27\right)$

上式可进一步表示为下列等效的关系式 

$\frac{d^2\mathrm{\&rho;}}{{\mathrm{dt}}^2}=\frac{\mathrm{\&mu;}}{{r_s}^3}[r_s-{\left(\frac{r_s}{r_c}\right)}^3{r}_c]+\mathrm{\&Delta;a}---\left(28\right)$

为了建立追踪星与目标星在动坐标系s-xyz中的相对运动方程，有 

$\frac{d^2\mathrm{\&rho;}}{{\mathrm{dt}}^2}=\frac{{\mathrm{\&delta;}}^2\mathrm{\&rho;}}{{\mathrm{\&delta;t}}^2}+2\stackrel{\&CenterDot;}{\mathrm{\&theta;}}\&times;v+\stackrel{\&CenterDot;}{\mathrm{\&theta;}}\&times;(\stackrel{\&CenterDot;}{\mathrm{\&theta;}}\&times;\mathrm{\&rho;})+\stackrel{\&CenterDot;\&CenterDot;}{\mathrm{\&theta;}}\&times;\mathrm{\&rho;}---\left(29\right)$

上式中的

和v分别为追踪星在目标星轨道坐标系中的相对加速度矢量和相对速度矢量，则 

$\frac{{\mathrm{\&delta;}}^2\mathrm{\&rho;}}{{\mathrm{\&delta;t}}^2}=\left[\begin{array}{c}\stackrel{\&CenterDot;\&CenterDot;}{x}\\ \stackrel{\&CenterDot;\&CenterDot;}{y}\\ \stackrel{\&CenterDot;\&CenterDot;}{z}\end{array}\right]---\left(30\right)$

$v=\left[\begin{array}{c}\stackrel{\&CenterDot;}{x}\\ \stackrel{\&CenterDot;}{y}\\ \stackrel{\&CenterDot;}{z}\end{array}\right]---\left(31\right)$

式(29)中的

分别为轨道参考坐标系旋转的角加速度矢量和角速度矢量。目标星的平均运动角速度为

由轨道要素描述的卫星位置和轨道角速度的形式如下 

$\left|r_s\right|=\frac{a(1-e^2)}{1+e\cos \mathrm{\&theta;}}---\left(32\right)$

$\stackrel{\&CenterDot;}{\mathrm{\&theta;}}=\frac{n{(1+e\cos \mathrm{\&theta;})}^2}{{(1-e^2)}^{\frac{3}{2}}}---\left(33\right)$

针对圆轨道e=0，由式(33)得并可得下列近似式 

$\left\{\begin{array}{c}\stackrel{\&CenterDot;}{n}=0\\ n=\left[\begin{array}{c}0\\ 0\\ n\end{array}\right]\end{array}\right.---\left(34\right)$

对于追踪星和目标星的近距离的相对运动情况，卫星间距ρ是小量，特别是相对高度不大（即r_s/r_c接近1）时，从式(28)可以看出，其右端第一项为中心引力加速度差表达式，通过取一次近似（即线性化）可以进行简化。下面为具体的简化过程，可以看出简化实质是对中心引力取一次近似（即线性化）。 

因为 

$r_c={{(x+r_s)}^2+y^2+z^2]}^{\frac{1}{2}}={({\mathrm{\&rho;}}^2+{r_s}^2+{2\mathrm{xr}}_s)}^{\frac{1}{2}}---\left(35\right)$

所以有 

${\left(\frac{r_s}{r_c}\right)}^3{=[1+{\left(\frac{\mathrm{\&rho;}}{r_s}\right)}^2+2\frac{x}{r_s}]}^{-\frac{3}{2}}---\left(36\right)$

在上式忽略

以及更高的幂次项，则近似可得一下公式 

${\left(\frac{r_s}{r_c}\right)}^3\&ap;1-3\frac{x}{r_s}---\left(37\right)$

将上式(37)和式(24)代入式(28)得 

$\frac{d^2\mathrm{\&rho;}}{{\mathrm{dt}}^2}=\frac{\mathrm{\&mu;}}{{r_s}^3}[r_s-(1-3\frac{x}{r_c}(r_s+\mathrm{\&rho;}))]+\mathrm{\&Delta;a}---\left(38\right)$

进一步简化上式，忽略小量xρ项，得 

将式(39)、式(30)、式(31)和式(34)代入式(29)可得 

$\left\{\begin{array}{c}\stackrel{\&CenterDot;\&CenterDot;}{x}-2n\stackrel{\&CenterDot;}{y}-3n^{2}x=a_x\\ \stackrel{\&CenterDot;\&CenterDot;}{y}+2n\stackrel{\&CenterDot;}{x}=a_y\\ \stackrel{\&CenterDot;\&CenterDot;}{z}+n^{2}z=a_z\end{array}\right.---\left(40\right)$

通过简化，将相对运动动力学方程化为一组常系数线性微分方程。式(40)称为hill方程，也称Clohessey-Whiltshire方程(简称C-W方程)。 

9.2.2关于xy方向和z方向的独立分析 

对于hill方程，很容易得出，其xy方向和z方向的运动是相互独立的。因此可以在初始时刻z方向速度

很小情况下，选定一个合适的z₀。又由于在圆锥台内的掠飞时间相对于目标星的GEO轨道周期而言比较小，因此可以认为在该过程中z方向不发生很大变化，所以可以把三维掠飞问题，近似按照二维问题去处理，即只考虑xy方向的相对掠飞轨迹。 

9.3关于掠飞过程中轨迹线性化近似与“主飘方向”分析法 

对于hill方程，若初始时刻的信息

给定，可写出某时刻的位置和速度信息，即 

$r\left(t\right)={\mathrm{\&Phi;}}_{11}r\left(0\right)+{\mathrm{\&Phi;}}_{12}\stackrel{\&CenterDot;}{r}\left(0\right)---\left(41\right)$

$\stackrel{\&CenterDot;}{r}={\mathrm{\&Phi;}}_{21}r\left(0\right)+{\mathrm{\&Phi;}}_{22}\stackrel{\&CenterDot;}{r}\left(0\right)---\left(42\right)$

其中： 

${\mathrm{\&Phi;}}_{11}=\left[\begin{array}{ccc}4-3\cos \left(\mathrm{nt}\right)& 0& 0\\ -6\mathrm{nt}+6\sin \left(\mathrm{nt}\right)& 1& 0\\ 0& 0& \cos \left(\mathrm{nt}\right)\end{array}\right]---\left(43\right)$

${\mathrm{\&Phi;}}_{12}=\left[\begin{array}{ccc}\frac{\sin \left(\mathrm{nt}\right)}{n}& \frac{2}{n}(1-\cos \left(\mathrm{nt}\right))& 0\\ -\frac{2}{n}(1-\cos \left(\mathrm{nt}\right))& -3t+\frac{4}{n}\sin \left(\mathrm{nt}\right)& 0\\ 0& 0& \frac{\sin \left(\mathrm{nt}\right)}{n}\end{array}\right]---\left(44\right)$

${\mathrm{\&Phi;}}_{21}=\left[\begin{array}{ccc}3n\sin \left(\mathrm{nt}\right)& 0& 0\\ -6n(1-\cos \left(\mathrm{nt}\right))& 0& 0\\ 0& 0& -n\sin \left(\mathrm{nt}\right)\end{array}\right]---\left(45\right)$

${\mathrm{\&Phi;}}_{22}=\left[\begin{array}{ccc}\cos \left(\mathrm{nt}\right)& 2\sin \left(\mathrm{nt}\right)& 0\\ -2\sin \left(\mathrm{nt}\right)& -3+4\cos \left(\mathrm{nt}\right)& 0\\ 0& 0& \cos \left(\mathrm{nt}\right)\end{array}\right]---\left(46\right)$

对于本发明申请，重点是分析相对轨迹r在允许区域内掠飞时间尽可能长，因此只对式(41)进行研究。注意到用上述表示方法分析，很容易得到r和

&以GEO轨道周期T进行周期变化的特性；又由于相对掠飞时间会远小于GEO轨道周期，因此可以认为nt对应的是小角度，所以可以对三角函数进行小角度近似化简处理，即 

sin(nt)≈nt 

                                  (47) 

cos(nt)=1 

对于式(43)和式(44)，带入式(47)得到： 

${\mathrm{\&Phi;}}_{11}\&ap;\left[\begin{array}{ccc}1& 0& 0\\ 0& 1& 0\\ 0& 0& 1\end{array}\right]---\left(48\right)$

${\mathrm{\&Phi;}}_{12}=\left[\begin{array}{ccc}t& 0& 0\\ 0& t& 0\\ 0& 0& t\end{array}\right]---\left(49\right)$

因此在小角度近似下，可以得到hill方程相对坐标r的三个分量的线性变化结论，即设r=[x;y;z]，r₀=[x_o;y₀;z₀]，

满足 

$\begin{array}{c}x=x_0+{\stackrel{\&CenterDot;}{x}}_{0}t\\ y=y_0+{\stackrel{\&CenterDot;}{y}}_{0}t\\ z=z_0+{\stackrel{\&CenterDot;}{z}}_{0}t\end{array}---\left(50\right)$

可见给定合适的进入位置r₀和离开位置r_f，并选取合适的

可以保证足够长的掠飞时间。 

注1：尽管从式(50)可以看出只要

的各分量选取得足够小，可以让掠飞时间达到很大，但是如果时间t过大，将破坏nt的小角度假设前提条件，因此这种基于线性化的分析方法将不再成立。所以建议合理选取

注2：用这种线性化分析方法设计出的掠飞轨迹，完全可以达到在规定的圆锥台范围内掠飞任务要求的时间长度。 

从式(50)可以看出，如果要求三轴分量同时满足保持相同掠飞时间的要求，对初始位置信息

的要求将过于严格。由于xy平面和z方向的运动是独立的，因此可以简单进行如下设计，即x或y一个方向的分量的掠飞时间满足相应的时间要求，另外两个方向的分量的掠飞时间只要长于该方向的时间即可。这种分析方法称之为“主飘方向”分析法。 

9.4关于从hill方程坐标系到两体模型轨道坐标系的转换 

9.4.1转换矩阵的求解 

在上述章节所设计的掠飞方案中，是按照目标星轨道为圆轨道的hill方程坐标系下进行设计的，并且没有考虑摄动的影响。因此将在考虑了摄动的轨道坐标系对算法进行分析。但由于GEO轨道是高轨（几万公里），所以摄动项相对较小，不会对掠飞过程造成太大的影响。 

由于高轨时摄动项影响较小，因此仍可认为目标星是处于一个圆轨道上运行，即偏心率e≈0。这样hill方程的三轴的主方向仍可以按照 

X：由地心指向目标星方向； 

Y：目标星速度方向； 

Z：右手定则。 

来确定。且认为三轴是垂直的。在本课题的所有仿真，都是基于上述三轴原则选取任意t时刻的hill方程三轴单位方向，同时也相当于按该原则确定了hill方程坐标系和惯性坐标系之间的余弦变换矩阵。 

假设在任意时刻t，目标星的位置和速度信息在惯性系表示分别为r_t，v_t，则有： 

$e_x=\frac{r_t}{\left|r_t\right|},e_y=\frac{v_t}{\left|v_t\right|},e_z=e_x\&times;e_y$

则从惯性系到hill方程坐标系的余弦矩阵为： 

$\mathrm{DCMTi}=\left[\begin{array}{c}{e}_x^T\\ e_y^T\\ e_z^T\end{array}\right]---\left(51\right)$

9.4.2关于追踪星位置和速度及速度脉冲的求解说明 

在已知某时刻的目标星的相对惯性系的位置和速度后，根据两星之间的相对位置和速度关系，并通过相关余弦矩阵，可以得到同一时刻追踪星的相对惯性系的位置和速度信息。 

9.4.2.1位置求解 

对于追踪星的位置信息求法，如图17所示，其中△r在hill方程坐标系下的表示就是hill方程中的xyz分量表示的向量△r_h，在惯性系下的表示为： 

△r=DCMTi^T△r_h    (52) 

9.4.2.2速度求解 

对于追踪星在惯性坐标系下的速度信息求法，根据hill方程的定义，相比于位置信息，要多考虑苛氏加速度的影响，即 

$\left[\begin{array}{c}\mathrm{\&Delta;}\frac{\mathrm{dx}}{\mathrm{dt}}\\ \mathrm{\&Delta;}\frac{\mathrm{dy}}{\mathrm{dt}}\\ \mathrm{\&Delta;}\frac{\mathrm{dz}}{\mathrm{dt}}\end{array}\right]={\mathrm{DCMTi}}^T\{\left[\begin{array}{c}\stackrel{\&CenterDot;}{x}\\ \stackrel{\&CenterDot;}{y}\\ \stackrel{\&CenterDot;}{z}\end{array}\right]+\left[\begin{array}{c}0\\ 0\\ n\end{array}\right]\&times;\left[\begin{array}{c}x\\ y\\ z\end{array}\right]\}---\left(53\right)$

其中

是追踪星在惯性坐标系下的相对目标星的绝对速度，[x,y,z]^T和

是hill方程下的分量，n为目标星轨道速度。 

9.4.2.3速度脉冲求解 

对于本任务的逼近过程，是采用若干次速度脉冲增量的方法实现的，因此有必要给出追踪星的瞬时速度增量在惯性系下的表示。 

由于是若干次速度脉冲增量实现的逼近过程，并且考虑了燃料最优的影响，因此每次的速度增量不是很大，所以本任务认为每次的速度脉冲增量是瞬时完成的。因此假设在每次速度脉冲前时间为t^-，加脉冲后时间为t⁺。有 

$\left[\begin{array}{c}\mathrm{\&Delta;}\frac{\mathrm{dx}}{{\mathrm{dt}}^{-}}\\ \mathrm{\&Delta;}\frac{\mathrm{dy}}{{\mathrm{dt}}^{-}}\\ \mathrm{\&Delta;}\frac{\mathrm{dz}}{{\mathrm{dt}}^{-}}\end{array}\right]={\mathrm{DCMTi}}^T\{{\left[\begin{array}{c}\stackrel{\&CenterDot;}{x}\\ \stackrel{\&CenterDot;}{y}\\ \stackrel{\&CenterDot;}{z}\end{array}\right]}_{t^{-}}+\left[\begin{array}{c}0\\ 0\\ n\end{array}\right]\&times;\left[\begin{array}{c}x\\ y\\ z\end{array}\right]\}---\left(54\right)$

和 

$\left[\begin{array}{c}\mathrm{\&Delta;}\frac{\mathrm{dx}}{{\mathrm{dt}}^{+}}\\ \mathrm{\&Delta;}\frac{\mathrm{dy}}{{\mathrm{dt}}^{+}}\\ \mathrm{\&Delta;}\frac{\mathrm{dz}}{{\mathrm{dt}}^{+}}\end{array}\right]={\mathrm{DCMTi}}^T\{{\left[\begin{array}{c}\stackrel{\&CenterDot;}{x}\\ \stackrel{\&CenterDot;}{y}\\ \stackrel{\&CenterDot;}{z}\end{array}\right]}_{t^{+}}+\left[\begin{array}{c}0\\ 0\\ n\end{array}\right]\&times;\left[\begin{array}{c}x\\ y\\ z\end{array}\right]\}---\left(55\right)$

则追踪星t时刻的变轨所需要的相对惯性系的速度脉冲为： 

$\begin{array}{c}\mathrm{\&Delta;v}=\left[\begin{array}{c}\mathrm{\&Delta;}\frac{\mathrm{dx}}{{\mathrm{dt}}^{+}}\\ \mathrm{\&Delta;}\frac{\mathrm{dy}}{{\mathrm{dt}}^{+}}\\ \mathrm{\&Delta;}\frac{\mathrm{dz}}{{\mathrm{dt}}^{+}}\end{array}\right]-\left[\begin{array}{c}\mathrm{\&Delta;}\frac{\mathrm{dx}}{{\mathrm{dt}}^{-}}\\ \mathrm{\&Delta;}\frac{\mathrm{dy}}{{\mathrm{dt}}^{-}}\\ \mathrm{\&Delta;}\frac{\mathrm{dz}}{{\mathrm{dt}}^{-}}\end{array}\right]={\mathrm{DCMTi}}^T\{{\left[\begin{array}{c}\stackrel{\&CenterDot;}{x}\\ \stackrel{\&CenterDot;}{y}\\ \stackrel{\&CenterDot;}{z}\end{array}\right]}_{t+}-{\left[\begin{array}{c}\stackrel{\&CenterDot;}{x}\\ \stackrel{\&CenterDot;}{y}\\ \stackrel{\&CenterDot;}{z}\end{array}\right]}_{t^{-}}\}\\ ={\mathrm{DCMTi}}^T{\mathrm{\&Delta;v}}_h\end{array}---\left(56\right)$

其中△v_h为hill方程系下的速度增量。 

9.5关于hill方程设计到两体问题设计思路的合理性说明 

对于该任务，目标星处于GEO轨道，平均轨道高度达到4万多公里，显然追踪星与目标星的距离的上下边界一般都远远小于该轨道高度；而一般轨迹转移的初始点与期望掠飞进入点的距离甚至远远小于掠飞区域内与目标星允许距离的下边界，因此更小于GEO轨道高度，因此保证了hill方程的成立条件和用于相对轨迹转移、逼近和掠飞的有效性。 

（十）与现有技术方案的比较 

如果要实现这一类任务的要求，除了本发明算法之外，从hill方程解的角度出发，结合传统方法的思想，可能还有一些悬停、伴飞和绕飞的方案，以下简单介绍几种并将其与本发明算法进行比较。 

10.1基于hill方程在无轨道控制下的方案 

hill方程在a_x=a_y=a_z=0时，由式(43)和式(44)所表示的r_t与r₀、

的关系，经过 适当的数学变换消去方程中的时间参数t，则可以将xy方向的运动表示成式(57)所示的椭圆方程： 

$\frac{{(x-x_{c0})}^2}{b^2}+\frac{{(y-y_{c0}+\frac{3}{2}{x}_{c0}\mathrm{nt})}^2}{{\left(2b\right)}^2}=1---\left(57\right)$

式中 

$x_{c0}={4x}_0+2\frac{{\stackrel{\&CenterDot;}{y}}_0}{n}---\left(58\right)$

$y_{c0}=y_0-2\frac{{\stackrel{\&CenterDot;}{x}}_0}{n}---\left(59\right)$

$b=\sqrt{{(\frac{2{\stackrel{\&CenterDot;}{y}}_0}{n}+{3x}_0)}^2+{\left(\frac{{\stackrel{\&CenterDot;}{x}}_0}{n}\right)}^2}---\left(60\right)$

式(57)描述的椭圆由上述3个参数确定。其中x_c0、y_c0与椭圆的中心有关，b决定了椭圆的长、短半轴，通过对这3个参数分析可以得到相应的运动规律。 

z方向的运动是由z₀、

所确定的振荡运动。 

10.1.1无轨道控制下的非闭合轨迹飞行 

当不满足x_c0=0时，xy方向运动的椭圆中心将沿着y方向飘移，形成非闭合轨迹的螺旋运动，且飘移速度与x_c0成正比。 

若用该方案来实现任务所描述的要求，计算复杂，难以确定x₀、y₀、z₀、

六个参数来满足在连续时间内保持在规定区域内飞行的要求。 

10.1.2无轨道控制下的伴随飞行 

当满足x_c0=0时，xy方向构成封闭的椭圆轨迹，且椭圆的长半轴为短半轴的2倍（即轨道平面内运动具有固定的偏心率

椭圆的中心位于y轴上，椭圆运动的周期为

此时可以导出伴随飞行的必要条件 

${\stackrel{\&CenterDot;}{y}}_0=-{2\mathrm{nx}}_0---\left(61\right)$

若用该方案来实现任务所描述的要求，希望能让整个伴随飞行的轨迹都在规定区域内，从而可以达到要求，但是z方向振荡运动的中心为z=0，无法满足任务中视线方向的任意性这一要求。 

10.1.3无轨道控制下的悬停 

在满足式(61)的基础上，若再满足b=0，则椭圆退化为y轴上的一个固定点y=y₀，此时可以导出相应的约束条件为x₀、

均为0，考虑z方向的运动可知，要想实现在某点的悬停，还要满足z₀、

为0，这种对初始状态苛刻的约束条件，同样无法满足任务中视线方向的任意性这一要求。 

10.1.4无轨道控制下的绕飞 

在满足式(61)的基础上，若进一步满足y_c0=0，则封闭椭圆的中心将位于相对运动坐标系的原点，这表明追踪星将环绕目标星形成环绕飞行，可以导出绕飞的必要条件 

$\left\{\begin{array}{c}{y}_0=2\frac{{\stackrel{\&CenterDot;}{x}}_0}{n}\\ {\stackrel{\&CenterDot;}{y}}_0={-2\mathrm{nx}}_0\end{array}\right.---\left(62\right)$

若用该方案来实现任务所描述的要求，希望绕飞轨迹穿过允许掠飞的区域，且时间满足要求。由于绕飞运动的周期就是目标星的飞行周期为

以目标星视场角范围±θ=±0.35°为例，可以粗略地计算绕飞轨迹在规定区域内的时间 

这个固定的时间不一定满足任务要求的几百秒需求。 

10.2脉冲控制下的“水滴”形轨迹飞行 

由hill方程得到的非闭合轨迹在xy面内投影的封闭部分像“水滴”而得名，如图18所示，可以通过周期性地在x坐标最大处施加xy面内的速度脉冲从而保持xy面投影的“水滴”形状，由于xy方向运动的周期

而z方向的运动周期

则最终闭合轨迹的周期T为T₀与T'的最小公倍数。 

若用该方案来实现任务所描述的要求，由于轨迹在xy面的投影始终为固定的一个“水滴”形状，所以可以根据规定的空间区域确定参数。由于任务要求的时间与目标星的运动周期相比很小，z方向的位置变化也很小，可以很好地使z方向的运动在任务要求的时间内保持在该空间区域。但由于任务所规定的空间区域较小，则T'较小，在任务要求的时间内需要多次周期性的xy面内速度脉冲，而任务要求在该段时间内会有类似实现 高精度的姿态指向等需求，每次加速度脉冲时尽管时间很短，但姿轨控耦合问题势必要影响姿态指向的精度，所以该方案存在一定的缺陷。 

10.3连续控制下的任意位置悬停 

考虑追踪星B在目标星A任意相对位置悬停，在轨道坐标系A-xyz中，设追踪星相对目标星的悬停方位角为α，高低角为β，相对距离为r，定义方法如图19所示。设目标星A在圆轨道上运行，速度为v₀，追踪星B的速度为v₁，显然追踪星也为圆轨道，只是圆心不一定在地心，显然追踪星B的轨道平面与目标星A的轨道平面之间的距离d保持不变，可以表示为 

d=AA'=OO'=rsinαcosβ=c nst    (64) 

根据图19，追踪星B的轨道中心为O’，有OO’平行与z轴，由于轨道平面O’A’B与轨道平面OAy平行，因此若想悬停，两航天器的速度必须满足 

$\left\{\begin{array}{c}({\stackrel{\&RightArrow;}{r}}_0\&times;{\stackrel{\&RightArrow;}{v}}_0)\&times;({\stackrel{\&RightArrow;}{\mathrm{\&rho;}}}_1\&times;{\stackrel{\&RightArrow;}{v}}_1)=0\\ \frac{v_0}{r_0}=\frac{v_1}{{\mathrm{\&rho;}}_1}={\mathrm{\&omega;}}_0=\sqrt{\frac{\mathrm{\&mu;}}{{r_0}^3}}\end{array}\right.---\left(65\right)$

式(65)表明两轨道平面法向量重合，两航天器角速度始终保持相等，ρ₁的定义见图19、图20，μ为地球引力常数。由于

因此有

矢量之间的夹角应该等于∠BO'A'，即 

$\mathrm{\&eta;}=\&angle;{\mathrm{BO}}^{\&prime;}{A}^{\&prime;}=\arccos \left(\frac{{\stackrel{\&RightArrow;}{v}}_0\&CenterDot;{\stackrel{\&RightArrow;}{v}}_1}{v_0{v}_1}\right)=\arctan \left(\frac{r\cos \mathrm{\&alpha;}\cos \mathrm{\&beta;}}{r_0+\sin \mathrm{\&beta;}}\right)---\left(66\right)$

追踪星B的轨道半径ρ₁可以表示为 

${\mathrm{\&rho;}}_1=\sqrt{r^2{\cos }^2\mathrm{\&beta;}{\cos }^2\mathrm{\&alpha;}+{(r_0+r\sin \mathrm{\&beta;})}^2}---\left(67\right)$

追踪星的速度大小为 

$v_1={\mathrm{\&omega;}}_0{\mathrm{\&rho;}}_1\sqrt{\frac{r^2{\sin }^2\mathrm{\&beta;}+r^2{\cos }^2\mathrm{\&alpha;}{\cos }^2\mathrm{\&beta;}+{r_0}^2+{2\mathrm{rr}}_0\sin \mathrm{\&beta;}}{{r_0}^3}\mathrm{\&mu;}}---\left(68\right)$

追踪星的向心加速度g^*如图20所示，由B指向O’，其大小可以表示为 

$g^{*}={{\mathrm{\&omega;}}_0}^2{\mathrm{\&rho;}}_1=\frac{\mathrm{\&mu;}}{{r_0}^3}\sqrt{r^2{\sin }^2\mathrm{\&beta;}+r^2{\cos }^2\mathrm{\&alpha;}{\cos }^2\mathrm{\&beta;}+{r_0}^2+{2\mathrm{rr}}_0\sin \mathrm{\&beta;}}---\left(69\right)$

追踪星本身受到地球万有引力作用，固有引力加速度g₁可以表示为 

$g_1=\frac{\mathrm{\&mu;}}{{r_1}^2}=\frac{\mathrm{\&mu;}}{{(r\sin \mathrm{\&beta;}+r_0)}^2+r^2{\cos }^2\mathrm{\&beta;}}---\left(70\right)$

则需要的控制加速度a可以表示为 

a=g^*-g₁    (71) 

同时考虑g₁与g^*之间的夹角σ为 

$\mathrm{\&sigma;}=\arctan \frac{d}{{\mathrm{\&rho;}}_1}=\arctan \frac{r\sin \mathrm{\&alpha;}\cos \mathrm{\&beta;}}{\sqrt{r^2{\cos }^2\mathrm{\&beta;}{\cos }^2\mathrm{\&alpha;}+{(r_0+r\sin \mathrm{\&beta;})}^2}}---\left(72\right)$

则a在追踪星径向和侧向需要施加的控制力加速度a_x、a_z可以分别表示为 

$\left\{\begin{array}{c}{a}_x=g^{*}\cos \mathrm{\&sigma;}-g_1\\ a_z=g^{*}\sin \mathrm{\&sigma;}\end{array}\right.---\left(73\right)$

若用该方案来实现任务所描述的要求，先将追踪星机动到悬停位置，然后施加速度脉冲，使追踪星的速度大小和方向满足式(66)、式(68)的要求，再对追踪星施加连续控制力，大小和方向满足式(71)、式(72)和式(73)。毫无疑问，这个悬停的方案一定能满足任务中规定区域和时间的要求，但由于要一直施加连续控制力，使姿轨控耦合的影响导致高精度姿态指向等任务需求难以实现，而且一直相对非合作目标停留在一点，非常容易引起对方的察觉，从而暴露身份，所以这个方案也有不足。 

10.4关于“主飘方向”分析法的优点说明 

主飘方向的本质就是在相对轨道坐标系x或y方向投影较短的方向，采用的是主要矛盾分析法。在找到主飘方向的期望起点与终点情况下，这种分析法将一个相对hill坐标系的三维掠飞问题，化成了一个一维初速度值求解问题，进而再通过一些简单的原则化成一个三维初速度值确定问题，分析思路简单。 

与无轨道控制的非闭合轨迹飞行方案相比，该发明算法具有参数设计简单的优点；与无轨道控制的伴随飞行和悬停方案相比，具有视场方向可以任意的优点；与无轨道控制的绕飞方案相比，具有掠飞时间可设计的优点；与脉冲控制的“水滴”形轨迹飞行方案相比，具有掠飞区域中完全不进行轨道控制，以实现足够长时间内完成高精度姿态指向控制等任务需求的优点；与连续控制的悬停方案相比，除了避免姿轨耦合控制，利于提高姿态指向控制等任务精度外，还具有不易暴露身份的优点。 

具体实施例：下面是一个算法仿真验证实验。 

（一）仿真参数 

1.1目标星轨道参数 

目标星处于GEO轨道，初始轨道六根数为：a₁=4.225×10⁷m，e₁=0，i₁=5°，Ω₁=31°， ω₁=0°，M0₁=0°。 

卫星摄动系数参数 

气动系数C_D=1，阻力系数C_d=2.2，反射系数C_r=0.8。 

1.2掠飞参数 

目标星视线方向矢量e_-hill=[-49647;-5194;2863]，追踪星在目标星视线方向最小距离R_min=50km，最远距离取R_max=99.2km，视线角θ=0.35°。预期掠飞时间按800s左右设计。 

（二）掠飞轨迹设计 

首先按0节设计思路设计给出主飘方向为y方向。关于初始位置和速度的确定，按照0、0节的原则，计算得到进入点Y0位置坐标为[-9.8437;-1.0908;0.5680]×104m，初始的速度为：[0.200;1.5222;-0.100]m/s。相对惯性系的绝对速度脉冲增量按照0节确定，假设通过轨道机动控制，使追踪星到达掠飞进入点时满足设计的初始速度，仿真观察掠飞段轨迹。 

（三）基于开普勒方程以及hill方程的仿真验证 

给出基于开普勒轨道方程以及hill方程的算法有效性仿真验证。两星在开普勒两体轨道模型基础上考虑多种轨道摄动的影响，仿真结果是将追踪星相对于目标星的位置矢量实时地转换回目标星的hill方程坐标系，并且与不考虑轨道摄动，且在hill动力学方程下进行仿真的情况进行了对比。 

相对轨迹在目标星轨道坐标系下的三维曲线如图21所示，表明两种坐标系下进行掠飞仿真得到的仿真曲线基本可以重合；三维曲线局部放大图如图22所示，表明由于摄动的存在，惯性系下仿真得到的轨线和hill系下不考虑摄动仿真得到的轨线存在一定的差别，但差别较小。根据图23看出，掠飞轨迹的y和z轴分量在掠飞过程中线性化程度保持较好。图24是在允许圆锥台区域内掠飞段的时间间隔示意图，在该组初始相对参数下，追踪星可以在目标星视场内掠飞775秒，与“主飘方向”分析法的分析结果基本吻合。 

Claims (6)

1.航天器相对轨道的控制方法，其特征在于它通过下述步骤实现：一、追踪航天器进入目标航天器的视线角范围内而且追踪航天器进入二者之间确定的距离范围内；二、计算并确定追踪航天器期望轨道的起点、末点和初始入轨速度，并确定主飘方向；三、追踪航天器在期望轨道的起点，以上述计算并确定的初始入轨速度进入轨道，并在期望轨道的末点脱离轨道，从而完成对目标航天器的掠飞动作。

2.根据权利要求1所述的航天器相对轨道的控制方法，其特征在于步骤二中主飘方向是从Y₀点到Y_f的轨迹；各点的hill方程下的坐标表示计算如下：

Y₀和Y_f的x坐标相同，为

x_f=X₀-Rsinβ    (1)

其中X₀为圆锥体底面中心x坐标，R为底面圆半径，β的计算满足：

$\mathrm{\&beta;}=\arctan \left(\frac{e_2}{e_1}\right)---\left(2\right)$

由于Y₀和Y_f两点都在圆锥台的表面上，因此其与原点的连线与目标航天器视线方向l成θ，通过求解如下方程即可确定两点的y方向坐标值，即

$\mathrm{ar}\cos ((x_f{e}_1+y_{o,f}{e}_2+z_0{z}_3)/\sqrt{x_f^2+y_{0,f}^2+z_0^2})=\mathrm{\&theta;}---\left(3\right)$

该方程有两个解，其中y_o,f为两点的y方向坐标。

3.根据权利要求2所述的航天器相对轨道的控制方法，其特征在于步骤二中设期望掠飞时间是t_lve，则主飘方向的期望速度分量计算为：

${\stackrel{.}{y}}_0=\frac{|Y_f-Y_0|}{t_{\mathrm{lve}}}---\left(4\right).$

4.根据权利要求1所述的航天器相对轨道的控制方法，其特征在于步骤二中从观测的目标航天器位置和速度数据出发，计算航天器期望轨道；轨道计算如下：

求半长轴a：

$a=\frac{\mathrm{\&mu;r}}{2\mathrm{\&mu;}-{\mathrm{rv}}^2}---\left(5\right)$

其中v为航天器轨道速度；

求偏心率e：

$e=\sqrt{1-p/a}---\left(6\right)$

其中p为椭圆半通径。

轨道倾角i与升交点赤经Ω计算如下：

i=arccos(h_z/h)    (7)

Ω=-arctan(2h_x/h_y)    (8)

求真近点角f：

cosf=(p/r-1)/e    (9)

求近地点幅角ω：

令i_Ω为地心指向升交点的单位矢量，在地心惯性系中表示为i_Ω=cosΩ·x_i+sinΩ·y_i。故可由下式确定轨道角u

$\cos u=\frac{1}{r}{r}^T{i}_{\mathrm{\&Omega;}}---\left(10\right)$

进而可得近地点幅角为

ω=u-f    (11)

当e<10^-6时，轨道称为圆轨道，有ω=0，f=u。

5.根据权利要求1所述的航天器相对轨道的控制方法，其特征在于在步骤二中还计算出任意时刻航天器的位置矢量和速度矢量，称为星历表计算；有

$\begin{array}{c}r=a(\cos E-e)\&CenterDot;P+a\sqrt{1-e^2}\sin E\&CenterDot;Q\\ v=-\sqrt{\mathrm{\&mu;}/p}\sin f\&CenterDot;P+\sqrt{\mathrm{\&mu;}/p}(e+\cos f)\&CenterDot;Q\end{array}---\left(12\right)$

其中，P和Q分别表示卫星近地点和半通径方向的单位矢量；通过旋转，P和Q的表达式如下所示：

$P=R_3(-\mathrm{\&Omega;})R_1(-i)R_3(-\mathrm{\&omega;})\left[\begin{array}{c}1\\ 0\\ 0\end{array}\right]\left[\begin{array}{c}\cos \mathrm{\&Omega;}\cos \mathrm{\&omega;}-\sin \mathrm{\&Omega;}\sin \mathrm{\&omega;}\cos i\\ \sin \mathrm{\&Omega;}\cos \mathrm{\&omega;}+\cos \mathrm{\&Omega;}\sin \mathrm{\&omega;}\cos i\\ \sin \mathrm{\&omega;}\sin i\end{array}\right]---\left(13\right)$

$Q=R_3(-\mathrm{\&Omega;})R_1(-i)R_3(-\mathrm{\&omega;})\left[\begin{array}{c}0\\ 1\\ 0\end{array}\right]\left[\begin{array}{c}-\cos \mathrm{\&Omega;}\sin \mathrm{\&omega;}-\sin \mathrm{\&Omega;}\cos \mathrm{\&omega;}\cos i\\ -\sin \mathrm{\&Omega;}\sin \mathrm{\&omega;}+\cos \mathrm{\&Omega;}\cos \mathrm{\&omega;}\cos i\\ \cos \mathrm{\&omega;}\sin i\end{array}\right]---\left(14\right)$

。

6.根据权利要求1所述的航天器相对轨道的控制方法，其特征在于步骤二中还包括下述步骤：在已知某时刻的目标星的相对惯性系的位置和速度后，根据两星之间的相对位置和速度关系，并通过相关余弦矩阵，可以得到同一时刻追踪星的相对惯性系的位置和速度信息；

对于追踪星的位置信息，△r在hill方程坐标系下的表示就是hill方程中的xyz分量表示的向量△r_h，在惯性系下的表示为：

△r=DCMTi^T△r_h    (15)

对于追踪星在惯性坐标系下的速度信息，根据hill方程的定义，相比于位置信息，要多考虑苛氏加速度的影响，即

$\left[\begin{array}{c}\mathrm{\&Delta;}\frac{\mathrm{dx}}{\mathrm{dt}}\\ \mathrm{\&Delta;}\frac{\mathrm{dy}}{\mathrm{dt}}\\ \mathrm{\&Delta;}\frac{\mathrm{dz}}{\mathrm{dt}}\end{array}\right]={\mathrm{DCMTi}}^T\{\left[\begin{array}{c}\stackrel{\&CenterDot;}{x}\\ \stackrel{\&CenterDot;}{y}\\ \stackrel{\&CenterDot;}{z}\end{array}\right]+\left[\begin{array}{c}0\\ 0\\ n\end{array}\right]\&times;\left[\begin{array}{c}x\\ y\\ z\end{array}\right]\}---\left(16\right)$

其中

是追踪星在惯性坐标系下的相对目标星的绝对速度，[x,y,z]^T和是hill方程下的分量，n为目标星轨道速度；

追踪星的瞬时速度增量在惯性系下的表示，由于是若干次速度脉冲增量实现的逼近过程，并且考虑了燃料最优的影响，因此每次的速度增量不是很大，每次的速度脉冲增量是瞬时完成的，假设在每次速度脉冲前时间为t^-，加脉冲后时间为t⁺。有

$\left[\begin{array}{c}\mathrm{\&Delta;}\frac{\mathrm{dx}}{{\mathrm{dt}}^{-}}\\ \mathrm{\&Delta;}\frac{\mathrm{dy}}{{\mathrm{dt}}^{-}}\\ \mathrm{\&Delta;}\frac{\mathrm{dz}}{{\mathrm{dt}}^{-}}\end{array}\right]={\mathrm{DCMTi}}^T\{{\left[\begin{array}{c}\stackrel{\&CenterDot;}{x}\\ \stackrel{\&CenterDot;}{y}\\ \stackrel{\&CenterDot;}{z}\end{array}\right]}_{t^{-}}+\left[\begin{array}{c}0\\ 0\\ n\end{array}\right]\&times;\left[\begin{array}{c}x\\ y\\ z\end{array}\right]\}---\left(17\right)$

和

$\left[\begin{array}{c}\mathrm{\&Delta;}\frac{\mathrm{dx}}{{\mathrm{dt}}^{+}}\\ \mathrm{\&Delta;}\frac{\mathrm{dy}}{{\mathrm{dt}}^{+}}\\ \mathrm{\&Delta;}\frac{\mathrm{dz}}{{\mathrm{dt}}^{+}}\end{array}\right]={\mathrm{DCMTi}}^T\{{\left[\begin{array}{c}\stackrel{\&CenterDot;}{x}\\ \stackrel{\&CenterDot;}{y}\\ \stackrel{\&CenterDot;}{z}\end{array}\right]}_{t^{+}}+\left[\begin{array}{c}0\\ 0\\ n\end{array}\right]\&times;\left[\begin{array}{c}x\\ y\\ z\end{array}\right]\}---\left(18\right)$

则追踪星t时刻的变轨所需要的相对惯性系的速度脉冲为：

$\begin{array}{c}\mathrm{\&Delta;v}=\left[\begin{array}{c}\mathrm{\&Delta;}\frac{\mathrm{dx}}{{\mathrm{dt}}^{+}}\\ \mathrm{\&Delta;}\frac{\mathrm{dy}}{{\mathrm{dt}}^{+}}\\ \mathrm{\&Delta;}\frac{\mathrm{dz}}{{\mathrm{dt}}^{+}}\end{array}\right]-\left[\begin{array}{c}\mathrm{\&Delta;}\frac{\mathrm{dx}}{{\mathrm{dt}}^{-}}\\ \mathrm{\&Delta;}\frac{\mathrm{dy}}{{\mathrm{dt}}^{-}}\\ \mathrm{\&Delta;}\frac{\mathrm{dz}}{{\mathrm{dt}}^{-}}\end{array}\right]={\mathrm{DCMTi}}^T\{{\left[\begin{array}{c}\stackrel{\&CenterDot;}{x}\\ \stackrel{\&CenterDot;}{y}\\ \stackrel{\&CenterDot;}{z}\end{array}\right]}_{t+}-{\left[\begin{array}{c}\stackrel{\&CenterDot;}{x}\\ \stackrel{\&CenterDot;}{y}\\ \stackrel{\&CenterDot;}{z}\end{array}\right]}_{t^{-}}\}\\ ={\mathrm{DCMTi}}^T{\mathrm{\&Delta;v}}_h\end{array}---\left(19\right)$

其中△v_h为hill方程系下的速度增量。

Images