2016, Vol. 31
2. 中国国土资源航空物探遥感中心, 北京 100083
2. China Aero Geophysical Survey and Remote Sensing Center for Land and Resources, Beijing 100083, China
19世纪末，匈牙利物理学家Baron Roland von Eötvös发明扭称后，重力梯度测量进入勘探领域，并在早期的油气勘探中发挥了重要作用.20世纪90年代末，航空重力梯度系统投入商业地球物理勘探，随着仪器精度的提高和处理解释技术的逐步完善，航空重力梯度测量技术以其高效快速、灵活机动和高空间分辨率的优势，在固体矿产勘查和油气勘探中发挥着越来越重要的作用.目前商用的航空重力梯度测量系统均为Lockheed Martin公司基于Bell Aerosapce公司GGI(Gravity Gradient Instrument)技术制造，分为部分张量航空重力梯度系统(Falcon)(张昌达，2005；舒晴等，2007； Christensen，2013；Moore et al.，2014)和全张量航空重力梯度系统(Air-FTG、FTGeX)两类.近十年来，设计精度优于1E的超导航空重力梯度系统研制取得了重要进展，英国ARKeX公司研制的EGG(舒晴等，2007)、澳大利亚西奥大学与力拓公司联合研制的VK-1(Anstie et al.，2010)、加拿大GEDEX公司与斯坦福大学联合研制的HD-AGG超导航空重力梯度仪均已完成样机研制，近几年来一直处于试验飞行阶段，有望成为下一代航空重力梯度勘探主力产品.
重力梯度张量有9个分量，其中有5个独立分量，对于重力梯度测量结果地质解释而言，Uzz分量更能直观反映地质目标体空间分布特征，而多分量组合反演和解释也更能得到合理的结果.商用的航空重力梯度张量测量系统(Falcon)仅测量Uxy和Uxx-Uyy分量，而一直处于试飞阶段的超导航空重力梯度测量系统EGG、VK-1、HD-AGG也仅能测量单一张量分量或者张量分量组合.因此，基于部分张量航空重力梯度系统测量结果换算得到其他重力梯度张量分量是后续处理和解释的基础，其对丰富测量成果、提供多元化重力梯度张量信息和精化地质解释结果具有重要意义.
当前，国内外重力梯度张量理论方法研究主要集中在重力梯度数据处理(Yuan et al.，2013； D’ Andrea and Grujic，2013)、正演数值模拟(陈召曦等，2012；陈涛等，2015)、三维反演(Li，2001；王浩然等，2013；陈闫等，2014；杨娇娇，2015)、地质应用方面(Lane and Peljo，2004；Cevallos，2012；于平等，2013； Cevallos et al.，2014)和精化大地水准面及恢复地球重力场研究(孟嘉春和刘尚余，1993；吴星等，2011；蔡林等，2012；郑伟等，2014)等方面，而重力梯度张量分量转换处理研究很少涉及.本文从重力位频率域微分通用公式出发，建立了针对Falcon部分张量航空重力梯度系统的转换处理公式，并通过模型试验检验了算法的正确性和有效性.
1 重力梯度张量频率域转换公式重力梯度是重力位的二阶导数，写成张量形式为
|
(1) |
根据重力梯度张量定义，利用傅立叶变换的微分性质，在频率域(波数域)进行转换可将其简化.实际测量中，尤其对于航空测量而言，测量结果多数为平面或近平面数据，我们可以直接利用重力梯度分量的二维频谱之间的关系.根据傅立叶变换定义，重力位关于x-y平面的二维傅氏变换可表示为
|
(2) |
重力位对x、y求导的频率域表达式为
|
(3) |
重力位对z的方向导数求解较为复杂，详细推导过程可参考富里叶变换和位场谱分析方法及其应用(吴宣志等，1987)，这里直接给出结果，定义z轴正向朝下，当观测平面位于质量体之上时，可得：
|
(4) |
综合(3)式和(4)式可得重力位U(x，y，z)对x、y、z的频率域微分通用公式为
|
(5) |
根据(1)式和(5)式，可得常用重力梯度张量分量与重力位在频率域(波数域)存在下列关系
|
(6) |
对于类似于Falcon的部分张量航空重力梯度系统，假设某一高度平面z=z0测量结果为Txy=Uxy、TUV=Uxx-Uyy，通过二维傅立叶变换可得其二维频谱分别为$\tilde{T}$xy(kx，ky，z0)和$\tilde{T}$UV(kx，ky，z0)，将其带入(6)式，经过简单变换可得基于Txy(kx，ky，z0)的重力梯度张量频率域转换公式为
|
(7) |
同理，基于$\tilde{T}$UV(kx，ky，z0)的转换公式为
|
(8) |
通过(7)式或(8)式可得到其他重力梯度张量分量的二维频谱，再进行二维傅立叶反变换，便可根据部分张量航空重力梯度系统测量结果Txy或TUV换算得到重力梯度张量其他分量的空间域场值.对于其他类型的部分张量重力梯度系统，只要能测量得到平面上某一个分量或分量组合，就能按照相同的原理进行重力张量转换处理.
2 重力梯度转换处理算法设计重力梯度张量分量转换处理是基于平面数据二维傅立叶变换来实现的，对于空间域连续的二维函数 f(x，y)与其对应的傅立叶变换为F(kx，ky)，其正反变换表达式为
|
(9) |
由于实际测量获得的数据为有限长离散序列，需要采用离散算法，二维傅立叶变换的离散形式为
|
(10) |
根据(7)式、(8)式和(10)式，利用二维离散傅立叶变换可以进行重力梯度张量分量转换处理，图 1给出了算法的流程.
|
图 1 重力梯度张量分量转换处理算法流程 Figure 1 The flowchart of gravity gradient tensor frequency domain transformation |
为检验频率域重力梯度张量转换处理算法的正确性和有效性，我们以理论模型数据和理论模型数据加白噪声数据作为输入进行了转换处理试验.
3.1 理论模型数据试验设计的矩形棱柱体模型如图 2所示，模型参数：体积50 m×10 m×6 m，x范围为-25~25 m，y范围为-5～5 m，z范围为47～53 m，模型的中心坐标(0，0，50 m)，密度为1.5 g/cm3.
|
图 2 矩形棱柱体模型示意图 Figure 2 The right rectangular prism 3D model |
地面测网大小4000 m×4000 m，网格间距10 m×10 m，x、y范围-2000～2000 m，节点数共401×401.
模型重力梯度张量分量异常正演数值模拟结果如图 3所示.
|
图 3 模型重力梯度张量异常平面等值线图 Figure 3 Gravity gradient tensor anomaly modeling result from the prism model |
模型Uxy、Uxx-Uyy分量转换处理结果分别如图 4、5所示，与理论计算结果相比，转换得到的张量分量虽然在边部有一定畸变，但仍完整保留了原始异常的绝大部分信息，转换误差主要由重力梯度张量转换公式中的频率域奇点引起(舒晴等，2016).重力梯度张量转换处理误差统计结果见表 1，相对而言Uzz分量转换效果最好，异常峰—峰值转换结果与理论值相对误差小于2%，其余分量相对误差差别不大，均能控制在5%以内，取得了较好的效果，从理论上验证了算法的正确性和有效性.
|
图 4 利用Uxy分量重力梯度张量转换处理结果 Figure 4 The gravity gradient transformation result based on Uxy component |
|
图 5 利用Uxx-Uyy分量重力梯度张量转换处理结果 Figure 5 The gravity gradient transformation result based on Uxx-Uyycomponent |
|
表 1 重力梯度张量分量异常频率域转换处理误差统计表 Table 1 The error statistics of the gravity gradient transformation in frequency domain |
为了进一步验证重力梯度张量转换处理算法的实用性，对图 2所示模型的正演结果添加白噪声模拟实测数据，添加白噪声后Uxy和Uxx-Uyy数据信噪比为20 dB，如图 6所示，上方为理论值，下方为添加白噪声后的结果.
|
图 6 添加白噪声前后Uxy和Uxx-Uyy对比图 Figure 6 Comparison between theoretical data and synthetic data contaminated by white noise |
图 7、图 8分别给出了添加白噪声后Uxy和Uxx-Uyy数据重力梯度张量转换处理结果.可以看出，理论数据受白噪声污染后转换结果中各重力梯度张量分量的等值线出现了不同程度的畸变，原始数据种的噪声通过转换处理传递到了其他张量分量之中，但各分量异常主要特征仍清晰可见，验证了转换算法的有效性和稳定性.根据添加白噪声前后重力梯度张量转换处理模型试验结果对比可见：转换处理效果很大程度上取决于输入数据的精度，如果原始测量数据信噪比较高，或在进行转换处理之前进行相应的去噪处理，将大大改善转换处理效果.
|
图 7 利用添加白噪声的Uxy分量重力梯度张量转换处理结果 Figure 7 The gravity gradient transformation result based on Uxy synthetic data |
|
图 8 利用添加白噪声的Uxx-Uyy分量重力梯度张量转换处理结果 Figure 8 The gravity gradient transformation result based on Uxx-Uyy synthetic data |
基于重力位微分的频率域通用表达式，我们建立了针对Falcon部分张量航空重力梯度系统测量结果Uxy和Uxx-Uyy分量的重力梯度张量分量频率域转换公式，设计了重力梯度张量分量转换处理算法.理论模型试验结果表明重力梯度张量分量转换结果虽然在边部发生了一定程度的畸变，但完整再现了张量分量异常的主要特征，各分量异常峰-峰值误差控制在5%以内，取得了较好的效果，证明了算法的正确性.模拟实测噪声模型试验表明重力梯度张量分量转换处理过程中数据噪声会传递到转换后的分量，从而引起各分量异常等值线畸变，但异常主要特征仍然得以保留，验证了算法的稳定性和有效性.重力梯度张量分量频率域转换处理方法适用于于其他部分张量重力梯度测量系统测量结果的转换处理.
致谢 感谢审稿专家提出的修改意见和编辑部的大力支持！| [] | Anstie J, Aravanis T, Johnston P, et al. 2010. Preparation for flight testing the VK1 gravity gradiometer[C].//Airborne Gravity 2010 Abstracts from the ASEG-PESA Airborne Gravity 2010 Workshop. Australia:ASEG-PESA, 2010. |
| [] | Cai L, Zhou Z B, Zhu Z, et al .2012. Spectral analysis for recovering the Earth's gravity potential by satellite gravity gradients[J]. Chinese J. Geophys. (in Chinese), 55 (5) : 1565–1571. DOI:10.6038/j.issn.0001-5733.2012.05.014 |
| [] | Capriotti J, Li Y G, Krahenbuhl R. 2015. Joint inversion of gravity and gravity gradient data using a binary formulation[C].//International Workshop on Gravity, Electrical and Magnetic Methods and their applications. Chengdu, China:SEG, 326-329. |
| [] | Cevallos C .2014. Automatic generation of 3D geophysical models using curvatures derived from airborne gravity gradient data[J]. Geophysics, 79 (5) : G49–G58. DOI:10.1190/geo2013-0436.1 |
| [] | Cevallos C, Kovac P, Lowe S J .2013. Application of curvatures to airborne gravity gradient data in oil exploration[J]. Geophysics, 78 (4) : G81–G88. DOI:10.1190/geo2012-0315.1 |
| [] | Chen T, Zhang G B, Suo K, et al .2015. A study of different wavelet basis functions in forward modeling of gravity gradient anomaly[J]. Geophysical and Geochemical Exploration (in Chinese), 39 (S) : 91–97. |
| [] | Chen Y, Li T L, Fan C S, et al .2014. The 3D focusing inversion of full tensor gravity gradient data based on conjugate gradient[J]. Progress in Geophysics (in Chinese), 29 (3) : 1133–1142. DOI:10.6038/pg20140317 |
| [] | Chen Z X, Meng X H, Guo L H, et al .2012. Three-dimensional fast forward modeling and the inversion strategy for large scale gravity and gravimetry data based on GPU[J]. Chinese J. Geophys. (in Chinese), 55 (12) : 4069–4077. DOI:10.6038/j.issn.0001-5733.2012.12.019 |
| [] | Christensen A N. 2013. Results from FALCON® airborne gravity gradiometer surveys over the Kauring AGG Test site[C].//23rd Geophysical Conference. ASEG Expanded Abstracts, 1-4. |
| [] | D'Andrea L, Grujic M. 2013. The effects of density contrast surfaces on Airborne Gravity Gradiometry (AGG) data interpretation[C].//23rd Geophysical Conference. ASEG Expanded Abstracts, 1-4. |
| [] | Dransfield M H, Buckingham M J, Edwards C, et al .1991. Gravity gradiometry for geophysical prospecting[J]. Exploration Geophysics, 22 (1) : 107–110. DOI:10.1071/EG991107 |
| [] | Johnston P. 2012. Self-gradient effects for airborne gravity gradiometry[C].//22nd Geophysical Conference. ASEG Expanded Abstracts, 1-3. |
| [] | Lane R, Peljo M. 2004. Estimating the pre-mining gravity and gravity gradient response of the Broken Hill Ag-Pb-Zn Deposit[C].//17th Geophysical Conference. ASEG Expanded Abstracts, 1-5. |
| [] | Li X, Li Y G, Meng X H, et al. 2011. Practical issues in the processing and inversion of airborne gravity gradiometry data[C].//International Workshop on Gravity, Electrical and Magnetic Methods and Their Applications. Beijing, China:SEG, 33. |
| [] | Li Y G. 2001a. 3-D Inversion of gravity gradiometer data[C].//2001 SEG Annual Meeting. SEG Technical Program Expanded Abstracts, 1470-1473. |
| [] | Li Y G. 2001b. Processing gravity gradiometer data using an equivalent source technique[C].//2001 SEG Annual Meeting. SEG Technical Program Expanded Abstracts, 1466-1469. |
| [] | Li Y G. 2015. Understanding curvatures in gravity gradiometry and their application to source estimation[C].//2015 SEG Annual Meeting. SEG Technical Program Expanded Abstracts, 1482-1489. |
| [] | 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[J]. Chinese J. Geophys. (in Chinese), 56 (7) : 2245–2256. DOI:10.6038/cjg20130712 |
| [] | Lumley J M, White J P, Barnes G, et al. 2004. A superconducting gravity gradiometer tool for exploration[C].//Airborne Gravity 2004 Geoscience Australia Record. Australia:ASEG-PESA, 2004. |
| [] | Meng J C, Liu S Y .1993. Research on accuracies for satellite gravity gradiometry and Earth's gravitational field[J]. Acta Geophysica Sinica (in Chinese), 36 (6) : 725–739. |
| [] | Moore D, Dransfield M, Christensen A, et al. 2014. High-resolution FALCON® airborne gravity gradient data adds values to seismic exploration in the Yakka Munga area of Canning Basin, West Australia[C].//Beijing 2014 International Geophysical Conference & Exposition. Beijing, China:SEG, 1273-1276. |
| [] | Shu Q, Zhou J X, Yin H .2007. Present research situtation and development trend of airborne gravity gradiometer[J]. Geophysical & Geochemical Exploration (in Chinese), 31 (6) : 485–488. |
| [] | Shu Q, Zhu X Y, Zhou J X, et al .2015. Forward modeling expressions for right rectangular prism with constant density[J]. Geophysical & Geochemical Exploration (in Chinese), 39 (6) : 1217–1222. |
| [] | Shu Q, Zhang W Z, Zhou J X, et al .2016. Error Study on gravity gradient tensor components transformation[J]. Geophysical & Geochemical Exploration (in Chinese), 40 (1) : 111–116. |
| [] | Wang H R, Chen C, Du J S .2013. 3-D inversion of gravity gradient tensor data and its applications[J]. Oil Geophysical Prospecting (in Chinese), 48 (3) : 474–481. |
| [] | Wetherley S, Moore D. 2015. Using FALCON® airborne gravity gradiometry for oil and gas exploration recent case studies[C].//Proceedings of the 12th SEGJ International Symposium. Tokyo, Japan:SEG, 60-63. http://cn.bing.com/academic/profile?id=b20f318d22895ea20784d02e8df20a8b&encoded=0&v=paper_preview&mkt=zh-cn |
| [] | Wu X, Wang K, Feng W, et al .2011. Method of tensor invariant based on non-full tensor satellite gravity gradients[J]. Chinese J. Geophys. (in Chinese), 54 (4) : 966–976. DOI:10.3969/j.issn.0001-5733.2011.04.011 |
| [] | Wu X Z, Liu G H, Xue G Q, et al .1987. Fourier transform and spectrum analysis method for potential field (in Chinese)[M]. Beijing: Surveying and Mapping Press . |
| [] | Yu P, Ding J H, Huang D N, et al .2013. Optimization of density and velocity model utilizing full tensor gravity gradient data[J]. Chinese J. Geophys. (in Chinese), 56 (9) : 3134–3144. DOI:10.6038/cjg20130923 |
| [] | Yuan Y, Huang D N, Yu Q L, et al .2013. Noise filtering of full-gravity gradient tensor data[J]. Applied Geophysics, 10 (3) : 241–250. DOI:10.1007/s11770-013-0391-3 |
| [] | Zeng H L .2005. Gravity field and gravity exploration (in Chinese)[M]. Beijing: Geological Publishing House . |
| [] | Zhang C D .2005. Several new types of airborne gravimetric systems and airborne gravity gradiometric systems[J]. Geophysical & Geochemical Exploration (in Chinese), 29 (6) : 471–476. |
| [] | Zheng W, Xu H Z, Zhong M, et al .2014. Precise recovery of the Earth's gravitational field by GRACE Follow-On satellite gravity gradiometry method[J]. Chinese J. Geophys. (in Chinese), 57 (5) : 1415–1423. DOI:10.6038/cjg20140506 |
| [] | 蔡林, 周泽兵, 祝竺, 等.2012. 卫星重力梯度恢复地球重力场的频谱分析[J]. 地球物理学报, 55 (5) : 1565–1571. DOI:10.6038/j.issn.0001-5733.2012.05.014 |
| [] | 陈涛, 张贵宾, 索奎, 等.2015. 不同小波基函数在重力梯度异常正演计算中的应用研究[J]. 物探与化探, 39 (S) : 91–97. |
| [] | 陈闫, 李桐林, 范翠松, 等.2014. 重力梯度全张量数据三维共轭梯度聚焦反演[J]. 地球物理学进展, 29 (3) : 1133–1142. DOI:10.6038/pg20140317 |
| [] | 陈召曦, 孟小红, 郭良辉, 等.2012. 基于GPU并行的重力、重力梯度三维正演快速计算及反演策略[J]. 地球物理学报, 55 (12) : 4069–4077. DOI:10.6038/j.issn.0001-5733.2012.12.019 |
| [] | 刘金钊, 柳林涛, 梁星辉, 等.2013. 基于实测重力异常和地形数据确定重力梯度的研究[J]. 地球物理学报, 56 (7) : 2245–2256. DOI:10.6038/cjg20130712 |
| [] | 孟嘉春, 刘尚余.1993. 卫星重力梯度测量与地球引力场的精度研究[J]. 地球物理学报, 36 (6) : 725–739. |
| [] | 舒晴, 周坚鑫, 尹航.2007. 航空重力梯度仪研究现状及发展趋势[J]. 物探与化探, 31 (6) : 485–488. |
| [] | 舒晴, 朱晓颖, 周坚鑫, 等.2015. 矩形棱柱体重力梯度张量异常正演计算公式[J]. 物探与化探, 39 (6) : 1217–1222. |
| [] | 舒晴, 张文志, 周坚鑫, 等.2016. 重力梯度张量分量转换处理误差研究[J]. 物探与化探, 40 (1) : 111–116. |
| [] | 王浩然, 陈超, 杜劲松.2013. 重力梯度张量数据的三维反演方法与应用[J]. 石油地球物理勘探, 48 (3) : 474–481. |
| [] | 吴星, 王凯, 冯炜, 等.2011. 基于非全张量卫星重力梯度数据的张量不变量法[J]. 地球物理学报, 54 (4) : 966–976. DOI:10.3969/j.issn.0001-5733.2011.04.011 |
| [] | 吴宣志, 刘光海, 薛光琦, 等.1987. 富里叶变换和位场谱分析方法及其应用[M]. 北京: 地质出版社 . |
| [] | 于平, 丁见焕, 黄大年, 等.2013. 应用全张量重力梯度数据优化密度和速度模型[J]. 地球物理学报, 56 (9) : 3134–3144. DOI:10.6038/cjg20130923 |
| [] | 曾华霖.2005. 重力场与重力勘探[M]. 北京: 地质出版社 . |
| [] | 张昌达.2005. 几种新型的航空重力测量系统和航空重力梯度测量系统[J]. 物探与化探, 29 (6) : 471–476. |
| [] | 郑伟, 许厚泽, 钟敏, 等.2014. 基于GRACE Follow-On卫星重力梯度法精确反演地球重力场[J]. 地球物理学报, 57 (5) : 1415–1423. DOI:10.6038/cjg20140506 |