基于自适应K-SVD的能量泄漏恢复研究

何真,曹思远,郝婳婕,段绪远,叶力夏提,张韵佩

(1.中国石油新疆油田分公司 勘探开发研究院，新疆 克拉玛依 834000;2.中国石油大学(北京) 地球物理学院，北京 102249;3.中国石油大学(北京) CNPC物探重点实验室，北京 102249;4.中国石油新疆油田分公司 数据公司信息研究所，新疆 克拉玛依 834000;5.中国石油新疆油田分公司 新港公司，新疆 克拉玛依 834000)

0 引言

在地震资料处理中，随机噪声的压制研究一直都是热点与难点，处理结果的好坏直接影响着后续解释与储层预测等工作。目前生产实际中常用的去噪手段有小波变换、奇异值分解(SVD)、F-X域去噪、Radon变换、三维块匹配算法(BM3D)等一系列方法。其中，小波变换[1]作为一种经典的数学变换，对于不同特征的带噪信号采用软阈值函数或硬阈值函数进行降噪处理，虽有一定效果但仍然有待提高，该方法在提出以后也不断有专家学者对其进行讨论与改进[2-3]。奇异值分解算法(SVD)[4]借助地震信号之间的相关性与随机噪声之间的不相关来进行信噪分离，根据能量的不同分为若干奇异值，大奇异值对应有效信号，小奇异值对应噪声，但地震同相轴斜率较大时则不适用该方法。后来刘志鹏、赵伟等[5]提出局部频率域SVD对地震波进行去噪，基于有效频带频谱存在较好的相关性来进行噪声压制。f-x域去噪[6-7]方法假定相干信号在f-x域可预测而随机噪声不可预测的原理对信号和噪声加以区分来达到去噪的目的。Radon变换去噪主要有线性、抛物线和双曲线拉东变换3种去噪方法[8-10]，根据同相轴的不同时差特征来选取合适的变换类型进行去噪。三维块匹配算法[11-12](BM3D)是近些年热度较高的图像去噪算法之一，它将噪声图像沿垂直与水平方向搜索相似块并分组，再进行协同滤波，该方法解决了图像去噪后模糊的问题，保留了图像的边缘细节，在一定程度上优于传统的空间域降噪技术和频率域降噪技术。

尽管去噪方法和理论在不断改进和完善，但在常规处理中由于参数选择的不合理或者基于一些不合理假设导致在对信号进行降噪的同时难免会“伤及”一部分有效信号，即残差剖面中存在有效信号的情形，称其为能量泄漏现象。对于噪声中含有的有效信号部分很少有人进行研究。本文提出基于优质样本的学习字典进行能量恢复研究，优质样本构建的字典包含较多的原地震信息以及较少的随机噪声，并且与残差剖面中残留的有效信号有较强的相关性，可以有效地将其提取并回收利用。

1 自适应K-SVD算法能量泄漏恢复原理

1.1 K-SVD算法

实际地震数据可表示为：

f=s+σ*n,

(1)

其中：f为含噪地震数据；s为有效信号；n为随机噪声；σ为噪声强度，一般不可知。压制随机噪声的过程就是将σ*n这一部分从地震数据f中去除，现用任意手段将信号f进行去噪处理得到：

D_noise=Method(f) ,

(2)

其中：D_noise表示去噪后信号，Method表示任一种去噪方法，可以是小波阈值去噪、SVD分解、中值滤波等。

一次去噪后信号D_noise即本文研究的部分，将其写成如下形式：

D_noise=D*X,

(3)

i=1,2,…,n

(4)

T0为稀疏表示中非零分量个数的最大值。字典D的更新是逐列进行，假设字典D和稀疏矩阵X都是固定的，当更新字典D中第h列原子dh时，则有：

(5)

(6)

(7)

(8)

(9)

Σ=diag(σ1,σ2,σ3,…,σr) 。

(10)

σr(r=1,2,3,…,r)表示矩阵的所有非零奇异值，r为矩阵的秩。

(11)

(12)

1.2 正交匹配追踪(OMP)

OMP算法是根据匹配度选择最优的原子进入原子集，它在更新信号残差时，减去的是当前残差在所有已被选择的原子上的正交投影，这样就使得该算法在迭代过程中避免对已选的原子的重复挑选。并且通过求解一个最小二乘问题找到原信号的最优近似解，再对残差进行更新迭代，最后用原子的线性组合来表征信号。具体算法如下：

1)输入m×n的字典D以及需要稀疏重构的信号y，给定预设误差ε,稀疏度K；

2)初始化残差ε0=y，索引集Λ0=Ø，迭代次数t=1；

3)将字典D每列与残差εt做内积，选择内积最大的一个即信号在原子垂直投影长度中最长的一个，并记录其索引φt；

4)令Λt=Λt-1∪{φt};

