CN103091710B - 一种逆时偏移成像方法及装置 - Google Patents

Abstract

本发明是关于一种逆时偏移成像方法及装置，所述的方法包括：获取震源波场信息和检波点波场信息；根据获取的震源波场信息和检波点波场信息，将震源波场分离为震源上行波场和震源下行波场，将检波点波场分离为检波点上行波场和检波点下行波场；对分离后的震源波场和检波点波场进行逆时偏移成像处理，得到逆时偏移成像剖面。本发明实施例提供的逆时偏移成像方法及装置，消除了逆时偏移过程中所产生的低频噪音，具有很高的保幅性，并保持了内部波场计算所具有的良好的并行特征，非常在适合在图形处理器GPU/中央处理器CPU异构平台上进行，极大的提高逆时偏移的计算效率。

Description

一种逆时偏移成像方法及装置

技术领域

本发明是关于地球物理勘探地震数据处理技术领域，尤其是关于石油天然气勘探地震数据处理技术领域，具体来说是关于一种逆时偏移成像方法及装置。

背景技术

随着油气勘探技术的不断发展，油气勘探的目标逐渐转向复杂地表条件及复杂地下地质构造的“双复杂”区域，而传统的成像方法已经无法满足“双复杂”区域的成像。目前国内大规模工业化应用的深度偏移成像方法主要是基于射线理论的Kirchhoff积分法和射线束类偏移方法，但这类方法是波动方程的近似解，它们对于高陡倾角并不能很好成像，而且复杂构造通常伴随着强烈的横向速度变化，因此结果会有多路径，从而会限制基于射线理论偏移方法的计算效率。为解决上述技术问题，现有技术提供一种基于单程波的波动方程偏移算法，利用波动方程的傍轴近似理论实现波场的外推，在特定的角度下可以很好的成像，但是当地层接近甚至超过90度的时候，这个方法同样会失效。

为解决上述技术问题，现有技术提出一种逆时偏移方法，逆时偏移方法最早是由Whitemore N D.等人在1982年SEG年会上就提出的，但受制于其巨大的计算量和高昂的价格，这个方法一直没有大规模应用于实际生产，而近几年，随着超大计算机群及GPU(Graphic Processing Unit，图形处理器)技术的发展，逆时偏移技术逐渐被大规模应用。叠前逆时偏移方法基本可以分为三步：1、在时间域进行炮点波场的正向外推；2、在时间域进行检波点波场反向外推；3、利用有效的成像条件在反射层出现的位置构建成像。目前，大规模工业化的逆时偏移波场外推，仍然是通过下述常密度声波方程(1)实现的：

$\frac{1}{V^2(x,y,z)}\frac{{\∂}^{2}P}{\∂t^2}=\frac{{\∂}^{2}P}{\∂x^2}+\frac{{\∂}^{2}P}{\∂y^2}+\frac{{\∂}^{2}P}{\∂z^2}---\left(1\right)$

其中t代表时间坐标，x，y，z代表空间坐标，P表示波场在空间位置(x，y，z)处的压力场，V表示纵横向可变的介质速度函数。在逆时偏移成像过程中，炮点波场和检波点波场都是通过方程(1)实现外推的，而实现成像时通常利用炮点和检波点外推波场的零延迟互相关来构建成像条件，如下方程(2)：

$I(x,y,z)=\underset{0}{\overset{T_{\max }}{\∫}}s(t,x,y,z)r(t,x,y,z)\mathrm{dt}---\left(2\right)$

方程(2)中，s是正向外推的震源波场，r是反向外推的检波点波场，I为成像结果，T_max为波场外推的最大时间；利用此成像条件进行逆时偏移需要实现三个环节，首先，将震源波场进行正向延拓，保存每一时间步的波场信息；然后，将检波点波场进行反向延拓，保存每一时间步的波场信息；最后，分别读取同一时刻的震源波场和检波点波场进行成像运算，累加入成像空间，便可得到最终成像结果。

