地球物理学报  2020, Vol. 63 Issue (1): 287-297   PDF    
黏声方程Q值反射波反演
邹鹏, 程玖兵     
同济大学海洋地质国家重点实验室, 上海 200092
摘要:地震波在非弹性介质中的衰减效应常用品质因子Q度量.相对准确的Q模型对提高强衰减介质中地震波成像的质量至关重要.本文提出了黏声介质反射波形反演(QRWI)方法来重建地下宏观Q模型.在缺乏大偏移距和低频地震数据时,该方法以黏声波方程为波场传播引擎,利用反射波核函数对模型中深部的敏感性去提取背景Q值.当速度高、低波数成分均已知时,基于波形拟合的QRWI可以获得较高分辨率的反演结果.由于地下介质速度的高波数扰动很难准确估计,本文通过引入峰值频移目标函数,极大地降低了QRWI对速度高波数成分的依赖.理论合成数据实验结果表明,本文方法反演得到的宏观Q模型可以满足衰减补偿逆时偏移成像的要求.
关键词: 衰减介质      品质因子      峰值频率      反射波形反演     
Visco-acoustic wave equation reflection inversion for the Q model
ZOU Peng, CHENG JiuBing     
State Key Laboratory of Marine Geology, Tongji University, Shanghai 200092, China
Abstract: Seismic wave attenuation in anelastic media is characterized by the quality factor (Q). An appropriate Q model is very important for improving the quality of seismic imaging in strong attenuation setting. In this paper, we develop a visco-acoustic wave equation reflection waveform inversion (QRWI) to update the background Q model. When the low frequencies and long-offset components of seismograms are unavailable, QRWI use the visco-acoustic wave equation as a wave propagator and exploit the reflection kernel to update the deep part of the Q model. If both the low-and high-wavenumber velocity components are known, QRWI can obtain relatively high resolution inverted model by matching the waveform directly. Because of the difficulty to accurately estimate the high-wavenumber velocity perturbations, we further introduce the misfit function of peak frequency shifts to reduce the dependency on the high-wavenumber components of the velocity model. The synthetic example demonstrates that the inverted Q model by QRWI can be confidently used for attenuation compensation in reverse time migration.
Keywords: Anelastic media    Quality factor    Peak frequency    Reflection waveform inversion    
0 引言

非弹性衰减广泛存在于真实地球介质中.地震尺度内的非弹性(本征)衰减主要与介质内部介观尺度(波长与颗粒尺度之间)的弹性非均匀性或流体饱和非均匀性诱导的波致流体流现象有关(Pride et al., 2003Müller et al., 2010).一方面,地震波在衰减介质中传播会耗散大量的能量,在地震记录上表现为振幅衰减和速度频散.此时,如果不采取针对性的处理会降低地震成像的质量,比如会出现衰减区域及其下方照明严重不足、分辨率低、成像位置偏差等(Hargreaves and Calvert, 1991Deng and McMechan, 2008Wang,2008刘伟和张剑锋,2018).另一方面,实验室观测和实际地震资料分析皆已表明,表征衰减性质的品质因子(Q)与岩石物性、流体类型及饱和度等因素有关(Winkler and Nur, 1982).可靠的Q模型可作为某些油气储层识别和烃类检测的一个判断标志(Frisillo and Stewart, 1980).因此,在油气勘探开发过程中,定量评估非弹性衰减效应,构建可靠本征衰减(Q)模型十分重要.

国内外学者对Q值估计问题做了大量研究.由于原理简单、计算效率高,射线层析Q反演方法(Brzostowski and McMechan, 1992Hu et al., 2011)被广泛应用于天然地震学(Chen and Clayton, 2012)与勘探地震学(Rossi et al., 2007)领域.除了速度和Q值会影响地震波振幅外,几何扩散、震源及检波器与地表的耦合特性以及震源的辐射模式均会影响地震波的振幅信息.在线性黏弹理论下,衰减系数与地震波的频率成正比关系(Kjartansson,19179),即高频成分的衰减比低频成分严重.为此,常以地震波的频率信息(尤其是质心频率)为数据拟合对象来区分衰减效应与弹性效应(Quan and Harris, 1997高静怀等,2008).然而,高频近似假设使得射线层析方法很难适应强横向非均匀介质.全波形反演(FWI)通过使模拟数据和观测数据达到最佳匹配来推测地下介质参数(Tarantola, 1984).从理论上讲,它可以考虑介质的诸多属性和各种波现象,同时具有层析和偏移的功能(Mora,1989),是目前利用地震数据推断地下介质性质最彻底的方法.近年来,由于计算能力和实际观测数据质量的提高,尤其是逆时偏移技术的成功应用,推动了FWI的迅猛发展(Virieux and Operto, 2009).近20年以来,考虑衰减效应的全波形反演(QFWI)方法也得到了发展,并逐渐走向实际应用(Song et al., 1995Liao and McMechan, 1996Smithyman et al., 2009Malinowski et al., 2011Zhou et al., 2018).最近,Yang等(2018)在贝叶斯框架下,用振幅随频率的线性变化关系反演品质因子,并将其用到了四川盆地实际资料中.

