带正则化校正的TTI介质qP波方程及其逆时偏移方法和应用

刘文卿，王西文，王宇超，雍学善，王小卫，张涛

1 成都理工大学油气藏地质及开发工程国家重点实验室，成都　610059　2　中国石油勘探开发研究院西北分院，兰州　730020



带正则化校正的TTI介质qP波方程及其逆时偏移方法和应用

刘文卿1,2，王西文2，王宇超2，雍学善2，王小卫2，张涛2

1 成都理工大学油气藏地质及开发工程国家重点实验室，成都6100592中国石油勘探开发研究院西北分院，兰州730020

摘要各向异性研究对地下介质精确成像有着重要的意义,在当前计算机硬件迅速发展及宽方位地震数据采集日益普遍的情况下，成像必须考虑介质的各向异性.逆时偏移是基于双程波动方程的较为精确的数值解的成像方法，所以相对于其他地震成像方法，它具有很大的优势,譬如不受反射界面的倾角限制、偏移速度结构合适时能够使回转波及多次波正确成像.在各向同性介质中，可使用标量波方程来模拟波场.而在各向异性介质中，P波和SV波是相互耦合的，即不存在单纯的标量波传播，通常利用能代表耦合波场中P波分量运动学特征的拟声波(qP波)进行偏移成像.本文中，我们推导出了TTI介质下qP波控制方程.该方程可采用显式有限差分格式进行求解.通过声学近似，若沿对称轴方向的剪切波速度为零，对于对称轴方向不变且ε≥δ的模型来说，可得到稳定的数值解.但对于TTI介质来说，由于沿对称轴方向各向异性参数是变化的，声学近似会引起波场传播及数值计算的不稳定.因此，我们提出了正则化有限横波的方法，很好地解决了这一问题.最后，给出了Foothill模型的测试结果及某探区实际资料试算结果，展示了采用这个方程进行复杂TTI模型正演和高质量逆时偏移成像结果，证实了该方法的正确性和实际资料应用中的有效性.

关键词TTI介质； 各向异性； 波场耦合； 有限横波； 正则化； 逆时偏移

1引言

随着宽方位地震勘探技术的发展及勘探精度的不断提高，各向异性地震成像受到业界更多的关注.逆时偏移方法基于精确的双程波动理论，在时间方向上对震源波场进行正向延拓、检波点波场反向延拓(Baysal et al.，1983；McMechan，1983；Whitmore，1983)，正反向延拓的波场零相位相关产生成像结果.该方法具有无反射倾角限制、能适应强横向变速,能够对回转波及多次反射波进行正确成像的特点.在各向同性介质中，可以对标量声波方程进行求解，得到逆时偏移成像结果.严格来说实际接收的地震波场更接近于弹性波方程所描述的波场特征.因此，一些学者考虑能否直接用弹性波方程外推成像(Sun et al.，1986；Chang and McMechan，1986).但利用弹性波方程进行波场外推成像存在很多困难，首先是多参数建模困难，其次是直接外推耦合的方程组计算量巨大，另外多分量数据波场分解及相关成像难度很大.因此，目前开展弹性波逆时偏移成像实用化的条件还不成熟.但并不能否认弹性波逆时偏移成像的重要性.

对于各向异性介质，严格意义上的标量波是不存在的.但对于一般的偏移成像，与各向同性介质相似，可通过声学近似的方式，能够找到一个能代表TI介质中P波运动学特征的控制方程.利用这样的方程就能进行P波波场的外推和偏移成像.Alkhalifah(1998)首先提出了拟声波近似思想，通过声学近似简化耦合的频散关系，引入了TI介质拟声波近似频散关系并导出方便数值求解的qP波控制方程.尽管这个标量波场的频散关系具有和矢量波场中P波相似的运动学特征，但它会产生伪横波.对于VTI介质，基于Alkhalifah的频散关系式，许多学者已经推导出了拟声波双程波动方程并用于VTI介质逆时深度偏移成像中，Alkhalifah(2000)，Zhou等(2006)，Du等(2008)，Fletcher等(2008)，Duveneck等(2008)利用不同辅助变量简化了拟声波方程的求解.同时在二维情况下通过解耦的P-SV关系导出‘纯P’波方程，但‘纯P波’方程的快速数值算法仍不成熟(Liu et al.，2009；Chu et al.,2011).

