CN112949049B - 一种利用带限思想计算重力异常低阶径向导数的方法 - Google Patents

Abstract

本发明涉及一种利用带限思想计算重力异常低阶径向导数的方法，其主要技术特点是：将地球外部重力异常Poisson积分式的解析核函数表示为球谐级数展开式，通过求导方法得到重力异常低阶径向导数球谐级数表达式；采用移去‑恢复技术将重力异常低阶径向导数球谐级数计算式截断为与重力观测值频谱范围相一致的带限求和式，并利用全球重力位模型补偿重力异常低阶径向导数远区截断误差，计算得到重力异常低阶径向导数。本发明设计合理，能够精确地计算重力异常低阶径向导数，提高了重力异常低阶径向导数的计算精度，解决了重力异常低阶径向导数计算结果的不稳定性问题，可广泛用于物理大地测量领域。

Description

一种利用带限思想计算重力异常低阶径向导数的方法

技术领域

本发明属于物理大地测量技术领域，尤其是一种利用带限思想计算重力异常低阶径向导数的方法。

背景技术

重力异常径向导数在大地测量学和地球物理勘探中具有非常重要的应用价值。高精度的重力异常低阶径向导数(即一阶径向导数、二阶径向导数和三阶径向导数)是反演近地表异常物体和确定矿产资源分布规律，同时也是推求地球内部重力和潜器水下重力辅助导航必不可少的基础资料。

计算重力异常低阶径向导数的解析核函数在球面边界存在奇异性，导致计算结果不稳定。位场理论的球谐表达式具有较好的稳定性，为克服计算重力异常低阶径向导数的解析核函数的奇异性问题提供了新途径。在实际应用中受观测数据覆盖范围限制，无法做到全球覆盖，实际计算过程中重力异常低阶径向导数的全球积分式需要进行适用观测数据保障条件的改化，以保证计算结果的可靠性。

综上所述，如何克服重力异常低阶径向导数计算结果的奇异性问题以提高重力异常低阶径向导数的计算精度是目前迫切需要解决的问题。

发明内容

本发明的目的在于克服现有技术的不足，提出一种设计合理、精度高且稳定性强的利用带限思想计算重力异常低阶径向导数的方法。

本发明解决其技术问题是采取以下技术方案实现的：

一种利用带限思想计算重力异常低阶径向导数的方法，包括以下步骤：

步骤1、将地球外部重力异常Poisson积分式的解析核函数表示为球谐级数展开式，通过求导方法得到重力异常低阶径向导数球谐级数表达式；

步骤2、采用移去-恢复技术将重力异常低阶径向导数球谐级数计算式截断为与重力观测值频谱范围相一致的带限求和式，并利用全球重力位模型补偿重力异常低阶径向导数远区截断误差，计算得到重力异常低阶径向导数。

而且，所述步骤1得到的重力异常低阶径向导数球谐级数表达式为：

式中，Δg为外部空间计算点

重力异常；Δg_R为球面上流动点

处的已知观测重力异常；R为地球椭球平均半径；r为计算点地心向径；

为计算点的纬度和经度；

为流动点的纬度和经度；σ为单位球面；dσ为单位球面的面积元；ψ为计算点至流动点之间的球面角距；

是计算点至流动点之间的空间距离；K(r,ψ)为积分核函数；P_n(cosψ)为n阶勒让德多项式级数。

而且，所述步骤2的具体实现方法为：利用移去-恢复技术，从重力异常观测值中移去参考重力异常，得到残差重力异常；从积分式核函数中移去与参考场对应阶次的核函数球谐表达式，得到截断核函数，使之与残差重力异常的频谱匹配；基于截断核函数和残差重力异常的局域积分得到残差重力异常低阶径向导数；利用全球重力场位模型高阶信息进行远区效应补偿；恢复参考重力异常低阶径向导数；对步骤1得到的重力异常低阶径向导数进行改化处理，得到重力异常低阶径向导数的计算公式并进行计算。

而且，所述重力异常低阶径向导数的计算公式为：

式中，δΔg_R为Δg_R的残差重力异常；g_1pref、g_2pref和g_3pref分别是参考重力异常一阶、二阶和三阶径向导数；

和

分别是重力异常一阶、二阶和三阶径向导数的改化核函数；

和

分别是一阶、二阶和三阶径向导数的远区补偿。

本发明的优点和积极效果是：

本发明设计合理，其将地球外部重力异常Poisson积分式的解析核函数表示为球谐级数展开式，并通过求导方法得到重力异常低阶径向导数球谐级数表达式；利用移去-恢复技术，将重力异常低阶径向导数球谐级数计算式截断为与重力观测值频谱范围相一致的带限求和式，同时利用全球重力位模型补偿重力异常低阶径向导数远区截断误差，精确地计算重力异常低阶径向导数，提高了重力异常低阶径向导数的计算精度，解决了重力异常低阶径向导数计算结果的不稳定性问题。

