CN108319566B - 基于向上延拓的航空重力点对点向下延拓解析方法 - Google Patents

Abstract

本发明涉及一种基于向上延拓的航空重力点对点向下延拓解析方法，其主要技术特点是包括：依据泰勒级数展开思想，建立航空重力向上延拓与向下延拓之间的点对点解析延拓模型；利用Poisson积分方程进行航空重力稳向上延拓，基于点对点解析延拓模型进行航空重力向下延拓反的解算。本发明设计合理，实现了航空重力向下延拓反问题的稳定解算，有效地解决了向下延拓解算固有的不适定性问题，同时利用超高阶位模型EGM2008建立的模拟标准场数据对最小二乘向下延拓解析模型解算结果的合理性和有效性进行了数值验证，试验证明本发明实用易行，具有较高的应用价值。

Description

基于向上延拓的航空重力点对点向下延拓解析方法

技术领域

本发明属于海洋重力技术领域，尤其是一种基于向上延拓的航空重力点对点向下延拓解析方法。

背景技术

航空重力位场向上和向下延拓是地球物理和大地测量数据分析应用中最重要的技术环节之一。航空重力位场向下延拓的目的是将数据观测面外推到更接近于位场源体，通过改善观测数据的信噪比，增强地球内部浅层质量分布异常的映射强度，来凸显地球重力场的局部变化特征，以提高位场数据解释推断的可靠性。但航空重力向下延拓在数学上属于典型的不适定反问题，其延拓算子对高频噪声具有明显的放大作用，很小的观测噪声也会引起延拓解算结果严重偏离真实的问题。

发明内容

本发明的目的在于克服现有技术的不足，提出一种基于向上延拓的航空重力点对点向下延拓解析方法，解决航空重力向下延拓的不适定问题。

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

一种基于向上延拓的航空重力点对点向下延拓解析方法，包括以下步骤：

步骤1、依据泰勒级数展开思想，建立航空重力向上延拓与向下延拓之间的点对点解析延拓模型；

步骤2、利用Poisson积分方程进行航空重力稳向上延拓，基于点对点解析延拓模型进行航空重力向下延拓反的解算；

所述步骤2采用如下Poisson积分方程进行航空重力稳向上延拓：

式中，r_p＝R+h_p，r_q＝R+h_q，R为地球椭球平均半径；Δg_pq代表与计算点Q相对应的飞行高度面上的重力异常；Δg_p为飞行高度面流动点的重力异常；

为计算点与流动点之间的空间距离，ψ为计算点与流动点之间的球面角距。

所述步骤1的实现方法为：

设海拔高度为h_o的待求地球表面O点的重力异常为Δg_o，已知对应飞行高度面h_p上的空中重力异常为Δg_p，则由泰勒级数展开模型可得到Δg_o和Δg_p之间的关系为：

式中，Δh_po＝h_p-h_o代表空间P点相对于地面O点的高度差；

代表空中重力异常Δg_p在P点的n阶垂向偏导数；δΔg_po代表重力异常Δg_p到Δg_o的向下延拓改正数；

设位于飞行高度面上方、海拔高度为h_q处的空中重力异常为Δg_q，令Δh_qp＝h_q-h_p，则根据重力场解析延拓理论，将Δg_q和Δg_p之间的关系表示为:

式中，δΔg_pq代表重力异常Δg_p到Δg_q的向上延拓改正数；

取Δh_qp＝Δh_po，将上述两公式相加并处理得到：

其中，R_2n(Δh_po)为向下延拓改正数的余项，设：

R_2n(Δh_po)＝0

则得到如下点对点解析延拓模型：

Δg_o＝2Δg_p-Δg_q

式中，航空重力异常Δg_p是已知的，Δg_q可依据Poisson积分方程由Δg_p向上延拓计算得到。

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

1、本发明依据泰勒级数展开模型，建立航空重力向上延拓与向下延拓之间的点对点解析延拓模型；利用Poisson积分方程进行航空重力稳向上延拓，基于点对点解析延拓模型实现航空重力向下延拓反问题的稳定解算，有效地解决了向下延拓解算固有的不适定性问题。

2、本发明利用超高阶位模型EGM2008建立的模拟标准场数据对最小二乘向下延拓解析模型解算结果的合理性和有效性进行了数值验证，试验证明本发明实用易行，具有较高的应用价值。

具体实施方式

以下通过实施例对本发明做进一步详述。

一种基于向上延拓的航空重力点对点向下延拓解析方法，包括以下步骤：

步骤1、依据泰勒级数展开思想，建立航空重力向上延拓与向下延拓之间的点对点解析延拓模型。

在本步骤中，假设海拔高度为h_o的待求地球表面O点的重力异常为Δg_o，已知对应飞行高度面h_p(P点)上的空中重力异常为Δg_p，则由泰勒级数展开模型可知，Δg_o和Δg_p之间的关系可表示为:

式中，Δh_po＝h_p-h_o代表空间P点相对于地面O点的高度差；

代表空中重力异常Δg_p在P点的n阶垂向偏导数；δΔg_po代表重力异常Δg_p到Δg_o的向下延拓改正数。

