CN102116870A - 弹性波高斯束叠前深度偏移技术 - Google Patents

Abstract

本发明提供弹性波高斯束叠前深度偏移技术，步骤为1)根据选择的高斯束初始宽度，通过高斯窗将单炮记录划分为一系列重叠的区域，每个区域的中心为一个beam center；2)将多分量地震记录利用局部倾斜叠加分解为不同波型的局部平面波分量；3)分别从震源和beam center处试射弹性动力学高斯束，计算出每条高斯束有效范围内的走时、振幅信息；4)根据推导PP波和PS波的成像公式，利用高斯束的走时、振幅、极性信息以及地表的局部倾角信息计算对应的某个beam center的成像值；5)重复3)、4)步的运算，得到单炮的成像结果，最终的偏移结果是所有单炮成像结果的叠加；该方法可以对多次波至进行成像，具有较高的精度，并且还兼顾计算效率和灵活性，是一种具有较高应用价值的多分量地震成像方法。

Description

弹性波高斯束叠前深度偏移技术

所属技术领域

本发明属于多分量地震资料处理领域，是一种适用于多波多分量地震数据的深度域成像方法。

背景技术

随着物探技术装备水平的提高，多波(多分量)地震勘探日益受到工业界的重视，并已逐渐进入工业化生产。多分量地震技术在非均质性储层的含油气预测、油气藏的精细描述和动态监测中已显示出独特的优势和巨大的应用潜力。然而，目前已有的适用于多分量地震数据的成像方法，还不够完善，往往存在成像精度以及计算效率之间的矛盾。

在传统的多分量地震数据处理过程中，往往首先通过波场分离方法将地震数据分为纵(P)波和横(S)波，然后利用传统的声波偏移成像方法，对分离后的P波和S波分别进行成像。然而，上述处理方法没有充分利用弹性波场的矢量特性；此外，在波场分离的过程中往往不能将不同波型的能量完全分离，残余的非本型波能量，往往导致成像结果中的大量噪声串扰，进而严重影响成像效果。

基于矢量波场延拓的弹性波地震成像方法，是解决上述问题的有效途径。由于常规单程波成像方法无法正确描述耦合的弹性波传播过程，现有的弹性波成像方法都是在时间域进行的，其大致可以分为两类：一类是弹性波逆时偏移，其基于弹性波波动方程，将多分量地震数据作为边界条件进行逆时外推，然后利用激发时成像条件或者互相关成像条件求取成像值。Sun和McMechan(1986)首先将逆时偏移应用到弹性波偏移中，此后，Chang和McMechan(1994)，Botelho和Stoffa(1991)，Yan和Sava(2008)以及杜启振和秦童(2009)等人将其进行推广并应用到三维以及各向异性介质的偏移成像中。另一类是弹性波Kirchhoff偏移，其通过计算本型波以及转换波的走时、振幅以及极性信息，随后沿着所选择波型的走时轨迹对波场能量进行叠加并求取成像值。Kuo和Dai(1984)基于均匀各向同性介质中位移格林函数的远场近似，首先提出了该方法。此后，Sena和Toksoz(1993)，Hokstad(2000)以及Druzhinin(2003)等人对其进行了发展。上述两类方法，具有各自的优势及不足，其中弹性波逆时偏移虽然具有极高的成像精度，但是其计算量非常巨大，而且其成像结果往往伴有大量的低频干扰以及转换波成像时难以解决的极性反转问题。弹性波Kirchhoff偏移，虽然具有灵活高效的特点，但是其难以对多次波至进行成像，致使其对地下复杂地质构造的成像精度不高。

发明内容

本发明的目的在于提供一种弹性波高斯束叠前深度偏移技术，是一种兼顾精度与效率的弹性波地震成像方法，以多分量地震数据为输入，利用弹性波动力学高斯束的振幅、走时以及极性信息来进行矢量波场的延拓，并且利用波场的矢量投影来提取纵横波的成像值。

本发明所采用技术方案是：弹性波高斯束叠前深度偏移技术，实现步骤为：

(1)根据选择的高斯束初始宽度，通过高斯窗将单炮记录划分为一系列重叠的区域，每个区域的中心为一个beam center；

(2)将多分量地震记录利用局部倾斜叠加分解为不同波型的局部平面波分量；

(3)分别从震源和beam center处试射弹性动力学高斯束，计算出每条高斯束有效范围内的走时、振幅信息；

(4)根据推导PP波和PS波的成像公式，利用高斯束的走时、振幅、极性信息以及地表的局部倾角信息计算对应的某个beam center的成像值；

(5)重复(3)、(4)步的运算，得到单炮的成像结果，最终的偏移结果是所有单炮成像结果的叠加。

其中步骤(3)中的计算方法是：

在S所描述的观测面上，假设由震源x_s激发，接收点x_r接收到的两分量弹性波地震记录为u_i(x_r；ω)，则反向延拓的弹性波位移场u_m(x；x_r；ω)可以通过Kirchhoff-Helmholtz积分来表示：

$u_m(x;x_r;\mathrm{\ω})={\∫}_S{\mathrm{dx}}_r[t_i(x_r;\mathrm{\ω})G_{\mathrm{im}}^{*}(x;x_r;\mathrm{\ω})-u_i(x_r;\mathrm{\ω}){\mathrm{\Σ}}_{\mathrm{im}}^{*}(x;x_r;\mathrm{\ω})]---\left(1\right)$

其中*代表复共轭；G_lm(x；x_r；ω)为位移格林张量，其代表x_r处l方向单位体力所造成的x处位移的m方向的分量；t_i(x_r)为x_r处应力；∑_im(x；x_r)为应力格林张量；

