2013, Vol. 56
2. 中海油研究总院, 北京 100027
2. Research Institute, CNOOC, Beijing 100027, China
Kirchhoff偏移法是目前最灵活、高效的三维叠前偏移方法.与其它的偏移方法相比，Kirchhoff偏移法具有以下几个特点：1.对速度模型的要求不高，不需要高精度的速度模型即可得到较为精确的成像效果；2.对野外采集数据的适应性很好，能很好的处理各类采集数据；3.计算效率较高，尤其在使用MPI并行的情况下能够提高效率.Kirchhoff偏移法发展至今，对其的研究主要集中在两个方面：1.走时计算：早期的走时计算主要使用直射线理论，该方法简单、直观且对于不复杂的地层成像有着不错的效果，因此一直沿用至今.很多的地震数据处理软件，如Omega、Paradigm等，也都是采用这种方法.但是随着对成像精度要求的不断提升，直射线理论的成像效果远不能满足需要.1969年，Tance & Koehler[1]对Slotnick[2]提出的时距曲线方程进行Taylor展开，形成弯曲射线理论，在不过多增加计算时间的情况下，有效地提高成像效果.随后数十年，不断有学者对该理论方程进行优化完善，提出了多种优化方程，如Castle[3]提出的时移曲线公式，Tsvankin和Thomsen[4]提出的优化动校正方程，Sun等[5]提出的分式优化六阶项方程.这些方程都在一定程度上提高了走时计算的精度.而针对横向变速介质，Liu和Wang [6]提出了非对称走时Kirchhoff叠前时间偏移格林函数；保幅权函数主要是消除几何扩散对地震数据的影响，使反射波振幅得到很好的效果.Kirchhoff的3D保幅叠前偏移最早由Bleistein[7]提出，而后张宇等[8]在此基础上将其完善.本文基于弯曲射线理论和保幅权函数公式，以SEG三维盐丘模型和某区域三维实际资料为研究对象，加入了MPI并行技术，去假频技术，实现三维保幅弯曲射线Kirchhoff叠前时间偏移，并最终应用于三维实际数据.
2 Kirchhoff偏移原理Kirchhoff偏移法是一种绕射求和偏移方法.基于绕射求和原理的偏移方法是在输入空间（x，t）上搜索所有的能量，只要绕射源（惠更斯二次震源）存在，那么在输出空间（x，z）中就能标出它的确切位置.对（x，z）空间中的每一点，在相应的（x，z）空间中沿着它的惠更斯二次震源绕射曲线轨迹进行搜索，把搜索到的各点振幅相加，然后放到（x，z）空间中的这个点上.也就是说，绕射求和就是沿着双曲线轨迹，直接作振幅叠加，这里的叠加双曲线（即求和轨迹）是受速度函数控制的.假设速度-深度模型为水平层状介质，则速度函数即为双曲线顶点的均方根速度.这里的（x，z）是对于深度偏移而言的，在时间偏移前后，就像剖面中所见那样，振幅的和值实际上被放在（x，τ）空间中了，这里的τ是偏移后该点所在位置的时间.
在进行绕射求和之前，必须考虑如下三个因素：
（1） 表明振幅随角度变化的倾斜因子或方向因子，它表示为传播方向与垂直轴Z之间夹角的余弦；
（2） 球面扩散因子，它在2D波动空间中用
（3） 子波整形因子，对2D偏移，设计一个45°常相位谱，振幅谱正比于频率平方根；对3D偏移，这个因子的相移为90°，振幅谱与频率成正比.
而结合了倾斜因子、球面扩散因子和子波整形因子的绕射求和偏移方法称为Kirchhoff偏移法.
Yilmaz[9]的《地震勘探原理》一书中，从Kirchhoff得到的解从求格林函数出发，对Kirchhoff偏移法的数学表达式进行了详细的推导.得到：
|
第一项决定于波场的垂向梯度.第二项叫作近场源项，因为它以1/r2衰减.这两项在地震偏移中都被忽略.剩下的第三项叫作远场源项，它是Kirchhoff偏移的基础.以离散形式写出：
|
式中，Δx和Δy为纵测线和横测线的道间距，Pout为在区域窗口A内用输入波场Pin得到的偏移输出.
Kirchhoff求和法要求：
（1）计算通过三维空间变速介质的非零偏移距走时；
（2）基于标量波动方程的Kirchhoff积分解，沿着计算的走时轨迹的振幅缩放比例并求和.
3 Kirchhoff叠前时间偏移积分法走时计算一般来说，叠前时间偏移基于均匀介质和直射线假设，采用双平方根算子（即DSR公式）计算地震波走时.但是由于DSR方程未考虑射线弯曲对走时的影响，走时难免会有较明显的误差，特别是以与陡倾角或大偏移距相应的广角入射时尤为显著，最终影响偏移成像的质量.而基于射线追踪的Kirchhoff叠前时间偏移是一种高频近似方法，利用射线追踪计算衍射走时需对速度场进行平滑处理，故实际生产中对复杂介质情况鲜见成效.同时其庞大的计算工作量也是制约其发展的主要原因.而本文采取的弯曲射线计算走时，则能在一定程度上弥补DSR公式的长偏移距误差，同时也不会增加大量的计算时间.
该方法的提出源于1969年，Tance&Koehler基于水平层状介质，将Slotnick给出的时距曲线参数方程：
|
在这种方法中，最直观，简单的计算走时，就是不做任何处理的保留高阶项，但是由于该方程式近似得来，并不能保证其在任意偏移距的收敛性，故而不宜直接使用.
刘洋等[11]考虑介质非均质性影响，在高阶时距方程的基础上导出了分式二阶动校正方程：
|
其中VQ=A31/4
本文即利用该动校正方程进行Kirchhoff叠前时间偏移走时计算.
4 保幅权函数随着勘探技术和生产要求的提高，地震勘探除了追求高精度的地下成像，还追求偏移剖面能够反映地下的岩性信息.这些信息在地震波的振幅中有所体现，但是由于地震波在传播过程中受到多方面因素的影响，其振幅往往发生畸变.常规的偏移方法并没有对这些畸变因素进行校正，其结果不能正确地反映地下的岩性信息.保幅偏移正是基于这个问题提出的.
Bleistein[12]首先给出了3D Kirchhoff叠前时间偏移的保幅权函数表达式：
|
其中，ξ描述了炮点和检波点之间的参数坐标，A（x，xs）、A（x，xr）分别表示炮点到成像点、成像点到检波点的振幅；|h（x，ξ）|称作Beylkin行列式，描述了地面坐标到地下成像坐标的转换关系，可以用下式表达：
|
其中，
而后Schleicher等[13]、Winbow等[14]、Zhang等[15]、Lee等[16]分别对该保幅权函数进行简化. Zhang等推导的近似权函数
首先用三维SEG盐丘模型对本文的Kirchhoff叠前时间偏移进行了检验，而后又将该方法应用到某区域的一套三维地震资料中.在运用过程中，为了提高计算速度和增强成像效果，运用了如下技术：
5.1 去假频技术Kirchhoff积分叠前时间偏移可能会受到三种假频的影响，分别为（1）由于时间采样不足所导致的输入数据假频；（2）由于偏移计算方法不当所产生的算子假频；（3）由于偏移输出网格选取不当所产生的成像假频.在这里，我们主要讨论算子假频问题.
早在1994年，Lumley等[17]就给出了二维Kirchhoff偏移的算子假频，并且得到理论界的广泛认同.但是对于如何正确地推广到三维则一度成为争论的焦点，Abma等[18]、Wang [19-20]分别提出了Kirchhoff时间偏移三维算子假频不同的公式表述.本文采用的是Lumley提出的三角滤波去假频技术.
5.2 MPI并行技术随着并行计算机技术的不断发展，高效、大尺度、大数据量的偏移成像变得可能.Kirchhoff叠前时间偏移计算有着两大特点：一是单道处理实现偏移，没有道与道之间的交叉运算；二是各成像点单独成像，各成像点之间没有交叉运算，这样的特点使得其实现并行计算变得很简单.只需要五个最基本的并行函数（MPI_INIT，MPI_COMM_SIZE，MPI_ COMM_RANK，MPI_BARRIER，MPI_FINALIZE），将输入的地震道数据按照一定的顺序（按炮或者按道）放入各自节点进行计算，将各节点的运算结果输出并相加就得到了最终的偏移结果.这样的MPI并行计算技术，能够大大的提高运算效率，减少运算时间.
5.3 SEG盐丘模型结果利用三维保幅弯曲射线Kirchhoff叠前时间偏移方法对SEG三维模型进行偏移成像，得到整个盐丘的偏移结果，并从中提取出两个方向的偏移剖面，与理论模型进行了对比.盐丘模型大小为676×676× 200（dx=dy=dz=20m），成像网格大小为20m× 20m，检波距为20m.共有45炮，每炮有201×201个检波点接收（即40401道），625个采样点，采样间隔8ms.
图 1显示的是SEG三维盐丘模型Inline方向的其中一个剖面，其为深度域剖面.该剖面由一些比较平缓，倾角不大的地层构成，右侧边缘有一陡倾角地层，同时在下部地层出现速度的反转.图 2显示的是本文所用方法所得到的该模型对应的偏移剖面，其为时间域剖面.该剖面能够比较清晰的看到前两个地层，但是有速度反转的下部地层十分模糊，即使成像也与模型不吻合（最底部一层两侧向上弯曲），断层和陡倾角也没有成像出现.
|
图 1 Inline模型剖面 Fig. 1 Model section of Inline |
|
图 2 Inline偏移剖面 Fig. 2 Migration section of Inline |
图 3显示的是SEG三维盐丘模型Crossline方向的其中一个剖面，其为深度域剖面.该剖面与Inline相比，较为复杂.不仅有陡倾角、断层，还有速度的反转和突变.图 4显示的是本文所用方法所得到的该模型对应的偏移剖面，其为时间域剖面.
|
图 3 Crossline模型剖面 Fig. 3 Model section of Crossline |
|
图 4 Crossline偏移剖面 Fig. 4 Migration section of Crossline |
从该剖面能够比较清晰地看出前四个地层，但是有速度反转的下部地层比较模糊（底部）.剖面两侧未能成像，原因在于模型的道集资料未能覆盖该区域.
分析上述剖面可以看出，对于一般地层，Kirchhoff叠前时间偏移方法可以较好的成像，但是对于有速度反转的情况，则很难得到较好的成像效果.在有速度反转的情况下，不但是在速度反转区域不能很好的成像效果，而且会影响到其上部地层的成像效果.因此，在选择使用Kirchhoff叠前时间偏移方法前，一定要先预测该区域有无速度反转区存在.
5.4 三维实际数据利用三维保幅弯曲射线Kirchhoff叠前时间偏移方法对某区域三维模型进行了偏移成像，得到了Inline和Crossline方向的剖面，表明本文提出的方法能够适用于三维真实地震数据的偏移成像.
图 5是常规方法的Inline方向偏移结果，图 6是本文方法的Inline方向偏移结果.通过与常规方法的偏移结果对比可以看出，本文提出的方法能够有效的提高成像结果的精确性和信噪比，能够展示常规偏移剖面不具备的、更为丰富的地层和构造信息.图 7是图 5的局部放大，图 8是图 5的局部放大.图 9是两种方法的细节对比.
|
图 5 常规方法偏移剖面 Fig. 5 Migration section of conventional method |
|
图 6 本文方法偏移剖面 Fig. 6 Migration section of presented method |
|
图 7 常规方法偏移剖面局部放大图 Fig. 7 Migration section of conventional method (zoom) |
|
图 8 本文方法偏移剖面局部放大图 Fig. 8 Migration section of presented method (zoom) |
|
图 9 两种方法对比图 (a)常规方法; (b)本文方法. Fig. 9 Comparison of two method (a) Conventional method; (b) Presented method. |
本文所提出的三维保幅弯曲射线Kirchhoff叠前时间偏移，首先是弯曲射线法的走时计算确保了偏移成像的准确性，其次保幅权函数的加入有效地提高了性噪比，同时结合去假频技术和MPI并行技术，提高了成像的分辨率和计算的效率.该方法与常规Kirchhoff叠前时间偏移，有较好的成像结果.
该方法应用于SEG盐丘模型数据可以看到，在没有速度反转的区域，成像效果较好，出现反转则无法精确成像.因此我们建议在使用该方法前，对目标区域有提前的判断，在没有速度反转的情况下使用该方法.同时该方法应用于某区域三维实际资料得到了令人满意的结果，可见该方法可以应用于工业生产.综合而言，该方法是处理三维数据成像中较简洁，快速，高效的算法.
| [1] | Taner M T, Koehler F. Velocity spectra-digital computer derivation and applications of velocity functions. Geophysics , 1969, 34(6): 859-881. DOI:10.1190/1.1440058 |
| [2] | Slotnick M. Lessons in Seismic Computing: A Memorial to the Author. Lawrence: Society of Exploration Geophysicists, 1959 . |
| [3] | Castle R J. A theory of normal moveout. Geophysics , 1994, 59(6): 983-999. DOI:10.1190/1.1443658 |
| [4] | Tsvankin I, Thomsen L. Nonhyperbolic reflection moveout in anisotropic media. Geophysics , 1994, 59(8): 1290-1304. DOI:10.1190/1.1443686 |
| [5] | Sun C W, Wang H W, Nartinez R D. Optimized 6th order NMO correction for long-offset seismic data. 72nd Annual International Meeting, SEG, Expanded Abstracts , 2002. |
| [6] | 刘洪, 王秀闽, 曾锐, 等. 单程波算子积分解的象征表示. 地球物理学进展 , 2007, 22(2): 463–471. Liu H, Wang X M, Zeng R, et al. Symbol description to integral solution of one-way wave operator. Progress in Geophysics (in Chinese) , 2007, 22(2): 463-471. |
| [7] | Bleistein N, Cohen J, Jaramillo H. True-amplitude transformation to zero offset of data from curved reflectors. Geophysics , 1999, 64(1): 112-129. DOI:10.1190/1.1444509 |
| [8] | 张宇. 振幅保真的单程波方程偏移理论. 地球物理学报 , 2006, 49(5): 1410–1430. Zhang Y. The theory of true amplitude one-way wave equation migration. Chinese J. Geophys. (in Chinese) , 2006, 49(5): 1410-1430. |
| [9] | Yilmaz O Z. Seismic Data Analysis. Tulsa: Society of Exploration Geophysicists, 2001 . |
| [10] | Hake H, Helbig K, Mesdag C S. Three-term Taylor series for t2 -x2-curves of P-and S-waves ober layered transversely isotropic ground. Geophysical Prospecting , 1984, 32(5): 828-850. DOI:10.1111/gpr.1984.32.issue-5 |
| [11] | 刘洋. 反射波分式展开时距方程及其精度分析. 石油物探 , 2003, 12(4): 441–447. Liu Y. Fraction expansion time-distance equation of reflection wave and its accuracy analysis. Geophysical Progress for Petroleum (in Chinese) , 2003, 12(4): 441-447. |
| [12] | Bleistein N. On the imaging of reflectors in the earth. Geophysics , 1987, 52(7): 931-942. DOI:10.1190/1.1442363 |
| [13] | Schleicher J, Tygel M, Hubral P. 3-D true-amplitude finite-offset migration. Geophysics , 2005, 58(8): 1112-1126. |
| [14] | Winbow G A, Schneider W A. Weights for 3-D controlled amplitude prestack time migration. 69th Annual International Meeting, SEG, Expanded Abstracts , 1999: 1110-1113. |
| [15] | Zhang Y, Gray S, Yong J. Exact and approximate weights for Kirchhoff migration. 70th Annual International Meeting, SEG, Expanded Abstracts , 2000: 1036-1039. |
| [16] | Lee S, King D, Lin S. Efficient true-amplitude weights in Kirchhoff time migration. 74th Annual International Meeting. SEG, Expanded Abstracts , 2004: 1089-1092. |
| [17] | Lumley D E, Claerbout J F, Bevc D. Anti-aliased Kirchhoff 3-D migration. 64th Annnual Internat. Meeting, SEG, Expanded Abstracts , 1994: 1282-1285. |
| [18] | Abma R, Sun J, Bernitsas N. Antialiasing methods in Kirchhoff migration. Geophysics , 1999, 64(6): 1783-1792. DOI:10.1190/1.1444684 |
| [19] | Wang C S. Kirchhoff antialiasing methods in the imaging domain. SEG Annual Meeting , 2000. |
| [20] | Wang C S. Accurate spatial directional derivative: Antialiasing Kirchhoff migration. Geophysics , 2004, 69(3): 811-816. DOI:10.1190/1.1759467 |