CN102169189B - 深水层间多次波消除方法 - Google Patents

Abstract

本发明涉及深水层间多次波消除方法，应用于油气勘探地震资料处理中，它包括1、消除自由表面多次波；2、对原始数据抽道集；3、将抽道集得到的数据(x_s，x_g；t)进行二维傅立叶变换得到频率-波数域数据(k_s，k_g；ω)；4、利用逆散射级数和频率-波数域数据D(k_s，k_g；ω)推导出n阶散射波场D_n；5、根据D_n得到消除与自由表面相关多次波后的波场D₀；6、依据散射波场形成层间多次波的两个约束条件，及下行散射点的深度Z_g比相邻的两个散射点Z_g，和Z_g，迭代计算出V_n，得到有关层间多次波的级数项；7、将逆散射级数项进行反傅立叶变换，得到d₁(x_s，x_g；t)；8、将d(x_s，x_g；t)与d₁(x_s，x_g；t)数据匹配相减，得到消除层间多次波以后的CMP道集。本发明基于量子力学的散射理论，提出深水层间多次波消除技术，填补了地震资料深水层间多次波消除的空白。

Description

深水层间多次波消除方法

技术领域

本发明涉及油气勘探领域，特别是关于海洋油气勘探地震采集中深水层间多次波的消除方法。

背景技术

在油气勘探地震采集中，采集的地震纪录含有反射波和多次波。目前的地震处理技术都是基于反射波的成像理论，多次波被都要被看成噪声消除掉，而消除多次波，同时不能破坏有效波是多次波消除技术的困难所在。多次波可以分成与表面相关的多次波和层间多次波，现在的油气勘探地震资料处理理论的多次波消除技术只能消除与表面相关的多次波，无法消除层间多次波，尤其无法消除深水层间多次波。

层间多次波的定义是在自由表面以下的任何一个界面，至少发生一次下行反射的波，下行反射的次数称为层间多次波的阶数。在图1中，一阶层间多次波发生了一次下行反射，二阶层间多次波发生了两次下行反射。此外我们还看出产生下行反射的条件是：下行反射的散射点深度要低于相邻两个散射点的深度。此外，从图2中可以推断出：反射波和折射波的形成同散射点的个数没有关系，即从一个到多个散射点均可以形成反射波和折射波；而层间多次波同散射点的个数有关，至少需要三个散射点才能形成层间多次波。图3中给出三个散射点不同空间位置组合形成的波，来分析空间位置组合对层间多次波形成的重要性，可以得出生成层间多次波的两个限定条件：第一，在形成层间多次波的所有散射点中，层间多次波发生在下行反射的散射点的深度至少要比相邻两个散射点浅；第二，层间多次波只要要三阶散射才能形成；最终成功将散射体的介质属性和层间多次波的形成联系起来，从而能够准测预测出散射体扰动波场产生的层间多次波，如图5～8所示，图5是消除层间多次波前模型道集剖面，该CMP道集是在一个三层模型条件生成，其中中间层为高速层，可以较为明显的分辨出层间多次波；图6是消除层间多次波后模型道集剖面，该剖面是在利用深水层间多次波消除方法以后得到的CMP剖面，可以清楚的发现层间多次波得到的明显的消除，较好的提升了CMP道集的分辨率；图7是消除层间多次波前地震叠加剖面，该剖面是实际地层中发生速度突变产生了较为明显的层间多次波典型剖面；图8是消除地震层间多次波后地震剖面，针对图7剖面中较为明显的层间多次波，利用深水层间多次波消除方法，对每一个CMP道集消除层间多次波以后叠加后得到的消除层间多次波以后的地震剖面，通过剖面可以发现，该剖面中去浅层地层相似的层间多次波得到了有效压制，地震资料的真实可靠性得到了明显提升。

从方法研究来看，层间多次波的特殊性在于它同反射波的不可分离性，在不同的域中其可分离性均很低，采用常规的滤波类方法，无法实现层间多次波的消除。目前，层间多次波的预测和衰减是多次波研究领域的一个重点、难点和热点。

