地球物理学报  2013, Vol. 56 Issue (9): 3134-3144   PDF    
应用全张量重力梯度数据优化密度和速度模型
于平 , 丁见焕 , 黄大年 , 马国庆 , 赵震宇 , 李丽丽 , 肖锋     
吉林大学地球探测科学与技术学院, 长春 130021
摘要: 重力梯度全张量数据具有高灵敏度反映探测目标局部微弱变化的特点, 可以对单一地震数据在建立速度模型过程中的不确定性, 尤其是各向异性速度模型建立过程中的精细结构进行补充修正.本文主要研究重震两类数据与各自所反映的物理模型关系以及密度模型与速度模型之间的相互联系和制约关系, 从而建立重力全张量数据与地震速度模型间的信息互补关系和模型参数间的解析表达式.在此基础上, 提出多参数加权正反演算法流程, 逐层异常分离技术和约束条件下迭代反演模型修正方法, 并实现模块编程和验证.实验数据和结果表明, 重力全张量数据和地震数据多参数融合能够明显提高地震速度建模精度, 减小模型解释的不确定性.
关键词: 重力梯度全张量      约束反演      密度模型      速度模型     
Optimization of density and velocity model utilizing full tensor gravity gradient data
YU Ping, DING Jian-Huan, HUANG Da-Nian, MA Guo-Qing, ZHAO Zhen-Yu, LI Li-Li, XIAO Feng     
College of Geo-Exploration Science and Technology, Jilin University, Changchun 130021, China
Abstract: Full tensor gravity gradient (FTG) data have relatively high sensitivity for representing small changing of exploration objects. It can compensate the drawback of velocity model created by seismic data, especially when anisotropic velocity model is used. In this study, gravity and seismic data are dealt with for obtaining the relationship between their causative sources as well as interactive constraints. The attempt can assist in establishing effective relation between the FTG data and the seismic velocity model. Due to the relation, a multiple parameter forward calculation and inversion algorithm and workflow are proposed by combination with effective weighting. An inversion method of layer-by-layer iterative stripping from shallow to deep is created in the study in order to improve density and velocity model. Programming and implementation has been carried out to test the feasibility of the method. Synthetic and field data test show that the high resolution velocity model can be obtained successfully by means of combination of FTG data and seismic presentation, which indicates that ambiguity of model interpretation can be greatly eliminated..
Key words: Full tensor gravity gradient      Constrained inversion      Density model      Velocity model     
1 引言

采用重力与地震数据融合解决速度建模问题已经进展多年,重力数据的精度以及由此反演出的密度模型的准确度成为数据融合的关键问题[1-8].由于传统重力场数据质量较差,制约了反演模型的质量和精度,限制了融合方法向解决复杂问题的深入推进.近十年来,西方发达国家将一些敏感性军事工程探测技术转向民用领域,包括高灵敏度的重力探测技术和数据,扩展了重力勘探技术和数据的应用范围,推进相关方法研究向高精度解决方案发展[9-10].在高灵敏度重力探测技术中,重力梯度全张量(FTG)数据,在能源、资源探测领域的成功应用改变了地震勘探领域研究人员对非震数据应用所持的保守看法[11-12].理论上,FTG数据由九个分量组成,其中5个分量为独立分量,分别记录了重力异常在三维空间的梯度微弱变化特性.较之传统重力场数据,FTG数据具有更高的分辨率,可以更加细致地反映场源的空间构造形态和密度属性变化规律[13].因此,对FTG数据进行多分量联合反演,可以充分提取张量数据的丰富信息,有效弥补地震各向异性速度建模过程中的一些主观假设和算法自身带来的不确定性[14].

近年来,由于新的探测技术带来的重力数据精度的提高,地震与非震勘探方法间融合技术不断改进,主要体现在对数据精细结构和模型局部变化之间的联系以及信息间的制约关系[15-19].其中,在研究密度体和速度体之间联系的过程中,尤其是突变边界条件局部变化区域,这些与梯度变化有关的高频端信息起到至关重要的作用.在重震数据融合研究方面,主要技术难点在于,基于周期变化和时间变化的地震数据如何与基于深度变化和频段变化的重力密度数据实现融合.在模型研究的基础上,数学表达形式的研究将体现彼此关系和特征、体现相关参数的参与和贡献加权、体现数学模型的联合表达形式和计算控制流程方式.在地震勘探领域长期研究过程中,基于地震波场建立各向异性速度模型的努力已经取得了卓有成效的进展,为相关非震数据信息应用对接提供了技术支持,为面向速度模型的校正分析而采取引进多参数、多属性、非震数据内容奠定了基础.本文通过利用高精度重力梯度全张量数据集特性[20-28],落实与地震数据相融合的具体技术,包括采用约束控制和局部聚焦反演方法技术、加权控制技术,对张量数据进行多分量联合反演,充分提取张量数据带来的丰富信息,改善速度建模过程中对局部精细结构的优化调整,达到完善地震速度建模的目标.

