1 引 言

地震数据的偏移成像是为油气勘探开发提供地下构造图像的主要方法技术，也是地震数据处理中的核心处理步骤^[1]，地震数据的偏移成像方法研究同样也是当前地震数据处理领域中的核心理论问题，尤其是针对地下高陡复杂构造区的地震数据偏移成像.

高陡复杂构造地震数据的偏移成像是当前国内外地震数据偏移成像方法理论研究的热点，主要有三个研究方向^[2-3]:

（1）基于高频近似射线理论的偏移成像方法研究^[4-5]；

（2）基于行波分解的单程波波动方程偏移成像方法研究^[6]；

（3）基于全波方程的逆时偏移成像方法研究^[7-9].对于射线偏移成像方法，由于其所基于的射线理论对于地下复杂构造区速度强横向变化的不适应性，使之在复杂构造区地震数据的偏移成像中的应用受到了一定的限制.近年发展起来的利用高斯束的射线偏移成像方法虽然可在一定程度改善复杂构造区地震数据的偏移成像^[10]，但由于计算的复杂性，使之还处于进一步的完善之中.单程波波动方程偏移方法是近三十年来主要研究发展的复杂构造区地震数据的偏移成像方法^[11-18]，虽然相对于基于射线理论的偏移成像方法存在计算量大的问题，但通过采用微机集群的并行算法已解决了计算效率方面存在的问题.由于当前偏移成像方法中的单程波方程都是通过深度方向的上下行波分解达到的^[19-20]，因此单程波波动方程偏移方法对于复杂构造区地震数据的偏移成像，尤其是高陡复杂构造区地震数据的偏移成像，主要存在大角度构造偏移成像不准确的问题.全波方程的逆时偏移成像方法是被认为是当前解决复杂构造区地震数据偏移成像问题（尤其是高陡复杂构造区地震数据的偏移成像问题）的最有效方法，但其面临的计算量，特别是三维地震数据偏移成像的计算量^[21]，即使采用微机集群的大规模并行计算也会存在计算效率低而导致计算成本极高的问题，这样的问题限制了逆时偏移成像方法在复杂构造区地震数据偏移成像中的应用.近年由Sava^[22]提出的Riemannian坐标系的波场外推与成像方法和由Shan^[23]提出倾斜坐标系的波场外推与成像方法为利用单程波波动方程进行高陡构造成像提供一种策略，即在偏移成像的波场传播中考虑不同的波场优势传播方向，但是他们的这两种方法在偏移成像的具体计算中特别麻烦，涉及到坐标变换和插值.Zhang和McMechan^[24]在回转波的偏移成像中提出了垂向外推与横向外推相结合的偏移成像方法，但他们只研究了二维叠后数据的偏移成像，而没有研究三维叠前数据的偏移成像.Sava 提出的Riemannian坐标系的波场外推与成像方法和Shan提出倾斜坐标系的波场外推与成像方法以及Zhang与McMechan提出的方法在实质上是一致的.

本文针对当前利用深度方向上下行波分解得到的单程波波动方程偏移方法在复杂构造区地震数据的偏移成像中的优劣，再结合逆时偏移成像方法中因为利用了全波方程而能同时考虑深度方向的上下行波和水平方向的前后行波与左右行波以提高高陡复杂构造区地震数据偏移成像效果的思想，提出在单程波波动方程偏移成像中利用波场垂向外推与水平方向外推相结合的手段解决高陡复杂构造地震数据的偏移成像问题，即利用常规的上下行波单程波波动方程偏移成像方法理论解决复杂构造的中低角度构造的偏移成像问题，再利用水平方向行波单程波波动方程偏移成像方法理论解决复杂构造的中高角度构造的偏移成像问题.把本文提出的基于波场垂向和水平方向外推的高陡构造偏移方法应用于模型数据的偏移成像，取得了满意的结果.

2 方法原理 2.1 波场的方向分解

震源激发的地震波和散射体产生的散射波是向着空间的各个方向传播的，只是对于反射地震勘探，由于观测系统中检波器的排列长度有限，通常把深度方向（垂向）视为地震波的优势传播方向，这就导致了我们常用的上下行波分解.地震数据中水平和近水平方向传播的地震波比较少，但它们对改善高陡复杂构造的偏移成像十分重要，因为它们通常是高陡复杂构造产生的散射地震波.我们知道对于式（1）所表述的三维波动方程，x、y和z三个方向的坐标是对等的，是没有主次之分的.

(1)

采用垂直向下为深度z方向的右手螺旋坐标系，对式（1）进行x、y、z和t的四维Fourier变换，得到波动方程的频散关系，即

(2)

假定垂向为地震波场的优势传播方向，则得到常见的上下行波分解，即

(3)

如果假定水平方向为地震波场的优势传播方向，则可得到x方向的前后行波分解和y方向的左右行波分解，即

(4)

(5)

由分解式（3）可得到上下行波的频散关系，即

(6)

由分解式（4）可得到前后行波的频散关系，即

(7)

由分解式（5）可得到左右行波的频散关系，即

(8)

2.2 波场的垂向外推和水平方向外推

由上面得到的上下、前后和左右行波频散关系，我们可以写出上下行波的单程波传播方程，

(9)

前后行波的单程波传播方程，

(10)

左右行波的单程波传播方程，

(11)

利用上述行波的单程波传播方程（9）、（10）和（11），我们就可以进行波场的垂向外推和水平方向外推.根据文献[25]提出的广义屏波场传播算子，我们可写出波场在z方向上的递进外推计算式，

(12)

式中，p（x，y，z_i；ω）和p（x，y，z_i+Δz；ω）分别表示z_i上的已知波场和z_i+Δz上的外推波场；Δz为外推步长；F_x，y和F_Kx，Ky^-1分别表示x和y方向上的正反Fourier变换；v₀（z_i）代表z_i上的背景速度；ε（x，y，z_i）为z_i上相对于v₀（z_i）的速度扰动；k_z（z_i）表示z_i上z方向的波数；α 代表下面表示式，即

仿照式（12），我们可写出波场在x方向上的递进外推计算式，即

(13)

式中，p（y，z，x_i；ω）和p（y，z，x_i+Δx；ω）分别表示x_i上的已知波场和x_i+Δx上的外推波场；Δx为外推步长；F_y，z和F_Ky，Kz^-1分别表示y和z方向上的正反Fourier变换；v₀（x_i）代表x_i上的背景速度；ε（y，z，x_i）为x_i上相对于v₀（x_i）的速度扰动；K_x（x_i）表示x_i上x方向的波数；α 代表下面表示式，即

同样也可写出波场在y方向上的递进外推计算式，即

(14)

式中，p（z，x，y_i；ω）和p（z，x，y_i+Δy；ω）分别表示y_i上的已知波场和y_i+Δy上的外推波场；Δy为外推步长；F_z，x和F_Kz，Kx^-1分别表示z和x方向上的正反Fourier变换；v₀（y_i）代表y_i上的背景速度；ε（y，z，y_i）为y_i上相对于v₀（y_i）的速度扰动；K_y（y_i）表示y_i上y方向的波数；α 代表下面表示式，即

2.3 波场垂向外推与水平方向外推的联合成像

利用上述的x、y和z三个方向的波场递进外推算子，我们就可以进行共炮道集地震数据的基于波场垂向外推与水平方向外推的联合成像.z方向的波场外推与成像也就是常规的单程波波动方程叠前深度偏移成像，x方向的波场外推与成像计算过程和y方向的波场外推与成像计算过程与深度方向外推的叠前深度偏移成像是十分相似的.但在波场的x方向和y方向外推过程中，对于共炮道集中各记录道波场的外推可采用逐步逐道累加的方式进行，这是与z方向的波场外推与成像所不同的.

假定共炮道集地震数据经z方向的波场外推和成像得到的偏移成像结果为I_z（x，y，z），经x方向的波场外推和成像得到的偏移成像结果为I_x（x，y，z），经y方向的波场外推和成像得到的偏移成像结果为I_y（x，y，z），则波场垂向外推与水平方向外推相结合的偏移成像结果I（x，y，z），有

(15)

式中，β 和λ 为加权因子.加权因子的取值与地震数据的信噪比以及地下构造倾角有关，但是地下构造倾角在一般的偏移成像中是难以知道的，因此加权因子的取值主要考虑数据的信噪比.对于高信噪比的地震数据，β 和λ 可取为1，即各个方向偏移成像结果的等权相加.由于本文方法的主要目的是解决高陡构造的构造成像问题，因此没有考虑加权因子对偏移成像结果保幅性的影响.