假设位于飞行高度面上方、海拔高度为h_q(Q点)处的空中重力异常为Δg_q，令Δh_qp＝h_q-h_p，那么，根据重力场解析延拓理论，同样可以将Δg_q和Δg_p之间的关系表示为:

式中，δΔg_pq代表重力异常Δg_p到Δg_q的向上延拓改正数，其它符号意义同前。取Δh_qp＝Δh_po，那么将式(1)和(2)两式相加并略加整理得：

这里称R_2n(Δh_po)为向下延拓改正数的余项，由于垂向梯度随着阶次的增大而急剧减小，故可忽略该余项的影响，即令：

R_2n(Δh_po)＝0 (5)

那么就得到如下近似的点对点向下解析延拓计算式(点对点解析延拓模型)：

Δg_o＝2Δg_p-Δg_q (6)

式中，航空重力异常Δg_p是已知的，Δg_q可依据Poisson积分方程由Δg_p向上延拓计算得到。

步骤2、利用Poisson积分方程进行航空重力稳向上延拓，基于点对点解析延拓模型进行航空重力向下延拓反的求解。

本步骤采用的Poisson积分方程为：

式中，r_p＝R+h_p，r_q＝R+h_q，R为地球椭球平均半径；Δg_pq代表与计算点Q相对应的飞行高度面上的重力异常；Δg_p为飞行高度面流动点的重力异常；

为计算点与流动点之间的空间距离，ψ为计算点与流动点之间的球面角距。使用式(7)进行向上延拓计算可避开常见的积分奇异性问题。

下面采用超高阶位模型作为标准场开展数值计算检验及分析比较研究。以美国本土西部一个3°×3°区块(

37°N～40°N；λ:248°E～251°E)作为试验区，选用EGM2008位模型模拟产生不同高度面2′×2′网格重力异常及其垂向偏导数的“真值”。选择美国本土作为试验区，是考虑到EGM2008位模型在美国地区具有更高的逼近度。该区块属于地形变化比较剧烈的大山区，试验效果具有一定的代表性。由EGM2008位模型计算空中网格重力异常的公式为：

式中，GM为地球引力常数，R为地球椭球平均半径，r＝R+h，h为计算点高度。

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

和

为完全规格化位系数。连续对式(8)求r的偏导数，可得到不同阶次的垂向偏导数计算模型。

这里取R＝6371km,h_i＝ikm(i＝0,1，…，10)，使用EGM2008模型(以360阶次作为参考场)分别计算11个对应于r_i＝R+h_i球面上的2′×2′网格重力异常“真值”Δg_i(tru)，同时计算5km高度面(r₅＝R+h₅)上的垂向偏导数。表1列出了4个不同高度面上的EGM2008位模型(361～2160阶次)残差重力异常的统计结果，表2列出了5km高度面垂向偏导数计算值的统计结果。本发明采用的试验方案具体设计为：以5km高度面作为航空重力测量的观测面，即假设Δg₅(tru)是已知的观测量，依次采用不同的向下延拓模型，由Δg₅(tru)分别计算4km、3km、2km、1km和0km高度面上的重力异常

将计算值与相对应的“真值”Δg_i(tru)(i＝0,…，4)作比较，从而获得对应延拓模型解算精度的评估参数。

表1 不同高度面位模型残差重力异常统计结果(单位：10^-5m﹒s^-2)

高度(km)	最小值	最大值	平均值	均方根
					0	-107.09	189.19	-0.44	34.28
1	-88.22	162.05	-0.39	30.26
					3	-68.34	120.33	-0.33	24.07
5	-55.22	90.66	-0.28	19.56

表2 5km高度面垂向偏导数统计结果

偏导数	最小值	最大值	平均值	均方根
					一阶(mGal/km)	-3.45985	2.09489	0.00589	0.94238
二阶(mGal/km<sup>2</sup>)	-0.28745	0.49396	0.00049	0.09966
					三阶(mGal/km<sup>3</sup>)	-0.08566	0.04196	-0.00022	0.01449
四阶(mGal/km<sup>4</sup>)	-0.01119	0.01833	0.00006	0.00287

为了检验点对点延拓模型解式(6)的计算效果，首先使用向上延拓积分式(7)，由Δg₅(tru)依据公式(8)分别计算6km、7km、8km、9km和10km高度面上的重力异常

求计算值与相对应“真值”Δg_i(tru)(i＝6，…，10)的互差，可得到向上延拓计算模型的精度评价，具体结果见表3。其中，积分半径统一取为ψ₀＝30′。为了减小积分边缘效应对评估结果的影响，计算区域边缘30′范围内的数据不参加对比分析，这个区域的数据也不再参加下一步的向下延拓计算检验。

表3 向上延拓模型计算精度检核(单位：10^-5m˙s^-2)

由表3结果看出，当观测数据不存在噪声干扰时，向上延拓模型计算结果可以达到很高的精度水平，即使是1km的延拓高度差，计算误差也不超过0.3mGal。分别将前面计算得到的6个高度面的向上延拓重力异常