从实际应用的角度看，在对地震波叠加或者偏移处理后，层间多次波不仅会形成虚假界面同相轴，而且也会同真实界面同相轴重合，影响到真实反射率的恢复，这样会导致错误的地震解释，甚至会导致勘探的失败。特别是在勘探成本比较高的区域，层间多次波的消除具有实际意义且是提高勘探成功率的重要技术手段之一。

在国外，一些地球物理学家在从事这方面的研究工作，最有影响力的是荷兰Delft大学的DELPHI小组和Houston大学的MOSRP小组。目前这两个小组在层间多次波预测衰减方面做了前沿性和开创性的研究工作，其中，以MOSRP小组提出的点散射物理模型方法较好，有不依赖速度模型的优点。特别是Houston大学MOSRP小组采用的点散射模型，利用量子力学中逆散射级数来进行反演研究，提出了用其来构建层间多次波，并取得了很好的效果。MOSRP小组提出的算法能够很好地预测层间多次波，该算法不仅能一次预测一个层位相关的层间多次波，还能预测所有层位相关的层间多次波，较其它预测算法优越。但是该算法计算量非常大，目前难以在实际地震资料处理中得以应用。

在国内，基于散射理论消除层间多次波这方面的技术还是空白。

因而无论从理论研究角度，还是从实际应用角度，无论是陆地地震勘探，还是海洋地震勘探，基于散射理论来研究层间多次波均具有十分重要的意义。

发明内容

针对上述问题，本发明的目的是提供一种在油气勘探地震采集中消除深水层间多次波的方法。

为实现上述目的，本发明采取以下技术方案：一种深水层间多次波消除方法，其特征在于，它包括如下步骤：

1)运用常规方法将原始地震数据解编、去噪、消除自由表面多次波；

2)运用常规方法建立观测系统，对原始数据抽道集得到共中心点道集；

3)将共中心点道集数据进行二维傅立叶变换得到频率-波数域数据 $\stackrel{\‾}{D}(k_s,k_g;\mathrm{\ω});$

4)建立微分方程：

${\tilde{L}}_0{G}_0=-\mathrm{\δ}(\stackrel{\→}{r}-{\stackrel{\→}{r}}_s)---\left(1\right),$

$\tilde{L}G=-\mathrm{\δ}(\stackrel{\→}{r}-{\stackrel{\→}{r}}_s)$

其中是参考及真实波场的微分算子，是波场点和源点的位置，设扰动算子为散射场算子为则有

$\tilde{V}=\tilde{L}-{\tilde{L}}_0,$

上式写为散射基本方程为

(2)式进一步展成逆散射级数形式，

如果把写成级数，即为

$\tilde{V}=V_1+V_2+V_3+\·\·\·={\mathrm{\Σ}}_{n=1}^{\∞}{V}_n---\left(4\right),$

将(4)式带到(3)式中并展开，取同阶项相等，则

一阶项为

二阶项为

三阶项为

n阶项为

G₀可分解为与自由表面有关的项和与自由表面无关的项自由表面多次波由项产生，而一次反射波和层间多次波由项产生，将代入(3)式并提取与有关的项，则得到

一阶项为

$D_1=G_0^d{V}_1{G}_0^d---\left(9\right),$

二阶项为

$D_2=G_0^d{V}_2{G}_0^d-G_0^d{V}_1{G}_0^f{V}_1{G}_0^d---\left(10\right),$

三阶项为

$D_3=G_0^d{V}_3{G}_0^d-G_0^d{V}_1{G}_0^f{V}_2{G}_0^d---\left(11\right),$

N阶项为

$D_n=G_0^d{V}_n{G}_0^d-G_0^d{V}_1{G}_0^f{V}_{n-1}{G}_0^d---\left(12\right);$

5)从(9)式解出V₁，代入(12)式则可解出V₂和D₂，再将V₁和V₂代入(11)式解出V₃和D₃，如此可递归直至解出D_n，则得到消除与自由表面相关多次波后的波场D₀为

$D_0={\mathrm{\Σ}}_{n=1}^{\∞}{G}_0^d{V}_n{G}_0^d={\mathrm{\Σ}}_{n=1}^{\∞}{D}_n---\left(13\right);$

