CN113466940B - 一种深层叠前地震数据噪声压制方法及系统 - Google Patents
一种深层叠前地震数据噪声压制方法及系统 Download PDFInfo
- Publication number
- CN113466940B CN113466940B CN202110674853.3A CN202110674853A CN113466940B CN 113466940 B CN113466940 B CN 113466940B CN 202110674853 A CN202110674853 A CN 202110674853A CN 113466940 B CN113466940 B CN 113466940B
- Authority
- CN
- China
- Prior art keywords
- noise
- data
- sub
- seismic data
- module
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 50
- 230000001629 suppression Effects 0.000 title claims abstract description 39
- 238000012545 processing Methods 0.000 claims abstract description 50
- 230000009466 transformation Effects 0.000 claims abstract description 28
- 230000001427 coherent effect Effects 0.000 claims abstract description 18
- 238000000638 solvent extraction Methods 0.000 claims abstract description 7
- 230000000903 blocking effect Effects 0.000 claims description 13
- 238000001914 filtration Methods 0.000 claims description 12
- 239000000126 substance Substances 0.000 claims description 5
- 238000003825 pressing Methods 0.000 claims description 3
- 238000005192 partition Methods 0.000 claims description 2
- 230000010287 polarization Effects 0.000 claims description 2
- 230000011218 segmentation Effects 0.000 claims 1
- 238000003860 storage Methods 0.000 description 14
- 238000004590 computer program Methods 0.000 description 9
- 230000006870 function Effects 0.000 description 7
- 238000010586 diagram Methods 0.000 description 6
- 230000008569 process Effects 0.000 description 6
- AYFVYJQAPQTCCC-GBXIJSLDSA-N L-threonine Chemical compound C[C@@H](O)[C@H](N)C(O)=O AYFVYJQAPQTCCC-GBXIJSLDSA-N 0.000 description 5
- 238000004364 calculation method Methods 0.000 description 5
- 230000009286 beneficial effect Effects 0.000 description 3
- 238000004422 calculation algorithm Methods 0.000 description 3
- 230000015556 catabolic process Effects 0.000 description 3
- 238000006731 degradation reaction Methods 0.000 description 3
- 239000004215 Carbon black (E152) Substances 0.000 description 2
- 229930195733 hydrocarbon Natural products 0.000 description 2
- 150000002430 hydrocarbons Chemical class 0.000 description 2
- 238000004519 manufacturing process Methods 0.000 description 2
- 238000009966 trimming Methods 0.000 description 2
- 230000001133 acceleration Effects 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 230000015572 biosynthetic process Effects 0.000 description 1
- 230000021615 conjugation Effects 0.000 description 1
- 238000007796 conventional method Methods 0.000 description 1
- 238000005520 cutting process Methods 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 238000009826 distribution Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 230000000717 retained effect Effects 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/36—Effecting static or dynamic corrections on records, e.g. correcting spread; Correlating seismic signals; Eliminating effects of unwanted energy
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/32—Transforming one recording into another or one representation into another
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/36—Effecting static or dynamic corrections on records, e.g. correcting spread; Correlating seismic signals; Eliminating effects of unwanted energy
- G01V1/364—Seismic filtering
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/20—Trace signal pre-filtering to select, remove or transform specific events or signal components, i.e. trace-in/trace-out
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/30—Noise handling
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/40—Transforming data representation
- G01V2210/48—Other transforms
Landscapes
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Acoustics & Sound (AREA)
- Environmental & Geological Engineering (AREA)
- Geology (AREA)
- General Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Geophysics (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
本发明公开了一种深层叠前地震数据噪声压制方法及系统，对原始叠前数据进行分块；对分块后的其中一块原始叠前地震数据利用中值估计算子估算该块数据的噪声方差等级；对估计过噪声方差等级的原始叠前地震数据分块进行高维连续小波变换；依据相干干扰的角度确定分块需要保留的连续小波变换变换结果角度范围；根据分块的噪声标准差确定其需要保留的子带内的阈值；利用确定好的每个保留子带的阈值对需要保留子带中的小波系数进行阈值处理；对阈值处理后的小波变换子带进行高维小波逆变换得到重建的高信噪比信号；重复以上步骤直至所有分块处理完成，将得到的所有高信噪比分块重新拼接为完整的高信噪比数据作为完成噪声压制的深层叠前地震数据。
Description
技术领域
本发明属于勘探地球物理技术领域，具体涉及一种深层叠前地震数据噪声压制方法及系统。
背景技术
近年来，随着浅层油气资源被大量发掘，油气资源勘探逐渐转向地下深层区域。地震波在地层中衰减较快，深层反射能量弱，致使深层地震资料信噪比较低。另外，面波在深层会产生更大影响，进一步降低地震资料信噪比，虽然可以通过叠加多次采集的地震资料(叠前资料)得到信噪比更高的叠后资料，但众多叠前处理以及分析仍然需要高信噪比的叠前数据。
现有技术1：基于预测的滤波方法(如FX，TX，多项式拟合等)。此类方法利用相邻道间的线性预测关系来从含噪地震信号中估计出反射信号。此类方法在目标地层信号为线性同相轴时表现良好，而为非线性同相轴时去噪结果较差；此外，此类方法并不具有区分同相轴角度的分辨能力。
现有技术2：基于高维离散小波域的滤波方法。离散小波变换适宜对点奇异目标进行稀疏表示，通过阈值处理等方法即可实现对随机噪声较好的压制。但离散小波原子具有各项同性的特点，这使得其无法表示信号的方向信息，因此很难处理具备一定方向的结构性噪声。
发明内容
本发明所要解决的技术问题在于针对上述现有技术中的不足，提供一种深层叠前地震数据噪声压制方法及系统，在高维小波域中对深层叠前地震信号进行方向与阈值滤波，从而得到高信噪比的地震资料。
本发明采用以下技术方案：
一种深层叠前地震数据噪声压制方法，包括以下步骤：
S1、获取原始叠前数据并进行分块；
S2、利用中值估计算子估算步骤S1分块后的任意一块原始叠前地震数据的噪声方差等级；
S3、对步骤S2估算的原始叠前地震数据进行高维连续小波变换，得到一系列高维连续小波变换的系数子带；
S4、依据相干干扰的角度确定步骤S3得到的系数子带中需要保留的高维连续小波变换系数子带；
S5、根据步骤S2中得到的噪声方差等级确定步骤S4中需要保留的高维连续小波变换系数子带内的阈值
S6、利用步骤S5确定的阈值
对需要保留的高维连续小波变换系数子带中的小波系数进行阈值处理；
S7、对步骤S6阈值处理后的小波变换子带进行高维小波逆变换得到重建的高信噪比信号
S8、重复步骤S2～S7直到步骤S1获取的所有分块全部处理完毕，然后将得到的所有高信噪比信号
重新拼接为完整的高信噪比叠前地震数据。
具体的，步骤S1中，将原始叠前地震数据分为多个N×N×N的数据块，且在数据块之间设置有重叠，分块数为：
其中，n为重叠量，
表示向上取整，Nx×Ny×Nz为三维叠前数据体。
具体的，步骤S2中，噪声标准差等级
具体为：
其中，MEDIAN为取中值的算子，F(k)为对分块做三维傅里叶变换得到的结果，[kmin,kmax]为F(k)的高频范围。
具体的，步骤S3中，对任意选取的一个未完成处理的分块f(x)进行高维连续小波变换得到高频子带CWT(f；b,a,θ,φ)如下：
其中，
表示原数据块通过高维傅立叶变换得到的波数域数据，
为在参数组(a,θ,φ；k)下的小波原子通过高维傅立叶变换得到的波数域，
为
的共轭，IFFT表示逆傅里叶变换。
具体的，步骤S4中，通过相干干扰角度确定数据体中的噪声倾角θcoherent_noise；保留连续小波变换变换结果CWT(f；b,a,θ,φ)中角度θ小于θcoherent_noise的子带，即只保留[-θcoherent_noise+Δθ,θcoherent_noise-Δθ]角度范围内的子带。
具体的，步骤S5中，各高频子带的阈值thre(a)为：
其中，
是每个子带的阈值系数，其取值与1/a3成正比；c为一控制噪声压制程度的常数，通常情况下取为3，
为标准差。
具体的，步骤S6中，根据步骤S5中确定的阈值系数d进行硬限幅滤波，具体为：
其中，CWTthr(f；b,a,θ,φ)为经过阈值滤波的小波系数，CWT(f；b,a,θ,φ)为未经滤波的小波系数。
具体的，步骤S7中，高信噪比信号
其中，CWTthr(f；b,a,θ,φ)为阈值滤波后的小波系数，a为小波原子的尺度因子，φ为小波原子在方位角上的极化方向的角度，θ为小波原子在倾角上计划方向的角度。
具体的，步骤S8具体为：将步骤S1中分割出的所有数据块进行步骤S2至步骤S7的操作后，重新按照分割时的排列拼接起得到噪声压制后的数据；相邻的两个数据块在拼接前，分别去掉对应方向上
的点，n为重叠量。
本发明的另一个技术方案是，一种深层叠前地震数据噪声压制系统，包括：
分块模块，获取原始叠前数据并进行分块；
变换模块，利用中值估计算子估算分块模块分块后的任意一块原始叠前地震数据的噪声方差等级；
估算模块，对变换模块估算的原始叠前地震数据进行高维连续小波变换，得到一系列高维连续小波变换的系数子带；
选择模块，依据相干干扰的角度确定估算模块得到的系数子带中需要保留的高维连续小波变换系数子带；
阈值模块，根据变换模块中得到的噪声方差等级确定选择模块中需要保留的高维连续小波变换系数子带内的阈值
处理模块，利用阈值模块确定的阈值
对需要保留的高维连续小波变换系数子带中的小波系数进行阈值处理；
重建模块，对处理模块经阈值处理后的小波变换子带进行高维小波逆变换得到重建的高信噪比信号
压制模块，将分块模块获取的所有分块全部处理完毕，然后将得到的所有高信噪比信号
重新拼接为完整的高信噪比叠前地震数据。
与现有技术相比，本发明至少具有以下有益效果：
本发明一种深层叠前地震数据噪声压制方法，使用连续高维小波变换对深层叠前地震数据资料进行噪声压制处理；首先根据计算机内存大小等因素对数据进行分块以适合处理，然后根据各分块绝对值的中值来估计其噪声标准差，使用一种新的基于快速傅立叶变换的快速连续小波变换算法对原始高维地震数据分块进行高维连续小波变换得到一系列高维连续小波变换系数子带，再通过地震数据中相干噪声的方向确定需要保存的高维连续小波变换系数子带，由估计的噪声标准差确定各个子带的阈值，并对保留的子带进行阈值处理，对经过阈值处理的子带进行小波反变换得到压制噪声后的分块地震数据，全部分块处理完毕后对数据分块进行拼接得到完整压噪数据。相比较于常规压制叠前地震数据噪声的方法，本发明使用高维连续小波的方法，可以同时处理随机噪声和相干噪声，并能取得良好的噪声压制结果。
进一步的，通过对数据进行分块来解决过大数据导致的计算机内存溢出以及运算时间开销过大的问题。通过数据块间设置重叠使得后续处理可以通过裁切的方式解决离散傅里叶变换导致的边缘劣化问题。
进一步的，使用新的快速连续小波变换算法，有利于减少计算的时间复杂度，使高维连续小波变换在地震数据噪声压制中有了实用价值。
进一步的，估计出数据中的噪声标准差
标准差与一个控制滤波强度的常数c及一个正比于
的阈值系数的乘积既是该子带的阈值d(a,θ,φ)。
进一步的，通过去除子带滤除方向性噪声，有利于滤除与反射信号不同方向上的相干噪声。
进一步的，通过对需保留的子带的阈值处理，有利于在保留有效信号的同时去除绝大部分随机噪声，这是由于随机噪声在同一子带内分布较均匀，而有效信号则分布更加集中，表现为同一子带内仅有小部分小波系数绝对值较大，其中包含绝大部分有效信号，而其余小波系数较小的部分则包含绝大部分噪声，因此通过阈值处理即可压制大部分随机噪声而保留有效信号。
进一步的，在步骤S7中使用快速连续小波反变换，有利于加速连续小波反变换过程，有利于进一步提高处理效率。
进一步的，通过对处理过后的数据块进行边缘裁切以避免傅里叶变换导致的边缘劣化现象，且由于步骤S1中分块时已经留有重叠量，保证裁切不会导致数据丢失。再通过数据拼接得到完整的经过噪声压制的深层叠前地震数据。
进一步的，实际工区中采集的数据，其方向性噪声与地层数据分布往往相似或者有一定的规律，有利于使用相同的参数对整个工区的数据进行处理，有利于处理效率的提升。[xk1]
综上所述，本发明可有效、快速的实现对深层叠前地震资料的噪声压制，所采用的方法为快速高维连续小波变换、方向子带选择及硬限幅阈值处理，具有能较好处理方向性噪声的特点，同时也可以去除随机噪声。
下面通过附图和实施例，对本发明的技术方案做进一步的详细描述。
附图说明
图1为本发明流程图；
图2为某实际CMP道集信号图；
图3为图2中的信号经过本专利方法压制噪声的结果；
图4为图2与图3的差剖面；
图5为某实际偏移前单次覆盖体inline线图；
图6为图5信号经过本发明方法压制噪声的结果图；
图7为图5与图6的差剖面。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
应当理解，当在本说明书和所附权利要求书中使用时，术语“包括”和“包含”指示所描述特征、整体、步骤、操作、元素和/或组件的存在，但并不排除一个或多个其它特征、整体、步骤、操作、元素、组件和/或其集合的存在或添加。
还应当理解，在本发明说明书中所使用的术语仅仅是出于描述特定实施例的目的而并不意在限制本发明。如在本发明说明书和所附权利要求书中所使用的那样，除非上下文清楚地指明其它情况，否则单数形式的“一”、“一个”及“该”意在包括复数形式。
还应当进一步理解，在本发明说明书和所附权利要求书中使用的术语“和/或”是指相关联列出的项中的一个或多个的任何组合以及所有可能组合，并且包括这些组合。
在附图中示出了根据本发明公开实施例的各种结构示意图。这些图并非是按比例绘制的，其中为了清楚表达的目的，放大了某些细节，并且可能省略了某些细节。图中所示出的各种区域、层的形状及它们之间的相对大小、位置关系仅是示例性的，实际中可能由于制造公差或技术限制而有所偏差，并且本领域技术人员根据实际所需可以另外设计具有不同形状、大小、相对位置的区域/层。
本发明提供了一种深层叠前地震数据噪声压制方法，使用高维连续小波变换对深层叠前地震数据资料进行噪声压制处理。本发明首先根据计算机内存大小等因素对数据进行分块以适合处理，然后根据各分块绝对值的中值来估计其噪声标准差，使用一种新的基于快速傅立叶变换的快速连续小波变换算法对原始高维地震数据分块进行高维连续小波变换得到系数子带，再通过地震数据中相干噪声的方向确定需要保存的子带，由估计的噪声标准差确定阈值，并对保留的子带进行阈值处理，对经过阈值处理的子带进行小波反变换得到压制噪声后的分块地震数据，全部分块处理完毕后对数据分块进行拼接得到完整压噪数据。
请参阅图1，本发明一种深层叠前地震数据噪声压制方法，包括以下步骤：
S1、对原始叠前数据进行分块；
以一个三维叠前数据体为例，其规格为Nx×Ny×Nz。根据计算机内存及并行运算能力选取适当的分块大小，将原始叠前地震数据分为多个N×N×N的数据块。为解决由于后续不同数据块阈值选择以及数值计算误差导致的分块边界不平滑问题，数据块间需要有部分重叠，以保证在处理完成后拼接时可以通过相应的处理解决分块边界不连续性。
若x方向、y方向及z方向的重叠量均为n点，此时分块数将为
其中符号
表示向上取整。
S2、对步骤S1分块后的其中一块原始叠前地震数据利用中值估计算子估算该块数据的噪声方差等级；
对任意选取的一个未完成处理的分块f(x)，对其做三维傅里叶变换得到的结果为F(k)。选取其高频范围[kmin,kmax]，在此范围之内求得傅里叶变换系数绝对值|F(k)|的中值
其中，MEDIAN是取|F(k)|的中值估计算子。则数据中估计的噪声标准差等级
与|F(k)|的中值
有关：
S3、对步骤S2选定的一块原始叠前地震数据进行高维连续小波变换得到一系列高维连续小波变换的系数子带；
以步骤S2中的未完成处理的分块f(x)为例，其快速连续小波变换结果CWT(f；b,a,θ,φ)为：
其中，参数b表示小波原子的偏移向量(三维数据体下是一个三维向量)，ψ是提前选择的一个具有方向选择性的小波函数，a表示小波原子缩放的尺度因子，θ表示小波原子绕倾角旋转角度，φ表示小波原子绕方位角旋转角度。rθ,φ(x-b)表示对向量x-b绕倾角旋转θ，绕方位角旋转φ。一块三维原始叠前地震数据经过连续小波变换变换后的结果CWT(f；b,a,θ,φ)为六维，而计算时间复杂度为O(N6NaNθNφ)，N为数据体分块的棱长，Na为尺度的个数，Nθ为倾角的个数，Nφ是方位角的个数。
上述计算量极为巨大，导致的计算代价极大，因此使用逆傅里叶变换改写(1)式以进行加速：
其中，IFFT指逆傅里叶变换，
表示原数据块通过高维傅立叶变换得到的波数域数据，
指在参数组(b,a,θ,φ)下的小波原子通过高维傅立叶变换得到的波数域形式，
为
的共轭。此时其时间复杂度为O(N3NaNθNφ)；通过上述加速方法将使时间消耗下降至原先的
在固定尺度a、倾角θ和方位角φ的情况下，可利用(2)式快速计算连续小波变换，得到的连续小波变换结果被称为一个子带。
S4、依据相干干扰的角度确定步骤S3需要保留的高维连续小波变换系数子带；
若叠前地震数据中相干干扰的倾角为θcoherent_noise，考虑到相干干扰的视速度较小，则只保留连续小波变换变换结果CWT(f；b,a,θ,φ)中角度小于θcoherent_noise的子带，即只保留[-θcoherent_noise+Δθ,θcoherent_noise-Δθ]角度范围内的高维连续小波变换系数子带(其中Δθ为保证角度滤波效果的一个较小角度值，一般取值为5度)。
S5、根据步骤S2中得到的噪声标准差确定步骤S4中各个需要保留子带内的阈值；
以步骤S3式(2)中求出噪声标准差
为基础，确定每个需要保留的子带内的阈值thre(a)：
其中，
是每个子带的阈值系数，其取值与1/a3成正比；c为一控制噪声压制程度的常数，通常情况下取为3。
S6、利用步骤S5中确定的每个需要保留子带的阈值
对需要保留子带中的小波系数进行阈值处理；
在步骤S4确定的需要保留的子带中，根据S5中确定的阈值
进行硬限幅滤波：
S7、对步骤S6阈值处理后的小波变换子带进行高维小波逆变换得到重建的高信噪比信号
完成降噪处理；
对一组经过阈值滤波的小波系数CWTthr(f；b,a,θ,φ)进行快速逆小波变换得到重建的高信噪比信号
其中，Cψ是一个与小波函数相关的常数。
S8、重复步骤S2～S7直到所有的分块儿全部处理完毕，然后将得到的高信噪比信号
重新拼接为完整的高信噪比叠前地震数据。
当所有数据分块都按照上述流程处理完成后，对每个分块数据压制噪声后的数据体拼接起来即可得到噪声压制后的数据；由于傅里叶变换导致边缘劣化，需要裁切掉一部分边缘无效数据；由于不同的数据块由于不同的阈值选取及数值化计算误差等原因，其压制噪声效果会有一定的差异，在拼接时会使数据产生明显的分界线；因此在步骤S1分割时使数据间有部分重叠，此部分重叠需要通过在本步骤拼接过程中裁切掉。
以x方向为例，假设数据块f1与f2为x方向上相邻的数据块，有n个点的重叠量，则数据块f1在x方向去掉后
的点，数据块f2在x方向去掉前
的点，再进行拼接。n值需要通过小批量实验确定。
本发明再一个实施例中，提供一种深层叠前地震数据噪声压制系统，该系统能够用于实现上述深层叠前地震数据噪声压制方法，具体的，该深层叠前地震数据噪声压制系统包括分块模块、变换模块、估算模块、选择模块、阈值模块、处理模块、重建模块以及压制模块。
其中，分块模块，获取原始叠前数据并进行分块；
变换模块，利用中值估计算子估算分块模块分块后的任意一块原始叠前地震数据的噪声方差等级；
估算模块，对变换模块估算的原始叠前地震数据进行高维连续小波变换，得到一系列高维连续小波变换的系数子带；
选择模块，依据相干干扰的角度确定估算模块得到的系数子带中需要保留的高维连续小波变换系数子带；
阈值模块，根据变换模块中得到的噪声方差等级确定选择模块中需要保留的高维连续小波变换系数子带内的阈值
处理模块，利用阈值模块确定的阈值
对需要保留的高维连续小波变换系数子带中的小波系数进行阈值处理；
重建模块，对处理模块经阈值处理后的小波变换子带进行高维小波逆变换得到重建的高信噪比信号
压制模块，将分块模块获取的所有分块全部处理完毕，然后将得到的所有高信噪比信号
重新拼接为完整的高信噪比叠前地震数据。
本发明再一个实施例中，提供了一种终端设备，该终端设备包括处理器以及存储器，所述存储器用于存储计算机程序，所述计算机程序包括程序指令，所述处理器用于执行所述计算机存储介质存储的程序指令。处理器可能是中央处理单元(Central ProcessingUnit，CPU)，还可以是其他通用处理器、数字信号处理器(Digital Signal Processor、DSP)、专用集成电路(Application Specific Integrated Circuit，ASIC)、现成可编程门阵列(Field-Programmable GateArray，FPGA)或者其他可编程逻辑器件、分立门或者晶体管逻辑器件、分立硬件组件等，其是终端的计算核心以及控制核心，其适于实现一条或一条以上指令，具体适于加载并执行一条或一条以上指令从而实现相应方法流程或相应功能；本发明实施例所述的处理器可以用于深层叠前地震数据噪声压制方法的操作，包括：
获取原始叠前数据并进行分块；利用中值估计算子估算分块后的任意一块原始叠前地震数据的噪声方差等级；对估算的原始叠前地震数据进行高维连续小波变换，得到一系列高维连续小波变换的系数子带；依据相干干扰的角度确定一系列高维连续小波变换的系数子带中需要保留的高维连续小波变换系数子带；根据噪声方差等级确定需要保留的高维连续小波变换系数子带内的阈值
利用需要保留的高维连续小波变换系数子带内的阈值
对需要保留子带中的小波系数进行阈值处理；对阈值处理后的小波变换子带进行高维小波逆变换得到重建的高信噪比信号
重复以上步骤直到所有的分块儿全部处理完毕，然后将得到的高信噪比信号
重新拼接为完整的高信噪比叠前地震数据。
本发明再一个实施例中，本发明还提供了一种存储介质，具体为计算机可读存储介质(Memory)，所述计算机可读存储介质是终端设备中的记忆设备，用于存放程序和数据。可以理解的是，此处的计算机可读存储介质既可以包括终端设备中的内置存储介质，当然也可以包括终端设备所支持的扩展存储介质。计算机可读存储介质提供存储空间，该存储空间存储了终端的操作系统。并且，在该存储空间中还存放了适于被处理器加载并执行的一条或一条以上的指令，这些指令可以是一个或一个以上的计算机程序(包括程序代码)。需要说明的是，此处的计算机可读存储介质可以是高速RAM存储器，也可以是非不稳定的存储器(non-volatile memory)，例如至少一个磁盘存储器。
可由处理器加载并执行计算机可读存储介质中存放的一条或一条以上指令，以实现上述实施例中有关深层叠前地震数据噪声压制方法的相应步骤；计算机可读存储介质中的一条或一条以上指令由处理器加载并执行如下步骤：
获取原始叠前数据并进行分块；利用中值估计算子估算分块后的任意一块原始叠前地震数据的噪声方差等级；对估算的原始叠前地震数据进行高维连续小波变换，得到一系列高维连续小波变换的系数子带；依据相干干扰的角度确定一系列高维连续小波变换的系数子带中需要保留的高维连续小波变换系数子带；根据噪声方差等级确定需要保留的高维连续小波变换系数子带内的阈值
利用需要保留的高维连续小波变换系数子带内的阈值
对需要保留子带中的小波系数进行阈值处理；对阈值处理后的小波变换子带进行高维小波逆变换得到重建的高信噪比信号
重复以上步骤直到所有的分块儿全部处理完毕，然后将得到的高信噪比信号
重新拼接为完整的高信噪比叠前地震数据。
为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。通常在此处附图中的描述和所示的本发明实施例的组件可以通过各种不同的配置来布置和设计。因此，以下对在附图中提供的本发明的实施例的详细描述并非旨在限制要求保护的本发明的范围，而是仅仅表示本发明的选定实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
以某实际地震数据资料为例。
请参阅图2，图2为某实际地区的叠前CMP道集数据，其中噪声掩埋了大部分有效信号，噪声包括大量的随机噪声以及大量有方向的结构性噪声。请参阅图3、图4，图3为经过本专利方法压制噪声的数据，图中展示了本方法去除了绝大部分随机噪声与结构性噪声，有效信号得以显露，说明本方法对深层叠前地震数据的噪声压制效率较高，图4为图3与图2的差剖面，图中没有肉眼看见的有效数据结构，说明本方法不会损伤有效信号。请参阅图5，为某实际偏移前单次覆盖体inline线，噪声包括大量随机噪声及少量结构性噪声。请参阅图6、图7，图6为本发明方法压制噪声后的数据，绝大部分噪声被去除，可以看清原本被掩埋的信号结构。图7为图5与图6的差剖面，其中不包含肉眼可观察到的有效信号结构。
综上所述，本发明一种深层叠前地震数据噪声压制方法及系统，对深层叠前地震数据中的随机噪声与带有方向的结构性噪声的压制效果都较好，且不会对有效信号造成损伤。相比现有方法，本发明方法能兼顾对随机噪声和结构性噪声的去除。
本领域内的技术人员应明白，本申请的实施例可提供为方法、系统、或计算机程序产品。因此，本申请可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且，本申请可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。
本申请是参照根据本申请实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。
这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。
这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。
以上内容仅为说明本发明的技术思想，不能以此限定本发明的保护范围，凡是按照本发明提出的技术思想，在技术方案基础上所做的任何改动，均落入本发明权利要求书的保护范围之内。
Claims (9)
1.一种深层叠前地震数据噪声压制方法，其特征在于，包括以下步骤：
S1、获取原始叠前数据并进行分块；
S2、利用中值估计算子估算步骤S1分块后的任意一块原始叠前地震数据的噪声标准差等级；
S3、对步骤S2估算的原始叠前地震数据进行高维连续小波变换，得到一系列高维连续小波变换的系数子带；
S4、依据相干干扰的角度确定步骤S3得到的系数子带中需要保留的高维连续小波变换系数子带，通过相干干扰角度确定数据体中的噪声倾角θcoherent_noise；保留连续小波变换变换结果CWT(f；b,a,θ,φ)中角度θ小于θcoherent_noise的子带，即只保留[-θcoherent_noise+Δθ,θcoherent_noise-Δθ]角度范围内的子带；
S5、根据步骤S2中得到的噪声方差等级确定步骤S4中需要保留的高维连续小波变换系数子带内的阈值
S6、利用步骤S5确定的阈值
对需要保留的高维连续小波变换系数子带中的小波系数进行阈值处理；
S7、对步骤S6阈值处理后的小波变换子带进行高维小波逆变换得到重建的高信噪比信号
S8、重复步骤S2～S7直到步骤S1获取的所有分块全部处理完毕，然后将得到的所有高信噪比信号
重新拼接为完整的高信噪比叠前地震数据。
2.根据权利要求1所述的方法，其特征在于，步骤S1中，将原始叠前地震数据分为多个N×N×N的数据块，且在数据块之间设置有重叠，分块数为：
其中，n为重叠量，
表示向上取整，Nx×Ny×Nz为三维叠前数据体。
3.根据权利要求1所述的方法，其特征在于，步骤S2中，噪声标准差等级
具体为：
其中，MEDIAN为取中值的算子，F(k)为对分块做三维傅里叶变换得到的结果，[kmin,kmax]为F(k)的高频范围。
4.根据权利要求1所述的方法，其特征在于，步骤S3中，对任意选取的一个未完成处理的分块f(x)进行高维连续小波变换得到高频子带CWT(f；b,a,θ,φ)如下：
其中，
表示原数据块通过高维傅立叶变换得到的波数域数据，
为在参数组(a,θ,φ；k)下的小波原子通过高维傅立叶变换得到的波数域，
为
的共轭，IFFT表示逆傅里叶变换。
5.根据权利要求1所述的方法，其特征在于，步骤S5中，各高频子带的阈值
为：
其中，
是每个子带的阈值系数，其取值与1/a3成正比；c为一控制噪声压制程度的常数，通常情况下取为3，
为噪声标准差等级。
6.根据权利要求1所述的方法，其特征在于，步骤S6中，根据步骤S5中确定的各高频子带的阈值
进行硬限幅滤波，具体为：
其中，CWTthr(f；b,a,θ,φ)为阈值滤波后的小波系数，CWT(f；b,a,θ,φ)为未经滤波的小波系数。
7.根据权利要求1所述的方法，其特征在于，步骤S7中，高信噪比信号
其中，CWTthr(f；b,a,θ,φ)为阈值滤波后的小波系数，a为小波原子的尺度因子，φ为小波原子在方位角上的极化方向的角度，θ为小波原子在倾角上计划方向的角度。
8.根据权利要求1所述的方法，其特征在于，步骤S8具体为：将步骤S1中分割出的所有数据块进行步骤S2至步骤S7的操作后，重新按照分割时的排列拼接起得到噪声压制后的数据；相邻的两个数据块在拼接前，分别去掉对应方向上
的点，n为重叠量。
9.一种深层叠前地震数据噪声压制系统，其特征在于，包括：
分块模块，获取原始叠前数据并进行分块；
变换模块，利用中值估计算子估算分块模块分块后的任意一块原始叠前地震数据的噪声方差等级；
估算模块，对变换模块估算的原始叠前地震数据进行高维连续小波变换，得到一系列高维连续小波变换的系数子带；
选择模块，依据相干干扰的角度确定估算模块得到的系数子带中需要保留的高维连续小波变换系数子带，通过相干干扰角度确定数据体中的噪声倾角θcoherent_noise；保留连续小波变换变换结果CWT(f；b,a,θ,φ)中角度θ小于θcoherent_noise的子带，即只保留[-θcoherent_noise+Δθ,θcoherent_noise-Δθ]角度范围内的子带；
阈值模块，根据变换模块中得到的噪声方差等级确定选择模块中需要保留的高维连续小波变换系数子带内的阈值
处理模块，利用阈值模块确定的阈值
对需要保留的高维连续小波变换系数子带中的小波系数进行阈值处理；
重建模块，对处理模块经阈值处理后的小波变换子带进行高维小波逆变换得到重建的高信噪比信号
压制模块，将分块模块获取的所有分块全部处理完毕，然后将得到的所有高信噪比信号
重新拼接为完整的高信噪比叠前地震数据。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110674853.3A CN113466940B (zh) | 2021-06-17 | 2021-06-17 | 一种深层叠前地震数据噪声压制方法及系统 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110674853.3A CN113466940B (zh) | 2021-06-17 | 2021-06-17 | 一种深层叠前地震数据噪声压制方法及系统 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN113466940A CN113466940A (zh) | 2021-10-01 |
| CN113466940B true CN113466940B (zh) | 2023-01-03 |
Family
ID=77870475
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110674853.3A Active CN113466940B (zh) | 2021-06-17 | 2021-06-17 | 一种深层叠前地震数据噪声压制方法及系统 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113466940B (zh) |
Families Citing this family (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN117132715B (zh) * | 2023-10-24 | 2024-02-02 | 之江实验室 | 基于物理驱动噪声鲁棒的飞行时间图像重建方法和装置 |
Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105182418A (zh) * | 2015-09-11 | 2015-12-23 | 合肥工业大学 | 一种基于双树复小波域的地震信号降噪方法及系统 |
Family Cites Families (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2008112462A2 (en) * | 2007-03-09 | 2008-09-18 | Fairfield Industries Incorporated | Geophone noise attenuation and wavefield separation using a multi-dimensional decomposition technique |
| US8352192B2 (en) * | 2008-03-28 | 2013-01-08 | Exxonmobil Upstream Research Comapny | Method for performing constrained polarization filtering |
| US9625593B2 (en) * | 2011-04-26 | 2017-04-18 | Exxonmobil Upstream Research Company | Seismic data processing |
| CN103543469B (zh) * | 2012-07-17 | 2016-12-21 | 中国石油化工股份有限公司 | 一种基于小波变换的小尺度阈值去噪方法 |
| EP3193193B1 (en) * | 2016-01-12 | 2021-09-29 | CGG Services SAS | Ava-compliant enhancement of pre-stack frequency spectrum of seismic data |
| CN106597539B (zh) * | 2016-12-28 | 2019-07-12 | 中国石油化工股份有限公司 | 针对黄土塬地区的曲波域Radon变换噪声压制方法 |
| CN107942376A (zh) * | 2018-01-02 | 2018-04-20 | 郑州云海信息技术有限公司 | 一种基于改进的FastICA算法的地震信号处理方法 |
| CN112633225B (zh) * | 2020-12-31 | 2023-07-18 | 矿冶科技集团有限公司 | 矿用微震信号滤波方法 |
-
2021
- 2021-06-17 CN CN202110674853.3A patent/CN113466940B/zh active Active
Patent Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105182418A (zh) * | 2015-09-11 | 2015-12-23 | 合肥工业大学 | 一种基于双树复小波域的地震信号降噪方法及系统 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN113466940A (zh) | 2021-10-01 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Li et al. | A method for low-frequency noise suppression based on mathematical morphology in microseismic monitoring | |
| Chen et al. | Random noise reduction using a hybrid method based on ensemble empirical mode decomposition | |
| CN113466940B (zh) | 一种深层叠前地震数据噪声压制方法及系统 | |
| US20090122061A1 (en) | Seismic data processing | |
| US10393899B2 (en) | Automatic tracking of faults by slope decomposition | |
| CN111736224B (zh) | 一种压制叠前地震资料线性干扰方法、存储介质及设备 | |
| CN107656312B (zh) | 基于分方位角叠加的河道砂体预测方法及装置 | |
| CN108957553B (zh) | 动校正量递推修正的无拉伸畸变动校正方法及装置 | |
| CN109782344B (zh) | 沉积层序边界识别方法及装置 | |
| Wang et al. | Robust singular value decomposition filtering for low signal-to-noise ratio seismic data | |
| CN111352158B (zh) | 地震信号增强方法及装置 | |
| CN111929726B (zh) | 地震相干数据体处理方法及装置 | |
| CN113743193B (zh) | 一种叠前地震资料线性干扰压制方法及系统 | |
| CN112198547A (zh) | 深层或超深层地震资料处理方法及装置 | |
| CN112835095B (zh) | 地震资料低幅度构造成图方法及装置 | |
| CN110596756B (zh) | 基于自适应混合复扩散模型的沙漠地震勘探噪声压制方法 | |
| CN111781642B (zh) | 地震数据层间多次波衰减方法及装置 | |
| CN111781644B (zh) | 浅中地层地震数据的线性干扰衰减方法及装置 | |
| CN111610559B (zh) | 目的层叠前深度偏移成像方法及装置 | |
| CN111538077A (zh) | 基于倾角约束的叠前深度偏移方法及装置 | |
| CN112764108A (zh) | 一种基于改进经验小波变换的新型地震资料噪声压制算法 | |
| CN115908458B (zh) | 一种深海区域干涉条纹提取方法、装置及存储介质 | |
| CN113050187B (zh) | 滤波方法及装置、计算机设备及计算机可读存储介质 | |
| Hadjadj et al. | VSP wavefield separation using the gray-scale Hough transform: synthetic data | |
| CN113466938B (zh) | 积分法叠前深度偏移方法及装置 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |