F-K变换与预测反褶积压制多次波效果对比

雍 凡 ， 罗水余 ， 李颜贵 ， 刘建生 ， 刘子龙 ， 蒋正中

(中国地质科学院 地球物理与地球化学勘查研究所，廊坊 065000)

0 前言

多次波的存在影响地震成像的真实性和可靠性，干扰地震资料的解释，是地震资料处理中不可回避的一个重要问题。用于压制多次波的方法可分为两大类:①基于有效波和多次波之间差异的滤波方法；②基于波动方程的预测减去法，通过波动方程模拟波场或反演地震数据来预测多次波，然后把他从原始地震数据中减去的压制方法[1-3]。

根据多次波发生下行反射的位置，可以把多次波分为两类:自由界面多次波和层间多次波。自由界面多次波(Surface-related multiple)，是指至少在自由界面发生一次下行反射所形成的波；层间多次波( Internal-multiple)，是指所有下行反射发生在除自由界面以外的其它反射界面的波。

对于海上数据，与自由界面相关的多次波通常是多次波的主要形式，通常情况下层间多次波很难识别。而对于陆上数据，由于近地表存在未固结层和地表条件的多变性，使得陆上数据的表面多次波不像海上资料那么清楚，层间多次波就成为多次波的重要形式[4-7]。

作者通过对地质模型进行声波时差方程正演的方法取得两组含不同类型多次波的正演数据，然后选取和应用了F-K域压制多次波和预测反褶积两种压制多次波的方法对含有与自由表面相关的多次波和层间多次波的地震资料进行了处理，并对比结果比较各方法的优劣及总结各自的适用条件。

1 正演模型数据

我们建立一个低速层中包含高速层的地质模型，并通过声波有限差分的数值计算方法生成了95炮的模拟地震记录。

建立一个四层水平层状速度模型(图1)，模型大小为 2 000 m×7 000 m，对应的速度和层厚分别为：第一层 1 500 m/s，层厚为500 m；第二层 3 000 m/s，层厚 为500 m；第三层 2 000 m/s，层厚为500 m；第四层 4 000 m/s，层厚为500 m。如此形成了三个波阻抗界面，从浅到深分别记为T1、T2和T3。

炮点从 1 500 m处开始布置，炮间距为40 m，共95炮。检波器炮点开始布设，最小偏移距为0 m，道间距为10 m，共201道接收，单边激发，接收排列随着炮点移动而移动。震源采用零相位的25 Hz雷克子波，记录采样率为1ms，记录长度为 2 200 ms。通过模拟放炮共得到95炮地震记录。

正演时如果将模型表面和边界设置为吸收边界，炮记录中只含有一次波和层间多次波。如果将模型表面设置为自由表面，则炮记录将包含与自由面相关的多次波记录。

如图2(a)所示，在最终的单炮地震记录上除了直达波和三个分界面上产生的三个一次反射波(T1、T2和T3)外，还有一系列的层间多次波。层间多次波由于比一次波至少多经过了二次反射，所以振幅比同时间到达的一次波要小得多，并且以一定的时间间隔重复出现。

图1 正演所用速度模型Fig．1 Velocity model

图2 含有多次波的合成地震记录Fig．2 Synthetized seismic shot gather with multiple(a)仅含层间多次波；(b)包含了与自由界面相关的多次波

如果考虑与自由界面相关的多次波，把地表设定为自由界面。由于自由界面的反射系数为“1”，与自由界面相关的多次波的振幅比层间多次波要大得多[8]，即使经过多次反射，依然在单炮记录上可见其影响(图2(b))。

2 F-K域多次波的压制

2.1 F-K变换的基本原理

在地震信号处理中，傅立叶变换是应用最为广泛的变换之一。它将信号分解为不同频率的指数函数(即正弦和余弦)。对于具有两个变量(如CMP道集的时间和炮检距)的函数来讲，F-K变换就是两个变量的双重傅立叶变换。首先，把地震数据d(x,t)由时间域变换到频率域，即

(1)

然后对空间方向再进行傅立叶变换，由空间域变换到波数域，即

(2)

在方程(1)和方程(2)中，f是频率；kx是水平波数，也可以认为是空间频率。傅立叶变换在实际地震数据处理中总是采用离散形式。一个由宽带子波组成的双曲线同相轴可以看作是由许多不同角度的平面波构成，其在F-K域中扩展为一个三角形区域。

利用两个与上面类似的变换，可以首先将波数-频率域数据变换到频率-空间域

(3)

再由频率-空间域变回时-空域

