2010, Vol. 53
2. 中国石油化工股份有限公司石油物探技术研究院, 南京 210014
2. SINOPEC Geophysical Research Institute, Nanjing 210014, China
在1978年Kirchhoff积分法叠前时间偏移就已被提出[1]，由于其原理简单，计算效率高并能适应不规则观测系统，目前在实际生产中得到广泛应用[2].Kirchhoff积分法偏移在理论方面比较简单，主要有三个核心问题：反假频问题、偏移孔径自适应计算问题[3~5]、走时计算方法问题.反假频方法方面有空变滤波反假频的方法[6, 7]，有三角滤波反假频方法[8]等，本研究采用前者.走时计算方法有常规的直射线方法、弯曲射线方法[9, 10]、非对称走时方法[11]等.当然，更精细的成像还要考虑采集脚印引起的加权算子[12, 13]的问题.为了后续的AVA分析，有输出角度道集的Kirchhoff积分法叠前时间偏移[14, 15].
一般情况下，Kirchhoff偏移的实现是按输出道的观点把一个输入的地震道的振幅按Kirchhoff积分偏移公式发散到不同时刻的等时面上，这是一个单道处理的方法，因此其具体的实现非常灵活.对于一个上千平方公里探区的巨大规模数据体的积分法叠前偏移，要开发出效率较高的、具有很强适用性的积分法叠前偏移程序，这其中还有很多问题需要解决.为了实现巨大规模数据体的高效率Kirchhoff积分法叠前偏移成像，本文着力于解决内存需求量、I/O量及单道成像的精度和效率这三方面的问题.三方面问题的解决达到很好平衡的方法就是一个较好的Kirchhoff积分叠前偏移的实现方案.
本文提出了面向大规模数据处理的Kirchhoff积分法叠前时间偏移实现方案，数值实验验证了本文提出的实现方案的优越性.在此较完善的实现方案的基础上，可以容易地把更优越的积分类偏移方法迅速推向实用化.
2 Kirchhoff积分法叠前体偏移的实现方案一般地，进行Kirchhoff积分叠前偏移的思路是按成像线输出成像结果.这种方案适用于小数据体.但是，在进行大规模数据的积分法叠前成像处理时，这种方案的缺陷就非常明显：（1）由于偏移孔径很大，输入数据的重复次数很多，可以用偏移孔径除成像线间距来评估输入数据的成像次数.这大大增加了I/O量.（2）Corss-line方向的偏移孔径是事先人为确定的，不是自适应计算得到的.
根据Kirchhoff积分偏移是单道进行，并产生等偏移距成像结果的特点，可以模拟实现共偏移距数据体的偏移.这样叠加前的全局成像数据体大小等于叠加后的全局成像数据体，放在每个计算节点附带的局部盘上.待所有等偏移距的成像数据体收集在一起，放到全局盘上，通过重排序，就可以产生叠前成像道集.成像道集加工处理后叠加，就形成最终的叠前偏移成像结果.
本文设计的实现方案为：把数据体按共偏移距道集排放，每一个共偏移距数据体单独成像，把3D Kirchhoff积分叠前偏移问题变成一个类似叠后3D Kirchhoff积分偏移问题.这样可以使每个共偏移距道集的成像结果在内存和局部盘（每个计算节点挂载的硬盘）之间通讯.当前共偏移距道集偏移完成后，把成像结果移到全局盘上.这样的策略可以用当前的大多数PC-Cluster机群系统实现巨大规模数据体的Kirchhoff积分叠前偏移成像处理.接下来，每一输入道的每一时刻的振幅按Kirchhoff积分偏移公式发散到整个成像体空间中.这要求每个共偏移距道集成像结果必须驻留在内存中，但是巨大数据体成像可能的内存申请将超出计算机允许的最大内存.我们把拟深度样点从浅到深分成若干份，每一段的全部横向成像空间可以驻留在内存中，当前段的体偏移完成后写到本地盘上.清空内存，进行下一拟深度段的Kirchhoff积分叠前体偏移.每一输入道按Kirchhoff积分偏移公式的处理是在炮检连线规定的旋转坐标系中进行的，炮检连线在旋转坐标系中与横轴平行，在这样的坐标系中进行一个地震道的Kirchhoff积分叠前时间偏移比较容易计算偏移孔径、处理反假频.然后，把在旋转坐标系中的成像结果投影回到原始成像坐标系中.这种一次输出成像数据体的方案称为体偏移方案，按线输出成像结果的方案称为线偏移方案.
Kirchhoff积分叠前偏移成像的并行实现方案也是提高效率的重要方面，这主要考虑PC-Cluster机群的内存、局部盘、全局盘和I/O通讯效率与Kirchhoff积分偏移的实现特点之间的最佳匹配.李建江等[16]从一般意义上分析了基于共享存储的叠前深度偏移的并行算法，这对当前以PC-Cluster机群为主的Kirchhoff积分叠前偏移的实现具有普遍适用的意义，本文提出的并行方案也是基于主从模式.实际上，作业粒度大小的划分及其与计算机结构与性能的匹配是更重要的问题.刘国锋等[11]提出了Kirchhoff积分法的两种并行实现方案，其中的成像点并行方案是每一个处理器面向所有的地震道，这一方案在面对巨大规模数据体及计算节点较多时，会因每个计算节点读入地震道的次数太多而导致I/O量太大，计算性能达不到最优，因此这并非一个最佳方案.这里提出的方案是把共偏移距地震数据的共偏移距号作为一级索引，每个共偏移距道集内的线号作为二级索引来给主从模式中的每个计算节点分配任务.这样基本兼顾了计算机性能与积分法偏移本身的特点，达到计算与I/O的基本平衡.
因此，该方案的流程为：（1）把数据按偏移距进行排序.（2）对每个共偏移距数据体进行并行偏移.按共偏移距数据体中的In-line线分配作业.每条线一个进程.每个进程的成像输出空间大小等于叠加后的全局成像数据体，可以放在局部盘上.（3）一个共偏移距数据体偏移完成后，归约收集每个处理器对应的局部盘上的成像结果.把归约叠加后的结果移到全局盘上.清除局部盘上的空间，为下一个共偏移距数据的偏移做准备.（4）所有偏移距的数据偏移完成，全局盘上就有了所有偏移距数据对应的成像结果.重排序可以得到每点的成像道集.成像道集修饰叠加可以得到最后的PSTM成像结果.
该方案的优点包括：（1）三维数据体中的每一道仅仅被输入一次，形成该道对应的由等时面构成的三维偏移响应时，计算量不会重复.（2）Cross-line方向的偏移孔径可以计算得出，不必人为填写.得到的成像质量更高.
3 单个输入道的Kirchhoff积分叠前体偏移 3.1 基本原则单个输入道的Kirchhoff积分叠前体偏移，原则上可以在原始坐标系中进行，如图 1中由（x，y）定义的坐标系.但是，由于成像剖面与旋转椭球体表示的等时面斜交，自适应偏移孔径的计算和反假频处理都不方便.为此，我们在由输入道炮检点连线定义的旋转坐标系中进行单个输入道的Kirchhoff积分叠前体偏移.计算每个地震道对应的由等时面构成偏移响应的具体做法如下.
|
图 1 旋转坐标系下三维旋转椭球体形状的成像响应 Fig. 1 3D imaging response in rotated coordinates system |
由于在（x'，y'）坐标系中，炮检连线平行于x'坐标，使得计算便利.因此先在（x'，y'）坐标系中计算三维旋转椭球体形状的成像响应，然后将该成像响应椭球体映射到（x，y）坐标系中，这是真实成像坐标系.反假频、自适应偏移孔径计算都在（x'，y'）坐标系中进行.
由于3D响应是一个旋转椭球体，仅仅需要计算该椭球体的1/4就可以了.其余3/4利用对称性充填，这样可以大大提高计算效率.这样的做法与均方根速度的概念吻合，也与时间偏移仅仅适用于水平层状介质的理论一致.在实际处理中，也可以用成像点处的均方根速度来计算每一成像点对应的走时，这样做可以提高成像精度.
为了完整起见，列出所用的3D Kirchhoff积分偏移公式：
|
(1) |
其中，
|
(2) |
其中，hx和hy分别是In-line和Cross-line方向的半偏移距.式（1）表示对一个CMP道集进行绕射叠加.另外，在偏移成像时一般忽略式（1）中的第二项.
3.2 自适应偏移孔径的选择在（x'，y'）规定的旋转坐标系中，对于每一个成像时间τ，在地表处存在如下的椭圆关系：
|
(3) |
其中，h是一个输入道的半炮检距.τ是垂向双程走时.显然，长半轴
相对于输入道的中点坐标，可以得到如下的y'方向偏移成像范围：
|
(4a) |
|
(4b) |
对于上述y'方向偏移成像范围内的每一个固定y'，对应的x'方向偏移孔径计算方法为
|
(5) |
相对于输入道的中点坐标，可以得到如下的x'方向的偏移孔径：
|
(6a) |
|
(6b) |
y'方向偏移成像范围内的所有固定的y'对应的成像线计算完毕，即完成了3D响应椭球体（等时面椭球体）的计算.
3.3 偏移响应算子的截断在（x'，y'）规定的旋转坐标系中，对于由输入炮检点连线为横坐标轴的旋转椭球体，y'=0时，它退化为二维的椭圆；x'=0时，它退化为一个圆.
y'=0时二维椭圆方程为
|
(7) |
x'=0时的圆方程为
|
(8) |
关于算子截断，可以从两个角度来考虑：一个是从绕射时距曲线的角度（等价于用方程（7））考虑，一个是从圆方程（方程（8））的角度考虑.
如从绕射时距曲线的角度考虑，相对于二维椭圆方程（7），有对应的绕射波时距关系：
|
(9) |
当偏移距为零时，绕射波时距关系为
|
(10) |
可以导出零偏移距剖面绕射波收敛的最大偏移孔径：
|
(11) |
当偏移距不为零时，假设对共偏移距道集做动校正，从（9）式可以得到
|
(12) |
其中，t0max是动校正后的最大双程走时.
方程（11）和（12）的使用很不方便.分析（11）和（12）式中τ/t0max的物理意义知道，
|
(13) |
其中，α为反射界面的倾角.
把（13）式代入方程（11）和（12），得到
|
(14) |
从圆方程（8）知道，
|
(15) |
利用（14）和（15）式可以进行算子截断处理，本质上就是利用反射界面的最大角度来控制In-line和Cross-line方向的偏移孔径.可以选取不同的角度来控制In-line和Cross-line两个方向偏移孔径的大小.在一定的地下反射结构情况下，这样做既提高成像质量，又提高计算效率.一般情况下，X方向给定的最大成像角度αxmax大于Y方向的最大成像角度αymax.
3.4 对X方向及Y方向偏移孔径的限制当成像深度较大时（τ较大），由（14）和（15）式计算出的孔径会比较大，而处于较大孔径处的成像点所对应的走时可能大于原始数据的记录时间Tmax，因为偏移时取不到振幅值而浪费了计算，这样的最大孔径就不是有效孔径.因此，可以根据最大记录时间来进一步限制最大孔径.
在旋转坐标系中，仅计算产生1/4个成像椭球体.从原点出发，当成像点离原点越来越远时，偏移孔径越来越大，当成像深度较大时，算出的走时可能会大于观测数据的走时.若某一个深度上计算出的走时大于观测数据的走时，不再进行更深部分的成像，进入下一道（具有更大偏移孔径的成像道）进行成像处理，走时的判断是类似的.这样处理可以节省大量的计算时间.
另外，可以利用成像范围对旋转坐标系中的成像孔径进行限制.一般地，炮检连线方位角不很大，在旋转坐标系中，大部分旋转椭球体等时面都不落在成像输出范围内，这会造成计算量的浪费，尤其是成像输出范围比较小的时候.这也说明该体偏移方案适合大规模数据体的积分法偏移成像.因为In-line方向的成像范围一般比较大，仅对随成像深度增加而增大的Y孔径进行判断，超出成像范围的大部分成像点不进行计算，这可明显提高计算效率.
3.5 旋转坐标系中3D响应椭球体向原始坐标系的投影旋转坐标系（x'，y'）与原始坐标系（x，y）之间的坐标变换关系为
|
(16) |
根据（X'min，X'max）和Y'min（，Y'max）确定的范围，把前面形成的3D响应椭球体通过（16）式映射到（x，y）坐标系中.
旋转角度的计算公式如下：
|
(17) |
必须注意到，从旋转坐标中均匀的成像网格点向原始坐标映射时会导致不同点落到同一个网格点的问题，这样会导致实际的输出网格点振幅分配不合理，同时会导致单道成像结果有空道.这两个问题可以通过以下方法解决.
设在旋转坐标系下计算出的样点值p（x'，y'），在原始坐标系下对应的样点值为p（x，y），现要将其‘发散’到周围4个网格点上.如图 2所示.
|
图 2 振幅分配示意图 Fig. 2 The amplitude is scatted to the surround four points |
设4个网格点为f（xi，yj），f（xi，yj+1），f（xi+1，yj），f（xi+1，yj+1），则应该有：
|
(18) |
其中，
|
coe1，coe2，coe3，coe4分别是将当前已计算得到的样点值‘发散’到周围4个样点值的权系数.诸权系数按照面积比计算：
|
满足：coe1 +coe2 +coe3 +coe4=1.
4 反假频处理每固定一个Y，由计算得到Xmax，由下式给出Pxmax和Pymax：
|
(19) |
|
(20) |
由下式算出射线参数Pmax=Pxmax+Pymax.最后得到反假频的具体公式：
|
(21) |
其中，
|
(22) |
然后用（21）式计算fantialiasing，及实施空变滤波反假频.
5 偏移速度问题的解释理论意义下，一个输入的某个时刻反射子波像在成像空间中对应一个等时面.该等时面的形成可以用该道CMP处的RMS速度计算走时，这样的计算是基于水平层状或缓横向变速假设的.
如图 3所示，在横向变速不大或者反射界面水平的情况下，用P2点的反射时间来近似代替P1和P3的反射时间是没问题的.但在横向变速较大以及成像深度较深时，P1及P3与P2相距甚远，相应的走时差异较大，此时用P2点的反射时间来代替P1和P3的反射时间会使得成像质量下降，这时叠前时间偏移的问题会凸显出来.
|
图 3 Kirchhoff积分PSTM均方根速度应用的解释 Fig. 3 The explanation of Vrms used in Kirchhoff integral PSTM |
实际编写程序时，我们并不能沿着等时面放成像子波.因为等时面对应的成像时间τ是横向变化的，而且很难用（23）式算出：
|
(23) |
因此，程序实现过程中先对成像时间τ循环，然后对横向成像位置循环.这样相当于把该道对应的中点处的RMS速度移到成像点处来计算走时.即：对于每个固定的τ，每个成像点位置的成像均用该道对应的中点处的RMS速度来计算走时.显然，这样的计算是基于水平层状假设的，或缓横向变速假设的.这对陡倾角成像会有影响，但在实际应用上却获得了成功.
改进的方法是在每个成像点处用该点的RMS速度计算走时（也可以考虑引入非对称路径走时计算或各向异性介质情况下走时的计算等），这样做可能会提高走时计算的精度和成像的精度.但是增加了很多的计算量.
6 三维Kirchhoff积分法体偏移的并行I/O方案3D Kirchhoff积分叠前体偏移的并行化计算把共偏移距道集中的偏移距作为一级并行计算索引；把每个共偏移距道集中的线号作为二级并行计算索引.用两级索引给每个处理器分配任务.这样做的目的是让并行粒度比较适中，适应当前计算机集群处理器个数很多的情况.
具体的并行计算采用主从模式分配计算作业.一级并行索引用来分组处理器；二级并行索引用来给处理器分配完成某一条线的数据成像处理任务.
当前实现的并行方案是二级并行计算索引按共偏移距数据体中的地震道号，而不是线号.这样做的优点是并行方案简单，有利于断点恢复.但付出的代价是单个共偏移距数据的单道成像用很多处理器并行计算，每个处理器上的成像结果在归约叠加时需要较大的通讯量，这就要求通讯效率较高，且通讯讯号要稳定.
7 非整样点上的振幅取样问题在Kirchhoff积分偏移成像过程中，由（23）式计算走时会遇到非整样点上的振幅取样问题.可以用SINC插值来计算不等于整样点处的走时对应的振幅.但这样会付出较大的计算代价.因为这步插值是在最内核循环进行的，因此付出的计算代价巨大.
可以用折中的方法，事先把地震数据插值成比较密的时间采样.计算走时如果不等于插值后的整样点处的走时就取最接近的整样点处的振幅值，误差不会很大.事先插值可以用SINC插值来解决，此时在外循环进行，计算量会很小.这是一个比较好的解决方案.
8 Kirchhoff3D体偏移的断点保护断点保护的实现分两种情况：（1）只输出偏移剖面情况下的断点保护.此时设定两个类型相同的临时文件轮换存放中间偏移剖面，临时文件写完后，在日志文件（*.log）中记录当前完成的偏移距和偏移用时等信息.若遇到断点重提作业时，主程序先从日志文件中读取最近完成的偏移距标号等信息，从而根据该偏移距标号把对应的中间偏移剖面数据读入内存中，然后进行下一个偏移距的偏移计算.（2）输出成像道集和偏移剖面情况下的断点保护.此时成像道集文件作为临时文件，当一个偏移距数据体中的所有道做完之后，主节点对各个节点的偏移结果进行归约，并将归约后的数据按照成像道集的格式写到道集文件中，道集写完后在日志文件中记录当前完成的偏移距标号和偏移用时等信息.遇到断点重提作业时，主程序先从日志文件中读取最近完成的偏移距标号等信息，从而进行下一个偏移距的偏移计算.完成所有共偏移距数据体，就得到最终的成像道集文件，对成像道集修饰处理后叠加求和，得到叠前偏移剖面.
以上两种断点保护方案的最小保护单位是一个偏移距，这是结合了具体的计算量和存储量决定的.当然，按照体偏移的实现方案，可以进行更加细致的断点保护.
9 偏移效果试验根据前面所讨论的方法，我们在iCluster叠前地震成像软件系统中实现了3D Kirchhoff积分叠前时间偏移的体偏移方案.主要目的是为各种3D Kirchhoff积分叠前偏移的实用化奠定软件方面的基础.
首先测试的数据是3DSEG岩丘模型.该模型是一个国际通用的测试复杂构造成像和速度模型建立方法的理论地质模型.该模型数据有三套.通常称为数据A、数据B和数据C.数据A为沿AA'线采集的一束线.数据B是沿AA'的一条二维线数据.数据C共有两套：一套宽方位角的，共25束线；一套是窄方位角的，共50束线.
数据A的观测参数为：共138炮，每炮6线排列，每线65道.炮间距80 m，线间距80 m，道间距40m.记录长度5s，采样点数625个，时间采样率8ms.
宽方位角数据的观测参数为：25束线，炮线间距320m.炮间距和线间距都是80m，道间距40m.每束线96炮，每炮65道；记录长度5s，采样点数625个，时间采样率8ms.
窄方位角数据的观测参数为：50束线，炮线间距160m.炮间距80m，线间距40m，道间距40m.每束线96炮，每炮65道；记录长度5s，采样点数625个，时间采样率8ms.
图 4展示的是对宽方角数据的3DKirchhoff积分叠前时间偏移结果.其中，图 4a显示的是体偏移的3D Kirchhoff积分叠前时间偏移的第130线的结果；图 4b是按线偏移方式进行3D Kirchhoff积分叠前时间偏移的第130线的结果.该模型的偏移参数包括：两种偏移方法的成像网格都是20×40m，线偏移In-line方向的最大偏移孔径（全孔径）为7000m，体偏移In-line方向最大偏移孔径为8000m；In-line方向的最大成像角度两者均为75°（相当于算子截断）；最大炮检距两者都是1391.5m；线偏移的Cross-line方向孔径为480m（即每条成像线的输入线为7条线），体偏移的Cross-line孔径是计算出来的，自适应的，其最大成像角度给定为45°；对原始输入道SINC插值加密倍数为4，即将原始的8ms采样加密成2ms.图 5a展示的是体偏移的3D Kirchhoff积分叠前时间偏移的第90线的结果；图 5b是按线偏移方式进行3D Kirchhoff积分叠前时间偏移的第90线的结果，参数同前面所述.红线标识的区域中明显可以看到体偏移的结果好于线偏移的方式.事实上，体偏移和线偏移从理论上是没有差异的，但实现方式差异比较大，线偏移的Cross-line方向的偏移孔径不是按速度和偏移参数自适应选择的，而是事先人为决定的.这是导致线偏移成像质量达不到较佳的主要原因.另外，反假频滤波器设计也难于达到最优.
|
图 4 线130的三维SEG岩丘模型体偏移（a）及线偏移（b）结果 Fig. 4 The volume Kirchhoff PSTM result (a) and the line target Kirchhoff PSTM result of 3D SEG salt model (line 130) |
|
图 5 线90的三维SEG岩丘模型体偏移结果（a）及线偏移结果（b） Fig. 5 The volume Kirchhoff PSTM result of 3D SEG salt model (line 90) (a) and the line target Kirchhoff PSTM result of 3D SEG salt model (line 90) (b) |
图 6和图 7分别展示了3DKirchhoff积分叠前时间偏移成像结果中的两条线.该工区的满覆盖面积54.8km2，总道数594万道，记录时间6s，采样2ms.该工区的特点是有负向背斜构造，包含比较大的陡倾角构造.b图是Omega处理系统的叠前成像结果，a图是iCluster系统的叠前成像结果.该实际数据的偏移参数包括：两者的成像网格都是25× 25m，In-line方向的最大偏移孔径（全孔径）两者均为8000m；In-line方向的最大成像角度均为60°；最大炮检距都是12000m；体偏移的Cross-line孔径是计算出来的，自适应的，其最大成像角度给定为40°（相当于算子截断），Omega处理系统没有此参数.对原始输入道SINC插值加密倍数为2，即将原始的2ms采样加密成1ms.总体而言，两种软件系统的成像效果差不多.红框标识出的部分是二者的差异明显之处.图 6a中红框的背斜界面同相轴更连续. 图 7a中断点更干脆.经处理人员的对比分析，给出了测试结论：叠前时间体偏移程序运行稳定，结果正确.其偏移效果明显好于线偏移效果，与商业软件Omega偏移结果基本一致.测试结果说明按本文所述的方法理论和实现方案给出的iCluster系统中的叠前时间偏移成像模块是正确的和实用的.
|
图 6 某三维实际数据体偏移结果（线600）（a）及三维实际数据Omega偏移某条线结果（线600）（b） Fig. 6 The volume Kirchhoff PSTM result of field data (a) and the Kirchhoff PSTM result of field data by seismic data processing software Omega (b) |
|
图 7 某三维实际数据体偏移结果（线610）（a）及某三维实际数据Omega偏移结果（线610）（b） Fig. 7 The volume Kirchhoff PSTM result of field data (a) and the Kirchhoff PSTM result of field data by seismic data processing software Omega (b) |
关于计算效率，该实际数据大小为66.4G，总共594万道，两种软件的偏移同样利用曙光PC-Cluster机群共160个进程并行偏移，iCluster系统的体偏移耗时12h，处理效率为每进程35.4 M/h；Omega处理系统的Kirchhoff积分偏移耗时22.6h，平均处理效率为每进程18.8M/h，相比之下体偏移的效率要比Omega处理系统高.所用机群环境软硬件参数如表 1，2所示.
|
表 1 计算机群系统软硬件环境 Table 1 Computer environment |
|
|
表 2 两个系统处理相同数据的效率对比 Table 2 The comparison of migration time on the same field data set produced by iCluster and Omega software system |
本文重点探讨大规模地震数据的Kirchhoff积分叠前偏移的计算机实现中存在的问题.为了实现巨大规模数据体的高效率Kirchhoff积分法叠前偏移成像，解决了内存需求量、I/O量及单个输入道Kirchhoff积分偏移成像的精度和效率这三方面的问题.
针对这三个问题，本文设计出的实现方案为：（1）根据计算机最大内存允许量把拟成像深度样点从浅到深分成若干份，每一段所用的成像空间可以驻留在内存中，当前段的体偏移完成后写到本地盘上.清空内存，进行下一拟深度段的体偏移.（2）内存足够大时，拟成像深度段数为1，原始数据的输入次数仅为1次，I/O量达到最小；而段数越多，需要的I/O量越大.（3）每一输入道的按Kirchhoff积分偏移公式的处理是在按炮检连线规定的旋转坐标系中进行的，炮检连线在旋转坐标系中与横轴平行.在这样的坐标系中进行一个地震道的Kirchhoff积分叠前时间偏移比较容易处理偏移孔径的计算、反假频的处理.而且可以只计算1/4个成像椭球体，其余3/4利用对称性充填，提高计算效率.然后，把在旋转坐标系中的成像结果投影回到原始的成像坐标系中.
数值实验验证了本文提出的实现方案的优越性.在此较完善的实现方案的基础上，可以容易地把更优越的积分法偏移方法迅速推向实用化.
单个地震道单独Kirchhoff积分偏移显然不是最高效的办法，在数据域Beam-Ray [17]处理的基础上，我们正把基于Beam-Ray的成像域中[18]的Kirchhoff积分叠前时间偏移，同时把压制采集脚印的方法融合到该实现方案中.相信这会进一步提高Kirchhoff积分偏移的成像质量.
| [1] | Schneider W A. Integral formulation in two and three dimensions. Geophysics , 1978, 43: 49-76. DOI:10.1190/1.1440828 |
| [2] | Bevc D. Imaging complex structure with semi recursive Kirchhoff migration. Geophysics , 1997, 62: 577-588. DOI:10.1190/1.1444167 |
| [3] | Sun J. Limited-aperture migration. Geophysics , 2000, 65: 584-595. DOI:10.1190/1.1444754 |
| [4] | Sun J. On the limited aperture migration in two dimensions. Geophysics , 1998, 63: 984-994. DOI:10.1190/1.1444409 |
| [5] | Sun S, Bancroft J. How much does the migration aperture actually contribute to the migration result? 71st Annual International Meet, SEG, Expanded Abstracts, 2001 |
| [6] | Lumley D, Claerbout J, Bevc D. Anti-aliased Kirchhoff 3D migration. 64th Ann. Internat. Mtg., Soc. Expl. Geophys., Expanded Abstracts, 1994, 1282~1285 |
| [7] | Sun J, Abma R, Bernitsas N. Anti-aliasing in Kirchhoff migration. 69th Ann. Internat. Mtg., Soc. Expl. Geophys., Expanded Abstracts, 1999, 1134~1137 |
| [8] | Yilmaz Ö. Seismic data analysis:Processing, Inversion, and interpretation of seismic data. SEG Investigations in Geophysics, No.10. 2001 |
| [9] | Walter K. Curved-ray time migration can improve seismic imaging. Oil & Gas Journal International Petroleum News and Technology , 2002: 1-8. |
| [10] | Sun C, Martinez R D.Amplitude preserving 3D prestack Kirchhoff time migration for V (Z) and VTI media.72 SEG Internat Mtg Expanded Abstracts MIG4.6.2002 |
| [11] | 刘国峰, 刘洪, 王秀闽, 等. Kirchhoff积分叠前时间偏移的两种走时计算及并行算法. 地球物理学进展 , 2009, 24(1): 131–136. Liu G F, Liu H, Wang X M. Two kinds of traveling time computation and parallel computing method of Kirchhoff migration. Progress in Geophys. (in Chinese) , 2009, 24(1): 131-136. |
| [12] | Canning A, Gerald H F. Gardner. Reducing 3-D acquisition footprint for 3-D DMO and 3-D prestack migration. Geophysics , 1998, 63: 1177-1183. DOI:10.1190/1.1444417 |
| [13] | Canning A, Gardner G H F. Reducing 3-D acquisition footprint for 3-D DMO and 3-D prestack migration, 3-D Seismic Exploration, 71st Ann. Internat. Mtg:Soc. of Expl. Geophys. Expanded Abstracts, 2001 |
| [14] | 程玖兵, 王楠, 马在田. 表驱三维角度域Kirchhoff叠前时间偏移成像方法. 地球物理学报 , 2009, 52(3): 792–800. Cheng J B, Wang N, Ma Z T. Table-driven 3-D angle-domain imaging approach for Kirchhoff prestack time migration. Chinese J. Geophy. (in Chinese) , 2009, 52(3): 792-800. |
| [15] | 邹振, 刘洪, 刘红伟. Kirchhoff叠前时间偏移角度道集. 地球物理学报 , 2010, 53(5): 1207–1214. Zou Z, Liu H, Liu H W. Common angle gathers based on Kirchhoff prestack time migration. Chinese J. Geophy. (in Chinese) , 2010, 53(5): 1207-1214. |
| [16] | 李建江, 舒继武, 王有新, 等. 一种基于共享存贮的叠前深度偏移并行算法. 软件学报 , 2002, 13(12): 2231–2237. Li J J, Shu J W, Wang Y X. A parallel algorithm for prestack depth migration based on shared memory. Journal of Software (in Chinese) , 2002, 13(12): 2231-2237. |
| [17] | 王华忠, 张元巧, 任浩然. 叠前地震数据射线束道集叠加压制噪音. 石油物探 , 2007, 46(6): 356–360. Wang H Z, Zhang Y Q, Ren H R. Prestack noise suppression with ray-beam gather stacking. Geophysical Prospecting for Petroleum (in Chinese) , 2007, 46(6): 356-360. |
| [18] | Sun Yonghe, Qin Fuhao. 3-D prestack Kirchhoff beam migration for depth imaging. Geophysics , 2000, 65: 1592-1603. DOI:10.1190/1.1444847 |