在黏声波形反演中,速度和衰减参数有很强的耦合性或trade-off效应(Song et al., 1995Hak and Mulder, 2011).Pratt等(2002)认为反演前期的振幅变化主要是由于速度误差引起的,在反演过程中应先构建比较准确的速度模型,然后估计Q模型.Kamei和Pratt(2008)用合成跨井数据比对了同时反演和顺序反演的结果,发现顺序反演是减弱参数耦合效应的有效途径.在前期处理已获得较好速度模型的前提下,Dutta和Schuster(2016)通过匹配峰值频率信息降低非衰减效应对反演的影响,从目标函数出发降低了QFWI的非线性程度.

在主动源观测条件下,FWI往往受观测孔径的制约.当缺少大孔径的回转波或折射波数据时,深部模型往往只被反射波照明.因此,利用反射波进行中深部的速度建模已成为共识.当反射数据占主导时,常规FWI目标函数的梯度(即核函数)中反射波的偏移响应与利用反射波反演背景模型的核函数相比,数值上往往高出一个数量级.此时,FWI主要体现偏移功能,侧重于对模型高波数成分的更新(Chi et al., 2015Alkhalifah and Wu, 2016).由于模型中低波数成分没有得到合理的重构,高波数的偏移成像结果也是不对的.为了解决反射波数据FWI的上述问题,Xu等(2012)提出了反射波形反演(RWI)方法,即将速度模型分解为背景模型和扰动模型两部分,利用反射波数据侧重更新背景速度.除了RWI这种数据域反演方法外,波动方程偏移速度分析(WEMVA)在成像域通过惩罚地下非零炮检距的能量去更新宏观速度模型(Symes,2008).最近,Shen等(2018)提出了一种在WEMVA框架下的Q模型反演方法.如何在RWI框架下达到类似目标尚无相关文献报道.

基于上述认识,假定在前处理已获得速度参数模型的前提下,本文聚焦基于反射波数据的Q值反演.首先在RWI框架下建立面向中深层背景Q模型的黏声反射波形反演(QRWI)方法.然后,为了降低对速度模型高波数扰动的依赖性,以峰值频移作为目标函数改进QRWI方法.最后用数值实验来验证方法的可行性,并给出结论.

1 波形匹配QRWI方法 1.1 RWI基本原理

常规FWI是通过使模拟数据和观测数据的波形达到最佳匹配来推断地下介质参数.在时间域,常以如下目标函数最小化为准则:

(1)

式中dobs是观测数据,dcal(m)是模拟数据,m代表介质的模型参数.在现有计算条件下,该优化问题一般采用梯度类局部寻优算法求解(Tarantola, 1984Virieux and Operto, 2009).根据伴随状态法(Plessix,2006),目标函数的梯度可用正传波场u与反传残差波场(伴随场)u零延迟互相关求取,即

(2)

其中⊗表示零延迟互相关运算符.在Born近似下,正传波场和伴随波场均可分解为背景场和散射场的叠加.于是方程(2)可以分解为(Chi et al., 2015):

(3)

FWI同时包含波动方程层析和偏移两种功能成分.方程(3)中第一项同时包含了初至波层析分量(透射波核函数)和反射波偏移响应分量,第二、三项分别是反射波震源端、检波器端层析分量(合称反射波核函数或反射波波路径),第四项则是界面处的高阶响应,一般不用于更新模型.图 1展示了单炮检对情况下FWI的核函数及其分量,可见不同类型的震相对目标函数梯度(核函数)的贡献和对模型更新的影响范围不一样,在反演的不同阶段,针对性地选取核函数中合适的成分是成功的关键(董良国等,2015).在有限偏移距和低频缺失的情况下,中深层的模型主要被反射波照明.然而在FWI梯度中,高波数的反射波偏移响应压制了低波数的反射波波路径,使得FWI很难用反射波来更新中深层背景参数.

