CN115166827B - 基于反褶积成像条件的最小二乘偏移成像方法、设备及存储介质 - Google Patents

Abstract

本发明涉及地球物理勘探技术领域，为基于反褶积成像条件的最小二乘偏移成像方法、设备及存储介质，本发明通过在最小二乘偏移成像方法中构建频率域波场传播算子，根据反褶积成像条件构建互为共轭的正演算子和偏移算子，通过对观测地震记录进行改造，从而实现等价于应用反褶积成像条件进行成像结果更新，使得到的地下反射率模型更加精确，可以有效消除多次波成像结果中的串扰噪声，提高成像结果的分辨率和信噪比。同时保证了最小二乘偏移成像方法的稳定性和收敛性，加快算法收敛速度，使最小二乘偏移方法能够得到更广泛的应用。

Description

基于反褶积成像条件的最小二乘偏移成像方法、设备及存储介质

技术领域

本发明涉及地球物理勘探技术领域，具体涉及一种基于反褶积成像条件的最小二乘偏移成像方法、设备及存储介质。

背景技术

随着地震波油气勘探目标复杂程度的不断提高和地震勘探技术的发展，工业界对地下照明、地下反射率估计等问题的精度也提出了更高要求。地震偏移成像方法也逐渐由叠后偏移走向叠前偏移，从基于射线理论走向基于波动方程理论，为地下高精度结构成像和反射率信息获取提供了更加坚实的理论基础和实施方案。基于波动方程的叠前深度偏移成像方法正得到越来越多的关注。

在众多偏移成像方法和技术中，由于基于互相关成像条件的最小二乘偏移能够在反映地下结构信息的同时，提供相对准确的地下反射率信息，因此得到了目前工业界的广泛关注及初步应用，其原理是将偏移成像纳入反演框架进行计算，主要包括一个用于生成模拟记录的正演算子、一个偏移算子和一个线性反演求解器。

对于传统最小二乘偏移而言，其正演算子通常采用波恩近似(散射近似)或WKBJ近似(反射近似)将波动方程线性化来得到的，其偏移算子可通过伴随状态法求取，其线性反演求解器即是通过优化算法，如最速下降法、共轭梯度法、LSQR算法等，迭代求解地下反射率模型。

然而，在目前的最小二乘偏移成像方法中，还存在以下问题：

(1)传统最小二乘偏移实现过程中采用的是互相关成像条件，其会导致算法对于反射率估计准确性差，同时无法有效消除多次波成像中的串扰噪声。

(2)性能更高的成像条件无法直接在最小二乘偏移成像方法中应用，因为这会破坏正演算子和偏移算子之间的共轭性，导致算法无法收敛。

发明内容

有鉴于此，本发明的第一个目的在于提出一种基于反褶积成像条件的最小二乘偏移成像方法，通过反褶积成像条件构建共轭算子，并通过对观测记录进行改造，实现等价于应用反褶积成像条件进行成像结果更新，实现加快算法收敛速度、更准确估计地下反射率模型、能够有效消除多次波成像中串扰噪声。

基于同样的发明构思，本发明的第二个目的在于提出一种计算机设备。

基于同样的发明构思，本发明的第三个目的在于提出一种存储介质。

本发明的第一个目的可以通过采取如下技术方案实现：

一种基于反褶积成像条件的最小二乘偏移成像方法，包括以下步骤：

构建波场传播算子；

根据反褶积成像条件，构建正演算子和偏移算子；

观测地震记录结果，获取第一地震记录；

根据第一地震记录，获得初始成像结果；

对第一地震记录进行改造，得到第二地震记录，使根据第二地震记录每一次迭代生成的成像结果更新量等价于应用了基于反褶积成像条件的偏移成像；

根据成像结果，使用构建的正演算子，计算获得模拟地震记录；

计算模拟地震记录与第二地震记录之间的记录残差；

根据残差，应用偏移算子计算目标函数对应的成像结果更新量，利用线性反演求解器获取新的成像结果；

判断是否满足预设终止条件，若满足则输出成像结果；否则继续计算残差和成像结果，进行迭代。

进一步的，所述构建波场传播算子，是在频率域中构建波场传播算子。

