CN111198402A - 基于轨道掩膜微分算子的地球重力场模型建模方法 - Google Patents

Abstract

本发明公开了基于轨道掩膜微分算子的地球重力场模型建模方法，涉及卫星重力探测技术领域，包括步骤：对轨道状态信息、星载加速度计和卫星姿态数据进行数据预处理；利用轨道掩膜微分算子计算卫星运动速度继而获得卫星动能和旋转位；通过卫星动能、旋转位、保守力位和非保守力位建立满足能量守恒定律的观测方程，在最小二乘准则下估算位系数改正值，再加上参考重力场模型位系数即可确定最终的地球重力场模型。相较于传统利用轨道牛顿多项式微分速度确定地球重力场模型的方法，本发明提出的方法有效抑制高频误差放大实现降噪，同时灵活性和计算效率高，获得高精度卫星速度观测值，进一步改进地球重力场模型精度水平。

Description

基于轨道掩膜微分算子的地球重力场模型建模方法

技术领域

本发明涉及卫星重力探测技术领域，特别是涉及基于轨道掩膜微分算子的地球重力场模型建模方法。

背景技术

地球重力场是地球最基本的一个物理场，反映了地球邻近空间及其内部物质迁移运动状态，因此确定高精度地球重力场模型是研究地球相关学科的热点问题。利用卫星重力测量技术改善重力场中长波信号，以及基于卫星跟踪卫星技术确定地球重力场模型的能量守恒法简洁高效被广泛采用。

能量守恒法是根据能量守恒定律建立卫星运动速度与地球重力场模型位系数间的线性函数关系，其关键问题是获取高精度的卫星运动速度值。计算卫星运动速度常用方法有经典的牛顿微分法，即通过牛顿插值的一次微分公式计算卫星速度，导致高频误差放大，降低速度观测数据的信噪比，进而影响重力场模型的建模精度。

发明内容

本发明实施例提供了基于轨道掩膜微分算子的地球重力场模型建模方法，可以解决现有技术中存在的问题。

本发明提供了基于轨道掩膜微分算子的地球重力场模型建模方法，包括以下步骤：

对重力专用卫星的轨道状态信息、星载加速度计数据和原始姿态数据进行数据预处理，获得连续精密轨道和非保守力加速度；

构建以计算时刻为中心的窗口函数，利用掩膜微分算子计算所述窗口函数中心时刻的卫星运动速度，并根据卫星运动速度计算卫星动能，同时根据卫星的位置和地球自转平均速度计算旋转位；

将所述连续精密轨道结合参考重力场模型计算正常重力位，通过背景模型确定潮汐位，根据所述正常重力位和潮汐位获得保守力位，对所述非保守力加速度通过沿轨积分得到非保守力位；在所述卫星动能中扣除旋转位、保守力位和非保守力位，根据能量守恒定律建立观测方程；针对所述观测方程利用最小二乘准则估算球谐位系数改正值，将所述球谐位系数改正值加上参考重力场模型的位系数确定最终的地球重力场模型。

本发明中的基于轨道掩膜微分算子的地球重力场模型建模方法，基于轨道掩膜微分算子确定卫星运动速度，继而计算卫星动能，再在卫星动能中扣除旋转位、保守力位和非保守力位，最后根据能量守恒定律确定地球重力场位系数。相较于传统利用轨道牛顿多项式微分速度确定地球重力场模型的方法，本发明提出的方法能够有效抑制高频误差放大实现降噪，同时灵活性和计算效率高，获得高精度卫星速度观测值，进一步改进地球重力场模型精度水平。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为本发明中方法的简化流程示意图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

参照图1，本发明提供了基于轨道掩膜微分算子的地球重力场模型建模方法，该方法包括以下步骤：

步骤1，数据预处理：对重力专用卫星的轨道状态信息、星载加速度计数据和原始姿态数据进行数据预处理，获得连续精密轨道和非保守力加速度。

具体地，对轨道状态信息进行内插和粗差探测处理后，得到连续精密轨道；

对连续精密轨道进行时空基准转换，将连续精密轨道转换至惯性系；

对星载加速度计数据和连续精密轨道进行时间同步，同时对原始姿态数据进行线性内插，计算惯性系下的非保守力加速度。

在进行时空基准转换时，采用经典的春分点法，即通过岁差、章动、地球自转和地球极旋四个矩阵实现。

步骤2，利用轨道掩膜微分算子计算卫星动能和旋转位：将连续精密轨道转换至惯性系后，利用掩膜微分算子计算窗口中心时刻的卫星运动速度观测值，进一步得到对应的卫星动能，同时结合卫星的位置和地球自转平均速度计算旋转位。

基于轨道掩膜微分算子确定卫星运动速度观测值，其原因在于有效抑制高频误差放大实现降噪，同时灵活性和计算效率高，获得高精度卫星速度观测值，进一步改进地球重力场模型精度水平。由步骤1获得的连续精密轨道采用轨道掩膜微分算子计算卫星运动速度，其推导过程如下。

