CN102105813A - 3d地震随机噪声的结构无关分析 - Google Patents

Abstract

一种利用3D算子在三维(3D)地震数据中识别和显示随机噪声以降低地震结构对噪声识别的影响的系统和方法。所述3D算子使用3D中所要求的性能的陈述导出。所述3D算子逐像素地应用到所述3D叠后数据中的每个像素，以在3D显示中显示图像或输出基本上与所述图像结构无关噪声的估计结果。结果显示以彩色产生，以表明噪声振幅，以方便在原始显示中的噪声区域定位。

Description

3D地震随机噪声的结构无关分析

技术领域

一般来说，本发明涉及图像处理领域，确切地说，涉及图像数据的抑制，以便在包含地质结构，比如断层的叠后(post-stacked)三维(3D)地震数据中估计和识别随机噪声。

背景技术

数据的可靠性是地震数据解释中的主要关注，因为钻探位置的选择以及钻探操作将成功的可能性都是基于做出钻探相关决策所用数据的准确性和可靠性。所以，随机噪声的估计和显示是数据可靠性评价中最重要的因素之一。这种评价又能够用在油藏工程研究中。

在地震勘探领域中，已经提出了许多噪声分析技术。在Dash，B和K.A.Obaidullah的“Determination of signal and noise statisticsusing correlation theory”，Geophysics，1970，35，24-32(“Dash和Obaidullah”)中介绍了在3D叠后地震数据的显示中提供噪声估计的一种方法。Dash和Obaidullah假设各道的图像信号相关而噪声不相关，应用相关理论从互相关推断出信号功率，从自相关推断出总功率。然后从差中提取噪声功率。缺点在于，当跨越断层或不整合时地震道不相关。所以这种噪声估计方法受地质结构影响而偏差。

Potter和Roden所介绍的方法使用水平分量记录估计垂直方向的噪声。这个方法仅仅适于强方向性噪声而不适用于3D叠后数据。在图像噪声方差估计中通常有三个基本步骤。它们是：(1)图像结构抑制，(2)局部方差估计以及(3)全局方差估计。它们之中最重要的步骤是图像结构抑制。见Potter，T.F.和R.B.Roden的“Seismicnoise estimation using horizontal components”，Geophysics，1999，No.4，pp.617-632。Rank等人公开了在水平和垂直方向上都使用简单的双抽头(two-tap)差分滤波器以抑制图像结构。已经发现，为了得到更好的结果应当串联滤波器。不过，这种简单的双抽头FIR滤波器也在滤波后的图像中遗留了许多边缘信息，使得在全局方差估计阶段需要特殊的后处理以校正噪声方差。见Rank，K.，M.Lendi和R.Unbehauen1999年的，“Estimation of Image Noise Variance，”IEEProc.-Vis.Image Signal Process.，Vol.146，No.2，pp.80-84。

为了在二维中更好地抑制图像结构，在J.Immerkaer的“FastNoise Variance Estimation”，Computer Vision and ImageUnderstanding，1996，Vol.64，No.2，pp.300-302(“Immerkaer)”中介绍了使用两个Laplacian滤波器的差作为掩码以滤波图像。在Immerkaer的文章中对于3D掩码或用于平滑以去除Immerkaer掩码的各向异性效应的装置没有提出建议。

以下专利公开了在包含图像和与数据可靠性无关的噪声的图像中降低噪声的方法。例如，USP 5461655公开了医学造影(radiography)成像背景下降低噪声的方法和装置。已公开的申请US 20040066978公开了结合医学成像使用的图像处理方法和图像处理装置。USP7085426和USP 7130484的公开内容介绍了如在医学应用中和武器运输的x射线扫描中使用Volterra滤波器增强噪声图像的轮廓。已公开的申请US20020012472公开了3D光学荧光显微图像时间序列的可视化方法，确切地说压缩4D图像数据以加速其可视化的方法(见图3至图6)。

USP6847738中公开了采用指定滤波器增强TV图像显示中清晰度的方法。最后，USP 5844627公开了在处理数字视频信号中降低空间噪声的结构和方法。

这些专利的现有技术方法都没有提高在随机噪声存在时识别有效结构的能力。所以，本发明的目的是基于随机噪声估计提供所得到的地震属性的可靠性评价。

发明内容