进一步的，构建的波场传播算子，用公式表示为：

其中，p表示传播的波场，ω是角频率，c是参考速度，v是地震波在地下传播的速度，a和b为差分系数，x和y为坐标，

表示拉普拉斯算子。

进一步的，所述根据第一地震记录，获得初始成像结果，具体包括：

利用波场传播算子传播震源函数或第一地震记录，应用反褶积成像条件进行单程波偏移成像，获得初始成像结果。

进一步的，所述根据反褶积成像条件，构建正演算子和偏移算子，包括：

对最小二乘偏移最优化问题进行加权修改，并在最小二乘偏移方法中引入反褶积成像条件；

根据最小二乘偏移问题中的正演算子和基于反褶积成像条件的偏移算子，构建两个相互共轭的新的正演和偏移算子。

进一步的，对最小二乘偏移最优化问题进行加权修改，并在最小二乘偏移方法中引入反褶积成像条件，用公式表示为：

其中，m为地下反射率模型，d_obs为第一地震记录，L表示最小二乘偏移的正演算子，L_d*表示基于反褶积成像条件的偏移算子，w为加权系数。

进一步的，所述根据反褶积成像条件，构建正演算子和偏移算子，包括：

构造一对相互共轭的新的正演算子和偏移算子，使其等于最小二乘偏移问题中的正演算子和基于反褶积成像条件的偏移算子的乘积：

其中，L为最小二乘偏移的正演算子，L_d ^*为基于反褶积成像条件的偏移算子，Λ和Λ^*分别为所构建的正演算子和偏移算子；

推导得到对应的正演和偏移过程的表达式：

其中，m为地下反射率模型，ω为角频率，d′_obs为第二地震记录；G_s为震源格林函数，S为正传地震波场，G_r为接收格林函数，上标*表示伴随算子；x_r表示地震波检波器的空间位置。

进一步的，对第一地震记录进行改造，得到第二地震记录，具体包括：

在第一地震记录上应用一个变换矩阵T，得到第二地震记录，其表达式为：

d′_obs＝Td_obs

其中，d_obs为第一地震记录，d′_obs为第二地震记录，变换矩阵T的具体形式为：

其中，G_r为接收格林函数，G_s为震源格林函数，S为正传地震波场，ω为角频率，上标*表示伴随算子。

本发明的第二个目的可以通过采取如下技术方案实现：

一种计算机设备，包括处理器以及用于存储处理器可执行程序的存储器，其特征在于，所述处理器执行存储器存储的程序时，实现上述的基于反褶积成像条件的最小二乘偏移成像方法。

本发明的第三个目的可以通过采取如下技术方案实现：

一种存储介质，存储有程序，其特征在于，所述程序被处理器执行时，实现上述的基于反褶积成像条件的最小二乘偏移成像方法。

本发明相对于现有技术具有如下的有益效果：

(1)本发明提供的基于反褶积成像条件的波动方程最小二乘偏移方法中构建的新的相互共轭的正演和偏移算子保证了算法的稳定性和收敛性，同时为快速线性反演求解器的使用提供了便利条件；

(2)本发明使用改造后的地震观测记录进行成像和残差计算，保证了在新的正演及偏移算子的作用下，实现等价于应用反褶积成像条件进行成像结果更新，有助于加快算法收敛速度。

(3)本发明提供的方法能够有效提高算法收敛速度，能够为地震勘探提供更加精确的地下反射率模型，同时可以有效消除多次波成像结果中的串扰噪声，提高成像结果的分辨率和信噪比。

附图说明

图1是本发明实施例1提供的基于反褶积成像条件的波动方程最小二乘偏移方法的流程示意图；

图2是本发明实施例1提供的传统基于互相关成像条件和反褶积成像条件的最小二乘多次波偏移方法中正演与偏移算子的点积测试结果示意图；

图3是本发明实施例1所采用的地下真实反射率模型示意图；

图4是本发明实施例1中传统基于互相关成像条件的常规多次波偏移成像结果(初始成像结果)示意图；

图5是本发明实施例1中传统基于互相关成像条件的最小二乘多次波偏移成像结果(迭代次数为5次)示意图；

图6是本发明实施例1中基于反褶积成像条件的常规多次波偏移成像结果(初始成像结果)示意图；

