[发明专利]一种压制地震资料异常噪声的方法、存储介质及设备

权利要求书

1.一种压制地震资料异常噪声的方法，其特征在于，包括以下步骤：

S1、读取含有高振幅异常噪声的原始地震资料，根据地震资料中噪声的复杂程度选取混合高斯分布的高斯成分个数和迭代次数；

S2、针对原始地震资料中有效信号的特点选取固定稀疏变换基；

S3、利用线性化Bregman方法得到地震资料中有效信号的稀疏表示，根据稀疏表示x^*与固定稀疏变换基的乘积和原始地震资料得到关于有效信号的估计和异常噪声估计，将异常噪声建模为混合高斯分布，利用EM算法更新每个采样点上异常噪声的权重矩阵w；

S4、重复步骤S3直到原始地震资料中异常噪声压制完成。

2.根据权利要求1所述的方法，其特征在于，步骤S1中，将原始地震资料中二维地震剖面y表示为有效信号s与异常噪声e的线性叠加；利用混合高斯分布拟合异常噪声e，根据异常噪声e的复杂程度确定迭代次数。

3.根据权利要求1所述的方法，其特征在于，步骤S2中，基于wrapping方法的离散Curvelet，其正变换步骤如下：

S201、对二维信号做二维FFT，得到信号的二维傅里叶矩阵如下：

其中，(n₁,n₂)表示数据中的采样点；

S202、对每个尺度和角度方向的频率矩阵做加窗处理如下：

其中，为窗函数，j表示第j个尺度，l表示第l个角度；

S203、对步骤S202中得到的每个矩阵做环绕处理，标准遵循Candes准则得到如下：

其中，W表示对加窗后的频率矩阵做环绕处理；

S204、对步骤S203中得到的做二维FFT处理获得Curvelet变换系数矩阵C^D(j,l,k)，其中，k表示位移。

4.根据权利要求3所述的方法，其特征在于，Curvelet反变换的实现步骤如下：

S205、将每个尺度和角度方向的C^D(j,l,k)序列做二维FFT变换得到如下傅里叶序列如下：

S206、将步骤S205中得到的每个尺度和方向对的傅里叶矩阵与对应窗的环绕矩阵相乘得到数据如下：

S207、将步骤S206中的数据进行解环绕处理，将得到的数据相加得到原始信号的傅里叶矩阵如下：

S208、对步骤S207中得到的数据进行IFFT处理得到重构的原始信号。

5.根据权利要求1所述的方法，其特征在于，步骤S3中，迭代过程具体为：

混合高斯分布中高斯成分个数K，总迭代次数I，曲波变换字典D；每步迭代j增加1，并计算x^j；当小于预设的值或迭代次数j＝I时，迭代终止；输出稀疏表示系数x^*＝x^j和有效信号

6.根据权利要求5所述的方法，其特征在于，稀疏系数的目标函数如下：

其中，D∈R^m×n为固定字典变换基，x∈R^n×1为对应于字典D的稀疏系数，⊙为Hadamard乘积，y为原始地震资料中二维地震剖面。

7.根据权利要求5所述的方法，其特征在于，异常噪声e建模为混合高斯分布如下：

其中，表示第k个以0为均值、为方差的高斯分布，共K个高斯分布；π_k为噪声混合比例。

8.根据权利要求5所述的方法，其特征在于，混合高斯分布的权重矩阵w更新公式为：

其中，i＝1,2,...,m，k＝1,2,...K，m为含噪数据y中的元素个数，γ_ik为自定义变量，为方差的高斯分布。