2 基本原理

通过重力与地震数据融合解决速度建模问题,首先建立FTG数据与场源属性和几何形态之间的参数关系,然后建立场源参数与地震波传播速度和走时关系,最终建立FTG数据与所修改的速度模型关系,并在相互约束条件下实现彼此间的模型修正.FTG数据与密度模型间的关系式,可以通过对重力位任意形体积分表达式进行求导展开获得.在笛卡尔坐标系下,针对三维层状密度模型体(图 1),重力场垂直分量的正演公式可以表达为傅氏变换后的波数域形式(简称:“层状介质Parker公式”[29]):

(1)

图 1 由上下起伏界面(hthb)定义的单层介质密度模型,密度(ρ)可沿水平方向变化,观测平面重力场和重力梯度可由正演计算和实际观测获得 Fig. 1 The density model defined by top interface ht and base interface hb, ρ variates along horizontal. The oberservation plane gravity field and gradient gravity can be derived from forward calculation and observation data

其中,Gz(xyz0)为在笛卡尔坐标系下(xyz0)观测点的垂向重力场z0>0,F[]为空间域向波数域转换的二维傅氏正变换算子;k分别为2维波数域中波数(或频率)和径向频率,其数值与数据网格化间距有关;C为引力参数;z0为观测平面高度,大于地球水准面高度;htn(xy)和hbn(xy),分别为层状介质中上下界面几何位置参数并作为第n阶求和项因子;ρ(xy)为该层的重力密度. (1)式主要反映波数域中重力场垂直分量与层状介质场源密度和分布形态关系.

对(1)式中F[Gz(xyz0)]项进行代换处理,可得到波数域重力位表达形式:

(2)

对(2)式中第一求和项整理,可得到相应的反演计算公式:

(3)

针对FTG数据集合形式,由(2)式可以整理出重力位沿三个方向的二阶空间导数响应并在此基础上形成矩阵表达式:

(4)

其中,等式左侧的表示波数域中9个重力梯度响应,仅5个响应为相对独立量,故表示为上三角矩阵形式,上三角阵中的第6个量满足位场拉普拉斯方程条件,即Gzz=-(Gxx+Gyy).

在波数域中,(4)式表示重力梯度张量数据集与层状介质密度和分布形态之间的解析关系,建立了针对模型体参数实现快速正演计算的基础[30-34](图 1).由于张量数据的空间特性,该矩阵形式既反映了分量与模型参数的对应关系,也反映了分量间和模型间的相互联系,以及通过矩阵表达式所反映的定量联系程度.这些关系和表达方式既实现正演计算,同时也构成了重力数据多分量和场源多参数反演分析的基础.从此基础出发,可以研究各分量对密度模型体正反演过程中所起的作用,从而确定各分量可利用的程度并设计加权因子,控制各分量对模型修正过程和贡献程度,达到有效利用高精度数据实现对密度和速度模型进行相互约束和修正目的,提高重震数据联合反演的效率和精度[35-36](图 2).

图 2 约束反演流程图表示FTG数据从不同角度精细反映地质目标分布,在局部约束条件下,通过FTG数据联合反演可以建立高精度的密度模型体 Fig. 2 FTG data set reflect geological object from different respects, so that high resolution density model can be built by joint inversion with the data set under constraint

通过对(1)-(4)式进行简单推导整理,最终可建立在波数域中进行张量数据集联合反演层状介质场源密度或几何参数的计算公式:

(5)

其中,F-1{}为二维离散傅氏逆变换算子并对应于上述公式(1)中正变换算子F[],wz为针对重力场垂直分量贡献的加权因子,wij(i=xyzj=xyz)为针对重力梯度响应贡献的加权因子,它们都是波数域中uv的函数,也可表示为wij(uv)和wz(uv). (5)式中的左侧表示平面分布场源的单层密度和厚度组合,右侧组合表示与FTG数据和迭代参数有关的参数组合,由此构成密度和厚度同时反演的迭代表达式.尽管利用(5)式不能同时确定密度和厚度参数,却为后续讨论中修正地震波传播的速度参数提供了一个可行的调整量.