然而，方程(2)所定义的成像条件通常会产生大量强振幅的低频噪音，有时会在一些反射界面处出现偏移假象，特别是在强反射界面上方，这些低频噪音甚至会污染有效信号而使浅层的地质构造被彻底淹没；在实际应用过程中，为有效消除这些低频噪音，一般会采用以下三种主要去噪思路，即：在波场传播过程中去噪、在应用成像条件时去噪或成像后滤波法去噪。然而，三种去噪方案中的任何一种都会或多或少伤害到有效信号，从而破坏资料的保幅性；通过分析研究，发现产生这些低频噪音的主要原因是震源波场和检波点波场的互相关不仅会在反射界面处产生振幅值，同时在波传播的整个路径的非反射点处，也会产生振幅，如果将这些振幅沿着时间求和，就会产生强振幅的低频噪音。为了消除这种低频噪音，必须要采用各种去噪方法，但不论何种去噪方法，都会破坏资料的有效信号，而破坏逆时偏移的保幅性。

发明内容

为克服现有技术中存在的问题，本发明提供一种逆时偏移成像方法及装置，在实际应用过程中，并不会产生低频噪音，具有很高的保幅性。

本发明提供一种逆时偏移成像方法，所述的方法包括：获取震源波场信息和检波点波场信息；根据获取的震源波场信息和检波点波场信息，将震源波场分离为震源上行波场和震源下行波场，将检波点波场分离为检波点上行波场和检波点下行波场；对分离后的震源波场和检波点波场进行逆时偏移成像处理，得到逆时偏移成像剖面；

所述的获取震源波场信息和检波点波场信息包括：

获取单炮数据，并对单炮数据进行叠前预处理；

对叠前预处理输出的单炮数据的正向波场外推及逆时外推，获得各个时间步的震源波场信息和检波点波场信息。

本发明还提供一种逆时偏移成像装置，所述的装置包括：波场信息获取单元，用于获取震源波场信息和检波点波场信息；波场分离单元，用于根据获取的震源波场信息和检波点波场信息，将震源波场分离为震源上行波场和震源下行波场，将检波点波场分离为检波点上行波场和检波点下行波场；成像单元，用于对分离后的震源波场和检波点波场进行逆时偏移成像处理，得到逆时偏移成像剖面；

所述的波场信息获取单元包括：

预处理模块，用于获取单炮数据，并对单炮数据进行叠前预处理；

外推处理模块，用于对叠前预处理输出的单炮数据的正向波场外推及逆时外推，获得各个时间步的震源波场信息和检波点波场信息。

本发明实施例提供的逆时偏移成像方法及装置，消除了逆时偏移过程中所产生的低频噪音，具有很高的保幅性，并保持了内部波场计算所具有的良好的并行特征，非常适合在GPU/CPU(Central Processing Unit，中央处理器)异构平台上进行，极大的提高逆时偏移的计算效率。

附图说明

此处所说明的附图用来提供对本发明的进一步理解，构成本申请的一部分，并不构成对本发明的限定。在附图中：

图1是本发明实施例提供的一种逆时偏移成像方法的流程图。

图2是本发明实施例提供的同一震源波场的不同传播路径的水平层状模型示意图。

图3是本发明实施例提供的一种逆时偏移成像方法的流程图。

图4是本发明实施例提供的一种逆时偏移成像装置模块图。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚明白，下面结合实施方式和附图，对本发明做进一步详细说明。在此，本发明的示意性实施方式及其说明用于解释本发明，但并不作为对本发明的限定。

本发明实施例提供一种逆时偏移成像方法及装置，以下结合附图对本发明进行详细说明。

实施例一

图1是本发明实施例提供的一种逆时偏移成像方法的流程图，如图1所示，逆时偏移成像方法包括如下步骤：

S101，获取震源波场信息和检波点波场信息。

在本发明实施例中，步骤S101可以包括如下子步骤：

获取单炮数据，并对单炮数据进行叠前预处理；

