CN111680354B - 近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法 - Google Patents

Abstract

本发明公开了一种近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法，基于地球扁率(J₂项)的轨道摄动模型，利用地球扁率对升交点赤经，近地点幅角和平近点角的长期摄动的平均影响，使得卫星轨道平面相对地球的平均漂移周期和轨道交点周期共振，形成轨道相对地球表面星下点和摄影点轨迹的重访回归。通过给定的轨道倾角、偏心率、回归天数和周期数后，可以求解轨道半长轴，继而解析求解星下点轨迹自交点的经度，再数值筛选出对应精度下自交点的纬度，即可得到一个回归周期内星下点轨迹自交点的分布情况。本发明给出了一种实用的近地回归轨道星下点轨迹自交点计算方法，在近地测绘、遥感、侦察、InSAR等回归轨道卫星测量任务中有重要应用价值。

Description

近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法

技术领域

本发明属于遥感卫星技术领域，尤其涉及一种近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法。

背景技术

根据任务需求，遥感卫星常常需要在一定时间内对特定区域进行重复观测来完成对地成像和3D地形建模等任务，对重点关注地区需要定时重访回归探测。实现这种功能一般采取近地回归轨道，具有成本低，成像效果好，重访时间短，星下点密集，可探测范围广等优点。

近地轨道由于受地球非球形引力摄动(尤其是J₂项)的影响比较明显，长期来看会引起卫星轨道的较大漂移，对于遥感卫星的重访探测有巨大影响，因此，对于设计星下点重访的卫星轨道而言，必需考虑地球非球形引力的摄动影响。

传统的近地卫星对地观测常采用太阳同步轨道，虽然实现了对地面指定点的定时重访观测，然而也存在观测区域和卫星观测角相对固定，相邻观测时间间隔较长(一日一次)、计算量大，易丢失等缺点。

发明内容

本发明的技术解决问题：克服现有技术的不足，提供一种近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法，克服了传统方法中依靠数值搜索星下点轨迹自交点位置计算量大，易丢失的缺点，提高了计算效率，达到星载计算机的处理能力要求，实现了星下点轨迹自交点位置的在线计算。

为了解决上述技术问题，本发明公开了一种近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法，包括：

根据轨道倾角、轨道偏心率、回归周期天数和回归圈数，迭代求解回归轨道半长轴；

根据迭代求解得到的回归轨道半长轴和地球扁率摄动系数J₂下的轨道运行规律，对近地回归轨道卫星星下点经纬度和摄影点经纬度分别进行计算；

根据计算得到的星下点经纬度和摄影点经纬度，分别得到一个回归周期内星下点轨迹自交点的分布和摄影点轨迹自交点的分布。

在上述近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法中，根据轨道倾角、轨道偏心率、回归周期天数和回归圈数，迭代求解回归轨道半长轴，包括：

选定轨道偏心率e、轨道倾角i、回归周期天数D和回归圈数N；

确定回归轨道半长轴a的迭代解算公式：

其中，μ表示地球引力常数，W_E表示地球东向角速度，J₂表示地球扁率摄动系数，R_E表示地球平均半径；

根据迭代初始值

利用牛顿迭代法，解算得到回归轨道半长轴a的值。

在上述近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法中，选定轨道偏心率e、轨道倾角i、回归周期天数D和回归圈数N，包括：

选定轨道偏心率e和轨道倾角i；

根据选定的轨道偏心率e和轨道倾角i，计算得到地球扁率摄动系数J₂下的轨道参数：

其中，Ω表示升交点赤经，ω表示近地点幅角，M表示平近点角，W_Ω表示Ω在地球扁率摄动系数J₂下的平均变化率，W_ω表示ω在地球扁率摄动系数J₂下的平均变化率，W_M表示M在地球扁率摄动系数J₂下的平均变化率；

根据W_E和W_Ω，确定地球与轨道平面的相对关系的周期P_E-O：

根据W_ω和W_M，确定航天器绕地球的周期P_nod：

对于互为质数的正整数D和N，满足：

DP_E-O＝NP_nod。