在(5)式中,加权因子的确定可由维纳滤波原理实现.在确定滤波因子w(xy)过程中,希望输出Gs(xy)代表所反演的目标层异常,Gobs(xy)代表所有层的总异常.在范数2条件下,可要求希望输出和实际滤波输出,Gobs(xy)w(xy)满足优化条件:

目标函数最小化:

(6)

其中,〈〉是统计平均处理算符.(6)式表明,在波数域中,维纳滤波转换函数W(uv)满足维纳最小均方差原理[37]

(7)

其中,Gobs(uv)为空间域观测值Gobs(xy)的傅氏变换结果,Gobs*(uv)为Gobs(uv)的共轭函数;Gs(uv)为希望输出能谱函数.在位场接近对称分布的特殊情况下,(7)式可近似简化为:

(8)

其中,分别表示希望输出和总观测值的平均径向能谱.

因此,基于(7)和(8)式,可确定(5)式中针对反演目标层的非常数项加权因子,wz(uv)和wij(uv) (i=xyzj=xyz),它们类似于异常分解的滤波响应,将各自目标异常从总观测值中分解出来,实现由异常数据反演场源的计算原则.值得一提的是,Gs(uv)的计算可由控制断面中的目标层正演异常确定,即目标层正演异常的平均径向能谱.在此基础上,还可加入用于获取可靠反演结果的约束条件和参数加权因子.例如,地震和测井解释约束条件和相关的局部锁定控制条件、不同张量数据对反演目标贡献的加权控制、以及修正步长控制和对反演结果的可信度评估因子.

解决上述问题的关键是,建立部分地震解释数据约束条件下FTG数据联合反演过程,充分利用参数控制过程中的先验信息,实现相互间对未知量的迭代计算过程控制,使局部参数约束问题可以在解析表达式中得以体现并实现针对性调整.其基本思想和具体步骤是,通过两类参数融合矩阵和联合反演目标函数的建立,实现参数间的相互约束或信息补充.重震模型体多参数融合可表示为矩阵形式:

(9)

其中m1m2分别是地震速度模型和重力密度模型参数向量.在向量空间范数2条件下,针对压制白噪音干扰目标,可建立重震参数联合反演目标函数.利用梯度张量对模型参数在多维空间的挠动进行调整.

3 模型结构和参数计算

基于上述基本原理和计算流程,可对重力场和FTG数据进行联合反演,确定模型体密度和几何参数,并由此对速度模型进行参数融合.首先,通过传统分析手段,由地震数据确定初始速度模型,了解模型建立过程中的不确定性分布.然后,针对性引入FTG数据实施约束反演,由反演出的密度信息对速度模型中的不确定性部分进行修正.密度模型和速度模型参数间的转换关系的确定,可在上述公式基础上获得.对体数据的进一步的调整,还可以通过局部钻井数据、地震剖面数据和位场数据的统计关系微调密度和速度模型修正量,完善(9)式中所表示的相互制约形式.在对速度模型进行精细调整的实际应用中,不同的构造地质模型通常是对应不同的速度和密度模型转换关系,形成不同的参数转换组合和计算流程.

(1)盐丘盖层构造速度-密度分布模型(图 3),随深度增加,速度和密度模型参数间的关系可由下述递推过程确定,实现彼此间的约束和调整:

图 3 盐丘盖层构造速度-密度分布模型示意图,v1(xy)和ρ1(xy)分别表示第一层分布的横向不均匀的速度和密度参数,以此类推到其他层 Fig. 3 Velocity and density distribution theoretical model of salt-cover

其中,Mask1由密度模型反演结果F[ρ1(xy)(h1(xy)-h2(xy))]确定.同理,依次类推,可确定Maski,以及hiviTi(i=1,2…,5)等参数.