6)假设源点位置为(x_s，z_s)，接收点位置为(x_g，z_g)，则对于参考介质任意点(x，z)的波场可通过(1)式解出，即(1)式可写为

$(\frac{{\▿}^2}{{\mathrm{\ρ}}_0}+\frac{{\mathrm{\ω}}^2}{{\mathrm{\κ}}_0})G_0(x,z,x_s,z_s;\mathrm{\ω})=-\mathrm{\δ}(x-x_s)[\mathrm{\δ}(z-z_s)-\mathrm{\δ}(z+z_s)]---\left(14\right),$

其中，ρ₀为参考介质的密度，κ₀为参考介质的体积模量，

(14)式对x作傅立叶变换，则有

$(\frac{1}{{\mathrm{\ρ}}_0}\frac{d^2}{{\mathrm{dz}}^2}+\frac{q^2}{{\mathrm{\ρ}}_0})G_0(k_x,z,x_s,z_s;\mathrm{\ω})=-\frac{1}{\sqrt{2\mathrm{\π}}}{e}^{-{\mathrm{ik}}_x{z}_x}[\mathrm{\δ}(z-z_s)-\mathrm{\δ}(z+z_s)]---\left(15\right)$

(15)式的通解为

$G_0(k_x,z,x_s,z_s;\mathrm{\ω})=\frac{{\mathrm{\ρ}}_0}{\sqrt{2\mathrm{\π}}}\frac{e^{-{\mathrm{ik}}_x{z}_x}}{-2\mathrm{iq}}[e^{\mathrm{iq}|z-z_s|}-e^{\mathrm{iq}|z+z_s|}]---\left(16\right)$

其中， $q=\mathrm{sgn}\left(\mathrm{\ω}\right)\sqrt{{(\mathrm{\ω}/c_0)}^2-k_x^2},$ 为垂直波数， $c_0=\sqrt{{\mathrm{\κ}}_0/c_0},$ 为海水的纵波速度，因为 $G_0=G_0^d+G_0^f,$ 则

$G_0^d(k_x,z,x_s,z_s;\mathrm{\ω})=\frac{{\mathrm{\ρ}}_0}{\sqrt{2\mathrm{\π}}}\frac{e^{-{\mathrm{ik}}_x{z}_x}}{-2\mathrm{iq}}{e}^{\mathrm{iq}|z-z_s|}---\left(17\right)$

将(17)式代入(9)式解出V1：

$V_1(k_s,k_g;\mathrm{\ω})={\left[e^{i{q}_g{z}_g}{e}^{{\mathrm{iq}}_s{z}_s}\right]}^{-1}{q}_s{q}_g{D}_1(k_s,q_s,k_g,q_g;\mathrm{\ω})---\left(18\right)$

其中， $q_g=\mathrm{sgn}\left(\mathrm{\ω}\right)\sqrt{{(\mathrm{\ω}/c_0)}^2-k_g^2},$ $q_s=\mathrm{sgn}\left(\mathrm{\ω}\right)\sqrt{{(\mathrm{\ω}/c_0)}^2-k_s^2}$

将(17)式代入(13)式解出

${\stackrel{\‾}{D}}_n(k_s,k_g;\mathrm{\ω})=\frac{1}{\mathrm{i\π}{\mathrm{\ρ}}_{0}S\left(\mathrm{\ω}\right)}{\∫}_{-\∞}^{\∞}{e}^{\mathrm{iq}(z_g+z_s)}{\stackrel{\‾}{D}}_1(k_s,k;\mathrm{\ω}){\stackrel{\‾}{D}}_{n-1}(k,k_g;\mathrm{\ω})\mathrm{dk},(n=\mathrm{2,3,4},\·\·\·)---\left(19\right)$

对于(19)式，当n＝2时有

${\stackrel{\‾}{D}}_2(k_s,k_g;\mathrm{\ω})=\frac{1}{\mathrm{i\π}{\mathrm{\ρ}}_{0}S\left(\mathrm{\ω}\right)}{\∫}_{-\∞}^{\∞}{e}^{\mathrm{iq}(z_g+z_s)}{\stackrel{\‾}{D}}_1(k_s,k;\mathrm{\ω}){\stackrel{\‾}{D}}_1(k,k_g;\mathrm{\ω})\mathrm{dk}---\left(20\right)$