7)判断是否满足εt<ε，若满足，停止迭代，否则返回3)。

1.3 相关系数匹配准则

(13)

(14)

此时停止OMP迭代，输出结果。

1.4 本文算法实现流程

图1 本文算法流程图Fig.1 Algorithm flow chart

2 实验结果与分析

实验分为对模型资料的处理以及实际资料的处理，首先通过楔形模型和倾斜模型的处理证明本文方法的正确性和有效性，再对实际资料进行分析、去噪、能量恢复获得高信噪比高保真度的地震剖面。

实验中为统一标准，设置地震数据块大小为8×8，字典中的原子数为256，初始化字典采用DCT字典。在评价去噪效果时采用能量叠加法[15]计算信噪比，公式为：

(15)

(16)

式中：M、N分别为采样点数和地震道数，SNR为信噪比，E为地震道总能量，EN为噪声总能量，ES为有效信号能量。

2.1 楔形模型去噪分析

图2a是楔形模型，道数30道，采样点数500，子波主频35 Hz；图2b是加入了50%的随机噪声，其信噪比为2.18 dB。现在利用KSVD算法分别对无噪和带噪情况下的模型求取特征字典，得到图2c(带噪字典)，图2d(无噪字典)。不难发现在无噪情况下，特征字典提取的都是有效信号的成分，而噪声情况下的字典则较为“模糊”。现对加噪信号分别采用带噪和无噪的字典进行稀疏求解，得到去噪结果图2e和图2f。信噪比分别为18.34 dB和20.05 dB，两种方法对噪声均有较好的压制作用，但对比残差剖面图2g、h和单道15道信号图2i、j、k、l来看，采用无噪字典去噪的剖面在信号的去噪效果和保幅性上更胜一筹。

2.2 倾斜模型能量泄漏恢复效果分析

上文2.1中通过楔形模型证明利用一个优质的学习样本(即无噪声或噪声尽可能少的情况下的剖面)进行字典学习对噪声压制以及信号重构的效果要优于传统对带噪信号进行学习(即普通样本)的效果。接下来采用倾斜模型，对其加入一定比例随机噪声，先利用传统去噪手段对其进行去噪处理，再分别对去噪后剖面和原始带噪剖面作为样本进行字典学习，将学习到的字典对残差剖面进行稀疏分解，实现信号和噪声的二次分离，最后将分离出的有效信号补回去噪后的剖面并分析讨论两种方法的恢复结果。

a—楔形模型;b—加噪楔形模型;c—噪声情况下构建的字典;d—无噪声情况下构建的字典;e—噪声字典去噪剖面;f—无噪声字典去噪剖面;g—噪声字典去噪后残差;h—无噪声字典去噪后残差;i—楔形模型第15道信号;j—加噪楔形模型第15道;k—噪声字典去噪后第15道;l—无噪字典去噪后第15道a—wedge model;b—noise wedge model;c—dictionary constructed under noise;d—noise dictionary denoising profile;e—noise dictionary denoising profile;f—noiseless dictionary denoising profile;g—noise dictionary denoising residuals;h—noiseless dictionary denoising residuals;i—wedge model 15th signal;j—noise wedge model 15th;k—noise dictionary denoising after 15th;l—No15 after noiseless dictionary denoising图2 楔形模型试验结果Fig.2 Wedge model test results

图3a为倾斜模型，道数81道，采样点数801，子波主频55 Hz；图3b加入65%的随机噪声，信噪比为3.31 dB。利用传统小波变换对其进行去噪处理，得到图3c，信噪比为19.62dB，去噪后剖面信噪比有了大幅度的提升，但从残差剖面观察(图3d)，去除的噪声中存在大量有效信号。以普通样本和优质样本为基础进行字典学习对残差剖面进行处理，得到图3e和图3f的信噪分离结果。从分离结果可以看出，以优质样本学习的字典在残差剖面中信噪分离的效果更加彻底，将其分别补回去噪后剖面得到能量恢复结果(图3g、h)，信噪比分别为18.13dB和21.25 dB。可见本文提出的方法对有效信号的提取更加精准和可靠，下文的单道信号分析中更能清晰地体现这一结论。

a—倾斜模型;b—加噪倾斜模型;c—小波阈值去噪剖面;d—小波阈值去噪残差;e—普通样本分离的有效信号;f—优质样本分离的有效信号;g—普通样本能量恢复结果;h—优势样本能量恢复结果a—tilt model;b—charging tilt model;c—wavelet threshold denoising profile;d—wavelet threshold denoising residual;e—effective signal for normal sample separation;f—effective signal for good sample separation;g—ordinary sample energy recovery results;h—high quality sample energy recovery results图3 倾斜模型剖面试验结果Fig.3 Tilt model section test results