对叠前预处理输出的单炮数据的正向波场外推及逆时外推，获得各个时间步的震源波场信息和检波点波场信息。

具体的，可以加载地震采集的单炮数据；对加载单炮数据进行静校正、去噪等叠前预处理工作，以便获得高信噪比、高保幅的单炮数据；提供相对准确的速度场在GPU/CPU异构平台上实现地震单炮数据的正向波场外推及逆时外推，获得各个时间步的震源波场信息和检波点波场信息并存储。

S102，根据获取的震源波场信息和检波点波场信息，将震源波场分离为震源上行波场和震源下行波场，将检波点波场分离为检波点上行波场和检波点下行波场。

图2是本发明实施例提供的同一震源波场的不同传播路径的水平层状模型示意图，如图2所示，n代表法向方向；I代表入射波；T代表透射波；R代表反射波，沿着波的传播路径，炮点波场和检波点波场通常有着相同的方向，同时，入射波场和反射波场在反射界面的法向方向的投影有着相反的方向。按照震源波场和检波点波场的传播方向，将两个波场分离成其各自的上行波场和下行波场，进行波场分离以后的震源和检波点波场的相关性与其传播方向在法向的投影有关，当他们有着相反的投影方向时，两个波场分离以后，至少有一个会是零，仅会在反射点处产生成像，而在非反射点处就不会产生成像。

在本发明实施例中，提出了一种有效的将震源波场和检波点波场分离成其上、下行分量的方法，主要是将上、下行波在F-K域对其进行分离。

在地表某一点(x，y，z＝z0)处，震源波场s(t，x，y，z)和检波点波场r(t，x，y，z)可以通过傅里叶变换实现震源和检波点波场上、下行波的分离。在傅里叶域，上行波场和下行波场能够被成功分离，并且可以分别表示为如下方程：