具体实施方式

一种利用带限思想计算重力异常低阶径向导数的方法，包括以下步骤：

步骤1、将地球外部重力异常Poisson积分式的解析核函数表示为球谐级数展开式，通过求导方法得到重力异常低阶径向导数球谐级数表达式。

在本步骤中，地球外部重力异常Poisson积分式的解析计算式为：

式中，Δg为外部空间计算点

重力异常；Δg_R为球面上流动点

处的已知观测重力异常；R为地球椭球平均半径；r为计算点地心向径；

为计算点的纬度和经度；

为流动点的纬度和经度；σ为单位球面；dσ为单位球面的面积元；ψ为计算点至流动点之间的球面角距；

是计算点至流动点之间的空间距离；K(r,ψ)为积分核函数。

重力异常低阶径向导数可表示为：

将解析核函数K(r,ψ)用球谐级数展开式表示为：

式中P_n(cosψ)为n阶勒让德(Legendre)多项式级数。

对式(4)求低阶径向导数，并令r＝R，可得球面上的核函数径向导数为：

将式(5)代入式(3)，可得到球谐级数展开表达的重力异常低阶径向导数积分计算式：

步骤2、依据各类重力异常观测经滤波处理后均表现为一类有限频谱带宽信号的特性，并顾及远区效应影响，采用移去-恢复技术，将重力异常低阶径向导数球谐级数计算式截断为与重力观测值频谱范围相一致的带限求和式，同时利用全球重力位模型补偿重力异常低阶径向导数远区截断误差，计算得到重力异常低阶径向导数。具体实现方法如下：

由于公式(6)表达的重力异常低阶径向导数计算模型是理论计算式，计算重力异常低阶径向导数要求全球积分，但在实际应用中受观测数据覆盖范围限制，无法做到全球覆盖，实际计算过程中重力异常低阶径向导数的全球积分式需要进行适用观测数据保障条件的改化，以保证计算结果的可靠性。

在本步骤中，顾及实测数据局域保障条件，引入全球重力场位模型，利用移去-恢复技术，从重力异常观测值中移去参考重力异常，以得到残差重力异常；从积分式核函数中移去与参考场对应阶次的核函数球谐表达式，以得到截断核函数，使之与残差重力异常的频谱匹配；基于截断核函数和残差重力异常的局域积分得到残差重力异常低阶径向导数；利用全球重力场位模型高阶信息进行远区效应补偿，以削弱远区截断误差的影响；恢复参考重力异常低阶径向导数。对式(6)改化处理，得到利用带限思想计算重力异常低阶径向导数的公式，具体表达式为：

式中δΔg_R为Δg_R的残差重力异常；g_1pref、g_2pref和g_3pref分别是参考重力异常一阶、二阶和三阶径向导数；

和

分别是重力异常一阶、二阶和三阶径向导数的改化核函数；

和

分别是一阶、二阶和三阶径向导数的远区补偿。

式(7)中，残差重力异常δΔg_R的计算式为：

δΔg_R＝Δg_R-Δg_Rref (8)

式中，Δg_Rref是球面上的参考重力异常，具体计算式为：

式中L是参考场位模型阶次；Δg_n(θ,λ)是重力异常n阶拉普拉斯面球谐函数，具体表达式为：

式中(θ,λ)是计算点的余纬和经度；GM为地球引力常数；

为完全规格化缔合勒让德函数；

和

为完全规格化地球位系数，其余符号意义同前。

式(7)中参考重力异常一阶、二阶和三阶径向导数g_1pref、g_2pref和g_3pref的具体表达式为：

式中符号意义同前。

式(7)中重力异常一阶、二阶和三阶径向导数的改化核函数

和

的具体表达式为：

式中N为与重力异常观测值滤波尺度对应的最高频谱阶次，其他符号意义同前。

式(7)中一阶、二阶和三阶径向导数的远区补偿

和

的具体表达式为：

其中：

下面采用全球重力场位模型EGM2008作为标准场开展数值计算检验及分析比较。

选取重力场变化比较剧烈的马里亚纳海沟作为试验区，具体覆盖范围为：(

11°N～14°N；λ:143°E～146°E)。选取r＝R+h，R＝6371km，由EGM2008模型分别计算3组分别对应于h₀＝0km、h₆＝6km、h₁₀＝10km高度面上的1′×1′网格剩余重力异常“真值”