图7是本发明实施例1中基于反褶积成像条件的最小二乘多次波偏移成像结果(迭代次数为5次)示意图。

图8为本发明实施例2的计算机设备示意图。

具体实施方式

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明的一部分实施例，而不是全部的实施例，基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

实施例1：

本实施例提供了一种基于反褶积成像条件的最小二乘偏移成像方法，包括以下步骤：

S100、构建波场传播算子；

本实施例中，使用傅里叶有限差分法在频率域中构建单程波场传播算子，如式(1)所示：

其中，p表示传播的波场，ω是角频率，c是参考速度，v是地震波在地下传播的速度，a和b为差分系数，x和y为坐标。

由于该算子是在频率域中实现计算的，有利于反褶积成像条件的实施。其对应的波场正传和反传过程可用格林函数来表示，如式(2)所示：

其中，p_s和p_r分别表示传播的波场，G_s和G_r分别表示震源和接收格林函数，S是震源函数，d_obs是第一地震记录。

S200、根据反褶积成像条件，构建正演算子和偏移算子，包括以下步骤：

S201、对表示传统基于互相关成像条件的最小二乘偏移方法的最优化问题，即式(3)，进行加权改造(也叫加权修改)，使得反褶积成像条件可以应用于最小二乘偏移方法中，如式(4)所示：

其中，m为地下反射率模型，d_obs是第一地震记录，L表示传统最小二乘偏移方法中的正演算子，

表示基于反褶积成像条件的偏移算子，w为加权系数；

S202、由于L和

之间非共轭，直接求解所述改造后最优化问题会导致算法收敛性差甚至不收敛，为了保证算法的收敛性，对改造后最优化问题应用一个变换矩阵T，对最优化问题进行进一步改造，如式(5)所示：

S203、定义两个互为共轭的算子Λ和Λ^*，并使其满足式(6)：

LL_d ^*＝ΛΛ^*         (6)

即：

式(7)中，d′_obs＝Td_obs，为第二地震记录；推导整理后，得到互为共轭的新正演算子Λ和偏移算子Λ^*对应的最小二乘偏移中的正演过程如式(8)所示，偏移过程如式(9)所示：

其中，m为地下反射率模型，ω为角频率，G_s为震源格林函数，S为正传地震波场，G_r为接收格林函数，上标*表示伴随算子；x_r表示地震波检波器的空间位置。

S300、观测地震记录结果，获取第一地震记录；

本实施例中，观测地震记录为人工爆炸源激发(陆上工区为炸药、海上工区为气枪)，高密度地震检波器阵列接收的主动源地震数据；目标工区地下速度模型由观测地震数据应用层析成像或全波形反演获得；另外，为提高地震成像质量，可对观测地震数据进行去噪(去除随机噪音、面波噪音、无线电干扰、邻炮干扰、线性噪音等)、切除(切除地震记录中的折射波成分)、反褶积等操作。

S400、根据第一地震记录，获得初始成像结果；

本实施例中，将目标工区观测地震记录及地下速度模型作为输入，利用S100所构建的波场传播算子传播震源函数或观测地震记录，应用反褶积成像条件，获得初始成像结果，具体计算如式(10)所示：

其中，m为地下反射率模型，ω为角频率，G_s为震源格林函数，S为正传地震波场，G_r为接收格林函数，d_obs是第一地震记录，上标*表示伴随算子。

S500、对第一地震记录进行改造，得到第二地震记录，使根据第二地震记录每一次迭代生成的成像结果更新量等价于应用了基于反褶积成像条件的偏移成像；

本实施例中，对第一地震记录进行改造的方式是在第一地震记录上应用一个变换矩阵T，得到第二地震记录，即d′_obs＝Td_obs，其中变换矩阵T如式(11)所示：

可以发现，改造后的地震数据符合S203构建得出的新正演算子，将其应用于S203所述构建得出的新偏移算子，其形式如式(12)所示：

这表示可以使用互为共轭的算子Λ和Λ^*进行最小二乘偏移方法计算，在保证算法收敛性的同时，可以获得基于反褶积成像条件的成像结果。

S600、根据成像结果，使用构建的正演算子，计算获得模拟地震记录；

