CN109143341A - 基于Hampel范数的降秩滤波方法 - Google Patents
基于Hampel范数的降秩滤波方法 Download PDFInfo
- Publication number
- CN109143341A CN109143341A CN201810946764.8A CN201810946764A CN109143341A CN 109143341 A CN109143341 A CN 109143341A CN 201810946764 A CN201810946764 A CN 201810946764A CN 109143341 A CN109143341 A CN 109143341A
- Authority
- CN
- China
- Prior art keywords
- hampel
- norm
- matrix
- noise
- reduced
- 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.)
- Pending
Links
- 238000000034 method Methods 0.000 title claims abstract description 47
- 238000001914 filtration Methods 0.000 title claims abstract description 32
- 239000011159 matrix material Substances 0.000 claims abstract description 52
- 230000003247 decreasing effect Effects 0.000 claims abstract description 14
- 238000000354 decomposition reaction Methods 0.000 claims description 11
- BTCSSZJGUNDROE-UHFFFAOYSA-N gamma-aminobutyric acid Chemical compound NCCCC(O)=O BTCSSZJGUNDROE-UHFFFAOYSA-N 0.000 claims description 7
- 238000004364 calculation method Methods 0.000 claims description 6
- 108010094028 Prothrombin Proteins 0.000 claims description 5
- AGVAZMGAQJOSFJ-WZHZPDAFSA-M cobalt(2+);[(2r,3s,4r,5s)-5-(5,6-dimethylbenzimidazol-1-yl)-4-hydroxy-2-(hydroxymethyl)oxolan-3-yl] [(2r)-1-[3-[(1r,2r,3r,4z,7s,9z,12s,13s,14z,17s,18s,19r)-2,13,18-tris(2-amino-2-oxoethyl)-7,12,17-tris(3-amino-3-oxopropyl)-3,5,8,8,13,15,18,19-octamethyl-2 Chemical compound [Co+2].N#[C-].[N-]([C@@H]1[C@H](CC(N)=O)[C@@]2(C)CCC(=O)NC[C@@H](C)OP(O)(=O)O[C@H]3[C@H]([C@H](O[C@@H]3CO)N3C4=CC(C)=C(C)C=C4N=C3)O)\C2=C(C)/C([C@H](C\2(C)C)CCC(N)=O)=N/C/2=C\C([C@H]([C@@]/2(CC(N)=O)C)CCC(N)=O)=N\C\2=C(C)/C2=N[C@]1(C)[C@@](C)(CC(N)=O)[C@@H]2CCC(N)=O AGVAZMGAQJOSFJ-WZHZPDAFSA-M 0.000 claims description 5
- 230000008569 process Effects 0.000 description 8
- 230000009466 transformation Effects 0.000 description 8
- 230000009467 reduction Effects 0.000 description 7
- 230000000694 effects Effects 0.000 description 6
- 230000001427 coherent effect Effects 0.000 description 5
- 230000008602 contraction Effects 0.000 description 5
- 238000005457 optimization Methods 0.000 description 4
- 238000013459 approach Methods 0.000 description 3
- 230000008859 change Effects 0.000 description 3
- 238000010586 diagram Methods 0.000 description 3
- 238000005516 engineering process Methods 0.000 description 3
- 230000001360 synchronised effect Effects 0.000 description 3
- 238000012935 Averaging Methods 0.000 description 2
- 230000007547 defect Effects 0.000 description 2
- 238000011156 evaluation Methods 0.000 description 2
- 238000012545 processing Methods 0.000 description 2
- 239000000654 additive Substances 0.000 description 1
- 230000000996 additive effect Effects 0.000 description 1
- 230000003321 amplification Effects 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000015572 biosynthetic process Effects 0.000 description 1
- 238000013016 damping Methods 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 235000013399 edible fruits Nutrition 0.000 description 1
- 239000011521 glass Substances 0.000 description 1
- 238000009499 grossing Methods 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 230000002401 inhibitory effect Effects 0.000 description 1
- 230000005764 inhibitory process Effects 0.000 description 1
- 238000003780 insertion Methods 0.000 description 1
- 230000001788 irregular Effects 0.000 description 1
- 238000012886 linear function Methods 0.000 description 1
- 229940050561 matrix product Drugs 0.000 description 1
- 239000000203 mixture Substances 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 238000003199 nucleic acid amplification method Methods 0.000 description 1
- 238000012887 quadratic function 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
- 238000011084 recovery Methods 0.000 description 1
- 238000010183 spectrum analysis Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 230000002123 temporal effect Effects 0.000 description 1
- 230000017105 transposition Effects 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
-
- 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
- G01V2210/324—Filtering
- G01V2210/3248—Incoherent noise, e.g. white noise
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
本发明公开一种基于Hampel范数的降秩滤波方法，针对现有技术不能有效地在高斯白噪同时又能同时去除规则噪声的问题，本发明的方法包括：S1、采用基于Hampel三段下降函数的稳健范数构建目标函数；S2、采用迭代重加权最小M估计算法求解步骤S1中的目标函数；S3、根据步骤S2的求解结果，得到Hankel矩阵的稳健估计；能够很好地除高斯和非高斯噪声。
Description
技术领域
本发明属于地震数据处理领域，特别涉及一种降秩滤波技术。
背景技术
地震勘探时所采集到的野外地震资料中伴随着大量的噪声，如何从地震资料中分离噪声在地球物理勘探中一直是一个具有挑战性的课题。从野外采集到的地震噪声主要包含了随机和规则的噪声这两类。所谓随机噪声指的是道与道之间的能量没有规律，具有随机性，本质上为有色随机噪声；相干噪声是指地震剖面之间噪声能量具有规律性，那么通过采集和叠加的简单办法来压制它就变得非常困难，强相干噪声常常使得地震资料质量变得非常差。如何能够有效压制相干噪声，对提高地震资料的质量是非常重要的。
由于强能量的规则和随机噪声会对各种资料处理的负面影响是明显的，比如它们会直接影响反褶积的效果，降低速度分析的精度等。对于不规则/相干噪声衰减，人们主要通过常用的变换域方法tau-p变换、拉东变换、f-x域滤波、eigenimage滤波器、基于小波变换分解,radial trace域和curvelet变换进行去除，但是这些方法有时会产生明显的副作用(如滤波后的道集记录上会出现明显的蚯蚓化现象等)。此外这些方法也有较苛刻的应用前提，比如道集上道数不足以及规则噪声的沿空间方向的倾角变化均会降其应用效果。对于随机噪声，去噪方法很多，比如经典结构导向滤波，中值滤波、保边滤波等。但是没有一种方法既能够去除高斯白噪同时又能同时去除规则噪声。所以在实际地震资料处理时，先需要去除高斯白噪声，然后再去除规则噪声，这样做一方面非常费时费力；另外一方面随机噪声的平滑效果会影响规则噪声的去除。
发明内容
为解决上述技术问题，本发明提供一种基于Hampel范数的降秩滤波方法，采用基于Hampel范数的稳健矩阵分解方法，能够有效去除各种类型和各种强度能量的高斯噪声和非高斯噪声。
本发明采用的技术方案为：一种基于Hampel范数的降秩滤波方法，包括：
S1、采用基于Hampel三段下降函数的稳健范数构建目标函数；
S2、采用迭代重加权最小M估计算法求解步骤S1中的目标函数；
S3、根据步骤S2的求解结果，得到Hankel矩阵的稳健估计。
进一步地，步骤S1的目标函数表达式为：ρ(ei)表示基于Hampel三段下降函数；且ρ(ei)表达式如下：
其中，N表示样本数目，ξ、Δ1和Δ2是Hampel三段下降函数的阈值参数，ei表示噪声，|·|表示取绝对值。
更进一步地，根据地震数据噪声先验信息，确定阈值参数ξ、Δ1和Δ2的取值。
进一步地，步骤S2具体为：
S21、根据步骤S1的目标函数，得到基于Hampel范数下的降秩滤波结果：
其中，U表示第一因子矩阵，V表示第二因子矩阵，M表示地震数据矩阵，m表示M的列向量，arg min(·)表示函数取值最小时的自变量取值，上标H表示复共轭；
S22、分别采用IRLM分解算法更新U和V。
更进一步地，采用IRLM分解算法更新U具体包括以下步骤：
A1、初始化第一因子矩阵U0；
A2、根据前一次迭代结果，计算当前迭代误差：
r′i＝VUi-1-M
其中，i表示迭代次数，且i为正整数；
A3、根据步骤A2计算的当前迭代误差计算加权矩阵：
A4、利用共轭梯度算法求解加权最小二乘问题：
Ui＝(VTW′iV)-1VTW′iM，
A5、达到迭代收敛条件或最大迭代次数，则停止迭代输出估计的第一因子矩阵U。
进一步地，采用IRLM分解算法更新V，具体包括以下步骤：
B1、初始化因子矩阵V0；
B2、根据前一次迭代结果，计算当前迭代误差：
ri＝UVi-1-M
其中，i表示迭代次数，且i为正整数；
B3、根据步骤B2计算的当前迭代误差计算加权矩阵：
B4、利用共轭梯度算法求解加权最小二乘问题：
Vi＝(UTWiU)-1UTWiM，
B5、达到迭代收敛条件或最大迭代次数，则停止迭代输出估计的第二因子矩阵V。
进一步地，步骤S3所述Hankel矩阵的稳健估计的表达式为：
本发明的有益效果：本发明的方法首先采用基于Hampel范数的目标函数代替原来截断SVD所采用的L2范数目标函数，然后采用两个低秩的因子矩阵乘积去逼近Hankel矩阵，最后迭代权值M估计法用于近似求解基于Hampel范数的目标函数，从而获得Hankel矩阵稳健估计；可以很好的去除高斯和非高斯噪声。
附图说明
图1为Hankel矩阵低秩放大特性示意图；
图2为本发明实施例提供的方案流程图；
图3为几种常用的鲁棒性范数及其相应的加权函数；
其中，图3(a)L2范数(L2norm)，图3(b)L2范数的加权函数(L2weight)，图3(c)L1范数(L1norm)，图3(d)L1范数的加权函数(L1weight)，图3(i)Biweight范数(Biweightnorm)，图3(j)Biweight范数的加权函数(Biweight weight)，图3(k)Hampel三段下降函数(Hampel norm)，图3(l)Hampel三段下降函数的加权函数(Hampel weight)；
图4为本发明实施例提供的原始叠后地震剖面图一；
图5为本发明实施例提供的降噪后的叠后剖面图一；
图6为本发明实施例提供的经本发明方法剖面一去除的干扰；
图7为本发明实施例提供的原始叠后地震剖面图二；
图8为本发明实施例提供的降噪后的叠后剖面二；
图9为本发明实施例提供的经本发明方法剖面二去除的干扰。
具体实施方式
为便于本领域技术人员理解本发明的技术内容，下面结合附图对本发明内容进一步阐释。
本发明目的在于提出一种方法在去随机噪声的同时又能有效去除规则噪声；主要考虑：(1)将随机噪声视为高斯分布的噪声，将规则噪声看做非高斯分布的噪声；(2)理想的地震数据具有低秩特性，而随机噪声和规则噪声都会破坏地震数据的低秩结构；从而将随机和规则噪声同步去除方法的问题转变为低秩约束下的高斯和非高斯噪声去除问题。
因此，本发明提出基于Hampel范数的降秩方法来求解以上技术问题，具体思路如下：一、既然随机噪声和规则噪声会破坏地震数据的低秩结构，也就是含噪的地震数据的秩是高的，那么如果将含噪的数据进行降秩，则相当于恢复得到理想地震数据；二是采用结构化矩阵方式，比如Hankel矩阵，可以进一步放大地震数据的低秩特征和含噪数据的升秩特性，相当于是秩的放大镜；三、降秩过程中，普通范数只会去除高斯白噪声，如要去除非高斯噪声，则需要引入Hampel范数去除它。
本实施例中简要介绍降秩滤波方法的基本原理，基于降秩滤波的地震噪声衰减和地震数据重建的详细内容可见文献“Oropeza,V.,and M.D.Sacchi,2011,Simultaneousseismic data denoising and reconstruction via multichannel singular spectrumanalysis:Geophysics,76,no.3,V25–V32,doi:10.1190/1.3552706”。降秩滤波具体主要可以分成两个步骤：
(1)基于滑动时窗的地震数据Hankel变换
本实施例主要以降秩滤波的二维(t-x)地震数据中的实现为例进行说明，但是降秩滤波也可以在广泛应用到三维和五维的地震数据中。降秩滤波主要利用地震数据的低秩性，通过Hankel变换可以进一步降低理想地震信号的低秩性，同时升高含噪地震资料的秩。Hankel变换主要包括以下两个分步骤：
步骤一：对地震数据进行傅里叶变换：小窗口内中的地震数据可以通过频率-空间域内平面波的叠加获得：
其中，j＝1,2,...,N是空间轴上的地震道集索引值，ω表示时间频率。假设地震数据由具有不同射线参数Pk的K个线性事件组成。这里，Ak(ω)表示第k平面波的复振幅，Δx表示地震剖面之间的空间间隔。
步骤二：Hankel变换，降秩滤波将单频D(ω)＝(D1(ω)，D2(ω)，...，DN(ω))T嵌入信号来构建如下Hankel矩阵：
其中，表示Hankel算子。
为了便于计算，本实施例选择使Hankel矩阵逼近为方阵，，省略符号ω，并且对所有频率进行分析。对于K平面波的叠加，可得rank(M)＝K。
(2)基于Hankel域的降秩滤波
由于外加噪声D将会增加矩阵M的秩，降秩滤波器通过降低秩，从而达到降噪的目的。降秩滤波可以表示成如下：
其中，是反对角平均算子，是降秩操作符，它通过秩K矩阵逼近M，是Hankel算子。运算符通过平均对角线将Hankel形式转换成向量。对于多维信号则采用块Hankel矩阵和阻尼对角平均算子。降秩操作可以通过截断的SVD，随机SVD或采用Lanczos双对角化方法和高效矩阵向量乘法的快速算法来实现，快速算法采用快速傅里叶变换。
如图2所示，基于上述原理介绍，本发明的技术方案为：一种基于Hampel范数的降秩滤波方法，包括：
S1、采用基于Hampel三段下降函数的稳健范数构建目标函数；具体过程如下：
(1)高斯噪声条件下的地震数据矩阵低秩分解
传统意义上对于地震数据矩阵M的秩K近似逼近过程是假定高斯条件下的低秩分解，具体是通过求解以下基于L2范数的最优化问题得到：
s.t.rank(M)＝K
其中，||·||F为Frobenius范数，为矩阵
公式(1)中最优化问题本质上是基于Frobenius范数的凸优化，只能平滑高斯噪声。由于采用了Frobenius范数，公式(1)的解具有解析表达式，具体由截断奇异值分解(TSVD)给出：
其中，(第一因子矩阵)和(第二因子矩阵)是包含特征向量的矩阵与K个最大特征值σq相关的向量，q＝1,2,...,K，矩阵是特征值为对角元素的对角矩阵。Frobenius范数对噪声平滑本质上对噪声进行均值滤波，对于规则噪声等所谓离群点或者野值是无法进行处理的。
(2)非高斯噪声条件下的地震数据矩阵低秩分解
为了压制非高斯的规则噪声，本发明通过引入稳健范数将公式(1)的问题进行广义化，得到了随机和规则噪声条件下都稳健的降秩滤波器，其具体定义式为：
其中，ρ(·)为随机和规则噪声条件下都鲁棒的范数。
如前所述，实际地震数据资料除了高斯噪声以外，还同时存在非高斯噪声，比如有线性干扰等形成的相干噪声。本发明最核心的问题就转变为如何设计出一个稳健范数，使得同步压制以上两种类型的范数。如公式(1)所示，基于Frobenius范数的SVD分解可以有效去除高斯噪声，那么新的稳健范数应该包含有Frobenius范数的部分，同时应该包含有非高斯降噪部分的内容。
对于同步压制实际地震资料中的相干和非相干噪声的稳健范数，本发明具体设计思路为：①如图2(a)和(b)所示L2范数消除高斯噪声，如图2(c)和(d)所示L1范数由于其加权函数本质上对数据做一个中值滤波，所以能够处理野值和离群值。但是对野值或离群值较大或者较多的情况，L1范数也无法进行消除，这时候一般会引入Biweight范数对大的野值或离群值进行截止；②但是从如图2(i)和(j)可得，Biweight范数缺陷也是非常明显，及时对大的野值或离群值可以截止，对小的高斯噪声可以平滑，但是无法去除中等大小的规则噪声；③基于以上各类范数缺陷，本发明提出了基于Hampel三段下降函数的稳健范数，对小的高斯随机噪声，采用L2范数；对中等大小的规则噪声，采用L1范数；对于大的规则噪声，直接进行截止操作。
本发明提出的基于Hampel函数的稳健范数来同步压制高斯和非高斯噪声，其Hampel三段下降函数具体定义为：
其中，ξ、Δ1和Δ2是Hampel三段下降函数的阈值，在实际地震资料处理过程中，阈值根据信噪比进行取值。在本实施例中，将称作Hampel函数，图2(k)和(l)给出了Hampel函数的具体示意图，从图1和公式(4)可得：
①当噪声ei小于阈值ξ时，Hampel函数为二次函数；
②当噪声ξ≤|ei|＜Δ2时，Hampel函数为一次函数；
③当噪声|ei|≥Δ2时，Hampel函数为一常数。
图1中a slice data表示切片数据，Block Hankel matrix表示Hankel矩阵。
Hampel函数结合了L1范数、L2范数和Biweight范数的优点：当地震噪声较小时，则使用L2范数平滑噪声；当存在中等强度的规则噪声时，则使用L1范数压制噪声；当存在大量高能量的规则和随机噪声时，直接对其进行截断操作，直接忽略了大的规则和随机噪声对降秩滤波器的影响。
阈值参数ξ、Δ1和Δ2用来控制Hampel范数对非高斯规则噪声和高斯噪声的抑制能力，阈值越小，Hampel范数抗规则噪声能力越强。Hampel范数参数估计主要依赖于规则噪声的先验分布概率。但是通常情况下，随机和规则噪声误差信号e的统计分布是未知的，本实施例中为了简便计算，先假设噪声误差信号e服从高斯分布，同时包含有非高斯规则噪声。通过估计没有规则噪声误差信号e的分布，实现检测到并且能够抑制规则噪声的目的。这样，误差大于某一给定阈值T的概率可表示为：
其中，为误差函数，表示不含非高斯噪声的误差信号e的标准差σe的估计。令θξ＝Pr{|e|＞ξ}，分别表示|e|大于ξ、Δ1和Δ2的概率。选择合适的θξ、和就可以确定出阈值参数ξ、Δ1和Δ2。例如，选取ξ、Δ1和Δ2分别为地震数据的95％、97.5％和99％置信区间，可以得到相应的阈值应取为：
接下来便是如何估计常用的方法便是利用递归公式对瞬时方差估计进行滤波：
其中，0＜λδ＜1为遗忘因子，k表示当前迭代次数，σe(k)表示第k次迭代误差ek的瞬时方差；但是σe(k)容易受到脉冲噪声的影响，所以本发明采用另外一种稳健性的估计方法：
其中，mad(·)表示平均绝对中值偏差(Median Absolute Deviation from themedian)，用来抑制非高斯规则噪声对方差估计的影响。
S2、采用迭代重加权最小M估计算法求解步骤S1中的目标函数；具体过程为：
对式(4)定义的目标函数进行优化求解，便可以得到基于Hampel范数下的降秩滤波结果。具体优化方法如下：
本发明选择使用一种新的迭代重加权最小M估计(Iterative re-weighted leastM-estimate,IRLM)算法(详细可参考文献：Y.Liu,J.S.Zhang,.Iterative Re-weightedLeast M-estimate algorithm and its application to AVO inversion[J],Exploration Geophysics)，M估计方法相比其它最小二乘法算法更稳健。根据式(4)，IRLM算法的加权函数可表示为：
其中，sign(·)表示符号函数。
图2(l)给出了Hampel函数的加权函数示意图，从图中可以看出，当噪声大于Δ2时，加权函数等于零，在计算误差的时候就完全排除了大的脉冲噪声的影响，从而具有更强的稳健性。由此，本发明采用的IRLM算法来分解交替求解公式(6)更新因子矩阵U与因子矩阵V；以下以因子矩阵V的更新过程为例进行说明(k表示迭代次数)：
Step1:根据上一次迭代结果(V0表示迭代初始模型)，计算当前时刻迭代误差：
rk＝UVk-1-M，
Step2:根据当前迭代误差计算加权矩阵：
Step3:利用共轭梯度算法求解加权最小二乘问题：
Vk＝(UTWkU)-1UTWkM，
上标T表示矩阵转置运算。
Step4:返回迭代步骤step1，直到满足设定的迭代收敛条件，输出估计结果。本发明中设定的迭代条件可以采用迭代误差小于预设门限，一般当迭代误差小于0.001时，则停止迭代，在具体的应用中预设门限可以根据实际需要修改。
S3、根据步骤S2的求解结果，得到Hankel矩阵的稳健估计。根据步骤S2估计得到的U和V，获得Hankel矩阵的估计 表示矩阵M的估计。
图4至图6为本实施例提供的原始叠后地震剖面图一的处理过程；图4为原始叠后地震剖面图一，图5为降噪后的叠后剖面一，图6为经本发明方法剖面一去除的干扰；图7只图9为本实施例提供的原始叠后地震剖面图二的处理过程；图7为原始叠后地震剖面图二，图8为降噪后的叠后剖面二，图9为经本发明方法剖面二去除的干扰；可以看出，实际去噪结果符合理论假设，并且实际效果非常好，降噪效果很明显。图4至图9中的Rows表示行，Columns表示列。
本领域的普通技术人员将会意识到，这里所述的实施例是为了帮助读者理解本发明的原理，应被理解为本发明的保护范围并不局限于这样的特别陈述和实施例。对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的权利要求范围之内。
Claims (7)
1.一种基于Hampel范数的降秩滤波方法，其特征在于，包括：
S1、采用基于Hampel三段下降函数的稳健范数构建目标函数；
S2、采用迭代重加权最小M估计算法求解步骤S1中的目标函数；
S3、根据步骤S2的求解结果，得到Hankel矩阵的稳健估计。
2.根据权利要求1所述的一种基于Hampel范数的降秩滤波方法，其特征在于，步骤S1的目标函数表达式为：ρ(ei)表示基于Hampel三段下降函数；且ρ(ei)表达式如下：
其中，N表示样本数目，ξ、Δ1和Δ2是Hampel三段下降函数的阈值参数，ei表示噪声，|·|表示取绝对值。
3.根据权利要求2所述的一种基于Hampel范数的降秩滤波方法，其特征在于，根据地震数据噪声先验信息，确定阈值参数ξ、Δ1和Δ2的取值。
4.根据权利要求2所述的一种基于Hampel范数的降秩滤波方法，其特征在于，步骤S2具体为：
S21、根据步骤S1的目标函数，得到基于Hampel范数下的降秩滤波结果：
其中，U表示第一因子矩阵，V表示第二因子矩阵，M表示地震数据矩阵，m表示M的列向量，argmin(·)表示函数取值最小时的自变量取值，上标H表示复共轭；
S22、分别采用IRLM分解算法更新U和V。
5.根据权利要求4所述的一种基于Hampel范数的降秩滤波方法，其特征在于，采用IRLM分解算法更新U具体包括以下步骤：
A1、初始化第一因子矩阵U0；
A2、根据前一次迭代结果，计算当前迭代误差：
r′i＝VUi-1-M
其中，i表示迭代次数，且i为正整数；
A3、根据步骤A2计算的当前迭代误差计算加权矩阵：
A4、利用共轭梯度算法求解加权最小二乘问题：
Ui＝(VTW′iV)-1VTW′iM，
A5、达到迭代收敛条件或最大迭代次数，则停止迭代输出估计的第一因子矩阵U。
6.根据权利要求4所述的一种基于Hampel范数的降秩滤波方法，其特征在于，采用IRLM分解算法更新V，具体包括以下步骤：
B1、初始化因子矩阵V0；
B2、根据前一次迭代结果，计算当前迭代误差：
ri＝UVi-1-M
其中，i表示迭代次数，且i为正整数；
B3、根据步骤B2计算的当前迭代误差计算加权矩阵：
B4、利用共轭梯度算法求解加权最小二乘问题：
Vi＝(UTWiU)-1UTWiM，
B5、达到迭代收敛条件或最大迭代次数，则停止迭代输出估计的第二因子矩阵V。
7.根据权利要求5或6所述的一种基于Hampel范数的降秩滤波方法，其特征在于，步骤S3所述Hankel矩阵的稳健估计的表达式为：
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201810946764.8A CN109143341A (zh) | 2018-08-20 | 2018-08-20 | 基于Hampel范数的降秩滤波方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201810946764.8A CN109143341A (zh) | 2018-08-20 | 2018-08-20 | 基于Hampel范数的降秩滤波方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN109143341A true CN109143341A (zh) | 2019-01-04 |
Family
ID=64790298
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201810946764.8A Pending CN109143341A (zh) | 2018-08-20 | 2018-08-20 | 基于Hampel范数的降秩滤波方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN109143341A (zh) |
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111010145A (zh) * | 2019-12-10 | 2020-04-14 | 西南大学 | 基于范数正则化离散线性系统的滤波方法、离散线性系统 |
| CN111257929A (zh) * | 2020-02-17 | 2020-06-09 | 成都理工大学 | 一种奇异值衰减的降秩去噪方法 |
| CN114966861A (zh) * | 2022-05-17 | 2022-08-30 | 成都理工大学 | 基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法 |
Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2010088190A2 (en) * | 2009-01-30 | 2010-08-05 | Geco Technology B.V. | Deriving a function that represents data points |
| CN106646612A (zh) * | 2016-12-20 | 2017-05-10 | 中国地质大学（北京） | 基于矩阵降秩的地震数据重建方法 |
-
2018
- 2018-08-20 CN CN201810946764.8A patent/CN109143341A/zh active Pending
Patent Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2010088190A2 (en) * | 2009-01-30 | 2010-08-05 | Geco Technology B.V. | Deriving a function that represents data points |
| CN106646612A (zh) * | 2016-12-20 | 2017-05-10 | 中国地质大学（北京） | 基于矩阵降秩的地震数据重建方法 |
Non-Patent Citations (4)
| Title |
|---|
| KE CHEN 等: "Robust reduced-rank seismic denoising", 《SEG TECHNICAL PROGRAM EXPANDED ABSTRACTS 2013》 * |
| 刘洋: "非高斯AVO三参数反演算法及其应用研究", 《中国博士学位论文全文数据库·基础科学辑》 * |
| 李志勇 等: "基于Hampel三截尾函数的储层弹性和物性参数同步反演", 《石油物探》 * |
| 贾超: "稳健回归分析方法在变形监测中的应用", 《中国优秀硕士学位论文全文数据库》 * |
Cited By (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111010145A (zh) * | 2019-12-10 | 2020-04-14 | 西南大学 | 基于范数正则化离散线性系统的滤波方法、离散线性系统 |
| CN111257929A (zh) * | 2020-02-17 | 2020-06-09 | 成都理工大学 | 一种奇异值衰减的降秩去噪方法 |
| CN111257929B (zh) * | 2020-02-17 | 2021-03-19 | 成都理工大学 | 一种奇异值衰减的降秩去噪方法 |
| CN114966861A (zh) * | 2022-05-17 | 2022-08-30 | 成都理工大学 | 基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法 |
| CN114966861B (zh) * | 2022-05-17 | 2024-03-26 | 成都理工大学 | 基于Lp伪范数和γ范数稀疏低秩约束的地震去噪方法 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Xue et al. | Amplitude-preserving iterative deblending of simultaneous source seismic data using high-order radon transform | |
| Huang et al. | Signal extraction using randomized-order multichannel singular spectrum analysis | |
| Liu et al. | A 1D time-varying median filter for seismic random, spike-like noise elimination | |
| CN106646612B (zh) | 基于矩阵降秩的地震数据重建方法 | |
| CN108828670B (zh) | 一种地震数据降噪方法 | |
| CN109143341A (zh) | 基于Hampel范数的降秩滤波方法 | |
| Zhou | A POCS method for iterative deblending constrained by a blending mask | |
| CN107179550B (zh) | 一种数据驱动的地震信号零相位反褶积方法 | |
| CN107144879A (zh) | 一种基于自适应滤波与小波变换结合的地震波降噪方法 | |
| CN103364832A (zh) | 一种基于自适应最优核时频分布的地震衰减定性估计方法 | |
| Wang et al. | Seismic data denoising for complex structure using BM3D and local similarity | |
| CN103364826A (zh) | 基于独立分量分析的地震盲源反褶积方法 | |
| CN113887398A (zh) | 一种基于变分模态分解和奇异谱分析的gpr信号去噪方法 | |
| CN105301648A (zh) | 一种获取共反射面元叠加参数的方法 | |
| CN107621654A (zh) | 一种基于最大相关熵的地震叠后波阻抗反演方法 | |
| Huang et al. | Randomized-order multichannel singular spectrum analysis for simultaneously attenuating random and coherent noise | |
| Jiang et al. | Seismic wavefield information extraction method based on adaptive local singular value decomposition | |
| Wang et al. | Multicomponent seismic noise attenuation with multivariate order statistic filters | |
| US11634175B2 (en) | Dip angle-steering median filtering method based on a niche differential evolution algorithm | |
| Bing et al. | A robust random noise suppression method for seismic data using sparse low-rank estimation in the time-frequency domain | |
| Aharchaou et al. | Prestack Q compensation with sparse tau-p operators | |
| CN107678065B (zh) | 提高地震分辨率的保构造井控空间反褶积方法和装置 | |
| Porsani et al. | An adaptive local-slope SVD filtering approach to enhance events on seismic sections | |
| Han et al. | Seismic signal detection algorithm based on gs transform filtering and emd denoising | |
| CN105954799B (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 | ||
| RJ01 | Rejection of invention patent application after publication | ||
| RJ01 | Rejection of invention patent application after publication |
Application publication date: 20190104 |