CN102854533A - 一种基于波场分离原理提高地震资料信噪比的去噪方法 - Google Patents

Abstract

本发明提供了一种提高地震资料信噪比的去噪方法，属于石油地球物理勘探领域。本发明利用数据分块技术与波场分离技术相结合来压制随机噪声进而提高地震资料的信噪比；所述方法首先对地震数据进行分块处理，然后对每个数据块进行复数域SVD滤波，利用波场分离技术进行特征值分解，利用优选特征值重构数据，并替换原始频率切片进而实现压制随机噪声。与传统的假设同相轴水平的SVD技术相比，该方法无需对同相轴进行任何相对时差校正，且利用本方法的去噪效果比较稳健，本发明采用的分块处理技术能有效地克服实际地震资料中倾角变化的非稳态，且可以大大缩短SVD和重构占用的时间。

Description

一种基于波场分离原理提高地震资料信噪比的去噪方法

技术领域

本发明属于石油地球物理勘探领域，具体涉及一种基于波场分离原理提高地震资料信噪比的去噪方法。

背景技术

随机噪声是地震资料处理中遇到的一类主要干扰波，在地震记录上主要表现为波形杂乱无章，频带很宽，视速度不确定，无一定的传播方向，在一定程度上能够满足统计规律。随机噪声的存在会降低地震资料的信噪比，影响资料处理的质量。压制随机噪声的方法主要分为3类：

1，基于预测技术的方法(f-x反褶积和投影反褶积等)；

2，基于统计原理的方法(中值滤波和多项式拟合等)；

3，基于数学变换的方法(小波包[1]、奇异值分解技术(SVD/K-L)[2]、Radon变换[3]和模式分解等(此处的[1]、[2]、[3]指参考文献，其中，[1]王振国，汪恩华.小波包相关阀值去噪[J].石油物探，2002，41(4)：400-405；[2]吕景贵，刘振彪，管叶君.压制叠前相干噪音的速度变换域滤波方法[J].石油物探，2001，40(4)：94-99；[3]张军华，吕宁，雷凌，等.抛物线拉冬变换消除多次波的应用要素分析[J].石油地球物理勘探，2004，39(4)：398-405))。

其中f-x(y)域反褶积是地震资料处理中压制随机噪声的主流技术。

从20世纪70年代开始，SVD技术被广泛应用于信号处理，包括SVD(基于奇异值分解)技术压制噪声及有损压缩数字图像；噪声压制、倾角滤波、垂直地震剖面数据波场分离和剩余静校正。这些技术计算量大，需要复杂的修正才能适用于倾角数据。而且基于SVD的噪声压制技术是在t-x域中实现的，一般要求信号为水平，并不适合同相轴交错的情况。当有效信号同相轴呈倾斜和弯曲时，SVD会损害有效信号

发明内容

本发明的目的在于解决上述现有技术中存在的难题，提供一种基于波场分离原理提高地震资料信噪比的去噪方法，将数据分块技术与波场分离技术相结合来压制随机噪声，无需对同相轴进行任何相对时差校正，在假设线性同相轴条件下能适应同相轴交错的复杂情况。

本发明是通过以下技术方案实现的：

一种基于波场分离原理提高地震资料信噪比的去噪方法，所述方法利用数据分块技术与波场分离技术相结合来压制随机噪声进而提高地震资料的信噪比；所述方法首先对地震数据进行分块处理，然后对每个数据块进行复数域SVD滤波，利用波场分离技术进行特征值分解，利用优选特征值重构数据，并替换原始频率切片进而实现压制随机噪声。

所述方法包括以下步骤：

(1)按时间-空间方向对地震数据进行分块处理，获得一组数据块；

(2)复数域SVD滤波步骤：将步骤(1)中分解得到的每一个数据块变换到频率域中，构建Hankel矩阵，然后利用波场分离技术对Hankel矩阵进行特征值分解，并从大到小选取1或2个奇异值，最后利用优选特征值重构原始数据，并替换原始频率切片进而实现压制随机噪声。

所述步骤(1)包括以下步骤：

(11)计算时空窗的个数：

N＝(T-S-D-OW)/(1-O)W            (1-1)

(1-1)式中：T为纵、横向最大长度向量；S为纵、横向起始向量；D为纵、横向延迟向量；O为纵、横向重叠百分比；W为时空窗纵、横向跨度；

(12)计算起始位置：