图 1 两层介质单炮、单检波器FWI的核函数及其分量 (a) FWI核函数;(b)透射波核函数;(c)反射波偏移响应;(d)源端反射波核函数;(e)检波器端反射波核函数.为了便于观察,显示分解后的核函数时按照各自的能量进行了归一化处理. Fig. 1 The kernel and sub-kernels of FWI for a single shot and receiver pair (a) FWI kernel; (b) Diving wave kernel; (c) Reflection wave migration-ellipse; (d) Reflection wave kernel for source side; (e) Reflection wave kernel for receiver side. Note the energy is normalized for display.

RWI本质上是利用低波数的反射波波路径(图 1de),即公式(3)中第二、三项来更新中深层背景参数.实现上主要有两种方式:一种是在光滑的背景模型中引入参数扰动(图 2),用Born正演来获取上行反射波(Xu et al., 2012Ma and Hale, 2013);另一种是通过方向分解来区分上、下行波(Wang et al., 2016),然后挑选对应方向的波场计算反射波核函数.另外,单程波方程自然地考虑了波的传播方向并有较高的计算效率,也已被学者们用于RWI之中(Dong et al., 2018).本文基于Born正演来推导QRWI的梯度表达式.

图 2 RWI波场分解及波路径示意图 模型参数的高波数成分(mr)扰动背景场u00产生背向散射场δuδ,可以看作是“透射波”从震源出发沿波路径到达反射位置mr然后再返回接收点. Fig. 2 Schematic of RWI The high-wavenumber component mr perturbs wavefields u0 and 0 in the low-wavenumber background ms and contributes to back-scattering wavefields δu and δ, which can be considered "transmission" waves propagating along the wave paths between mr and the source and receiver.
1.2 QRWI方法

RWI方法需要将模型m分解为光滑的低波数背景(ms)和高波数扰动(mr),其中前者控制地震波传播的运动学信息,后者作为二次源产生散射(反射)波.QRWI需要考虑的模型参数包括速度(v)和品质因子(Q).在顺序反演的思想下,为了简化问题,先假设速度已知.为了便于推导,采用(非弛豫)弹性模量K(K=ρv2,密度ρ为常数)和弛豫时间τ作为模型参数.一般认为,Q对地震波的影响主要发生在传播路径上,即传播效应占主导(赵静等,2013),所以本文只关注背景Q模型的反演,不考虑高波数Q扰动的估计问题.为了用Born正演产生反射数据,将体积模量K分解为光滑背景K0和高波数扰动δK,即

(4)

选用Blanch等(1995)给出的含有一个力学弛豫模型的标准线性固体(SLS)模型来表征黏声介质中的衰减特性.在Born近似下,QRWI的正演(物理状态)方程可描述为

(5)

其中

(6)

以及w=[P, vx, vz, rp, δP, δvx, δvz, δrp]T为状态变量,其中Pvxvzrp分别表示声压、质点振动速度的xz分量和记忆变量的背景波场分量,以及它们对应的反射场为δP、δvxδvz分别表示对时间、xz方向偏导运算;s=[s, 0, 0, 0, 0, 0, 0, 0]T表示加载在声压分量上的点震源,T为转置算子.应力、应变的弛豫时间τσ, τε, τ与品质因子Q以及参考频率fω有关.在地震勘探应用中,通常选取震源子波主频为参考频率.于是上述弛豫时间满足:

(7)

以及

(8)

反射波形反演目标函数只涉及反射信号.根据伴随状态法,QRWI目标函数J(m)对弛豫参数τ的梯度表示为

(9)

式中是物理状态变量w的伴随状态变量,满足如下伴随方程:

(10)

其中Δd表示观测数据与模拟数据的反射波形残差,是波形残差最小二乘目标泛函的伴随源;算子L*是黏声介质波动方程正演算子L的伴随算子,具有如下表达式:

(11)

方程(9)中的梯度分别对应方程(3)中的第二、三、四项.如前文所说,第四项是界面处的高阶响应,不用于更新模型.计算出目标函数的梯度之后,可用最速下降法来迭代更新弛豫参数τ(x):

(12)

式中k是迭代次数,Pc(x)是预条件算子,更新步长α可通过任何反向追踪的线性搜索方法求解(Nocedal and Wright, 1999).在每一个迭代步内,当τ(x)更新之后用方程(8)将其转换为Q(x)模型.

2 峰值频率匹配QRWI方法 2.1 峰值频移目标函数

在实际地震数据处理中,可以借助层析成像或偏移速度分析构建宏观(低波数)速度模型,高波数的速度扰动通常是借助最小平方偏移或全波形反演获取.在黏声介质中,这两种实现除了需要低波数的背景速度外,还需要输入Q参数.在没有Q值输入的情况下,重建高波数的速度扰动是很有挑战的.高波数速度成分对地震波的振幅影响很大,这对基于波形匹配的QRWI方法应用带来了困扰.

