CN108072892B - 一种自动化的地质构造约束层析反演方法 - Google Patents

Abstract

本发明公开了一种一种自动化的地质构造约束层析反演方法，包括以下步骤：首先提取观测数据；利用观测数据，在初始速度模型中反传，进行叠前深度偏移，得到初始的偏移剖面；从初始的偏移剖面中自动提取地层构造属性；根据观测数据、初始速度模型以及地层构造属性，利用射线追踪或敏感核函数构建构造约束层析方程组；迭代反演求解速度更新量，更新速度模型，最终更新后获得最终的速度模型。本发明中的整个层析反演迭代过程全自动化，不需要进行人工干预，这提高了整个速度估计过程的效率，减少了人工工作量，缩短了处理周期；由于引入了地下地质倾角来约束速度模型，反演得到的速度模型具有较高的精度，满足后续偏移成像的要求。

Description

一种自动化的地质构造约束层析反演方法

技术领域

本发明涉及油气勘探开发中的地震成像与反演领域，具体涉及一种自动化的地质构造约束层析反演方法。

背景技术

地震勘探技术的目标是利用地震波成像技术，实现对地下结构的定位、识别和描述，为地下油气藏的勘探提供直观、可靠的依据。而地震波成像技术主要包括偏移成像和反演成像两个方面。偏移成像的本质是利用观测到的地震波场记录进行正、反向传播，同时消除地震波的传播效应，最后获取地下地质结构图像的过程；反演成像的本质是根据观测数据和地球物理模型参数之间的函数关系，逆向映射求取地球物理模型的过程。因此，从本质上说，反演成像的应用比偏移成像的应用更为广泛。常规的地震反演成像主要包括地震层析反演、最小二乘叠前深度偏移和AVO/AVA反演这三项核心技术，而地震层析反演技术更是后两项技术得以成功实现的基础和前提。

从理论上，地震层析反演技术主要有两个方向：基于射线理论的射线类层析和基于波动理论的波动方程层析。而射线束层析是近年研究的热点，射线束层析是介于射线类层析和波动方程层析之间的折中方法，兼顾了层析效率高和相对稳定的优点，其类型主要有高斯束层析、菲涅尔体层析、胖射线层析和高斯波包层析等。

常规射线层析具有灵活、高效的优势，然而其对应的反问题往往是稀疏、病态的，所以在反演求解层析方程的过程中需要加入对反演模型的先验约束。在反问题的求解过程中加入模型的构造约束是目前最有效的层析预条件技术，目前工业界常用的构造约束技术能较好地控制地下速度模型的趋势，但往往需要人工拾取地下地质层位，造成工作流程繁琐、耗时耗力。

地震层析反演技术能估计出较为准确的地下宏观速度模型，是一种兼具效率与实用性的速度反演技术，然而由于数据采集不理想以及射线理论的限制，反演过程往往是病态的。常规解决方案是在偏移剖面上人工拾取地下地质层位，并以此约束地下速度模型，该方法在实际操作过程中十分耗时耗力，甚至达不到勘探开发的效率需求。

因此，研发出一种构造约束层析反演方法，避免人工拾取地下地质层位的繁重工作，从而对于快速、准确构建地下速度模型，使得勘探开发效率更高。

发明内容

针对地震层析反演技术中存在的偏移剖面上人工拾取地下地质层位，并以此约束地下速度模型，该方法在实际操作过程中十分耗时耗力，甚至达不到勘探开发的效率需求的问题，本发明提出一种自动化的构造约束层析反演方法，旨在避免人工拾取地下地质层位的繁重工作，从而对于快速、准确构建地下速度模型，使得勘探开发效率更高。

本发明提供的一种自动化的地质构造约束层析反演方法，包括以下步骤：

S100：首先提取观测数据；

S200：利用观测数据，在初始速度模型中反传，进行叠前深度偏移，得到初始的偏移剖面；

S300：从初始的偏移剖面中自动提取地层构造属性；

S400：根据观测数据、初始速度模型以及地层构造属性，利用射线追踪计算层析敏感核函数，构建构造约束层析方程组；

S500：迭代反演求解速度更新量，更新速度模型，最终更新后获得最终的速度模型。