$S_u(\mathrm{\&omega;},k_z)=\left\{\begin{array}{cc}s(\mathrm{\&omega;},k_z)& \mathrm{\&omega;}{k}_z\&GreaterEqual;0\\ 0& \mathrm{\&omega;}{k}_z<0\end{array}\right.---\left(3\right)$

$S_d(\mathrm{\&omega;},k_z)=\left\{\begin{array}{cc}s(\mathrm{\&omega;},k_z)& \mathrm{\&omega;}{k}_z<0\\ 0& \mathrm{\&omega;}{k}_z\&GreaterEqual;0\end{array}\right.---\left(4\right)$

$R_u(\mathrm{\&omega;},k_z)=\left\{\begin{array}{cc}R(\mathrm{\&omega;},k_z)& \mathrm{\&omega;}{k}_z\&GreaterEqual;0\\ 0& \mathrm{\&omega;}{k}_z<0\end{array}\right.---\left(5\right)$

$R_d(\mathrm{\&omega;},k_z)=\left\{\begin{array}{cc}R(\mathrm{\&omega;},k_z)& \mathrm{\&omega;}{k}_z<0\\ 0& \mathrm{\&omega;}{k}_z\&GreaterEqual;0\end{array}\right.---\left(6\right)$

式中，S_u、S_d分别代表震源上、下行波场的傅里叶变换，R_u、R_d分别代表检波点上、下行波场的傅里叶变换。

S103，对分离后的震源波场和检波点波场进行逆时偏移成像处理，得到逆时偏移成像剖面。

在本发明实施例中，将分离后的震源波场和检波点波场转换至时间域，便可得到基于波场分离的逆时偏移的时间域成像条件，转换至时间域的表达式分别为：

$s_u(t,x,y,z)=\underset{-\&infin;}{\overset{+\&infin;}{\&Integral;}}\underset{-\&infin;}{\overset{+\&infin;}{\&Integral;}}{S}_u(\mathrm{\&omega;},k_z)e^{i(\mathrm{tw}+z{k}_z)}\mathrm{d\&omega;d}{k}_z---\left(7\right)$

$s_d(t,x,y,z)=\underset{-\&infin;}{\overset{+\&infin;}{\&Integral;}}\underset{-\&infin;}{\overset{+\&infin;}{\&Integral;}}{S}_d(\mathrm{\&omega;},k_z)e^{i(\mathrm{tw}+z{k}_z)}\mathrm{d\&omega;d}{k}_z---\left(8\right)$

$r_u(t,x,y,z)=\underset{-\&infin;}{\overset{+\&infin;}{\&Integral;}}\underset{-\&infin;}{\overset{+\&infin;}{\&Integral;}}{R}_u(\mathrm{\&omega;},k_z)e^{i(\mathrm{tw}+z{k}_z)}\mathrm{d\&omega;d}{k}_z---\left(9\right)$

$r_d(t,x,y,z)=\underset{-\&infin;}{\overset{+\&infin;}{\&Integral;}}\underset{-\&infin;}{\overset{+\&infin;}{\&Integral;}}{R}_d(\mathrm{\&omega;},k_z)e^{i(\mathrm{tw}+z{k}_z)}\mathrm{d\&omega;d}{k}_z---\left(10\right)$

在逆时偏移中，震源和检波点波场外推是利用全程声波波动方程展开的，它包含震源波场和检波点波场两个组成部分，可以分别表示为如下表达式：

s(t,x,y,z)＝s_d(t,x,y,z)+s_u(t,x,y,z)         (11)

r(t,x,y,z)＝r_d(t,x,y,z)+r_u(t,x,y,z)            (12)

式中，s代表震源波场，s^d代表震源的下行波场，s_u代表震源的上行波场，r代表检波点波场，r^d代表检波点下行波场，r_u代表检波点上行波场。分别将其带入方程(2)便可得到四项互相关结果:

$\begin{array}{c}I(x,y,z)=\underset{0}{\overset{T\max }{\&Integral;}}{s}_d(t,x,y,z)r_u(t,x,y,z)\mathrm{dt}+\underset{0}{\overset{T\max }{\&Integral;}}{s}_u(t,x,y,z)r_d(t,x,y,z)\mathrm{dt}\\ +\underset{0}{\overset{T\max }{\&Integral;}}{s}_d(t,x,y,z)r_d(t,x,y,z)\mathrm{dt}+\underset{0}{\overset{T\max }{\&Integral;}}{s}_u(t,x,y,z)r_u(t,x,y,z)\mathrm{dt}\\ =I_{z1}(x,y,z)+I_{z2}(x,y,z)+I_{z3}(x,y,z)+I_{z4}(x,y,z)\end{array}---\left(13\right)$

方程(13)中，I_z1表示震源下行波场和检波点上行波场的互相关(表示为I_z1)，其实质是单程波的波动方程偏移；I_z2表示震源上行波场和检波点下行波场的互相关；I_z3表示震源下行波场和检波点下行波场的互相关；I_z4表示震源上行波场和检波点上行波场的互相关。依据射线理论，有些反射点处，入射波为下行波，反射波为上行波；而有些反射点处，入射波为上行波，反射波为下行波，因此将两个反射都参与进来，才能够利用射线追踪进行成像，也就是说，只需要I_z1和I_z2两项，便可实现逆时偏移成像。剩余的两个部分产生的结果便是低频噪音，也就是说利用传播方向相同的两个波场，不管是下行波I_z3还是上行波I_z4，将它们沿着整个射线路径进行互相关，其互相关结果便是成像空间的低频噪音。此时，可以得到一个能有效消除低频噪音的互相关成像条件:

I(x,y,z)＝I_z1(x,y,z)+I_z2(x,y,z)       (14)

方程(14)的成像条件，在时间域下最终表达为：

$I(x,y,z)=\underset{0}{\overset{T_{\max }}{\&Integral;}}[s_d(t,x,y,z)r_u(t,x,y,z)+r_d(t,x,y,z)s_u(t,x,y,z)]\mathrm{dt}---\left(15\right)$

根据方程(15)对分离后的震源和检波点波场进行逆时偏移成像，输出最终逆时偏移成像剖面。

本发明实施例提供的逆时偏移成像方法，沿着波的传播路径将反射点和非反射点区别开来，然后只在反射点处利用成像条件，消除了逆时偏移过程中所产生的低频噪音，具有很高的保幅性。

本发明实施例提供的逆时偏移成像方法，保持了内部波场计算所具有的良好的并行特征，非常在适合在GPU/CPU异构平台上进行，因此会极大的提高逆时偏移的计算效率。

实施例二

为了测试本发明的准确性和稳定性，本实施例对地球物理界内最典型的BP(BritishPetroleum，英国石油)数学模型进行了测试。图3是本实施例提供的一种逆时偏移成像方法的流程图，如图3所示，逆时偏移成像方法包括如下步骤：

S301，将BP模型正演得到的地震炮记录加载到工区。

S302，利用有效的去噪方法，分步分域去除地震记录中的各种噪音，然后重新抽取到炮域。

S303，提供准确的用于正演单炮记录的BP数学模型作为逆时偏移的速度—深度模型。

S304，在GPU/CPU异构平台上实现步骤S302输出的地震单炮数据的正向波场外推及逆时外推，并存储各个时间步的波场信息。

S305，按照本发明提供的波场分离方法，对步骤S304所存储的波场进行上、下行波的分离。

S306，利用本发明提供的逆时偏移成像条件，对分离后的震源和检波点波场进行逆时偏移成像，并输出最终逆时偏移成像剖面。

将应用本发明提到的逆时偏移成像条件得到的成像结果与传统互相关条件得到的逆时偏移成像结果进行比较，可以得知，利用本发明实施例提供的成像条件得到的成像结果明显优于传统成像条件得到的成像结果，特别是对强反射界面上方的构造识别能力更加强大。

在本发明实施例中，为了测试本发明的普适性，特别是针对野外实际采集资料的有效性，选取了某复杂地表、复杂地下地质构造区域的150平方公里的实测数据进行了逆时偏移测试。按照上述步骤，采用本发明所提供的逆时偏移成像条件，得到的成像结果，将其与采用传统逆时偏移成像条件，得到的去噪后的偏移结果进行对比。可以看到，采用本发明提到的成像条件，对于资料的浅层构造成像更加准确，并且增加了资料的保幅性。对于陡倾角构造的识别能力，也明显优于采用传统逆时偏移成像条件得到的成像结果。

本发明实施例提供的逆时偏移成像方法，与传统的逆时偏移成像条件相比，最大的优势在于它克服了传统逆时偏移成像条件会产生大量强振幅的低频噪音的缺陷，直接在利用成像条件的过程中就略去了低频噪音项，并不损失有效的地震信号，因此具有较高的资料保幅性；与此同时，本发明提供的逆时偏移成像条件，保持了内部波场计算所具有的良好的并行特征，非常在适合在GPU/CPU异构平台上进行，因此会极大的提高逆时偏移的计算效率。

实施例三

图4是本发明实施例提供的一种逆时偏移成像装置模块图，逆时偏移成像装置采用如实施例一所述的方法生成逆时偏移成像剖面，如图4所示，逆时偏移成像装置包括：

波场信息获取单元401，用于获取震源波场信息和检波点波场信息。

在本发明实施例中，波场信息获取单元401还可以包括：

预处理模块，用于获取单炮数据，并对单炮数据进行叠前预处理；

外推处理模块，用于对叠前预处理输出的单炮数据的正向波场外推及逆时外推，获得各个时间步的震源波场信息和检波点波场信息。

波场分离单元402，用于根据获取的震源波场信息和检波点波场信息，将震源波场分离为震源上行波场和震源下行波场，将检波点波场分离为检波点上行波场和检波点下行波场。

在本发明实施例中，波场分离单元402通过傅里叶变换将震源波场分离为震源上行波场和震源下行波场，将检波点波场分离为检波点上行波场和检波点下行波场。

成像单元403，用于对分离后的震源波场和检波点波场进行逆时偏移成像处理，得到逆时偏移成像剖面。

在本发明实施例中，所述成像单元403包括：

时间域转换模块，用于将分离后的震源波场和检波点波场转换至时间域，得到逆时偏移的时间域成像条件；

逆时偏移成像模块，用于根据所述的逆时偏移的时间域成像条件对分离后的震源和检波点波场进行逆时偏移成像，生成逆时偏移成像剖面。

本发明实施例提供的逆时偏移成像装置，沿着波的传播路径将反射点和非反射点区别开来，然后只在反射点处利用成像条件，消除了逆时偏移过程中所产生的低频噪音，具有很高的保幅性。

本发明实施例提供的逆时偏移成像装置，保持了内部波场计算所具有的良好的并行特征，非常在适合在GPU/CPU异构平台上进行，因此会极大的提高逆时偏移的计算效率。

以上所述的具体实施方式，对本发明的目的、技术方案和有益效果进行了进一步详细说明，所应理解的是，以上所述仅为本发明的具体实施方式而已，并不用于限定本发明的保护范围，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (6)

1.一种逆时偏移成像方法，其特征在于，所述的方法包括：

获取震源波场信息和检波点波场信息；

根据获取的震源波场信息和检波点波场信息，将震源波场分离为震源上行波场和震源下行波场，将检波点波场分离为检波点上行波场和检波点下行波场；

对分离后的震源波场和检波点波场进行逆时偏移成像处理，得到逆时偏移成像剖面；

所述的获取震源波场信息和检波点波场信息包括：

获取单炮数据，并对单炮数据进行叠前预处理；

对叠前预处理输出的单炮数据的正向波场外推及逆时外推，获得各个时间步的震源波场信息和检波点波场信息。

2.根据权利要求1所述的逆时偏移成像方法，其特征在于，所述的将震源波场分离为震源上行波场和震源下行波场，将检波点波场分离为检波点上行波场和检波点下行波场包括：

通过傅里叶变换将震源波场分离为震源上行波场和震源下行波场，将检波点波场分离为检波点上行波场和检波点下行波场。

3.根据权利要求1所述的逆时偏移成像方法，其特征在于，所述对分离后的震源波场和检波点波场进行逆时偏移成像处理，得到逆时偏移成像剖面包括：

将分离后的震源波场和检波点波场转换至时间域，得到逆时偏移的时间域成像条件；

根据所述的逆时偏移的时间域成像条件对分离后的震源和检波点波场进行逆时偏移成像，生成逆时偏移成像剖面。

4.一种逆时偏移成像装置，其特征在于，所述的装置包括：

波场信息获取单元，用于获取震源波场信息和检波点波场信息；

波场分离单元，用于根据获取的震源波场信息和检波点波场信息，将震源波场分离为震源上行波场和震源下行波场，将检波点波场分离为检波点上行波场和检波点下行波场；

成像单元，用于对分离后的震源波场和检波点波场进行逆时偏移成像处理，得到逆时偏移成像剖面；

所述的波场信息获取单元包括：

预处理模块，用于获取单炮数据，并对单炮数据进行叠前预处理；

外推处理模块，用于对叠前预处理输出的单炮数据的正向波场外推及逆时外推，获得各个时间步的震源波场信息和检波点波场信息。

5.根据权利要求4所述的逆时偏移成像装置，其特征在于，所述的波场分离单元通过傅里叶变换将震源波场分离为震源上行波场和震源下行波场，将检波点波场分离为检波点上行波场和检波点下行波场。

6.根据权利要求4所述的逆时偏移成像装置，其特征在于，所述成像单元包括：

时间域转换模块，用于将分离后的震源波场和检波点波场转换至时间域，得到逆时偏移的时间域成像条件；

逆时偏移成像模块，用于根据所述的逆时偏移的时间域成像条件对分离后的震源和检波点波场进行逆时偏移成像，生成逆时偏移成像剖面。