根据本发明，提供了利用新的3D算子在3D地震数据中识别和显示随机噪声以降低地震结构对噪声识别的影响的系统和方法。所述新算子使用3D中所要求的性能的陈述导出。所述算子然后逐像素地应用到3D叠后数据中的每个像素，以在3D显示中提供或估计基本上与图像结构无关的噪声。在优选实施例中，结果显示以彩色产生，以表明噪声振幅，从而方便在原始显示中的噪声区域定位。

作为替代，所述显示能够以灰度色调呈现。不过，正如本领域的技术人员将理解，通过使用彩色将方便视觉显示的解释任务。

如本文所用，术语“像素”意指所述3D叠后数据或所述数据的滤波后版本中最小的可寻址点，它给出了所述点的振幅。

本发明基于3×3的掩码，其成员由性能需求的陈述而不是由如在现有技术的Immerkaer掩码中的两个Laplacian滤波器的相减而形成。本发明的方法提供了Immerkaer掩码的特征从2D到3D的扩展。

因此，本发明通过使用3D掩码以抑制包含结构和随机噪声的3D图像中的结构而克服了现有技术的限制，从而有助于决定哪个图像数据可能是结构，哪个可能是噪声。

正如以下更详细的介绍，本发明利用基于性能的方法导出Immerkaer 2D掩码，然后将所述方法扩展到3D。另外，本发明提供了平滑以除去Immerkaer掩码的各向异性效应。本发明将在存在结构的情况下正确地执行随机噪声估计，并且将不受地震结构的影响。

附图说明

本专利或申请文件包含以彩色绘制的至少一幅图。美国专利和商标局在收到请求和必要费用支付后将提供具有彩色图的本专利或专利申请出版的副本。

以下参考附图进一步介绍本发明，其中：

图1是适于本发明实施的基于计算机的系统的示意图；

图2是具有垂直边缘的图像上的片段；

图3是3D掩码的图片表达；

图3A是本发明算法的一个实施例的流程图；

图4a是地震数据的2D显示；

图4b是使用噪声分析的现有技术方法时图4a中噪声的2D显示；

图4c是使用根据本发明的方法开发的2D掩码时图4a中噪声的2D显示；

图4d显示了表明图4b和图4c中随机噪声等级所用的彩色范围；

图5a是2D地震数据的另一个实例；

图5b是使用根据本发明的方法开发的2D掩码时图5a的2D数据中存在的噪声的估计结果；

图5c显示了表明图5b中随机噪声等级所用的彩色范围；

图6a是3D叠后地震数据的实例；

图6b显示了表明图6a中地震振幅所用的彩色范围；

图6c是将本发明的掩码应用到图6a的地震数据所达到的结果的呈现；

图6d显示了表明图6c中随机噪声等级所用的彩色范围；

图7显示了在3D叠后地震数据集中存在的噪声区域的时间切片；

图8(a)至图8(c)分别是垂直、水平和对角线台阶数据的矩阵表达；

图9是恒定斜率数据的矩阵表达。

具体实施方式

参考图1，图中显示了适于本发明实施的基于计算机的系统100，系统100包括非易失性存储器118，非易失性存储器118接收并存储来自数据采集系统(未显示)的叠后数据120。随机存取存储器120为图像分析及其噪声分量隔离提供工作空间，而非易失性存储器118为操作系统、输入/输出和3D显示软件提供存储。非易失性存储器118也为根据以下结合图3A讨论的掩码和算法抑制图像结构的指令提供存储。在优选实施例中，处理器116执行全部的内务处理和分析计算，而在替代实施例中分析计算在分开的数字信号处理器(DSP)(未显示)中执行以加速数值计算。叠后地震数据可以经由网络连接112、替代串行连接(未显示)或经由比如闪存或光盘的可移动介质输入到基于计算机的系统100中。优选Linux操作系统，但是在替代实施例中也可以使用Unix操作系统或Windows操作系统。

基于计算机的系统100的其他元件包括：键盘122，显示器126，优选情况下显示器126是可以结合打印设备(未显示)一起使用的彩色电子显示器，以及定点设备124。

将参考现有技术的方法和技术介绍本发明的3D掩码的开发过程。

Immerkaer(1996)公开的快速噪声方差估计算法基于两个拉普拉斯(Laplacian)滤波器的差而抑制图像结构的影响。拉普拉斯滤波器的一般形式可以表示为：

