CN117111155B - 一种基于集成框架的微地震数据去噪方法 - Google Patents
一种基于集成框架的微地震数据去噪方法 Download PDFInfo
- Publication number
- CN117111155B CN117111155B CN202311384530.6A CN202311384530A CN117111155B CN 117111155 B CN117111155 B CN 117111155B CN 202311384530 A CN202311384530 A CN 202311384530A CN 117111155 B CN117111155 B CN 117111155B
- Authority
- CN
- China
- Prior art keywords
- data
- noise
- signal
- denoising
- harmonic
- 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 59
- 230000000737 periodic effect Effects 0.000 claims abstract description 18
- 230000001629 suppression Effects 0.000 claims abstract description 16
- 238000010183 spectrum analysis Methods 0.000 claims abstract description 13
- 238000004422 calculation algorithm Methods 0.000 claims description 49
- 239000011159 matrix material Substances 0.000 claims description 21
- 238000000354 decomposition reaction Methods 0.000 claims description 19
- 239000013598 vector Substances 0.000 claims description 14
- 230000008569 process Effects 0.000 claims description 11
- 230000000694 effects Effects 0.000 claims description 8
- 230000009466 transformation Effects 0.000 claims description 6
- 230000003044 adaptive effect Effects 0.000 claims description 5
- 238000004364 calculation method Methods 0.000 claims description 5
- 238000012935 Averaging Methods 0.000 claims description 4
- 230000010355 oscillation Effects 0.000 claims description 4
- 238000005070 sampling Methods 0.000 claims description 4
- 241000282414 Homo sapiens Species 0.000 claims description 3
- 239000000654 additive Substances 0.000 claims description 3
- 230000000996 additive effect Effects 0.000 claims description 3
- 238000006243 chemical reaction Methods 0.000 claims description 3
- 238000013507 mapping Methods 0.000 claims description 3
- 238000012360 testing method Methods 0.000 description 16
- 238000010586 diagram Methods 0.000 description 14
- 230000008901 benefit Effects 0.000 description 5
- 238000012545 processing Methods 0.000 description 4
- 238000013135 deep learning Methods 0.000 description 2
- 238000011161 development Methods 0.000 description 2
- 230000018109 developmental process Effects 0.000 description 2
- 238000001914 filtration Methods 0.000 description 2
- 230000036039 immunity Effects 0.000 description 2
- 230000009467 reduction Effects 0.000 description 2
- 238000004458 analytical method Methods 0.000 description 1
- 230000002146 bilateral effect Effects 0.000 description 1
- 230000015572 biosynthetic process Effects 0.000 description 1
- 229960001948 caffeine Drugs 0.000 description 1
- 239000002131 composite material Substances 0.000 description 1
- 238000013527 convolutional neural network Methods 0.000 description 1
- 230000003631 expected effect Effects 0.000 description 1
- 239000000284 extract Substances 0.000 description 1
- 238000005755 formation reaction Methods 0.000 description 1
- 238000009499 grossing Methods 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 238000012544 monitoring process Methods 0.000 description 1
- 238000007781 pre-processing Methods 0.000 description 1
- 229910052704 radon Inorganic materials 0.000 description 1
- SYUHGPGVQRZVTB-UHFFFAOYSA-N radon atom Chemical compound [Rn] SYUHGPGVQRZVTB-UHFFFAOYSA-N 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 230000011218 segmentation Effects 0.000 description 1
- 238000000926 separation method Methods 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 239000007787 solid Substances 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- RYYVLZVUVIJVGH-UHFFFAOYSA-N trimethylxanthine Natural products CN1C(=O)N(C)C(=O)C2=C1N=CN2C RYYVLZVUVIJVGH-UHFFFAOYSA-N 0.000 description 1
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
- G01V1/362—Effecting static or dynamic corrections; Stacking
-
- 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/282—Application of seismic models, synthetic seismograms
-
- 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
- G01V1/325—Transforming 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
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/14—Fourier, Walsh or analogous domain transformations, e.g. Laplace, Hilbert, Karhunen-Loeve, transforms
- G06F17/148—Wavelet transforms
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/10—Pre-processing; Data cleansing
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
- G06F18/213—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods
- G06F18/2135—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods based on approximation criteria, e.g. principal component analysis
-
- 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
- G01V2210/32—Noise reduction
-
- 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
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
- G06F2218/02—Preprocessing
- G06F2218/04—Denoising
- G06F2218/06—Denoising by applying a scale-space analysis, e.g. using wavelet analysis
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
- G06F2218/08—Feature extraction
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/30—Assessment of water resources
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Remote Sensing (AREA)
- Theoretical Computer Science (AREA)
- Data Mining & Analysis (AREA)
- General Life Sciences & Earth Sciences (AREA)
- Geology (AREA)
- Geophysics (AREA)
- Acoustics & Sound (AREA)
- Environmental & Geological Engineering (AREA)
- Computer Vision & Pattern Recognition (AREA)
- General Engineering & Computer Science (AREA)
- Mathematical Physics (AREA)
- Evolutionary Computation (AREA)
- Bioinformatics & Computational Biology (AREA)
- Evolutionary Biology (AREA)
- Artificial Intelligence (AREA)
- Computational Mathematics (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Pure & Applied Mathematics (AREA)
- Algebra (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
本发明公开了一种基于集成框架的微地震数据去噪方法，包括，基于周期性的谐波噪声，将有效信号视为非周期，通过奇异谱分析，提取数据的周期状谐波信号，获得去除谐波的信号；基于去除谐波的信号，根据频率的差异，在时频域划分有效信号和随机噪声；然后利用非局部波形相似性，计算具有波形相似结构的信号均值，实现微地震数据中的随机噪声压制。本发明将非局部波形相似性引入微地震去噪中，利用随机噪声的随机性和有效信号的非局部相似性完成数据随机噪声的衰减，同时使用奇异谱分析对谐波噪声进行压制，并与非局部波形相似性结合，可以压制微地震数据的大部分噪声，具有较大实际应用价值。
Description
技术领域
本发明属于微地震数据处理领域，尤其涉及一种基于集成框架的微地震数据去噪方法。
背景技术
在油气田开发阶段，微地震接收地震压裂产生的振动，由此判断地下结构、实时监测压裂情况，及时调整压裂参数，从而降低开发成本。由于震源能量低，传播路径复杂，通常微地震采集得到的数据中包含大量噪声，包括呈随机分布的蓝色噪声、谐波噪声等。国内外大量学者对随机噪声的压制方法进行大量研究，基本包括三大类：单道信号的噪声压制方法、具有多道结构特征的随机噪声压制算法和深度学习类算法。单道信号的噪声压制方法（Huang et al.， 2020）需要考虑信号自身的特征（如频率、振幅等），如，基于滤波器的双边滤波去噪（Zhang and Guo， 2017）、局部投影去噪（Chen et al.， 2018）、基于递归滤波器的随机噪声压制算法（Mousavi and Langston， 2016）等，也包括利用自适应信号频率分解算法进行噪声衰减，包括：EMD（Empirical Mode Decomposition， Gómez and Velis，2016; Li et al.， 2020）、EEMD（Ensemble Empirical Mode Decomposition， Han andBaan， 2015; Jia et al.， 2016; Zuo et al.， 2019）和VMD（Variational modedecomposition， Chen et al.， 2021）等以及时频域去噪算法（Vera Rodriguez et al.，2012;Mousavi et al.， 2016; Akram et al.， 2016）。具有多道结构特征的随机噪声压制算法考虑数据的横向连续性，选择合适的稀疏变换，在变换域中完成噪声衰减，包括Radon变换（Sabbione et al.， 2015）、 变换（Forghani-Arani et al.， 2013）、Shearlet变换（Liang et al.， 2018）等。深度学习类包括：自编码器卷积神经网络（Dong et al.，2019）、自适应稀疏字典学习及去噪（Wang et al.， 2020）等。总的来说，这些算法基本针对数据的随机噪声，并且大部分根据数据的可预测性和频率上的差异完成降噪。
而在图像处理领域，以非局部块相似性为正则化算法，可完成数据的去噪、修补、去模糊等操作，由于叠后地震数据的非局部块相似性特征，众多学者将BM3D (Dabov etal. 2007)、NLM (Non-local mean， Buades et al. 2008)算法引入地震数据中，同压缩感知（Donoho， 2006）结合，实现地震数据的去噪（Buades et al.， 2005; Wang et al.，2019）、重建（Chen et al.， 2021）、高分辨率（Chen et al.， 2022）等工作。地震信号作为反射系数的低通滤波，其波形存在大量的相似，这也是反演多解性的原因。因此地下地层结构不能完全相似，就像不可能存在两片完全相同的树叶。因此，基于块匹配的非局部相似性或多或少会存在一些问题。除随机噪声外，微地震中还可能含有振幅不稳定的规则谐波噪声，目前少有文献对微地震这种噪声的产生机理深入研究。
发明内容
为解决上述问题，本发明提供了如下方案：一种基于集成框架的微地震数据去噪方法，包括：
基于周期性的谐波噪声，将有效信号视为非周期，通过奇异谱分析，提取数据的周期状谐波信号，获得去除谐波的信号；
基于所述去除谐波的信号，根据频率的差异，在时频域划分有效信号和随机噪声；然后利用非局部波形相似性，计算具有波形相似结构的信号均值，实现微地震数据中的随机噪声压制。
优选地，通过奇异谱分析，提取数据的周期状谐波信号，获得去除谐波的信号的过程包括，基于构造在时间序列上的特定矩阵进行奇异值分解，从所述时间序列中分解获得趋势、振荡分量和周期规则的谐波噪声。
优选地，基于构造在时间序列上的特定矩阵进行奇异值分解，从所述时间序列中分解获得趋势、振荡分量和周期规则的谐波噪声的过程包括，
对于长度为的信号/>，将所述信号/>映射为长度为/>的向量序列，获得个长度为/>的向量；
基于所述向量组建Hankel矩阵，对所述Hankel矩阵进行奇异值分解，然后进行分组操作，根据下标集合划分获得互不相交的子集；其中，所述子集代表不同的趋势成分；
根据信号特征选择部分子集参与重建，进行对角线平均，得到分解后的序列。
优选地，所述向量的表达式为：
所述Hankel矩阵的表达式为：
。
优选地，对所述Hankel矩阵进行奇异值分解，然后进行分组操作的过程包括，
令，其中/>为矩阵/>的特征值，且/>，而是矩阵/>对应于特征值的标准正交向量；
令，将轨迹矩阵/>的SVD写成：
；其中，/>；
然后进行分组操作，将下标集合划分成/>个互不相交的子集，每个子集代表不同的趋势成分。
优选地，所述根据信号特征选择部分子集参与重建，进行对角线平均，得到分解后的序列的过程包括，
令，则对应于/>的合成矩阵/>，则；
将每个矩阵变换为一个长度为/>的新序列，得到分解后的序列；
令的元素为/>，令，/>，如果/>，则，否则/>，利用对角线平均公式，将矩阵/>转换为序列/>。
优选地，所述对角线平均公式的表达式为：
。
优选地，基于所述去除谐波的信号，根据频率的差异，在时频域划分有效信号和随机噪声的过程包括，
使用所有尺度计算小波变换系数的叠加振幅，利用小波变换将数据转换到时频域后，所述数据在低频时间分辨率低，频率分辨率高，高频时间分辨率高，频率分辨率低；
基于在数据时频系数中的部分频带内，微地震数据的随机噪声频带较低，且分布连续，能量均衡；有效信号呈间断式能量分布，分布在较高频带的特性将高频从微地震事件的频率范围中分离出来；
所述叠加振幅的计算公式为：
其中，是刻度数，/>表示尺度/>和时间采样点/>处的小波变换系数，/>是总时间点数。
优选地，利用非局部波形相似性，计算具有波形相似结构的信号均值的过程包括，
根据信号的非局部均值算法，单道含加性噪声的微地震数据可以建模为无噪数据/>和噪声/>的和；公式表达式为：
对于给定的位置，通过部分信号加权相加的方式计算获得估计值/>，定义这部分信号的定义域/>，则/>表达式为，
其中，表示归一化系数，使用波形相似性大小计算权重，描述为，
式中，为人为设置的带宽参数，/>表示/>周围局部波形集合，包含/>个样本；
在纯噪声位置计算自适应参数，即，
其中，为纯噪声信号的瞬时振幅在/>和，/>表示虚数单位，/>和均表示常数，用于将与信噪比相关的参数/>放缩到指定范围。
优选地，所述微地震数据去噪方法还包括基于相对误差和余弦相似度对去噪能力进行检测；
其中，所述相对误差表示去噪后数据和干净数据的差异，相对误差越小，去噪效果越好；所述相对误差定义为：
所述余弦相似度表示去噪后数据和干净数据的相似程度，相似度越接近于100%，表示去噪后数据和干净数据越相似，定义为：
。
与现有技术相比，本发明具有如下优点和技术效果：
本发明将奇异谱分析方法与非局部波形相似性原理结合，这两种方法的联合使用可以压制微地震数据的大部分噪声，具有较大实际应用价值。
本发明利用奇异谱分析算法压制数据中的谐波噪声，基于非局部波形相似性提取有效信号。通过奇异谱分析算法考虑了数据中的谐波影响，并且充分利用有效信号和随机噪声的频率差异及两者的波形相似性差异完成降噪。
附图说明
构成本申请的一部分的附图用来提供对本申请的进一步理解，本申请的示意性实施例及其说明用于解释本申请，并不构成对本申请的不当限定。在附图中：
图1为本发明实施例的基于集成框架的微地震数据去噪处理流程图；
图2为本发明实施例的时频域分割数据示意图；
图3为本发明实施例的非局部波形相似性的计算原理图；
图4为本发明实施例的合成地震记录及加噪后数据图；
图5为本发明实施例的算法去噪后数据及误差图；
其中，（a）为算法去噪后的数据图；
（b）为去噪后数据和无噪数据的差图；
图6为本发明实施例的去噪前后时间域对比图；、
其中，（a）表示TF-Thre方法去噪前的时间域对比图；
（b）表示TF-Thre方法去噪后的时间域对比图；
（c）表示SSA-TF-Thre方法去噪前的时间域对比图；
（d）表示SSA-TF-Thre方法去噪后的时间域对比图；
（e）表示VMD方法去噪前的时间域对比图；
（f）表示VMD方法去噪后的时间域对比图；
（g）表示SSA-NLWM方法去噪前的时间域对比图；
（h）表示SSA-NLWM方法去噪后的时间域对比图；
图7为本发明实施例的算法鲁棒性测试图；
图8为本发明实施例的算法抗噪性测试图；
图9为本发明实施例的TF-Thre算法和NLWM算法去噪结果图；
其中，（a）表示TF-Thre算法去噪结果图；
（b）表示NLWM算法去噪结果图；
图10为本发明实施例的含噪多道压裂数据图；
图11为本发明实施例的微震数据去噪后剖面图；
其中，（a）表示采用TF-Thre算法的微震数据去噪后剖面图；
（b）表示采用VMD算法的微震数据去噪后剖面图；
（c）表示采用SSA-NLWM算法的微震数据去噪后剖面图。
具体实施方式
需要说明的是，在不冲突的情况下，本申请中的实施例及实施例中的特征可以相互组合。下面将参考附图并结合实施例来详细说明本申请。
需要说明的是，在附图的流程图示出的步骤可以在诸如一组计算机可执行指令的计算机系统中执行，并且，虽然在流程图中示出了逻辑顺序，但是在某些情况下，可以以不同于此处的顺序执行所示出或描述的步骤。
如图1所示，本发明所提供的一种基于集成框架的微地震数据去噪方法，包括，奇异谱分析（SSA）和非局部波形均值（NLWM）。对于周期性的谐波噪声，将有效信号视为非周期，通过奇异谱分析，提取数据的周期状谐波信号。此外，对于去除谐波的信号，将只包含随机噪声和有效信号；然后根据频率的差异，首先在时频域大致划分信号和随机噪声，利用非局部波形相似性，计算具有波形相似结构的信号均值，达到压制随机噪声的目的。这两种方法的串联可以压制绝大部分微地震数据的噪声。
利用非局部均值算法压制随机噪声会忽略频率域特征。为此，本实施例通过时频域划分方案分离部分噪声和有效信号。最后，由于非局部波形的相似度寻找相似的信号结构，而对周期性较强的谐波无效，本实施例遵循“先规则后随机”的去噪规律，引入奇异谱分析压制谐波。
进一步地，本实施例的奇异谱分析是一种处理非线性时间序列数据的方法，可以对时间序列进行分析和预测。基于构造在时间序列上的特定矩阵的奇异值分解（SVD），从一个时间序列中分解出趋势、振荡分量和噪声。SSA具有非常广泛的适用性，对于时间序列，既不需要假设参数模型，也不需要假设平稳性条件。
具体地，对于长度为的信号/>，将所述信号/>映射为长度为/>的向量序列，获得/>个长度为/>的向量；
这些向量组成Hankel矩阵：
对Hankel矩阵进行奇异值分解，令，其中/>为矩阵/>的特征值，且/>，而/>是矩阵/>对应于特征值的标准正交向量；
令，这种情况下，轨迹矩阵/>的SVD可以写成：
；其中，/>；
然后进行分组操作，将下标集合划分成/>个互不相交的子集，每个子集代表不同的趋势成分，因此可以根据信号特征选择部分子集参与重建。令/>，则对应于/>的合成矩阵/>，则有；
为了重建信号，需要进行对角线平均，将上式的每个矩阵变换为一个长度为的新序列，即得到分解后的序列。令/>的元素为，令/>，，如果/>，则/>，否则，利用对角线平均公式，将矩阵/>转换为序列/>。其中，对角线平均公式的表达式为：
通过上述方法即可实现含谐波微地震数据的去噪。
进一步地，针对时频域信号和噪声的划分，本实施例基于时频域信号噪声分割策略进行信号的预处理，算法假设：信号的能量集中在时频系数的高频部分，而噪声则由值较小，且频率较低的系数表示，且分布面积广。算法根据时频系数的峰值分布，大致区分有效信号和随机噪声的分布频率范围。
利用小波变换（Continuous wavelet transformation，CWT）将数据转换到时频域后，其低频时间分辨率低，频率分辨率高；高频时间分辨率高，频率分辨率低。在数据时频系数中的部分频带内，微震数据的随机噪声频带较低（较高尺度），且分布连续，能量均衡；而有效信号呈间断式能量分布，分布在较高频带（较低尺度）。利用这一特性将高频从微地震事件的频率范围中分离出来。为此，使用所有尺度计算CWT系数的叠加振幅：
其中，是刻度数，/>表示尺度/>和时间采样点/>处的小波变换系数，/>是总时间点数。如图2所示，微地震数据的时频系数上可分为4个区域，区域1是低频随机噪声，本部分的算法主要针对该噪声进行处理；区域2是和有效信号同频带的随机噪声；区域3是有效信号；区域4是较规则的谐波噪声。高低频的差异使CF具有两个不同的峰值。如果从CF的双峰中发现了能量较强的低频，则通过找到一个最佳阈值尺度因子/>的方法来分离波形的低频和高频成分。/>可以将小波系数分离为两个部分，上段包含低频分量的系数，下段包含高频分量的系数，从而完成随机噪声的预分离处理。
进一步地，根据信号的非局部均值算法，单道含加性噪声的微地震数据可以建模为无噪数据/>和噪声/>的和；公式表达式为：
对于给定的位置，通过部分信号加权相加的方式计算获得估计值/>，定义这部分信号的定义域/>，则/>表达式为，
其中，表示归一化系数，使用波形相似性大小计算权重，类似地震反演中的波形指示建模，波形间的相似性越高，认为权重越大，贡献越多。具体描述为，
式中，为人为设置的带宽参数，/>表示/>周围局部波形集合，包含/>个样本；值越大，表示高斯类型权重的方差越大，即该值越大，平均能力越强，去噪能力越高，但也会因此损失部分有效信号。而因为微地震数据通常具有多道数据，且各道信噪比存在差异，且在无微震有效信号位置基本均被噪声覆盖，因此，在纯噪声位置计算自适应参数/>，即，
其中，为纯噪声信号的瞬时振幅在/>和，/>表示虚数单位，/>和均表示常数，用于将与信噪比相关的参数/>放缩到指定范围。
本实施例的非局部波形均值的权重取决于两段波形的相似度，而不是/>和/>的物理距离。具体的相似度计算方式如图3所示，图中含噪信号具备较高相似度，平均后可以获得几个信号的公共部分，即完成了去噪。且这种平均是非局部的，不限制计算成本的情况下，可以完成信号的任意两个波形的匹配。在实际应用中，可以设定相似波形的搜索范围/>，以降低计算量。
测试一：
进一步地，本实施例的微地震数据去噪方法还包括基于相对误差和余弦相似度对去噪能力进行检测；
其中，所述相对误差表示去噪后数据和干净数据的差异，相对误差越小，去噪效果越好；所述相对误差定义为：
所述余弦相似度表示去噪后数据和干净数据的相似程度，相似度越接近于100%，表示去噪后数据和干净数据越相似，定义为：
首先使用一个简单模型测试非局部波形均值（NLWM）的去噪能力，基于随机稀疏反射系数同Ricker子波褶积合成地震记录如图4所示，其中，线条①为添加全频带随机噪声后的数据，线条②为干净数据，在去噪测试中，选择的对比算法为：小波阈值去噪（WTD）、经验小波变换（EWT）和变分模态分解（VMD）算法。
测试结果如图5所示，可以发现NLWM和VMD的去噪结果类似，去噪后信号和无噪信号的相似度差距较小（图5（a），NLWM和VMD的相似度分别为97%和96%），且去噪后信号和无噪信号的差剖面上也证明了这个结论（图5a，NLWM和VMD的相对误差分别为23%和27%）。总体上，NLWM和VMD去噪效果优于WTD和EWT算法。但这种数值上的微弱优势并不能证明NLWM在随机噪声压制方面的优势。但是需注意，在如箭头所示的位置，无波形显示，NLWM在这部分的去噪能力较强，而微震信号绝大部分为稀疏分布，这使得NLWM在微震数据的随机噪声压制中具有得天独厚的优势。
进一步地，使用稀疏模型继续测试算法的去噪能力，本模型测试部分通过合成采样间隔为1ms的地震数据，有效信号使用主频50Hz的Ricker子波和随机反射系数褶积形成（图6（a）、6（c）、6（e）和6（g）黑色实线），同时加入正弦波噪声和蓝色噪声模拟微震数据中的谐波和随机噪声。合成模型如图6（b）、6（d）、6（f）和6（h）黑色实线所示。各算法的临时简称及算法来源如表1所示。
表1
对比方案选择时频域随机噪声阈值压制算法（下称TF-Thre）、时频域随机噪声压制算法同SSA的结合（下称SSA-TF-Thre），同时具备压制谐波和随机噪声能力的变分模态分解算法（VMD）和我们建议的SSA同非局部波形均值的结合（SSA-NLWM），上述四种方法的预期效果及算法来源见表1。去噪结果如图6所示，建议的SSA-NLWM表现出良好的去噪性能，其它三种方法存在不同程度的残余噪声（图6（g）和6（h））；同时，比较图6（a）和图6（c）可以发现，SSA具有相当可靠压制周期噪声能力，具体的，其算法性能在Sacchi（2009）中有详细说明，这体现了算法在本模型的适用性。
为了测试算法的鲁棒性和抗噪性，本实施例取消了随机数种子的设置，这使得每次执行算法生成的噪声、随机反射系数各不相同，基于该模型，使用TF-Thre、SSA-TF-Thre和SSA-NLWM进行算法鲁棒性测试（VMD依赖分解的参数，需要人工选择固有模态，未参与大量测试）。共测试1000次，相对误差曲线（高斯平滑后）如图7所示，建议的SSA-NLWM误差一直最小（线条①），同时SSA-TF-Thre（线条②）和SSA-NLWM（线条①）相差“随机噪声”。
然后固定随机数种子，改变信噪比，测试算法的去噪能力，如图8所示，可以发现，随信噪比的降低，几种方法的误差均随之增大，但是建议的SSA-NLWM误差一直最小（线条①）。并且，SSA-TF的稳定性要优于其它两种算法，这是因为算法对谐波无压制作用。
测试二：
基于S. Mostafa Mousavi (2015)文章中应用的开源数据进行测试，如图9所示，建议的NLWM去噪后的赋值和原始信号的匹配度更高，特别对于弱信号（箭头）；并且如图中虚线框所示，在强振幅处，算法也有着准确匹配振幅的优势，去噪效果略优于TF-Thre算法。
另外，使用一个多道的实际压裂监测微震数据进行去噪测试，原始数据如图10所示，数据共26道，其中第1道的信噪比最低，第3-24道有不同程度的谐波噪声干扰，分别应用TF-Thre、VMD和SSA-NLWM进行噪声压制测试，测试结果如图11所示。TF-Thre对谐波能量强的数据处理效果较差，且存在部分残余噪声（如图11（a）第1道）；VMD算法剩余噪声较多，且损伤有效信号（图11（b））。建议的方法虽然也损失部分有效信号，但是仍体现出较高的去噪能力（图11（c））。
以上，仅为本申请较佳的具体实施方式，但本申请的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本申请揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本申请的保护范围之内。因此，本申请的保护范围应该以权利要求的保护范围为准。
Claims (8)
1.一种基于集成框架的微地震数据去噪方法，其特征在于，包括：
基于周期性的谐波噪声，将有效信号视为非周期，通过奇异谱分析，提取数据的周期状谐波信号，获得去除谐波的信号；
基于所述去除谐波的信号，根据频率的差异，在时频域划分有效信号和随机噪声；
基于所述有效信号，利用非局部波形相似性，计算具有波形相似结构的信号均值，实现微地震数据中的随机噪声压制；
基于所述去除谐波的信号，根据频率的差异，在时频域划分有效信号和随机噪声的过程包括，
使用所有尺度计算小波变换系数的叠加振幅，利用小波变换将数据转换到时频域后，所述数据在低频时间分辨率低，频率分辨率高，高频时间分辨率高，频率分辨率低；
基于在数据时频系数中的部分频带内，微地震数据的随机噪声频带较低，且分布连续，能量均衡；有效信号呈间断式能量分布，分布在较高频带的特性将高频从微地震事件的频率范围中分离出来；
所述叠加振幅的计算公式为：
其中，/>是刻度数，/>表示尺度和时间采样点/>处的小波变换系数，/>是总时间点数；
利用非局部波形相似性，计算具有波形相似结构的信号均值的过程包括，
根据信号的非局部均值算法，单道含加性噪声的微地震数据Y建模为无噪数据X和噪声D的和；公式表达式为：
Y=X+D
对于给定的位置，通过部分信号加权相加的方式计算获得估计值/>，定义这部分信号的定义域/>，则/>表达式为，
其中，/>表示归一化系数，使用波形相似性大小计算权重/>，描述为，
式中，/>为人为设置的带宽参数，/>表示/>周围局部波形集合，包含/>个样本；
在纯噪声位置计算自适应参数，即，
其中，/>为纯噪声信号的瞬时振幅在/>和，/>表示虚数单位，/>和/>均表示常数，用于将与信噪比相关的参数/>放缩到指定范围。
2.根据权利要求1所述的基于集成框架的微地震数据去噪方法，其特征在于，通过奇异谱分析，提取数据的周期状谐波信号，获得去除谐波的信号的过程包括，基于构造在时间序列上的Hankel矩阵进行奇异值分解，从所述时间序列中分解获得趋势、振荡分量和周期规则的谐波噪声。
3.根据权利要求2所述的基于集成框架的微地震数据去噪方法，其特征在于，基于构造在时间序列上的Hankel矩阵进行奇异值分解，从所述时间序列中分解获得趋势、振荡分量和周期规则的谐波噪声的过程包括，
对于长度为的信号/>，将所述信号/>映射为长度为/>的向量序列，获得个长度为/>的向量；
基于所述向量组建Hankel矩阵，对所述Hankel矩阵进行奇异值分解，然后进行分组操作，根据下标集合划分获得互不相交的子集；其中，所述子集代表不同的趋势成分；
根据信号特征选择部分子集参与重建，进行对角线平均，得到分解后的序列。
4.根据权利要求3所述的基于集成框架的微地震数据去噪方法，其特征在于，所述向量的表达式为：
所述Hankel矩阵的表达式为：
。
5.根据权利要求3所述的基于集成框架的微地震数据去噪方法，其特征在于，对所述Hankel矩阵进行奇异值分解，然后进行分组操作的过程包括，
令，其中/>为矩阵/>的特征值，且/>，而是矩阵/>对应于特征值的标准正交向量；
令，将轨迹矩阵/>的SVD写成：
；其中，/>；
然后进行分组操作，将下标集合划分成/>个互不相交的子集/>，每个子集代表不同的趋势成分。
6.根据权利要求3所述的基于集成框架的微地震数据去噪方法，其特征在于，所述根据信号特征选择部分子集参与重建，进行对角线平均，得到分解后的序列的过程包括，
令，则对应于/>的合成矩阵/>，则；
将每个矩阵变换为一个长度为/>的新序列，得到分解后的序列；
令的元素为/>，/>，令，如果/>，则，否则/>，利用对角线平均公式，将矩阵/>转换为序列/>。
7.根据权利要求6所述的基于集成框架的微地震数据去噪方法，其特征在于，所述对角线平均公式的表达式为：
。
8.根据权利要求1所述的基于集成框架的微地震数据去噪方法，其特征在于，所述微地震数据去噪方法还包括基于相对误差和余弦相似度对去噪能力进行检测；
其中，所述相对误差表示去噪后数据和干净数据的差异，相对误差越小，去噪效果越好；所述相对误差定义为：
所述余弦相似度表示去噪后数据和干净数据的相似程度，相似度越接近于100%，表示去噪后数据和干净数据越相似，定义为：
。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202311384530.6A CN117111155B (zh) | 2023-10-25 | 2023-10-25 | 一种基于集成框架的微地震数据去噪方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202311384530.6A CN117111155B (zh) | 2023-10-25 | 2023-10-25 | 一种基于集成框架的微地震数据去噪方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN117111155A CN117111155A (zh) | 2023-11-24 |
| CN117111155B true CN117111155B (zh) | 2023-12-26 |
Family
ID=88795191
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202311384530.6A Active CN117111155B (zh) | 2023-10-25 | 2023-10-25 | 一种基于集成框架的微地震数据去噪方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN117111155B (zh) |
Families Citing this family (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN117290640B (zh) * | 2023-11-27 | 2024-01-26 | 天津领语未来智能科技有限公司 | 用于非线性信号处理的奇异谱谐波分解方法 |
Citations (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| FR2845483A1 (fr) * | 2002-10-04 | 2004-04-09 | Geophysique Cie Gle | Procede de reduction du bruit harmonique dans des signaux vibrosismiques. |
| US7260021B1 (en) * | 2003-01-08 | 2007-08-21 | Westerngeco L.L.C. | Method of harmonic noise attenuation in correlated sweep data |
| WO2011144215A2 (fr) * | 2010-05-20 | 2011-11-24 | Entreprise Nationale De Geophysique - Enageo- Filiale Du Groupe Sonatrach | Méthode d'atténuation du bruit harmonique en vibrosismique par filtrage temps-variant avec référence |
| CN106680874A (zh) * | 2016-12-08 | 2017-05-17 | 西安交通大学 | 基于波形形态特征稀疏化建模的谐波噪声压制方法 |
| CN107957566A (zh) * | 2017-11-17 | 2018-04-24 | 吉林大学 | 基于频率选择奇异谱分析的磁共振测深信号提取方法 |
| CN108254789A (zh) * | 2018-01-26 | 2018-07-06 | 东北石油大学 | 一种压制地震数据随机噪音的去噪方法、设备及介质 |
| CN110837121A (zh) * | 2018-08-16 | 2020-02-25 | 中国石油化工股份有限公司 | 地震资料中的随机噪音消除方法及系统 |
| CN114428346A (zh) * | 2020-10-14 | 2022-05-03 | 中国石油化工股份有限公司 | 一种谐波压制方法及系统 |
| CN115877461A (zh) * | 2022-12-23 | 2023-03-31 | 吉林大学 | 一种基于多尺度注意力交互网络的沙漠地震噪声抑制方法 |
-
2023
- 2023-10-25 CN CN202311384530.6A patent/CN117111155B/zh active Active
Patent Citations (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| FR2845483A1 (fr) * | 2002-10-04 | 2004-04-09 | Geophysique Cie Gle | Procede de reduction du bruit harmonique dans des signaux vibrosismiques. |
| US7260021B1 (en) * | 2003-01-08 | 2007-08-21 | Westerngeco L.L.C. | Method of harmonic noise attenuation in correlated sweep data |
| WO2011144215A2 (fr) * | 2010-05-20 | 2011-11-24 | Entreprise Nationale De Geophysique - Enageo- Filiale Du Groupe Sonatrach | Méthode d'atténuation du bruit harmonique en vibrosismique par filtrage temps-variant avec référence |
| CN106680874A (zh) * | 2016-12-08 | 2017-05-17 | 西安交通大学 | 基于波形形态特征稀疏化建模的谐波噪声压制方法 |
| CN107957566A (zh) * | 2017-11-17 | 2018-04-24 | 吉林大学 | 基于频率选择奇异谱分析的磁共振测深信号提取方法 |
| CN108254789A (zh) * | 2018-01-26 | 2018-07-06 | 东北石油大学 | 一种压制地震数据随机噪音的去噪方法、设备及介质 |
| CN110837121A (zh) * | 2018-08-16 | 2020-02-25 | 中国石油化工股份有限公司 | 地震资料中的随机噪音消除方法及系统 |
| CN114428346A (zh) * | 2020-10-14 | 2022-05-03 | 中国石油化工股份有限公司 | 一种谐波压制方法及系统 |
| CN115877461A (zh) * | 2022-12-23 | 2023-03-31 | 吉林大学 | 一种基于多尺度注意力交互网络的沙漠地震噪声抑制方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN117111155A (zh) | 2023-11-24 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Jian et al. | On the denoising method of prestack seismic data in wavelet domain | |
| CN117111155B (zh) | 一种基于集成框架的微地震数据去噪方法 | |
| CN110031899B (zh) | 基于压缩感知的弱信号提取算法 | |
| Zhou et al. | Robust noise attenuation based on nuclear norm minimization and a trace prediction strategy | |
| Xu et al. | Monochromatic noise removal via sparsity-enabled signal decomposition method | |
| Chen et al. | Plane-wave orthogonal polynomial transform for amplitude-preserving noise attenuation | |
| Huang et al. | Erratic noise suppression using iterative structure‐oriented space‐varying median filtering with sparsity constraint | |
| Zheng et al. | The surface wave suppression using the second generation curvelet transform | |
| Lin et al. | A branch construction-based CNN denoiser for desert seismic data | |
| Zhou et al. | A hybrid method for noise suppression using variational mode decomposition and singular spectrum analysis | |
| CN113887398A (zh) | 一种基于变分模态分解和奇异谱分析的gpr信号去噪方法 | |
| CN110780349A (zh) | 一种基于增强块匹配精度的加权核范数最小化算法及沙漠地震中低频噪声抑制方法和应用 | |
| CN109901224A (zh) | 一种地震资料低频信号保护压制噪声方法 | |
| Du et al. | Study on optical fiber gas-holdup meter signal denoising using improved threshold wavelet transform | |
| CN109959964B (zh) | 一种高铁震源地震信号的宽频背景噪声压制方法 | |
| Lin et al. | Research on microseismic denoising method based on CBDNet | |
| Feng | Seismic random noise attenuation using effective and efficient dictionary learning | |
| Jia et al. | Separating ground-roll from land seismic record via convolutional neural network | |
| Pan et al. | Symplectic transformation based variational Bayesian learning and its applications to gear fault diagnosis | |
| CN113109873B (zh) | 一种基于秩残差约束的沙漠地区地震信号噪声抑制方法 | |
| Liu et al. | Adaptive time-reassigned synchrosqueezing transform for seismic random noise suppression | |
| CN113567129A (zh) | 一种列车轴承振动信号基于ceemd的降噪方法 | |
| Lin et al. | Structure-oriented CUR low-rank approximation for random noise attenuation of seismic data | |
| CN112363217A (zh) | 一种地震数据随机噪声压制方法及系统 | |
| Zhou et al. | Coherent noise attenuation by kurtosis-guided adaptive dictionary learning based on variational sparse representation |
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 |