即为地震数据的频率-波数域形式，

由公式(14)～公式(16)，依据散射波场形成层间多次波的两个约束条件，及下行散射点的深度Z_g比相邻的两个散射点Z_g’和Z_g”小，采用公式(12)进行迭代计算出V_n，即得到有关层间多次波的级数项；

7)将得到的有关层间多次波的逆散射级数项进行反傅立叶变换，得到空间时间域的层间多次波并将多次波成像；

8)运用常规方法将原始地震数据与预测出来的层间多次波数据匹配相减，从而得到消除层间多次波以后的CMP道集。

本发明由于采取以上技术方案，其具有以下优点：1、现在的油气勘探地震资料处理中，多次波消除技术只能消除与表面相关的多次波，无法消除层间多次波，尤其无法消除深水层间多次波。本发明基于量子力学的散射理论，提出深水层间多次波消除技术，填补了地震资料深水层间多次波消除的空白。2、我们根据Weglein提出的算法，对其算法做了进一步的推导，在水平层状介质的假设前提下(即在公式(20)中假设z_g，z_s相等)，将其积分路径做了修改，并将其频率域计算方法变换到时间域进行计算，大大提高了计算效率，其应用的效果非常明显，不仅仅能够提高地震资料速度分析的精度，也能够对其进行衰减，得到更好的叠加剖面。通过对南海深部地震数据资料的实际处理，得到了较好的效果。

附图说明

图1是地震数据预处理图；

图2是正演反射波和层间多次波散射模型，形成层间多次波的图形；

图3是三个散射点不同位置形成的波形图；

图4是将抽道集得到的CMP数据d(x_s，x_g；t)进行二维傅立叶变换得到频率-波数域数据D(k_s，k_g；ω)；

图5是消除层间多次波前模型道集剖面；

图6是消除层间多次波后模型道集剖面；

图7是消除层间多次波前地震叠加剖面；

图8是消除地震层间多次波后地震剖面。

具体实施方式

下面结合附图和实施例对本发明的进行详细的描述。

散射理论是扰动分析的一种。广义上，它描述了介质的扰动同波场扰动之间的关系。散射正演也被称为正向散射，是指根据参考介质和参考波场以及扰动算子，最后计算出实际介质的波场。散射反演也被称为逆散射，是指根据参考介质、参考波场和实际波场求取实际介质相对于参考介质的扰动，即散射体的介质属性。

我们根据地震波传播的弹性波动方程，将地震波传播的介质视为均匀介质和扰动介质的叠加，根据量子力学的散射理论，导出基于均匀介质的背景场和基于扰动介质的散射场，而总的波场是背景场和散射场之和，散射场可以用扰动介质和背景场的级数来表示。而散射级数本身也即是与散射体介质属性有关的数学表达式。

对于海洋地震勘探，比较方便的方法就是假设边界在自由表面的半空间的水为参考介质，在自由表面下的水中的源产生的对于参考介质的参考波场的格林函数为G₀，而实际接收到的真实波场的格林函数为G，对于压力波场在弹性介质中传播，它们满足下面的微分方程：

${\tilde{L}}_0{G}_0=-\mathrm{\δ}(\stackrel{\→}{r}-{\stackrel{\→}{r}}_s),---\left(1\right)$

$\tilde{L}G=-\mathrm{\δ}(\stackrel{\→}{r}-{\stackrel{\→}{r}}_s)$

其中是参考及真实波场的微分算子，是波场点和源点的位置。设扰动算子为散射场算子为则有

$\tilde{V}=\tilde{L}-{\tilde{L}}_0,$

上式可以写为散射基本方程(李普曼-施温格方程)为

上式可以进一步展成逆散射级数形式，

如果把    也写成级数，即

$\tilde{V}=V_1+V_2+V_3+\·\·\·={\mathrm{\Σ}}_{n=1}^{\∞}{V}_n---\left(4\right)$

