CN111190226B - 一种三维地震数据面波噪声压制方法 - Google Patents

Abstract

本发明公开了一种三维地震数据面波噪声压制方法，包括以下步骤：输入含面波噪声的三维地震数据、观测系统文件及Radon变换参数；构建三维圆锥Radon变换，将含面波噪声的三维地震数据体变换到圆锥Radon域；构建目标函数，采用最小平方方法，求解得到去除面波噪声的三维数据体地震道；判断所有数据体是否已经都已经去噪完成，如果未完成，则返回步骤2；如果已完成，输出去除面波噪音的三维地震数据体。本发明三维地震数据面波噪声压制方法，采用三维最小二乘圆锥Radon变换方法压制三维数据体中线性面波噪声，将三维地震数据体中的线性的面波噪声在每一时刻聚焦到了固定点上，使得三维数据体中的面波噪声得到了更准确地压制。

Description

一种三维地震数据面波噪声压制方法

技术领域

本发明属于石油勘探技术领域，具体涉及一种三维地震数据面波噪声压制方法。

背景技术

面波作为地震勘探中的一种相干噪声，广泛存在于地震记录中，属于一种具有能量强，频率低，震动延续时间长的规则干扰波，在地震记录上的表现比较明显，能够严重影响地震采集得到数据的分辨率和信噪比，妨碍解释人员对其进行正确地分析，因此，对地震信号中的面波进行压制就成为了地震数据处理阶段的一个重要环节。

其中，对于面波噪声的压制，大都采用Radon变换方法。在地震资料面波噪声压制处理中，较为常用的为积分路径为线性的Radon变换，该变换称为τ-p变换，或称之为倾斜叠加。为满足复杂介质数据处理的需求，相继出现了抛物Radon变换和双曲Radon变换。这些技术依赖于地震同相轴与Radon基函数之间剩余时差的相似性，例如线性、抛物线和双曲线Radon基函数等。这些相似性决定了地震同相轴在拉东域中的聚焦性，直接影响了Radon变换的地震数据处理质量。但是，三维数据体中的面波噪声呈圆锥形分布，这与最小二乘三维τ-p变换的基函数完全不匹配。这种不匹配使得传统三维τ-p Radon变换无法准确地压制三维地震数据体中的面波噪声。

发明内容

本发明的目的是为克服上述现有技术的不足，提供一种三维地震数据面波噪声压制方法，本发明方法相比于传统三维τ-p变换，三维地震数据体中的面波噪声得到了更准确地压制，设计合理，具有良好的效果。

为实现上述目的，本发明采用如下技术方案：

一种三维地震数据面波噪声压制方法，包括以下步骤：

步骤1：输入含面波噪声的三维地震数据、观测系统文件及Radon变换参数；

步骤2：构建三维圆锥Radon变换，将含面波噪声的三维地震数据体变换到圆锥Radon域；

步骤3：构建目标函数，采用最小平方方法，求解得到去除面波噪声的三维数据体地震道；

步骤4：判断所有数据体是否已经都已经去噪完成，如果未完成，则返回步骤2；如果已完成，输出去除面波噪音的三维地震数据体，去除面波噪声结束。

优选的，所述步骤2中，三维圆锥Radon变换的表达式为：

式(1)中，D(t,x,y)为三维地震数据，t为记录时间，(x,y)为空间坐标值，M_δ(τ,p_x,p_y)为圆锥拉东正变换系数，τ为截距时间，p_x为x测线方向的慢度，p_y为y测线方向的慢度，δ表示狄拉克函数，

为圆锥Radon变换的基函数。

优选的，所述步骤3包括以下步骤：

步骤31：对方程(1)的两边关于t进行一维傅里叶变换：

式(2)中，d为频率域的三维地震数据，m_δ为频率域的圆锥Radon变换的基函数，ω为角频率，i为虚数单位；

步骤32：将式(2)采用矩阵表示法表示，式(2)简写为

d(ω)＝L(ω)m_δ(ω) (3)

式(3)中，d(ω)为与角频率有关的频率域的三维地震数据的矩阵向量，m_δ(ω)为与角频率有关的频率域圆锥Radon变换的基函数；L(ω)为逆三维圆锥Radon变换矩阵，L(ω)由公式(4)表示：