进一步的，步骤S200中所述叠前深度偏移为高斯束叠前深度偏移。

进一步的，步骤S300中所述地层构造属性包括地层倾角信息和地层位置信息。

进一步的，步骤S300中所述的自动提取的方法为采用结构张量算法处理偏移剖面，提取偏移剖面中的地层倾角信息。

进一步的，所述结构张量算法如下：

其中g_x为地震图像沿水平方向的梯度，

g_y为地震图像沿垂直方向的梯度，

<·>为二维高斯光滑滤波，

G为结构张量算子。

进一步的，步骤S400中所述的构造约束层析方程组为：

S^TL^TLSu＝S^TL^Tτ

其中L为射线层析核函数，

S为预条件算子，

τ是正演计算地震波传播走时与接收数据走时的差，

u为层析反演的速度模型更新量；

其中预算子S的表达式如下：

S＝(I)+D^TGD^-1

其中，I是单位阵，

D是梯度算子

D^T是D的转置矩阵

G是结构张量算子。

进一步的，采用最小二乘QR方法(LSQR)求解构造约束层析方程组S^TL^TLSu＝S^TL^Tτ。

本发明提供的一种自动化的地质构造约束层析反演方法。有现有技术相比具有以下优势：第一，本发明中的整个层析反演迭代过程全自动化，不需要进行人工干预，这提高了整个速度估计过程的效率，减少了人工工作量，缩短了处理周期；第二，由于引入了地下地质倾角来约束速度模型，反演得到的速度模型具有较高的精度，满足后续偏移成像的要求。

本发明的其他优点、目标和特征在某种程度上将在随后的说明书中进行阐述，并且在某种程度上，基于对下文的考察研究对本领域技术人员而言将是显而易见的，或者可以从本发明的实践中得到教导。本发明的目标和其他优点可以通过下面的说明书，权利要求书，以及附图中所特别指出的结构来实现和获得。

附图说明

附图用来提供对本发明的进一步理解，并且构成说明书的一部分，与本发明的实施例共同用于解释本发明，并不构成对本发明的限制。在附图中：

图1是本发明实施例一中的自动化的地质构造约束层析反演方法的流程图；

图2是本发明实施例二中采用的多层背斜理论模型；

图3是本发明实施例二中的多层背斜理论模型的平滑速度模型；

图4是本发明实施例二中的多层背斜理论模型的正演地震记录；

图5是本发明实施例二中的初始速度模型；

图6是本发明实施例二中的初始速度叠前深度偏移剖面图；

图7是本发明实施例二中的在初始偏移剖面上提取的地层倾角信息剖面图；

图8是本发明实施例二中的地质构造约束层析反演迭代50次更新后的速度模型；

图9是本发明比较例中无地质构造约束层析反演迭代50次更新后的速度模型；

图10是初始速度模型、正确速度模型、本发明实施例二地质构造约束层析反演速度模型和比较例无地质构造约束层析反演速度模型在CIP＝488(其中CIP是指commonimage point，即共成像点)处的偏移后生成的角道集的对比图；

其中，图10a为初始速度模型偏移后的角道集；图10b为正确速度模型偏移后的角道集；图10c为本发明实例二中的地质构造约束层析反演速度模型偏移后的角道集；图10d为本发明比较例中无地质构造约束层析反演速度模型偏移后的角道集。

另外，在图2、图3、图5、图8和图9四幅图中，左侧为色标，右侧为图件，其中图件的横坐标和纵坐标的单位均为km，即千米；色标的单位为m/s，即米/秒；图6和图7的横坐标和纵坐标的单位均为km，即千米。

具体实施方式

下面将结合附图对本发明作进一步说明。

实施例一

如图1所示为本实施例中的自动化的地质构造约束层析反演方法的流程图，本实施例中的自动化的地质构造约束层析反演方法，包括以下步骤：

S100：首先提取观测数据；

S200：利用观测数据，在初始速度模型中返传，进行叠前深度偏移，得到初始的偏移剖面；

S300：从初始的偏移剖面中自动提取地层构造属性；

S400：根据观测数据、初始速度模型以及地层构造属性，利用射线追踪计算敏感核函数，构建构造约束层析方程组；

S500：迭代反演求解速度更新量，更新速度模型，最终更新后获得最终的速度模型。