本实施例中，根据步骤S400所得到的初始成像结果，或步骤S800所得到的新的成像结果，应用步骤S203构建的正演算子反偏移计算相对应的模拟地震记录。

S700、计算模拟地震记录与第二地震记录之间的记录残差；

本实施例中，将步骤S600中获得的模拟地震记录与S500中获得的第二地震记录作差，得到记录残差。

S800、根据残差，应用偏移算子计算目标函数对应的成像结果更新量，利用线性反演求解器获取新的成像结果；

S900、判断是否满足预设终止条件，若满足则输出成像结果；否则继续计算残差和成像结果，进行迭代。

本实施例中，进行迭代的方法是：将步骤S800获取的新的成像结果作为步骤S600的成像结果输入，并再次执行步骤S600-S900。

本实施例中，预设终止条件为：满足预设条件中的任意一条，包括：

第一预设条件，与残差相关的目标函数小于设定阈值；

第二预设条件，算法迭代次数大于设定阈值。

图2为传统基于互相关成像条件和反褶积成像条件的最小二乘偏移方法中正演与偏移算子的点积测试结果示意图。可见，本发明中构建得出的正演与偏移算子几乎严格共轭，这保证了算法的收敛性。

图3为本实施例中采用的地下真实反射率模型。图4和图5分别为本发明实施例中使用传统基于互相关成像条件获得的常规偏移成像结果和最小二乘偏移成像结果(迭代次数为5次)，这些结果对应输入速度模型大小为3201*1201个网格点，网格间距为25英尺，输入地震数据为包含表面多次波的500个长度为20秒的共炮点道集；可以看出图4中的成像结果存在强烈的串扰噪声，在矩形方框内和箭头标记部分尤其明显，图5经过5次最小二乘偏移的成像结果中除箭头所表示的最强的串扰外的大部分串扰几乎都得到了有效抑制。

图6和图7分别为本发明实施例中使用基于反褶积成像条件获得的常规偏移成像结果和最小二乘偏移成像结果(迭代次数为5次)；可以看出图6中成像结果的串扰噪声明显少于图4中成像结果中的，尤其是对应图4矩形方框内的部分，展示了反褶积成像条件的优越性，图7经过5次基于反褶积成像条件的最小二乘偏移成像结果几乎完全解决了串扰噪声的问题，通过图7和图5的对比可以发现，使用新构建的正演和偏移算子，可以更加稳健地解决最小二乘偏移这个线性反演问题，加快算法收敛效率，同时有效解决传统多次波偏移方法中串扰噪声的问题，提高成像结果的分辨率和信噪比。

因此，本实施例中构建的新的相互共轭的正演和偏移算子保证了算法的稳定性和收敛性，同时为快速线性反演求解器的使用提供了便利条件；本实施例使用改造后的地震观测记录进行成像和残差计算，保证了在新的正演及偏移算子的作用下，实现等价于应用反褶积成像条件进行成像结果更新，有助于加快算法收敛速度；本实施例提供的方法能够有效提高算法收敛速度，能够为地震勘探提供更加精确的地下反射率模型，同时可以有效消除多次波成像结果中的串扰噪声，提高成像结果的分辨率和信噪比。

实施例2：

如图8所示，本实施例提供了一种计算机设备，该计算机设备包括通过系统总线801连接的处理器802、存储器、输入装置803、显示装置804和网络接口805，该处理器用于提供计算和控制能力，该存储器包括非易失性存储介质806和内存储器807，该非易失性存储介质806存储有操作系统、计算机程序和数据库，该内存储器807为非易失性存储介质中的操作系统和计算机程序的运行提供环境，处理器802执行存储器存储的计算机程序时，实现上述实施例1的基于反褶积成像条件的最小二乘偏移成像方法，具体如下：

构建波场传播算子；

根据反褶积成像条件，构建正演算子和偏移算子；

观测地震记录结果，获取第一地震记录；

根据第一地震记录，获得初始成像结果；

对第一地震记录进行改造，得到第二地震记录，使根据第二地震记录每一次迭代生成的成像结果更新量等价于应用了基于反褶积成像条件的偏移成像；

根据成像结果，使用构建的正演算子，计算获得模拟地震记录；