重力专用卫星的轨道时间序列采用多项式逼近：

式中，x(t₀+δ_m)_N表示观测时刻卫星的轨道位置矢量；t₀表示计算历元；δ_m表示计算历元与实际观测历元的差值；N表示多项式阶数；c_k表示待定的多项式系数；M表示观测个数；Δt表示观测数据采样间隔；[·]表示向下取整符号。

为降低多项式逼近有限阶次引起的截断误差，观测历元往往关于中心时刻对称分布，要求多项式阶数偶数以及观测个数奇数。此时，对应的逐历元方程组式表述成：

上式的矩阵形式记作：

x_M×1＝B_M×(N+1)c_(N+1)×1

其中x_M×1表示轨道位置向量；B_M×(N+1)表示范德蒙特系数阵；c_(N+1)×1表示待估的多项式系数。

需要说明的是，当观测个数恰好等于多项式阶数，称作插值过程；当观测个数远多于多项式阶数，称作拟合逼近过程。可见，多项式阶数和观测个数的不同组合形式使得算法灵活性高。

利用经典最小二乘平差法解算上式得到：

式中

表示多项式系数估值；P表示观测量的权阵；W_(N+1)×M表示滤波核函数。

对轨道时间序列多项式微分可得：

仅考虑中心时刻，此时δ＝0，则滤波核函数W中的元素γ＝[0 1 0 … 0]^T。

由此可得：

式中c₁表示多项式系数的第二项；υ表示掩膜微分算子。由此可知，卫星运动速度观测值实质上仅与多项式系数c₁有关，换句话说，卫星运动速度值仅是滤波核函数W第2行各元素与观测值的线性组合。

需要说明的是，当多项式阶数N、观测个数M、观测数据采样间隔Δt给定时，此时滤波核函数W为常数阵，使得掩膜微分算子υ各元素亦为常数，因此微分滤波系数仅需计算一次即可，有利于提高计算效率。

因此卫星运动速度矢量的各元素表示为

对应的矩阵形式表述：

式中

表示卫星运动速度值；x表示卫星轨道位置矢量；F表示掩膜微分滤波系数：

数值微分不可避免地导致高频误差放大，该误差与观测数据频率成正比，因此可以增加采样间隔来提高信噪比。也就是说，在上述计算过程中顾及微分间隔，通过上述掩膜微分算子计算卫星运动速度矢量的各元素

：

式中υ_j表示掩膜微分算子元素；δt表示掩膜微分算子微分间隔；Δt表示观测数据采样间隔；τ表示掩膜微分滤波系数非零元素间隔；

表示微分warm-up(预热)数；M表示观测个数。

本步骤中，改进的解算策略采用轨道掩膜微分算子计算卫星运动速度：首先是在计算时刻左右两边各选取相同个数的观测值，构成以计算时刻为中心的窗口函数，然后对连续精密轨道进行最小二乘多项式拟合，得到轨道时间序列的多项式形式，通过线性组合滤波核函数与观测值计算窗口中心时刻的卫星运动速度；随后移动窗口，计算下一时刻，以此类推，组成平滑波谱，从而计算卫星运动速度的整个数据集。

步骤3，基于能量守恒定律估计地球重力场位系数：将连续精密轨道结合参考重力场模型计算正常重力位，同时通过背景模型获得潮汐位，将正常重力位和潮汐位求和得打保守力位，对非保守力加速度通过沿轨积分得到非保守力位。在卫星动能中扣除旋转位、保守力位和非保守力位后，根据能量守恒定律建立观测方程，最后，针对观测方程利用最小二乘准则估算球谐位系数改正值，将该球谐位系数改正值加上参考重力场模型的位系数即可确定最终的地球重力场模型。推导过程如下：

根据能量守恒定律建立的观测方程为：

引力位及其梯度分量：

其中T表示扰动位，E₀表示能量积分常数，U₀表示正常重力位，

表示卫星动能，

表示卫星运动速度矢量，

表示由地球自转引起的旋转位，ω表示地球自转角速度，V_t表示三体问题全部的直接和间接潮汐位，包括日月引力位、固体/海洋/大气、海洋/大气负荷等潮汐影响因素，V_c表示非保守力位，(r,θ,λ)分别表示地固球坐标下的地心向径、地心余纬和地心经度，GM和R分别表示地心引力常数和地球平均半径，l和m分别表示球谐级数展开式的阶和次，L表示确定地球重力场模型的最大阶数，

和

分别表示完全规格化的l阶m次球谐位系数改正值，

表示完全规格化的l阶m次缔合勒让德函数。

需要说明的是，基于轨道掩膜微分算子确定地球重力场模型的方法是先由轨道掩膜微分算子计算卫星运动速度，继而获得卫星动能，再在卫星动能中扣除旋转位、保守力位和非保守力位后在惯性系下根据能量守恒定律建立观测方程，对该观测方程使用最小二乘准则估算球谐位系数改正值，最后在球谐位系数改正值的基础上加上参考重力场模型的位系数即可确定最终的地球重力场模型。