在上述近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法中，根据求解得到的回归轨道半长轴和地球扁率摄动系数J₂下的轨道运行规律，对近地回归轨道卫星星下点经纬度进行计算，包括：

计算格林威治赤经α_G：

其中，MJDC表示简化儒略日，

MJD表示日期码；

计算卫星赤纬δ和赤经α：

α＝arcsin(sin(i)×sin(u))

其中，u表示纬度幅角；

计算星下点大地经度λ和星下点大地纬度φ：

λ＝α

其中，

表示地球扁平率；

取星下点轨迹上交叉两点的经纬度为(λ₁,φ₁)和(λ₂,φ₂)，则有：λ₁＝λ₂，φ₂＝φ₁；

由三角函数性质，知：u₁＝u₂+2kπ，u₁+u₂＝(2k+1)π，则：

其中，u₁表示(λ₁,φ₁)的纬度幅角，u₂表示(λ₂,φ₂)的纬度幅角，u₁的模与u₂的模相等，均为u₀；

由球面三角性质知：

其中，α_G1表示纬度幅角u₁对应的格林威治赤经，α_G2表示纬度幅角u₂对应的格林威治赤经；

由三角函数性质知：

(λ₁+α_G1-Ω)-(λ₂+α_G2-Ω)＝2nπ

(λ₁+α_G1-Ω)-(λ₂+α_G2-Ω)＝(2n+1)π

相减得：

(λ₁+α_G1-Ω)+(λ₂+α_G2-Ω)＝nπ

假定初始时刻t₀＝0，则：

λ₁+λ₂＝nπ+2Ω-2α_G0-W_E(t₁+t₂)

其中，α_G0表示初始时刻t₀对应的格林威治赤经，t₁和t₂分别表示第一次和第二次经过自交点时刻；

又：

u₁+u₂＝2u₀+W_u(t₁+t₂)＝(2k+1)π

W_u＝W_M+W_ω

则有：

求解可得星下点轨迹交叉点经度λ_交＝λ₁＝λ₂；

以纬度幅角为自变量，仿真求解和数值搜索，得到一个回归周期内星下点轨迹在经度λ_交下的纬度φ_交。

在上述近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法中，根据求解得到的回归轨道半长轴和地球扁率摄动系数J₂下的轨道运行规律，对近地回归轨道卫星摄影点经纬度进行计算，包括：

确定卫星飞行方向上某位置摄影点(λ′,φ′)与星下点对地心的张角为γ：

其中,β表示摄影点与卫星连线与卫星星下点之间的夹角，H表示星下点与卫星之间的距离，λ′表示摄影点经度，φ′表示摄影点纬度；

则，摄影点的经纬度按照如下坐标转换实现：

L＝I·R_z(Ω)·R_x(i)·R_z(u)·R_y(-γ)

其中，L表示以摄影点为原点的“天-东-北”坐标系，I表示惯性坐标系，x轴指向春分点，R_z(Ω)、R_x(i)、R_z(u)、R_y(-γ)表示坐标转换矩阵；

L＝[1 0 0]^T

I＝[cosλ′cos(φ′+α_G)cosλ′sin(φ′+α_G)sinλ′]^T

展开可得：

求解得到摄影点(λ′,φ′)的经纬度值。

本发明具有以下优点：

(1)本发明公开了一种近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法，可实现卫星探测器在星下点和摄影点轨迹的自交点处的短期高频次覆盖，获取高质量的成像数据，利于成像和建模，且提高了卫星单次发射的商业价值。

(2)本发明公开了一种近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法，采用半解析的方法计算回归轨道星下点和摄影点轨迹的自交点位置，满足了星载计算机的计算能力，实现了星上自主判断成像时刻，避免了星上误判拍照时机，减少了地面工作人员工作量。

(3)本发明公开了一种近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法，通过半解析求解回归轨道星下点和摄影点轨迹的自交点位置，可系统的得到星下点轨迹的自交叉点的分布与演化规律，可为商业卫星的发射轨道以及卫星长时间在轨服务提供有利参考。

