一种基于二维曲波变换的三维地震数据去噪方法
【技术领域】
[0001] 本发明涉及的是地震数据去噪的方法，具体是一种基于二维曲波变换的三维地震 数据去噪方法。 技术背景
[0002] 随着我国地震勘探进程的不断向前发展，山区，沙漠，厚黄土、砂土等复杂地区的 地震勘探项目逐渐增多，而这些地区表层激发条件都不够理想，再加上野外各种外界环境 的干扰，所采集到的地震数据包含各种严重的噪声干扰，掩盖了有效波信息，使得有效信号 同相轴模糊不清，能量相对较弱，信噪比低。尽管在野外可以采取相应的一些抗噪措施来提 高信噪比。但是，在多风季节或受其它复杂条件限制时，仅靠野外采集阶段完全压制随机噪 声是不可能，这就需要在室内进行有效地去噪工作，提高叠前地震资料的信噪比，以便后续 资料的处理。
[0003] 目前，基于数学变换的去噪方法较多，这些方法充分利用地震数据在数学变换域 内稀疏特点，通过对变换域中的系数进行阈值处理，从而达到去噪的目的。而为了达到所需 要的去噪效果，这就要求稀疏表示所采用的基函数能够捕捉到地震波前，只要保留少数较 大的稀疏系数就可以表示出原始数据的主要特征，而滤除的大量小系数不影响数据的主要 特征。
[0004] 假设采用多尺度多方向的二维曲波变换进行去噪，由于二维曲波变换由各向异性 的曲线状基元所构成，可以更加稀疏地表示地震波前特征，从而能够弥补了其他数学变换 方法的不足。然而在去噪领域中，常规方法常被用于压制二维地震数据的随机噪声，但随机 噪声存在于三维地震数据空间中，导致以往二维去噪方法无法获得满意的去噪结果。而对 于曲波变换如果想真正压制三维地震数据随机噪声，理论上讲应该需要使用三维曲波变 换，但三维曲波变换运算速度较慢，处理时间过长，远满足不了海量数据处理的要求。
[0005] 因此，本发明进行折中考虑，将二维曲波变换应用于三维地震资料去噪处理，提出 在每次去噪过程中逐次对时间切片进行多尺度多方向二维曲波正变换，从而得到不同尺度 层的曲波系数，同时根据有效信号系数和噪声系数的分布不一样，改变以往采用单一全局 阈值会损伤部分有效波的做法，提出采用局部阈值方法，也即每一尺度都采用一个阈值参 数，这样可以分别提取出每一尺度的有效波系数，最后再进行二维曲波反变换，从而完成整 个处理流程。
【发明内容】
[0006] 本发明的目的是为了能够快速高精度去除地震勘探数据中的噪声干扰，而提供了 一种基于二维曲波变换的三维地震数据去噪方法。
[0007] 本发明提出一种基于二维曲波变换的三维地震数据去噪方法，首先抽取含噪三维 地震数据的时间切片，对其进行多尺度多方向二维曲波变换得到曲波域系数，然后在曲波 域采用局部阈值法去噪方法，对二维曲波变换后的每一个尺度都选取一个软阈值算子，通 过阈值处理，得到各尺度下的有效波曲波系数，最后将提取出来的有效波曲波系数进行逆 变换重构出地震信号，从而达到去噪目的。
[0008] 一个含噪声的三维地震信号的模型可以表示成如下的形式：
[0009] f(i, j,t) = s(i, j,t)+k · e(i, j,t) i = 1,2, . . . ,m, j = 1,2,. . · ,n,t = 1,2,..., k
[0010] 式中f(ij，t)为三维地震数据含噪信号，s(ij，t)为三维不含噪地震信号，e(i， j，t)为噪声信号，k表示噪声水平值。去噪的过程就是从含噪信号f(i，j，t)中，提取真实信 号s( i，j，t)，去除噪声干扰信号e( i，j，t)。
[0011] 在数据去噪过程中，采用局部阈值法，具体去噪步骤如下：
[0012] (1)首先逐次抽取三维含噪地震数据的时间切片，然后选择分解尺度记为N，将含 噪时间切片进行N尺度二维曲波变换，得到含噪的曲波系数，这些系数中少量的较大曲波系 数能够代表信号本身，而大部分较小值的曲波系数则表示高频噪声干扰信号。
[0013] (2)根据各分解尺度上曲波系数的分布情况，选择与分解尺度有关的局部阈值参 数，以反映曲波系数在不同尺度上的不同特征，然后对该尺度曲波系数分量进行软阈值处 理，保留有效波曲波系数。
[0014] (3)根据经过软阈值处理后的各尺度有效波曲波系数分量，进行地震信号的二维 曲波反变换，对反变换后的时间切片进行组合，所得到的数据体即为本发明去噪后的三维 地震数据。
[0015] 讲一击.所怵一维曲浊亦拖的宙^/为.
[0016]
[0017]式中：ΦιΜ表示曲波函数，j，l，k分别表示尺度，方向和位置参数，f(x)为
[0018] 地震数据，其频率域定义式为：
[0019]
[0020] 经过变换后得到的曲波系数，用(：{」}{1}仏1氺2)表示其结构，其中」表示尺度，1表 示方向，（lu，k 2)表示j尺度1方向上的矩阵系数。
[0021] 进一步，局部阈值参数，其表达式如下：
[0022] \c = a〇j
[0023] 式中λ。表示局部阈值参数，a表示与尺度方向有关的经验参数，〇表示最高尺度层 曲波系数标准差的估计值，τ代表曲波域内所有曲波系数的标准差估计值。
[0024] 进一步，所述软阈值算子，其表达式如下：
[0025]
[0026] 式中F表示软阈值算子，λ。表示与分解尺度有关的局部阈值参数，不同分解尺度的 阈值参数不一样。
[0027] 本发明的优点:本发明采用了具有多尺度和多方向性的二维曲波变换进行三维地 震资料去噪，通过逐次对时间切片进行去噪处理，实现了基于二维曲波变换的三维地震资 料去噪处理方法，从而避免了计算速度慢，效率低的缺点，大幅度地提高了计算效率，节省 了运算时间。同时在去噪过程中，提出局部阈值参数的去噪方法，尽可能地保护了微弱的有 效波信号，从而使反射波同相轴更加连续、清晰，提高了地震数据的信噪比，降低了对计算 机内存的要求。
【附图说明】
[0028] 图1是本发明实施例中三维地震数据去噪流程图。
[0029] 图2是原始地震数据及其加噪对比图。
[0030] 图3是为含噪数据二维曲波6尺度分解图。
[0031] 图4是局部阈值去噪结果图。
[0032] 图5是去除的噪声剖面图。
【具体实施方式】
[0033]以下实施案例用于说明本发明，但不用来限制本发明的范围。
[0034] 实施例1
[0035]实现该方法的步骤主要包括，去噪方程的构建，二维曲波变换，局部阈值法，阈值 算子处理等。具体步骤如下：
[0036]步骤1:去噪方程的构建。一个含噪声的三维地震信号的模型可以表示成如下的形 式：
[0037] f(i,j,t)=s(i,j,t)+k?e(i,j,t) i = l,2,...,m,j = l,2,...,n,t = l,2,...,k
[0038] 式中f(ij，t)为三维地震数据含噪信号，s(ij，t)为三维不含噪地震信号，e(i， j，t)为噪声信号，k表示噪声水平值。去噪的过程就是从含噪信号f(i，j，t)中，提取真实信 号s( i，j，t)，去除噪声干扰信号e( i，j，t)。
[0039] 步骤2:首先逐次抽取三维含噪地震数据的时间切片，然后选择合适的分解尺度 (记为N)，将含噪时间切片进行N尺度二维曲波变换，得到含噪的曲波系数，这些系数中少量 的较大曲波系数能够代表信号本身，而大部分较小值的曲波系数则表示高频噪声干扰信 号。
[0040] 所述二维曲波变换的定义为：
[0041]
[0042] 式中：表示曲波函数，j，l，k分别表示尺度，方向和位置参数，f(x)为地震数 据，其频率域定义式为：
[0043
[0044] 经过变换后得到的曲波系数，可用(：{」}{1}仏1氺2)表示其结构，其中」表示尺度，1 表示方向，（lu，k 2)表示j尺度1方向上的矩阵系数。
[0045] 步骤3:根据各分解尺度上曲波系数的分布情况，选择与分解尺度有关的局部阈值 参数，以反映曲波系数在不同尺度上的不同特征，然后对该尺度曲波系数分量进行软阈值 量化处理，保留有效波曲波系数。
[0046] 去噪后的有效波信号可以下述方法估算出，即
[0047] S = C-HF(CS))
[0048] 式中C表示二维曲波变换，0-1表示逆二维曲波变换，F表示软阈值算子，其表达式如 下：
[0049]
[0050]式中λ。表示与分解尺度有关的阈值参数，不同分解尺度的阈值参数不一样，其表 达式如下：
[0051] Xc = a〇T
[0052] 式中λ。表示局部阈值参数，a表示与尺度方向有关的经验参数，〇表示最高尺度层 曲波系数标准差的估计值，τ代表曲波域内所有曲波系数的标准差估计值。