式(4)中，m和n为计数变量，x_R为R处的x坐标，x_S为S处的x坐标，y_T为T处的y坐标，y_U为U处的y坐标，

和

分别表示x_R和y_T处的慢度，M和N为计算点数，分别由式(5)和(6)计算得到：

M＝n_xn_y (5)

式(5)中，n_x为纵测线方向上的道数，n_y为横测线方向上的道数；式(6)中，

为纵测线方向上的慢度采样数，

为横测线方向上的慢度采样数；

式(4)中，R、S、T和U均为中间变量，分别由公式(7)～公式(10)计算得到：

U＝floor(m/n_x) (10)

式(7)～(8)中，％为取余数符号；式(9)～(10)中，floor()为向下舍入函数；

步骤33：对于方程(3)，构建如下所示的目标函数：

式(11)中，|| ||₂表示L2范数值；

步骤34：采用最小二乘法对式(11)进行求解；等式(11)的三维圆锥Radon变换L2范数的最小二乘解为

m_δ＝L^H(LL^H)^-1d (12)

式(12)中，上标H表示共轭；令中间矩阵z由下式表示：

z＝(LL^H)^-1d (13)

令中间变量A由下式表示：

A＝LL^H (14)

则方程(12)可改写为：

m_δ＝L^Hz (15)

方程(13)可改写为：

Az＝d (16)

矩阵A的元素可表示为：

式(17)中，

表示不同点的x坐标，

表示不同点的y坐标，a_0,0,...,a_M-1,M-1为矩阵A的元素值；

步骤35：采用共轭梯度算法求解式(16)；

步骤36：将式(16)的解与L^H相乘，得到最终的三维最小二乘圆锥Radon变换结果。

本发明的有益效果是：

本发明三维地震数据面波噪声压制方法，采用三维最小二乘圆锥Radon变换方法压制三维数据体中线性面波噪声；与常规三维最小二乘τ-p变换不同，本发明中的三维最小二乘圆锥Radon变换方法将三维地震数据体中的线性的面波噪声在每一时刻聚焦到了固定点上，使得三维数据体中的面波噪声得到了更准确地压制。

附图说明

构成本申请的一部分的说明书附图用来提供对本申请的进一步理解，本申请的示意性实施例及其说明用于解释本申请，并不构成对本申请的不当限定。

图1是本发明三维地震数据面波噪声压制方法的流程示意图。

具体实施方式

应该指出，以下详细说明都是例示性的，旨在对本申请提供进一步的说明。除非另有指明，本文使用的所有技术和科学术语具有与本申请所属技术领域的普通技术人员通常理解的相同含义。

需要注意的是，这里所使用的术语仅是为了描述具体实施方式，而非意图限制根据本申请的示例性实施方式。如在这里所使用的，除非上下文另外明确指出，否则单数形式也意图包括复数形式，此外，还应当理解的是，当在本说明书中使用术语“包含”和/或“包括”时，其指明存在特征、步骤、操作、器件、组件和/或它们的组合。

下面结合附图和实施例对本发明进一步说明。

如图1所示，一种三维地震数据面波噪声压制方法，包括以下步骤：

步骤1：输入含面波噪声的三维地震数据、观测系统文件及Radon变换参数；

步骤2：构建三维圆锥Radon变换，将含面波噪声的三维地震数据体变换到圆锥Radon域；

三维地震数据体的面波噪声呈圆锥形分布，这与传统最小二乘三维τ-p变换的基函数(线性Radon基函数、抛物Radon基函数、双曲Radon基函数)不匹配；构建三维圆锥Radon变换，将含面波噪声的三维地震数据体变换到圆锥Radon域；

具体地，所述步骤2中，三维圆锥Radon变换的表达式为：

式(1)中，D(t,x,y)为三维地震数据，t为记录时间，(x,y)为空间坐标值，M_δ(τ,p_x,p_y)为圆锥拉东正变换系数，τ为截距时间，p_x为x测线方向的慢度，p_y为y测线方向的慢度，δ表示狄拉克函数，

