CN105093305A - 一种确定地震数据反射界面位置的方法 - Google Patents

Abstract

本发明提供了一种确定地震数据反射界面位置的方法，属于石油地球物理勘探领域。本方法包括：(1)设地震波为f(t)，对其进行Hilbert变换，求取实地震信号的解析信号，得到地震波f(t)的瞬时振幅、瞬时相位剖面θ(t)和瞬时频率；(2)基于地震波瞬时相位剖面θ(t)进行一阶差分运算，得到瞬时相位一阶差分G(t)；(3)对地震数据f(t)进行极性反转处理得到F(t)＝-f(t)，对F(t)进行Hilbert变换，得到瞬时相位剖面θ₁(t)，再对瞬时相位剖面θ₁(t)进行一阶差分运算，得到极性反转后的地震波瞬时相位一阶差分j(t)，对其求相反数得到J(t)＝-j(t)。

Description

一种确定地震数据反射界面位置的方法

技术领域

本发明属于石油地球物理勘探领域，具体涉及一种确定地震数据反射界面位置的方法。

背景技术

地球物理勘探中，反射波的传播在地壳介质的分界面处会产生速度和方向的变化，这些都是地球物理信号的奇异性，判断出奇异性的大小和位置就可以对异常现象作出解释，所以反射界面位置的确定对地球物理勘探研究具有特殊意义。当地下岩层存在分界面时，返回地面的反射波就会引起地面振动，这时，一个排列上所有地震道也就都会跳动，因而会在地震记录上出现波峰、波谷构成的同相轴。一般地层含油气后其与上下岩层差别较大，这时在地震记录上会出现振幅比较强的同相轴。在地震记录上不仅能找到可能的含油气地层(当然还要参考其他信息)，还能计算出含油气层的倾斜角度及埋藏深度。所以，研究反射界面的位置有着很重要的价值。

为了有效提取地震剖面中的同相轴，地球物理学界提出了多种有效的计算方法，如神经网络、小波分析、模式识别、混沌理论、图像边缘检测等等。但是利用神经网络进行同相轴追踪，效率低以及大量的迭代运算满足不了实时性处理。为了有效解决上述方法存在的不足，很多学者便开始采用简单的图像处理方法，如选择具有良好特性的、改进的Sobel算子对地震剖面图像进行边缘检测，但得到的同相轴比较模糊、分辨率不高，且由于噪声因素引起的边缘信息也被提取出来，造成同相轴的错判。而利用边缘检测的Canny算子的方法，则先对图像边缘提取，然后提取边缘检测的结果中的灰度突变区域的包络线的中轴作为同相轴的检测追踪结果，当分辨率较高时可以得到比较准确的同相轴检测结果，但Canny算子对噪声比较敏感，且参数的选取关系着是否精确提取，而一般地震数据先验知识都是不知道的，所以参数的选择上存在很大的误差。

发明内容

本发明的目的在于解决上述现有技术中存在的难题，提供一种确定地震数据反射界面位置的方法，针对利用神经网络进行同相轴追踪，效率低以及大量的迭代运算满足不了实时性处理的缺点，通过对地震数据进行Hilbert变换，得到虚地震记录，进而求得瞬时振幅和瞬时相位，并将瞬时振幅与极性反转前后地震记录的瞬时相位一阶差分之和进行加权，得到复合属性，利用这种复合属性可以比较准确、可靠地识别地震反射界面的位置，也可得到准确的同相轴检测结果。

本发明是通过以下技术方案实现的：

一种确定地震数据反射界面位置的方法，包括：

一种确定地震数据反射界面位置的方法，包括：

(1)设地震波为f(t)，对其进行Hilbert变换，求取实地震信号的解析信号，得到地震波f(t)的瞬时振幅、瞬时相位剖面θ(t)和瞬时频率；

(2)基于地震波瞬时相位剖面θ(t)进行一阶差分运算，得到瞬时相位一阶差分G(t)，则

G(t)＝θ(t+1)-θ(t)；

(3)对地震数据f(t)进行极性反转处理得到F(t)＝-f(t)，对F(t)进行Hilbert变换，得到瞬时相位剖面θ₁(t)，再对瞬时相位剖面θ₁(t)进行一阶差分运算，得到极性反转后的地震波瞬时相位一阶差分j(t)，对其求相反数得到J(t)＝-j(t)；

