一种基于透射方程的各向异性参数反演方法

摘要

本发明提供了一种基于透射方程的各向异性参数反演方法，属于地球物理反演领域。所述方法包括：(1)选取任意三个不同方位角对应的叠后地震数据体，通过稀疏脉冲的反演方法计算得到对应该三个不同方位角的反射系数(2)分别计算三个不同方位角的反射系数对应的透射系数(3)分析叠前方位道集获得入射角信息，所述入射角信息包括最大入射角i2、最小入射角i1和间隔；(4)利用方位角信息透射系数信息计算得到各向异性参数Δε(V)和Δδ(V)。

说明书

技术领域

本发明属于地球物理反演领域，具体涉及一种基于透射方程的各向异性参 数反演方法。

背景技术

各向异性参数在裂缝定量预测中起着举足轻重的作用，而由于各向异性参 数数量级远远小于地震数据数量级，为其计算带来了极大的困难。目前计算各 向异性参数主要是基于叠前方位道集，主要分为以下三类渠道：1.根据随地震 方位变化的速度信息，速度分析精度难以满足方法需要；2.根据随地震方位变 化的剩余时差信息，拾取剩余时差时存在较多误差；3.以简化Ruger(1998)方 程为基础(未经简化或者近似的Ruger方程的解是计算不出来的，必须要通过 简化方程或者一些近似关系将Ruger方程简化，才能得到Ruger方程的解。并 且现有方法中使用的Ruger方程都是Ruger反射系数方程，而不是Ruger透射 方程)，基于各种假设条件的弹性参数反演，虽然研究表明简化后的Ruger方程 对计算反射系数影响不大，但是由此产生的传递误差对于计算数量级极小的各 向异性参数影响还是十分巨大的。

发明内容

本发明的目的在于解决上述现有技术中存在的难题，提供一种基于透射方 程的各向异性参数反演方法，基于完整的Ruger透射系数方程，以叠后数据体 方位之间的差异为基础，得到各向异性参数Δδ^(V)，Δε^(V)的解。

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

一种基于透射方程的各向异性参数反演方法，包括：

(1)选取任意三个不同方位角对应的叠后地震数据体，通过稀疏脉 冲的反演方法计算得到对应该三个不同方位角的反射系数

(2)分别计算三个不同方位角的反射系数对应的透射 系数

(3)分析叠前方位道集获得入射角信息，所述入射角信息包括最大入射角i₂、 最小入射角i₁和间隔；

(4)利用方位角信息透射系数信息计算得 到各向异性参数Δε^(V)和Δδ^(V)。

所述步骤(2)是这样实现的：

分别对三个不同方位角的反射系数采用公式(11) 得到对应的透射系数

R_PP+T_p＝1(11)。

所述步骤(4)是这样实现的：

利用公式(6)和(7)计算得到各向异性参数Δε^(V)和Δδ^(V)：