从Hooke定律和运动方程出发，通过声学近似也可推导出VTI介质拟声波方程.从平面波传播的角度，VTI介质向TTI介质转化可以看成在旋转后的坐标系下的平面波传播.坐标系旋转后，仅增加了数学运算的复杂度，并不增加物理量纲.通过设置对称轴方向SV波速度为零，可推导出TTI介质下的拟声波方程(Fletcher et al.，2009；Duveneck and Bakker,2011).

本文将着重介绍一个新的TTI介质拟声波方程推导及实现策略，该方程从精确P-SV波频散关系推导而来，并且使用Thomsen(1986)定义的各向异性参数ε和δ进行参数化.对于对称轴方向不变且ε≥δ的介质(VTI或HTI介质)来说，通过声学近似，可得到稳定的数值解，可以直接将该方程应用于正演和偏移中去.然而，对于非均匀TTI介质中的正演和逆时偏移，各向异性对称轴方向的变化会引起数值计算的不稳定问题.为了得到更加稳定的TTI介质中qP波方程，Flecter等(2009)通过在声学近似方程中加入非零横波速度得到有限横波方程来增加方程求解稳定性.更自然的解决办法由Duveneck和Bakker(2011)首先提出，即对弹性波方程直接做声学近似，避开从频散关系推导qP波方程时对称轴方向缓变的假设，导出了TTI介质中与弹性波方程对应的qP波方程，该方程可以稳定模拟P波分量的传播.Zhang Y和Zhang H Z(2009)从构造自共轭微分算子角度导出了一种TTI介质中的qP波方程.该方程可以看作是Duveneck和Bakker(2011)导出的qP波方程的近似方程，其忽略了部分对称轴方向空间导数项以提高计算效率，但保留了Duveneck方程中微分算子的自共轭性来保证稳定性(Zhang Y and Zhang H Z,2009；Zhang H Z and Zhang Y，2011).本文采用一种新的方式增强qP波传播时的稳定性，采用的是有限横波正则化的方法，这样就可以消除由倾斜对称轴变化以及由ε<δ引起的传播不稳定问题，同时最大程度地减弱波场传播中的伪横波影响.在数值试算及实际资料试算部分，展示了采用所导出方程的正演和逆时偏移结果，证明了方法理论的正确性.

2TTI介质qP波控制方程推导

从平面波传播的角度，从VTI介质向TTI介质转化可看作是在旋转坐标系下的平面波传播.首先，基于Christoffel方程从相速度/频散关系出发，可导出VTI介质耦合频散关系：

(1)

(2)

则(1)式可以重新写为：

(3)

联立(2)、(3)式,并将其反变换到时空域，得到VTI介质中的有限横波方程(Zhou et al.,2006)：

×(p-q).

(4)

同理,通过求解非均匀TTI介质的Christoffel方程,可推导出旋转坐标系下TTI介质P-SV波耦合频散关系:

(5)

由观测坐标与旋转坐标系下波数矢量的关系可得：

(6)

(7)

将(4)式中的二阶微分算子换成(7)式旋转坐标系下的形式，即可得TTI介质中的有限横波方程：

(8)

将(7)式对应的空间微分算子用H1、H2重写可得：

(9)

则有限横波方程(8)可表示为：

(10)

式中，θ是倾角，φ是对称轴的方位角.注意对于VTI介质，微分算子H1只作用于垂直方向，而H2只作用于水平方向.两个算子都不包含空间交叉导数项.在各向同性介质中，式(10)中的两个二阶偏微分方程是等价的，波场p等价于波场q.

方程(10)可以采用伪谱法求解，求解中要用到波数域空间导数计算(Fletcher,2009).然而，由于伪谱法需要多次傅里叶变换，不利于大型数据的波传播并行计算.因此，采用时空域高阶有限差分法求解(10)式.尽管波场p和q都是四阶偏微分方程组的解，但是我们将波场p看作压力场，将q看作辅助波场时,式(10)可看作是两个二阶偏微分方程,与各向同性声波方程类似，可通过高阶差分方法求解(董良国，2000；刘红伟，2010；李博，2010；刘文卿，2012,2013)，本文不再讨论其解法.