(4)

F-K域常被用来分析不同的地震波的时间特性和速度特性。

2.2 在F-K域对多次波压制的效果

由于一次波与多次波在视速度上存在差异，可以选择对炮集或CMP道集进行NMO校正，再对其进行二维傅立叶变换，对于视速度小于一次波的与自由界面相关的多次波由于校正不足将投影在波数为负的区域，可以通过对正值波数里的数据进行充零，在反变化到T-X域从而达到压制多次波的目的[9-10]。

图3展示了对图1模型正演得到的包含与自由界面相关多次波的地震数据(图2)进行NMO校正后在F-K域滤波的过程，可以看出在大炮检距处，与自由表面相关的多次波和一次波具有明显的时差，通过F-K滤波后可以很容易消除这部分多次波。然而在小炮检距处，两者之间的时差较小，很难分开。从压制多次波前后的速度谱上可以看到低速的浅层低速的多次波的能量基本被消除。

如图4所示，只含有层间多次波的炮集在进行NMO校正后，层间多次波也被拉平了，甚至在大炮检距端还有动校不足的现象。说明层间多次波的均方根速度与一次波的速度接近，甚至大于深层一次波的速度。所以即使应用滤波器将F-K谱的一半象限充零(图4(d))，也不能起到衰减层间多次波的目的，层间多次波在滤波后的炮集上依然保留了绝大多数的能量(图4(c))。

综上所述，在大炮检距与自由界面相关的多次波与一次波有明显的时差，应用该方法能后有效的去除多次波。而在小炮检距，两者的时差很小，很难区分。另外层间多次波的速度和一次波的速度非常接近，在浅层存在高速层的情况下，层间多次波的均方根速度甚至可能大于同深度的一次波。针对上述两个问题，F-K域滤波的方法很难将一次波和多次波分离。

3 反褶积压制多次波

3.1 预测反褶积的基本原理

在地震勘探中，由于震源爆炸时岩石破坏圈和岩石塑性圈的作用，使得震源发出的尖脉冲到达弹性形变区时变成一个具有一定延续时间的稳定波形b(t)(地震子波(wavelet))。地层对震源脉冲的这种改造作用就相当于一个滤波器，通常称为大地滤波器。通过这个滤波器的作用，子波的高频成分损失，脉冲的频谱变窄，从而使激发时产生的尖脉冲经大地滤波后其延续时间加大。反褶积的作用就是为了去掉大地滤波器的影响，将地震子波b(t)压缩成原来的震源脉冲形式。预测反褶积在某种意义上可以说是一种广义的最小平方反褶积，脉冲反褶积可以看作是它的一个特例，并能用于研究一般的反褶积问题。在地震资料数字处理中，预测反褶积主要是用来压制多次波、海上鸣震等规则干扰波。

图3 F-K域压制自由表面多次波Fig．3 Surface-related multiple attenuation in F-K domain(a)包含多次波的一个经过NMO校正后的CMP道集，部分多次波校正不足；(b) 图(a)的F—K谱；(c) 图(a)经过F—K域滤波后的CMP道集；(d)进行多次波压制后的F-K谱

图4 F-K域压制层间多次波Fig．4 Internal multiple attenuation in F-K domain(a)包含层间多次波的一个经过NMO校正后的CMP道集；(b) 图(a)的F—K谱；(c) 图(a)经过F—K域滤波后的CMP道集；(d)进行多次波压制后的F-K谱

设地震子波b(t)满足最小相位条件，反射系数为白噪声，褶积模型为

(5)

则有t+l时刻的输出值

(6)

对照预测误差公式

(7)

e(t+l)=b′(t)*ξ(t+l)

(8)

即将一个子波的前部与反射系数的褶积就得到了预测误差(一次反射波)。反过来讲，用这种方法可压缩子波长度，提高地震资料的分辨率。

由于预测反褶积后，子波被切成l长，因此预测反褶积实际上是一种子波波形切除反褶积。特别当l=1时，子波变成了δ脉冲，以上预测反褶积实际上变成了脉冲反褶积。

设计一个预测因子

c(t)=(c(0),c(1),…,c(m))

(9)

使

(10)

这一问题显然符合维纳滤波的原理, 易得求解预测因子的矩阵方程组

(11)

方程组中左边rxx为地震记录的自相关，右边是期望输出x(t+l)即一个带时移的地震记录与原记录的自相关。求得c(t)后，预测出多次波来，再从含有多次波的实际地震记录中减去预测出的多次波，就可得到只含有一次反射波的地震记录。