将(4)式带到(3)式中并展开，取同阶项相等，则

一阶项为

二阶项为

三阶项为

n阶项为

而事实上G₀可分解为与自由表面有关的项和与自由表面无关的项即则自由表面多次波由项产生，而一次反射波和层间多次波由项产生。将代入(3)式并提取与有关的项，则得到

一阶项为

$D_1=G_0^d{V}_1{G}_0^d---\left(9\right)$

二阶项为

$D_2=G_0^d{V}_2{G}_0^d-G_0^d{V}_1{G}_0^f{V}_1{G}_0^d---\left(10\right)$

三阶项为

$D_3=G_0^d{V}_3{G}_0^d-G_0^d{V}_1{G}_0^f{V}_3{G}_0^d---\left(11\right)$

n阶项为

$D_n=G_0^d{V}_n{G}_0^d-G_0^d{V}_1{G}_0^f{V}_{n-1}{G}_0^d---\left(12\right)$

设散射场在接收面上记录到的波场为D，则波场D经过直达波(从震源点不经过反射和折射直接传播到检波器的波)切除和鬼波(经由海水空气界面反射而被检波器接收到的波)消除后得到波场即为D₁，而相应的D₂为一阶自由表面相关多次波，D₃为二阶自由表面相关多次波，D_n为n-1阶自由表面相关多次波，则消除与自由表面相关多次波后的波场D₀为

$D_0={\mathrm{\Σ}}_{n=1}^{\∞}{G}_0^d{V}_n{G}_0^d={\mathrm{\Σ}}_{n=1}^{\∞}{D}_n---\left(13\right)$

从(9)式解出V₁，代入(10)式则可解出V₂和D₂，在将V₁和V₂代入(11)式解出V₃和D₃，如此可递归直至解出D_n，这样就可通过(9)～(12)式的计算达到消除与自由表面相关多次波的目的。

通过(1)～(13)式我们就可以得到消除与自由表面相关多次波后的波场D₀，为了求解D₀，现假设源点位置为(x_s，z_s)，接收点位置为(x_g，z_g)，则对于参考介质任意点(x，z)的波场可通过(1)式解出，即(1)式可写为

$(\frac{{\▿}^2}{{\mathrm{\ρ}}_0}+\frac{{\mathrm{\ω}}^2}{{\mathrm{\κ}}_0})G_0(x,z,x_s,z_s;\mathrm{\ω})=-\mathrm{\δ}(x-x_s)[\mathrm{\δ}(z-z_s)-\mathrm{\δ}(z+z_s)]---\left(14\right)$

其中，ρ₀为参考介质的密度，κ₀为参考介质的体积模量。

(14)式对x作傅立叶变换，则

$(\frac{1}{{\mathrm{\ρ}}_0}\frac{d^2}{{\mathrm{dz}}^2}+\frac{q^2}{{\mathrm{\ρ}}_0})G_0(k_x,z,x_s,z_s;\mathrm{\ω})=-\frac{1}{\sqrt{2\mathrm{\π}}}{e}^{-{\mathrm{ik}}_x{z}_x}[\mathrm{\δ}(z-z_s)-\mathrm{\δ}(z+z_s)]---\left(15\right)$

(15)式的通解为

$G_0(k_x,z,x_s,z_s;\mathrm{\ω})=\frac{{\mathrm{\ρ}}_0}{\sqrt{2\mathrm{\π}}}\frac{e^{-{\mathrm{ik}}_x{z}_x}}{-2\mathrm{iq}}[e^{\mathrm{iq}|z-z_s|}-e^{\mathrm{iq}|z+z_s|}]---\left(16\right)$

其中， $q=\mathrm{sgn}\left(\mathrm{\ω}\right)\sqrt{{(\mathrm{\ω}/c_0)}^2-k_x^2},$ 为垂直波数， $c_0=\sqrt{{\mathrm{\κ}}_0/c_0},$ 为海水的纵波速度。而 $G_0=G_0^d+G_0^f,$ 则