(4)本发明公开了一种近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法，可通过对全球重点目标区域设计有针对的近地低轨回归轨道，实现广地域，可变观测角和小时间间隔的多次重访，基于回归轨道的近地观测卫星可以改善传统太阳同步轨道观测卫星的不足，更适用于测绘、遥感、侦察、InSAR等对地成像和3D地形建模等任务。

附图说明

图1是本发明实施例中一种近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法的步骤流程图；

图2是本发明实施例中一种坐标转换示意图；

图3是本发明实施例中一种星下点的轨迹示意图；

图4是本发明实施例中一种摄影点的轨迹示意图；

图5是本发明实施例中一种地心惯性坐标系下的轨道示意图；

图6是本发明实施例中一种星下点轨迹自交点经度示意图；

图7是本发明实施例中一种星下点轨迹自交点示意图；

图8是本发明实施例中一种摄影点轨迹部分自交点示意图。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明公开的实施方式作进一步详细描述。

实施例1

如图1，在本实施例中，该近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法，包括：

步骤101，根据轨道倾角、轨道偏心率、回归周期天数和回归圈数，迭代求解回归轨道半长轴。

在本实施例中，回归轨道半长轴的具体求解过程可以如下：

a)选定轨道偏心率e、轨道倾角i、回归周期天数D和回归圈数N。

首先，选定轨道偏心率e和轨道倾角i。

然后，根据选定的轨道偏心率e和轨道倾角i，计算得到地球扁率摄动系数J₂下的轨道参数：

其中，Ω表示升交点赤经，ω表示近地点幅角，M表示平近点角，W_Ω表示Ω在地球扁率摄动系数J₂下的平均变化率，W_ω表示ω在地球扁率摄动系数J₂下的平均变化率，W_M表示M在地球扁率摄动系数J₂下的平均变化率。

由于地球具有东向角速度W_E；同时，轨道平面具有角速度W_Ω(向东为正)，是由J₂造成的升交点赤经的平均变化率；因此，地球相对于轨道平面的角速度是W_E-W_Ω，可得地球与轨道平面的相对关系的周期P_E-O：

另一方面，航天器绕地球的周期P_nod，即交点周期是：

最后，对于互为质数的正整数D和N，满足：

DP_E-O＝NP_nod

即如果在整数个轨道日内航天器运行整数个圈，则航天器与地球的关系能重复出现，因而卫星星下点轨迹满足回归条件。

b)确定回归轨道半长轴a的迭代解算公式：

其中，μ表示地球引力常数，W_E表示地球东向角速度，J₂表示地球扁率摄动系数，R_E表示地球平均半径。

c)根据迭代初始值

利用牛顿迭代法，解算得到回归轨道半长轴a的值。

步骤102，根据迭代求解得到的回归轨道半长轴和地球扁率摄动系数J₂下的轨道运行规律，对近地回归轨道卫星星下点经纬度和摄影点经纬度分别进行计算。

在本实施例中，近地回归轨道卫星星下点经纬度的计算过程可以如下：

计算格林威治赤经α_G：

其中，MJDC表示简化儒略日，

MJD表示日期码。

计算卫星赤纬δ和赤经α：

α＝arcsin(sin(i)×sin(u))

其中，u表示纬度幅角。

计算星下点大地经度λ和星下点大地纬度φ：

λ＝α

其中，

表示地球扁平率。

取星下点轨迹上交叉两点的经纬度为(λ₁,φ₁)和(λ₂,φ₂)，则有：λ₁＝λ₂，φ₂＝φ₁。

由三角函数性质，知：u₁＝u₂+2kπ，u₁+u₂＝(2k+1)π，则：

其中，u₁表示(λ₁,φ₁)的纬度幅角，u₂表示(λ₂,φ₂)的纬度幅角，u₁的模与u₂的模相等，均为u₀。

由球面三角性质知：

其中，α_G1表示纬度幅角u₁对应的格林威治赤经，α_G2表示纬度幅角u₂对应的格林威治赤经。

由三角函数性质知：

(λ₁+α_G1-Ω)-(λ₂+α_G2-Ω)＝2nπ