3.2 应用预测反褶积压制多次波

不同于一次波不会以重复的模式出现，多次波具有明显的周期性。即使在单道的情况下，多次波也可以根据其周期性来识别。

图5是含有与自由表面相关多次波的单炮记录及其自相关，从处理前的地震数据(图5(a))中可以看出与自由表面相关多次出现的周期在零炮检距处约为0.6 s，自相关前0.05 s的能量代表了地震子波，小炮检距端0.6 s～0.7 s之间的同相轴代表了多次波的能量。但是随着偏移距的增加，多次波的周期逐渐减小，到偏移距为2 000 m的时候周期变为0.4 s。对该炮集应用预测距离为0.6 s，算子长度为0.8 s的预测反褶积，得到反褶积后的炮集及其对应的自相关(图5(b))。可以看到在小炮检距端，具有较好周期性的多次波得到了很好的压制。而随着炮检距的增加，多次波不再满足该周期性，预测反褶积的压制效果就不再明显。0.6 s～0.7 s之间的周期性能量不再存在，而在大炮检距端依然存在着较强的周期信号。

相对于前面讨论的与自由界面相关的多次波，层间多次波射线路径更为复杂，并且至少经过了三次反射，能量衰减很快。观察单炮记录和其自相关(图6(a)、(c))，可以看到地震子波的能量集中在自相关的前0.05 s内，0.3 s处存在较强的周期信号，并且该信号也有随着偏移距增加周期变小的趋势，但比起与自由表面相关的多次波这个变化趋势要小。图6((b)、(d))为应用0.3 s预测步长和0.8 s算子长度的预测反褶积后的单炮记录及其自相关，可以看到小炮检距的多次波被很好地压制了，而大炮检距依然残留着多次波的能量。

预测反褶积同样可以运用到叠加剖面上，但实际效果可能不能令人满意，主要是因为叠加处理使得各种多次波之间的振幅关系发生了很大的改变。

图5 含有与自由表面相关多次波的单炮记录及其自相关Fig．5 Surface-related multiple attenuation using predictive deconvolution(a)单炮记录；(b)自相关

图6 只含有层间多次波的单炮记录及其自相关处理前后对比图Fig．6 Internal multiple attenuation using predictive deconvolution(a)处理前单炮记录；(b)处理后单炮记录；(c)处理前自相关；(d)处理后自相关

图7 应用预测反褶积压制多次波前、后叠加剖面对比Fig．7 Stacked seismic sections before and after multiple attenuation using predictive deconvolution(a)压制多次波前叠加剖面；(b)压制多次波后叠加剖面

4 结论

1)F-K域多次波压制手段适用于多次波与一次波时差差异较大的情况，对与自由表面相关的多次波的压制效果要好于层间多次波。

2)预测反褶积对小炮检距附近严格满足周期性的多次波去除效果较好，并且可以应用到叠加剖面上。对于本文中采用的数据来说，预测反褶积处理层间多次波的效果要好于与自由表面相关的多次波。

3)F-K域多次波压制和预测反褶积两种方法具有互补性，可同时应用以改善处理效果。

参考文献：

[1] 李鹏.多次波问题的研究进展[J]. 地球物理学进展, 2006(3): 888-897.

[2] 牛滨华, 沈操,黄新武. 波动方程多次波压制技术的进展[J]. 地球物理学进展, 2002(3): 480-485.

[3] 王汝珍.多次波识别与衰减[J]. 勘探地球物理进展, 2003(Z1): 423-432.

[4] 崔凤林, 松辽盆地北部东部深层火山岩处理解释技术及效果[J]. 天然气工业, 2007(S1): 125-127.

[5] 刘洪林.火山机构的反射特征及其地震识别方法[J]. 大庆石油学院学报, 2008(4): 20-22,133.

[6] 韩站一. 火成岩下伏地层地震反射波场特征研究[J]. 石油物探, 2009(3): 226-231,15.

[7] 刘伟明. 基于CFP技术层间多次波压制方法研究[D].长春:吉林大学,2010.

[8] 顾汉明, 王纬,陈国俊. 复杂介质中地震多次反射波快速正演模拟[J]. 地球科学, 2001(5): 541-544.

[9] 梁尚勇, 胡斌与季红军, 苏北盆地复杂断块地区资料处理解释及效果[J]. 石油物探, 2002(3):367-371.

[10] 徐文君, 於文辉,卞爱飞. 地震资料多次波处理[J]. 工程地球物理学报, 2005(6): 417-424.