(4)将步骤(2)和步骤(3)的结果作为输入，利用下式对两个结果进行加权相加运算：

M(t)＝a₁*G(t)+a₂*J(t)

其中a₁+a₂＝1，a₁代表地震数据中波谷占的权重，a₂代表地震数据中波峰占的权重；

(5)将步骤(4)的结果M(t)与原始地震波f(t)相乘得到新的属性Q(t)：Q(t)＝M(t)*f(t)，并选取阈值滤去新的属性Q(t)中的较低值，获得准确的反射界面位置。

所述新的属性Q(t)是由不同振幅的脉冲组成，其中小振幅脉冲为干扰信号，而大振幅脉冲为真实反射界面的位置。

所述选取阈值滤去新的属性Q(t)中的较低值，获得准确的反射界面位置是这样实现的：

通过实验选取阈值；

将新的属性Q(t)中的振幅小于或等于阈值的脉冲滤去，得到的结果即为准确的反射界面位置。

与现有技术相比，本发明的有益效果是：

(1)本发明基于复合属性的固有性质进行二维地震数据反射界面识别，由于瞬时相位是一种独立于地震振幅，并且不受其影响的地震属性参数，可以清晰地反映一些弱反射波及深层反射波，更多地获取到深层构造区域信息，因此利用基于瞬时相位的复合振幅法进行反射界面位置的识别，能更好地反映构造几何形态的细节。

(2)本发明算法逻辑简单，易于实现；且计算速度快、费用少，并可检测地下构造形态及地下地层单元的形状，可在石油地球物理勘探领域推广应用。

附图说明

图1本发明方法的步骤框图。

图2反射系数序列模型及其合成地震记录，其中实线为合成记录，虚线为复地震道虚部。

图3反射系数序列模型合成记录的瞬时相位。

图4合成记录瞬时相位的的差分属性。

图5反射系数序列模型合成记录极性反转后记录，实线为合成地震记录极性反转后的记录，虚线为极性反转后复地震道的虚部。

图6极性反转后合成记录的瞬时相位。

图7极性反转后合成记录瞬时相位的的差分属性。

图8反射系数序列模型合成记录瞬时相位的差分属性与极性反转后合成记录瞬时相位差分属性的相反数加权相加的结果。

图9反射系数序列模型合成记录瞬时相位的差分属性与极性反转后合成记录瞬时相位差分属性的相反数加权相加后加阈值的结果。

图10反射系数序列模型反射界面识别结果。

图11楔形体地质模型。

图12楔形体地质模型的合成地震记录波形图。

图13楔形体地质模型的合成地震记录(变密度显示)。

图14楔形体地质模型合成地震记录的复地震道虚部。

图15楔形体地质模型合成地震记录的瞬时相位。

图16楔形体地质模型合成地震记录瞬时相位一阶差分属性。

图17楔形体地质模型合成地震记录极性反转后剖面。

图18楔形体地质模型合成地震记录极性反转后复地震道虚部。

图19合成地震记录极性反转后的瞬时相位。

图20合成地震记录极性反转后瞬时相位的一阶差分。

图21本发明识别楔形体地质模型的合成地震记录反射界面的结果。

具体实施方式

下面结合附图对本发明作进一步详细描述：

本发明利用瞬时相位属性的固定优势，通过特殊运算得到复合属性，通过复合属性剖面来准确、可靠地识别地震反射界面的位置。

本发明采用Hilbert变换，通过求解实地震信号的解析信号，将可反映地震信号局部变化情况的地震波的瞬时振幅、瞬时相位和瞬时频率等信息分离开，将瞬时振幅与极性反转前后地震记录的瞬时相位一阶差分之和进行加权，得到复合属性，并利用该属性结果进行反射界面位置的识别。

如图1所示，所述方法包括以下步骤：

(1)设地震数据为f(t)，对其进行Hilbert变换，求取实地震信号的解析信号，得到地震波的瞬时振幅、瞬时相位、瞬时频率等信息。

Hilbert变换的基本原理：某一实函数f(t)的Hilbert变换定义为：

$\tilde{f}\left(t\right)=H\left[f\left(t\right)\right]=\frac{1}{\mathrm{\π}}{\∫}_{-\∞}^{+\∞}\frac{f\left(\mathrm{\τ}\right)}{t-\mathrm{\τ}}\mathrm{d\τ}---\left(1\right)$