为圆锥Radon变换的基函数。

步骤3：构建目标函数，采用最小平方方法，求解得到去除面波噪声的三维数据体地震道；

具体地，所述步骤3包括以下步骤：

步骤31：对方程(1)的两边关于t进行一维傅里叶变换：

式(2)中，d为频率域的三维地震数据，m_δ为频率域的圆锥Radon变换的基函数，ω为角频率，i为虚数单位；

步骤32：对于一个给定的频率ω，方程(2)是一个复杂的线性系统；因此，将式(2)采用矩阵表示法表示，式(2)简写为

d(ω)＝L(ω)m_δ(ω) (3)

式(3)中，d(ω)为与角频率有关的频率域的三维地震数据的矩阵向量，m_δ(ω)为与角频率有关的频率域圆锥Radon变换的基函数；L(ω)为逆三维圆锥Radon变换矩阵，L(ω)由公式(4)表示：

式(4)中，m和n为计数变量，x_R为R处的x坐标，x_S为S处的x坐标，y_T为T处的y坐标，y_U为U处的y坐标，p_xR和p_yT分别表示x_R和y_T处的慢度，M和N为计算点数，分别由式(5)和(6)计算得到：

M＝n_xn_y (5)

式(5)中，n_x为纵测线方向上的道数，n_y为横测线方向上的道数；式(6)中，

为纵测线方向上的慢度采样数，

为横测线方向上的慢度采样数；

式(4)中，R、S、T和U均为中间变量，分别由公式(7)～公式(10)计算得到：

式(7)～(8)中，％为取余数符号；式(9)～(10)中，floor()为向下舍入函数；

步骤33：对于方程(3)，构建如下所示的目标函数：

式(11)中，|| ||₂表示L2范数值；

步骤34：采用最小二乘法对式(11)进行求解；

等式(11)的三维圆锥Radon变换L2范数的最小二乘解为

m_δ＝L^H(LL^H)^-1d (12)

式(12)中，上标H表示共轭；令中间矩阵z由下式表示：

z＝(LL^H)^-1d (13)

令中间变量A由下式表示：

A＝LL^H (14)

则方程(12)可改写为：

m_δ＝L^Hz (15)

方程(13)可改写为：

Az＝d (16)

矩阵A的元素可表示为：

式(17)中，

表示不同点的x坐标，

表示不同点的y坐标，a_0,0,...,a_M-1,M-1为矩阵A的元素值；

步骤35：由式(17)可知，矩阵A为一厄米矩阵，因此采用共轭梯度算法求解式(16)；

步骤36：将式(16)的解与L^H相乘，得到最终的三维最小二乘圆锥Radon变换结果。

步骤4：判断所有数据体是否已经都已经去噪完成，如果未完成，则返回步骤2；如果已完成，输出去除面波噪音的三维地震数据体，去除面波噪声结束。

实施例：

本发明三维地震数据面波噪声压制方法，应用于三维地震数据体，取得了理想的计算效果。

采用本发明三维地震数据面波噪声压制方法对含面波噪声的三维地震数据体进行面波噪声压制，与采用传统三维最小二乘τ-p变换方法相比，具有以下优势：

采用本发明三维地震数据面波噪声压制方法得到的零截距时间片与采用传统三维最小二乘τ-p变换方法得到的零截距时间片进行对比发现，本发明方法可以将三维面波噪声和反射波聚焦到每个时间片上的一个点，而传统三维最小二乘τ-p变换将其分散在椭圆上；因此，本申请三维地震数据面波噪声压制方法在慢度方向上具有比三维最小二乘τ-pRadon变换更高的分辨率。

采用本发明三维地震数据面波噪声压制方法得到的压制面波噪声的三维数据体、得到的分离的面波噪声与采用传统三维最小二乘τ-p变换方法得到的压制面波噪声的三维数据体、得到的分离的面波噪声进行对比发现，本发明方法在消除三维面波噪声和减少反射波损耗方面取得了更好的效果。