3带正则化校正的TTI介质qP波方程推导

方程(10)即为Fletcher等(2009)提出的所谓的有限横波方程计算(即剪切波速度不为0).从方程(10)出发，根据Alkhalifah的思想，将对称轴方向剪切波速度设为0，可得到一个更为简化的TTI介质拟声波方程，这个方程更容易求解，但数值求解也更加不稳定(Fletcher et al.,2009)，Zhang Y和Zhang H Z(2009)采用伪谱法计算时，也发现了类似的不稳定现象.当然，我们可以对模型做平滑处理能使波场稳定传播，但这种方式同时会显著地改变波的运动学特征，影响成像精度.值得注意的是，有限横波方程也并非完全稳定，这种不稳定性通常开始于模型中对称轴倾角存在剧烈变化的地方.本文提出了一种带正则化项的有限横波方法，数值试验表明相较于有限横波方程，这种方程能使波场更加精确稳定地传播下去.下面讨论带正则化项的有限横波方程实现过程.

所谓带正则化项方程，即是设计一个滤波器，对波场传播中不稳定的高波数成分做适当衰减，以达到减少频散或使波场计算稳定的目的.首先来具体推导各向同性介质中带正则化项平面波解的形式，以形象说明正则化项的作用.对于构建的带正则化项的方程：

(11)

方程(11)两边同时对时间/空间做傅氏变化，经过整理并带入平面波解通解形式，可得到各向同性介质下带正则化项声波方程平面波解形式：

(12)

其中k=(kx,ky,kz)，x=(x,y,z)分别为波数以及空间坐标向量,σ为正则化系数(通常取值(0～1)).与不带正则化项的声波方程平面波解相比，(12)式多出了一项振幅衰减项，其物理意义为，对平面波波场中传播时间较长的高波数成分做衰减.利用这一思想，假设波场中不稳定成分主要为高波数情况下，构造TTI介质中带正则化项的有限横波方程，使得波场数值计算更加稳定.对于方程(10)可以构建成带正则化项的有限横波方程：

(13)

4带正则化校正的TTI介质qP波方程的逆时偏移成像方法

带正则化项的方程(13)看起来似乎要比原始的有限横波方程(10)复杂许多，但在实际计算中，其增加的计算量和储存量并不多.原因在于正则化项中旋转坐标系下波场微分无需重复计算，仅需多储存两个波场，记录上一时刻所有正则化项波场之和对应的两个波场，然后对时间求微分即可.实验表明选择合理的σ值，即使倾斜同相轴的方位及倾角变化较大，依然能够保证波场稳定地传播.这也是这种方式实际资料计算中比较可行的原因.

(14)

其初始条件为：

(15)

其中s代表震源位置，f(t)为模拟震源函数.在此RTM实现方案中，首先将边界层最里侧的波场记录下来.二维情况下是四个边界面处的波场记录，三维情况下是六个边界面处的三维波场.视机器内存大小，可存储在内存中，也可存储在硬盘上.

其次，在求解炮波场正向外推基础上再求解波场沿时间逆向外推定解问题，如式(16)—(17).

(16)

其边值条件为：

(17)

其中Bp(x,t)、Bq(x,t)分别是两个波场在边界处的值.利用每个时间层上的边界波场，解上述边值问题，沿时间逆时外推可实现炮波场逆向外推.此时，与检波点波场的逆向外推相同步.因此在每个时间步，利用成像条件进行炮波场与检波点波场零延迟相关成像.

5数值试算

