2018, Vol. 61
2. 中国石油大学地球科学与信息学院, 山东青岛 266555
2. School of Geosciences, China University of Petroleum, Shandong Qingdao 266555, China
随着陆地以及海上OBC多分量地震勘探采集的广泛应用，针对于转换波地震数据的处理方法正在得到越来越多的重视(刘国峰等，2010；杜启振等，2011；岳玉波等，2012；王华忠等，2013；孙章庆等，2016).由于转换波共转换点是随时间和空间变化的，常规的基于NMO+DMO+叠后偏移的处理流程，难以获得理想效果.因此，叠前时间偏移对转换波地震数据的成像处理来说尤为重要.
虽然转换波Kirchhoff时间偏移同常规的纵波时间偏移具有相同的理论基础(Kuo and Dai, 1984; Schleicher et al., 1993)，但具体的实现过程却存在较大的差异，这是因为转换波的下行P波和上行S波对应着不同类型的传播速度，而且转换波的反射路径并不对称.因此，需要根据转换波数据的特点对转换波偏移的实现过程进行特殊设计.本文对VTI介质转换波Kirchhoff叠前时间偏移实现过程中的若干关键技术环节进行了研究，目的在于在保证计算效率的同时，提升偏移的成像效果.首先，根据Zhang等(2000)、Dellinger等(2000)提出的纵波偏移加权函数的简化思路，推导了简化的转换波时间偏移振幅加权函数，避免了处于偏移最内层循环的加权函数复杂的计算过程，从而兼顾了偏移计算的高效和成像振幅的准确性；其次，根据直射线近似假设和采样定理(Abma et al., 1999; Zhang et al., 2001)，推导了转换波偏移的最大无假频频率计算公式，确保了算子假频的有效滤除；接下来，根据转换波传播路径的特点(Lüth et al., 2005)，选择使用了依据共转换点(CCP)对称的偏移孔径，从而保证了偏移孔径的优化选取.
在Kirchhoff偏移的实现过程中，为了保证计算效率，往往预先将成像道网格点抽稀为粗网格，并且在粗网格上进行走时等必要参数的计算，然后利用计算效率较高的线性插值来计算细成像网格上的参数值，进而进行偏移信号的计算.虽然线性插值在偏移过程中用于插值缓慢变化的走时等参数场时精度已经足够，然而线性插值只是一种零阶的插值算法，其无法满足Kirchhoff偏移所必须的参数场至少一阶的连续性，从而导致偏移信号在插值节点处的不平滑，所造成的高频噪声能量，足以对偏移信号中的高频成分产生不良影响，尤其是在偏移叠加剖面中，由于叠加过程中上述高频噪声能量得到增强，使得叠加剖面信号中的高频成分被完全破坏.为解决上述问题，本文提出利用三次卷积插值算法进行上述插值计算过程(Keys, 1981)，该方法可以保证插值参数场至少一阶的连续性，而且同线性插值相比，该方法并没有增加太多的计算量.
1 方法原理Kirchhoff偏移的实质是对地震反射记录沿绕射走时面的积分叠加(Schleicher et al., 1993；Bleistein et al., 2001).因此，地震波走时的计算是Kirchhoff偏移的核心，目前通用的转换波走时计算公式为(Li et al., 2007):
|
(1) |
式(1)适用于二维、三维VTI介质转换波走时的计算.其中，
Kirchhoff偏移通过加权函数来确保成像结果的相对振幅保持.Bleistein等(2001)给出的振幅加权函数为
|
(2) |
式中|∇(tp+ts)|为走时梯度项，Ap、As分别为震源和接收点到成像点的射线理论振幅，h为Beylkin行列式.由于加权函数的计算位于偏移过程的最内层循环，直接对上式进行计算，计算量巨大，严重影响计算效率.
针对该问题，Zhang等(2000)、Dellinger等(2000)提出通过假设地下介质为常速来简化加权函数的计算，在保证计算效率的同时，尽可能的保证成像振幅的准确性.在此应用该思路，进行简化的转换波偏移加权函数的推导.假设纵、横波速度分别为Vp和Vs，速度比为γ，则式(2)中的各个参量可以得到极大的简化，公式为
|
(3) |
|
(4) |
|
(5) |
其中，θ是反射张角，z是成像点的深度，z= 
|
(6) |
算子反假频是Kirchhoff偏移过程中需要特别注意的环节(Abma et al., 1999).正确的转换波偏移反假频公式需要考虑不同波型速度对算子空间导数的影响.根据尼奎斯特采样定理，在三维情况下，转换波Kirchhoff偏移的最大无假频频率可以表示为
|
(7) |
其中，Δx和Δy为inline和crossline方向的道间距，∂t/∂x和∂t/∂y为偏移算子的空间偏导数，其计算为反假频的关键，其计算过程需要考虑上行转换波的影响.在此利用直射线近似进行计算，得：
|
(8) |
其中，tp、ts分别是下行P波和上行S波的单程走时.rpx、rpy、rsx、rsy分别为连接震源、接收点到成像点的空间矢量在inline和crossline方向的分量.在求取fmax后，便可以通过设计使用合适的滤波器进行算子反假频处理，在Kirchhoff偏移中最常用的是三角滤波器(Lumley et al., 2001).
1.3 偏移孔径偏移孔径的选取关乎偏移计算效率和成像效果(Sun, 2000).根据转换波的传播特点，优化的转换波的偏移孔径应当沿着CCP来对称选取(图 1).
|
图 1 转换波偏移孔径 Fig. 1 Aperture definition of PS-wave migration |
Thomsen (1999)给出了的转换波共转换点计算公式.在层状介质假设条件下，转换波共转换点到震源的水平距离xccp可以表示为
|
(9) |
其中，h为炮检距，tps0为转换波垂向双程走时，Vps为转换波均方根速度，
在Kirchhoff偏移的实现过程中，往往首先在预先定义的粗成像网格上计算走时等参数信息，然后通过插值算法将上述参数插值到细成像网格上，以提高偏移计算效率.线性插值是Kirchhoff偏移中最常用的插值算法，然而线性插值只是一种低阶的插值方法，无法保证插值后参数场一阶的连续性，会导致偏移信号在插值节点处的不平滑，所造成的高频噪声能量足以对信号中的较弱的高频成分产生不良影响，尤其是在最终的偏移叠加剖面中，由于偏移叠加过程中插值节点的高频噪声得到增强，使得叠加剖面中的高频成分被完全破坏.
针对上述问题，本文选择三次卷积插值(Keys, 1981)代替线性插值进行细网格上参数场的插值运算.三次卷积插值是一种局部的分段逼近sinc函数的三次函数，具备一阶的连续性，简化后的插值公式为
|
(10) |
其中，x代表细网格成像点，g(x)为插值结果，xk－1、xk、xk+1、xk+2分别为x点临近的四个粗网格成像点(图 2)，ck－1、ck、ck+1、ck+2为对应的参数值，s=(x－xk)/h，h粗网格的采样间隔.在实际计算中，可以预先将上式中的插值系数表计算出来，这样每次细网格上的插值运算，只需要四次乘法和三次加法，同线性插值(两次加法，两次乘法)相比所增加的计算量，并不会对偏移整体的运行时间产生太大影响.
|
图 2 三次卷积插值原理 Fig. 2 Sketch of cubic convolution interpolation |
加权函数的作用是在偏移过程中消除地震波几何扩散对成像振幅的影响，从而恢复地下界面正确的反射率，本例通过一个简单的三维褶积合成数据进行加权函数的测试.假设地下纵横波速度为常数，并且存在三个因密度差异而产生的反射界面，每个界面的反射系数都是相等的常数，地震波的能量衰减仅由几何扩散导致.图 3a展示了合成转换波记录的二维剖面，图 3b为沿三层反射记录提取的归一化振幅，可以看到由于几何扩散的影响，反射信号振幅存在不同程度的衰减.应用本文所给出的加权函数对合成记录进行偏移运算，得到如图 3c所示的成像结果，同样沿反射层位提取归一化振幅，得到的结果如图 3d所示，可以看到偏移结果在一定的照明范围内正确的恢复了地下界面的反射系数，证明了本文所给出的加权函数的正确性.
|
图 3 加权函数测试 (a)转换波共炮检距记录; (b)沿图(a)中各反射层提取的归一化振幅; (c)偏移成像结果; (d)沿图(c)中各成像层位提取的归一化振幅. Fig. 3 Migration tests of weight functions (a) PS-wave common offset gather; (b)Normalized amplitudes extracted along the three reflection events in (a); (c) Migrated image; (d) Normalized amplitudes extracted along the three migrated events in (c). |
分别利用一个实际数据的脉冲响应和叠加成像结果测试反假频的应用效果.图 4a和图 4b分别为应用常规纵波反假频准则以及应用本文转换波反假频准则的算子脉冲响应结果，可以看到图 4a中脉冲响应的陡倾角部分存在明显的假频能量，而在图 4b中上述假频能量则得到了有效的滤除.由此可以证明，在转换波Kirchhoff偏移过程中，需要考虑S波对偏移算子的影响，由于S波速度明显低于P波速度(尤其在浅层)，公式(7)和公式(8)所计算得到的转换波偏移的最大无假频频率要明显低于纵波计算结果.
|
图 4 脉冲响应对比 (a)应用纵波反假频的脉冲响应; (b)应用本文转换波反假频的脉冲响应. Fig. 4 Comparison of impulse response (a) PP-wave anti-aliasing formula; (b) PS-wave anti-aliasing formula. |
图 5a和图 5b分别为应用纵波反假频以及应用本文转换波反假频的偏移叠加结果，可以看到图 5a中存在大量的算子假频造成的噪声能量，严重影响了浅层的成像质量，而图 5b则有效的消除了上述噪声干扰，成像剖面的连续性和信噪比有了明显提高.
|
图 5 反假频叠加结果对比 (a)应用纵波反假频的偏移叠加结果; (b)应用本文转换波反假频的偏移叠加结果. Fig. 5 Comparison of anti-aliasing stacked images (a) PP-wave anti-aliasing formula; (b) PS-wave anti-aliasing formula. |
接下来进行偏移孔径应用效果的测试.图 6a和图 6b分别为应用基于ACP(渐进转换点)对称孔径和基于CCP对称孔径的偏移叠加结果，在对比的过程中使用了相同的孔径函数，可以看到同CCP孔径偏移结果相比，ACP孔径偏移结果丢失了很多浅层成像构造，证明了基于CCP对称偏移孔径的优势.虽然可以通过扩大浅层孔径来改善基于ACP孔径的偏移结果，但是势必会增加不少的计算量，而且会引入更多的偏移画弧噪声.
|
图 6 叠加结果对比 (a)应用纵波反假频的偏移叠加结果; (b)应用本文转换波反假频的偏移叠加结果. Fig. 6 Comparison of Stacked images (a) ACP-based aperture; (b) CCP-based aperture. |
最后分别采用脉冲响应和偏移叠加结果来验证不同插值算法对偏移剖面中高频成分的影响，由于两种插值算法的全频带成像结果没有可以分辨的区别，其差异只体现在信号的高频端，因而在此只显示应用高通滤波(＞40 Hz)后的成像结果.图 7a和图 7b分别为应用线性插值的偏移脉冲响应和偏移叠加结果，可以看到插值节点上产生的平层状的高频能量，明显破坏了剖面中的有效信号成分，在叠加剖面中更为明显.图 8a和图 8b分别显示了应用三次卷积插值的偏移脉冲响应和偏移叠加结果，可以看到成像结果完全不存在图 7中所示的高频干扰，高频成分得到了较好的保留.
|
图 7 应用线性插值的转换波偏移结果 (a)偏移脉冲响应; (b)偏移叠加结果. Fig. 7 Filtered images obtained with linear interpolation (a) Impulse response; (b) Stacked image. |
|
图 8 应用三次卷积插值的转换波偏移结果 (a)偏移脉冲响应; (b)偏移叠加结果. Fig. 8 Filtered images obtained by cubic convolution interpolation (a) Impulse response; (b) Stacked image. |
由于具有较高的计算效率和灵活性，Kirchhoff偏移依然是当前业界内最常用的偏移成像方法.虽然理论基础较为简单，但是其计算效率和应用效果却很大程度上依赖于算法的设计和关键环节的实现技巧.本文对转换波Kirchhoff叠前时间偏移过程中的几个关键技术环节进行了讨论分析，给出了相应的解决方案和实现思路，并且通过数值算例的测试证明了上述算法的正确性和有效性.同时，本文还提出利用三次卷积插值代替线性插值进行细成像网格上走时等参数场的计算，有效的避免了插值节点处的高频噪声干扰，保证了成像结果中高频成分的正确性.
Abma R, Sun J, Bernitsas N.
1999. Antialiasing methods in Kirchhoff migration. Geophysics, 64(6): 1783-1792.
DOI:10.1190/1.1444684 |
|
Bleistein N, Cohen J K, Stockwell J W. 2001.
Mathematics of Multidimensional Seismic Imaging, Migration, and Inversion. New York: Springer.
|
|
Dellinger J A, Gray S H, Murphy G E, et al.
2000. Efficient 2.5-D true-amplitude migration. Geophysics, 65(3): 943-950.
DOI:10.1190/1.1444790 |
|
Du Q Z, Li F, Hou B, et al.
2011. Angle-domain migration velocity analysis based on elastic-wave Kirchhoff prestack depth migration. Chinese Journal of Geophysics, 54(5): 1327-1339.
DOI:10.3969/j.issn.0001-5733.2011.05.022 |
|
Keys R G.
1981. Cubic convolution interpolation for digital image processing. IEEE Transactions on Acoustics, Speech, and Signal Processing, 29(6): 1153-1160.
DOI:10.1109/TASSP.1981.1163711 |
|
Kuo J T, Dai T F.
1984. Kirchhoff elastic wave migration for the case noncoincident source and receiver. Geophysics, 49(2): 1223-1238.
|
|
Li X Y, Dai H C, Mancini F.
2007. Converted-wave imaging in anisotropic media:theory and case studies. Geophysical Prospecting, 55(3): 345-363.
DOI:10.1111/gpr.2007.55.issue-3 |
|
Liu G F, Liu H, Li B, et al.
2010. Direct pre-stack time migration on rugged topography. Oil Geophysical Prospecting, 45(2): 196-200.
|
|
Lumley D E, Claerbout J F, Bevc D. 2001. Anti-aliased Kirchhoff 3-D migration. //71th Ann. Internat Mtg., Soc. Expi. Geophys. . Expanded Abstracts.
|
|
Lüth S, Buske S, Giese R, et al.
2005. Fresnel volume migration of multicomponent data. Geophysics, 70(6): S12-S129.
|
|
Schleicher J, Tygel M, Hubral P.
1993. 3-D true amplitude finite-offset migration. Geophysics, 58(8): 1112-1126.
DOI:10.1190/1.1443495 |
|
Sun J.
2000. Limited-aperture migration. Geophysics, 65(2): 584-595.
DOI:10.1190/1.1444754 |
|
Sun Z Q, Sun J G, Yue Y B, et al.
2016. Ray tracing using GMM-ULTI method in the random media with complex topography. Chinese Journal of Geophysics, 59(7): 2615-2627.
DOI:10.6038/cjg20160725 |
|
Thomsen L.
1999. Converted-wave reflection seismology over inhomogeneous, anisotropic media. Geophysics, 64(3): 678-690.
DOI:10.1190/1.1444577 |
|
Wang H Z, Liu S Y, Yang Q Y, et al.
2013. Exploration strategy and imaging method in mountain area. Oil Geophysical Prospecting, 48(1): 151-159.
|
|
Yue Y B, Li Z C, Qian Z Q, et al.
2012. Amplitude-preserved Gaussian beam migration under complex topographic conditions. Chinese Journal of Geophysics, 55(4): 1376-1383.
DOI:10.6038/j.issn.0001-5733.2012.04.033 |
|
Zhang Y, Gray S, Sun J, et al. 2001. Theory of migration anti-aliasing. //71th Ann. Internat. Mtg., Soc. Expi. Geophys. Expanded Abstracts, 997-1000.
|
|
Zhang Y, Gray S, Young J. 2000. Exact and approximate weights for Kirchhoff migration. //70th Ann. Internat. Mtg., Soc. Expi. Geophys. Expanded Abstracts, 1036-1039.
|
|
杜启振, 李芳, 侯波, 等.
2011. 角度域弹性波Kirchhoff叠前深度偏移速度分析方法. 地球物理学报, 54(5): 1327–1339.
DOI:10.3969/j.issn.0001-5733.2011.05.022 |
|
刘国峰, 刘洪, 李博, 等.
2010. 起伏地表直接叠前时间偏移. 石油地球物理勘探, 45(2): 196–200.
|
|
孙章庆, 孙建国, 岳玉波, 等.
2016. 复杂山地随机介质GMM-ULTI法射线追踪. 地球物理学报, 59(7): 2615–2627.
DOI:10.6038/cjg20160725 |
|
王华忠, 刘少勇, 杨勤勇, 等.
2013. 山前带地震勘探策略与成像处理方法. 石油地球物理勘探, 48(1): 151–159.
|
|
岳玉波, 李振春, 钱忠平, 等.
2012. 复杂地表条件下保幅高斯束偏移. 地球物理学报, 55(4): 1376–1383.
DOI:10.6038/j.issn.0001-5733.2012.04.033 |
|