和飞行高度面上的观测重力异常Δg₅(tru)代入向下延拓计算式(6)，可求得相对应延拓高度差下的向下延拓重力异常：

求上述计算值

相对应“真值”Δg_i(tru)的互差，可得到向下延拓计算式(6)的精度评价，具体结果见表4。

表4 点对点向下延拓模型计算精度检核(单位：10^-5m˙s^-2)

高度差(km)	最小值	最大值	平均值	均方根
					1	-0.34	0.28	-0.01	0.10
2	-1.86	1.01	-0.01	0.36
					3	-4.45	2.50	-0.01	0.88
4	-8.19	4.61	-0.01	1.61
					5	-13.21	7.40	-0.02	2.59

表4中延拓高度差为ikm对应于计算面的高度为(5-i)km。由表4结果看出，虽然向上延拓模型的计算误差随计算高度的增大而减小，但点对点向下延拓模型的计算精度随延拓高度差的增大而降低，说明向下延拓模型的截断误差对延拓计算结果有较大影响。当延拓高度差小于3km时，该模型的计算精度优于1mGal。

为了进一步考察数据观测噪声对向下延拓计算结果的影响，这里在模拟观测量Δg₅(tru)中分别加入±1mGal、±3mGal和±5mGal的白噪声，对应生成三组带噪声的观测量，并重复前面的计算和比对检核过程。加入观测噪声后，与表3相对应的向上延拓模型计算精度检核结果(这里只列出互差均方根值(rms)，下同)如表5所示，与表4相对应的向下延拓模型计算精度检核结果如表6所示。

表5 噪声影响下的向上延拓误差均方根值(单位：10^-5m˙s^-2)

高度(km)	6	7	8	9	10
						噪声＝1mGal	0.26	0.21	0.17	0.14	0.12
噪声＝3mGal	0.70	0.56	0.45	0.37	0.31
						噪声＝5mGal	1.18	0.95	0.78	0.65	0.54

表6 噪声影响下的点对点向下延拓误差均方根值(单位：10^-5m˙s^-2)

高度差(km)	1	2	3	4	5
						噪声＝1mGal	0.25	0.41	0.89	1.61	2.58
噪声＝3mGal	0.70	0.67	1.00	1.67	2.62
						噪声＝5mGal	1.16	0.97	1.10	1.66	2.58

对比表3和表5数值结果可以看出，数据观测噪声对向上延拓计算结果的影响很小，几乎不改变计算参数的变化形态，说明向上延拓模型作为低通滤波器的作用相当明显。进一步对比表4和表6结果可以看出，得益于向上延拓模型抑制高频噪声的优良特性，点对点向下延拓模型计算结果也几乎不受数据观测噪声的影响，在延拓高度差小于3km的高度内，仍可获得1mGal左右的计算精度，证明本发明给出的技术方案实用易行，具有较高的应用价值。

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

Claims (2)

1.一种基于向上延拓的航空重力点对点向下延拓解析方法，其特征在于包括以下步骤：

步骤1、依据泰勒级数展开思想，建立航空重力向上延拓与向下延拓之间的点对点解析延拓模型；

步骤2、利用Poisson积分方程进行航空重力稳向上延拓，基于点对点解析延拓模型进行航空重力向下延拓反的解算；

所述步骤2采用如下Poisson积分方程进行航空重力稳向上延拓：

式中，r_p＝R+h_p，r_q＝R+h_q，R为地球椭球平均半径；Δg_pq代表与计算点Q相对应的飞行高度面上的重力异常；Δg_p为飞行高度面流动点的重力异常；

为计算点与流动点之间的空间距离，ψ为计算点与流动点之间的球面角距。

2.根据权利要求1所述基于向上延拓的航空重力点对点向下延拓解析方法，其特征在于：所述步骤1的实现方法为：

设海拔高度为h_o的待求地球表面O点的重力异常为Δg_o，已知对应飞行高度面h_p上的空中重力异常为Δg_p，则由泰勒级数展开模型可得到Δg_o和Δg_p之间的关系为：

式中，Δh_po＝h_p-h_o代表空间P点相对于地面O点的高度差；

代表空中重力异常Δg_p在P点的n阶垂向偏导数；δΔg_po代表重力异常Δg_p到Δg_o的向下延拓改正数；

设位于飞行高度面上方、海拔高度为h_q处的空中重力异常为Δg_q，令Δh_qp＝h_q-h_p，则根据重力场解析延拓理论，将Δg_q和Δg_p之间的关系表示为:

式中，δΔg_pq代表重力异常Δg_p到Δg_q的向上延拓改正数；

取Δh_qp＝Δh_po，将上述两公式相加并处理得到：

其中，R_2n(Δh_po)为向下延拓改正数的余项，设：

R_2n(Δh_po)＝0

则得到如下点对点解析延拓模型：

Δg_o＝2Δg_p-Δg_q

式中，航空重力异常Δg_p是已知的，Δg_q可依据Poisson积分方程由Δg_p向上延拓计算得到。