为了测试采用带正则化校正的TTI介质qP波方程波场传播的稳定性,我们选择了FoothillTTI模型中参数变化最为剧烈的部分作为测试模型.图1展示了该模型的参数数据，模型中包含了速度、倾角及各向异性参数.图2展示的是不同控制方程不同时刻的波场快照,图2a在正演方程中采取不带正则化项的有限横波方程，在模型中存在对称轴倾角剧烈变化的地方，波场传播不稳定.模型正演中采用25Hz的雷克子波作为震源，震源放在(3550，0)m处,计算网格在两个方向上均为10m.利用空间10阶、时间2阶有限差分去逼近控制方程中相应导数项.图2c、2d分别为采用带正则化校正的TTI介质qP波方程2000ms、3000ms的波场快照，与图2a、2b同一时刻的波场快照对比可以看出，此时波的传播是稳定的.通过数值试算可看到，直接由频散关系导出的TTI介质中的有限横波qP波方程在忽略参数的空间的导数时，由于受对称轴倾角变化的影响仍然会产生数值求解不稳定的现象，我们通过正则化的方式衰减高波数成分，最终使得整个波场计算稳定.

我们将带正则化校正的TTI介质qP波方程方法应用到TTI逆时偏移中,Foothill模型作为该测试所用模型，模型参数如图1所示，波场传播的计算均是采用高阶有限差分法来求解,我们采用的是时间方向2阶差分，空间方向10阶差分的有限差分算法.图3展示了逆时偏移成像结果，反射层都正确归位，陡倾角地层都能成像.但在地层倾角较大、各向异性较强的区域，其下覆地层影响较大，各向同性及VTI各向异性逆时偏移成像效果都不理想，而TTI各向异性逆时偏移能使陡倾角地层及下覆地层准确成像.因此，对于各向异性介质的地震数据,采用TTI各向异性算法进行成像可得到更为精确的结果.

6实际资料试算与分析

本文所采用的实际数据是我国西部某油田的资料.该区发育大规模岩溶缝洞型储层，非均质性及各向异性都比较明显.以往由于各向异性的影响，人们认为窄方位角资料效果好，如今方位各向异性现象逐步被人们认识，它可以提供更多的裂缝信息和储层信息.因此，常规各向同性地震成像方法由于没有考虑地层各向异性特性，很难准确刻画缝洞的空间 位置.传统的偏移距道集是反射点所有信息的综合，并没有准确包含倾角信息和方位角信息，无法将角度信息从道集中提取出来、也无法精确描述与方位有关的各向异性(王西文，2010).因此，只有全方位的角度道集才能获得准确的全方位速度信息、倾角信息及与方位有关的各向异性参数信息，可以基于全方位共反射角道集进行各向异性速度建模,建模流程如图4.

图1　Foothill模型参数(a) vP速度； (b) 倾角模型； (c) epsilon模型； (d) delta模型.

图2　不同方程传播的波场快照(a) 有限横波方程2 s波场； (b) 有限横波方程3 s波场； (c) 正则化有限横波方程2 s波场； (d) 正则化有限横波方程3 s波场.

图3　Foothill 模型不同方法逆时偏移结果对比(a) 各向同性逆时偏移结果； (b) VTI各向异性逆时偏移结果； (c) TTI各向异性逆时偏移结果.

图5为传统Kirchhoff方法获得的共反射点(CRP)道集与全方位反射角道集对比.从道集特征可看出，常规的CRP道集是反射点所有信息的综合，掩盖了不同方位的速度变化特征及各向异性特征.全方位反射角道集定义了一个新的数据结构(KorenandRavve,2011)，其角度为5°、10°、15°、20°、25°、30°、35°，每一个角度包含全方位的信息.从道集特征可看出，随反射角度增大，相同反射角,不同方位道集不平直，具有不同的速度，其各向异性特征更加明显.利用全方位反射角道集不同方位的剩余时差(图5中绿线)可获得有效的各向异性参数，通过网格层析成像方法建立精确的各向异性速度模型,如图6.

图7a为各向同性逆时偏移成像结果，图7c为TTI各向异性逆时偏移成像结果，由于TTI各向异性逆时偏移成像考虑整个波场信息、地层的各向异性参数及倾角和方位角等空间属性参数，所以串珠成像更聚焦、深度误差更小，成像精度更高.从成像结果及钻井误差对比分析得出，各向同性逆时偏移结果深度误差为190m(图7a)，而各向异性逆时偏移结果深度误差为30m(图7b,图7c),对于TTI各向异性逆时偏移成像由于上覆地层倾角不大，垂向消除深度误差的效果与VTI介质逆时偏移基本相 同，但在横向上串珠位置相差12.5m，与钻井更加吻合.同时深层成像TTI介质各向异性逆时偏移成像相对VTI介质及各向同性介质的成像有较大的改善.