T₀＝S+D+(i-l)W(1-O)    1≤i≤N       (1-2)

(13)采用镶边函数对边界进行平滑过渡，所述镶边函数采用Hanning函数，即：

$W\left(x\right)=\left\{\begin{array}{cc}\frac{1}{2}[1+\cos \left(\frac{\mathrm{\πx}}{a}\right)]& \left|x\right|\≤a\\ 0& \left|x\right|>a\end{array}\right.---(1-3)$

(1-3)式中：a为镶边带长度。

(14)设分时窗操作算子为P，Cadzow滤波算子为F，则整个数据分块处理过程表示为：

$\stackrel{\‾}{d}=\left[W_1{P}^{\′}{W}_2\mathrm{FP}\right]d---(1-4)$

(1-4)式中：W₁和W₂表示分块正、反过程的加权系数算子，由(1-3)式计算得到的，d指的是原始地震数据。

步骤(2)包括以下步骤：

(21)对每一个数据块进行傅里叶正变换；

(22)对选定频率范围的，构建n×n的复数Hankel矩阵A，其元素为每道的时间频率变换值；选择合适的奇异值个数k，利用SVD计算获得选取有限个奇异值后构建的矩阵F_k(A)；沿反对角线平均，恢复矩阵Hankel结构；利用F_k(A)中的值替换每道Fourier变换值；

(23)进行傅里叶反变换，获得计算结果。

与现有技术相比，本发明的有益效果是：

1)与其它传统的SVD特征值技术不同，本发明的方法对同相轴倾斜或交错情况没有任何特殊要求，并且不需要对同相轴进行相对时差校正；

2)假定特定时空范围内剖面含有有限个倾角，模型数据和实际数据处理结果表明，只取第一个特征值会对有效信号造成显著损害，而本发明取前3个显著特征值进行去噪效果比较稳健。

3)在小时空窗内地震数据的倾角变换可以看作是稳态的，实际资料处理结果表明，本发明采用分块处理的方式可以克服数据倾角的非稳态变化。

4)SVD和重构需要的时间随着数据体的增大而增大，本发明采用分块处理的方法可以大大缩短SVD和重构占用的时间。

附图说明

图1是本发明方法实施例中的分块处理的示意图。

图2是本发明方法实施例中的复数域波场分离的效果图。

图3是采用现有技术中的常规f-x域反褶积滤波结果图。

具体实施方式

下面结合附图对本发明作进一步详细描述：

一种提高地震资料信噪比的去噪方法，所述方法利用数据分块技术与波场分离技术相结合来压制随机噪声进而提高地震资料的信噪比；所述方法首先对地震数据进行分块处理，然后对每个数据块进行复数域SVD滤波，利用波场分离技术(又称为复数域奇异值分解技术、奇异值滤波、奇异值分解及重构，简称SVD)进行特征值分解，利用优选特征值重构数据，并替换原始频率切片进而实现压制随机噪声。

所述方法包括以下步骤：

(1)数据分块处理步骤：将地震数据按时间-空间方向对地震数据进行分块处理(图1给出的是分块的示意图)

地震信号同相轴在不同的时空位置上倾角不同，从统计规律上来说，同相轴的倾角变化为非稳态的。应用奇异值滤波一般选取较小的秩数字(1-3)，即假定倾角个数不超过3，如果将其应用在整个剖面上，难免会损害有效信号。在数据处理过程中，本发明采用分块处理的方法来克服上述非稳态的矛盾。另外奇异值分解及重构(即SVD技术)所占用的时间随数据体阶数的增大呈级数增加，而本发明采用数据分块处理处理(又称为数据分时空窗处理、时窗分割处理等)可以大大缩短SVD所占用的时间。

分块个数的计算公式为

N＝(T-S-D-OW)/(1-O)W    (1-1)

式中：T为纵、横向最大长度向量；S为纵、横向起始向量；D为纵、横向延迟向量；O为纵、横向重叠百分比；W为时空窗纵、横向跨度。

起始位置计算公式为

T₀＝S+D+(i-l)W(1-O)   1≤i≤N    (1-2)

为了减弱截断效应，在数据分块处理时必须对边界进行平滑过渡，采用的镶边函数为Hanning函数，即

