应用阈值控制的弹性波高斯束偏移方法

2023-11-26孟繁琨李振春付继有张凯徐学成刘强

石油地球物理勘探 2023年5期

孟繁琨，李振春*，付继有，张凯，徐学成，刘强

（1. 中国石油大学（华东），山东青岛 266580；2. 山东省深层油气重点实验室，山东青岛 266580；3. 中国石油长庆油田页岩油开发分公司，甘肃庆阳 745000）

0 引言

随着勘探地球物理技术的进步，地震勘探由传统的单分量声波勘探向多波、多分量勘探发展，勘探对象也逐渐转向复杂目标。弹性波地震勘探可有效利用纵、横波信息，降低勘探结果的多解性，且有助于提高地震勘探分辨率。为此，国内众多学者从弹性波正演模拟、速度建模和偏移成像等方面展开大量研究，深化了对弹性波传播规律的认识［1-7］，尤其在弹性波偏移成像领域取得了很多成果。然而，针对弹性波偏移成像方法的研究大多是基于常规油气藏，而对于地震资料信噪比较低的探区，常规偏移成像方法的应用效果具有一定局限性。因此，探究一种适用于低信噪比数据的弹性波高斯束偏移方法尤为重要。

弹性波偏移方法包括波动方程类和射线类两种。Sun 等［8］通过研究弹性波波动方程，于1986 年首先提出一种波动方程类弹性波偏移方法；经过Jia等［9］的发展，形成对观测系统适应性更强的各向同性介质角度域弹性波逆时偏移（RTM）方法。RTM 等偏移成像类方法虽然具有较高成像精度，但在计算效率方面欠理想。与此同时，Kuo 等［10］与Pao 等［11］研究了一种弹性波Kirchhoff 积分法偏移成像方法。Kirchhoff 类弹性波偏移具有较高计算效率，然而常规射线类方法受到多值走时等问题影响，且无法解决焦散区和阴影区的成像问题，在应用中存在一定局限性。为此，Červen ý 等［12］将高斯束方法引入地球物理学领域，在Hill［13-14］、Gray 等［15］努力下，实现了高斯束偏移成像方法。Nowack 等［16］将高斯束偏移方法引入共炮域，提高了其对观测系统的适应性。作为一种改进的Kirchhoff偏移方法，高斯束偏移方法不仅具有较高计算效率，且能有效克服常规射线类偏移方法的局限性，解决了焦散区和阴影区的成像问题。后来，Vinje 等［17］、吴建文［18］在常规声波高斯束偏移基础上，提出基于Tau-p域阈值滤波的声波控制束偏移成像方法；杨晶等［19］进一步研究了时间域声波高斯束偏移成像方法。

然而以上研究均基于声波假设，无法有效解决弹性波偏移成像问题。在声波高斯束偏移发展的同时，Katchalov等［20］、岳玉波［21］等进一步研究了多波多分量弹性波偏移成像方法；秦宁［22］通过研究各向异性介质射线追踪算法，发展了弹性波各向异性高斯束逆时偏移方法。只是以上研究均采用常规弹性波高斯束偏移成像思路，未提及复杂油气藏勘探中低信噪比地震资料成像问题。

本文基于前人研究成果，将阈值控制压制噪声方法引入常规弹性波高斯束偏移，提出一种基于阈值控制的弹性波高斯束偏移成像优化方法，通过滤除低波场值噪声，达到提高剖面成像质量的目的。文中采用洼陷模型和断块模型进行测试，试算结果验证了本文方法的正确性和有效性。相较于常规弹性波高斯束偏移成像方法，本文方法可有效压制噪声，增大偏移剖面信噪比，提高低信噪比地震资料成像质量。

1 阈值控制的弹性波高斯束偏移基本原理

运动学射线追踪与动力学射线追踪是实现高斯束偏移的关键。通过运动学射线追踪计算中心射线路径和走时，以动力学射线追踪方程计算动力学射线追踪参数。在射线追踪基础上计算得到高斯束，进而表征格林函数，推导出波场正向延拓和反向延拓公式，利用互相关条件实现弹性波高斯束偏移成像。

1.1 运动学射线追踪方程

本文采用岳玉波［21］给出的各向同性介质弹性波运动学射线追踪方程，震源处采用纵波射线追踪，而检波点处束中心分别采用纵波和横波做射线追踪，得到中心射线走时和路径信息

式中：xi为直角坐标系（x，z）中的位置；pPi、pSi分别为纵（P）、横（S）波慢度分量；τ为沿射线方向传播的走时；VP、VS分别为地下介质的纵、横波速度；s为二维射线中心坐标系（图1）中射线Ω上某点到选定参考点S的弧长。

1.2 动力学射线追踪方程

动力学射线追踪广泛应用于高频近似下的波场计算与地震反演中［23］。本文利用动力学射线追踪计算动力学射线追踪参数，用于中心射线的振幅及波前计算

上述式中：P（s）和Q（s）为动力学射线追踪参数；n代表射线Ω附近一点到S点的距离。

1.3 阈值控制函数

本文将阈值函数引入Tau-p变换，修改了倾斜叠加变换公式，并将其应用于偏移成像过程，实现基于阈值控制的弹性波高斯束偏移成像优化方法。

已知不同波型的多分量地震记录的加窗倾斜叠加为