式中Ti为第i反射层的双程走时,Hi为第i反射层的深度,Maski为重力约束反演参数,可完成参数锁定和动态调整两项功能,并根据上述迭代反演过程确定贡献权限,实现目标函数极小化迭代过程.上述速度和密度参数间的递推关系由模型构造和地震波走时确定.重力反演过程的调整控制部分(如盐丘盖层起伏)应用了地震先验信息约束控制条件(如水平层控制),同时结合浅层到深层波段信息的递推关系(如利用(5)式反演时应用了wz(uv)和wij(uv) (i=xyzj=xyz)异常分离因子),最后还要对约束量(如Maski)进行系统性调整,由此形成相互约束条件下的修正参数的递推计算过程,以下给出计算原理和具体实现过程:考虑第一层时(即i=1),.通常情况下,上覆水平层时间T1(xy)可以通过自激自收双程走时确定,此时Maski为第一层密度ρ0(xy),同时对(5)式实施控制点深度约束迭代反演,从而有

(10)

同时注意到(10)式中的观测平面重力场和FTG异常分离因子wij(uv)(i=xyzj=xyz)由(8)式通过利用剖面或径向平均功率谱数据确定,希望输出为第一层的正演重力和FTG数据.

同理,针对第二层有,

等式右侧三项因子的组合由(5)式确定,相关异常分离因子wij(uv)(i=xyzj=xyz)由(8)式确定,希望输出为第二层的正演重力和梯度场,以此类推,用于第三和更深层的计算过程.

值得注意的是,上述递推公式中的T4为盐丘底面的双程走时,采用的展开和计算原理相同.但是可利用Mask控制参数调整实现部分地震数据(如位于盐丘边缘的缓起伏部分)约束位场反演.反演结果可对地震波能量传输偏弱的盐丘中心区实施修正.

(2)高角度推覆构造速度-密度分布模型(图 4),基于上述同样原理,随深度增加和角度变化,速度和密度模型参数间的关系可由下述递推过程完成:

图 4 高角度推覆构造速度密度理论分布模型示意图,v1(xy)和ρ1(xy)分别表示第一层分布的横行不均匀的速度和密度参数,以此类推到其他层 Fig. 4 Velocity and density distribution theoretical model of wide angle structure complexity

综上所述,图 2和两个速度-密度分布模型假设给出了模型体由浅入深逐层约束反演计算流程.具体做法是,由模型结构建立递推计算过程,通过利用(5)式和先验信息确定Maski值的约束区和反演区、从而固定约束区集中计算反演区目标层的密度和厚度组合参数,ρ(xy)(ht(xy)-hb(xy));通过利用(7)或(8)式和控制剖面,计算加权控制量wij(uv),实现针对反演对象的重力场和FTG异常分离.在计算过程中,还可以通过调整计算参数,例如,阻尼因子和迭代步长,实现(9)式中目标函数呈现递减趋势变化,进而控制速度模型局部参数调整量变化趋势,最终实现重震数据相互约束条件下的密度-速度模型反演计算.计算结果可为反射层标定提供一组密度层约束参数并计算出相关转换速度,从而达到对密度-速度模型的局部修正目的.

4 理论模型试算

上述理论和计算流程可通过编程形成计算模块和插件,插入到相关地震数据解释软件,形成了具有重震数据联合反演和融合解释功能的软件集成系统,并测试上述分析结果.首先采用盐丘理论模型进行测试(见图 5).其中假定密度参数为已知参数,重点完成对模型几何参数的修正.图 5a表示在钻孔和地震解释数据实现垂向和边部区域约束条件下,对重力场数据进行迭代反演,对中部隆起盐丘初始模型进行修正,即采用Mask掩膜对不同深度约束区和反演区进行划分和计算.计算流程包括正演计算待修正模型的响应并利用数据拟合差转换成对模型的改正量(见图 5b).依次完成迭代计算过程,并逐次完成对单层模型体的修正,然后至上而下完成对多层模型体的修正(见图 5c);计算结果还同时包括深度估计中可能出现的误差范围(见图 5c中的粉红线),有助于评估模型修正结果的可信度;图 5d给出了由待修正模型正演数据和理论数据拟合差确定的迭代计算收敛过程,反映了整个迭代过程向着既定目标平稳递减收敛,也说明了计算流程和实现方法的有效性.理论上,对重力场数据进行分析的方法同样适用于重力梯度数据,尤其在数据含噪音情况下,梯度数据由于灵敏度高应该表现更为出色.