(λ₁+α_G1-Ω)-(λ₂+α_G2-Ω)＝(2n+1)π

相减得：

(λ₁+α_G1-Ω)+(λ₂+α_G2-Ω)＝nπ

假定初始时刻t₀＝0，则：

λ₁+λ₂＝nπ+2Ω-2α_G0-W_E(t₁+t₂)

其中，α_G0表示初始时刻t₀对应的格林威治赤经，t₁和t₂分别表示第一次和第二次经过自交点时刻；

又：

u₁+u₂＝2u₀+W_u(t₁+t₂)＝(2k+1)π

W_u＝W_M+W_ω

则有：

上述方程为一个丢番图方程，可通过解析计算求解得星下点轨迹交叉点经度λ_交＝λ₁＝λ₂；然后，以纬度幅角为自变量，仿真求解和数值搜索，得到一个回归周期内星下点轨迹在经度λ_交下的纬度φ_交。

在本实施例中，近地回归轨道卫星摄影点经纬度的计算过程可以如下：

确定卫星飞行方向上某位置摄影点(λ′,φ′)与星下点对地心的张角为γ：

其中，β表示摄影点与卫星连线与卫星星下点之间的夹角，H表示星下点与卫星之间的距离，λ′表示摄影点经度，φ′表示摄影点纬度。

摄影点的经纬度可按照图2所示的坐标转换实现，即，摄影点的经纬度按照如下坐标转换实现：

L＝I·R_z(Ω)·R_x(i)·R_z(u)·R_y(-γ)

其中，L表示以摄影点为原点的“天-东-北”坐标系，I表示惯性坐标系，x轴指向春分点，R_z(Ω)、R_x(i)、R_z(u)、R_y(-γ)表示坐标转换矩阵。

L＝[1 0 0]^T

I＝[cosλ′cos(φ′+α_G)cosλ′sin(φ′+α_G)sinλ′]^T

展开可得：

求解得到摄影点(λ′,φ′)的经纬度值。

由于摄影点轨迹的自交叉点分布与星下点轨迹的自交叉点分布有一一对应关系，即一个摄影点轨迹的自交点对应一个星下点轨迹的自交叉点，所以摄影点轨迹的自交点可以参照星下点轨迹的自交点方法计算。

步骤103，根据计算得到的星下点经纬度和摄影点经纬度，分别得到一个回归周期内星下点轨迹自交点的分布和摄影点轨迹自交点的分布。

在本实施例中，基于任务需求和设计方案设计卫星轨道高度，轨道倾角和轨道偏心率，使星下点轨迹的自交点落在目标点的中心区域，同时确保摄影点轨迹的边界能够覆盖整个目标。考察轨道高度、轨道倾角和回归周期和回归天数之间的关系，轨道高度高或者回归周期长，会使得重访时间间隔变大，有利于提高成像质量，但是不能实现对目标区域的短期高频次覆盖；轨道倾角大，则星下点与摄影点组成的探测区域分布变大，可探测的纬度层更高，但是覆盖区域会变得稀疏。根据任务需求和星上载荷限制，选取最佳的轨道设计方案和载荷。基于所优选的星上载荷，设计卫星平台的构型、能源系统、姿轨控系统，遥控遥测、数传等分系统的工程设计细节。最后进行快速获取高质量的成像数据，数据分析以及图像处理等后续工作。

实施例2

在上述实施例的基础上，下面结合一个实例进行说明。

在本实施例中，星下点和摄影点轨迹如图3和图4所示，中心部分轨迹为星下点轨迹，侧面轨迹为摄影点轨迹，通过星下点轨迹和摄影点轨迹分别形成的自交点可以确定覆盖范围的中心区域和边界范围。

具体计算步骤如下：

(1)选定如下参数：e＝0，i＝94°；根据轨道高度初值估计D＝3，N＝29；将J₂摄动下升交点赤经，近地点幅角和真近点角的漂移速率表示为回归轨道半长轴的函数。

为了简化系统，引入归一化单位，其中，特征长度为地球平均半径R_E＝6378140m，单位角速度为