地震衰减不仅影响振幅而且改造相位.由线性黏弹理论可知,衰减系数与地震波的频率成正比.高频成分的衰减强于低频成分的衰减导致地震波峰值频率往低频移动.几何扩散和透射损失等弹性效应虽影响地震波的振幅但不改变峰值频率信息.因此,衰减介质中峰值频移量主要受Q值的影响.Quan和Harris(1997)首次在射线Q层析中将峰值频率信息作为匹配对象.近期,Dutta和Schuster(2016)将峰值频率移动目标函数扩展到了波动方程层析中,在理论合成数据和实际数据中都获得了很好的效果.这种匹配地震数据次生属性(secondary attribute)或骨架数据(skeletonized data)的反演方法(Luo and Schuster, 1991)虽会牺牲反演的分辨率但会提供一个相对稳健的反演结果.本文将把峰值频移目标函数引入到QRWI中,以降低它对高波数速度成分的依赖.

QRWI的峰值频移目标函数定义为

(13)

式中分别是观测与模拟反射震相的峰值频率.在QRWI中,峰值频移与波形匹配目标函数的梯度有相同的表达式和计算流程,不同的是伴随源.此时伴随源可写成ΔfPobs,具体的推导过程可以参考Dutta和Schuster(2016)附录A.

2.2 峰值频率对速度高波数成分的依赖性

用Born正演来近似产生反射数据时,不仅需要光滑的背景速度模型,而且还需要扰动速度模型.在理论模型正演中,通常用高斯平滑来分解速度模型的高、低波数成分.本文将一维两层速度模型减去其高斯平滑模型得到的扰动模型定义为平滑差模型.在实际情况中,背景速度相对容易获取,速度扰动则难以准确获取,常用RTM图像作为“伪扰动速度”模型以保证反射数据走时信息的正确.用峰值频移作为目标函数时,要求地下介质与Q无关的弹性效应不改变地震子波形状,进而不改变峰值频率信息.下面将调查几种扰动速度模型对峰值频率的影响.

图 3a所示,对上层速度为2.1 km·s-1,下层速度为2.3 km·s-1的两层速度模型,用平滑半径为80 m的高斯平滑结果作为Born模拟的背景速度模型,分别用RTM像、界面处的单脉冲、界面上下极性相反的双值脉冲以及平滑差模型作为速度扰动,得到的模拟反射波见图 3b.相较于声波方程全波场模拟结果,无论是用RTM像还是用界面处的单脉冲作为扰动速度模型,其模拟反射波子波波形均有较大的改造.用界面上下极性相反的双值脉冲作为扰动速度模型产生的模拟结果较好地保持了子波的形状.用平滑差模型作为扰动速度模型,Born模拟产生的结果与全波场模拟结果最吻合.此实验说明Born模拟的扰动速度模型在反射界面上下应具有极性相反的结构,且最好从正负最大值连续变化到零.

图 3 Born正演的不同扰动速度模型(a)及其正演结果(b) (a)第1道红线为真实速度模型,黑线为背景速度模型;第2、3、4、5道分别为RTM像、界面位置单值脉冲、界面上下极性相反的双值脉冲以及平滑半径为80 m的平滑差模型;(b)红线为全波场模拟结果,其余是对应于(a)中各速度扰动的Born模拟结果. Fig. 3 Different perturbation velocity models (a) and results (b) of Born modeling In the first plot, the red curve of the first trace is true velocity model and the black curve is the smoothness background velocity model; the rest traces are the image of RTM, a pulse, a double value pulses with different signs, and the difference between the red curve of the first trace and the black curve, respectively. In the second plot, the first trace is the result from full wave equation modeling, and the rest traces corresponding to the Born modeling results of first plot.

进一步地,图 4a用平滑半径为160 m的平滑模型作为层状模型Born模拟的背景速度模型,分别用不同平滑尺度的平滑差模型作为扰动速度模型,其模拟结果如图 4b所示.从图中可以看出,Born正演的反射波形及其振幅谱对扰动速度的光滑性不敏感.此外,当下层速度从1.7 km·s-1变化到3.5 km·s-1时,用80 m平滑半径分解高、低波数速度成分的Born正演结果(图 5a)虽具有不同的波形振幅,数据极性也可能发生变化,它们都具有相同的峰值频率信息(图 5b).