$L_{\mathrm{\α}}=\frac{1}{1+\mathrm{\α}}\×\left(\begin{array}{ccc}\mathrm{\α}& 1-\mathrm{\α}& \mathrm{\α}\\ 1-\mathrm{\α}& -4& 1-\mathrm{\α}\\ \mathrm{\α}& 1-\mathrm{\α}& \mathrm{\α}\end{array}\right),---\left(1\right)$

其中α是0与1之间的参数。不难证明，Immerkaer的掩码是L₁与L₀之间差的二倍。

$2(L_1-L_0)=\left(\begin{array}{ccc}1& 0& 1\\ 0& -4& 0\\ 1& 0& 1\end{array}\right)-2\left(\begin{array}{ccc}0& 1& 0\\ 1& -4& 1\\ 0& 1& 0\end{array}\right)=\left(\begin{array}{ccc}1& -2& 1\\ -2& 4& -2\\ 1& -2& 1\end{array}\right).---\left(2\right)$

Immerkaer掩码中固有的性质可以描述为：

(1)以Immerkaer掩码滤波恒定灰度级时，输出为零。Immerkaer掩码项的零和性质确保了这一点。

(2)以Immerkaer掩码滤波恒定斜率斜坡时，输出为零。Immerkaer掩码项的零和性质也保证了这个性质。图9中使用矩阵表达显示了恒定斜率斜坡。

(3)当图像具有垂直或水平台阶边缘时，如图8(a)和图8(b)中的矩阵表达所示，用Immerkaer掩码滤波的输出也为零。

(4)当图像具有对角线台阶边缘时，如由图8(c)中的矩阵表达所示，用Immerkaer掩码滤波的输出不为零。

(5)Immerkaer掩码在三个主轴方向上是对称的。

考虑这些特性，前三个是好的特性，并可用于抑制诸如恒定灰度级、斜坡边缘和台阶边缘之类的图像结构。已经发现第四个性质是Immerkaer模式的主要缺陷。这是因为所述算子的各向异性并能够用平滑模式部分地减轻。

本发明利用的掩码包括以下所期望的设计特性：

(1)所述掩码系数应当为零和。这个特性将保证所述掩码将抑制恒定灰度级和恒定斜率斜坡。

(2)当掩码的范围包含水平边缘或垂直边缘时，滤波的输出应当为零。

(3)掩码在三个主轴方向上应当是对称的。

(4)噪声方差守恒增加为随后的测试情况和3D掩码的要求。图像噪声被假设为独立的高斯(Gaussian)噪声。当将高斯噪声通过所设计的掩码时，该噪声便不再是空间独立的，而是空间相关的噪声。所以，噪声方差将不改变。

以上要求能够以公式表示为2D掩码的最优化问题，以便显示出基于掩码性能的最优化过程将导致2D的Immerkaer掩码。

根据以上的要求(1)，所述掩码(mask)应当具有以下结构：

$\mathrm{mask}=\left(\begin{array}{ccc}{x}_2& x_3& x_2\\ x_3& x_1& x_3\\ x_2& x_3& x_2\end{array}\right)---\left(3\right)$

在3式的掩码模板中，必须确定x₁、x₂和x₃三个系数。对称要求(1)已经由掩码系数x₁、x₂和x₃的位置所满足。零和要求(1)暗示以下关系：

x₁+4x₂+4x₃＝0.    (4)

根据要求(2)并假设如图2中所示的垂直边缘，导出下式：

2x₂+x₃＝0.    (5)

归一化要求(4)产生

$x_1^2+4x_2^2+4x_3^2=1.---\left(6\right)$

式4、式5和式6组合如下：