即归一化系统中，一个周期的时间为2π，对应实际系统中未单位化的时间为

基于此，以下主要轨道数据均为归一化之后的数据。

考虑卫星轨道高度为600km，纬度幅角u＝0.80356(rad)，升交点赤经Ω＝0rad；根据卫星的动力学方程，以及初值[0；1.0940734；0；1.6398；0；0.80356]^T，卫星在地心惯性坐标系下的运动如图5所示。

(2)迭代求解回归轨道半长轴。

根据选定的参数可求解得到：

W_E＝7.2722×10^-5，W_Ω＝1.08462×10^-7，W_ω＝7.58518×10^-7

结合回归轨道半长轴a的迭代解算公式：

求得回归轨道半长轴a＝6.86635×10⁶m，满足当前回归周期天数和回归圈数条件下摄动轨道的精确值。

(3)考虑J2项对平近点角漂移速度影响：W_M＝-7.660836×10^-7，W_u＝W_M+W_ω，解丢番图方程得星下点轨迹自交点的经度值，得到自交点经度位置为步长0.36961°的均匀分布值，其分布图如图6所示，图6中星下点轨迹自交点计算结果与STK仿真分析一致，验证了结果的正确性，与传统方法相比解析求解出了经度值，提高了精度的同时大大简化了计算量。

在求解自交点纬度值时，以时间为自变量，数值搜索自交点经度下的纬度值，可以采用二分法等方法逐次逼近，至达到精度要求为止，采用二分法求得的自交点分布图如图7所示，以时间为变量，通过经度对纬度进行搜索，大大提高了计算效率；由图7知相邻两列自交点关于两列自交点经度中心线与赤道的交点呈中心对称分布，每隔一列自交点分布在同一纬度层上，每一纬度层上自交点纬度相差不大。

(4)假设卫星摆角为45°，即卫星摄影点与卫星连线和卫星星下点之间的夹角β为45°，卫星摄影点示意图如图8所示，图8中每个自交点附近的东南西北四个方向上分别存在一个摄影点，四个摄影点沿卫星摄影点轨迹所围区域为短时间内卫星能再覆盖的区域，可实现对该区域内目标的短期再次访问，绘制图像和对比分析。

本发明虽然已以较佳实施例公开如上，但其并不是用来限定本发明，任何本领域技术人员在不脱离本发明的精神和范围内，都可以利用上述揭示的方法和技术内容对本发明技术方案做出可能的变动和修改，因此，凡是未脱离本发明技术方案的内容，依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化及修饰，均属于本发明技术方案的保护范围。

本发明说明书中未作详细描述的内容属于本领域技术人员的公知技术。

Claims (4)

1.一种近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法，其特征在于，包括：

根据轨道倾角、轨道偏心率、回归周期天数和回归圈数，迭代求解回归轨道半长轴；

根据迭代求解得到的回归轨道半长轴和地球扁率摄动系数J₂下的轨道运行规律，对近地回归轨道卫星星下点经纬度和摄影点经纬度分别进行计算；

根据计算得到的星下点经纬度和摄影点经纬度，分别得到一个回归周期内星下点轨迹自交点的分布和摄影点轨迹自交点的分布；

其中：

根据迭代求解得到的回归轨道半长轴和地球扁率摄动系数J₂下的轨道运行规律，对近地回归轨道卫星星下点经纬度和摄影点经纬度分别进行计算，包括：

计算格林威治赤经α_G：

其中，MJDC表示简化儒略日，

MJD表示日期码；

计算卫星赤纬δ和赤经α：

α＝arcsin(sin(i)×sin(u))

其中，u表示纬度幅角，i表示轨道倾角，μ表示地球引力常数，Ω表示升交点赤经；

计算星下点大地经度λ和星下点大地纬度φ：

λ＝α

其中，

表示地球扁平率；R_E表示地球平均半径；a表示回归轨道半长轴；

取星下点轨迹上交叉两点的经纬度为(λ₁,φ₁)和(λ₂,φ₂)，则有：λ₁＝λ₂，φ₂＝φ₁；