计算模拟地震记录与第二地震记录之间的记录残差；

根据残差，应用偏移算子计算目标函数对应的成像结果更新量，利用线性反演求解器获取新的成像结果；

判断是否满足预设终止条件，若满足则输出成像结果；否则继续计算残差和成像结果，进行迭代。

实施例3：

本实施例提供了一种存储介质，该存储介质为计算机可读存储介质，其存储有计算机程序，所述计算机程序被处理器执行时，实现上述实施例1的基于反褶积成像条件的最小二乘偏移成像方法，具体如下：

构建波场传播算子；

根据反褶积成像条件，构建正演算子和偏移算子；

观测地震记录结果，获取第一地震记录；

根据第一地震记录，获得初始成像结果；

对第一地震记录进行改造，得到第二地震记录，使根据第二地震记录每一次迭代生成的成像结果更新量等价于应用了基于反褶积成像条件的偏移成像；

根据成像结果，使用构建的正演算子，计算获得模拟地震记录；

计算模拟地震记录与第二地震记录之间的记录残差；

根据残差，应用偏移算子计算目标函数对应的成像结果更新量，利用线性反演求解器获取新的成像结果；

判断是否满足预设终止条件，若满足则输出成像结果；否则继续计算残差和成像结果，进行迭代。

需要说明的是，本实施例的计算机可读存储介质可以是计算机可读信号介质或者计算机可读存储介质或者是上述两者的任意组合。计算机可读存储介质例如可以是但不限于电、磁、光、电磁、红外线、或半导体的系统、装置或器件，或者任意以上的组合。计算机可读存储介质的更具体的例子可以包括但不限于：具有一个或多个导线的电连接、便携式计算机磁盘、硬盘、随机访问存储器(RAM)、只读存储器(ROM)、可擦式可编程只读存储器(EPROM或闪存)、光纤、便携式紧凑磁盘只读存储器(CD-ROM)、光存储器件、磁存储器件、或者上述的任意合适的组合。

在本实施例中，计算机可读存储介质可以是任何包含或存储程序的有形介质，该程序可以被指令执行系统、装置或者器件使用或者与其结合使用。而在本实施例中，计算机可读信号介质可以包括在基带中或者作为载波一部分传播的数据信号，其中承载了计算机可读的程序。这种传播的数据信号可以采用多种形式，包括但不限于电磁信号、光信号或上述的任意合适的组合。计算机可读信号介质还可以是计算机可读存储介质以外的任何计算机可读存储介质，该计算机可读信号介质可以发送、传播或者传输用于由指令执行系统、装置或者器件使用或者与其结合使用的程序。计算机可读存储介质上包含的计算机程序可以用任何适当的介质传输，包括但不限于：电线、光缆、RF(射频)等等，或者上述的任意合适的组合。

上述计算机可读存储介质可以以一种或多种程序设计语言或其组合来编写用于执行本实施例的计算机程序，上述程序设计语言包括面向对象的程序设计语言—诸如Java、Python、C++，还包括常规的过程式程序设计语言—诸如C语言或类似的程序设计语言。程序可以完全地在用户计算机上执行、部分地在用户计算机上执行、作为一个独立的软件包执行、部分在用户计算机上部分在远程计算机上执行、或者完全在远程计算机或服务器上执行。在涉及远程计算机的情形中，远程计算机可以通过任意种类的网络，包括局域网(LAN)或广域网(WAN)连接到用户计算机，或者，可以连接到外部计算机(例如利用因特网服务提供商来通过因特网连接)。

显然，上述所述的实施例只是本发明的一部分实施例，而不是全部实施例，本发明不限于上述实施例的细节，任何所属技术领域的普通技术人员对其所做的适当变化或修饰，皆视为不脱离本发明的专利范畴。

Claims (10)

1.一种基于反褶积成像条件的最小二乘偏移成像方法，其特征在于，包括以下步骤：

构建波场传播算子；

对表示基于互相关成像条件的最小二乘偏移方法的最优化问题，进行改造，得到适用于最小二乘偏移方法的反褶积成像条件；根据反褶积成像条件，构建正演算子和偏移算子；

观测地震记录结果，获取第一地震记录；