图 5 井孔及边部先验信息约束条件下的盐丘理论模型反演结果 (a)密度体反演初始模型,上图表示针对模型体的重力异常响应,下图为密度模型;红线表示理论模型和异常,蓝线表示有待反演修正的初始模型和对应正演异常,数据拟合差为2.2mGal;(b)数据第一次迭代反演及模型修正结果,数据拟合差为1.7mGal;为确保收敛分别使用阻尼和松弛因子0.001和0.2;(c)数据第30次迭代反演及模型修正结果,数据拟合精度0.047,含深度计算误差范围(粉红色);(d)数据迭代反演及模型修正30次理论模型数据和模型正演数据拟合误差结果.S1,S2,S3表示地层,D1D2D3D4为密度. Fig. 5 Constraint inversion result using noise-free data set due to a salt dome controlled by two boreholes and the neighboring parts (a) Initial model and forward calculation (in blue) respect to synthetical model and data (in red) ResSTD value is 2.2(mGal); (b) Results at iteration 1, using damping and relaxation value of 0.001 and 0.2 respectively. ResSTD is 1.7 (mGal); (c) Inversion result of iteration 30 with ResSTD=0.047(mGal), with the depth estimation errorbars; (d) Plot of the data fit standard deviation of the 30 iteration steps.
5 实际数据试算

实际数据来源于美国地调局USGS-Openfile测试数据,主要解决盐丘隆起部分速度建模过程中由于能力传输造成的模糊区问题.初步验证表明,采用Mask约束盐丘缓起伏区,并对盐丘上涌造成的地震解释模糊区进行重力和重力梯度约束反演,效果理想达到了预期目标(见图 6).

图 6 FTG和地震数据约束反演及速度模型修正验证结果 各图上方为数据响应,下方为模型参数分布,蓝线为重力场,红线为水平梯度异常.(a)由地震SEGY数据建立初始密度模型,盐丘起伏部分成像模糊;(b)对重力和梯度数据进行联合反演得到的密度模型;(c)由单一地震数据确定的速度模型;(d)用密度数据进行局部修正后的速度模型.速度模型局部修正主要针对剖面图中部盐丘起伏部分进行. Fig. 6 Test results of constraint inversion of FTG and seismic data sets inconjunction with velocity model constraint and correction The upper map is FTG data, the blue line is gravity field and red line is horizontal gradient anomaly. The lower one are models of velocity and density. (a) The initial density model created with seismic data; (b) The density model obtained by the joint inversion with gravity and FTG data; (c) Shows the velocity model built with seismic data and (d) shows velocity module is modified by FTG data inversion for central part of velocity parameters.

试算结果表明,由于引入精度较高的FTG数据,从而获取对复杂盐丘地层的层间变化尤其是高角度变化反映明显的非震数据信息.该信息有助于对速度模型形成的后期阶段进行密度参数辅助性局部修正,可提高速度模型构建和调整的效率和精度.本方法可利用地震数据部分可靠信息,如近似水平层部分和盐丘上部和薄层边缘层,作为先验信息,对(密度体)反演目标实施反演约束,进而专注于反演盐丘起伏、厚层和高角度变化部分.该部分往往由于地震波能量传输在盐丘构造条件下衰减率高以及高角度分析难度大等问题,而造成地震解释中相对模糊部分.另外,利用地质构造和地震数据局部特点构建先验信息,形成局部约束条件下的重力梯度全张量数据反演,还能够有效压制地球物理数据反演过程中的非惟一性和不确定性影响因素.文中采用的Maski控制参量相当于平面控制掩膜,可在平面上根据不同层位先验信息分配约束区和反演区,并确定反演区中反演对象由密度和层位两类参数构成,而非传统反演中的单一的属性和几何因子.由Maski、不同模型结构和不同层位的双程走时形成对应于各层的递推计算序列,对模型体各层位由浅入深实施计算和调整,最终确定针对速度模型的局部调整方案.值得注意的是,相对重力场而言,FTG数据从不同方向和角度反映了密度模型空间变化精细特点,由此构成的联合反演结果除了可以用于对地震速度模型进行局部修正外,还具有减小多解性问题影响的特点.可商榷的理由是,观测到的重力场和FTG数据各量间相对独立,从不同角度反映了场源密度和形体参数变化更为精细的结构,其数据集联合反演,尤其在约束条件下反演,理论上可以极大改善单一数据反演带来的问题.

6 讨论与认识

计算表明,由于采用了本文给出的计算表达式,计算速度明显加快.例如,计算7个面积性网格分布的正演分量,即含6个梯度分量和重力场,不再是计算单一量的7倍关系.由表达式(4)可以推测,当完成一个重力分量的计算也就基本上完成了对其余6个量的计算.这对由(5)式构成的多批次迭代反演尤为重要,节省了大量计算时间.