图4　全方位角度域网格层析各向异性速度建模流程

7结论

对于复杂构造成像，叠前逆时偏移是一个有力的工具.宽方位数据采集的普及使得各向异性逆时偏移变得必不可少.本文着重讨论了TI介质下的声学近似、基于频散关系的控制方程推导及高阶有限差分数值解法.基于Christoffel方程从相速度/频散关系出发推导出的qP方程，通过强制方程中沿对称轴方向的剪切波速度为0，可以得到一个简 化的波动方程.这个退化的方程可以提高求解的效率，并且在ε-δ≥0的VTI介质中是稳定的.若将其推广到TTI介质，在速度模型复杂时存在数值计 算的不稳定问题.我们采用一个新的TTI介质拟声波方程及其解决方案，在有限横波基础上提出了一种带正则化项的有限横波方法，保证波场能精确稳定地传播下去，同时最大程度地减弱SV反射系数的影响.

图5　全方位反射角道集(a)与传统CRP道集(b)对比

图6　各向异性参数剖面(a) epslion参数剖面； (b) delta参数剖面； (c)层速度参数剖面； (d) 地层倾角参数剖面； (e)方位角参数剖面.

图7　各向同性逆时偏移(a)、VTI各向异性逆时偏移(b)与TTI各向异性逆时偏移(c)对比

FoothillTTI模型的正演和逆时偏移试验表明，在波场传播步数较多时，常规的有限横波方程也开始出现计算不稳定的现象，通过正则化的方式衰减高波数成分，最终使得整个波场计算稳定，同时成像效果也最好.我们提出的方案在3DTTI介质实际资料中也能得到非常好的效果.在模型参数剧烈变化的区域，上述方法仍然会产生SV波.为了进行P波偏移，这些SV波被看作是噪声.但在实际资料偏移中，我们并没有在最终偏移的剖面上发现来自这些SV能量的污染，可能是SV能量相对于P波能量比较弱的缘故，这些SV波能量是否在RTM角道集上有反映，还有待于进一步研究.从具有实际物理含义及数值计算稳定角度出发，我们应在弹性波框架下做近似，即TI介质中qP控制方程的推导应在时空域进行，这些是下一步研究的方向.

致谢感谢同济大学王华忠老师、周阳博士在研究中给予的指导与帮助.

References

Alkhalifah T.1998.Acoustic approximations for processing in transversely isotropic media.Geophysics,63(2): 623-631.

Alkhalifah T.2000.An acoustic wave equation for anisotropic media.Geophysics,65(4): 1239-1250.

Baysal E,Kosloff D D,Sherwood J W C.1983.Reverse time migration.Geophysics,48(11): 1514-1524.

Chang W F,McMechan G A.1986.Reverse-time migration of offset vertical seismic profiling data using the excitation-time imaging condition.Geophysics,51(1): 67-84.

Chu C L,Nacy B K,Anno P D.2011.An accurate and stable wave equation for pure acoustic TTI modeling.∥ 81st Ann.Internat.Mtg.,Soc.Expl.Geophys.,Expanded Abstracts,179-184.

Dong L G,Ma Z T,Cao J Z,et al.2000.A staggered-grid high-order difference method of one-order elastic wave equation.Chinese J.Geophys.(in Chinese),43(3): 411-419.

Du X,Fletcher R P,Fowler P J.2008.A new pseudo-acoustic wave equation for VTI media.∥ 70thEAGE Conference & Exhibition.Rome,Italy.Duveneck E,Milcik P,Bakker P M,et al.2008.Acoustic VTI wave equations and their application for anisotropic reverse-time migration.∥ 78th Ann.Internat.Mtg.,Soc.Expl.Geophys.,Expanded Abstracts,2186-2190.

Duveneck E,Bakker P M.2011.Stable P-wave modeling for reverse-time migration in tilted TI media.Geophysics,76(2): S65-S75.