由三角函数性质，知：u₁＝u₂+2kπ，u₁+u₂＝(2k+1)π，则：

其中，u₁表示(λ₁,φ₁)的纬度幅角，u₂表示(λ₂,φ₂)的纬度幅角，u₁的模与u₂的模相等，均为u₀；

由球面三角性质知：

其中，α_G1表示纬度幅角u₁对应的格林威治赤经，α_G2表示纬度幅角u₂对应的格林威治赤经；

由三角函数性质知：

(λ₁+α_G1-Ω)-(λ₂+α_G2-Ω)＝2nπ

(λ₁+α_G1-Ω)-(λ₂+α_G2-Ω)＝(2n+1)π

相减得：

(λ₁+α_G1-Ω)+(λ₂+α_G2-Ω)＝nπ

假定初始时刻t₀＝0，则：

λ₁+λ₂＝nπ+2Ω-2α_G0-W_E(t₁+t₂)

其中，α_G0表示初始时刻t₀对应的格林威治赤经，t₁和t₂分别表示第一次和第二次经过自交点时刻，W_E表示地球东向角速度；

又：

u₁+u₂＝2u₀+W_u(t₁+t₂)＝(2k+1)π

W_u＝W_M+W_ω

则有：

求解可得星下点轨迹交叉点经度λ_交＝λ₁＝λ₂；其中，W_ω表示ω在地球扁率摄动系数J₂下的平均变化率，W_M表示M在地球扁率摄动系数J₂下的平均变化率；

以纬度幅角为自变量，仿真求解和数值搜索，得到一个回归周期内星下点轨迹在经度λ_交下的纬度φ_交。

2.根据权利要求1所述的近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法，其特征在于，根据轨道倾角、轨道偏心率、回归周期天数和回归圈数，迭代求解回归轨道半长轴，包括：

选定轨道偏心率e、轨道倾角i、回归周期天数D和回归圈数N；

确定回归轨道半长轴a的迭代解算公式：

其中，J₂表示地球扁率摄动系数；

根据迭代初始值

利用牛顿迭代法，解算得到回归轨道半长轴a的值。

3.根据权利要求2所述的近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法，其特征在于，选定轨道偏心率e、轨道倾角i、回归周期天数D和回归圈数N，包括：

选定轨道偏心率e和轨道倾角i；

根据选定的轨道偏心率e和轨道倾角i，计算得到地球扁率摄动系数J₂下的轨道参数：

其中，ω表示近地点幅角，M表示平近点角，W_Ω表示Ω在地球扁率摄动系数J₂下的平均变化率；

根据W_E和W_Ω，确定地球与轨道平面的相对关系的周期P_E-O：

根据W_ω和W_M，确定航天器绕地球的周期P_nod：

对于互为质数的正整数D和N，满足：

DP_E-O＝NP_nod。

4.根据权利要求3所述的近地回归轨道卫星星下点和摄影点轨迹自交点的计算方法，其特征在于，根据求解得到的回归轨道半长轴和地球扁率摄动系数J₂下的轨道运行规律，对近地回归轨道卫星摄影点经纬度进行计算，包括：

确定卫星飞行方向上某位置摄影点(λ′,φ′)与星下点对地心的张角为γ：

其中,β表示摄影点与卫星连线与卫星星下点之间的夹角，H表示星下点与卫星之间的距离，λ′表示摄影点经度，φ′表示摄影点纬度；

则，摄影点的经纬度按照如下坐标转换实现：

L＝I·R_z(Ω)·R_x(i)·R_z(u)·R_y(-γ)

其中，L表示以摄影点为原点的“天-东-北”坐标系，I表示惯性坐标系，x轴指向春分点，R_z(Ω)、R_x(i)、R_z(u)、R_y(-γ)表示坐标转换矩阵；

L＝[1 0 0]^T

I＝[cosλ′cos(φ′+α_G) cosλ′sin(φ′+α_G) sinλ′]^T

展开可得：

求解得到摄影点(λ′,φ′)的经纬度值。

Images