就FTG联合反演而言,尽管采用了一系列约束条件和针对各层位异常分离加权滤波改进方法,使反演结果在非惟一性和计算过程稳定性方面得到了改善,减小了人为操作影响.但是,反演结果仍然无法彻底摆脱人为操作带来的随机影响传统问题.主要体现在控制参数的选择上,不同的迭代精度控制,阻尼因子和迭代步长,加权因子的确定以及实施约束的条件和范围,都会不同程度造成计算结果的变化.另外,对波数域零频响应处理,(7)和(8)式中维纳滤波转换函数处理以及网格化和扩边处理等均对反演结果有不同程度影响,属于反演计算经典问题,并得到了较好解决,相关文献对此做了介绍,故不在此文讨论中.

反演目标由该层的密度和形体参数组合构成ρ(xy)·(ht(xy)-hb(xy)),处理中可以利用整个模型体层与层间的递推关系加以确定.同时,还可以加入地震解释数据对该量进行事后约束,使所需要的密度参数变化专注于对速度调整时需要的部分.理想的解决方案应该是利用钻孔数据实施约束.另外,采用剖面控制下的维纳滤波方法确定异常分离响应,可获得与反演对象相对应的异常数据,对应的逼近程度取决于剖面结构模型对整体结构模型所具有的代表程度.通常情况下,由于实施了约束反演,所分离出的异常逼近程度对局部反演结果影响不大,即反演结果的相对起伏变化被限制在所约束的变换范围内.文中建立的模型相对简单,旨在表明重力梯度数据应用于速度模型校正的有效性.在更为复杂情况下,实施局部约束反演仍需做进一步研究,对此将另行讨论.

7 结论

本文利用FTG数据具有的高敏感度变化特性和相对独立性变化特点,利用多参数融合可以压制干扰、改善模型参数评估的基本原理,着重研究局部约束条件下重力场和FTG数据联合反演解释方法,研究在给定地质构造模型条件下,密度模型与速度模型之间的相互联系和制约关系,以及递推公式,较好解决了重力与地震数据融合及速度建模问题.提出并实现对地震速度模型局部不确定因素进行局部调整的方法和处理流程.同时完成了编程和插件模块,实现与其他综合分析方法的联合.初步试验效果表明,由于引入FTG数据,可以对速度模型形成后的模糊区继续进行密度参数辅助条件下的局部修正调整,高效率优化模型构建参数.另外,重力场和FTG数据的联合反演结果除了可以用于对地震速度模型进行局部修正外,还具有减小多解性影响问题.

本文主要针对可行性研究所需的理论基础进行研究,并初步完成模型数据验证和部分实测数据验证.在实际应用中,仍需要针对大量的实测数据和不同探测目标所形成的具体情况,开展更深入的实验研究,不断积累经验,推动多参数数据融合技术向更深层次发展.

致谢

感谢英国Bridgeporth地球物理公司首席科学家Mark Davies博士在本文研究过程中给予的有益讨论和建议.