优选的，步骤S200中所述叠前深度偏移为高斯束叠前深度偏移。

优选的，步骤S300中所述地层构造属性包括地层倾角信息和地层位置信息。

优选的，步骤S300中所述的自动提取的方法为采用结构张量算法处理偏移剖面，提取偏移剖面中的地层倾角信息。

本实施例中的结构张量法拾取地下地质倾角的原理及方法如下：

本实施例采用结构张量算法提取偏移图像的局部倾角信息。设A为二维地震图像，二维图像A中表示空间方向信息的结构张量由图像梯度值定义，结构张量表示区域的变化方向和沿变化方向的变化量大小，地震地层纹理和断层纹理由局部各点方位信息变化关系确定。引入Gaussian函数模糊局部细节，使得结构张量突出显示区域内信号的复杂性。对二维图像，结构张量是一个2*2的矩阵：

其中g_x与g_y代表地震图像沿水平方向和垂直方向的梯度，<·>代表二维高斯光滑滤波。

针对半正定矩阵G，特征值与特征向量可由求解|G-λI|＝0得到：

λ₁：最大特征值，张量能量在第一个特征张量方向v₁的能量，

λ₂：最小特征值，张量能量在第二个特征张量方向v₂的能量，

(λ₁-λ₂)/λ₁：线性度，反映局部方向的一致性。

特征向量描述了图像局部线性结构的方向性，针对图像的每个点，特征向量v₁正交于图像的主结构方向，特征向量v₂平行于图像的主结构方向，该点的局部倾角方向即为向量v₂的方向。结构张量算子G理想地包含了地下地质结构的局部构造特征。可用作下一步构造层析预条件算子。

本实施例中构建构造约束层析方程组方法如下：

常规层析方程组可表达如下：

L△m＝τ (3)

其中L是射线层析核函数，在地震层析反演领域多有介绍，本发明不就此展开。△m是速度模型更新量，τ是正演计算地震波传播走时与接收数据走时的差。

考虑模型预条件，即△m＝Su，则构造约束层析方程组可以表示为：

S^TL^TLSu＝S^TL^Tτ (4)

式中，预条件算子S是包含地质构造信息的光滑算子，该方程构建的即为地质构造约束预条件的层析方程，其对应的解为预条件后的光滑解。

所以，将地质构造信息加入构建的光滑矩阵中是地质构造约束预条件的关键点。地下介质在模型参数化后，其基本的地质规律没有改变，所以参数之间必然存在着一定的联系。层析中数据测量的精度由偏移剖面上的反射面倾角和在共成像到集中拾取成像深度的精度决定。所以，地质构造中反射面的倾角信息和散射点在深度上的分布规律为地质意义上的平滑提供了一个可行的方式，而且不依赖于先验信息。

本发明所选取的预条件算子S表达如下：

S＝(I)+D^TGD^-1 (5)

其中，I是单位阵，D是梯度算子

D^T是D的转置矩阵

G是结构张量算子。

本发明利用最小二乘QR方法(LSQR)求解矩阵方程组(4)，该方法是一种迭代的方法，可以在最小二乘意义下高效地求解大规模稀疏矩阵。

实施例二

本实施例通过多层背斜理论模型测试检验地质构造约束的可行性和有效性。首先选取理论模型，如图2为理论速度模型，该模型背斜发育，背斜顶点处具有反射层数较多、纵向速度变化较大的特点。该模型的横向和纵向网格数分别为1201和601，横向和纵向网格间距分别为10m和5m。图3为平滑后的理论速度模型，结合图3与高斯束正演得到反射数据，即如图4所示的正演地震记录剖面图，将图4视为观测数据。

首先利用图4所示的观测数据在初始速度模型中反传得到如图5所示的始速度模型剖面图，进行高斯束叠前深度偏移，得到如图6所示的初始速度叠前深度偏移剖面剖面图。图6中的初始速度模型视浅层为海水，即速度已知为常数，第一反射界面作为海底面并且深度已知，背斜部分的速度由浅到深按常梯度变化。然后，在初始偏移剖面图6中利用自动拾取技术提取地层倾角和位置信息，图7是在初始偏移剖面上提取的倾角场，作为地质构造约束的主要依据。输入观测数据、初始速度模型以及地层倾角信息，构建高斯光滑矩阵S和线性化矩阵L，建立预条件层析方程，反演求解速度更新量，更新速度模型。最后，更新迭代50次，最终更新后获得如图8所示的新的速度模型。