图 4 不同平滑尺度模型分解的Born模拟结果 (a)速度模型,其中第1道红线为真实层状速度,黑线为平滑半径为160 m的光滑背景,第2、3、4、5道分别为平滑半径为40 m、80 m、120 m、160 m的平滑差模型;(b)地震响应,其中红线是黏声波方程模拟结果,其余对应于(a)图中各速度扰动模型的Born正演结果;(c)振幅谱. Fig. 4 Different scales decomposition of velocity for Born modeling In the plot (a), the red curve of the first trace is true velocity model and the black curve is the background velocity model with 160m smooth radius; the second, third, fourth, fifth traces denote the perturbation velocity model with 40 m, 80 m, 120 m, 160 m smooth radius, respectively. Plot (b) is the seismic response corresponding plot (a), and plot (c) is the corresponding amplitude spectrum.
图 5 两层模型不同反差情况下的Born正演上层速度固定为2.1 km·s-1,下层速度从1.7 km·s-1变化到3.5 km·s-1,其中分离高、低波数速度结构的高斯平滑半径为80 m; (a) Born正演模拟结果;(b)归一化振幅谱. Fig. 5 Two-layer model with 2.1 km·s-1 in the top layer and the velocity of the bottom layer varies from 1.7 km·s-1 to 3.5 km·s-1, the Gauss smooth radius is 80 m for decomposition the background and perturbation velocity for the Born modeling: (a) the Born modeling seismic traces; (b) the normalized amplitude spectra.

上述数值实验表明,在Born正演中,用平滑差模型代表高波数速度扰动可以很好地保持地震反射子波峰值频率信息,且基本不受平滑尺度、速度反差幅度的影响.因为背景Q反演不追求过高的分辨率,所以在偏移图像中拾取一些主要界面,置入平滑差模型,或在主界面内填入块状速度模型然后平滑做差,均可构建出合理的伪扰动速度模型.因此,基于峰值频移目标函数的QRWI算法可按如下四步进行:第一步用背景速度获取成像剖面并提取主要反射界面;第二步利用反射界面信息构建伪扰动速度模型;第三步用背景速度、伪扰动速度和初始Q模型实施Born正演,计算QRWI梯度与步长;第四步迭代更新Q模型直至目标函数收敛.

3 数值实验

本文用理论模型合成数据验证QRWI方法的有效性.如图 6,模型宽4 km,深2.5 km,中深部包含一个强的Q异常体,最小Q值为30.在其表面均匀布设20炮,炮间距为200 m;检波器对称分布于震源的两边,间距为10 m且最大偏移距为1 km.用峰值频率为30 Hz的雷克子波作为爆炸震源.在反演中用均匀模型(Q=200)作为初始模型.在勘探地震频段,对地震波有较强吸收作用的Q的变化范围为20~200(Quan and Harris, 1997Sams et al., 1997Pride et al., 2003),因此在反演过程中将Q值限制在这个范围内.当用波形匹配作为目标函数时,认为速度的高、低波数成分均已知.当峰值频率匹配作为目标函数时,仅已知速度的低波数成分,高波数速度模型用平滑差伪扰动速度模型代替.

图 6 理论模型波形匹配QRWI反演 (a)速度模型;(b)真实Q模型;(c)反演Q模型. Fig. 6 The inversion results of QRWI based on waveform difference misfit function (a) The velocity model; (b) The true Q model; (c) Inverted Q model.

在已知高、低波数速度模型的情况下,波形匹配QRWI可以很好地重建中深部Q模型(图 6).当高波数速度模型未知时,基于RTM像(图 7b)提取主要反射界面,然后置入平滑差模型构建伪扰动速度模型(图 7c),峰值频率匹配QRWI最终也较好地恢复了中深部的Q值(见图 7d图 8).由于反射数据在地表附近覆盖不够,QRWI在浅层只能得到一个等效解.要降低浅表层反演的多解性,需加入先验信息的约束或者多种反演方法相结合(如由QFWI先反演浅层Q模型).在本文的数值实验中,我们通过给定第一层的Q值来约束QRWI由于浅层数据覆盖不足产生的多解性.接下来,用QRTM来验证反演结果的可靠性.图 9分别展示了用初始Q模型、波形匹配反演Q模型和峰值频率匹配反演Q模型补偿的逆时偏移结果.初始模型QRTM成像结果在强衰减区域的能量损失并未得到有效补偿,且由于速度频散的影响使得成像位置有偏差.两种反演模型QRTM均很好地补偿了衰减区的能量损失,合理校正了速度频散效应从而实现反射波的正确归位.

