留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

逆时偏移波场延拓中的几何扩散问题

方刚 杜启振 栾锡武

方刚, 杜启振, 栾锡武. 逆时偏移波场延拓中的几何扩散问题[J]. 海洋地质与第四纪地质, 2017, 37(1): 162-167. doi: 10.16562/j.cnki.0256-1492.2017.01.019
引用本文: 方刚, 杜启振, 栾锡武. 逆时偏移波场延拓中的几何扩散问题[J]. 海洋地质与第四纪地质, 2017, 37(1): 162-167. doi: 10.16562/j.cnki.0256-1492.2017.01.019
FANG Gang, DU Qizhen, LUAN Xiwu. THE GEOMETRY SPREADING OF WAVEFIELD EXTRAPOLATION FOR REVERSE TIME MIGRATION[J]. Marine Geology & Quaternary Geology, 2017, 37(1): 162-167. doi: 10.16562/j.cnki.0256-1492.2017.01.019
Citation: FANG Gang, DU Qizhen, LUAN Xiwu. THE GEOMETRY SPREADING OF WAVEFIELD EXTRAPOLATION FOR REVERSE TIME MIGRATION[J]. Marine Geology & Quaternary Geology, 2017, 37(1): 162-167. doi: 10.16562/j.cnki.0256-1492.2017.01.019

逆时偏移波场延拓中的几何扩散问题


doi: 10.16562/j.cnki.0256-1492.2017.01.019
详细信息
    作者简介:

    方刚(1984—), 男, 博士, 助理研究员,研究方向为地震波传播理论与方法、地震成像,E-mail: fanggeo@163.com

  • 基金项目:

    国家自然科学基金项目 41504109

    中央高校基本科研业务费专项项目 11CX06002A

    山东省自然科学基金项目 BS2015HZ008

  • 中图分类号: P315

THE GEOMETRY SPREADING OF WAVEFIELD EXTRAPOLATION FOR REVERSE TIME MIGRATION

More Information
图(4)
计量
  • 文章访问数:  392
  • HTML全文浏览量:  66
  • PDF下载量:  2
  • 被引次数: 0
出版历程
  • 收稿日期:  2016-03-14
  • 修回日期:  2016-07-18
  • 刊出日期:  2017-02-28

逆时偏移波场延拓中的几何扩散问题

doi: 10.16562/j.cnki.0256-1492.2017.01.019
    作者简介:

    方刚(1984—), 男, 博士, 助理研究员,研究方向为地震波传播理论与方法、地震成像,E-mail: fanggeo@163.com

基金项目:

国家自然科学基金项目 41504109

中央高校基本科研业务费专项项目 11CX06002A

山东省自然科学基金项目 BS2015HZ008

  • 中图分类号: P315

摘要: 从理论上分析逆时偏移波场延拓过程中的振幅变化有助于发展真振幅的逆时偏移方法,提高成像结果的精度和保真性。本文研究逆时偏移波场延拓中的几何扩散对振幅的影响,首先基于高频渐进理论,给出正向和逆时延拓波场的WKBJ近似表示;然后应用稳相原理,证明逆时偏移中逆向波场延拓能够补偿正向波场延拓中的几何扩散效应;最后通过数值算例验证理论推导的结果,分析逆时偏移自动补偿几何扩散后的成像振幅。该研究为真振幅逆时偏移方法研究提供理论基础。

