应急遥感在国计民生、国防安全等领域发挥着重要作用，而传统单一卫星无法满足应急遥感全球覆盖、快速重访、多源信息获取的迫切需求。微纳卫星造价低廉，可批量生产，发射形式多样，因而基于微纳航天器的分布式遥感技术成为解决上述瓶颈的关键。对于分布式系统而言，其功能实现常依赖于系统的空间几何构型，构型不同功能也不同^[1]。如何在现有微纳航天器的技术基础上开展凝视、推扫、视频、多星组网成像等多种遥感模式的关联性研究，是中国“高性能空天一体化组网监测系统技术”的重要组成部分。

目前用于描述编队航天器之间的相互运动方法有2种^[2]：1)基于相对位置速度等运动状态的动力学方法，2)基于轨道要素差的运动学方法。前者主要有Hill模型^[3]及其改进型包括Tschauner-Hempel模型^[4]、Schweighar-Sedwick模型^[5]、Vadali-Alfriend-Vaddi模型^[6]等；Hill模型及其改进型在队形长期维持过程中需要消耗燃料，因而不适用于微纳卫星系统编队。后者主要有一阶轨道要素差法^[7]、J₂摄动不变求解法^[8]、半长轴修正法^[9]、相对漂移率最小优化求解法^[10]等方法；上述方法的缺点是优化求解过程复杂，难以满足应急遥感编队航天器数量众多、灵活多变，且工作模式不同则编队构型不同的需求。

针对多模式分布式遥感任务典型想定和燃料最省、求解高效等目标，本文研究了分布式遥感微纳航天器集群的自然编队构型设计。提出了一种直接参数求解方法，可以在2种典型任务模式下直接获取编队中各伴飞航天器的轨道参数，并能轻易地将编队航天器的数量拓展到N；STK仿真验证结果表明，编队在一个工作周期内队形保持良好，队形扩散满足任务需求。

1 多航天器编队自然演化模型 1.1 坐标系定义

为了描述参考航天器与伴飞航天器之间的相对运动，对编队坐标系进行定义如图 1所示。

图中OXYZ为地心坐标系，Sxyz为参考航天器轨道坐标系，其定义如下：坐标系原点S为参考航天器质心；x轴在参考航天器轨道面内，指向S点与地心连线的反方向；y轴在参考航天器轨道面内与x轴垂直，指向航天器运动方向；z轴垂直于轨道面，并与x、y轴成右手系。图中矢量r_r由地心O指向参考航天器质心S，矢量r_fk由地心O指向伴飞航天器质心F，矢量ρ由参考航天器质心S指向伴飞航天器质心F。下标r对应参考航天器，下标f对应伴飞航天器。

1.2 相对运动方程

令σ_r={a_r，e_r，i_r，Ω_r，ω_r，M_r}分别表示参考航天器的半长轴、偏心率、轨道倾角、升交点赤经、近地点幅角、平近点角；σ_fk={a_fk，e_fk，i_fk，Ω_fk，ω_fk，M_fk}则表示第k个伴飞航天器的轨道要素，k=1，2，…，N。由轨道根数差描述的一阶近似相对运动^[11]方程为

$ \begin{array}{l} {x_{{\rm{f}}k}} \approx \left\| {{\mathit{\boldsymbol{r}}_{{\rm{f}}k}}} \right\| - \left\| {{\mathit{\boldsymbol{r}}_{\rm{r}}}} \right\|,\\ {y_{{\rm{f}}k}} \approx \left\| {{\mathit{\boldsymbol{r}}_{{\rm{f}}k}}} \right\|\left( {\Delta {u_{{\rm{f}}k}} + \Delta {\mathit{\Omega }_{{\rm{f}}k}}\cos {i_{\rm{r}}}} \right),\\ {z_{{\rm{f}}k}} \approx \left\| {{\mathit{\boldsymbol{r}}_{{\rm{f}}k}}} \right\|\left( { - \Delta {\mathit{\Omega }_{{\rm{f}}k}}\sin {i_{\rm{r}}}\cos \Delta {u_{\rm{r}}} + \Delta {i_{{\rm{f}}k}}\cos {u_{\rm{r}}}} \right). \end{array} $