图 7 理论模型峰值频率匹配QRWI反演 (a)光滑背景速度模型;(b) RTM图像(红线为拾取的反射界面);(c) Born正演采用的伪扰动速度模型;(d)反演Q模型. Fig. 7 The inversion results of QRWI based on frequency shift misfit function (a) The smooth background velocity model; (b) The image of RTM; (c) The pseudo high wavenumber velocity perturbation; (d) The inverted Q model.
图 8 伪井曲线对比 (a)x=1 km;(b) x=2 km;(c) x=3 km.红线是真实Q模型,蓝线是波形匹配反演结果,绿线是峰值频率匹配反演结果. Fig. 8 Pseudo-wells comparison at (a) x=1 km, (b) x=2 km and (c) x=3 km The red, blue, and green curves represent the true Q model, waveform difference based QRWI inversion result and peak frequency shift based QRWI inversion result, respectively.
图 9 黏声补偿逆时偏移(QRTM) (a)初始模型QRTM;(b)波形匹配QRWI反演模型QRTM;(c)峰值频率匹配QRWI反演模型QRTM. Fig. 9 Q-compensated RTM images using initial Q model (a) and waveform difference misfit based QRWI inverted Q model (b) and peak frequency shift based QRWI inverted Q model (c)

由于波形匹配目标函数拟合的是全部波形信息,其分辨率会高于峰值频移目标函数,但前提是要同时已知速度高、低波数成分.如图 10,从分辨率和目标泛函收敛情况来看,虽然峰值频率匹配不如波形匹配,但它不需要知道准确的高波数速度成分,这对QRWI走向实用至关重要.

图 10 目标函数收敛曲线(蓝线表示波形残差目标函数,红线表示峰值频移目标函数) Fig. 10 The history of convergence (the blue and red curves represent waveform difference and peak frequency shift misfit function, respectively)

最后通过在观测记录中加入高斯白噪来检测峰值频率匹配QRWI方法的抗噪性.图 11分别展示了不含噪声以及信噪比为16 dB、10 dB的第10炮合成地震记录.当信噪比为16 dB时,从观测记录中提取各震相的峰值频率信息相对容易,其反演结果(图 12a)比较好地刻画了背景Q值模型.当信噪比低至10 dB时,噪声掩盖了中深层的反射波能量,很难直接从观测记录中提取各震相的峰值频率信息.此时,对含噪数据做中值滤波处理后,峰值频率匹配QRWI反演结果(图 12b)也能大致刻画深部Q值异常体,但分辨率有所降低.

图 11 第10炮合成记录(a) 无噪声; (b) 含信噪比为16 dB高斯噪声; (c) 含信噪比为10 dB高斯噪声 Fig. 11 The 10th shot synthetic seismogram:(a) no noise, (b) with 16 dB Gauss noise, and (c) 10 dB Gauss noise
图 12 含高斯噪声数据峰值频率匹配QRWI反演结果 (a) 信噪比16 dB; (b) 信噪比10 dB Fig. 12 The peak frequency shift based QRWI inversion result for contaminated data, with 16 dB (a) and 10 dB (b) Gauss noise
4 结论

本文在RWI框架下,提出了分别基于反射波形匹配与峰值频率匹配的Q值反演(即QRWI)方法.在缺乏大偏移距和低频地震数据时,它可以利用反射波数据对模型中深部照明,重建相应区域的宏观Q模型.数值分析表明,在Born正演合成反射数据时,峰值频率远不如波形本身对速度模型的高波数(扰动)成分敏感.在反射界面处置入层状速度平滑差模型,构造出的伪扰动速度就可以保证Born正演反射波峰值频率信息的可靠性.相对于波形匹配QRWI,虽然峰值频率匹配QRWI方法反演分辨率略有降低,但它极大地降低了对速度模型高波数成分的依赖.理论模型数值实验表明,在满足各自条件的情况下,两种方法的反演结果均能为带衰减补偿的逆时偏移提供可靠的Q模型.数值实验表明,峰值频率匹配QRWI方法对高斯噪声有一定的抗噪性,反演结果的分辨率随着信噪比的减小而降低.

本文在顺序反演的思路下,假设宏观速度模型已知情况下反演Q参数.在获得相对准确的Q模型之后又能促进速度模型的优化与更新.在RWI框架下如何联合反演速度和衰减参数值得进一步研究.