$G_0^d(k_x,z,x_s,z_s;\mathrm{\ω})=\frac{{\mathrm{\ρ}}_0}{\sqrt{2\mathrm{\π}}}\frac{e^{-{\mathrm{ik}}_x{z}_x}}{-2\mathrm{iq}}{e}^{\mathrm{iq}|z-z_s|}---\left(17\right)$

将(17)式代入(9)式就可解出V₁

$V_1(k_s,k_g;\mathrm{\ω})={\left[e^{i{q}_g{z}_g}{e}^{{\mathrm{iq}}_s{z}_s}\right]}^{-1}{q}_s{q}_g{D}_1(k_s,q_s,k_g,q_g;\mathrm{\ω})---\left(18\right)$

其中， $q_g=\mathrm{sgn}\left(\mathrm{\ω}\right)\sqrt{{(\mathrm{\ω}/c_0)}^2-k_g^2},$ $q_s=\mathrm{sgn}\left(\mathrm{\ω}\right)\sqrt{{(\mathrm{\ω}/c_0)}^2-k_s^2}$

将(17)式代入(13)式解出

${\stackrel{\‾}{D}}_n(k_s,k_g;\mathrm{\ω})=\frac{1}{\mathrm{i\π}{\mathrm{\ρ}}_{0}S\left(\mathrm{\ω}\right)}{\∫}_{-\∞}^{\∞}{e}^{\mathrm{iq}(z_g+z_s)}{\stackrel{\‾}{D}}_1(k_s,k;\mathrm{\ω}){\stackrel{\‾}{D}}_{n-1}(k,k_g;\mathrm{\ω})\mathrm{dk},(n=\mathrm{2,3,4},\·\·\·)---\left(19\right)$

对于(19)式，当n＝2时有

${\stackrel{\‾}{D}}_2(k_s,k_g;\mathrm{\ω})=\frac{1}{\mathrm{i\π}{\mathrm{\ρ}}_{0}S\left(\mathrm{\ω}\right)}{\∫}_{-\∞}^{\∞}{e}^{\mathrm{iq}(z_g+z_s)}{\stackrel{\‾}{D}}_1(k_s,k;\mathrm{\ω}){\stackrel{\‾}{D}}_1(k,k_g;\mathrm{\ω})\mathrm{dk}---\left(20\right)$

而即为原始地震数据的频率-波数域形式。

最后通过反Fourier变换得到原始地震数据，从而达到将散射体的介质属性与层间多次波和真实地震数据联系起来的目的，也就是我们能够用原始地震数据资料，不依靠速度模型、不借助其他速度相关资料能够准确预测层间多次波的目的。

上述讲述的实质是准确预测层间多次波的过程，当层间多次波被预测到后，还要进一步消除掉它。因此作为消除层间多次波的一个完整的实施方案，它具体应该包括如下步骤：

1、运用常规方法将原始地震数据解编、去噪、消除自由表面多次波；

2、运用常规方法建立观测系统，对原始数据抽道集得到共中心点(CMP)道集；

3、将抽道集得到的CMP数据进行二维傅立叶变换得到频率-波数域数据 $\stackrel{\‾}{D}(k_s,k_g;\mathrm{\ω}),$ 见图4；

4、依据公式(3)和得到的频率波数域数据D(k_s，k_g；ω)，即可推导出n阶散射波场D_n；

5、从(9)式解出V₁，代入(12)式则可解出V₂和D₂，再将V₁和V₂代入(11)式解出V₃和D₃，如此可递归直至解出D_n，则得到消除与自由表面相关多次波后的波场D₀为

$D_0={\mathrm{\Σ}}_{n=1}^{\∞}{G}_0^d{V}_n{G}_0^d={\mathrm{\Σ}}_{n=1}^{\∞}{D}_n---\left(13\right);$

6、由公式(14)～公式(16)，依据散射波场形成层间多次波的两个约束条件，及下行散射点的深度Z_g比相邻的两个散射点Z_g’和Z_g”小，采用公式(12)进行迭代(由图2、3可以见，地震波扰动可以同时生成反射波和层间多次波，并且只有经过三阶以上的散射才可能生成层间多次波，从而得出了散射形成层间多次波的第一个条件)，即可计算出V_n，即得到有关层间多次波的级数项；