若假设S为自由地表，则根据自由应力边界条件t(x；ω)＝0，从而式(1)变为：

$u_m(x;x_r;\mathrm{\ω})=-{\∫}_S{\mathrm{dx}}_r{u}_i(x_r;\mathrm{\ω}){\mathrm{\Σ}}_{\mathrm{im}}^{*}(x;x_r;\mathrm{\ω})---\left(2\right)$

并取G_lm(x；x_r；ω)的高频近似解，代入式(2)，得

$u_m(x;x_r;\mathrm{\ω})=u_m^p(x;x_r;\mathrm{\ω})+u_m^s(x;x_r;\mathrm{\ω})$

$=-\mathrm{i\ω}\underset{v}{\mathrm{\Σ}}{\∫}_S{\mathrm{dx}}_r\mathrm{\ρ}\left(x_r\right)U_m^{v*}(x;x_r;\mathrm{\ω})[u_1(x_r;\mathrm{\ω})W_1^v\left(x_r\right)+u_2(x_r;\mathrm{\ω})W_2^v\left(x_r\right)]$

(3)

式(3)为解耦的弹性波波场延拓公式，其中，

分别为位移场u_m(x；x_r；ω)中的P波和S波成分；

权值具有以下形式

$\left\{\begin{array}{c}{W}_1^p\left(x_r\right)=2v_s^2\left(x_r\right)p_2^p\left(x_r\right)e_1^p\left(x_r\right)\\ W_2^p\left(x_r\right)=2v_s^2\left(x_r\right)p_2^p\left(x_r\right)e_2^p\left(x_r\right)+\left(\frac{v_p^2\left(x_r\right)-2v_s^2\left(x_r\right)}{v_p\left(x_r\right)}\right)\\ W_1^s\left(x_r\right)=v_s^2\left(x_r\right)p_2^s\left(x_r\right)e_1^s\left(x_r\right)+v_s^2\left(x_r\right)p_1^s\left(x_r\right)e_2^s\left(x_r\right)\\ W_2^s\left(x_r\right)=-2v_s^2\left(x_r\right)p_1^s\left(x_r\right)e_1^s\left(x_r\right)\end{array}\right.---\left(4\right)$

利用弹性动力学高斯束的叠加积分来计算式(3)中的位移矢量

得到基于高斯束表示的弹性波场延拓公式；

高斯束的初始波前为平面的，将多分量地震记录分解为不同波型不同初始方向的局部平面波，然后利用相对应高斯束进行波场延拓，其具体实现过程如下：

对地震记录加入一系列高斯窗，将其分为一系列重叠的局部区域，高斯窗函数具有如下性质

$\frac{\mathrm{\ΔL}}{\sqrt{2\mathrm{\π}}{w}_0}\sqrt{\left|\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right|}\underset{L}{\mathrm{\Σ}}\exp [-|\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\left|\frac{(x_r-L)}{2w_0^2}\right]\≈1---\left(5\right)$

其中，x＝L为高斯窗的中线，也是波束中心L的水平坐标；ΔL为其水平间隔，将式(5)代入式(3)得

$u_m(x;x_r;\mathrm{\ω})=-\frac{\mathrm{i\ω\ΔL}}{\sqrt{2\mathrm{\π}}{w}_0}\sqrt{\left|\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right|}\underset{v}{\mathrm{\Σ}}\underset{L}{\mathrm{\Σ}}{\∫}_S{\mathrm{dx}}_r\mathrm{\ρ}\left(x_r\right)U_m^{v*}(x;x_r;\mathrm{\ω})[u_1(x_r;\mathrm{\ω})W_1^v\left(x_r\right)$ (6)

$+u_2(x_r;\mathrm{\ω})W_2^v\left(x_r\right)\left]\exp \right[-\left|\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right|\frac{(x_r-L)}{2w_0^2}]$

引入相移校正因子，将

通过L处出射的高斯束来表示

$U_m^v(x;x_r;\mathrm{\ω})\≈{\mathrm{\Ψ}}^v\∫{\hat{u}}_m^v(x;L;\mathrm{\ω})\exp [-\mathrm{i\ω}{p}_1^v\left(L\right)(x_r-L)]\frac{{\mathrm{dp}}_1^v\left(L\right)}{p_2^v\left(L\right)}---\left(7\right)$

将上式代入式(6)，令ρ(x_r)≈ρ(L)，

并交换积分次序，便可以得到反向延拓的P波位移以及S波位移

$u_m^p(x;x_r;\mathrm{\ω})=-\frac{\mathrm{\ΔL\ω}}{4\mathrm{\π}}\underset{L}{\mathrm{\Σ}}\∫\frac{d{p}_1^p\left(L\right)}{p_2^p\left(L\right)}\sqrt{\mathrm{\ρ}\left(L\right)}{\hat{u}}_m^{p*}(x;L;\mathrm{\ω})---\left(8\right)$

$\×[W_1^p\left(L\right)D_1^p(L;p_1^p;\mathrm{\ω})+W_2^p\left(L\right)D_2^p(L;p_1^p;\mathrm{\ω})]$

$u_m^s(x;x_r;\mathrm{\ω})=-\frac{\mathrm{\ΔL\ω}}{4\mathrm{\π}}\underset{L}{\mathrm{\Σ}}\∫\frac{d{p}_1^s\left(L\right)}{p_2^s\left(L\right)}\sqrt{\mathrm{\ρ}\left(L\right)}{\hat{u}}_m^{s*}(x;L;\mathrm{\ω})---\left(9\right)$

$\×[W_1^s\left(L\right)D_1^s(L;p_1^s;\mathrm{\ω})+W_2^s\left(L\right)D_2^s(L;p_1^s;\mathrm{\ω})]$

其中，为对不同波型多分量地震记录的加窗局部倾斜叠加

$D_n^v(L;P_1^v;\mathrm{\ω})=\sqrt{\frac{\left|\mathrm{\ω}\right|}{2\mathrm{\π}}}{\∫}_S{\mathrm{dx}}_r{u}_n(x_r;\mathrm{\ω})\exp [\mathrm{i\ω}{p}_1^v\left(L\right)(x_r-L)-|\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\left|\frac{{(x_r-L)}^2}{2w_0^2}\right]---\left(10\right)$

权值此时为

$W_1^p\left(L\right)=2{\mathrm{\γ}}^2\left(L\right)p_2^p\left(L\right)e_1^p\left(L\right),$ $W_2^p\left(L\right)=2{\mathrm{\γ}}^2\left(L\right)p_2^p\left(L\right)e_2^p\left(L\right)+\left(\frac{1-2{\mathrm{\γ}}^2\left(L\right)}{v_p\left(L\right)}\right)$ $W_1^s\left(L\right)=p_2^s\left(L\right)e_1^s\left(L\right)+p_1^s\left(L\right)e_2^s\left(L\right),$ $W_2^s\left(L\right)=-2p_1^s\left(L\right)e_1^s\left(L\right),$ $\mathrm{\γ}\left(L\right)=\frac{v_s\left(L\right)}{v_p\left(L\right)}$

(11)

其中步骤(4)中，根据Claerbout成像原理，首先将震源位移波场

通过弹性动力学高斯束来表示

$U_m^p(x;x_s;\mathrm{\ω})\≈\frac{i}{4\mathrm{\π}{v}_p^2\left(x_s\right)}\sqrt{\frac{{\mathrm{\ω}}_r{w}_0^2}{\mathrm{\ρ}\left(x_s\right)}}\∫\frac{d{p}_1^p\left(x_s\right)}{p_2^p\left(x_s\right)}{\hat{u}}_m^p(x;x_s;\mathrm{\ω})---\left(12\right)$

接下来根据P波以及S波的传播特点，结合式(8)，式(9)所表示的反向延拓的接收波场，定义如下成像公式

$I^{\mathrm{pp}}\left(x\right)={\∫U}_2^{p*}(x;x_s;\mathrm{\ω})u_2^p(x;x_r;\mathrm{\ω})\mathrm{d\ω}$

$=\frac{\mathrm{\ΔL}\sqrt{{\mathrm{\ω}}_r{w}_0^2}}{16{\mathrm{\π}}^2}\underset{L}{\mathrm{\Σ}}\∫\mathrm{d\ω}\frac{\mathrm{i\ω}}{v_p^2\left(x_s\right)}\sqrt{\frac{\mathrm{\ρ}\left(L\right)}{\mathrm{\ρ}\left(x_s\right)}}\∫\∫\frac{d{p}_1^p\left(x_s\right)d{p}_1^p\left(L\right)}{p_2^p\left(x_s\right)p_2^p\left(L\right)}---\left(13\right)$

$\×{\hat{u}}_2^{p*}(x;x_s;\mathrm{\ω}){\hat{u}}_1^{p*}(x;L;\mathrm{\ω})[W_1^p\left(L\right)D_1^p(L;p_1^p;\mathrm{\ω})+W_2^p\left(L\right)D_2^p(L;p_1^p;\mathrm{\ω})]$

$I^{\mathrm{ps}}\left(x\right)=\frac{\mathrm{\ΔL}\sqrt{{\mathrm{\ω}}_r{w}_0^2}}{16{\mathrm{\π}}^2}\underset{L}{\mathrm{\Σ}}\∫\mathrm{d\ω}\frac{\mathrm{i\ω}}{v_p^2\left(x_s\right)}\sqrt{\frac{\mathrm{\ρ}\left(L\right)}{\mathrm{\ρ}\left(x_s\right)}}\∫\∫\frac{d{p}_1^p\left(x_s\right)d{p}_1^s\left(L\right)}{p_2^p\left(x_s\right)p_2^s\left(L\right)}\mathrm{sgn}\left(\mathrm{\α}\right)---\left(14\right)$

$\×{\hat{u}}_2^{p*}(x;x_s;\mathrm{\ω}){\hat{u}}_1^{s*}(x;L;\mathrm{\ω})[W_1^s\left(L\right)D_1^s(L;p_1^s;\mathrm{\ω})+W_2^s\left(L\right)D_2^s(L;p_1^s;\mathrm{\ω})]$

其中，I^pp(x)为P-P单炮成像值，I^ps(x)为P-S单炮成像值，将所有炮记录进行偏移并叠加，便得到最终的弹性波成像结果。

本发明的有益效果是：该方法不但可以对多次波至进行成像，具有较高的精度，并且还兼顾计算效率和灵活性，是一种具有较高应用价值的多分量地震成像方法。

附图说明

图1是具体实施方式的平层模型偏移试算效果图。图中(a)水平分量单炮记录；(b)垂直分量单炮记录；(c)垂直分量记录P-P波单炮成像结果；(d)两分量记录的P-P单炮成像结果；(e)水平分量记录P-S波单炮成像结果；(f)经过极性校正后的两分量记录的P-S成像结果；(g)P-P波成像叠加结果；(h)P-S波成像叠加结果。

图2是P-S转换波在反射界面处的偏振图。

图3是实施例的复杂模型偏移试算结果图。其中(a)P波速度模型；(b)S波速度模型；(c)X分量单炮记录；(d)Z分量单炮记录；(e)P-P波成像叠加结果；(f)P-S波成像叠加结果。

具体实施方式

弹性波高斯束叠前深度偏移技术，实现步骤为：

(1)根据选择的高斯束初始宽度，通过高斯窗将单炮记录划分为一系列重叠的区域，每个区域的中心为一个beam center；

(2)将多分量地震记录利用局部倾斜叠加分解为不同波型的局部平面波分量；

(3)分别从震源和beam center处试射弹性动力学高斯束，计算出每条高斯束有效范围内的走时、振幅信息；

(4)根据推导PP波和PS波的成像公式，利用高斯束的走时、振幅、极性信息以及地表的局部倾角信息计算对应的某个beam center的成像值；

(5)重复(3)、(4)步的运算，得到单炮的成像结果，最终的偏移结果是所有单炮成像结果的叠加。

其中步骤(3)中的计算方法是：

在S所描述的观测面上，假设由震源x_s激发，接收点x_r接收到的两分量弹性波地震记录为u_i(x_r；ω)，则反向延拓的弹性波位移场u_m(x；x_r；ω)可以通过Kirchhoff-Helmholtz积分来表示：

$u_m(x;x_r;\mathrm{\ω})={\∫}_S{\mathrm{dx}}_r[t_i(x_r;\mathrm{\ω})G_{\mathrm{im}}^{*}(x;x_r;\mathrm{\ω})-u_i(x_r;\mathrm{\ω}){\mathrm{\Σ}}_{\mathrm{im}}^{*}(x;x_r;\mathrm{\ω})]---\left(1\right)$

其中*代表复共轭；G_lm(x；x_r；ω)为位移格林张量，其代表x_r处l方向单位体力所造成的x处位移的m方向的分量；t_i(x_r)为x_r处应力；∑_im(x；x_r)为应力格林张量；

若假设S为自由地表，则根据自由应力边界条件t(x；ω)＝0，从而式(1)变为：

$u_m(x;x_r;\mathrm{\ω})=-{\∫}_S{\mathrm{dx}}_r{u}_i(x_r;\mathrm{\ω}){\mathrm{\Σ}}_{\mathrm{im}}^{*}(x;x_r;\mathrm{\ω})---\left(2\right)$

并取G_lm(x；x_r；ω)的高频近似解，代入式(2)，得

$u_m(x;x_r;\mathrm{\ω})=u_m^p(x;x_r;\mathrm{\ω})+u_m^s(x;x_r;\mathrm{\ω})$

$=-\mathrm{i\ω}\underset{v}{\mathrm{\Σ}}{\∫}_S{\mathrm{dx}}_r\mathrm{\ρ}\left(x_r\right)U_m^{v*}(x;x_r;\mathrm{\ω})[u_1(x_r;\mathrm{\ω})W_1^v\left(x_r\right)+u_2(x_r;\mathrm{\ω})W_2^v\left(x_r\right)]$

(3)

式(3)为解耦的弹性波波场延拓公式，其中，

分别为位移场u_m(x；x_r；ω)中的P波和S波成分；

权值

具有以下形式

$\left\{\begin{array}{c}{W}_1^p\left(x_r\right)=2v_s^2\left(x_r\right)p_2^p\left(x_r\right)e_1^p\left(x_r\right)\\ W_2^p\left(x_r\right)=2v_s^2\left(x_r\right)p_2^p\left(x_r\right)e_2^p\left(x_r\right)+\left(\frac{v_p^2\left(x_r\right)-2v_s^2\left(x_r\right)}{v_p\left(x_r\right)}\right)\\ W_1^s\left(x_r\right)=v_s^2\left(x_r\right)p_2^s\left(x_r\right)e_1^s\left(x_r\right)+v_s^2\left(x_r\right)p_1^s\left(x_r\right)e_2^s\left(x_r\right)\\ W_2^s\left(x_r\right)=-2v_s^2\left(x_r\right)p_1^s\left(x_r\right)e_1^s\left(x_r\right)\end{array}\right.---\left(4\right)$

利用弹性动力学高斯束的叠加积分来计算式(3)中的位移矢量

得到基于高斯束表示的弹性波场延拓公式；

高斯束的初始波前为平面的，将多分量地震记录分解为不同波型不同初始方向的局部平面波，然后利用相对应高斯束进行波场延拓，其具体实现过程如下：

对地震记录加入一系列高斯窗，将其分为一系列重叠的局部区域，高斯窗函数具有如下性质

$\frac{\mathrm{\ΔL}}{\sqrt{2\mathrm{\π}}{w}_0}\sqrt{\left|\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right|}\underset{L}{\mathrm{\Σ}}\exp [-|\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\left|\frac{(x_r-L)}{2w_0^2}\right]\≈1---\left(5\right)$

其中，x＝L为高斯窗的中线，也是波束中心L的水平坐标；ΔL为其水平间隔，将式(5)代入式(3)得

$u_m(x;x_r;\mathrm{\ω})=-\frac{\mathrm{i\ω\ΔL}}{\sqrt{2\mathrm{\π}}{w}_0}\sqrt{\left|\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right|}\underset{v}{\mathrm{\Σ}}\underset{L}{\mathrm{\Σ}}{\∫}_S{\mathrm{dx}}_r\mathrm{\ρ}\left(x_r\right)U_m^{v*}(x;x_r;\mathrm{\ω})[u_1(x_r;\mathrm{\ω})W_1^v\left(x_r\right)$ (6)

$+u_2(x_r;\mathrm{\ω})W_2^v\left(x_r\right)\left]\exp \right[-\left|\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right|\frac{(x_r-L)}{2w_0^2}]$

引入相移校正因子，将通过L处出射的高斯束来表示

$U_m^v(x;x_r;\mathrm{\ω})\≈{\mathrm{\Ψ}}^v\∫{\hat{u}}_m^v(x;L;\mathrm{\ω})\exp [-\mathrm{i\ω}{p}_1^v\left(L\right)(x_r-L)]\frac{{\mathrm{dp}}_1^v\left(L\right)}{p_2^v\left(L\right)}---\left(7\right)$

将上式代入式(6)，令ρ(x_r)≈ρ(L)，

并交换积分次序，便可以得到反向延拓的P波位移

以及S波位移

$u_m^p(x;x_r;\mathrm{\ω})=-\frac{\mathrm{\ΔL\ω}}{4\mathrm{\π}}\underset{L}{\mathrm{\Σ}}\∫\frac{d{p}_1^p\left(L\right)}{p_2^p\left(L\right)}\sqrt{\mathrm{\ρ}\left(L\right)}{\hat{u}}_m^{p*}(x;L;\mathrm{\ω})---\left(8\right)$

$\×[W_1^p\left(L\right)D_1^p(L;p_1^p;\mathrm{\ω})+W_2^p\left(L\right)D_2^p(L;p_1^p;\mathrm{\ω})]$

$u_m^s(x;x_r;\mathrm{\ω})=-\frac{\mathrm{\ΔL\ω}}{4\mathrm{\π}}\underset{L}{\mathrm{\Σ}}\∫\frac{d{p}_1^s\left(L\right)}{p_2^s\left(L\right)}\sqrt{\mathrm{\ρ}\left(L\right)}{\hat{u}}_m^{s*}(x;L;\mathrm{\ω})---\left(9\right)$

$\×[W_1^s\left(L\right)D_1^s(L;p_1^s;\mathrm{\ω})+W_2^s\left(L\right)D_2^s(L;p_1^s;\mathrm{\ω})]$

其中，

为对不同波型多分量地震记录的加窗局部倾斜叠加

$D_n^v(L;P_1^v;\mathrm{\ω})=\sqrt{\frac{\left|\mathrm{\ω}\right|}{2\mathrm{\π}}}{\∫}_S{\mathrm{dx}}_r{u}_n(x_r;\mathrm{\ω})\exp [\mathrm{i\ω}{p}_1^v\left(L\right)(x_r-L)-|\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\left|\frac{{(x_r-L)}^2}{2w_0^2}\right]---\left(10\right)$

权值

此时为

$W_1^p\left(L\right)=2{\mathrm{\γ}}^2\left(L\right)p_2^p\left(L\right)e_1^p\left(L\right),$ $W_2^p\left(L\right)=2{\mathrm{\γ}}^2\left(L\right)p_2^p\left(L\right)e_2^p\left(L\right)+\left(\frac{1-2{\mathrm{\γ}}^2\left(L\right)}{v_p\left(L\right)}\right)$ $W_1^s\left(L\right)=p_2^s\left(L\right)e_1^s\left(L\right)+p_1^s\left(L\right)e_2^s\left(L\right),$ $W_2^s\left(L\right)=-2p_1^s\left(L\right)e_1^s\left(L\right),$ $\mathrm{\γ}\left(L\right)=\frac{v_s\left(L\right)}{v_p\left(L\right)}$

(11)

其中步骤(4)中，根据Claerbout成像原理，首先将震源位移波场

通过弹性动力学高斯束来表示

$U_m^p(x;x_s;\mathrm{\ω})\≈\frac{i}{4\mathrm{\π}{v}_p^2\left(x_s\right)}\sqrt{\frac{{\mathrm{\ω}}_r{w}_0^2}{\mathrm{\ρ}\left(x_s\right)}}\∫\frac{d{p}_1^p\left(x_s\right)}{p_2^p\left(x_s\right)}{\hat{u}}_m^p(x;x_s;\mathrm{\ω})---\left(12\right)$

接下来根据P波以及S波的传播特点，结合式(8)，式(9)所表示的反向延拓的接收波场，定义如下成像公式

$I^{\mathrm{pp}}\left(x\right)={\∫U}_2^{p*}(x;x_s;\mathrm{\ω})u_2^p(x;x_r;\mathrm{\ω})\mathrm{d\ω}$

$=\frac{\mathrm{\ΔL}\sqrt{{\mathrm{\ω}}_r{w}_0^2}}{16{\mathrm{\π}}^2}\underset{L}{\mathrm{\Σ}}\∫\mathrm{d\ω}\frac{\mathrm{i\ω}}{v_p^2\left(x_s\right)}\sqrt{\frac{\mathrm{\ρ}\left(L\right)}{\mathrm{\ρ}\left(x_s\right)}}\∫\∫\frac{d{p}_1^p\left(x_s\right)d{p}_1^p\left(L\right)}{p_2^p\left(x_s\right)p_2^p\left(L\right)}---\left(13\right)$

$\×{\hat{u}}_2^{p*}(x;x_s;\mathrm{\ω}){\hat{u}}_1^{p*}(x;L;\mathrm{\ω})[W_1^p\left(L\right)D_1^p(L;p_1^p;\mathrm{\ω})+W_2^p\left(L\right)D_2^p(L;p_1^p;\mathrm{\ω})]$

$I^{\mathrm{ps}}\left(x\right)=\frac{\mathrm{\ΔL}\sqrt{{\mathrm{\ω}}_r{w}_0^2}}{16{\mathrm{\π}}^2}\underset{L}{\mathrm{\Σ}}\∫\mathrm{d\ω}\frac{\mathrm{i\ω}}{v_p^2\left(x_s\right)}\sqrt{\frac{\mathrm{\ρ}\left(L\right)}{\mathrm{\ρ}\left(x_s\right)}}\∫\∫\frac{d{p}_1^p\left(x_s\right)d{p}_1^s\left(L\right)}{p_2^p\left(x_s\right)p_2^s\left(L\right)}\mathrm{sgn}\left(\mathrm{\α}\right)---\left(14\right)$

$\×{\hat{u}}_2^{p*}(x;x_s;\mathrm{\ω}){\hat{u}}_1^{s*}(x;L;\mathrm{\ω})[W_1^s\left(L\right)D_1^s(L;p_1^s;\mathrm{\ω})+W_2^s\left(L\right)D_2^s(L;p_1^s;\mathrm{\ω})]$

其中，I^pp(x)为P-P单炮成像值，I^ps(x)为P-S单炮成像值，将所有炮记录进行偏移并叠加，便得到最终的弹性波成像结果。

首先通过一个三层水平层状介质模型来测试本发明的正确性。模型网格为301×401，纵横向采样间隔都为10m，两个水平反射层分别位于1400m以及2700m处。随着深度的增大，P波速度分别为2500m/s，2900m/s，3200m/s，S波速度分别为1500m/s，1800m/s，2600m/s，密度设定为常数。模拟弹性波记录主频为30Hz，共有51炮，每炮151道，道间距为20m，时间采样点数为1500，采样间隔为2ms。图1a，图1b分别为第26炮单炮记录的X，Z分量。图1c为仅对Z分量记录的P-P单炮成像结果，可以明显看到Z分量记录中的S波成分所导致的串扰噪声(箭头所指处)，图1d为采用式(15)、(23)对两分量记录加权后得到的P-P单炮成像结果，此时串扰噪声得到了有效的压制。图1e为X分量记录的P-S单炮成像结果，图1f为采用式(15)，(23)对两分量记录加权后的P-S单炮成像结果，可以看到，相比于图1e，图1f中的串扰噪声得到了有效的压制(箭头所指处)，且成像结果中的极性反转现象得到了有效的校正(椭圆标识处)。图1g，1h分别为对所有炮偏移叠加后的P-P，P-S成像结果，可以两者均清晰的反映了界面的真实形态，且由于S波传播速度较慢，使得P-S成像结果具有相对更高的分辨率。

在对转换波进行成像时，往往出现转换波偏振方向不同而导致的成像剖面中的极性反转现象，严重影响叠加成像的效果。针对上述问题，本发明提出一种校正方法。图2显示了P-S转换波的传播过程，可以看到由于入射P波具有不同符号的入射角(α₁＞0，α₂＜0，以逆时针为正)，使得O₁点产生的P-S转换波具有正水平方向(以向右为正)的位移分量，O₂产生的P-S转换波具有负水平方向的位移分量，最终使R₁，R₂点接收到的X分量地震记录具有相反的极性。也就是说X分量记录中转换S波的极性同P波入射角的正负有关，因此可以根据P波入射角α的正负直接对成像结果进行校正，为此引入符号函数，得到式(14)所示的P-S转换波成像公式。

利用一复杂构造模型来测试本发明对复杂地质模型的成像效果。速度模型如图3a，3b所示，其网格大小为1000×550，纵横向采样间隔都为10m，5m，同样将密度设定为常数。模型共分为7层，其中包含断层、洼陷等复杂构造。模拟弹性波记录主频为25Hz，共有200炮，每炮361道，道间距为10m，时间采样点数为1200，采样间隔为2.1ms。图3c，3d分别为第1炮单炮记录的X，Z分量，可以看到由于地下构造复杂导致记录中的波场信息也非常复杂。图3e，3f分别为最终的P-P，P-S成像结果，可以看到两者均良好的反映了速度模型的构造情况，模型中的断层、洼陷以及起伏的底层均得到了较好的成像，证明了本文方法的有效性。

Claims (3)

1.弹性波高斯束叠前深度偏移技术，其特征在于实现步骤为：

(1)根据选择的高斯束初始宽度，通过高斯窗将单炮记录划分为一系列重叠的区域，每个区域的中心为一个beam center；

(2)将多分量地震记录利用局部倾斜叠加分解为不同波型的局部平面波分量；

(3)分别从震源和beam center处试射弹性动力学高斯束，计算出每条高斯束有效范围内的走时、振幅信息；

(4)根据推导PP波和PS波的成像公式，利用高斯束的走时、振幅、极性信息以及地表的局部倾角信息计算对应的某个beam center的成像值；

(5)重复(3)、(4)步的运算，得到单炮的成像结果，最终的偏移结果是所有单炮成像结果的叠加。

2.根据权利要求1所述的弹性波高斯束叠前深度偏移技术，其特征在于步骤(3)中的计算方法是：

在S所描述的观测面上，假设由震源x_s激发，接收点x_r接收到的两分量弹性波地震记录为u_i(x_r；ω)，则反向延拓的弹性波位移场u_m(x；x_r；ω)可以通过Kirchhoff-Helmholtz积分来表示：

$u_m(x;x_r;\mathrm{\ω})={\∫}_S{\mathrm{dx}}_r[t_i(x_r;\mathrm{\ω})G_{\mathrm{im}}^{*}(x;x_r;\mathrm{\ω})-u_i(x_r;\mathrm{\ω}){\mathrm{\Σ}}_{\mathrm{im}}^{*}(x;x_r;\mathrm{\ω})]---\left(1\right)$

其中*代表复共轭；G_lm(x；x_r；ω)为位移格林张量，其代表x_r处l方向单位体力所造成的x处位移的m方向的分量；t_i(x_r)为x_r处应力；∑_im(x；x_r)为应力格林张量；

若假设S为自由地表，则根据自由应力边界条件t(x；ω)＝0，从而式(1)变为：

$u_m(x;x_r;\mathrm{\ω})=-{\∫}_S{\mathrm{dx}}_r{u}_i(x_r;\mathrm{\ω}){\mathrm{\Σ}}_{\mathrm{im}}^{*}(x;x_r;\mathrm{\ω})---\left(2\right)$

并取G_lm(x；x_r；ω)的高频近似解，代入式(2)，得

$u_m(x;x_r;\mathrm{\ω})=u_m^p(x;x_r;\mathrm{\ω})+u_m^s(x;x_r;\mathrm{\ω})$

$=-\mathrm{i\ω}\underset{v}{\mathrm{\Σ}}{\∫}_S{\mathrm{dx}}_r\mathrm{\ρ}\left(x_r\right)U_m^{v*}(x;x_r;\mathrm{\ω})[u_1(x_r;\mathrm{\ω})W_1^v\left(x_r\right)+u_2(x_r;\mathrm{\ω})W_2^v\left(x_r\right)]$

(3)

式(3)为解耦的弹性波波场延拓公式，其中，

分别为位移场u_m(x；x_r；ω)中的P波和S波成分；

权值

具有以下形式

$\left\{\begin{array}{c}{W}_1^p\left(x_r\right)=2v_s^2\left(x_r\right)p_2^p\left(x_r\right)e_1^p\left(x_r\right)\\ W_2^p\left(x_r\right)=2v_s^2\left(x_r\right)p_2^p\left(x_r\right)e_2^p\left(x_r\right)+\left(\frac{v_p^2\left(x_r\right)-2v_s^2\left(x_r\right)}{v_p\left(x_r\right)}\right)\\ W_1^s\left(x_r\right)=v_s^2\left(x_r\right)p_2^s\left(x_r\right)e_1^s\left(x_r\right)+v_s^2\left(x_r\right)p_1^s\left(x_r\right)e_2^s\left(x_r\right)\\ W_2^s\left(x_r\right)=-2v_s^2\left(x_r\right)p_1^s\left(x_r\right)e_1^s\left(x_r\right)\end{array}\right.---\left(4\right)$

利用弹性动力学高斯束的叠加积分来计算式(3)中的位移矢量

得到基于高斯束表示的弹性波场延拓公式；

高斯束的初始波前为平面的，将多分量地震记录分解为不同波型不同初始方向的局部平面波，然后利用相对应高斯束进行波场延拓，其具体实现过程如下：

对地震记录加入一系列高斯窗，将其分为一系列重叠的局部区域，高斯窗函数具有如下性质

$\frac{\mathrm{\ΔL}}{\sqrt{2\mathrm{\π}}{w}_0}\sqrt{\left|\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right|}\underset{L}{\mathrm{\Σ}}\exp [-|\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\left|\frac{(x_r-L)}{2w_0^2}\right]\≈1---\left(5\right)$

其中，x＝L为高斯窗的中线，也是波束中心L的水平坐标；ΔL为其水平间隔，将式(5)代入式(3)得

$u_m(x;x_r;\mathrm{\ω})=-\frac{\mathrm{i\ω\ΔL}}{\sqrt{2\mathrm{\π}}{w}_0}\sqrt{\left|\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right|}\underset{v}{\mathrm{\Σ}}\underset{L}{\mathrm{\Σ}}{\∫}_S{\mathrm{dx}}_r\mathrm{\ρ}\left(x_r\right)U_m^{v*}(x;x_r;\mathrm{\ω})[u_1(x_r;\mathrm{\ω})W_1^v\left(x_r\right)$ (6)

$+u_2(x_r;\mathrm{\ω})W_2^v\left(x_r\right)\left]\exp \right[-\left|\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\right|\frac{(x_r-L)}{2w_0^2}]$

引入相移校正因子，将

通过L处出射的高斯束来表示

$U_m^v(x;x_r;\mathrm{\ω})\≈{\mathrm{\Ψ}}^v\∫{\hat{u}}_m^v(x;L;\mathrm{\ω})\exp [-\mathrm{i\ω}{p}_1^v\left(L\right)(x_r-L)]\frac{{\mathrm{dp}}_1^v\left(L\right)}{p_2^v\left(L\right)}---\left(7\right)$

将上式代入式(6)，令ρ(x_r)≈ρ(L)，

并交换积分次序，便可以得到反向延拓的P波位移

以及S波位移

$u_m^p(x;x_r;\mathrm{\ω})=-\frac{\mathrm{\ΔL\ω}}{4\mathrm{\π}}\underset{L}{\mathrm{\Σ}}\∫\frac{d{p}_1^p\left(L\right)}{p_2^p\left(L\right)}\sqrt{\mathrm{\ρ}\left(L\right)}{\hat{u}}_m^{p*}(x;L;\mathrm{\ω})---\left(8\right)$

$\×[W_1^p\left(L\right)D_1^p(L;p_1^p;\mathrm{\ω})+W_2^p\left(L\right)D_2^p(L;p_1^p;\mathrm{\ω})]$

$u_m^s(x;x_r;\mathrm{\ω})=-\frac{\mathrm{\ΔL\ω}}{4\mathrm{\π}}\underset{L}{\mathrm{\Σ}}\∫\frac{d{p}_1^s\left(L\right)}{p_2^s\left(L\right)}\sqrt{\mathrm{\ρ}\left(L\right)}{\hat{u}}_m^{s*}(x;L;\mathrm{\ω})---\left(9\right)$

$\×[W_1^s\left(L\right)D_1^s(L;p_1^s;\mathrm{\ω})+W_2^s\left(L\right)D_2^s(L;p_1^s;\mathrm{\ω})]$

其中，

为对不同波型多分量地震记录的加窗局部倾斜叠加

$D_n^v(L;P_1^v;\mathrm{\ω})=\sqrt{\frac{\left|\mathrm{\ω}\right|}{2\mathrm{\π}}}{\∫}_S{\mathrm{dx}}_r{u}_n(x_r;\mathrm{\ω})\exp [\mathrm{i\ω}{p}_1^v\left(L\right)(x_r-L)-|\frac{\mathrm{\ω}}{{\mathrm{\ω}}_r}\left|\frac{{(x_r-L)}^2}{2w_0^2}\right]---\left(10\right)$

权值

此时为

$W_1^p\left(L\right)=2{\mathrm{\γ}}^2\left(L\right)p_2^p\left(L\right)e_1^p\left(L\right),$ $W_2^p\left(L\right)=2{\mathrm{\γ}}^2\left(L\right)p_2^p\left(L\right)e_2^p\left(L\right)+\left(\frac{1-2{\mathrm{\γ}}^2\left(L\right)}{v_p\left(L\right)}\right)$ $W_1^s\left(L\right)=p_2^s\left(L\right)e_1^s\left(L\right)+p_1^s\left(L\right)e_2^s\left(L\right),$ $W_2^s\left(L\right)=-2p_1^s\left(L\right)e_1^s\left(L\right),$ $\mathrm{\γ}\left(L\right)=\frac{v_s\left(L\right)}{v_p\left(L\right)}$

(11)。

3.根据权利要求1所述的弹性波高斯束叠前深度偏移技术，其特征在于步骤(4)中，根据Claerbout成像原理，首先将震源位移波场

通过弹性动力学高斯束来表示

$U_m^p(x;x_s;\mathrm{\ω})\≈\frac{i}{4\mathrm{\π}{v}_p^2\left(x_s\right)}\sqrt{\frac{{\mathrm{\ω}}_r{w}_0^2}{\mathrm{\ρ}\left(x_s\right)}}\∫\frac{d{p}_1^p\left(x_s\right)}{p_2^p\left(x_s\right)}{\hat{u}}_m^p(x;x_s;\mathrm{\ω})---\left(12\right)$

接下来根据P波以及S波的传播特点，结合式(8)，式(9)所表示的反向延拓的接收波场，定义如下成像公式

$I^{\mathrm{pp}}\left(x\right)={\∫U}_2^{p*}(x;x_s;\mathrm{\ω})u_2^p(x;x_r;\mathrm{\ω})\mathrm{d\ω}$

$=\frac{\mathrm{\ΔL}\sqrt{{\mathrm{\ω}}_r{w}_0^2}}{16{\mathrm{\π}}^2}\underset{L}{\mathrm{\Σ}}\∫\mathrm{d\ω}\frac{\mathrm{i\ω}}{v_p^2\left(x_s\right)}\sqrt{\frac{\mathrm{\ρ}\left(L\right)}{\mathrm{\ρ}\left(x_s\right)}}\∫\∫\frac{d{p}_1^p\left(x_s\right)d{p}_1^p\left(L\right)}{p_2^p\left(x_s\right)p_2^p\left(L\right)}---\left(13\right)$

$\×{\hat{u}}_2^{p*}(x;x_s;\mathrm{\ω}){\hat{u}}_1^{p*}(x;L;\mathrm{\ω})[W_1^p\left(L\right)D_1^p(L;p_1^p;\mathrm{\ω})+W_2^p\left(L\right)D_2^p(L;p_1^p;\mathrm{\ω})]$

$I^{\mathrm{ps}}\left(x\right)=\frac{\mathrm{\ΔL}\sqrt{{\mathrm{\ω}}_r{w}_0^2}}{16{\mathrm{\π}}^2}\underset{L}{\mathrm{\Σ}}\∫\mathrm{d\ω}\frac{\mathrm{i\ω}}{v_p^2\left(x_s\right)}\sqrt{\frac{\mathrm{\ρ}\left(L\right)}{\mathrm{\ρ}\left(x_s\right)}}\∫\∫\frac{d{p}_1^p\left(x_s\right)d{p}_1^s\left(L\right)}{p_2^p\left(x_s\right)p_2^s\left(L\right)}\mathrm{sgn}\left(\mathrm{\α}\right)---\left(14\right)$

$\×{\hat{u}}_2^{p*}(x;x_s;\mathrm{\ω}){\hat{u}}_1^{s*}(x;L;\mathrm{\ω})[W_1^s\left(L\right)D_1^s(L;p_1^s;\mathrm{\ω})+W_2^s\left(L\right)D_2^s(L;p_1^s;\mathrm{\ω})]$

其中，I^pp(x)为P-P单炮成像值，I^ps(x)为P-S单炮成像值，将所有炮记录进行偏移并叠加，便得到最终的弹性波成像结果。

Images