References
Alkhalifah T, Wu Z D. 2016. The natural combination of full and image-based waveform inversion. Geophysical Prospecting, 64(1): 19-30. DOI:10.1111/1365-2478.12264
Blanch J O, Robertsson J O A, Symes W W. 1995. Modeling of a constant Q:methodology and algorithm for an efficient and optimally inexpensive viscoelastic technique. Geophysics, 60(1): 176-184. DOI:10.1190/1.1443744
Brzostowski M A, McMechan G A. 1992. 3-D tomographic imaging of near-surface seismic velocity and attenuation. Geophysics, 57(3): 396-403. DOI:10.1190/1.1443254
Chen T, Clayton R W. 2012. Structure of central and southern Mexico from velocity and attenuation tomography. J. Geophys. Res., 117(B9): B09302. DOI:10.1029/2012JB009233
Chi B X, Dong L G, Liu Y Z. 2015. Correlation-based reflection full-waveform inversion. Geophysics, 80(4): R189-R202. DOI:10.1190/geo2014-0345.1
Deng F, McMechan G A. 2008. Viscoelastic true-amplitude prestack reverse-time depth migration. Geophysics, 73(4): S143-S155. DOI:10.1190/1.2938083
Dong L G, Huang C, Chi B X, et al. 2015. Strategy and application of waveform inversion based on seismic data subset. Chinese J. Geophys. (in Chinese), 58(10): 3735-3745. DOI:10.6038/cjg20151024
Dong L G, Fan Z Y, Wang H Z, et al. 2018. Correlation-based reflection waveform inversion by one-way wave equations. Geophysical Prospecting, 66(8): 1503-1520. DOI:10.1111/1365-2478.12668
Dutta G, Schuster G T. 2016. Wave-equation Q tomography. Geophysics, 81(6): R471-R484. DOI:10.1190/geo2016-0081.1
Frisillo A L, Stewart T J. 1980. Effect of partial gas/brine saturation on ultrasonic absorption in sandstone. J. Geophys. Res., 85(B10): 5209-5211. DOI:10.1029/JB085iB10p05209
Gao J H, Yang S L, Wang D X. 2008. Quality factor extraction using instantaneous frequency at envelope peak of direct waves of VSP data. Chinese J. Geophys. (in Chinese), 51(3): 853-861.
Hak B, Mulder W A. 2011. Seismic attenuation imaging with causality. Geophysical Journal International, 184(1): 439-451. DOI:10.1111/j.1365-246X.2010.04848.x
Hargreaves N D, Calvert A J. 1991. Inverse Q filtering by Fourier transform. Geophysics, 56(4): 519-527. DOI:10.1190/1.1443067
Hu W Y, Liu J, Bear L, et al. 2011. A robust and accurate seismic attenuation tomography algorithm.//81st Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 2727-2731.
Kamei R, Pratt R G. 2008. Waveform tomography strategies for imaging attenuation structure with cross-hole data.//In 70th EAGE Conference and Exhibition incorporating SPE EUROPEC 2008. Rome, Italy.
Liao Q B, McMechan G A. 1996. Multifrequency visco-acoustic modeling and inversion. Geophysics, 61(5): 1371-1378. DOI:10.1190/1.1444060
Liu W, Zhang J F. 2018. De-absorption prestack time migration:steep dip structure imaging and its application. Chinese J. Geophys. (in Chinese), 61(2): 707-715. DOI:10.6038/cjg2018L0092
Luo Y, Schuster G T. 1991. Wave-equation traveltime inversion. Geophysics, 56(5): 645-653. DOI:10.1190/1.1443081
Ma Y, Hale D. 2013. Wave-equation reflection traveltime inversion with dynamic warping and full-waveform inversion. Geophysics, 78(6): R223-R233. DOI:10.1190/geo2013-0004.1
Malinowski M, Operto S, Ribodetti A. 2011. High-resolution seismic attenuation imaging from wide-aperture onshore data by visco-acoustic frequency-domain full-waveform inversion. Geophysical Journal International, 186(3): 1179-1204. DOI:10.1111/j.1365-246X.2011.05098.x
Mora P. 1989. Inversion=migration+tomography. Geophysics, 54(12): 1575-1586. DOI:10.1190/1.1442625
Müller T M, Gurevich B, Lebedev M. 2010. Seismic wave attenuation and dispersion resulting from wave-induced flow in porous rocks-A review. Geophysics, 75(5): 75A147-75A164. DOI:10.1190/1.3463417
Nocedal J, Wright S J. 1999. Numerical Optimization. New York: Springer.
Plessix R E. 2006. A review of the adjoint-state method for computing the gradient of a functional with geophysical applications. Geophysical Journal International, 167(2): 495-503. DOI:10.1111/j.1365-246X.2006.02978.x
Pratt R G, Hou F, Bauer K, et al. 2002. Waveform tomography images of velocity and inelastic attenuation from the Mallik 2002 Crosshole Seismic Surveys.//Proceedings of Scientific Results from the Mallik 2002 Gas Hydrate Production Research Well Program, Mackenzie Delta. North Territories, Canada: Geological Survey of Canada, 1-14.
Pride S R, Harris J M, Johnson D L, et al. 2003. Permeability dependence of seismic amplitudes. The Leading Edge, 22(6): 518-525. DOI:10.1190/1.1587671
Quan Y L, Harris J M. 1997. Seismic attenuation tomography using the frequency shift method. Geophysics, 62(3): 895-905. DOI:10.1190/1.1444197
Rossi G, Gei D, B hm G, et al. 2007. Attenuation tomography:an application to gas-hydrate and free-gas detection. Geophysical Prospecting, 55(5): 655-669. DOI:10.1111/j.1365-2478.2007.00646.x
Sams M S, Neep J P, Worthington M H, et al. 1997. The measurement of velocity dispersion and frequency-dependent intrinsic attenuation in sedimentary rocks. Geophysics, 62(5): 1456-1464. DOI:10.1190/1.1444249
Shen Y, Biondi B, Clapp R. 2018. Q-model building using one-way wave-equation migration Q analysis-Part 1:theory and synthetic test. Geophysics, 83(2): S93-S109. DOI:10.1190/geo2016-0658.1
Smithyman B, Pratt R G, Hayles J, et al. 2009. Detecting near-surface objects with seismic waveform tomography. Geophysics, 74(6): WCC119-WCC127. DOI:10.1190/1.3223313
Song Z M, Williamson P R, Pratt R G. 1995. Frequency-domain acoustic-wave modeling and inversion of crosshole data:Part II-Inversion method, synthetic experiments and real-data results. Geophysics, 60(3): 796-809. DOI:10.1190/1.1443818
Symes W W. 2008. Migration velocity analysis and waveform inversion. Geophysical Prospecting, 56(6): 765-790. DOI:10.1111/j.1365-2478.2008.00698.x
Tarantola A. 1984. Inversion of seismic reflection data in the acoustic approximation. Geophysics, 49(8): 1259-1266. DOI:10.1190/1.1441754
Virieux J, Operto S. 2009. An overview of full-waveform inversion in exploration geophysics. Geophysics, 74(6): WCC1-WCC26. DOI:10.1190/1.3238367
Wang F, Donno D, Chauris H, et al. 2016. Waveform inversion based on wavefield decomposition. Geophysics, 81(6): R457-R470. DOI:10.1190/geo2015-0340.1
Wang Y H. 2008. Seismic Inverse Q Filtering. Malden, MA: Blackwell Publishing.
Winkler K W, Nur A. 1982. Seismic attenuation:effects of pore fluids and frictionalsliding. Geophysics, 47(1): 1-15.
Xu S, Wang D L, Chen F, et al. 2012. Inversion on reflected seismic wave.//2012 SEG Annual Meeting. Las Vegas, Nevada: Society of Exploration Geophysicists, 1-7.
Yang X C, Teng L, Li J N, et al. 2018. Bayesian linearized amplitude-versus-frequency inversion for quality factor and its application. Journal of Geophysics and Engineering, 15(3): 751-758. DOI:10.1088/1742-2140/aa91f0
Zhao J, Gao J H, Wang D X, et al. 2013. Estimation of quality factor Q from pre-stack CMP records. Chinese J. Geophys. (in Chinese), 56(7): 2413-2428. DOI:10.6038/cjg20130727
Zhou W, Brossier R, Operto S, et al. 2018. Velocity model building by waveform inversion of early arrivals and reflections:a 2D case study with gas-cloud effects. Geophysics, 83(2): R141-R157. DOI:10.1190/geo2017-0282.1
董良国, 黄超, 迟本鑫, 等. 2015. 基于地震数据子集的波形反演思路、方法与应用. 地球物理学报, 58(10): 3735-3745. DOI:10.6038/cjg20151024
高静怀, 杨深林, 王大兴. 2008. 利用VSP资料直达波的包络峰值处瞬时频率提取介质品质因子. 地球物理学报, 51(3): 853-861. DOI:10.3321/j.issn:0001-5733.2008.03.026
刘伟, 张剑锋. 2018. 黏弹性叠前时间偏移:陡倾角构造成像与实际应用. 地球物理学报, 61(2): 707-715. DOI:10.6038/cjg2018L0092
赵静, 高静怀, 王大兴, 等. 2013. 利用叠前CMP资料估计介质品质因子. 地球物理学报, 56(7): 2413-2428. DOI:10.6038/cjg20130727