本发明三维地震数据面波噪声压制方法，采用三维最小二乘圆锥Radon变换方法压制三维数据体中线性面波噪声；与常规三维最小二乘τ-p变换不同，本发明中的三维最小二乘圆锥Radon变换方法将三维地震数据体中的线性的面波噪声在每一时刻聚焦到了固定点上，使得三维数据体中的面波噪声得到了更准确地压制。

上述虽然结合附图对本发明的具体实施方式进行了描述，但并非对本发明的限制，所属领域技术人员应该明白，在本发明的技术方案的基础上，本领域技术人员不需要付出创造性劳动即可做出的各种修改或变形仍在本发明的保护范围以内。

Claims (3)

1.一种三维地震数据面波噪声压制方法，其特征是，包括以下步骤：

步骤1：输入含面波噪声的三维地震数据、观测系统文件及Radon变换参数；

步骤2：构建三维圆锥Radon变换，将含面波噪声的三维地震数据体变换到圆锥Radon域；

步骤3：构建目标函数，采用最小平方方法，求解得到去除面波噪声的三维数据体地震道；

步骤4：判断所有数据体是否已经都已经去噪完成，如果未完成，则返回步骤2；如果已完成，输出去除面波噪音的三维地震数据体，去除面波噪声结束。

2.如权利要求1所述的三维地震数据面波噪声压制方法，其特征是，所述步骤2中，三维圆锥Radon变换的表达式为：

式(1)中，D(t,x,y)为三维地震数据，t为记录时间，(x,y)为空间坐标值，M_δ(τ,p_x,p_y)为圆锥拉东正变换系数，τ为截距时间，p_x为x测线方向的慢度，p_y为y测线方向的慢度，δ表示狄拉克函数，

为圆锥Radon变换的基函数。

3.如权利要求2所述的三维地震数据面波噪声压制方法，其特征是，所述步骤3包括以下步骤：

步骤31：对方程(1)的两边关于t进行一维傅里叶变换：

式(2)中，d为频率域的三维地震数据，m_δ为频率域的圆锥Radon变换的基函数，ω为角频率，i为虚数单位；

步骤32：将式(2)采用矩阵表示法表示，式(2)简写为

d(ω)＝L(ω)m_δ(ω) (3)

式(3)中，d(ω)为与角频率有关的频率域的三维地震数据的矩阵向量，m_δ(ω)为与角频率有关的频率域圆锥Radon变换的基函数；L(ω)为逆三维圆锥Radon变换矩阵，L(ω)由公式(4)表示：

式(4)中，m和n为计数变量，x_R为R处的x坐标，x_S为S处的x坐标，y_T为T处的y坐标，y_U为U处的y坐标，

和

分别表示x_R和y_T处的慢度，M和N为计算点数，分别由式(5)和(6)计算得到：

M＝n_xn_y (5)

式(5)中，n_x为纵测线方向上的道数，n_y为横测线方向上的道数；式(6)中，

为纵测线方向上的慢度采样数，

为横测线方向上的慢度采样数；

式(4)中，R、S、T和U均为中间变量，分别由公式(7)～公式(10)计算得到：

S＝m％n_x (8)

U＝floor(m/n_x) (10)

式(7)～(8)中，％为取余数符号；式(9)～(10)中，floor()为向下舍入函数；

步骤33：对于方程(3)，构建如下所示的目标函数：

式(11)中，|| ||₂表示L2范数值；

步骤34：采用最小二乘法对式(11)进行求解；等式(11)的三维圆锥Radon变换L2范数的最小二乘解为

m_δ＝L^H(LL^H)^-1d (12)

式(12)中，上标H表示共轭；令中间矩阵z由下式表示：

z＝(LL^H)^-1d (13)

令中间变量A由下式表示：

A＝LL^H (14)

则方程(12)可改写为：

m_δ＝L^Hz (15)

方程(13)可改写为：

Az＝d (16)

矩阵A的元素可表示为：

式(17)中，

表示不同点的x坐标，

表示不同点的y坐标，a_0,0,...,a_M-1,M-1为矩阵A的元素值；

步骤35：采用共轭梯度算法求解式(16)；

步骤36：将式(16)的解与L^H相乘，得到最终的三维最小二乘圆锥Radon变换结果。

Images