对于二维地震数据的偏移成像，有

(16)

波场垂向外推与水平方向外推相结合的偏移成像方法的计算量相对常规的单程波叠前深度偏移成像方法有一定的增加，对于二维地震数据的偏移成像约增加2倍计算量，对三维地震数据的偏移成像约增加4倍计算量.

3 数值试验

为验证上节提出的高陡构造偏移成像方法的正确性和有效性，我们分别进行了二维脉冲响应试验、SEG-EAGE 二维盐丘模型试验和一个二维高陡构造成像试验.

3.1 脉冲响应试验

在脉冲响应试验中，我们分别作了均匀速度模型的脉冲响应试验和非均匀盐丘体速度模型的脉冲响应试验.图 1和图 2 分别为均匀速度模型和非均匀速度模型的脉冲响应试验结果，图 2 中的背景为含有盐丘的非均匀速度模型.在本试验中取β=1.

图 1a和2a为波场垂向外推的结果，图 1b 和2b为波场水平外推的结果，图 1c和2c为波场垂向外推与水平方向外推相结合的结果.由于受单程波方程不能描述大角度传播的波场的影响，图 1a、1b和图 2a、2b所示的脉冲响应结果都缺失了大角度，这和理论是一致的.图 1c和2c所示的脉冲响应结果没有出现大角度的缺失现象，也就是说明利用波场垂向外推与水平方向外推相结合是有可能实现对90°构造的成像.

3.2 SEG-EAGE 二维盐丘模型试验

SEG-EAGE 二维盐丘模型的观测系统为右边放炮，左边接收，共有325炮，每炮有176道，炮间距48m，道间距24m，最小偏移距为零，记录长度5s，时间采样率8 ms.SEG-EAGE 二维盐丘模型具有速度横向变化剧烈和盐体侧翼倾角较陡的特点.

图 3a为SEG-EAGE 二维盐丘模型的常规叠前深度偏移成像结果，即仅考虑波场垂向外推而得到的偏移成像结果.图 3b为应用本文提出的波场垂向外推与水平方向外推相结合的偏移成像方法得到的偏移成像结果.对比图 3a和3b可看出，应用本文提出的偏移成像方法可以使模型上部陡倾角的盐体侧翼得到清楚的成像，见图中的圆圈所示区域.由于受观测系统中检波器排列长度的限制，图 3b中深部陡倾角构造的偏移成像结果相对于图 3a的结果没有得到改进，这是因为这些中深部陡倾角构造的大角度反射波场需要有大偏移距的检波器才能观测到.

3.3 高陡构造模型试验

高陡构造模型的观测系统为中间放炮，两边接收，共有110炮，每炮有241道，炮间距100m，道间距50m，最小偏移距为-6000 m，记录长度3.5s，时间采样率1ms.图 4a所示为高陡构造模型图.左边高速岩体上部的倾角最高达90°.

图 4b为高陡构造模型的常规叠前深度偏移成像结果，即仅考虑波场垂向外推而得到的偏移成像结果.图 4c为应用本文提出的波场垂向外推与水平方向外推相结合的偏移成像方法得到的偏移成像结果.对比图 4b和4c可看出，应用本文提出的偏移成像方法可以使模型上部岩体侧翼的90°倾角构造得到清楚的成像，见图中的圆圈所示区域.由于受观测系统中检波器排列长度的限制和高速岩体中深部侧翼凹向构造的特殊性，使得图 4c岩体中深部陡倾角构造的偏移成像结果相对于图 4b的结果没有得到明显改进，因为中深部侧翼凹向陡倾角构造的大角度反射波场需要有更大偏移距的检波器才能观测到.

4 结 论

不同于能同时考虑地震数据中各个传播方向成分的全波方程逆时偏移成像方法，本文提出的垂向波场外推与水平方向波场外推相结合的高陡构造偏移成像方法，是一种分步（序贯）考虑地震数据中各个传播方向成分的单程波波动方程偏移成像方法，因此相对于逆时偏移成像方法具有极高的计算效率优势.本文提出的高陡构造偏移成像方法也是目前常规的波动方程叠前深度偏移方法的补充和完善.