即 $\tilde{f}\left(t\right)=H\left[f\left(t\right)\right]=f\left(t\right)*\frac{1}{\mathrm{\πt}}=f\left(t\right)*h\left(t\right)---\left(2\right)$

H表示Hilbert变，*表示卷积，因此Hilbert变换相当于通过一个冲激响应为1/π的线性网络，该网络也称为Hilbert滤波器。实值函数的Hilbert变换的计算也可以在频率实现，如果x(t)，y(t)，h(t)的Fourier变换分别是X(f)、Y(f)、H(f)，那么可以推导出Y(f)与X(f)、H(f)，那么可以推导出Y(f)＝X(f)H(f)，式中X(f)为实信号x(t)的频谱；H(f)为h(t)的频谱；Y(f)为y(t)的频谱。其中

$H\left(f\right)=\left\{\begin{array}{cc}2& f\&GreaterEqual;0\\ 0& f<0\end{array}\right.---\left(3\right)$

综合(1)-(3)式，复地震道可表示为：

z(t)＝A(t)e^iθ(t)(4)

式中A(t)为瞬时振幅；θ(t)为瞬时相位。A(t)和θ(t)是两个基本的瞬时属性。瞬时振幅也称振幅包络或反射强度，相当于地震道包络函数。它反映了地震波能量的瞬时变化情况，与地震相位无关，可用来判断与岩性有关的地质体。瞬时振幅可表示为：

$A\left(t\right)=\sqrt{x^2\left(t\right)+y^2\left(t\right)}---\left(5\right)$

瞬时相位是一种独立于地震振幅，并且不受其影响的地震属性参数，因此瞬时相位能更好地反映构造几何形态的细节，不仅其本身是一种瞬时属性，而且是求得瞬时频率的基础。瞬时相位可追踪连续性差的弱反射波及极性变化的反射波。瞬时相位可表示为式(6)。

θ(t)＝tan^-1[y(t)/x(t)](6)

(2)基于地震波瞬时相位剖面θ(t)进行一阶差分运算，得到瞬时相位一阶差分G(t)，则

G(t)＝θ(t+1)-θ(t)(7)

(3)对地震数据f(t)进行极性反转处理得到F(t)＝-f(t)，用步骤(1)(2)方法得到极性反转后的地震波瞬时相位一阶差分j(t)，对其求相反数得到J(t)＝-j(t)；

(4)将步骤(2)和步骤(3)的结果作为输入，将两个结果进行加权相加运算如式(8)，

M(t)＝a₁*G(t)+a₂*J(t)(8)

其中a₁+a₂＝1。根据G(t)和J(t)在数据整体评估中占的重要作用设定比例，确定a1和a2的值，a1和a2的值加起来应为1，a1代表地震数据中波谷占的权重，a2代表地震数据中波峰占的权重。

(5)将步骤(4)的结果M(t)与原始地震数据f(t)相乘，并选取合适的阈值滤去较低值，获得准确的反射界面位置。

图1示出本发明二维地震数据反射界面识别的方法的具体流程。下面参照附图并结合实施例对本发明做进一步说明。

设计一个三层的模型，第一层与第三层速度均为2500m/s，第二层速度为3000m/s，即第一界面为正反射，第二界面为负反射。用20HZ的子波模拟地震道如图2中实线，再通过Hilbert变换将实地震道变成复地震道，复地震道虚部如图2中虚线，实地震道的瞬时相位为图3，对瞬时相位进行差分运算后结果如图4。从图2到图4中我们可以得出复地震道虚部相对实地震道发生了90。相移，当反射振幅为正时，复地震道虚部由负到正连续变化，在波峰位置时为零，合成地震记录的相位从-π/2到π/2连续变化，且在波峰位置相位为零，对其相位求差分后的值很小；当反射振幅为正时，复地震道虚部由正到负连续变化，在波谷位置虚地震记录值为零，合成记录的相位由π/2连续变化到π，在波谷位置由π直接转变为-π，其后连续变化到-π/2，对其相位求差分后在波谷位置的值为最大值。这一变化说明，当地层厚度大于四分之一波长时，相位差分后的异常位置为波谷的位置。而为了得到波峰的位置，我们对合成地震记录进行极性反转处理如图5中实线所示，用相同的方法得到其瞬时相位如图6，并对瞬时相位进行差分，得到合成记录极性反转后的波谷位置如图7所示，让其相位差分结果的相反数代表原始合成记录的波峰位置后与原始地震记录的相位差分结果进行相加运算，这样就准确的找到波峰和波谷的位置如图8所示。将其结果与原始合成地震记录进行加权相乘，即可得到模型的射界面位置如图9所示。选择阈值滤去低值后可得到准确的反射界面如图10所示。

下面将上述方法用于楔形体模型顶底界面的识别中：

设计一个楔形体模型如图11所示，采用20HZ雷克子波，进行自激自收地震正演模拟，然后利用该发明进行反射界面的识别。图12为合成地震记录波形图，图13为合成地震记录的变密度显示，首先，根据步骤一说明的方法，对楔形体模型合成地震数据运用Hilbert变换得到复地震道虚部(如图14所示)、瞬时相位(如图15所示)，然后根据步骤二说明的方法得到瞬时相位的差分属性(如图16所示)，接下来根据步骤三得到为原始合成记录极性反转后的记录(如图17所示)，再根据步骤一、二分别得到极性反转后记录的复地震道虚部(如图18所示)、瞬时相位(如图19所示)和瞬时相位的差分属性(如图20所示)。得到上述结果后，根据步骤四方法，将极性反转后记录的瞬时相位差分属性值乘以-1后与原始合成地震记录瞬时相位的差分属性进行加权相加运算，最后，根据步骤五说明的方法，将步骤四得到的结果与原始合成地震记录加权相乘，滤去小值，即可得到楔形体模型反射界面识别的结果(如图21所示)。从图21可以看出通过这种方法当薄层厚度大于四分之一波长时可以准确的得到楔形体模型的顶、底反射界面位置。

本发明采用Hilbert变换，通过求解实地震信号的解析信号，将可反映地震信号局部变化情况的地震波的瞬时振幅、瞬时相位和瞬时频率等信息分离开，将瞬时振幅与极性反转前后地震记录的瞬时相位一阶差分之和进行加权，得到复合属性，并利用该属性结果进行反射界面位置的识别。主要用于地震反射界面的位置的识别，通过该方法可以检测地震反射同相轴，提高剖面信噪比和分辨率，使反射同相轴在检测后的剖面上很容易被识别和追踪。

上述技术方案只是本发明的一种实施方式，对于本领域内的技术人员而言，在本发明公开了应用方法和原理的基础上，很容易做出各种类型的改进或变形，而不仅限于本发明上述具体实施方式所描述的方法，因此前面描述的方式只是优选的，而并不具有限制性的意义。

Claims (3)

1.一种确定地震数据反射界面位置的方法，其特征在于：所述方法包括：

(1)设地震波为f(t)，对其进行Hilbert变换，求取实地震信号的解析信号，得到地震波f(t)的瞬时振幅、瞬时相位剖面θ(t)和瞬时频率；

(2)基于地震波瞬时相位剖面θ(t)进行一阶差分运算，得到瞬时相位一阶差分G(t)，则

G(t)＝θ(t+1)-θ(t)；

(3)对地震数据f(t)进行极性反转处理得到F(t)＝-f(t)，对F(t)进行Hilbert变换，得到瞬时相位剖面θ₁(t)，再对瞬时相位剖面θ₁(t)进行一阶差分运算，得到极性反转后的地震波瞬时相位一阶差分j(t)，对其求相反数得到J(t)＝-j(t)；

(4)将步骤(2)和步骤(3)的结果作为输入，利用下式对两个结果进行加权相加运算：

M(t)＝a₁*G(t)+a₂*J(t)

其中a₁+a₂＝1，a₁代表地震数据中波谷占的权重，a₂代表地震数据中波峰占的权重；

(5)将步骤(4)的结果M(t)与原始地震波f(t)相乘得到新的属性Q(t)：Q(t)＝M(t)*f(t)，并选取阈值滤去新的属性Q(t)中的较低值，获得准确的反射界面位置。

2.根据权利要求1所述的确定地震数据反射界面位置的方法，其特征在于：所述新的属性Q(t)是由不同振幅的脉冲组成，其中小振幅脉冲为干扰信号，而大振幅脉冲为真实反射界面的位置。

3.根据权利要求1所述的确定地震数据反射界面位置的方法，其特征在于：所述选取阈值滤去新的属性Q(t)中的较低值，获得准确的反射界面位置是这样实现的：

通过实验选取阈值；

将新的属性Q(t)中的振幅小于或等于阈值的脉冲滤去，得到的结果即为准确的反射界面位置。