对图3倾斜模型时间域第75道，观察普通样本的字典学习恢复效果(图4d)和优质样本的字典学习恢复效果(图4e)，优质样本学习下恢复的结果形态更加接近原始信号，而普通样本学习下恢复的结果中含有较多的随机噪声，信噪分离程度相对不够理想。

a—倾斜模型;b—倾斜模型加噪;c—小波阈值去噪结果;d—普通样本能量恢复结果;e—优质样本能量恢复结果a—tilt model;b—tilt model noise;c—wavelet threshold denoising result;d—normal sample energy recovery result;e—good sample energy recovery result图4 倾斜模型第75道信号试验结果Fig.4 Test results of the 75th signal of tilt model

为更好地体现本文方法的实用性，现将带噪倾斜剖面分别采用f-x域去噪、改进小波阈值[16]、ASVD[17]3种方法对其进行一次去噪，再分别对3个去噪剖面进行字典构造进而提取残差剖面的有效信号，并对结果进行分析。

由于补偿恢复后的结果从时间域剖面中很难看出差异，所以对以上3种恢复结果只从单道信号进行观察讨论。从去噪结果看，3种方法均有较好的效果，噪声基本被去除，但f-x域去噪相对于改进小波阈值和ASVD方法，在去除噪声的同时也伤及了大量有效信号。仍然以恢复后的第75道为例，见图5g、h、i，对f-x域去噪、改进小波阈值以及ASVD三种去噪结果进行能量恢复后随机噪声含量均比较小，但在有效信号振幅能量上后两者恢复得更好，更加接近原始模型资料。经对比也不难发现：原始剖面降噪效果越好，有效信号保留越多，在构建字典原子时信号的特征越能被更好地提取和表达，从而实现更加理想的恢复效果。

a—f-x域去噪剖面;b—f-x域去噪残差;c—改进小波阈值去噪剖面;d—改进小波阈值去噪残差;e—ASVD去噪剖面;f—ASVD去噪残差;g—f-x域去噪后能量恢复结果(第75道);h—改进小波阈值去噪后能量恢复结果(第75道);i—ASVD去噪后能量恢复结果(第75道)a—f-x domain denoising profile;b—f-x domain denoising residual;c—improved wavelet threshold denoising profile;d—improved wavelet threshold denoising residual;e—ASVD denoising proflie;f—ASVD denoising residual;g—f-x domain denoising energy recovery result(75th channel);h—improved wavelet threshold denoising energy recovery result(75th channel);i—ASVD denoising energy recovery result(75th channel)图5 本文方法实用性试验结果Fig.5 Result of the practical test of the method

2.3 实际资料处理及效果分析

选取某实际二维叠加数据进行处理分析。图6a为某区块实际地震资料，剖面上随机噪声较为发育，能量较弱的信号基本被湮没，给接下来的小层解释、储层预测等工作带来了较大难度。实际资料的一次去噪选取fx域和ASVD法，进行f-x域去噪，得到去噪结果图6b和残差剖面图6c，采用本文方法对残差剖面处理实现噪声信号二次分离，将分离的有效信号加回去噪后剖面得到能量恢复结果(图6d)，同样地，采用上述方法得到ASVD的处理结果(图6e、f、g)，对于一次去噪结果(图6b和图6e)，两种方法均去除大部分随机噪声，使得有效波变得清晰，但残差剖面中仍能发现较多的有效信号。而恢复后的剖面图6d、g中同相轴信息被较好地提取(见红圈及红框部分)，剖面变得更加自然、连续。不难发现，ASVD处理后再进行能量恢复的剖面在细节表现上要优于f-x去噪后再恢复的结果。

a—某地区二维实际资料;b—f-x域去噪剖面;c—f-x域去噪残差;d—f-x域去噪后恢复剖面;e—ASVD去噪剖面;f—ASVD去噪残差;g—ASVD去噪后恢复剖面a—2D actual data in a region;b—f-x domain denoising profile;c—f-x domain denoising residual;d—f-x domain denoising recovery profile;e—ASVD denoising proflie;f—ASVD denoising residual;g—ASVD domain denoising recovery profile图6 实际资料处理Fig.6 Actual data processing

3 结论与建议

本文基于K-SVD字典学习方法进行地震能量泄露恢复研究，通过楔形、倾斜模型以及实际资料处理分析得到以下结论:

1)以无噪或者噪声尽可能小的剖面作为样本进行字典学习、稀疏表示得到的去噪结果优于传统的字典学习方法。

2)利用一次去噪后的剖面作为样本进行学习，原始剖面的降噪性、保真性越好，在构建字典原子时信号的特征越能被更好的提取和表达，从而实现更加理想的恢复效果。

3)根据不同的地震信号可以构建不同的特征字典，具有较好的自适应性，并且在适用范围上有一定应用前景。