根据第一地震记录，获得初始成像结果；

对第一地震记录进行改造，得到第二地震记录，使根据第二地震记录每一次迭代生成的成像结果更新量等价于应用了基于反褶积成像条件的偏移成像；

根据初始成像结果或新的成像结果，使用构建的正演算子，计算获得模拟地震记录；

计算模拟地震记录与第二地震记录之间的记录残差；

根据残差，应用偏移算子计算目标函数对应的成像结果更新量，利用线性反演求解器获取新的成像结果；

判断是否满足预设终止条件，若满足则输出成像结果；否则继续计算残差和成像结果，进行迭代。

2.根据权利要求1所述的基于反褶积成像条件的最小二乘偏移成像方法，其特征在于，所述构建波场传播算子，是在频率域中构建波场传播算子。

3.根据权利要求2所述的基于反褶积成像条件的最小二乘偏移成像方法，其特征在于，构建的波场传播算子，用公式表示为：

其中，p表示传播的波场，ω是角频率，c是参考速度，v是地震波在地下传播的速度，a和b为差分系数，x和y为坐标，

表示拉普拉斯算子。

4.根据权利要求1所述的基于反褶积成像条件的最小二乘偏移成像方法，其特征在于，所述根据第一地震记录，获得初始成像结果，具体包括：

利用波场传播算子传播震源函数或第一地震记录，应用反褶积成像条件进行单程波偏移成像，获得初始成像结果。

5.根据权利要求1所述的基于反褶积成像条件的最小二乘偏移成像方法，其特征在于，所述根据反褶积成像条件，构建正演算子和偏移算子，包括：

对表示基于互相关成像条件的最小二乘偏移方法的最优化问题进行加权修改，使反褶积成像条件适用于最小二乘偏移方法；

对加权修改后的最优化问题，应用变换矩阵以作进一步改造；

根据最小二乘偏移问题中的正演算子和基于反褶积成像条件的偏移算子，构建两个相互共轭的新的正演算子和偏移算子。

6.根据权利要求5所述的基于反褶积成像条件的最小二乘偏移成像方法，其特征在于，对最小二乘偏移最优化问题进行加权修改，使反褶积成像条件适用于最小二乘偏移方法，用公式表示为：

其中，m为地下反射率模型，d_obs为第一地震记录，L表示最小二乘偏移的正演算子，L_d*表示基于反褶积成像条件的偏移算子，w为加权系数。

7.根据权利要求6所述的基于反褶积成像条件的最小二乘偏移成像方法，其特征在于，所述根据反褶积成像条件，构建正演算子和偏移算子，包括：

构造一对相互共轭的新的正演算子和偏移算子，使其等于最小二乘偏移问题中的正演算子和基于反褶积成像条件的偏移算子的乘积：

其中，L为最小二乘偏移的正演算子，L_d ^*为基于反褶积成像条件的偏移算子，Λ和Λ^*分别为所构建的正演算子和偏移算子；

推导得到对应的正演和偏移过程的表达式：

其中，m为地下反射率模型，ω为角频率，d′_obs为第二地震记录；G_s为震源格林函数，S为正传地震波场，G_r为接收格林函数，上标*表示伴随算子；x_r表示地震波检波器的空间位置。

8.根据权利要求1-7中任一项所述的基于反褶积成像条件的最小二乘偏移成像方法，其特征在于，对第一地震记录进行改造，得到第二地震记录，具体包括：

在第一地震记录上应用一个变换矩阵T，得到第二地震记录，其表达式为：

d′_obs＝Td_obs

其中，d_obs为第一地震记录，d′_obs为第二地震记录，变换矩阵T的具体形式为：

其中，G_r为接收格林函数，G_s为震源格林函数，S为正传地震波场，ω为角频率，上标*表示伴随算子。

9.一种计算机设备，包括处理器以及用于存储处理器可执行程序的存储器，其特征在于，所述处理器执行存储器存储的程序时，实现权利要求1-8任一项所述的基于反褶积成像条件的最小二乘偏移成像方法。

10.一种存储介质，存储有程序，其特征在于，所述程序被处理器执行时，实现权利要求1-8任一项所述的基于反褶积成像条件的最小二乘偏移成像方法。

Images