English Abstract

  • 逆时偏移(Reverse time migration,RTM)[1-3]采用双程波动方程进行时间域的波场延拓构建成像所需的地震波场,能够处理复杂介质中传播的地震波(如多次波、棱柱波、回转波等),在成像中没有倾角限制,是目前成像精度最高的方法之一[4]。逆时偏移中的波场延拓决定了逆时偏移的成像精度和保幅性,波场延拓包括震源波场的正向延拓和检波波场的逆时延拓。

    真振幅成像从地震记录中恢复出地下真实反射系数的像,是AVO(振幅随偏移距变化)或AVA(振幅随角度变化)分析的基础[5]。以真振幅成像为目标,要求逆时偏移中的波场延拓同时保持地震波场传播的运动学和动力学特征。几何扩散、透射损失和吸收衰减是导致地震波传播过程中振幅失真的主要因素。Deng和McMechan[6, 7]、Du等[8]考虑地下强反射界面对地震振幅的影响,研究了针对逆时偏移补偿透射损失的方法。真实的地下介质具有粘滞性,会导致地震波在传播过程中发生振幅衰减和相位畸变。许多学者在地震波场延拓[9, 10]和逆时偏移[11, 12]的研究中考虑介质的粘滞性,补偿由吸收衰减造成的振幅衰减并校正相位的畸变,发展振幅保真的偏移成像方法[13]。波场延拓的几何扩散与所采用的方程有关,Zhang等[14]证明对于单程波方程[15]由于采用了近似的伪微分算子,用其进行波场延拓是不保幅的,随后基于高频渐进理论[16],Zhang等[17]给出了真振幅的单程波方程以及基于此方程的偏移成像方法。随着计算机性能的提高,计算量较大但精度更高的双程波方程逐渐取代单程波方程被用于波场延拓,通常认为采用双程波方程的波场延拓是保幅的,但这种“保幅性”还缺少理论上的证明,尤其是逆时偏移中波场延拓的几何扩散对成像振幅影响的理论分析还少有研究。

    本文基于高频渐进理论,分析基于双程波动方程的逆时偏移波场延拓中的几何扩散问题。首先给出基于双程波动方程正向和逆时延拓波场的高频渐进表示,然后通过稳相原理证明逆时波场延拓中的几何扩散能自动的补偿,最后利用数值算例验证该结论。

    • 考虑基于双程波动方程的震源波场正向延拓,在速度和密度可变的情况下,震源波场延拓满足方程

      $ \begin{array}{l} \left[ {\rho ({\bf{x}})\nabla \cdot \frac{1}{{\rho ({\bf{x}})}}\nabla + \frac{{{\omega ^2}}}{{{v^2}({\bf{x}})}}} \right]{P_F}({\bf{x}},\omega ) = \\ - F(\omega )\delta \left( {{\bf{x}} - {{\bf{x}}_s}} \right) \end{array} $

      (1)

      其中, PF(x, ω)表示由xs处震源激发,在x点处得到的频率为ω的地震波场,F(ω)表示震源的信号,ρ(x)表示介质的密度, v(x)表示波速。基于该方程可构建波场正向延拓算子,实现震源波场的正向外推。在高频条件下,正传波场的WKBJ格林函数为:

      $ G\left( {{\bf{x}},{{\bf{x}}_s},\omega } \right) \sim B\left( {{\bf{x}},{{\bf{x}}_s}} \right){{\rm{e}}^{i\omega \tau \left( {{\bf{x}},{{\bf{x}}_s}} \right)}} $

      (2)

      其中,G(x, xs, ω)表示xs处的点源在x处的地震响应;τ(x, xs)表示xs处的波场传播至x所需的时间,满足程函方程;B(x,xs)是方程(1)所描述波场的WKBJ首阶振幅,满足变密度的输运方程。Bleistein等[16]指出B(x, xs)与常密度情况下波场的WKBJ首阶振幅A(x, xs)存在如下关系:

      $ B\left( {{\bf{x}},{{\bf{x}}_s}} \right) = A\left( {{\bf{x}},{{\bf{x}}_s}} \right)\sqrt {\frac{{\rho \left( {{{\bf{x}}_s}} \right)}}{{\rho ({\bf{x}})}}} $

      (3)

      对于逆时延拓,以接收波场为边界条件,检波波场的逆时延拓满足方程

      $ \left\{ {\begin{array}{*{20}{l}} {\left[ {\rho ({\bf{x}})\nabla \cdot \frac{1}{{\rho ({\bf{x}})}}\nabla + \frac{{{\omega ^2}}}{{v{{({\bf{x}})}^2}}}} \right]{P_B}({\bf{x}},\omega ) = 0}\\ {{P_B}\left( {{{\bf{x}}_0},\omega } \right) = Q\left( {{{\bf{x}}_0},\omega } \right),{{\bf{x}}_0} = ({\bf{x}},y,0)} \end{array}} \right. $

      (4)

      其中x0表示地表检波器的位置,Q(x0, ω)表示接收到的地震记录。基于该方程可构建逆时延拓算子,逆时延拓波场的WKBJ非因果格林函数为:

      $ {G^*}\left( {{\bf{x}},{{\bf{x}}_s},\omega } \right) \sim B\left( {{\bf{x}},{{\bf{x}}_s}} \right){{\rm{e}}^{ - i\omega \tau \left( {{\bf{x}},{{\bf{x}}_s}} \right)}} $

      (5)

      其中,*表示复共轭,其相位和振幅与(2)式类似,分别满足程函方程和变密度的输运方程。

      考虑图 1所示的逆时偏移波场延拓,S0为地表观测面,S1为反射界面,实线和虚线分别表示正传和逆传波场的传播路径。波场的正向延拓过程中,由xs处的震源激发的波场在xS1点的入射波场为PFI(xS1, ω),该点的反射波场为PFR(xS1, ω), 根据高频渐近理论,震源波场可以表示为[16]

      $ {{P}_{F}}\left( \text{x},\omega \right)\sim A\left( \text{x},{{\text{x}}_{S}} \right)\sqrt{\frac{\rho \left( \text{x} \right)}{\rho \left( {{\text{x}}_{S}} \right)}}{{\text{e}}^{i\omega \tau \left( \text{x},{{\text{x}}_{S}} \right)}} $

      (6)

      图  1  逆时偏移波场延拓示意图

      Figure 1.  Wave propagation in reverse time migration

      在高频近似下,经反射面S1反射的波场可以表示为WKBJ格林函数沿着反射界面的积分[18]

      $ \begin{array}{l} P_F^R({\bf{x}},\omega ) \sim 2i\omega \int R \left( {{{\bf{x}}_{{S_1}}}} \right)P_F^I\left( {{{\bf{x}}_{{S_1}}},\omega } \right)A\left( {{{\bf{x}}_{{S_1}}},{\bf{x}}} \right)\\ \sqrt {\frac{{\rho \left( {{{\bf{x}}_{{S_1}}}} \right)}}{{\rho ({\bf{x}})}}} n \cdot \nabla \tau \left( {{{\bf{x}}_{{S_1}}},{\bf{x}}} \right){{\rm{e}}^{i\omega \tau \left( {{{\bf{x}}_{{S_1}}},{\bf{x}}} \right)}}{\rm{d}}{{\bf{x}}_{{S_1}}} \end{array} $

      (7)

      其中,R(xS1)表示界面S1的反射系数。

      对于逆时偏移中的波场逆时延拓过程,检波器接收到的检波波场为Q(xS0, ω),根据非因果格林函数(5)式,逆向传播的检波波场由下式给出[18]

      $ \begin{array}{l} {P_B}({\bf{x}},\omega ) \sim 2i\omega \int Q \left( {{{\bf{x}}_{{S_0}}},\omega } \right)A\left( {{{\bf{x}}_{{S_0}}},{\bf{x}}} \right)\sqrt {\frac{{\rho ({\bf{x}})}}{{\rho \left( {{{\bf{x}}_{{S_0}}}} \right)}}} {n_0}\\ \cdot \nabla \tau \left( {{{\bf{x}}_{{S_0}}},{\bf{x}}} \right){{\rm{e}}^{ - i\omega \tau \left( {{{\bf{x}}_{{S_0}}},{\bf{x}}} \right)}}{\rm{d}}{{\bf{x}}_{{S_0}}} \end{array} $

      (8)

      式(7)和式(8)分别为震源波场和检波波场的高频渐近表达式。

    • 几何扩散是地震波场传播过程中的一种重要的传播效应,它是指地震波场传播中波前面不断扩大,单位面积上的地震能量逐渐减少,地震振幅相应逐渐减弱的现象。为了从理论上推导这种传播效应所导致的振幅变化,在推导过程中采用准确的介质模型。所谓“准确模型”是相对于地震偏移中所用的光滑模型(介质参数渐变而未知,介质交界面模糊不清)而言,指的是所需要的介质参数是已知的,介质交界面清晰。实际中,接收波场Q(xS0, ω)的振幅受到粘滞性、几何扩散、透射损失、采集效应等因素的综合影响,这里仅考虑几何扩散对振幅的影响,假设介质是均匀无粘性的单层介质,忽略介质、采集孔径等其他因素的影响。

      逆时偏移只有当震源波场和检波波场在成像点处的振幅、相位都准确时,应用成像条件才能得到真振幅的成像结果[15]。因为逆时偏移采用双程波动方程构建波场正传算子,当所用的模型是准确模型时,双程波动方程得到的解的相位和振幅都是准确的, 震源波场正向传播中存在的几何扩散正是对物理过程的准确描述。检波波场的逆向传播也是采用双程波动方程,但这个过程是一个物理不可实现的非因果过程。逆时偏移通常将接收到的波场Q(xS0, ω)沿时间逆向传播,得到界面S1处的检波波场PBI(xS1, ω)(如图 1所示),在成像条件中用该检波波场代替反射波场PFR(xS1, ω)以获得最终的成像结果。只有当检波波场PBI(xS1, ω)与反射波场PFR(xS1, ω)的振幅相等,才能保证成像结果的保幅性。下面从理论上分析二者振幅的关系,给出检波波场逆向传播能够自动补偿正向传播过程中几何扩散的理论证明。

      根据(7)式,地表观测面S0上的检波器接收到的波场Q(xS0, ω)可以表示为如下积分

      $ \begin{array}{l} Q\left( {{{\bf{x}}_{{S_0}}},\omega } \right) \sim 2i\omega \int_{{S_1}} {P_F^R} \left( {{{\bf{x}}_{{S_1}}},\omega } \right)A\left( {{{\bf{x}}_{{S_1}}},{{\bf{x}}_{{S_0}}}} \right)\\ \sqrt {\frac{{\rho \left( {{{\bf{x}}_{{S_0}}}} \right)}}{{\rho \left( {{{\bf{x}}_{{S_1}}}} \right)}}} {n_{{S_1}}} \cdot {\nabla _\tau }\left( {{{\bf{x}}_{{{\rm{S}}_1}}},{{\bf{x}}_{{S_0}}}} \right){{\rm{e}}^{i\omega \tau \left( {{{\bf{x}}_{{S_1}}},{{\bf{x}}_{{S_0}}}} \right)}}{\rm{d}}{{\bf{x}}_{{S_1}}} \end{array} $

      (9)

      (9) 式实际上表示了用因果格林函数(2)式将界面xS1上的反射波场PFR(xS1, ω)延拓至xS0上的过程。

      根据(8)式,可将接收波场Q(xS0, ω)逆向延拓得到任意点x处的检波波场

      $ \begin{array}{l} P_B^I({\bf{x}},\omega ) \sim 2i\omega \int_{{S_0}} Q \left( {{{\bf{x}}_{{S_0}}},\omega } \right)A\left( {{{\bf{x}}_{{S_0}}},{\bf{x}}} \right)\\ \sqrt {\frac{{\rho ({\bf{x}})}}{{\rho \left( {{{\bf{x}}_{{S_0}}}} \right)}}} {n_{{S_0}}} \cdot \nabla \tau \left( {{{\bf{x}}_{{S_0}}},{\bf{x}}} \right){{\rm{e}}^{ - i\omega \tau \left( {{{\bf{x}}_{{S_0}}},{\bf{x}}} \right)}}{\rm{d}}{{\bf{x}}_0} \end{array} $

      (10)

      x点处的三维空间坐标为(x, y, z)。将(9)式代入(10)得

      $ \begin{array}{l} P_B^I({\bf{x}},\omega ) \sim {(2i\omega )^2}\int_{{S_0}} {\int_{{S_1}} {P_F^R} } \left( {{{\bf{x}}_{{S_1}}},\omega } \right)\\ \sqrt {\frac{{\rho ({\bf{x}})}}{{\rho \left( {{{\bf{x}}_{{S_1}}}} \right)}}} A\left( {{{\bf{x}}_{{S_1}}},{{\bf{x}}_{{S_0}}}} \right)A\left( {{{\bf{x}}_{{S_0}}},{\bf{x}}} \right){n_{{S_1}}} \cdot {\nabla _\tau }\left( {{{\bf{x}}_{{S_1}}},{{\bf{x}}_{{S_0}}}} \right){n_{{S_0}}} \cdot \\ \nabla \tau \left( {{{\bf{x}}_{{S_0}}},{\bf{x}}} \right){{\rm{e}}^{i\omega \left[ {\tau \left( {{{\bf{x}}_{{S_1}}},{{\bf{x}}_{{S_0}}}} \right) - \tau \left( {{{\bf{x}}_{{S_0}}},{\bf{x}}} \right)} \right]}}{\rm{d}}{{\bf{x}}_{{S_1}}}{\rm{d}}{{\bf{x}}_{{S_0}}} \end{array} $

      (11)

      根据文献[16],式(11)中格林函数的振幅和相位可以由下列式给出

      $ \begin{array}{l} A\left( {{{\bf{x}}_{{s_0}}},{{\bf{x}}_{{s_1}}}} \right) = \frac{1}{{4{\rm{ \mathsf{ π} }}\left| {{{\bf{x}}_{{s_1}}} - {{\bf{x}}_{{s_0}}}} \right|}},\tau \left( {{{\bf{x}}_{{s_1}}},{{\bf{x}}_{{s_0}}}} \right) = \\ \frac{{\left| {{{\bf{x}}_{{s_1}}} - {{\bf{x}}_{{s_0}}}} \right|}}{v},{n_{{s_1}}} \cdot \nabla \tau \left( {{{\bf{x}}_{{s_0}}},{{\bf{x}}_{{s_1}}}} \right) = \frac{{\cos {\theta _{{s_1}}}}}{v} \end{array} $

      $ \begin{array}{l} A\left( {{\bf{x}},{{\bf{x}}_{{s_0}}}} \right) = \frac{1}{{4{\rm{ \mathsf{ π} }}\left| {{{\bf{x}}_{{s_0}}} - {\bf{x}}} \right|}},\tau \left( {{{\bf{x}}_{{s_0}}},{\bf{x}}} \right) = \frac{{\left| {{{\bf{x}}_{{s_0}}} - {\bf{x}}} \right|}}{v},\\ {n_{{s_0}}} \cdot \nabla \tau \left( {{{\bf{x}}_{{s_0}}},{{\bf{x}}_{{s_1}}}} \right) = \frac{{\cos {\theta _{{s_0}}}}}{v} \end{array} $

      其中,v为波在介质中的传播速度,θS1为正传播场的传播方向与界面S1之法线方向之间的夹角,θS0为逆传播场的传播方向与界面S0法线方向之间的夹角。根据格林定理,S0的法方向是指向界面下方,而S1的法方向指向界面上方。利用上述形式,(11)式可表示为:

      $ \begin{array}{l} P_B^I({\bf{x}},\omega ) \sim {\left( {\frac{{i\omega }}{{2{\rm{ \mathsf{ π} }}}}} \right)^2}\int_{{S_0}} {\int_{{S_1}} {P_F^R} } \left( {{{\bf{x}}_{{S_1}}},\omega } \right)\sqrt {\frac{{\rho ({\bf{x}})}}{{\rho \left( {{{\bf{x}}_{{S_{\rm{1}}}}}} \right)}}} \\ \frac{1}{{\left| {{{\bf{x}}_{{S_1}}} - {{\bf{x}}_{{S_0}}}} \right|}}\\ \frac{1}{{\left| {{{\bf{x}}_{{S_0}}} - {\bf{x}}} \right|}}\frac{{\cos {\theta _{{S_1}}}}}{v}\frac{{\cos {\theta _{{S_0}}}}}{v}{{\rm{e}}^{i\omega \left[ {\frac{{\left| {{{\bf{x}}_{{S_1}}} - {{\bf{x}}_{{S_0}}}} \right|}}{v} - \frac{{\left| {{{\bf{x}}_{{S_0}}} - {\bf{x}}} \right|}}{v}} \right]}}{\rm{d}}{{\bf{x}}_{{S_1}}}{\rm{d}}{{\bf{x}}_{{S_0}}} \end{array} $

      (12)

      由式(12)难以直接分析波场PBI(x,ω)与PFR(xS1, ω)的振幅间的关系,下面应用稳相法对式进行稳相点分析。

      稳相点分析[16]可以用于估算如下的n维震荡(Fourier型)积分,

      $ I(\omega ) = \int_D f ({\bf{x}}){{\rm{e}}^{i\omega \varphi ({\bf{x}})}}{\rm{d}}{\bf{x}},{\bf{x}} = \left( {{{\bf{x}}_1}, \cdots ,{{\bf{x}}_n}} \right) $

      (13)

      当频率ω较大时,该积分可以在稳相点附近作渐近展开,得到如下渐近近似

      $ I(\omega ) \sim {\left( {\frac{{2{\rm{ \mathsf{ π} }}}}{\omega }} \right)^{\frac{m}{2}}}\frac{{f({\bf{\bar x}})}}{{\sqrt {\left| {detA} \right|} }}{{\rm{e}}^{i\omega \varphi ({\bf{\bar x}}) + i\frac{{\rm{ \mathsf{ π} }}}{4}{\rm{sgn}}A}} $

      (14)

      其中,x是稳相点,满足稳相条件

      $ \nabla \varphi ({\bf{\bar x}}) = 0 $

      (15)

      并且

      $ detA \ne 0,{A_{jk}} = \left[ {\frac{{{\partial ^2}\varphi ({\bf{\bar x}})}}{{\partial {x_j}\partial {x_k}}}} \right],j,k = 1,2, \cdots ,n $

      (16)

      这里Aφ的海森矩阵;sgnA表示矩阵A的正特征值的个数与负特征值的个数之差。xS0xS1为界面S0S1上的稳相点,对(12)作稳相点渐近展开有

      $ \begin{array}{l} P_B^I({\bf{x}},\omega ) \sim {\left( {\frac{{i\omega }}{{2{\rm{ \mathsf{ π} }}}}} \right)^2}{\left( {\frac{{2{\rm{ \mathsf{ π} }}}}{\omega }} \right)^2}\frac{{P_F^R\left( {{{{\bf{\bar x}}}_{{S_1}}},\omega } \right)}}{{\sqrt {\left| {detA} \right|} }}\sqrt {\frac{{\rho ({\bf{x}})}}{{\rho \left( {{{{\bf{\bar x}}}_{{S_1}}}} \right)}}} \\ \frac{1}{{\left| {{{{\bf{\bar x}}}_{{S_1}}} - {{{\bf{\bar x}}}_{{S_0}}}} \right|}}\frac{1}{{\left| {{{{\bf{\bar x}}}_{{S_0}}} - {\bf{x}}} \right|}}\frac{{\cos {\theta _{{S_1}}}\cos {\theta _{{S_0}}}}}{{{v^2}}}\\ {{\rm{e}}^{i\omega \left[ {\frac{{\left| {{{\mathit{\boldsymbol{\bar x}}}_{{S_1}}} - {{\mathit{\boldsymbol{\bar x}}}_{{S_0}}}} \right|}}{v} - \frac{{\left| {{{\mathit{\boldsymbol{\bar x}}}_{{S_0}}} - \mathit{\boldsymbol{x}}} \right|}}{v}} \right] + i\frac{{\rm{ \mathsf{ π} }}}{4}sgnA}} \end{array} $

      (17)

      由式(17)可知相函数

      $ \varphi \left( {{{{\bf{\bar x}}}_{{S_1}}} - {{{\bf{\bar x}}}_{{S_0}}}} \right) = \frac{{\left| {{{{\bf{\bar x}}}_{{S_1}}} - {{{\bf{\bar x}}}_{{S_0}}}} \right|}}{v} - \frac{{\left| {{{{\bf{\bar x}}}_{{S_0}}} - {\bf{x}}} \right|}}{v} $

      (18)

      φ的一阶导数为

      $ \frac{{\partial \varphi \left( {{{{\bf{\bar x}}}_{{S_1}}},{{{\bf{\bar x}}}_{{S_0}}}} \right)}}{{\partial {{{\bf{\bar x}}}_{{S_0}}}}} = \frac{1}{v}\left( {\frac{{{{{\bf{\bar x}}}_{{S_0}}} - {\bf{x}}}}{{\left| {{{{\bf{\bar x}}}_{{S_0}}} - {\bf{x}}} \right|}} + \frac{{{{{\bf{\bar x}}}_{{S_0}}} - {{{\bf{\bar x}}}_{{S_1}}}}}{{\left| {{{{\bf{\bar x}}}_{{S_1}}} - {{{\bf{\bar x}}}_{{S_0}}}} \right|}}} \right) $

      (19)

      $ \frac{{\partial \varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{\partial {{{\bf{\bar x}}}_{{S_1}}}}} = \frac{1}{v}\frac{{{{{\bf{\bar x}}}_{{S_1}}} - {{{\bf{\bar x}}}_{{S_0}}}}}{{\left| {{{{\bf{\bar x}}}_{{S_1}}} - {{{\bf{\bar x}}}_{{S_0}}}} \right|}} $

      (20)

      $ \frac{{\partial \varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{\partial {{\bar y}_{{S_0}}}}} = \frac{1}{v}\left( {\frac{{{{\bar y}_{{S_0}}} - y}}{{\left| {{{{\bf{\bar x}}}_{{S_0}}} - {\bf{x}}} \right|}} + \frac{{{{\bar y}_{{S_0}}} - {{\bar y}_{{S_1}}}}}{{\left| {{{{\bf{\bar x}}}_{{S_1}}} - {{{\bf{\bar x}}}_{{S_0}}}} \right|}}} \right) $

      (21)

      $ \frac{{\partial \varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{\partial {{\bar y}_{{S_1}}}}} = \frac{1}{v}\frac{{{{\bar y}_{{S_1}}} - {{\bar y}_{{S_0}}}}}{{\left| {{{{\bf{\bar x}}}_{{S_1}}} - {{{\bf{\bar x}}}_{{S_0}}}} \right|}} $

      (22)

      根据稳相条件(15)式,由(19)—(22)式可知

      $ x = {{{\bf{\bar x}}}_{{S_0}}} = {{{\bf{\bar x}}}_{{S_1}}},y = {{\bar y}_{{s_0}}} = {{\bar y}_{{s_1}}} $

      (23)

      则有

      $ {\theta _{{S_1}}} = {\theta _{{S_0}}} = 0 $

      (24)

      相函数φ的海森矩阵有如下的非零分量,

      $ \frac{{{\partial ^2}\varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{\partial {\bf{\bar x}}_{{S_0}}^2}} = \frac{1}{v}\left( {\frac{1}{{\left| {{{{\bf{\bar x}}}_{{S_0}}} - {\bf{x}}} \right|}} + \frac{1}{{\left| {{{{\bf{\bar x}}}_{{S_1}}} - {{{\bf{\bar x}}}_{{S_0}}}} \right|}}} \right) $

      (25)

      $ \frac{{{\partial ^2}\varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{{\partial ^2}{{{\bf{\bar x}}}_{{S_1}}}}} = \frac{1}{{v\left| {{{{\bf{\bar x}}}_{{{\rm{S}}_1}}} - {{{\bf{\bar x}}}_{{S_0}}}} \right|}} $

      (26)

      $ \frac{{{\partial ^2}\varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{{\partial ^2}{{\bar y}_{{S_0}}}}} = \frac{1}{{v\left| {{{{\bf{\bar x}}}_{{S_1}}} - {{{\bf{\bar x}}}_{{S_0}}}} \right|}} + \frac{1}{{v\left| {{{{\bf{\bar x}}}_{{S_0}}} - {\bf{x}}} \right|}} $

      (27)

      $ \frac{{{\partial ^2}\varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{{\partial ^2}{{\bar y}_{{S_1}}}}} = \frac{1}{{v\left| {{{{\bf{\bar x}}}_{{S_1}}} - {{{\bf{\bar x}}}_{{S_0}}}} \right|}} $

      (28)

      $ \frac{{{\partial ^2}\varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{\partial {{{\bf{\bar x}}}_{{S_0}}}\partial {{{\bf{\bar x}}}_{{S_1}}}}} = - \frac{1}{{v\left| {{{{\bf{\bar x}}}_{{S_1}}} - {{{\bf{\bar x}}}_{{S_0}}}} \right|}} $

      (29)

      $ \frac{{{\partial ^2}\varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{\partial {{\bar y}_{{S_0}}}\partial {{\bar y}_{{S_1}}}}} = - \frac{1}{{v\left| {{{{\bf{\bar x}}}_{{S_1}}} - {{{\bf{\bar x}}}_{{S_0}}}} \right|}} $

      (30)

      因此,海森矩阵的行列式为

      $ \sqrt {\left| {{\rm{det}}A} \right|} = \frac{1}{{{v^2}}}\frac{1}{{\left| {{{{\bf{\bar x}}}_{{S_1}}} - {{{\bf{\bar x}}}_{{S_0}}}} \right|}}\frac{1}{{\left| {{{{\bf{\bar x}}}_{{S_0}}} - {\bf{x}}} \right|}} $

      (31)

      海森矩阵A是一个4×4的矩阵,可以分解为4个2×2的子矩阵

      $ A = \left( {\begin{array}{*{20}{l}} {{D_{11}}}&{{D_{12}}}\\ {{D_{21}}}&{{D_{22}}} \end{array}} \right) $

      (32)

      其中

      $ {D_{11}} = \left( {\begin{array}{*{20}{c}} {{A_{11}}}&{{A_{12}}}\\ {{A_{21}}}&{{A_{22}}} \end{array}} \right) = \left( {\begin{array}{*{20}{c}} {\frac{{{\partial ^2}\varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{\partial {\bf{\bar x}}_{{S_0}}^2}}}&{\frac{{{\partial ^2}\varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{\partial {{{\bf{\bar x}}}_{{S_0}}}\partial {{{\bf{\bar x}}}_{{S_1}}}}}}\\ {\frac{{{\partial ^2}\varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{\partial {{{\bf{\bar x}}}_{{S_1}}}\partial {{{\bf{\bar x}}}_{{S_0}}}}}}&{\frac{{{\partial ^2}\varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{\partial {\bf{\bar x}}_{{S_1}}^2}}} \end{array}} \right) $

      (33)

      $ {D_{22}} = \left( {\begin{array}{*{20}{c}} {{A_{33}}}&{{A_{34}}}\\ {{A_{43}}}&{{A_{44}}} \end{array}} \right) = \left( {\begin{array}{*{20}{c}} {\frac{{{\partial ^2}\varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{\partial \bar y_{{S_0}}^2}}}&{\frac{{{\partial ^2}\varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{\partial {{\bar y}_{{S_0}}}\partial {{\bar y}_{{S_1}}}}}}\\ {\frac{{{\partial ^2}\varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{\partial {{\bar y}_{{S_1}}}\partial {{\bar y}_{{S_0}}}}}}&{\frac{{{\partial ^2}\varphi \left( {{{{\bf{\bar x}}}_{{S_0}}},{{{\bf{\bar x}}}_{{S_1}}}} \right)}}{{\partial \bar y_{{S_1}}^2}}} \end{array}} \right) $

      (34)

      $ {D_{12}} = \left( {\begin{array}{*{20}{c}} {{A_{13}}}&{{A_{14}}}\\ {{A_{23}}}&{{A_{24}}} \end{array}} \right) = 0,{D_{21}} = \left( {\begin{array}{*{20}{c}} {{A_{31}}}&{{A_{32}}}\\ {{A_{41}}}&{{A_{42}}} \end{array}} \right) = 0 $

      (35)

      因为A所具有的特殊结构,求解矩阵A的行列式转化为求解其两个子矩阵(33)式和(34)式的行列式。从(25)~(30)式不难发现,D11=D22,因此D11D22的特征值是成对出现的。考虑到

      $ \left| {{D_{11}}} \right| = \frac{1}{{{v^2}}}\frac{1}{{\left| {{{{\bf{\bar x}}}_{{S_1}}} - {{{\bf{\bar x}}}_{{S_0}}}} \right|}}\frac{1}{{\left| {{{{\bf{\bar x}}}_{{S_0}}} - {\bf{x}}} \right|}} > 0 $

      (36)

      D11的两个特征值的符号都是同号的,因此A的特征值中,正、负特征值的个数之差为

      $ {\rm{sgn}}A = 4 $

      (37)

      将式(24)、(31)、(37)带入式(17)可得

      $ \begin{gathered} P_B^I\left( {{\mathbf{x}},\omega } \right)\tilde{\ }P_F^R\left( {{{{\mathbf{\bar x}}}_{{S_1}}},\omega } \right)\sqrt {\frac{{\rho ({\mathbf{x}})}}{{\rho \left( {{{{\mathbf{\bar x}}}_{{S_1}}}} \right)}}} \hfill \\ {{\text{e}}^{i\omega \left[ {\frac{{\left| {{{{\mathbf{\bar x}}}_{{S_1}}} - {{{\mathbf{\bar x}}}_{{S_0}}}} \right|}}{v} - \frac{{\left| {{{{\mathbf{\bar x}}}_{{S_0}}} - {\mathbf{\bar x}}} \right|}}{v}} \right]}} \hfill \\ \end{gathered} $

      (38)

      x是界面S1上的稳相点时,即

      $ {\rm{x}} = {\overline {\rm{x}} _S}_{_1} $

      (39)

      则由式(38)可得

      $ P_B^I\left( {{{\bf{x}}_{{S_1}}},\omega } \right) \sim P_F^R\left( {{{\bf{x}}_{{S_1}}},\omega } \right) $

      (40)

      (40) 式表明界面S1xS1点上的反射波场PFR的首阶振幅与该点检波波场PBI的首阶振幅是相等的,地震波场正向传播过程中(图 1中由波场PFR传播至Q的过程)的几何扩散,在波场逆向传播过程中(图 1中由波场Q传播至PBI的过程)得到了自动的补偿。

    • 为了在理论上推导波场的振幅关系,(40)式是基于WKBJ格林函数给出的,反映了波场振幅首阶项之间的关系,这是在高频近似下的结果。在实际的逆时偏移计算中,往往采用有限差分等数值方法求解波动方程,没有高频近似的假设。下面通过一个算例来说明在一定精度的要求下,采用没有高频近似假设的方法(有限差分方法)求解得到波场振幅也满足(40)式的结果。图 2(a)为有限差分方法正演得到的由1.8 km深度处的反射界面产生的反射波场,图 2(b)是将地表接收到的反射波场逆向延拓得到的检波波场,图 2(c)是上述两地震波场的误差剖面。图 2(a)(b)中反射波场和检波波场的最大振幅分别为25.3和24.8,图 2(c)中误差的最大值为1.1,反射波场和检波波场间的相对误差小于5%。因此,在一定的精度要求下,有限差分方法得到的检波波场和反射波场也是近似相等的,也就是说波场正传过程中的几何扩散在波场的逆向传播过程中自动得到了补偿。

      图  2  逆时偏移几何扩散自动补偿的数值验证

      Figure 2.  Numerical verification of geometric spreading in RTM

      最后,基于本文给出的理论推导,我们用一个水平层状模型测试逆时偏移中的几何扩散对成像振幅的影响。对横向7 km宽,纵向3 km深,界面位于深度2 km位置处,上、下层介质的速度分别为4.5 km/s和2 km/s的水平层状模型进行逆时偏移成像。炮点位置在水平3.5 km处,采用中间放炮两边接收的观测方式,获得单炮逆时偏移成像剖面(图 3)。图 4比较了沿着地层界面提取的成像振幅和真实的反射系数,除了两端因为边界效应引起的误差外,成像振幅很好的吻合了真实的反射系数。逆时偏移成像能够获得构造位置准确的成像结果,同时在成像过程中补偿了几何扩散引起的振幅衰减,得到的成像振幅也能很好的反映反射系数随入射角的变化,具有较高的保幅性。

      图  3  水平层状模型单炮逆时偏移成像

      Figure 3.  The single shot RTM of a layers model

      图  4  逆时偏移成像振幅(虚线)与理论反射系数(实线)对比

      Figure 4.  The comparison of the RTM amplitude (dashed lines) and the theoratical reflection coefficience(solid line)

    • (1) 波场延拓中的几何扩散是导致地震成像振幅失真的主要因素之一,本文通过对逆时偏移正向波场延拓算子和逆时波场延拓算子进行高频渐进展开、应用稳相分析,给出了逆时偏移正向波场延拓过程中的几何扩散能够在逆时波场延拓中进行补偿的理论证明,为发展真振幅的逆时偏移方法奠定理论基础。

      (2) 数值算例表明,实际计算中取消高频近似的假设,采用有限差分方法进行波场延拓得到的波场振幅也符合理论推导的结论,逆时偏移补偿了几何扩散对振幅的影响,在不考虑介质的粘滞性和透射损失的情况下,逆时偏移成像振幅能够反映出地层真实反射系数的变化。

参考文献 (18)

目录

    /

    返回文章
    返回