7、将得到的有关层间多次波的逆散射级数项进行反Fourier变换，最后得到空间时间域的层间多次波d₁(x_s，x_g；t)，并将多次波成像；

8、运用常规方法将原始地震数据d(x_s，x_g；t)与预测出来的层间多次波d₁(x_s，x_g；t)数据匹配相减，从而得到消除层间多次波以后的CMP道集。

Claims (1)

1.一种深水层间多次波消除方法，其特征在于，它包括如下步骤：

1）运用常规方法将原始地震数据解编、去噪、消除自由表面多次波；

2）运用常规方法建立观测系统，对原始数据抽道集得到共中心点道集；

3）将共中心点道集数据进行二维傅立叶变换得到频率-波数域数据 $\stackrel{\‾}{D}(k_s{k}_g;\mathrm{\ω});$

4）建立微分方程：

${\tilde{L}}_0{G}_0=-\mathrm{\δ}(\stackrel{\→}{r}-{\stackrel{\→}{r}}_s)---\left(1\right),$

$\tilde{L}G=-\mathrm{\δ}\left(\stackrel{\→}{r}{-\stackrel{\→}{r}}_s\right)$

其中是参考及真实波场的微分算子，是波场点和源点的位置，设扰动算子为散射场算子为则有

上式写为散射基本方程为

（2）式进一步展成逆散射级数形式，

如果把写成级数，即为

$\tilde{V}=V_1+V_2+V_3+\·\·\·={\mathrm{\Σ}}_{n=1}^{\∞}{V}_n---\left(4\right),$

将（4）式带到（3）式中并展开，取同阶项相等，则

一阶项为

二阶项为

三阶项为

ｎ阶项为

G₀可分解为与自由表面有关的项和与自由表面无关的项自由表面多次波由项产生，而一次反射波和层间多次波由项产生，将代入（3）式并提取与有关的项，则得到

一阶项为

$D_1=G_0^d{V}_1{G}_0^d---\left(9\right),$

二阶项为

$D_2=G_0^d{V}_2{G}_0^d-G_0^d{V}_1{G}_0^f{V}_1{G}_0^d---\left(10\right),$

三阶项为

$D_3=G_0^d{V}_3{G}_0^d-G_0^d{V}_1{G}_0^f{V}_2{G}_0^d---\left(11\right),$

N阶项为

$D_n=G_0^d{V}_n{G}_0^d-G_0^d{V}_1{G}_0^f{V}_{n-1}{G}_0^d---\left(12\right);$

5）设散射场在接收面上记录到的波场为D，则波场D经过直达波切除和鬼波消除后得到波场即为D₁，而相应的D₂为一阶自由表面相关多次波，D₃为二阶自由表面相关多次波，D_n为n-1阶自由表面相关多次波；从（9）式解出V₁，代入（12）式则可解出V₂和D₂，再将V₁和V₂代入（11）式解出V₃和D₃，如此递归直至解出D_n，则得到消除与自由表面相关多次波后的波场D₀为

$D_0={\mathrm{\Σ}}_{n=1}^{\∞}{G}_0^d{V}_n{G}_0^d={\mathrm{\Σ}}_{n=1}^{\∞}{D}_n---\left(13\right);$

6）假设源点位置为(x_s,z_s)，接收点位置为(x_g,z_g)，则对于参考介质任意点(x,z)的波场可通过（1）式解出，即（1）式可写为

$(\frac{{\▿}^2}{{\mathrm{\ρ}}_0}+\frac{{\mathrm{\ω}}^2}{{\mathrm{\κ}}_0})G_0(x,z,x_s,z_s;\mathrm{\ω})=-\mathrm{\δ}(x-x_s)[\mathrm{\δ}\left({z-z}_s\right)-\mathrm{\δ}\left({z+z}_s\right)]---\left(14\right),$

其中，ρ₀为参考介质的密度，κ₀为参考介质的体积模量，

（14）式对x作傅立叶变换，则有