尽管已描述了本发明的优选实施例，但本领域内的技术人员一旦得知了基本创造性概念，则可对这些实施例作出另外的变更和修改。所以，所附权利要求意欲解释为包括优选实施例以及落入本发明范围的所有变更和修改。

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

Claims (5)

1.基于轨道掩膜微分算子的地球重力场模型建模方法，其特征在于，包括以下步骤：

对重力专用卫星的轨道状态信息、星载加速度计数据和原始姿态数据进行数据预处理，获得连续精密轨道和非保守力加速度；

构建以计算时刻为中心的窗口函数，利用掩膜微分算子计算所述窗口函数中心时刻的卫星运动速度，并根据卫星运动速度计算卫星动能，同时根据卫星的位置和地球自转平均速度计算旋转位；

将所述连续精密轨道结合参考重力场模型计算正常重力位，通过背景模型确定潮汐位，根据所述正常重力位和潮汐位获得保守力位，对所述非保守力加速度通过沿轨积分得到非保守力位；在所述卫星动能中扣除旋转位、保守力位和非保守力位，根据能量守恒定律建立观测方程；针对所述观测方程利用最小二乘准则估算球谐位系数改正值，将所述球谐位系数改正值加上参考重力场模型的位系数确定最终的地球重力场模型。

2.如权利要求1所述的基于轨道掩膜微分算子的地球重力场模型建模方法，其特征在于，对所述轨道状态信息、星载加速度计数据和原始姿态数据进行数据预处理的具体步骤为：

对所述轨道状态信息进行内插和粗差探测处理，得到所述连续精密轨道；

对所述连续精密轨道进行时空基准转换，将所述连续精密轨道转换至惯性系；

对所述星载加速度计数据和连续精密轨道进行时间同步，同时对所述原始姿态数据进行线性内插，计算惯性系下的非保守力加速度。

3.如权利要求1所述的基于轨道掩膜微分算子的地球重力场模型建模方法，其特征在于，计算窗口函数中心时刻的卫星运动速度的具体方法为：

在计算时刻左右两边各选取相同个数的观测值，构成以计算时刻为中心的窗口函数，然后对连续精密轨道进行最小二乘多项式拟合，得到轨道时间序列的多项式形式，引入掩膜微分算子后计算窗口函数中心时刻的卫星运动速度；随后移动窗口，计算下一时刻的卫星运动速度，以此类推，计算卫星运动速度的整个数据集。

4.如权利要求3所述的基于轨道掩膜微分算子的地球重力场模型建模方法，其特征在于，利用掩膜微分算子计算卫星速度的方法为：

重力专用卫星的轨道时间序列多项式表示为：

式中，x(t₀+δ_m)_N表示观测时刻卫星的轨道位置矢量，t₀表示计算历元，δ_m表示计算历元与实际观测历元的差值，N表示多项式阶数，c_k表示待定的多项式系数，M表示观测个数，Δt表示观测数据采样间隔，[·]表示向下取整符号；

与上式对应的逐历元方程组式为：

上式的矩阵形式记作：

x_M×1＝B_M×(N+1)c_(N+1)×1

其中x_M×1表示轨道位置向量；B_M×(N+1)表示范德蒙特系数阵；c_(N+1)×1表示待估的多项式系数；

利用经典最小二乘平差法解算上式得到：

式中

表示多项式系数估值，P表示观测量的权阵；W_(N+1)×M表示滤波核函数；

对所述轨道时间序列多项式微分可得：

仅考虑中心时刻，此时δ＝0，则滤波核函数W中的元素γ＝[0 1 0 … 0]^T，由此可得：

式中c₁表示多项式系数的第二项，υ表示掩膜微分算子；

通过上述掩膜微分算子计算卫星运动速度矢量的各元素

式中，δt表示掩膜微分算子微分间隔，τ表示掩膜微分滤波系数非零元素间隔，

表示微分预热数。

5.如权利要求1所述的基于轨道掩膜微分算子的地球重力场模型建模方法，其特征在于，根据能量守恒定律建立的观测方程为：

利用下面引力位及其梯度分量表达式计算所述球谐位系数改正值：

其中T表示扰动位，E₀表示能量积分常数，

表示卫星动能，U₀表示正常重力位，

表示由地球自转引起的旋转位，ω表示地球自转角速度，V_t表示潮汐位，V_c表示非保守力位，(r,θ,λ)分别表示地固球坐标下的地心向径、地心余纬和地心经度，GM和R分别表示地心引力常数和地球平均半径，l和m分别表示球谐级数展开式的阶和次，L表示确定地球重力场模型的最大阶数，

和

分别表示完全规格化的l阶m次球谐位系数改正值，

表示完全规格化的l阶m次缔合勒让德函数。

Images