比较例

其他步骤与实施例二相同，在本比较例中，不加入地质构造约束，如图9所示，为本比较例的层析反演迭代50次更新后的速度模型。

图10是初始速度模型、正确速度模型、本发明实施例二地质构造约束层析反演速度模型和比较例无地质构造约束层析反演速度模型在CIP＝488(其中CIP是指commonimage point，即共成像点)处的偏移后生成的角道集的对比图；

其中，图10a为初始速度模型偏移后的角道集；图10b为正确速度模型偏移后的角道集；图10c为本发明实例二中的地质构造约束层析反演速度模型偏移后的角道集；图10d为本发明比较例中无地质构造约束层析反演速度模型偏移后的角道集。

对比图10c和图10d，可以看出两者同样都是更新50次，图10c和图10d分别与图7比较，二者在深部位置的速度，非常明显；在无地质构造约束时深部速度更新比较慢，与正确速度相差较大。

对比结果发现：初始速度偏移道集存在上翘的现象(图10a)，说明初始速度偏小；而经过本发明地质构造约束的层析反演的速度模型经偏移后，角道集被拉平(图10c)；且相较于无地质构造约束层析反演的速度模型经偏移后的结果(图10d)，更加趋近于真实速度模型偏移后的角道集(图10b)，这表明本发明的地质构造约束层析反演方法能有效提高层析反演的精度，且通过实施例看出本发明整体自动化完成，不需要额外的人工干预，是一种自动化程度很高的构造约束层析反演方法。

虽然已经参考优选实施例对本发明进行了描述，但在不脱离本发明的范围的情况下，可以对其进行各种改进并且可以用等效物替换其中的部件。尤其是，只要不存在结构冲突，各个实施例中所提到的各项技术特征均可以任意方式组合起来。本发明并不局限于文中公开的特定实施例，而是包括落入权利要求的范围内的所有技术方案。

Claims (6)

1.一种自动化的地质构造约束层析反演方法，包括以下步骤：

S100：首先提取观测数据；

S200：利用观测数据，在初始速度模型中反传，进行叠前深度偏移，得到初始的偏移剖面；

S300：从初始的偏移剖面中自动提取地层构造属性；

S400：根据观测数据、初始速度模型以及地层构造属性，利用射线追踪计算层析敏感核函数，构建构造约束层析方程组；

S500：迭代反演求解速度更新量，更新速度模型，最终更新后获得最终的速度模型；

步骤S400中，确定构造约束层析方程组所依据的公式为：

S^TL^TLSu＝S^TL^Tτ

式中，L为射线层析核函数，S为预条件算子，τ是正演计算地震波传播走时与接收数据走时的差，u为层析反演的速度模型更新量；

预条件算子S的表达式如下：

S＝(I+D^TGD)^-1

其中，I是单位阵，

D是梯度算子

D^T是D的转置矩阵

G是结构张量算子。

2.根据权利要求1所述的自动化的地质构造约束层析反演方法，其特征在于，步骤S200中所述叠前深度偏移为高斯束叠前深度偏移。

3.根据权利要求1所述的自动化的地质构造约束层析反演方法，其特征在于，步骤S300中所述地层构造属性包括地层倾角信息和地层位置信息。

4.根据权利要求1所述的自动化的地质构造约束层析反演方法，其特征在于，步骤S300中所述的自动提取的方法为采用结构张量算法处理偏移剖面，提取偏移剖面中的地层倾角信息。

5.根据权利要求4所述的自动化的地质构造约束层析反演方法，其特征在于，所述结构张量算法如下：

其中g_x为地震图像沿水平方向的梯度，

g_y为地震图像沿垂直方向的梯度，

＜·＞为二维高斯光滑滤波，

G为结构张量算子。

6.根据权利要求5所述的自动化的地质构造约束层析反演方法，其特征在于，采用最小二乘QR方法(LSQR)求解构造约束层析方程组S^TL^TLSu＝S^TL^Tτ。

Images