其中，$B={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{2}i,$$C={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{2}i{\tan }^{2}i,$$D={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{4}i,$$E=$${\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{2}i{\tan }^{2}i+{\sin }^{4}i.$

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

(1)采用的是完整的Ruger透射系数方程，未对方程进行简化，理论上不存 在误差；

(2)具有叠后地震反演稳定性高、信噪比高的优点；

(3)输入数据范围可控，实际应用中易于实现；

(4)计算效率高，可直接获得各向异性参数Δε^(V)和Δδ^(V)。

附图说明

图1本发明方法的步骤框图

图2a模型数据中的Vp图

图2b模型数据中的Vs图

图2c模型数据中的Den图

图2d模型数据中的δ图

图2e模型数据中的ε图

图2f模型数据中的γ图

图3a正演得到的0°方位的叠前数据

图3b演得到的30°方位的叠前数据

图3c正演得到的60°方位的叠前数据

图4a反演得到的0°方位对应的反射系数

图4b反演得到的30°方位对应的反射系数

图4c反演得到的60°方位对应的反射系数

图5a计算得到的0°方位对应的透射系数

图5b计算得到的30°方位对应的透射系数

图5c计算得到的60°方位对应的透射系数

图6a最终结果各向异性参数Δε^(V)

图6b最终结果各向异性参数Δδ^(V)。

具体实施方式

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

目前绝大多数的各向异性参数反演是基于简化的Ruger反射系数方程，其 近似或者简化的过程引入了误差，而Ruger透射方程与反射系数方程具有相同 的各向异性特征。因此本发明基于完整的Ruger透射系数方程得到各向异性参 数Δδ^(V)，Δε^(V)的解，从理论上不引入误差，计算结果精度更高，并且具有叠后反 演的优点，反演结果可用于裂缝预测和页岩气预测等热点领域。

本发明以完整的Ruger透射方程为基础，计算出各向异性参数ε^(V)，δ^(V)，方 法适用于HTI(横向各向异性介质)。采用的Ruger透射方程与反射系数方程具 有相同的各向异性特征，并且未对方程进行简化，理论上不存在误差，计算结 果精度更高，并且具有叠后反演的优点，反演结果可用于裂缝预测和页岩气预 测等热点领域

Ruger(2002)借助弱各向异性的概念，推导出HTI介质的纵波近似透射系 数关系式：

式中：和分别为HTI介质上下两层的纵波、横波速度平均值以及密度平均 值；Δδ^(V)，Δε^(V)分别为上下两层各向异性参数差值；i和分别代表入射角和方位 角。

对于方程1按照入射角i方向进行叠加，方程变为：

其中：$A={\Sigma }_{i_1}^{i_2}1-\frac{1}{2}\frac{\mathrm{\Delta Z}}{\bar{Z}}+\frac{1}{2}\frac{\mathrm{\Delta \alpha }}{\bar{\alpha }}{\sin }^{2}i+\frac{1}{2}\frac{\mathrm{\Delta \alpha }}{\bar{\alpha }}{\sin }^{2}i{\tan }^{2}i,$$B={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{2}i,$$C={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{2}i{\tan }^{2}i,$$D={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{4}i,$$E={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{2}i{\tan }^{2}i+{\sin }^{4}i$i和分别代表入射角和方位角。

分别令

解由方程(3)(4)(5)组成的方程组，得

其中$B={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{2}i,$$C={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{2}i{\tan }^{2}i,$$D={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{4}i,$$E={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{2}i{\tan }^{2}i+{\sin }^{4}i,$i和分别代表入射角和方位角，分别代表时的透射系数，上述为HTI介质中近似各向异性参 数计算公式，在此基础上进行各向异性参数Δε^(V)和Δδ^(V)的反演研究。

反射系数通过稀疏脉冲反演计算，常规的稀疏脉 冲流程可分为三步：

(1)反射系数反演

采用最大似然反褶积进行反射系数的反演，最大似然反褶积对地层的假设 认为：地层的反射系数是由较大的反射界面的反射和具有高斯背景的小反射叠 加组合而成，导出一个最小目标函数：

$J=\frac{1}{R^2}\underset{K=1}{\overset{L}{\Sigma }}{r}^2\left(K\right)+\frac{1}{N^2}\underset{K=1}{\overset{L}{\Sigma }}{n}^2\left(K\right)-2M\ln \left(\lambda \right)-2(L-M)\ln (1-\lambda )---\left(8\right)$

式中，R2和N2分别为反射系数和噪音的均方值，r(K)和n(K)表示第K个 采样点的反射系数和噪音，M表示反射层数，L表示采样总数，λ表示给定反射 系数的似然值。通过多次迭代，求取反射系数。

(2)根据反射系数的反演结果结合阻抗趋势计算一个初始的波阻抗

根据最大似然反褶积计算得到的反射系数，结合初始阻抗模型，采用递推 算法，反演得到初始的波阻抗模型：

$Z\left(i\right)=Z(i-1)\frac{1+R\left(i\right)}{1-R\left(i\right)}---\left(9\right)$

式中，Z(i)为第i层的波阻抗值，R(i)为第i层的反射系数。

(3)结合井的约束条件进行波阻抗反演

约束稀疏脉冲反演对每一道依据目标函数对计算出的初始波阻抗进行调整， 包括对反射系数的调整(即根据初始波阻抗(这时有一个对应初始波阻抗的反 射系数)对实际波阻抗进行逼近，使用过迭代寻优的方法，最终得到一个近似 于实际波阻抗的结果，这时对应的反射系数就是所需要的反射系数)。目标优化 函数为：

F＝L_p(r)+λL_q(s-d)+α^-1L₁ΔZ(10)

式中，r为反射系数序列，Δz为与阻抗趋势的差序列，d为地震道序列，s为 合成地震道序列，λ为残差权重因子，α为趋势权重因子，p、q为L模因子。 具体的，右式第一项反映了反射系数的绝对值和，第二项反映了合成声波记录 与原始地震数据的差值，第三项为趋势约束项。

本方法需要的就是(3)中调整后的反射系数

根据zeoppritz方程，在垂直入射情况下，反射系数与透射系数存在如下 关系：

R_PP+T_p＝1(11)

而本方法反演是以叠后数据为基础得到的反射系数R_PP，因此是近似地符合垂直 入射条件，可以通过方程(11)计算得到透射系数T_p。

本发明基于透射方程的各向异性参数反演方法的实现步骤为：

(1)选取任意三个不同方位角对应的叠后地震数据体，通过稀疏脉 冲的反演方法计算得到对应该三个不同方位角的反射系数

(2)通过方程(11)分别计算对应的透射系数

(3)分析叠前方位道集获得入射角信息，所述入射角信息包括最大入射角i₂、 最小入射角i₁和间隔，然后通过公式$B={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{2}i,$$C={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{2}i{\tan }^{2}i,$$D={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{4}i,$$E={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{2}i{\tan }^{2}i+{\sin }^{4}i$计算出B、C、D、E的值；

(4)利用方位角信息透射系数信息和参数B、 C、D、E带入方程(6)(7)中，计算得到各向异性参数Δε^(V)和Δδ^(V)。

图1示出本发明一种基于透射方程的各向异性参数反演方法的具体流程。 下面参照附图并结合实施例对本发明做进一步说明。图2a-图2f为使用的各向 异性模型数据。图3a-图3c为正演得到的0°、30°、60°方位的叠前数据。

首先根据步骤(1)对0°，30°，60°方位的叠后数据体分别进行稀疏脉 冲反演，计算得到对应的反射系数体R_PP(0°)，R_PP(30°)，R_PP(60°)，如图4a- 图4c所示，然后根据方程(11)分别计算对应的透射系数如 图5a-图5c所示。

然后根据步骤(3)分析叠前方位道集的入射角信息，本次叠前地震数据沿 入射角叠加的范围选取为i₁＝0°，i₂＝30°，入射角间隔为1°，因此根据 $B={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{2}i,$$C={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{2}i{\tan }^{2}i,$$D={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{4}i,$$E={\Sigma }_{i_1}^{i_2}\frac{1}{2}{\sin }^{2}i{\tan }^{2}i+{\sin }^{4}i$分别计 算出B＝1.3607，C＝0.2636，D＝2.7214，E＝2.985。对于实际资料应用过程，由于 目前国内绝大多数观测系统为树状观测系统，导致从根本上讲叠前道集会存在 近、远偏移距(近远入射角)能量分布较弱，目前尚未有有效的方法可以解决 叠前道集能量问题，比较常见的方法是对叠前道集做能量均衡，但是由于缺乏 理论的依据，往往容易使原本存在的AVO特征扭曲，为之后的叠前反演带来巨 大的影响。本发明在叠前道集叠加至叠后数据处对于入射角范围i₁、i₂可根据叠 前道集质量进行分选，虽然在一定程度上造成了能量的损失，但是可以获得更 稳定的叠后数据，并且在理论上获得验证，反演效果不受影响。

接下来根据步骤(4)将已知的方位角信息透射系数和B、C、D、E的值，代入公式(6)、(7)中， 计算得到各向异性参数Δε^(V)和Δδ^(V)，如图6a、图6b所示，反演结果可以用与裂 缝定量预测。

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