(1)

其中：ΔΩ_fk、Δu_fk、Δi_fk分别为参考航天器与第k个伴飞航天器之间的升交点赤经、纬度幅角、轨道倾角之差，‖r_fk‖为第k个伴飞航天器到地心的距离。

1.3 摄动影响下的相对运动方程

在自然编队条件下，系统的构型通常会受到J₂等摄动因素的影响。J₂摄动对航天器的轨道参数的影响可用如下方程描述^[8]：

$ \left\{ \begin{array}{l} {{\dot a}_{{\rm{f}}k}} = {{\dot e}_{{\rm{f}}k}} = {{\dot i}_{{\rm{f}}k}} = 0,\\ {{\mathit{\dot \Omega }}_{{\rm{f}}k}} = - \frac{3}{2}{J_2}\frac{{{\mu ^{0.5}}R_{\rm{e}}^2}}{{a_{{\rm{f}}k}^{3.5}{{\left( {1 - e_{{\rm{f}}k}^2} \right)}^2}}}\cos {i_{{\rm{f}}k}},\\ {{\dot \omega }_{{\rm{f}}k}} = \frac{3}{4}{J_2}\frac{{{\mu ^{0.5}}R_{\rm{e}}^2}}{{a_{{\rm{f}}k}^{3.5}{{\left( {1 - e_{{\rm{f}}k}^2} \right)}^2}}}\left( {5{{\cos }^2}{i_{{\rm{f}}k}} - 1} \right),\\ {{\dot M}_{{\rm{f}}k}} = \frac{{{\mu ^{0.5}}}}{{a_{{\rm{f}}k}^{1.5}}} + \frac{3}{4}{J_2}\frac{{{\mu ^{0.5}}R_{\rm{e}}^2}}{{a_{{\rm{f}}k}^{3.5}{{\left( {1 - e_{{\rm{f}}k}^2} \right)}^{1.5}}}}\left( {5{{\cos }^2}{i_{{\rm{f}}k}} - 1} \right). \end{array} \right. $

(2)

J₂长期项对应的是平轨道根数，由式(2)可知参考航天器与伴飞航天器轨道差随时间变化为

$ \left\{ \begin{array}{l} \Delta {a_{{\rm{f}}k}}\left( t \right) = \Delta {a_{{\rm{f}}k}},\\ \Delta {i_{{\rm{f}}k}}\left( t \right) = \Delta {i_{{\rm{f}}k}},\\ \Delta {e_{{\rm{f}}k}}\left( t \right) = \Delta {e_{{\rm{f}}k}},\\ \Delta {\mathit{\Omega }_{{\rm{f}}k}}\left( t \right) = \Delta {\mathit{\Omega }_{{\rm{f}}k}} + \Delta {{\mathit{\dot \Omega }}_{{\rm{f}}k}}t,\\ \Delta {\omega _{{\rm{f}}k}}\left( t \right) = \Delta {\omega _{{\rm{f}}k}} + \Delta {{\dot \omega }_{{\rm{f}}k}}t,\\ \Delta {M_{{\rm{f}}k}}\left( t \right) = \Delta {M_{{\rm{f}}k}} + \Delta {{\dot M}_{{\rm{f}}k}}t. \end{array} \right. $

(3)

其中：Δa_fk、Δi_fk、ΔΩ_fk、ΔM_fk分别表示第k(k=1，2，…，N)个伴飞航天器与参考航天器之间的半长轴、轨道倾角、升交点赤经及平近点角之差，$\Delta {\mathit{\dot \Omega }_{{\rm{f}}k}}$、$\Delta {\dot M_{{\rm{f}}k}}$、$\Delta {\dot \omega _{{\rm{f}}k}}$分别表示伴飞航天器与参考航天器之间由摄动因素引起的升交点赤经、平近点角、近点幅角的漂移差。

将式(3)代入式(1)中，根据航天器轨道动力学模型对其进行化简。对于圆形轨道或近圆轨道，其偏心率e_r、e_fk均为一阶小量，略去二阶以上小量，可得从航天器相对于主航天器的编队相对位置随时间变化方程：

$ \left\{ \begin{array}{l} {x_{{\rm{f}}k}}\left( t \right) = \Delta {a_{{\rm{f}}k}} - {A_{{\rm{f}}k}}\cos \left( {{M_{\rm{r}}} + {{\dot M}_{\rm{r}}}t + {\varphi _{{\rm{f}}k}}} \right),\\ {y_{{\rm{f}}k}}\left( t \right) = 2{B_{{\rm{f}}k}}\sin \left( {{M_{\rm{r}}} + {{\dot M}_{\rm{r}}}t + {\varphi _{{\rm{f}}k}}} \right) + \\ \;\;\;\;\;\;\;\;\;\;\;{D_{\rm{r}}}\left( {\Delta {\omega _{{\rm{f}}k}} + \Delta {M_{{\rm{f}}k}} + \Delta {\mathit{\Omega }_{{\rm{f}}k}}\cos {i_{\rm{r}}}} \right) + \\ \;\;\;\;\;\;\;\;\;\;\;{D_{\rm{r}}}\left( {\Delta {{\dot \omega }_{{\rm{f}}k}} + \Delta {{\dot M}_{{\rm{f}}k}} + \Delta {{\mathit{\dot \Omega }}_{{\rm{f}}k}}\cos {i_{\rm{r}}}} \right)t,\\ {z_{{\rm{f}}k}}\left( t \right) = - {C_{{\rm{f}}k}}\cos \left[ {{M_{\rm{r}}} + {{\dot M}_{\rm{r}}}t + {\omega _{\rm{r}}} + \Delta {{\dot \omega }_{\rm{r}}}t + } \right.\\ \;\;\;\;\;\;\;\;\;\;\;\left. {2{e_{\rm{r}}}\sin \left( {{M_{\rm{r}}} + {{\dot M}_{\rm{r}}}t} \right) + {\varphi _{{\rm{f}}k}}} \right]. \end{array} \right. $

(4)

其中：

$ {A_{{\rm{f}}k}} = \sqrt {a_{{\rm{f}}k}^2e_{{\rm{f}}k}^2 + a_{\rm{r}}^2e_{\rm{r}}^2 - 2{a_{{\rm{f}}k}}{e_{{\rm{f}}k}}{a_{\rm{r}}}{e_{\rm{r}}}\cos \left( {\Delta {M_{{\rm{f}}k}} + \Delta {{\dot M}_{{\rm{f}}k}}t} \right)} , $

$ {B_{{\rm{f}}k}} = {D_{\rm{r}}}\sqrt {e_{{\rm{f}}k}^2 + e_{\rm{r}}^2 - 2{e_{{\rm{f}}k}}{e_{\rm{r}}}\cos \left( {\Delta {M_{{\rm{f}}k}} + \Delta {{\dot M}_{{\rm{f}}k}}t} \right)} , $

$ {C_{{\rm{f}}k}} = {D_{\rm{r}}}\sqrt {{{\left( {\Delta {\mathit{\Omega }_{{\rm{f}}k}} + \Delta {{\mathit{\dot \Omega }}_{{\rm{f}}k}}t} \right)}^2}{{\sin }^2}{i_{\rm{r}}} + \Delta i_{{\rm{f}}k}^2} , $

$ {D_{\rm{r}}} = {a_{\rm{r}}}\left[ {1 - {e_{\rm{r}}}\cos \left( {{M_{\rm{r}}} + {e_{\rm{r}}}\sin {M_{\rm{r}}}} \right)} \right], $

$ \sin {\varphi _{{\rm{f}}k}} = {D_{\rm{r}}}{e_{{\rm{f}}k}}\sin \left( {\Delta {M_{{\rm{f}}k}} + \Delta {{\dot M}_{{\rm{f}}k}}t} \right)/{B_{{\rm{f}}k}}, $

$ \cos {\varphi _{{\rm{f}}k}} = {D_{\rm{r}}}\left[ {{e_{{\rm{f}}k}}\cos \left( {\Delta {M_{{\rm{f}}k}} + \Delta {{\dot M}_{{\rm{f}}k}}t} \right) - {e_{\rm{r}}}} \right]/{B_{{\rm{f}}k}}, $

${\rm sin}{\phi _{{\rm{f}}k}} = {D_{\rm{r}}}\Delta {i_{{\rm{f}}k}}/{C_{{\rm{f}}k}}, $

$ \cos {\phi _{{\rm{f}}k}} = {D_{\rm{r}}}\left( {\Delta {\mathit{\Omega }_{{\rm{f}}k}} + \Delta {{\mathit{\dot \Omega }}_{{\rm{f}}k}}t} \right)\sin {i_{\rm{r}}}/{B_{{\rm{f}}k}}. $

式(4)即为近距离编队在J₂摄动影响下的时间演化模型。由式(4)可知，在队形演化过程中伴飞星与参考星的相对位置既有相位随时间变化的周期运动项，也有随时间推移而逐渐漂移的扩散项。

2 面向典型任务模式的自然编队构型设计

以6颗微纳卫星组成的分布式遥感系统为例。系统包含信息处理卫星、普查卫星、详查卫星、多光谱卫星、微光成像卫星和视频卫星^[2]。其功能划分为：信息处理卫星负责系统任务规划、星地信息传输、载荷信息处理等；普查卫星用于目标大范围搜索和位置确定，并将结果发送给信息处理卫星进行任务规划；其他各星根据不同任务进行构型转换，协同工作获取不同维度的目标影像，结果发送给信息处理卫星，由信息处理卫星发送到地面用户。

根据设计的任务需求，编队的典型工作模式主要有2种：1)重点目标协同观察模式；2)详查引导普查工作模式。

编队星间通信拓扑为单层树图通信拓扑，如图 2所示，即各成像卫星通过与信息处理卫星的星间通信获取地面控制信息，并通过星间通信将获取的图像信息传送到信息处理卫星，再由信息处理卫星传送到地面，而各成像卫星之间、成像卫星与地面之间不直接进行信息交换。当系统收到地面发送的进入某一工作模式或进行工作模式切换的指令后，编队进入队形重构状态，各航天器根据信息处理卫星的计算结果调整自身位置进行队形重构，在达到预定位置后发动机关闭编队进入队形维持状态。在队形维持期间，编队采用自然编队构型模式。在自然编队条件下编队应维持一个工作周期内队形扩散应不超出20 km的星间数传覆盖范围。

2.1 重点目标协同观测工作模式

对于重点目标观测任务来说，需要微纳航天器集群内的各卫星使用不同的成像手段，从不同角度对重点目标进行协同观测。详查卫星应用较高的分辨率对地面重点目标进行高精度观测；多光谱卫星获取目标地域对应的光谱信息；视频卫星以“太空录像”的方式连续多帧序列成像，对目标区域进行凝视，监测目标区域的动态变化及目标的运动矢量；普查卫星和微光成像卫星辅助上述3颗卫星从不同角度进行协同观测。

在重点目标协同观测模式下，若用各成像卫星之间的连线所包围的几何面积来衡量编队的观测基线，该模式下采用平面内绕飞构型能够在单层树图通信拓扑约束下拓展不同成像的观测基线，如图 3所示。

在该构型模式下，各成像卫星以信息处理卫星为中心构成绕飞队形，绕飞轨迹在ySz平面内的投影为半径为l的圆，即：

$ y_{{\rm{f}}k}^2\left( t \right) + z_{{\rm{f}}k}^2\left( t \right) = {l^2}. $

(5)

对于近圆轨道编队由式(5)可知，编队内绕飞中心为原点，且不随时间变化。代入式(4)可得：

$ \left\{ \begin{array}{l} {D_{\rm{r}}}\left( {\Delta {\omega _{{\rm{f}}k}} + \Delta {M_{{\rm{f}}k}} + \Delta {\mathit{\Omega }_{{\rm{f}}k}}\cos {i_{\rm{r}}}} \right) = 0,\\ {D_{\rm{r}}}\left( {\Delta {{\dot \omega }_{{\rm{f}}k}} + \Delta {{\dot M}_{{\rm{f}}k}} + \Delta {{\mathit{\dot \Omega }}_{{\rm{f}}k}}\cos {i_{\rm{r}}}} \right)t = 0,\\ y_{{\rm{f}}k}^2\left( t \right) + z_{{\rm{f}}k}^2\left( t \right) = 4B_{{\rm{f}}k}^2{\sin ^2}\left( {{M_{\rm{r}}} + {{\dot M}_{\rm{r}}}t + {\varphi _{{\rm{f}}k}}} \right) + \\ C_{{\rm{f}}k}^2{\cos ^2}\left( {{M_{\rm{r}}} + {{\dot M}_{\rm{r}}}t + {\omega _{\rm{r}}} + \Delta {{\dot \omega }_{\rm{r}}}t + {\phi _{{\rm{f}}k}}} \right) = {l^2}. \end{array} \right. $

(6)

对式(6)进行化简可得：

$ \left\{ \begin{array}{l} \Delta {\omega _{{\rm{f}}k}} + \Delta {M_{{\rm{f}}k}} + \Delta {\mathit{\Omega }_{{\rm{f}}k}}\cos {i_{\rm{r}}} = 0,\\ \Delta {{\dot \omega }_{{\rm{f}}k}} + \Delta {{\dot M}_{{\rm{f}}k}} + \Delta {{\mathit{\dot \Omega }}_{{\rm{f}}k}}\cos {i_{\rm{r}}} = 0,\\ {C_{{\rm{f}}k}} = 2{B_{{\rm{f}}k}} = l,\\ {\varphi _{{\rm{f}}k}} = {\omega _{\rm{r}}} + \Delta {{\dot \omega }_{\rm{r}}}t + {\phi _{{\rm{f}}k}}. \end{array} \right. $

(7)

由于编队中伴飞航天器的数目较多，若令N个伴飞航天器在ySz平面内的相位均匀分布，即：

$ {\varphi _{{\rm{f}}k}} = \frac{{2{\rm{ \mathsf{ π} }}k}}{N},\;\;\;\;\;k = 1,2, \cdots ,N. $

(8)

在树图通信拓扑的平面内绕飞模式下，摄动因素影响的长期项表现为队形的空间扩散，短期项影响为星间相位的短期波动。在实际任务中可以通过各卫星的姿态控制来实现编队的协同观测，因而可以在松散控制的前提下，忽略内摄动引起短期相位变化，即忽略$\Delta {\dot \omega _{{\rm{f}}k}}$、$\Delta {\mathit{\dot \Omega }_{{\rm{f}}k}}$、$\Delta {\dot M_{{\rm{f}}k}}$等项的影响。联立式(7)、(8)和(4)可得约束方程组：

$ \left\{ \begin{array}{l} \Delta {M_{{\rm{f}}k}} = \arctan \left[ {\frac{{l\sin \left( {2{\rm{ \mathsf{ π} }}k/N + {\omega _{\rm{r}}}} \right)}}{{l\cos \left( {2{\rm{ \mathsf{ π} }}k/N + {\omega _{\rm{r}}}} \right) + 2{D_{\rm{r}}}{e_{\rm{r}}}}}} \right],\\ \Delta {\mathit{\Omega }_{{\rm{f}}k}} = \left. {\frac{{l\cos \left( {2{\rm{ \mathsf{ π} }}k/N} \right)}}{{{D_{\rm{r}}}\sin {i_{\rm{r}}}}}} \right),\\ \Delta {i_{{\rm{f}}k}} = \frac{{l\sin \left( {2{\rm{ \mathsf{ π} }}k/N} \right)}}{{{D_{\rm{r}}}}},\\ \Delta {e_{{\rm{f}}k}} = \sqrt {e_{\rm{r}}^2 + l{e_{\rm{r}}}\cos \left( {{\omega _{\rm{r}}} + \frac{{2{\rm{ \mathsf{ π} }}k}}{N}} \right)/{D_{\rm{r}}} + \frac{{{l^2}}}{{4D_{\rm{r}}^2}}} - {e_{\rm{r}}},\\ \Delta {{a}_{\text{f}k}}=-\frac{14{{J}_{2}}R_{\text{e}}^{2}{{a}_{\text{r}}}\sin 2{{i}_{\text{r}}}}{47a_{\text{r}}^{2}+7{{J}_{2}}R_{\text{e}}^{2}\left( 3\cos 2{{i}_{\text{r}}}+1 \right)}\Delta {{i}_{\text{f}k}},\\ \Delta {\omega _{{\rm{f}}k}} = - \Delta {M_{{\rm{f}}k}} - \Delta {\mathit{\Omega }_{{\rm{f}}k}}\cos {i_{\rm{r}}}. \end{array} \right. $

(9)

其中：轨道差ΔM_fk、ΔΩ_fk、Δi_fk、Δa_fk、Δω_fk仅与参考航天器的轨道参数有关，即根据参考航天器的轨道参数即可由方程(9)直接求出k(k=1，2，…，N)个伴飞航天器的轨道参数。该方法与吴宝林等^[10]提出的J₂摄动下的轨道相对漂移率最小优化方法相比，优点在于简单、直接、能够轻易地将编队航天器的数量拓展到N。

参考星的预定入轨参数为：半长轴a_r=6 911.98 km，偏心率e_r=0.005 593 7，轨道倾角i_r=97.459 4°，升交点赤经Ω_r=64.026 1°，近地点幅角ω_r=187.981°，平近点角M_r=37.918 6°。代入到式(9)，令k=1，2，…，5，可得到编队中每颗伴星的轨道参数如表 1所示。

以某微纳卫星为例，其包络尺寸为324 mm×360 mm×400 mm，迎风面面积为890.9 cm²，对天面积为1 049.76 cm²，平台质量15.5 kg，载荷质量5 kg，卫星总质量20.5 kg，质心坐标系的转动惯量为J_x=0.28 kg·m²，J_y=0.33 kg·m²，J_z=0.35 kg·m²。

将卫星面质比，质心坐标系转动惯量J_x、J_y、J_z及表 1中各航天器的轨道根数等参数代入到STK中进行仿真。图 4给出144 h内伴飞航天器1与参考航天器之间的三维位置关系图。由图 4可知，在一个工作周期(96 h)内编队保持了较好的近圆型绕飞构型，且整个绕飞轨迹均在参考航天器20 km距离范围内，即满足编队在一个工作周期内的扩散不超过20 km的要求。

随着飞行时间的逐渐增加，绕飞轨迹向椭圆逐渐过渡。原因为在式(9)化简过程中忽略了$\Delta {\mathit{\dot \Omega }_{{\rm{f}}k}}$、$\Delta {\dot M_{{\rm{f}}k}}$在短期内的影响，随着时间增加，该摄动影响逐渐增加并最终破坏了分布式遥感系统的圆形绕飞构型。伴星2、3、4、5与参考星之间的位置关系与图 4相近，本文不再给出。

图 5给出了200 min内各伴飞航天器在参考航天器坐标系y轴投影的相位关系图。

由图 5可知伴飞航天器1、2、3、4、5在ySz平面内相位差恒定，且按照圆构型沿逆时针方向缓慢旋转。由于篇幅限制无法画出分布式系统在自然编队情况下一个工作周期内的相位关系，但从STK仿真结果来看，系统在一个工作周期内各伴飞航天器之间的相位差基本保持良好。但随着时间增加各伴飞航天器之间的相位差也逐渐发生变化，原因是在式(9)的化简过程中忽略了$\Delta {\dot \omega _{{\rm{f}}k}}$在短期的影响。

由于本文主要考虑的是遥感任务，为了使获取图像的尺度变化较小，遥感卫星一般情况下均选择圆形和近圆形轨道，因此近圆轨道假设是合理的。在参考星轨道符合近圆轨道假设的情况下，变更参考星轨道能够获得与前文的类似仿真结果。

2.2 普查引导详查工作模式

在普查引导详查工作模式下，编队采用领航编队模式。以6颗微纳航天器组成的分布式遥感系统为例，普查卫星在前利用较大的视场进行大范围搜查，微光成像卫星协同普查卫星进行普查工作，然后将普查信息传送到信息处理卫星；信息处理卫星对信息进行处理，对重点目标进行辨识和信息提取，并反馈到详查卫星；详查卫星利用较高的分辨率对重点目标进行详查，详查结果反馈到信息处理卫星；详查过程中，多光谱卫星协同详查卫星完成目标详查，视频卫星对详查范围进行视频成像监测，检查详查目标有无异常及侦查漏洞。其在领航模式下的编队构型图，如图 6所示。

在该模式下各航天器应满足以下条件：

条件1：z方向相对位置恒为零，

$ {z_{{\rm{f}}k}}\left( t \right) = 0. $

(10)

条件2：y方向相对位置不随时间漂移，

$ {D_{\rm{r}}}\left( {\Delta {\omega _{ \cdot {\rm{f}}k}} + \Delta {{\dot M}_{{\rm{f}}k}} + \Delta {{\mathit{\dot \Omega }}_{{\rm{f}}k}}\cos {i_{\rm{r}}}} \right)t = 0. $

(11)

条件3：在y方向相邻2个航天器之间的间距相等，且伴飞航天器1(普查卫星)在参考航天器之前，伴飞航天器2、3、4、5在参考航天器之后，设2个相邻航天器之间的间距为h，则有

$ \left\{ \begin{array}{l} {D_{\rm{r}}}\left( {\Delta {\omega _{{\rm{f}}k}} + \Delta {M_{{\rm{f}}k}} + \Delta {\mathit{\Omega }_{{\rm{f}}k}}\cos {i_{\rm{r}}}} \right) = h,\;\;\;\;k = 1;\\ {D_{\rm{r}}}\left( {\Delta {\omega _{{\rm{f}}k}} + \Delta {M_{{\rm{f}}k}} + \Delta {\mathit{\Omega }_{{\rm{f}}k}}\cos {i_{\rm{r}}}} \right) = \\ \;\;\;\;\left( { - k + 1} \right)h,\;\;\;k = 2,3, \cdots ,N. \end{array} \right. $

(12)

根据领航模式下的约束条件，由式(10)—(12)等约束条件及式(2)可得，当k=1时参考航天器与伴飞航天器之间的轨道要素参方程为

$ \left\{ \begin{array}{l} \Delta {a_{{\rm{f}}k}} = 0,\\ \Delta {i_{{\rm{f}}k}} = 0,\\ \Delta {\mathit{\Omega }_{{\rm{f}}k}} = 0,\\ \Delta {e_{{\rm{f}}k}} = 0,\\ \Delta {M_{{\rm{f}}k}} = h/{D_{\rm{r}}} - \Delta {\omega _{{\rm{f}}k}},\;\;\;\;\;当\;k = 1. \end{array} \right. $

(13)

当k=2，3，…，N时参考航天器与伴飞航天器之间的轨道要素参方程为

$ \left\{ \begin{array}{l} \Delta {a_{{\rm{f}}k}} = 0,\\ \Delta {i_{{\rm{f}}k}} = 0,\\ \Delta {\mathit{\Omega }_{{\rm{f}}k}} = 0,\\ \Delta {e_{{\rm{f}}k}} = 0,\\ \Delta {M_{{\rm{f}}k}} = \left( { - k + 1} \right)h/{D_{\rm{r}}} - \Delta {\omega _{{\rm{f}}k}},\;\;\;\;\;当\;k = 1. \end{array} \right. $

(14)

若令式(13)、(14)中Δω_fk=0，将预定入轨参数代入，并令k=1, 2，…，5，可得各颗伴飞星轨道参数如表 3所示。

图 7给出了STK仿真领航模式x方向前2 000 min内编队构型随时间的演化关系图。由仿真结果可知伴飞星与参考星之间的位置在x方向呈现小幅的周期性的波动，y方向基本保持恒定，z方向位置基本保持不变。由仿真结果可知，编队构型在甚至数月的时间内维持稳定，原因在于与树图通信拓扑约束下的圆形绕飞模式相比，在方程式(13)、(14)的求解过程中没有忽略$\Delta {\mathit{\dot \Omega }_{{\rm{f}}k}}$、$\Delta {\dot M_{{\rm{f}}k}}$、$\Delta {\dot \omega _{{\rm{f}}k}}$等摄动项引起轨道漂移差的影响。

由方程(13)、(14)可知，当系统处于领航模式时Δa_fk=Δi_fk=Δe_fk=0，由方程(2)可知分布式系统中各航天器在J₂摄动影响的轨道漂移值相等，即漂移差$\Delta {\mathit{\dot \Omega }_{{\rm{f}}k}}$=0、$\Delta {\dot M_{{\rm{f}}k}}$=0、$\Delta {\dot \omega _{{\rm{f}}k}}$=0，满足零J₂摄动影响的条件^[12]$\Delta {\mathit{\dot \Omega }_{{\rm{f}}k}}$=0和$\Delta {\dot M_{{\rm{f}}k}}$+$\Delta {\dot \omega _{{\rm{f}}k}}$=0，因而在领航模式下分布式航天系统的构型能够维持长期稳定，即间接利用了零J₂摄动影响相关的研究成果。

该方法与文[12]所使用的“零J₂摄动条件”方法相比具有更大的适用范围，如在重点目标协同跟踪工作模式时，系统需要采用的平面内绕飞构型就无法满足零J₂摄动条件；而在普查引导详查工作模式时采用该方法的求解结果满足零J₂摄动条件，即间接利用了零J₂摄动影响相关的研究成果。与文[10]中的“相对漂移率最小优化求解法”相比，该方法求解6颗微纳航天器编队所需的计算量为400多次乘加运算，计算量不到前者的1/10。

3 结论

以微纳航天器技术为基础, 分布式航天器系统已经成为当前航天领域发展的一个重要方向。受所携带燃料较少等条件限制，系统通常采用松散编队控制。在自然编队条件下，系统的构型通常会受到J₂等摄动因素的影响。本文针对圆轨道或近圆轨道提出了一种考虑摄动时间因素影响的集群航天器自然编队演化模型, 并在此基础上根据分布遥感编队的任务特点提出了一种参数直接求解的方法，可以在平面内绕飞模式、领航模式等典型模式下直接求解k=1，2，…，N个伴飞航天器的轨道参数。由STK仿真结果可知，对于树图通信拓扑下的平面内绕飞模式，编队能够在一个工作周期内维持构型基本不变，对于领航模式编队构型能够基本维持长期稳定，满足基于微纳航天器的分布式可重构遥感技术的任务需求。