参考文献
[1] Thomsen L. Weak elastic anisotropy. Geophysics , 1986, 51(10): 1954-1966. DOI:10.1190/1.1442051
[2] Alkhalifah T, Tsvankin I. Velocity analysis for transversely isotropic media. Geophysics , 1995, 60(5): 1550-1566. DOI:10.1190/1.1443888
[3] Liu Z Y. An analytical approach to migration velocity analysis. Geophysics , 1997, 62(4): 1238-1249. DOI:10.1190/1.1444225
[4] 刘礼农, 高红伟, 刘洪, 等. 三维VTI介质中波动方程深度偏移的最优分裂Fourier方法. 地球物理学报 , 2005, 48(2): 406–414. Liu L N, Gao H W, Liu H, et al. Wave equation depth migration with optimum split-step Fourier method in 3-D VTI media. Chinese J. Geophys. (in Chinese) , 2005, 48(2): 406-414.
[5] 郝重涛, 姚陈. 任意空间取向TI介质中体波速度特征. 地球物理学报 , 2007, 50(2): 546–555. Hao C T, Yao C. A theoretical study on the effect of ozone below cloud on total ozone retrieved by TOMS and a new inversion algorithm. Chinese J. Geophys. (in Chinese) , 2007, 50(2): 546-555.
[6] 郑海山, 张中杰. 横向各向同性(VTI)介质中非线性地震波场模拟. 地球物理学报 , 2005, 48(3): 660–671. Zheng H S, Zhang Z J. Synthetic seismograms of nonlinear seismic waves in anisotropic (VTI) media. Chinese J. Geophys. (in Chinese) , 2005, 48(3): 660-671.
[7] 陈洁, 温宁, 陈邦彦. 重磁电震联合反演研究进展. 地球物理学进展 , 2007, 22(5): 1427–1438. Chen J, Wen N, Chen B Y. Joint inversion of gravity-magnetic-electrical-seismic combination survey: progress and prospect. Progress in Geophysics (in Chinese) , 2007, 22(5): 1427-1438.
[8] 于鹏, 王家林, 吴健生, 等. 重力与地震资料的模拟退火约束联合反演. 地球物理学报 , 2007, 50(2): 529–538. Yu P, Wang J L, Wu J S, et al. Constrained joint inversion of gravity and seismic data using the simulated annealing algorithm. Chinese J. Geophys. (in Chinese) , 2007, 50(2): 529-538.
[9] Klingelé E E, Marson I, Kahle H G. Automatic interpretation of gravity gradiometric data in two dimensions: vertical gradient. Geophys. Prosp. , 1991, 39(3): 407-434. DOI:10.1111/gpr.1991.39.issue-3
[10] Marson I, Klingele E E. Advantages of using the vertical gradient of gravity for 3-D interpretation. Geophysics , 1993, 58(11): 1588-1595. DOI:10.1190/1.1443374
[11] Zhang C Y, Mushayandebvu M F, Reid A B, et al. Euler deconvolution of gravity tensor gradient data. Geophysics , 2000, 65(2): 512-520. DOI:10.1190/1.1444745
[12] Hood P. Gradient measurements in aeromagnetic surveying. Geophysics , 1965, 30(5): 891-802. DOI:10.1190/1.1439666
[13] Jekeli C. On optimal estimation of gravity from gravity gradients at aircraft altitude. Rev. Geophys. , 1985, 23(3): 301-311. DOI:10.1029/RG023i003p00301
[14] Vasco D W, Peterson J E, Majer E L. A simultaneous inversion of seismic traveltimes and amplitudes for velocity and attenuation. Geophysics , 1996, 61(6): 1738-1757. DOI:10.1190/1.1444091
[15] 杨辉, 戴世坤, 牟永光. 三维重力地震剥层联合反演. 石油地球物理勘探 , 2004, 39(4): 468–471. Yang H, Dai S K, Mu Y G. Joint layer-stripped inversion of 3-D gravity and seismic data. Oil Geophysical Prospecting (in Chinese) , 2004, 39(4): 468-471.
[16] 杜向东, 翁斌, 刘军荣, 等. TI介质偏移速度建模研究. 地球物理学报 , 2008, 51(2): 538–545. Du X D, Weng B, Liu J R, et al. Migration velocity modeling strategies of TI media. Chinese J. Geophys. (in Chinese) , 2008, 51(2): 538-545.
[17] 张秉铭, 滕吉文, 万志超. 横向各向同性介质中地震波速度分析及其意义. 地球物理学进展 , 1997, 12(1): 53–60. Zhang B M, Teng J W, Wan Z C. Velocity analysis and its implication in transversely isotropic media. Progress in Geophysics (in Chinese) , 1997, 12(1): 53-60.
[18] 冯锐, 陶裕录. 地震-重力联合反演中的非块状一致性模型. 地球物理学报 , 1993, 36(4): 463–475. Feng R, Tao Y L. A non-block consistent model in seismo-gravity inversion. Chinese J. Geophys. (in Chinese) , 1993, 36(4): 463-475.
[19] 王西文, 孟昭泰, 张工会, 等. 重力、地震联合反演在研究含油气盆地构造中的应用. 石油地球物理勘探 , 1991, 26(2): 201–209. Wang X W, Meng Z T, Zhang G H, et al. Application of combined gravitational and seismic inversion to studying geological structure inside a hydrocarbon-bearing basin. Oil Geophysical Prospecting (in Chinese) , 1991, 26(2): 201-209.
[20] 马国庆, 黄大年, 于平, 等. 改进的均衡滤波器在位场数据边界识别中的应用. 地球物理学报 , 2012, 55(12): 4288–4295. Ma G Q, Huang D N, Yu P, et al. Application of improved balancing filters to edge identification of potential field data. Chinese J. Geophys. (in Chinese) , 2012, 55(12): 4288-4295.
[21] Huang D N, Yu P, Feng X, et al. Gravity magnetic seismic data constraint inversion applied for all-in-one deep geological target interpretation. 2011 The international symposium on deep exploration into the Lithosphere. Beijing, 2011.
[22] Routh P S, Jorgensen G J, Kisabeth J L, et al. Inversion of gravity, tensor gravity, and magnetic field data significant for exploration. Oil & Gas Journal , 2002, 100(9): 40-43.
[23] Butler D K. Generalized gravity gradient analysis for 2D inversion. Geophysics , 1995, 60(4): 1018-1028. DOI:10.1190/1.1443830
[24] Routh P S, Jorgensen G J, Kisabeth J L, et al. Base of salt imaging using gravity and tensor gravity data. SEG annual meeting 2001. San Antonio, 2001: 1482-1484. http://www.oalib.com/references/18991035
[25] DiFrancesco D, Grierson A, Kaputa D, et al. Gravity Gradiometer systems-advances and challenges. // EGM 2007 International Workshop. Capri, Italy. http://www.oalib.com/references/18991036
[26] DiFrancesco D, Meyer T, Christensen A, et al. Gravity gradiometry-today and tomorrow. //11th SAGA Biennial Technical Meeting and Exhibition Swaziland, 2009: 80-83. http://www.oalib.com/references/18991037
[27] Murphy C A. Interpreting FTG Gravity data using horizontal tensor components. GEM 2007 International workshop, Capri, Italy, 2007. http://www.eageseg.org/data/egm2007/Sessione%20A/Oral%20papers/A_OP_04.pdf
[28] Murphy C A, Mumaw G R, Zuidweg K. Regional target prospecting in the Faroe-Shetland basin area using 3D-FTG gravity data. // 67th Meeting, EAGE, Expanded Abstracts, 2005: 503. http://www.oalib.com/references/18991068
[29] Blakely R J. Potential Theory in Gravity & Magnetic Applications. California: Cambridge University Press, 1995 .
[30] Murphy C A, Brewster J. Target delineation using full tensor gravity gradiometry data. ASEG 2007, Perth, Western Australia, Extended Abstracts, 2007. https://www.researchgate.net/profile/Colm_Murphy2/publication/242288632_Target_Delineation_using_Full_Tensor_Gravity_Gradiometry_data/links/545727c60cf2cf5164804344.pdf?inViewer=true&disableCoverPage=true&origin=publication_detail
[31] Mataragio J, Kieley J. Application of full tensor gradient invariants in detection of intrusion-hosted sulphide mineralization: implications for deposition mechanisms. First Break , 2009, 27(7): 95-98.
[32] Murphy C A, Dickinson J L. Geological mapping and targeting using invariant tensor analysis on full tensor gravity data. GEM 2010International Workshop, Capri, Italy, 2010. http://www.oalib.com/references/18991041
[33] Montana C J, Mickus K L, Peeples W J. Program to calculate the gravitational field and gravity gradient tensor resulting from a system of right rectangular prisms. Computers & Geosciences , 1992, 18(5): 587-602.
[34] Wu L, Tian J W. Automated gravity gradient tensor inversion for underwater object detection. Journal of Geophysics and Engineering , 2010, 7(4): 410-416. DOI:10.1088/1742-2132/7/4/008
[35] 王浩然, 陈超.重力梯度张量数据的反演方法研究. //中国地球物理学会第二十七届年会论文集.北京:中国地球物理学会, 2011: 702. Wang H R, Chen C. Study on inversion method of gravity gradient tensor data. // Chinese Geophysical Annual Meeting (in Chinese). Beijing: Chinese Geophysical Annual Meeting, 2011: 702. http://cpfd.cnki.com.cn/Article/CPFDTOTAL-ZGDW201110001638.htm
[36] 朱自强, 曾思红, 鲁光银. 重力张量数据的目标体边缘检测方法探讨. 石油地球物理勘探 , 2011, 3(3): 482–488. Zhu Z Q, Zeng S H, Lu G Y. Discussion on the edge detection technique for gravity gradient tensor data. Oil Geophysical Prospecting (in Chinese) , 2011, 3(3): 482-488.
[37] Pawlowski R S, Hansen R O. Gravity anomaly separation by Wiener filtering. Geophysics , 1990, 55(5): 539-548. DOI:10.1190/1.1442865