(i对应0km，6km，10km)和剩余重力异常一阶、二阶和三阶径向导数“真值”

(j＝1，2，3)。表1列出了3组不同高度面1′×1′网格剩余重力异常“真值”的统计结果，表2列出了相对应3组高度面一阶、二阶和三阶径向导数“真值”的统计结果。

表1不同高度面EGM2008模型重力异常统计结果(mGal)

高度面(km)	最小值	最大值	平均值	均方根值
					0	-78.48	132.75	-0.05	26.36
6	-41.18	74.21	-0.04	16.22
					10	-30.45	52.29	-0.04	12.00

表2不同高度面EGM2008模型重力异常径向偏导数统计结果

为了检验本发明的计算效果，首先以前面选定的3个高度面(h₀＝0km、h₆＝6km、h₁₀＝10km)上的位模型剩余重力异常

作为观测量，依据本发明公式(7)计算相对应高度面上的1′×1′网格一阶、二阶和三阶径向偏导数

将计算值

与相对应的“真值”

作比较，可获得不同高度面不同阶次偏导数计算模型的精度评估信息，具体比对统计结果列于表3。这里积分半径统一取为ψ₀＝0.5°，故计算区域外围0.5°范围内的比对结果不参加精度评估统计计算(下同)。为了对比分析评价径向偏导数积分模型改化前后的计算效果，本试验同时给出了采用传统算法(公式3)完成相同参量计算获得的精度评估结果，具体见表4。

表3利用本发明计算不同高度面重力异常低阶径向导数精度检核

表4利用传统算法计算不同高度面重力异常低阶径向导数精度检核

从表3检核结果可以看出，依据本发明计算重力异常一阶至三阶径向导数，均可获得比较满意的符合精度。从表3并结合表2统计结果可以看出，本发明的绝对精度(互差均方根值)都随计算高度面的增高和偏导数阶数的增大而提升，其相对精度(互差均方根值/径向导数均方根值)的变化趋势则正好相反，均随计算高度面的增高和偏导数阶数的增大而降低。计算高度面越高，其相对精度的降低幅度越明显。这个结果显然跟高度面越高，高阶偏导数的绝对量值越小有关，同时跟偏导数阶数越高，积分模型离散化误差的影响越大有关，这是符合理论分析预期的结果。进一步对比表3和表4统计结果可知，本发明结果精度明显优于传统算法，表明本发明实用易行，具有较高的应用价值。

需要强调的是，本发明所述的实施例是说明性的，而不是限定性的，因此本发明包括并不限于具体实施方式中所述的实施例，凡是由本领域技术人员根据本发明的技术方案得出的其他实施方式，同样属于本发明保护的范围。

Claims (1)

1.一种利用带限思想计算重力异常低阶径向导数的方法，其特征在于：包括以下步骤：

步骤1、将地球外部重力异常Poisson积分式的解析核函数表示为球谐级数展开式，通过求导方法得到如下重力异常低阶径向导数球谐级数表达式；

式中，Δg为外部空间计算点

重力异常；Δg_R为球面上流动点

处的已知观测重力异常；R为地球椭球平均半径；r为计算点地心向径；

为计算点的纬度和经度；

为流动点的纬度和经度；σ为单位球面；dσ为单位球面的面积元；ψ为计算点至流动点之间的球面角距；

是计算点至流动点之间的空间距离；K(r,ψ)为积分核函数；P_n(cosψ)为n阶勒让德多项式级数；

步骤2、采用移去-恢复技术将重力异常低阶径向导数球谐级数计算式截断为与重力观测值频谱范围相一致的带限求和式，并利用全球重力位模型补偿重力异常低阶径向导数远区截断误差，计算得到重力异常低阶径向导数；

所述步骤2的具体实现方法为：利用移去-恢复技术，从重力异常观测值中移去参考重力异常，得到残差重力异常；从积分式核函数中移去与参考场对应阶次的核函数球谐表达式，得到截断核函数，使之与残差重力异常的频谱匹配；基于截断核函数和残差重力异常的局域积分得到残差重力异常低阶径向导数；利用全球重力场位模型高阶信息进行远区效应补偿；恢复参考重力异常低阶径向导数；对步骤1得到的重力异常低阶径向导数重力异常低阶径向导数球谐级数表达式进行处理，得到如下重力异常低阶径向导数的计算公式并进行计算；

式中，δΔg_R为Δg_R的残差重力异常；g_1pref、g_2pref和g_3pref分别是参考重力异常一阶、二阶和三阶径向导数；

和

分别是重力异常一阶、二阶和三阶径向导数的改化核函数；

和

分别是一阶、二阶和三阶径向导数的远区补偿。