$W\left(x\right)=\left\{\begin{array}{cc}\frac{1}{2}[1+\cos \left(\frac{\mathrm{\πx}}{a}\right)]& \left|x\right|\≤a\\ 0& \left|x\right|>a\end{array}\right.---(1-3)$

式中：a为镶边带长度。

设分时窗操作算子为P，Cadzow滤波算子为F，则整个数据分块处理过程可以表示为

$\stackrel{\‾}{d}=\left[W_1{P}^{\′}{W}_2\mathrm{FP}\right]d---(1-4)$

式中：W₁和W₂表示分块正、反过程的加权系数算子，由(1-3)式计算得到的，d指的是原始地震数据。(1-4)是用数学向量的形式描述“分块-＞处理-＞合并”的过程，描述的是整个过程。

综上所述，利用(1-1)式首先计算每个空间方向的块数，然后利用(1-2)式计算每个分块的起始位置，同时利用(1-3)式计算每个分块的加权系数及反加权系数，这样就完成了数据分块处理。

(2)复数域SVD滤波步骤：将步骤(1)中分解得到的每一个数据块变换到频率域中，构建Hankel矩阵，然后利用复数域奇异值分解(SVD)技术对Hankel矩阵进行特征值分解，从大到小选取1或2个奇异值(又称为优选特征值)，最后重构数据。

对于时间频率域中的数据，数据的有效秩数不超过数据的倾角个数。设有n×n空间网格，其空间坐标为(x_p，y_q)，p＝1…n，q＝1…n。设B为n×n维矩阵，其元素为每一道频率为ω的离散傅里叶变化值。如果只有k个倾角，矩阵B具有如下形式：

${\left[B\right]}_{p,q}=\underset{i=1}{\overset{k}{\mathrm{\Σ}}}{b}_i{e}^{j(c_i{x}_p+d_i{y}_q)}=\underset{i=1}{\overset{k}{\mathrm{\Σ}}}{b}_i{e}^{j{c}_i{x}_p}{e}^{j{d}_i{y}_q}---(2-1)$

式中，c_i和d_i为与倾角和频率有关的实数；b_i为复数值，其幅度为该频率处倾斜同相轴的振幅。定义n个向量f_i和g_i(i＝1…k)：

$\begin{array}{cc}{\left[f_i\right]}_p=b_i{e}^{j{c}_i{x}_p}& p=1...n\end{array}---(2-2)$

$\begin{array}{cc}{\left[g_i\right]}_q=e^{j{d}_i{y}_q}& q=1...n\end{array}---(2-3)$

则有

$B=\underset{i=1}{\overset{k}{\mathrm{\Σ}}}{f}_i{g}_i^T---(2-4)$

矩阵B的各列为f_i(i＝1…k)的线性组合，因此矩阵B的秩最大为k，前k个特征值不为零，可以用前k个特征图像加权和重构数据。

对于二维频率域地震数据，利用Hankel矩阵估算二维数据的有效秩数。设有N个地震道(此N与式(1-1)中的时空窗个数N不是同样的意思，此处的为标量，而上面的为矢量)(这些空间位置不必等间隔分布)，其空间坐标为x_p，p＝1…N。对于一个时间频率ω，其空间采样值分别为

t₁，t₂，....，t_N

假设地震信号具有k个倾角，根据f-x域地震信号具有预测性的特点可以将地震信号表示为k个正弦信号之和

$t_p=\underset{i=1}{\overset{k}{\mathrm{\Σ}}}{a}_i{e}^{\mathrm{ip}{b}_i}---(2-5)$

式中：a_i为复数；b_i为实数。

数据Hankel矩阵形式为

矩阵A的最大秩数等于信号的倾角个数。一般选取

如果无噪剖面只有k个倾角，则F_k(A)＝A。上面的计算说明，矩阵的秩数小于等于同相轴的个数。

本发明方法的软件实现流程如下：

分块

对每一数据块

{

Fourier正变换；

对选定频率范围

{

构建n×n的复数Hankel矩阵，其元素为每道的时间频率变换值；

选择合适的k值，利用SVD计算F_k(A)；

沿反对角线平均，恢复矩阵Hankel结构；

利用F_k(A)中的值替换每道Fourier变换值

}

Fourier反变换；

计算结果

}

合并分块。

根据上述流程，在Cygwin(一个在windows平台上运行的unix模拟环境，是cygnus solutions公司开发的自由软件)下利用fortran语言实现了本发明的方法。实施细节及参数如下：分块大小7道*200采样点，选取头三个最大奇异值。

本发明将数据分块技术与波场分离技术相结合来压制随机噪声，在f-x域首先构建Hankel矩阵，然后利用复数域奇异值分解(SVD)技术进行特征值分解，利用优选特征值的加权，并替换频率切片实现压制随机噪声，与传统的假设同相轴水平的SVD技术相比，本发明无需对同相轴进行任何相对时差校正，在假设线性同相轴条件下能适应同相轴交错的复杂情况。模型数据表明，如果地层的倾角个数与选取的特征值个数相同，则对无噪数据没有损害。图2中a区域为合成记录，中间b区域为采用本发明方法处理的结果，右边c区域显示的是a-b的结果图，即原始数据与采用本发明方法处理结果的差。图3中a区域为f-x反褶积处理结果，b区域为2a-a的结果图，即原始数据与反褶积处理结果的差，从图2的c区域和图3的b区域可以看出，常规f-x域反褶积对有效信号损害比较厉害，而本方法基本上看不出对有效信号的损害，因此在处理边界数据时，本发明与f-x域反褶积相比具有明显的优势。

分块处理能有效地克服实际地震资料中倾角变化非稳态的情况。实际资料处理结果表明，本发明是一项稳健的压制随机噪声的方法。

上述技术方案只是本发明的一种实施方式，对于本领域内的技术人员而言，在本发明公开了应用方法和原理的基础上，很容易做出各种类型的改进或变形，而不仅限于本发明上述具体实施方式所描述的方法，因此前面描述的方式只是优选的，而并不具有限制性的意义。

Claims (4)

1.一种基于波场分离原理提高地震资料信噪比的去噪方法，其特征在于：所述方法利用数据分块技术与波场分离技术相结合来压制随机噪声进而提高地震资料的信噪比；所述方法首先对地震数据进行分块处理，然后对每个数据块进行复数域SVD滤波，利用波场分离技术进行特征值分解，利用优选特征值重构数据，并替换原始频率切片进而实现压制随机噪声。

2.根据权利要求1所述的基于波场分离原理提高地震资料信噪比的去噪方法，其特征在于：所述方法包括以下步骤：

(1)按时间-空间方向对地震数据进行分块处理，获得一组数据块；

(2)复数域SVD滤波步骤：将步骤(1)中分解得到的每一个数据块变换到频率域中，构建Hankel矩阵，然后利用波场分离技术对Hankel矩阵进行特征值分解，并从大到小选取1或2个奇异值，最后利用优选特征值重构原始数据，并替换原始频率切片进而实现压制随机噪声。

3.根据权利要求1所述的基于波场分离原理提高地震资料信噪比的去噪方法，其特征在于：所述步骤(1)包括以下步骤：

(11)计算时空窗的个数：

N＝(T-S-D-OW)/(1-O)W    (1-1)

(1-1)式中：T为纵、横向最大长度向量；S为纵、横向起始向量；D为纵、横向延迟向量；O为纵、横向重叠百分比；W为时空窗纵、横向跨度；

(12)计算起始位置：

T₀＝S+D+(i-l)W(1-O)   1≤i≤N    (1-2)

(13)采用镶边函数对边界进行平滑过渡，所述镶边函数采用Hanning函数，即：

$W\left(x\right)=\left\{\begin{array}{cc}\frac{1}{2}[1+\cos \left(\frac{\mathrm{\πx}}{a}\right)]& \left|x\right|\≤a\\ 0& \left|x\right|>a\end{array}\right.---(1-3)$

(1-3)式中：a为镶边带长度。

(14)设分时窗操作算子为P，Cadzow滤波算子为F，则整个数据分块处理过程为：

$\stackrel{\‾}{d}=\left[W_1{P}^{\′}{W}_2\mathrm{FP}\right]d---(1-4)$

(1-4)式中：W₁和W₂表示分块正、反过程的加权系数算子，由(1-3)式计算得到的，d指的是原始地震数据。

4.根据权利要求1所述的基于波场分离原理提高地震资料信噪比的去噪方法，其特征在于：所述步骤(2)包括以下步骤：

(21)对每一个数据块进行傅里叶正变换；

(22)对选定频率范围的，构建n×n的复数Hankel矩阵A，其元素为每道的时间频率变换值；选择合适的奇异值个数k，利用SVD计算获得选取有限个奇异值后构建的矩阵F_k(A)；沿反对角线平均，恢复矩阵Hankel结构；利用F_k(A)中的值替换每道Fourier变换值；

(23)进行傅里叶反变换，获得计算结果。

Images