Fletcher R P,Du X,Fowler P J.2008.A new pseudo-acoustic wave equation for TI media.∥ 78thAnnual.Internationaln Meeting,SEG,Expanded Abstracts,2082-2086.

Fletcher R P,Du X,Fowler P J.2009.Reverse time migration in tilted transversely isotropic (TTI) media.Geophysics,74(6): WCA179-WCA187.

Koren Z,Ravve L.2011.Full-azimuth subsurface angle domain wavefield decomposition and imaging Part I: Directional and reflection image gathers.Geophysics,76(1): S1-S13.

Li B,Liu H W,Liu G F,et al.2010.Computational strategy of seismic pre-stack reverse time migration on CPU/GPU.Chinese J.Geophys.(in Chinese),53(12): 2938-2943,doi: 10.3969/j.issn.0001-5733.2010.12.017.

Liu F Q,Morton S A,Jiang S S,et al.2009.Decoupled wave equations for P and SV waves in an acoustic VTI media.∥ 79th Annual International Meeting,SEG,Expanded Abstracts,2844-2848.

Liu H W,Li B,Liu H,et al.2010.The algorithm of high order finite difference pre-stack reverse time migration and GPU implementation.Chinese J.Geophys.(in Chinese),53(7): 1725-1733,doi: 10.3969/j.issn.0001-5733.2010.07.024.

Liu W Q,Wang Y C,Yong X S,et al.2012.Prestack reverse time migration on GPU/CPU co-parallel computation.Oil Geophysical Prospecting (in Chinese),47(5): 712-716.

Liu W Q,Wang X W,Liu H,et al.2013.Application of velocity modeling and reverse time migration to subsalt structure.Chinese J.Geophys.(in Chinese),56(2): 612-625,doi: 10.6038/cjg20130225.McMechan G A.1983.Migration by extrapolation of time-dependent boundary values.Geophysical Prospecting,31(3): 413-420.

Sun R,McMechan G A,Lee C S,et al.2006.Prestack scalar reverse-time depth migration of 3D elastic seismic data.Geophysics,71(5): S199-S207.Thomsen L.1986.Weak elastic anisotropy.Geophysics,51(10): 1954-1966.

Wang X W,Liu W Q,Wang Y C,et al.2010.Research and application of common reflection angle domain pre-stack migration.Chinese J.Geophys.(in Chinese),53(11): 2732-2738,doi: 10.3969/j.issn.0001-5733.2010.11.021.

Whitmore N D.1983.Iterative depth migration by backward time propagation.∥53rdAnnual.International Meeting,SEG,Expanded Abstracts,382-385.

Zhang H Z,Zhang Y.2011.Reverse time migration in vertical and tilted orthorhombic media.∥ 81st Ann.Internat.Mtg.,Soc.Expl.Geophys.,Expanded Abstracts,185-189.

Zhang Y,Zhang H Z.2009.A stable TTI reverse time migration and its implementation.∥ 79th Ann.Internat.Mtg.,Soc.Expl.Geophys.,Expanded Abstracts,2794-2798.

Zhou H,Zhang G,Bloor R.2006.An anisotropic acoustic wave equation for VTI media.∥ 68th EAGE Conference & Exhibition.SPE,EAGE.

附中文参考文献

董良国,马在田,曹景忠等.2000.一阶弹性波方程交错网格高阶差分解法.地球物理学报,43(3): 411-419.

李博,刘红伟,刘国峰等.2010.地震叠前逆时偏移算法的CPU/GPU实施对策.地球物理学报,53(12): 2938-2943,doi: 10. 3969/j.issn.0001-5733.2010.12.017.

刘红伟,李博,刘洪等.2010.地震叠前逆时偏移高阶有限差分算 法及GPU实现.地球物理学报,53(7): 1725-1733,doi: 10.3969/j.issn.0001-5733.2010.07.024.

刘文卿,王宇超,雍学善等.2012.基于GPU/CPU叠前逆时偏移研究及应用.石油地球物理勘探,47(5): 712-716.

刘文卿,王西文,刘洪等.2013.盐下构造速度建模与逆时偏移成像研究及应用.地球物理学报,56(2): 612-625,doi: 10.6038/cjg20130225.