$(\frac{1}{{\mathrm{\ρ}}_0}\frac{d^2}{{\mathrm{dz}}^2}+\frac{q^2}{{\mathrm{\ρ}}_0})G_0(k_x,z,x_s,z_s;\mathrm{\ω})=-\frac{1}{\sqrt{2\mathrm{\π}}}{e}^{-i{k}_x{z}_x}[\mathrm{\δ}(z-z_s)-\mathrm{\δ}(z+z_s)]---\left(\text{15}\right)$

（15）式的通解为

$G_0(k_x,z,x_s,z_s;\mathrm{\ω})=\frac{{\mathrm{\ρ}}_0}{\sqrt{2\mathrm{\π}}}\frac{e^{-i{k}_x{z}_x}}{-2\mathrm{iq}}[e^{\mathrm{iq}|z-z_s|}-e^{\mathrm{iq}|z+z_s|}]---\left(16\right)$

其中， $q=\mathrm{sgn}\left(\mathrm{\ω}\right)\sqrt{{(\mathrm{\ω}/c_0)}^2-k_x^2},$ 为垂直波数， $c_0=\sqrt{{\mathrm{\κ}}_0/c_0},$ 为海水的纵波速度，因为 $G_0=G_0^d+G_0^f,$ 则

$G_0^d(k_x,z,x_s{,z}_s;\mathrm{\ω})=\frac{{\mathrm{\ρ}}_0}{\sqrt{2\mathrm{\π}}}\frac{e^{{-\mathrm{ik}}_x{z}_x}}{-2\mathrm{iq}}{e}^{\mathrm{iq}|z-z_s|}---\left(17\right)$

将（17）式代入（9）式解出V₁：

$V_1(k_s,k_g;\mathrm{\ω})={\left[e^{{\mathrm{iq}}_g{z}_g}{e}^{{\mathrm{iq}}_s{z}_s}\right]}^{-1}{q}_s{q}_g{D}_1(k_s,q_s,k_g,q_g;\mathrm{\ω})---\left(18\right)$

其中， $q_g=\mathrm{sgn}\left(\mathrm{\ω}\right)\sqrt{{(\mathrm{\ω}/c_0)}^2-k_g^2},$ $q_s=\mathrm{sgn}\left(\mathrm{\ω}\right)\sqrt{{(\mathrm{\ω}/c_0)}^2-k_s^2}$

将（17）式代入（13）式解出

${\stackrel{\‾}{D}}_n(k_s,k_g;\mathrm{\ω})=\frac{1}{\mathrm{i\π}{\mathrm{\ρ}}_{0}S\left(\mathrm{\ω}\right)}{\∫}_{-\∞}^{\∞}{e}^{\mathrm{iq}(z_g+z_s)}{\stackrel{\‾}{D}}_1(k_s,k;\mathrm{\ω}){\stackrel{\‾}{D}}_{n-1}(k,k_g;\mathrm{\ω})\mathrm{dk}(n=\mathrm{2,3,4},\·\·\·)---\left(19\right)$

对于（19）式，当n=2时有

${\stackrel{\‾}{D}}_2(k_s,k_g;\mathrm{\ω})=\frac{1}{\mathrm{i\π}{\mathrm{\ρ}}_{0}S\left(\mathrm{\ω}\right)}{\∫}_{-\∞}^{\∞}{e}^{\mathrm{iq}(z_g+z_s)}{\stackrel{\‾}{D}}_1(k_s,k;\mathrm{\ω}){\stackrel{\‾}{D}}_1(k,k_g;\mathrm{\ω})\mathrm{dk}---\left(20\right)$

即为地震数据的频率-波数域形式，

由公式（14）～公式（16），依据散射波场形成层间多次波的两个约束条件，及下行散射点的深度Z_g比相邻的两个散射点Z_g’和Z_g’’小，采用公式(12)进行迭代计算出V_n,即得到有关层间多次波的级数项；

7）将得到的有关层间多次波的逆散射级数项进行反傅立叶变换，得到空间时间域的层间多次波并将多次波成像；

8）运用常规方法将原始地震数据与预测出来的层间多次波数据匹配相减，从而得到消除层间多次波以后的CMP道集。