地球物理学报  2020, Vol. 63 Issue (8): 3131-3143   PDF    
融合多源数据构建区域航空重力梯度扰动全张量
刘金钊1,2, 梁星辉2, 叶周润3, 刘站科4, 郎骏健2, 王国成2, 柳林涛2     
1. 中国地震局第一监测中心, 天津 300180;
2. 大地测量与地球动力学国家重点实验室, 中国科学院精密测量科学与技术创新研究院, 武汉 430077;
3. 合肥工业大学, 合肥 230009;
4. 自然资源部第一大地测量队, 西安 710054
摘要:通过联合全球重力位模型(EGM2008)、航空重力扰动数据和剩余地形模型(RTM)数据,基于频谱域(二维FFT变换)和空间域(Stokes数值积分)算法对毛乌素测区GT-2A航空重力测量系统采集的空中测线后处理重力扰动数据进行解算,构建了该地区的航空重力梯度扰动全张量.(1)残余航空重力扰动延拓结果表明:残余航空重力扰动经向下延拓至大地水准面,再向上延拓至航空高度后与原数据差值的标准差为1.0078 mGal,考虑边缘效应后,内缩计算范围得到的差值标准差减小至0.1269 mGal.(2)基于残余重力扰动数据(原航空高度数据及向下延拓数据),通过不同方案解算得到的梯度扰动结果表明:两种方案得到的研究区域重力梯度扰动各分量之差的最大标准差为6.4798E(Γyz分量),最小标准差为2.6968E(Γxy分量),内缩计算范围后得到的差值标准差最大值为1.8307E(Γzz分量),最小值为0.7223E(Γyz分量).本文的思路和方法可为未来我国自主构建航空重力梯度标定场提供参考.
关键词: 航空重力扰动      区域重力梯度标定场      二维FFT变换      Stokes数值积分      EGM2008      剩余地形模型     
Combining multi-source data to construct full tensor of regional airborne gravity gradient disturbance
LIU JinZhao1,2, LIANG XingHui2, YE ZhouRun3, LIU ZhanKe4, LANG JunJian2, WANG GuoCheng2, LIU LinTao2     
1. The First Monitoring and Application Center, China Earthquake Administration, Tianjin 300180, China;
2. State Key Laboratory of Geodesy and Earth's Dynamics, Innovation Academy for Precision Measurement Science and Technology. Chinese Academy of Science, Wuhan 430077, China;
3. Hefei University of Technology, Hefei 230009, China;
4. The First Geodetic Survey Brigade Ministry of Natural Resources, Xi'an 710054, China
Abstract: In this paper, by combining EGM2008, airborne gravity disturbance and Residual Terrain Model (RTM) data, we calculated the full tensor of gravity gradients disturbance over Mu US surveying area in Inner Mongolia. The gravity disturbance data are post-processed from flight line gravity data collected by GT-2A airborne gravity measurement system. The spectral domain (2D FFT) and spatial domain (Stokes' numerical integration) algorithms have been applied to computing the full tensor gravity gradients disturbance respectively. (1) The continuation results of the residual airborne gravity disturbance show that: the standard deviation of the difference between the original residual airborne gravity disturbance and its counterpart which firstly continued downward to the geoid surface and then continued upward back to the aviation height is 1.0078 mGal. In order to alleviate the edge effects, the standard deviation of the difference obtained from the inner calculation range is reducing to 0.1269 mGal. (2) Based on the two type residual gravity disturbance data (the original and the downward continuation data): The maximum standard deviation of the differences of the gradients disturbance components calculated through different computational schemes is 6.4798E (difference of Γyz component), and the minimum standard deviation is 2.6968E (difference of Γxy component). The maximum standard deviation of the difference obtained from the inner calculation range is reducing to 1.8307E (difference of Γzz component), and the minimum standard deviation is reducing to 0.7223E (difference of Γyz component). The methods and results proposed in this paper could provide a reference for the future construction on airborne gravity gradient calibration field site in China independently.
Keywords: Airborne gravity disturbance    Gravity gradients calibration field site    2D FFT    Stokes' numerical integration    EGM2008    RTM    
0 引言

相对于传统的重力场数据,重力梯度数据含有更多的高频信息量,在油气矿产资源勘探、大地测量学、地球物理学以及国防建设中具有较强的优越性和应用价值(Bell, 1998; 曾华霖, 1999; Hatch, 2004; 张昌达, 2005; 边少锋和纪兵, 2006; 舒晴等, 2007; Dransfield, 2007; DiFrancesco et al, 2009; 熊盛青, 2009; 叶周润, 2010; 梁星辉, 2012; 袁园, 2015刘金钊等,2015; 叶周润等,2015).

1886年,匈牙利科学家Loránd Eötvös设计了最早的重力梯度测量系统——扭称,随后几十年,北美、加拿大、澳大利亚、德国、意大利和英国等多个国家对重力梯度测量仪器及其关键技术进行了研制与实验,并于20世纪70到80年代研制成功了基于旋转加速度计技术的航空重力梯度测量系统以及于2007年研制成功的具有更高测量精度的超导重力梯度测量系统,同时基于冷原子技术的重力梯度测量系统也在不断测试发展中(表 1)(Murphy, 2004; Lumley et al., 2004; 罗嗣成, 2007; Anstie et al., 2010; Dransfield et al, 2010; Kieran et al, 2000; 梁星辉, 2012);2009年3月欧洲航天局成功发射了基于重力梯度测量技术的GOCE卫星,使得利用卫星重力梯度数据优化全球重力场模型成为大地测量研究领域的前沿和热点问题之一(Klees et al., 2000; 罗志才等, 2009; 郑伟等, 2010; 刘晓刚, 2011; Pail et al., 2011; 游为等, 2012).

表 1 目前主流航空重力梯度测量系统相关技术指标(Lee., 2001; Dransfield et al., 2004; Hinks et al., 2004; Murphy, 2010; Annecchione et al., 2007) Table 1 Related technical indicators of mainstream airborne gravity gradient measurement systems

因此及时开展重力梯度方面的深入研究能够更好地促进我国重力及重力梯度测量关键技术的向前发展,而高分辨率区域重力梯度场建模研究不仅能够为未来我国独立自主设计的重力梯度测量系统的外部检核提供条件,同时对于相关学科(如恢复中高频段重力场信号、区域重力场模型优化及大地水准面精化等;浅层地壳结构精细反演;水下潜艇重力梯度匹配导航、导弹精密轨道控制等)的探索研究及国民经济发展等也具有积极的科学意义(刘金钊,2017).

目前,解算区域重力梯度场主要是基于地形高程数据(边少锋和张赤军, 1999; 张赤军, 1999; Jekeli and Zhu, 2006; Rózsa and Tóth, 2007; Makhloof and Ilk, 2008; 武凛等, 2009; Eshagh et al, 2009; 刘繁明等, 2013)、重力数据(Mickus and Hinojosa, 2001; Zhu and Jekeli, 2006; 蒋甫玉等, 2011; Wang et al., 2012; Jiang et al., 2012; 刘金钊等, 2015)或者二者的联合(Rózsa and Tóth, 2005; Zhu and Jekeli, 2007, 2009; 钱东等, 2011; 刘金钊等, 2013; Bucha et al., 2016),利用频谱域或者空间域方法来实现.本文通过联合全球重力位模型(EGM2008)、航空重力扰动数据和剩余地形模型(RTM)数据,采用不同的解算方案,分别基于频谱域(二维FFT变换)和空间域(Stokes数值积分)算法,对毛乌素测区GT-2A航空重力测量系统采集的空中测线后处理重力扰动数据进行解算,构建了该地区的航空重力梯度扰动全张量.并对两种方案的解算结果进行了对比分析.

1 数据和方法

2015年4月17日至6月11日,自然资源部第一大地测量队采用引进的国内首台GT-2A航空重力仪在毛乌素地区实施了航空重力测量(宛家宽等,2017).毛乌素测区位于陕西省与内蒙古交界处,地貌以沙漠丘陵为主,地形起伏如图 1,高程统计见表 2.航空重力测量载体飞机选用塞斯纳208机型,机上搭载GT-2A航空重力测量系统和1台Ashtech双频GPS接收机,地面架设了2个GPS基准站.选取鄂尔多斯机场作为飞机起降保障机场,机场距毛乌素测区的飞行距离约为80 km(蒋涛等,2018).图 1中的黑色方框为本文计算结果的输出区域,该区域大小约139 km×60 km;灰色线条为测区飞行航线格网,航线格网间隔约0.55 km,红色方框为考虑边缘效应内缩计算范围后的结果输出区域,大小约137 km×58 km.解算区域外地形质量影响结果分析见附录.

图 1 研究区域SRTM地形起伏图(黑色方框为结果输出区域,灰色线条为测区飞行航线格网,红色方框为考虑边缘效应内缩计算范围后的结果输出区域) Fig. 1 SRTM topographic relief map of the study area (the black box in the figure is the result output area, the gray line is the grid of the flight route within the survey area, The red box is the result output area after shrinking the calculation range after considering the edge effect)
表 2 研究区域SRTM数据地形高程统计(单位:m) Table 2 Statistics of SRTM terrain elevation data (Units: meters)

本文航空重力梯度场构建技术路线见图 2:首先在实测航空重力扰动中去除全球重力位模型EGM2008的重力扰动及剩余地形模型RTM的航空高度重力效应的影响,得到残余航空重力扰动数据,作为后续航空重力梯度场解算的输入数据之一;一方面,根据频谱域方法(二维FFT变换),沿方案一进行梯度解算;另一方面,将残余航空重力扰动向下延拓至大地水准面得到大地水准面上的残余重力扰动,由Stokes积分得到方案二梯度解算结果.具体解算过程及结果对比分析见后续章节.

图 2 本文航空重力梯度场构建技术路线图 Fig. 2 Technological flowchart of construction for airborne gravity gradiometry calibration field
1.1 航空重力扰动数据预处理

实测航空重力扰动数据经插值后得到规则格网数据,这里选择的格网化插值方法是反距离加权插值,由于航空重力测线数据本身分布密集合理,简单的反距离加权插值即可达到较好的格网化效果,另外该方法对数据的处理较为简单,对数据的干预少,更能保持数据本身的特性(蒋涛等,2018).对85019个航空重力扰动数据进行格网化后,得到空间分辨率为0.3'×0.3'的重力扰动格网,相关统计量见表 3.

表 3 研究区域航空重力扰动统计(单位: mGal, 1 mGal=1×10-5ms-2) Table 3 Statistics of airborne gravity disturbance(Units: mGal, 1 mGal=1×10-5ms-2)
1.1.1 全球重力位模型解算重力扰动

基于EGM2008模型来解算研究区域的重力扰动,可以根据公式(1)来进行计算(Barthelmes,2013):

(1)

式中:GM是地心引力常数,r是航空高度到地心的径向距离,θλ为球坐标下的地心余纬和经度,CnmTSnmT是完全正则化的EGM2008扰动位球谐系数,Nmax是最大展开阶数,Pnmmn次完全正则化勒让德函数.解算得到的EGM2008航空重力扰动相关统计量见表 4.

表 4 研究区域EGM2008重力扰动统计(单位: mGal) Table 4 Statistics of EGM2008 gravity disturbance calculated in the study area (Units: mGal)
1.1.2 剩余地形模型解算重力效应

剩余地形模型RTM(Residual Terrain Model, Forsberg, 1984)通过SRTM高程数据(Jarvis et al., 2008)和DTM2006.0全球地形球谐展开模型(Pavlis et al., 2007)的差值来计算,得到的RTM数据根据公式(2)来求解相应的航空高度处的地形重力效应(Wang, 1990):

(2)

其中:是航空高度处计算点(x, y, z)与地形积分点(x′, y′, z′)的水平距离.计算得到的RTM重力效应见表 5.

表 5 研究区域RTM重力统计(单位: mGal) Table 5 Statistics of gravity effect which derived from RTM (Units: mGal)
1.1.3 残余航空重力扰动

从格网化后的航空重力扰动数据中减去EGM2008得到航空高度处重力扰动及RTM的地形重力效应,就得到了残余航空重力扰动数据,见表 6.作为本文利用不同解算方案构建区域梯度场的输入数据之一.

表 6 研究区域残余航空重力扰动统计(单位: mGal) Table 6 Statistics of residual airborne gravity disturbance (Units: mGal)
1.2 航空重力梯度扰动解算方案一

根据1.1节得到的残余航空重力扰动数据,我们基于两种不同的解算方案(图 2),分别利用频谱域算法和空间域算法求解最后的梯度张量.

解算方案一选择频谱域方法,基于重力数据与重力梯度的频谱关系,通过二维傅里叶变换FFT来获得研究区域的残余航空重力梯度扰动,然后恢复EGM2008和RTM的梯度信号,最终获得完整的航空重力梯度场.

1.2.1 基于FFT方法计算航空高度残余梯度扰动

在地球外部,重力位ϕ满足拉普拉斯方程,▽2ϕ=0.因此,重力位的傅里叶变换满足(kx, ky, kz)Φ(k)=0,其中k=(kx, ky, kz)是波数向量,kxkykz分别是xyz方向上的波数. kzkxky的关系可以表示成-ikz=|k|,其中|k|=(kx2+ky2)1/2(参考Blakely, 1996).最后,在频谱域,重力梯度与重力扰动的关系可以表示为

(3)

其中,分别是傅里叶变换和傅里叶逆变换,Gz是重力扰动gz(x, y)的二维傅里叶变换,

(4)

基于残余航空重力数据,利用二维FFT方法解算得到残余航空梯度扰动见图 3.

图 3 研究区域二维FFT方法解算的残余航空梯度扰动 Fig. 3 Residual airborne gradient disturbance calculated by 2D-FFT
1.2.2 全球重力位模型解算重力梯度扰动

这里选择的全球重力位模型是EGM2008.在球坐标系(r, θ, λ)下, 重力梯度扰动δΓiji, j∈{x, y, z},由公式(5)计算得到(Barthelmes,2013):

(5)

其中,.式中:GM是地心引力常数,r是航空高度到地心的径向距离,θλ为球坐标下的地心余纬和经度,CnmTSnmT是完全正则化的EGM2008扰动位球谐系数,Pnmmn次完全正则化勒让德函数.解算得到的EGM2008航空重力梯度扰动相关统计量见表 7.

表 7 研究区域重力位模型重力梯度扰动(单位:E,1E=10-9s-2) Table 7 Statistics of calculated EGM2008 gradients disturbance
1.2.3 剩余地形模型解算重力梯度效应

剩余地形模型RTM的重力梯度效应可根据棱柱数值积分法来求解,将格网RTM点看作具有常密度的棱柱体(这里我们选择标准岩石密度ρ=2670 kg·m-3),通过公式(6)至公式(11)求解各梯度分量(参考Jekeli and Zhu, 2006; Forsberg, 1984):

(6)

(7)

(8)

(9)

(10)

(11)

每个RTM点棱柱的尺寸定义为Δx, ΔyzRTM,解算得到的剩余地形模型的重力梯度效应见表 8.

表 8 研究区域剩余地形模型(RTM)重力梯度(单位:E) Table 8 Statistics of gradients derived from RTM (Units: E)
1.3 航空重力梯度扰动解算方案二 1.3.1 残余航空重力扰动向下延拓至大地水准面

位场数据的向下延拓是一个数据噪声得到放大的不适定过程,无论采用何种正则化解法或平滑方法都只是在一定程度上削弱或控制向下延拓误差(蒋涛等,2018).同时,位场向下延拓也能够增强异常的短波信号,或者突出场源分布的细节,是位场数据解释中一项非常有用的技术.这里我们选择Van Tran和Nguyen(2020)的位场延拓方法,该算法的基本思想是基于泰勒级数展开方法,通过不同高度可靠的向上延拓的位场组合来表示不同阶位场梯度,进而求得位场的向下延拓值,该算法在含有噪声的情况下,比希尔伯特方法和拉普拉斯方法更加稳定.

位场的向上延拓可以用公式(12)来表示,通过表示不同高度向上延拓方程,将位场的各阶导数(f, f′, f″, f(3), …)表示成不同高度的幂(Δh, (Δh)2, (Δh)3, …)及其延拓值(f(x, y, zh), f(x, y, z+2Δh), f(x, y, z+3Δh), …)的组合,进而带入向下延拓公式(13)及拓展式(14),完成位场的不同高度向下延拓值:

(12)

(13)

一般形式,将zz-(m-1)Δh代替:

(14)

向下延拓得到的研究区域大地水准面上残余航空重力扰动相关统计量见表 9.利用该数据作为我们航空重力梯度扰动解算方案二的输入数据.同时,我们又将大地水准面上的残余航空重力扰动向上延拓回航空高度,与原残余航空重力扰动数据进行了对比,分析该算法的内符合精度,相关结果分析见第2节.

表 9 研究区域残余航空重力扰动向下延拓统计(单位: mGal) Table 9 Statistics of residual gravity disturbance after downward continuation (Units: mGal)
1.3.2 基于Stokes积分计算航空高度残余梯度扰动

因为本文研究的梯度扰动主要是局部信号,在局部大地水准面足够平滑的情况下,我们经常利用平面来代替大地水准面进行计算.在平面近似下,Stokes函数能够在区域笛卡尔直角坐标系下近似表示为(Zhu and Jekeli, 2007)

(15)

其中,该公式含包了通过向上延拓来评估航空高度处的重力梯度.

平面近似下的Stokes积分公式为

(16)

其中σ是积分单元,,梯度可以表示成

(17)

进而,梯度张量的对角线分量可以表示为

(18)

非对角线分量可以表示为

(19)

实际数值计算中,可以将积分转化成求和来进行计算:

(20)

(21)

基于向下延拓至大地水准面上的残余重力数据,利用Stokes数值积分方法解算得到残余航空梯度扰动见图 4.与图 3比较,我们发现,向下延拓过程放大了输入数据中的短波信号,突出了场源分布的细节特征.在恢复全球重力位模型的梯度扰动及剩余地形模型的重力梯度效应后,我们会对方案一、方案二得到的全张量梯度信号进行对比分析,相关结果见第2节.

图 4 研究区域Stokes积分方法解算的残余航空梯度扰动 Fig. 4 Residual airborne gradient disturbance calculated by Stokes numerical integration
2 区域重力梯度标定场结果讨论 2.1 残余航空重力扰动延拓误差统计

作为输入数据之一的残余航空重力扰动数据(原数据及其延拓值),其中解算方案一是用原残余航空重力扰动数据求解残余梯度扰动,解算方案二是用向下延拓至大地水准面上的残余重力扰动求解残余梯度扰动(见图 2);向下延拓会增强原数据中的短波信号,突出细节特征,把向下延拓值再向上延拓回原航空高度,与原残余航空重力扰动进行比较,可以为我们分析解算得到的不同方案梯度扰动的定性分析提供参考(本文未就延拓过程对不同方案解算梯度扰动结果的影响进行定量分析,这部分工作将作为后续研究内容来开展).原残余航空重力扰动及其延拓值(先向下延拓至大地水准面,再向上延拓回航空高度),见图 5,二者差值的统计量见表 10.

图 5 残余航空重力扰动与其先向下再向上延拓后的值. (a)原残余航空重力扰动;(b)延拓后的残余航空重力扰动 Fig. 5 Residual airborne gravity disturbance and its continuation values. (a) Original airborne gravity disturbance; (b) Continuation values from (a)
表 10 残余航空重力扰动与其延拓值的差值统计(单位: mGal) Table 10 Statistics of the differences between the residual airborne gravity disturbance and its continuation values (Units: mGal)

图 5表 10中可以看出:残余航空重力扰动经向下延拓至大地水准面,再向上延拓至航空高度后与原数据差值的标准差为1.0078 mGal.考虑边缘效应的影响,我们对原计算范围进行内缩,重新计算的差值统计量见表 11,内缩计算范围后的差值标准差减小至0.1269 mGal.

表 11 内缩计算范围后残余航空重力扰动与其延拓值的差值统计(单位: mGal) Table 11 Statistics of the differences between the residual airborne gravity disturbance and its continuation values within inner calculation area (Units: mGal)
2.2 方案一与方案二航空重力梯度张量结果差值统计

根据图 2的航空重力梯度场构建技术路线,经二维FFT变换得到残余梯度扰动,恢复EGM2008解算的航空高度梯度扰动和RTM解算得到的航空高度梯度效应,最终得到方案一的梯度全张量结果,见图 6.同理,由向下延拓至大地水准面上的重力扰动,经Stokes数值积分得到航空高度处的残余梯度扰动,再恢复EGM2008解算的航空高度梯度扰动和RTM解算得到的航空高度梯度效应,最终得到方案二的梯度全张量结果,见图 7.

图 6 方案一计算航空高度梯度扰动张量 Fig. 6 Airborne gradients disturbance calculated from computational scheme one
图 7 方案二计算航空高度梯度扰动张量 Fig. 7 Airborne gradients disturbance calculated from computational scheme two

对比图 6图 7,我们发现两种方案解算得到的梯度各分量具有相同的变化趋势特征,所不同的是,由延拓值经Stokes数值积分得到的方案二结果表现出了更为明显的细节特征,这与我们在分析残余重力扰动延拓得到的结论是一致的.

两种方案解算得到的研究区域航空高度梯度扰动张量的差值见表 12,从表中可以看出两种方案得到的研究区域重力梯度扰动各分量之差的最大标准差为6.4798E(Γyz分量),最小标准差为2.6968E(Γxy分量);这里的计算结果不符值主要来源于两个方面:(1)由原残余重力扰动经延拓过程引入的误差;(2)由不同方案梯度解算算法引入的误差.这两类误差的定量统计超出了本文的研究范围,将作为后续工作进行深入分析.

表 12 两种方案航空高度梯度扰动差值的统计量(单位:E) Table 12 Statistics of the differences of the airborne gradients disturbance between computational schemes (Units: E)

与残余重力扰动延拓对比分析一样,我们内缩了计算范围后发现(见表 13):内缩计算范围得到的差值标准差最大值为1.8307E(Γzz分量),最小值为0.7223E(Γyz分量).目前国内没有实测重力梯度数据进行对比分析,然而由本文解算结果得到的内符合精度已经达到了2E以内,小于部分主流航空重力梯度测量系统相关技术指标(见表 1),我们认为本文的解算结果是可靠的,这将为未来我国自主构建航空重力梯度标定场提供一定参考.

表 13 内缩计算范围后两种方案航空高度梯度扰动差值的统计量(单位:E) Table 13 Statistics of the differences of the airborne gradients disturbance derived from two computational schemes within inner calculation area (Units: E)
3 结论

(1) 本文通过联合全球重力位模型(EGM2008)、航空重力扰动数据和剩余地形模型(RTM)数据,基于频谱域(二维FFT变换)和空间域(Stokes数值积分)算法对毛乌素测区GT-2A航空重力测量系统采集的空中测线后处理重力扰动数据进行解算,构建了该地区的航空重力梯度扰动全张量.通过分析发现,全球重力位模型对区域梯度场构建的贡献较小,区域梯度信号主要由区域重力数据和地形数据来贡献.

(2) 基于残余航空重力扰动数据(原数据及延拓数据),通过不同方案解算得到的梯度扰动结果表明:两种方案得到的研究区域重力梯度扰动各分量之差(内符合)的最大标准差为6.4798E(Γyz分量),最小标准差为2.6968E(Γxy分量),内缩计算范围得到的差值标准差最大值为1.8307E(Γzz分量),最小值为0.7223E(Γyz分量).两种方案计算结果不符值主要来源于两个方面:残余航空重力扰动的延拓误差以及频谱域(二维FFT)和空间域(Stokes数值积分)算法的计算误差.对最终结果误差(以及数据本身误差和分辨率的影响)的进一步量化分析将作为我们后续工作进行深入研究.

(3) 本文得到的区域梯度场是由全球重力位模型、航空重力扰动数据及剩余地形模型数据蕴含的梯度效应三部分组合得到,采用的是大地测量学中经典的移去—计算—恢复的解算策略,然而该方法是否会存在不同源数据解算的信号频谱重叠及频谱空区等问题是值得研究的问题;另外地形高程数据的密度近似由常密度至变密度(根据先验地质、地震及其他地球物理数据)的转变也是我们今后的研究重点.

(4) 本文的研究思路和方法为我国自主构建重力梯度标定场及其他地球物理科学研究提供参考.

附录A

为了量化计算研究区域外地形质量对输出结果的影响,我们基于本文重力梯度场的解算高度(航空高度),选取了不同范围的RTM剩余高程模型解算了区域中心单点的全张量重力梯度值,数据及解算结果见附图 A1附表 A1.

附图A1 不同RTM数据范围 Appendix Fig. A1 Different RTM data ranges
附表A1 区域中心点不同RTM数据范围解算的重力梯度值(单位:E) Appendix Table A1 Gravity Gradients of the center point of the region calculated from different RTM ranges (Units: E)

附图 A1中,A1(a)A1(b)A1(c)A1(d)的地形范围分别为2 rad,1.5 rad,1.0 rad和0.5 rad,计算得到区域中心点重力梯度值见附表 A1.

附表 A1中可以看出,不同范围地形质量对本文航空高度处解算的梯度张量各分量影响较小,本文的解算结果输出范围相对于参与解算的地形质量范围是可靠的.

致谢  感谢两位匿名审稿专家及期刊编辑部对本文工作提出的意见和建议,这对本文工作质量的提高及今后研究思路的拓展大有裨益;文中相关图件由Generic Mapping Tools (GMT)及MALAB制图软件完成;在此一并表示感谢.
References
Annecchione M, Main B. 2007. Benefits of a high performance airborne gravity gradiometer (HD-AGGTM) for resource exploration. http://www.dmec.ca/ex07-dvd/E07/pdfs/64P.pdf.
Anstie J, Aravanis T, Johnston P, et al. 2010. Preparation for flight testing the VK1 gravity gradiometer.//Airborne Gravity 2010. Geoscience Australia Record, 23, 5-12.
Barthelmes F. 2013. Definition of functionals of the geopotential and their calculation from spherical harmonic models: theory and formulas used by the calculation service of the International Centre for Global Earth Models (ICGEM). http://icgem.gfz-potsdam.de/ICGEM/.
Bell R E. 1998. Gravity gradiometry. Scientific American, 278(6): 74-79.
Bian S F, Zhang C J. 1999. Computation of topographic effects on vertical gravity gradient. Computing Techniques for Geophysical and Geochemical Exploration (in Chinese), 21(2): 133-140.
Bian S F, Ji B. 2006. The development and application of the gravity gradiometer. Progress in Geophysics (in Chinese), 21(2): 660-664.
Blakely R J. 1996. Potential Theory in Gravity and Magnetic Applications. Cambridge: Cambridge University Press.
Bucha B, Janák J, Papčo J, et al. 2016. High-resolution regional gravity field modelling in a mountainous area from terrestrial gravity data. Geophysical Journal International, 207(2): 949-966. DOI:10.1093/gji/ggw311
DiFrancesco D, Meyer T, Christensen A, et al. 2009. Gravity gradiometry-today and tomorrow.//11th SAGA Biennial Technical Meeting and Exhibition. Swaziland, 80-83.
Dransfield M. 2007. Airborne gravity gradiometry in the search for mineral deposits.//Presented at the EGM 2007 International Workshop in Capri. Italy, 7: 341-354.
Dransfield M, Le Roux T, Burrows. 2010. Airborne gravimetry and gravity gradiometry at Fugro Airborne Surveys, The ASEG-PESA Airborne Gravity 2010 Workshop. Geoscience Australia Record, 23: 49-57.
Dransfield M H, Lee J B. 2004. The FALCON airborne gravity gradiometer survey system. The ASEG-PESA Airborne Gravity 2004 Workshop. Geoscience Australia Record, 18: 15-18.
Eshagh M, Sjöberg L E. 2009. Topographic and atmospheric effects on GOCE gradiometric data in a local north-oriented frame:a case study in fennoscandia and iran. Studia Geophysica et Geodaetica, 53(1): 61-80. DOI:10.1007/s11200-009-0004-z
Forsberg R. 1984. A study of terrain reductions, density anomalies and geophysical inversion methods in gravity field modelling. Columbus: Department of Geodetic Science and Surveying.
Hatch D. 2004. Evaluation of a full tensor gravity gradiometer for kimberlite exploration. The ASEG-PESA Airborne Gravity 2004 Workshop. Geoscience Australia Record, 18, 73-79.
Hinks D, McIntosh R, Lane R J L. 2004. A comparison of the Falcon and Air-FTG airborne gravity gradiometer systems at the Kokong Test Block, Botswana.//Lane R J L ed. Airborne Gravity 2004-Abstracts from The ASEG-PESA Airborne Gravity 2004 Workshop. Geoscience Australia Record, 18, 125-134.
Jarvis A, Reuter H I, Nelson A, et al. 2008. Hole-filled SRTM for the globe Version 4. CGIAR-CSI SRTM 90m Database http://srtm.csi.cgiar.org.
Jekeli C, Zhu L Z. 2006. Comparison of methods to model the gravitational gradients from topographic data bases. Geophysical Journal International, 166(3): 999-1014. DOI:10.1111/j.1365-246X.2006.03063.x
Jiang F Y, Gao L K, Huang L Y. 2011. Study on gravity gradient tensor of oil-gas model. Journal of Jilin University (Earth Science Edition) (in Chinese), 41(2): 545-551.
Jiang F Y, Huang Y, Yan K. 2012. Full gravity gradient tensors from vertical gravity by cosine transform. Applied Geophysics, 9(3): 247-260.
Jiang T, Xiao X N, Dang Y M, et al. 2018. Evaluation of the external accord accuracy of airborne gravity data with upward continuation. Acta Geodaetica et Cartographica Sinica (in Chinese), 47(5): 567-574.
Kieran A C, Hatch D, Main B. 2000. Performance of the Gedex High-Definition Airborne Gravity Gradiometer. The ASEG-PESA Airborne Gravity 2010 Workshop. Geoscience Australia Record, 23, 37-43.
Klees R, Koop R, Visser P, et al. 2000. Efficient gravity field recovery from GOCE gravity gradient observations. Journal of Geodesy, 74(7-8): 561-571. DOI:10.1007/s001900000118
Lee J B. 2001. Falcon gravity gradiometer technology. Exploration Geophysics, 32(3-4): 247-250. DOI:10.1071/EG01247
Liang X H. 2012. Method and experiment research on airborne gravimetry[Ph. D. thesis] (in Chinese). Beijing: University of Chinese Academy of Sciences.
Liu F M, Zhang Y F, Jing X, et al. 2013. Gravity gradient Parker's forward method and application Based on cosine transform. Acta Geodaetica et Cartographica Sinica (in Chinese), 42(2): 177-190.
Liu J Z, Liu L T, Liang X H, et al. 2013. Based on measured gravity anomaly and terrain data to determine the gravity gradients. Chinese Journal of Geophysics (in Chinese), 56(7): 2245-2256. DOI:10.6038/cjg20130712
Liu J Z, Liu L T, Liang X H, et al. 2015. Comparison of methods to model gravity gradient field using gravity anomaly data. Geomatics and Information Science of Wuhan University (in Chinese), 40(12): 1677-1682.
Liu J Z. 2017. High resolution regional full tensor gravity gradient field modeling based on SRBFs and RTM technique. Project Documentation of National Natural Science Foundation of China (in Chinese).
Liu X G. 2011. Theory and methods of the earth's gravity field model recovery from GOCE data. Acta Geodaetica et Cartographica Sinica (in Chinese), 41(2): 315-315.
Lumley J M, White J P, Barnes G, et al. 2004. A superconducting gravity gradiometer tool for exploration.//The ASEG-PESA Airborne Gravity 2004 Workshop. ASEG-PESA, 18: 21-39.
Luo S C. 2007. Rotating accelerometer gravity gradiometer[Master's thesis] (in Chinese). Wuhan: Huazhong University of Science and Technology.
Luo Z C, Wu Y L, Zhong B, et al. 2009. Pre-processing of the GOCE satellite gravity gradiometry data. Geomatics and Information Science of Wuhan University (in Chinese), 34(10): 1163-1167.
Makhloof A A, Ilk K H. 2008. Effects of topographic-isostatic masses on gravitational functionals at the Earth's surface and at airborne and satellite altitudes. Journal of Geodesy, 82(2): 93-111.
Mickus K L, Hinojosa J H. 2001. The complete gravity gradient tensor derived from the vertical component of gravity:a Fourier transform technique. Journal of Applied Geophysics, 46(3): 159-174.
Murphy C A. 2004. The Air-FTG airborne gravity gradiometer system.//The ASEG-PESA Airborne Gravity 2004 Workshop. ASEG-PESA, 18: 7-14.
Murphy C A. 2010. Recent developments with Air-FTG®.//Lane R J ed. The ASEG-PESA Airborne Gravity 2010 Workshop. ASEG-PESA, 23: 142-151.
Pail R, Bruinsma S, Migliaccio F, et al. 2011. First GOCE gravity field models derived by three different approaches. Journal of Geodesy, 85(11): 819-843. DOI:10.1007/s00190-011-0467-x
Pavlis N K, Factor J K, Holmes S A. 2007. Terrain-related gravimetric quantities computed for the next EGM.//Proceedings of the 1st International Symposium of the International Gravity Field Service (IGFS). Istanbul. 318-323.
Qian D, Liu F M, Li Y, et al. 2011. Comparison of gravity gradient reference map composition for navigation. Acta Geodaetica et Cartographica Sinica (in Chinese), 40(6): 736-744.
Rózsa S, Tóth G. 2005. Prediction of vertical gravity gradients using gravity and elevation data.//Sansò F ed. A Window on the Future of Geodesy. Berlin Heidelberg: Springer, 344-349.
Rózsa S, Tóth G. 2007. The determination of the effect of topographic masses on the second derivatives of gravity potential using various methods.//Tregoning P, Rizos C eds. Dynamic Planet. Berlin Heidelberg: Springer, 391-396.
Shu Q, Zhou J X, Y H. 2007. Present research situtation and development trend of airborne gravity gradiometer. Geophysical and Geochemical Exploration (in Chinese), 31(6): 485-488.
Van Tran K, Nguyen T N. 2020. A novel method for computing the vertical gradients of the potential field:application to downward continuation. Geophysical Journal International, 220(2): 1316-1329.
Wan J K, Luo Z C, Liu Z K. 2017. Evaluating the precision of GT-2A airborne gravimeter measurements. Journal of Geodesy and Geodynamics (in Chinese), 37(10): 1096-1100.
Wang B Z, Gao S S, Liu K H, et al. 2012. High-accuracy practical spline-based 3D and 2D integral transformations in potential-field geophysics. Geophysical Prospecting, 60(5): 1001-1016. DOI:10.1111/j.1365-2478.2011.01026.x
Wang Y M. 1990. Terrain effects on geoid undulation computation. Manuscripta Geodaetica, 15: 23-29.
Wu L, Ma J, Zhou Y, et al. 2009. Modelling full-tensor gravity gradient maps for gravity matching navigation. Journal of System Simulation (in Chinese), 21(22): 7037-7041.
Xiong S Q. 2009. The present situation and development of airborne gravity and magnetic survey techniques in China. Progress in Geophysics (in Chinese), 24(1): 113-117.
Ye Z R. 2010. The primary research of airborne gravity gradiometry[Ph. D. thesis](in Chinese). Beijing: University of Chinese Academy of Sciences.
Ye Z R, Liu L T, Liang X H. 2015. The signal extraction of gravity gradient disturbance on Moho. Acta Geodaeticaet Cartographica Sinica, 44(6): 609-615.
You W, Fan D M, He Q B. 2012. Recovering Earth's gravitational field model using GOCE satellite orbits. Geomatics and Information Science of Wuhan University (in Chinese), 37(3): 294-297.
Yuan Y. 2015. Comprehensive analysis, processing and interpretation of the full tensor gravity gradient data[Ph. D. thesis] (in Chinese). Changchun: Jilin University.
Zeng H L. 1999. Present state and revival of gravity gradiometry. Geophysical and Geochemical Exploration (in Chinese), 23(1): 1-6.
Zhang C D. 2005. On the subject of airborne gravimetry and airborne gravity gradiometry. Chinese Journal of Engineering Geophysics (in Chinese), 2(4): 282-291.
Zhang C J. 1999. Determination of vertical gradient of gravity anomaly with topographic data. Chinese Science Bulletin, 44(11): 1029-1034. DOI:10.1007/BF02886024
Zheng W, Xu H Z, Zhong M, et al. 2010. Study progress in international satellite gravity gradiometry programs. Science of Surveying and Mapping (in Chinese), 35(2): 57-61.
Zhu L, Jekeli C. 2006. Comparing methods to model the local gravity gradients from gravity anomalies.//Proceedings of the Symposium of the International Gravity Field Service. Istanbul, Turkey.
Zhu L Z, Jekeli C. 2007. Combining gravity and topographic data for local gradient modeling.//Tregoning P, Rizos C, eds. Dynamic Planet. Berlin, Heidelberg: Springer, 130: 288-295.
Zhu L Z, Jekeli C. 2009. Gravity gradient modeling using gravity and DEM. Journal of Geodesy, 83(6): 557-567.
边少锋, 张赤军. 1999. 地形起伏对重力垂直梯度影响的计算. 物探化探计算技术, 21(2): 133-140. DOI:10.3969/j.issn.1001-1749.1999.02.004
边少锋, 纪兵. 2006. 重力梯度仪的发展及其应用. 地球物理学进展, 21(2): 660-664. DOI:10.3969/j.issn.1004-2903.2006.02.050
蒋甫玉, 高丽坤, 黄麟云. 2011. 油气模型的重力梯度张量研究. 吉林大学学报(地球科学版), 41(2): 545-551.
蒋涛, 肖学年, 党亚民, 等. 2018. 航空重力向上延拓的外部精度评估. 测绘学报, 47(5): 567-574.
梁星辉. 2012.航空重力测量方法及试验研究[博士论文].北京: 中国科学院大学.
刘繁明, 张迎发, 荆心, 等. 2013. 基于余弦变换的重力梯度Parker正演方法及其应用. 测绘学报, 42(2): 177-190.
刘金钊, 柳林涛, 梁星辉, 等. 2013. 基于实测重力异常和地形数据确定重力梯度的研究. 地球物理学报, 56(7): 2245-2256. DOI:10.6038/cjg20130712
刘金钊, 柳林涛, 梁星辉, 等. 2015. 利用重力异常数据构建重力梯度场的方法比较. 武汉大学学报·信息科学版, 40(12): 1677-1682.
刘金钊. 2017.基于SRBFs和RTM技术的高分辨率区域全张量重力梯度场建模研究.国家自然科学基金青年基金申请书.
刘晓刚. 2011. GOCE卫星测量恢复地球重力场模型的理论与方法. 测绘学报, 41(2): 315-315.
罗嗣成. 2007.旋转加速度计重力梯度仪[硕士论文].武汉: 华中科技大学.
罗志才, 吴云龙, 钟波, 等. 2009. GOCE卫星重力梯度测量数据的预处理. 武汉大学学报·信息科学版, 34(10): 1163-1167.
钱东, 刘繁明, 李艳, 等. 2011. 导航用重力梯度基准图构建方法的比较研究. 测绘学报, 40(6): 736-744.
舒晴, 周坚鑫, 尹航. 2007. 航空重力梯度仪研究现状及发展趋势. 物探与化探, 31(6): 485-488. DOI:10.3969/j.issn.1000-8918.2007.06.002
宛家宽, 罗志才, 刘站科. 2017. GT-2A航空重力仪精度评定. 大地测量与地球动力学, 37(10): 1096-1100.
武凛, 马杰, 周瑶, 等. 2009. 重力场匹配导航的全张量重力梯度基准图模拟. 系统仿真学报, 21(22): 7037-7041.
熊盛青. 2009. 我国航空重磁勘探技术现状与发展趋势. 地球物理学进展, 24(1): 113-117.
叶周润. 2010. 航空重力梯度测量之初步研究. .
叶周润, 柳林涛, 梁星辉. 2015. Moho面扰动重力梯度信息的提取. 测绘学报, (06): 609-615.
游为, 范东明, 贺全兵. 2012. 利用GOCE卫星轨道反演地球重力场模型. 武汉大学学报·信息科学版, 37(3): 294-297.
袁园. 2015.全张量重力梯度数据的综合分析与处理解释[博士论文].长春: 吉林大学.
曾华霖. 1999. 重力梯度测量的现状及复兴. 物探与化探, 23(1): 1-6. DOI:10.3969/j.issn.1000-8918.1999.01.001
张昌达. 2005. 航空重力测量和航空重力梯度测量问题. 工程地球物理学报, 2(4): 282-291. DOI:10.3969/j.issn.1672-7940.2005.04.007
张赤军. 1999. 用地形数据确定重力异常垂直梯度. 科学通报, 44(6): 656-661. DOI:10.3321/j.issn:0023-074X.1999.06.019
郑伟, 许厚泽, 钟敏, 等. 2010. 国际卫星重力梯度测量计划研究进展. 测绘科学, 35(2): 57-61.