王西文,刘文卿,王宇超等.2010.共反射角叠前偏移成像研究及应用.地球物理学报,53(11): 2732-2738,doi: 10.3969/j.issn.0001-5733.2010.11.021.

(本文编辑何燕)

基金项目国家油气专项“天然气地球物理烃类检测评价技术及应用(2016ZX05007-006)资助.

作者简介刘文卿，男，1977年生，在读博士研究生，2000年于成都理工大学获应用地球物理学士学位.主要从事地震速度建模与地震成像方法研究工作.E-mail：liuwq@petrochina.com.cn

doi:10.6038/cjg20160326 中图分类号P631

收稿日期2014-11-30，2015-12-28收修定稿

A regularized qP-wave equation for TTI media and its application to reverse time migration

LIU Wen-Qing1,2,WANG Xi-Wen2,WANG Yu-Chao2,YONG Xue-Shan2,WANG Xiao-Wei2,ZHANG Tao2

1StateKeyLaboratoryofOilandGasReservoirGeologyandExploitation,ChengduUniversityofTechnology，Chengdu610059,China2ResearchInstituteofPetroleumExploration&Development-Northwest,Petrochina,Lanzhou730020,China

AbstractThe research of anisotropy is significant for subsurface precise imaging.With rapid development of computation capability and gradual generalization of seismic acquisition with wide azimuth,anisotropy should be taken into consideration.The reverse time migration is based on the accurate solution of two-way wave equation,therefore compared with other methods,it has many advantages,such as no dip limitation and the ability to image turning waves and multiples.Scalar wave equation could be used to calculate the wave field in isotropic media.In anisotropic media,however,P and SV waves are coupled and there is no pure scalar wave mode.Consequently,pseudo acoustic wave (qP-waves) which can represent the kinematic property of P-component in the coupled wave field are usually used for imaging in anisotropic media.

In this paper,starting with the P-SV coupled dispersion relation for VTI media,we derive qP wave equation for TTI media,which can be solved using an explicit finite difference format.With the help of acoustic approximation,stable numerical solutions can be obtained for a model with vertical symmetric axis,ε≥δ and,if the velocity of shear waves becomes zero along the symmetric axis.However,acoustic approximation will cause the instability of wave field propagation and numerical computation as anisotropic parameters vary in the direction of symmetric axis in TTI media.To solve this problem,we proposed to regularize the finite shear wave equation by adding a regularized term to the original equation.The idea of our method is similar to designing a filter which could attenuate the high wave number components appropriately.We give out the specific form of regularized TTI qP wave equation.Then,the TTI reverse time migration (RTM) method using the regularized wave equation was discussed.In numerical implementation,the computing complexity of our regularized wave equation only increases by one wave field differential and the storage needed only increases by two more wave field.Numerical tests show that by choosing appropriate σ,we can obtain stable wave field even with sharp variation of azimuth and dip of tilted axis.

For numerical examples,we give out the modelling and RTM results on the Foothill TTI model and the imaging result of field data from a carbonate rock area.When prorogating time of the wavefield is much longer,the traditional finite shear wave equation shows wavefield instability.In contrast,our method can obtain wavefield without any instabilities.In areas with abrupt changes of model parameters,our method may still produce weak SV wave which are viewed as artificial noise in reverse time migration.But for real seismic data,we didn’t find any of these SV wave energy in the final migration profile which may be caused by the quite weaker SV energy compared with P wave energy.The RTM imaging with high quality for field seismic data demonstrates that our method works well and can be applied to field data effectively.

KeywordsTTI media; Anisotropy; Wave field coupling; Finite shear wave; Regularization; RTM

刘文卿，王西文，王宇超等.2016.带正则化校正的TTI介质qP波方程及其逆时偏移方法和应用.地球物理学报，59(3):1059-1069,doi:10.6038/cjg20160326.

Liu W Q,Wang X W,Wang Y C,et al.2016.A regularized qP-wave equation for TTI media and its application to reverse time migration.Chinese J.Geophys.(in Chinese),59(3):1059-1069,doi:10.6038/cjg20160326.