$\left\{\begin{array}{c}{x}_1+4x_2+4x_3=0\\ 2x_2+x_3=0\\ x_1^2+4x_2^2+4x_3^2=1\end{array}\right.---\left(7\right)$

解该方程组得到以下结果：

$\left(\begin{array}{c}{x}_1\\ x_2\\ x_3\end{array}\right)=\frac{1}{6}\left(\begin{array}{c}4\\ 1\\ -2\end{array}\right).---\left(8\right)$

将这个结果与Immerkaer掩码对比表明，它是Immerkaer掩码的归一化版本。这说明所述Immerkaer掩码也是在以上建立的准则下的最佳掩码。这个结果用作将掩码设计技术扩展到更高维的基础。

对于3D掩码，如同2D掩码的情况，对称要求由图3中展示的单元值的位置建立，在其配置下，3D对称掩码具有4个值要确定。

零和要求产生：

x₁+6x₂+12x₃+8x₄＝0.   (9)

边缘和拐角边缘结构的零输出要求产生：

x₂+4x₃+4x₄＝0.        (10)

x₃+2x₄＝0.            (11)

最后，归一化要求能够以公式表示为：

$x_1^2+6x_2^2+12x_3^2+8x_4^2=1.---\left(12\right)$

解式9至式12产生3D掩码单元为：

$\left(\begin{array}{c}{x}_1\\ x_2\\ x_3\\ x_4\end{array}\right)=\frac{1}{6\sqrt{6}}\left(\begin{array}{c}-8\\ 4\\ -2\\ 1\end{array}\right).---\left(12\right)$

作为替代，3D算子可以写为非归一化的整数格式：

$\mathrm{mask}=\left(\begin{array}{ccc}1& -2& 1\\ -2& 4& -2\\ 1& -2& 1\end{array}\right),\left(\begin{array}{ccc}-2& 4& -2\\ 4& -8& 4\\ -2& 4& -2\end{array}\right),\left(\begin{array}{ccc}1& -2& 1\\ -2& 4& -2\\ 1& -2& 1\end{array}\right).---\left(13\right)$

这个解具有类似于2D情况的某些特征。中心系数与对角线系数的符号相同，而与非对角线系数的符号不同。越接近中心的系数绝对值就越大。所以从x₁到x₄绝对值在减少。x₄系数的绝对值最小，因为从x₄到中心的距离最大。

根据本发明，完成抑制图像结构过程的方式为叠后图像中的每个像素值都以使用图3A中介绍的算法，通过将图像掩码值应用到图3中示意地显示的立方体中对应的像素所确定的新值替换。在优选实施例中，每个像素处的估计的噪声值然后用优选情况下具有11个像素窗口的矩形波串(box-car)滤波器平滑，以产生局部平均并降低新掩码的非各向同性特性。

实例

以下若干实例采用实际数据展示了将本发明的方法和系统应用到3D叠后地震数据的噪声估计的结果，连同以上介绍的矩形波串平滑后的结果。注意，在优选实施例中，每个像素幅度的范围从0到255，但是在其他实施例中，幅度的范围可以是其他值。

为了说明所述掩码的结构无关特性，将其应用到了图4a中的数据集，其结果在图4c中给出。同样，将Dash和Obaidullah的现有技术方法应用到相同的数据集，其结果在图4b中给出。对比图4b与图4c中所示的结果，可见在导致断层处道间相关性损失的300处的断层在图4b中显示为噪声区域310，但是在图4c中却没有出现。如图4b展示，Dash和Obaidullah的相关方法在断层位置产生了高噪声水平，因为相关方法不区分噪声和地质结构。如图4c所示，本发明的掩码没有显示出由于断层导致的任何偏差，并且产生了更精确噪声显示和估计。根据本领域公知的技术，图4c的高噪声区域以显著的颜色如红色显示在320处。

在典型的地震剖面中，高噪声区域往往伴随着复杂地质的不良数据采集或处理。另一方面，低噪声位置与强的反射层相关联。这展示在图5a中，其中在强反射层附近地震同相轴(event)清晰，而在“噪声泄漏”的垂直趋势中同相轴模糊。这些特征清楚地显示在图5b的噪声估计显示中。图5a与图5b的对比揭示了图5a中噪声的视觉观察结果与应用本发明的噪声估计方法的结果显示匹配得很好。同样，使用图6c中基本上仅有噪声的3D图像，不难识别图6a所示3D地震数据的高噪声/低可靠性区域。图6a的高噪声/低可靠性区域在图6c中的噪声估计体中以黑色显示。在优选实施例中，极高噪声区域在图6d的彩色图中以黑色显示和/或打印。图6b表示的显示使用了表明图6a的地震振幅显示的彩色范围。

使用例如图7中所示数据时间切片完成数据可靠性和/或噪声存在的附加确认。图7的高噪声区域以显著的颜色如红色指明。

通过目前优选的实施例和实例已经详细介绍了本发明。不过，根据本说明书，对本方法的其他实施例和修改对本领域的普通技术人员将是显而易见的。所以，本发明的范围由随后的权利要求书确定。

Claims (20)

1.一种用于从包含原始像素集以表示多维空间中结构和噪声的叠后多维地震图像中抑制结构图像分量的方法，所述方法包括以下步骤：

a.在处理器处接收地震图像数据；

b.通过输入设备接收以对在三维中所述地震图像数据的噪声抑制所要求的性能为特征的多维掩码的参数；

c.使用执行预定的图像处理程序的所述处理器，对第一多维空间中所述地震图像的所述原始像素集中的每个像素应用所述多维掩码，以在第二多维空间中产生新像素集；

d.在显示设备上显示所述第一多维空间；

e.在所述显示设备上显示所述第二多维空间；以及

f.对比显示以识别所述第一多维空间中的噪声，从而方便从所述叠后多维地震图像中抑制结构图像分量，包括所识别的噪声。

2.根据权利要求1的方法，其中，选择所述多维掩码参数以在所述第二多维空间中基本上消除所述结构图像分量的影响。

3.根据权利要求2的方法，其中，所述多维空间具有两个维度。

4.根据权利要求2的方法，其中，所述多维空间具有三个维度。

5.根据权利要求1的方法，其中，显示步骤(d)和(e)中每一个都包括：

f.在显示器上以彩色显示相应多维空间，从而指明所述地震图像数据中噪声区域的噪声振幅以预定的显著颜色或黑色出现。

6.根据权利要求5的方法，其中，所述步骤(f)包括以灰度显示所述相应多维空间。

7.根据权利要求1的方法，其中，所接收的参数指定所述多维掩码中的多个掩码系数具有零和，从而所述多维掩码抑制恒定灰度级别和恒定斜率斜坡。

8.根据权利要求1的方法，其中，所接收的参数指定，当使用所述多维掩码处理与水平边缘或垂直边缘对应的地震图像数据时，应用所述多维掩码的输出为零。

9.根据权利要求1的方法，其中，所接收的参数指定所述多维掩码在三个主轴方向上是对称的。

10.根据权利要求1的方法，其中，所接收的参数指定，当处理所述地震图像数据时，所述多维掩码执行噪声方差守恒，从而使所述地震图像数据中的所述噪声相关。

11.一种用于从包含原始像素集以表示多维空间中结构和噪声的叠后多维地震图像中抑制结构图像分量的系统，所述系统包括：

输入设备，用于输入以对在三维中所述地震图像数据的噪声抑制所要求的性能为特征的多维掩码的参数；

显示设备；以及

执行预定的图像处理程序的处理器，其中，所述处理器处理所述地震图像数据的方式为对第一多维空间中所述地震图像的所述原始像素集中的每个像素应用所述多维掩码，以在第二多维空间中产生新像素集；

其中，所述显示设备显示第一和第二多维空间，以允许用户对比所显示的多维空间以识别所述第一多维空间中的噪声，从而方便从所述叠后多维地震图像中抑制结构图像分量，包括所述识别的噪声。

12.根据权利要求11的系统，其中，选择所述多维掩码参数以在所述第二多维空间中基本上消除所述结构图像分量的影响。

13.根据权利要求12的系统，其中，所述多维空间具有两个维度。

14.根据权利要求12的系统，其中，所述多维空间具有三个维度。

15.根据权利要求11的系统，其中，显示器以彩色显示相应多维空间，从而表明指明所述地震图像数据中噪声区域的噪声振幅。

16.根据权利要求15的系统，其中，所述显示设备以灰度显示所述相应多维空间。

17.根据权利要求11的系统，其中，所接收的参数指定所述多维掩码中的多个掩码系数具有零和，从而所述多维掩码抑制恒定灰度级别和恒定斜率斜坡。

18.根据权利要求11的系统，其中，所接收的参数指定，当所述处理器使用所述多维掩码使用水平边缘或垂直边缘对应的所述地震图像数据进行处理时，应用所述多维掩码的输出为零。

19.根据权利要求11的系统，其中，所接收的参数指定所述多维掩码在三个主轴方向上是对称的。

20.根据权利要求11的系统，其中，所接收的参数指定当处理所述地震图像数据时，所述多维掩码执行噪声方差守恒，从而使所述地震图像数据中的所述噪声相关。

Images