式中：xs、xr对应为震源和接收点位置；um（xr，ω）为xr处接收到的弹性波地震记录，m=1、2，对应为x、z方向的分量；上标ν指P 波和S 波两种不同波型；ωref、w0分别为高斯束参考频率和初始宽度；为在束中心xL处ν型波慢度矢量的x方向分量。

通过在倾斜叠加过程中引入阈值参数ε，得到基于阈值控制的倾斜叠加公式

其中，阈值参数ε的选取需满足以下条件

式中：W为阈值控制系数，可据不同地震数据做相应赋值；为波场值算术平均值。

1.4 弹性波波场正、反向延拓公式

在二维射线中心坐标系（图1）中，可构建由参考点x0处出射，且经过计算空间任一点x处的高斯束位移公式

式中：φν为复值常数；ρ（s）为介质的密度数值；vν（s）为不同波型的相速度；eP为P 波在x处的高斯束极化矢量，表达式为其中Ω的切向矢量t为高斯束极化矢量的主分量，Ω的法向矢量n为次分量；同理，eS为S 波高斯束极化矢量，表达式为但与P 波不同，S 波中n为主分量，－t为次分量。

1.4.1 波场的正向延拓

根据式（8），将由xs处以不同角度出射，且对x位移有作用的高斯束叠加起来，便可得震源处正向延拓公式

1.4.2 阈值控制的波场反向延拓

在xs处激发、xr处接收到的记录用ui（xr，ω）表示，便可将弹性波反向波场的位移um（x，xs，ω）用Kirchhoff-Helmholtz积分表示为

式中：上标*表示复共轭，在xr处对其沿i方向施加单位体力，造成点x处沿m方向产生一定位移量，用格林张量表示；为xr处张量；nj为xr处沿外法线方向的单位矢量，应力格林张量其中Cijkl为四阶刚度系数，在二维各向同性介质弹性波中表示为

式中：λ（x）、μ（x）为拉梅弹性参数，满足λ（x）+2μ（x）为Kronecker Delta函数。

据文献［21］，可得基于阈值控制的反向延拓位移公式

1.5 阈值控制的弹性波高斯束偏移成像公式

根据Clearbout［24］成像法则，利用震源波场与接收点反向延拓的波场之间的互相关计算，可得到单炮的成像值。由式（9）、式（12）可得到阈值控制的PP 波与PS波成像公式

1.6 技术实现流程

阈值控制的弹性波高斯束偏移方法的具体实现流程如图2所示。

图2 阈值控制的弹性波高斯束偏移流程

2 模型试算

为了验证该方法在低信噪比地震数据中的适用性，本文利用洼陷模型和断块模型进行测试。

2.1 洼陷模型试算

设计的洼陷模型（图3）大小为1801×301，纵向网格间距为10 m，横向网格间距为10 m，炮点距为7 m，第一炮位于10 m 处，总计250 炮激发，每炮接收道数为1801，道间距为10 m，时间采样间隔为1 ms，采集时长为3001 ms。

图3 洼陷模型P 波（a）和S 波（b）速度场

为验证低信噪比地震记录对偏移结果的影响及基于阈值控制的弹性波高斯束偏移方法的有效性，在SU 系统中利用suaddnoise 指令对正演得到的炮记录做加噪处理，其指令原理如下［25］

式中：O（x，t）为输出地震记录；I（x，t）为输入地震记录；GN为高斯噪声；I（x，t）*为输入地震记录绝对值的平均值；|signal|max为输入地震记录最大绝对值；s（x，t）为每个记录采样点值；σ为输入记录的标准差。

信噪比（SNR）与噪声量（noise）的关系为［26］

式中：s（x，t）为地震记录；rand（x，t）为（x，t）点的随机数；σrand为所有采样点随机数的标准差。

从图4 所示地震记录可知：原始地震记录上的反射波信息较清晰，同相轴较收敛；当SNR 为5时，地震记录中的有效反射波信息受到噪声干扰（图4c、图4d）。针对这些低信噪比地震记录用常规方法做偏移成像（图5）时，会导致成像剖面上出现大量噪声，影响偏移成像质量（图5c、图5e）。而采用本文方法，所得成像剖面（图5d、图5f）较清晰，有效消除了噪声干扰，提高了地层的分辨率，尤其是浅层的效果更显著。

图4 洼陷模型原始炮记录的X（a）、Z（b）分量及SNR 为5 的X（c）、Z（d）分量炮记录

图5 采用不同方法对洼陷模型进行偏移成像的结果

2.2 复杂断块模型试算

设计的复杂断块模型（图6）大小为1000×550，纵向网格间距为10 m，横向网格间距为10 m，炮间距为5 m，第一炮位于2 m 处，总计200炮激发，每炮接收道数为1000，道间距为10 m，时间采样间隔为1 ms，采集时长为5501 ms。模型正演及加噪后的地震记录分别如图7所示。对比所得成像结果（图8）可见：采用常规弹性波高斯束偏移对低信噪比数据的成像剖面（图8c、图8e）上有大量噪声，有效同相轴能量被噪声干扰，地震剖面信噪比较低，整体成像质量较差；用本文方法的偏移剖面（图8d、图8f）上噪声明显减少，同相轴能量也更聚焦，因此该方法适用于低信噪比数据的高质量成像。