CN111239833B - 一种基于泊松比分解的k值鲁棒ypd叠前同时反演方法 - Google Patents

Abstract

本发明涉及一种基于泊松比分解的K值鲁棒YPD叠前同时反演方法，包括：一、测井曲线地震尺度匹配处理的步骤；二、实际叠前道集时深标定的步骤；三、实际叠前道集保AVO处理的步骤：四、K值扰动AVA鲁棒性分析的步骤；五、基于泊松比分解的K值鲁棒YPD叠前同时反演的步骤。本发明提供了一种适用于纵、横波速度变化显著的复杂地质条件的地下岩石杨氏模量、泊松比和密度的反演方法。相比Zong‑YPD反射系数近似反演，本发明通过泊松比分解降低了反射近似方程对K值的依赖度，构建了一套K值鲁棒的YPD直接反演方法及流程，有效提高了杨氏模量、泊松比和密度参数的反演精度。

Description

一种基于泊松比分解的K值鲁棒YPD叠前同时反演方法

技术领域

本发明涉及一种基于泊松比分解的K值鲁棒YPD叠前同时反演方法，属于勘探地球物理领域。

背景技术

YPD(Young’s modulus-Poisson’ratio-Density，简称YPD)叠前同时反演是指利用叠前地震、测井、地质信息直接反演杨氏模量、泊松比和密度。杨氏模量和泊松比用于表征介质的力学性质，是桥梁、隧道建设和致密油气储层压裂开发等应用中的关键工程参数。因此，研发预测地下介质杨氏模量和泊松比的技术具有显著现实意义。考虑方法效率和应用性，利用Zoeppritz方程各类近似的线性AVA反演是目前主流的叠前地震反演技术。其中，基于Fatti近似的纵、横波阻抗和密度反演方法工业应用范围最广，一般使用反演的纵、横波阻抗和密度间接计算杨氏模量和泊松比参数。

间接计算杨氏模量和泊松比包含乘、除与平方运算，放大了反演结果误差，增大了结果不确定度。为此，宗兆云等(杨氏模量和泊松比反射系数近似方程及叠前地震反演.地球物理学报,2013.55(11):3786-3794.)推导了由杨氏模量、泊松比和密度表示的YPD近似方程，实现了杨氏模量和泊松比的直接反演。基于精确Zoeppritz方程，通过杨氏模量、泊松比系数和密度替换Zoeppritz方程中的纵、横波速度和密度，Zhou等(Prestack amplitudeversus angle inversion for Young's modulus and Poisson's ratio based on theexact Zoeppritz equations.Geophysical Prospecting,2017,65(6):1462-1476.)，通过迭代重加权最小平方算法反演求取杨氏模量、泊松比和密度。此外，宋建国等(基于比值均方根的杨氏模量反演方法.地球物理学报,2018.61(4):1508-1518.)推导了杨氏模量、泊松比系数和密度的比值均方根的一阶泰勒，由广义线性间接反演得到杨氏模量和泊松比。

但是，上述现有技术中至少存在如下问题：(1)精确Zoeppritz方程在反演过程中，由于方程的非线性引起正、反演运算复杂，要对参数进行一阶近似，通过非线性或迭代线性化反演方法求取，运算量和效率限制了工业应用推广；(2)更重要的，无论是近似方程还是精确方程，都需要预知背景横纵波速度比的平方(即K值)，在具体反演过程中将K设定为常值。然而，对于非均质较强的地层，K值变化将引起YPD同时反演的不稳定及偏差。

发明内容

针对上述问题，本发明的目的是提供一种基于泊松比分解的K值鲁棒YPD叠前同时反演方法，该方法通过泊松比分解实现YPD反射系数的新参数化近似方程，最大限度降低反射近似方程中K值的影响，提高反演过程对K值的容错度，实现鲁棒性的YPD参数同时反演。

为实现上述目的，本发明采取以下技术方案：一种基于泊松比分解的K值鲁棒YPD叠前同时反演方法，包括：

一、测井曲线地震尺度匹配处理的步骤，包括：

①针对井上目的层段的测井数据进行预处理，具体为对测井曲线进行环境校正，剔除泥岩或含煤地层井眼垮塌引起的测量异常值，其中测井曲线包括地层密度曲线、纵波时差曲线和横波时差曲线；②当缺少横波测井资料时，使用环境校正后的测井曲线进行预测横波速度；③使用Backus平均及滑动窗函数法地层块化平均，进行测井曲线尺度向地震尺度的匹配处理，确保Backus处理后的地层密度、声波资料能够有效反映地震尺度响应，且消除K值中的高频异常对后续计算分析的影响；④在尺度匹配处理后测井曲线的数据基础上，通过公式(1)-(3)分别求取出杨氏模量、泊松比和K值，并生成相应的测井曲线，用于后续分析和反演：

式中，E为杨氏模量；υ为泊松比；ρ为地层密度；V_S为横波速度；V_P为纵波速度。

二、实际叠前道集时深标定的步骤：

结合校验炮或零偏VSP测井的时深关系，进行叠后井震标定，以此时深关系对实际叠前道集进行井震标定，以避免二类AVO引起叠后剖面中“暗点”的漏标，并确定是否需要进行剩余时差校正：当实际叠前道集存在显著未拉平或抖动时，则需要进行剩余时差校正。

三、实际叠前道集保AVO处理的步骤：

以井合成叠前道集(根据声波测井曲线、密度曲线、地震子波和褶积得到井合成叠前道集，属于现有技术)及井上目的层段的AVO信息为标准，对实际叠前道集进行剩余时差校正、随机噪声压制以及残余多次波压制。由此，在保护叠前道集AVO品质的基础上，可以提高数据品质。

四、K值扰动AVA鲁棒性分析的步骤，包括：

首先，选取井上目的层段的弹性参数构建双层介质模型，模拟对比Aki-Richards近似方程、Zong-YPD近似方程和基于泊松比分解的YPD近似方程的AVA曲线；然后，扰动K值模拟AVA曲线的变化，以此分析各个近似方程对K值的确定度和鲁棒性，具体包括以下：

①与纵、横波速度和地层密度变化相关的Aki-Richard反射系数近似方程为：

其中，θ为纵波入射角；纵波速度反射系数项为

横波速度反射系数项为

地层密度反射系数项为

ΔV_P为相邻地层的纵波速度的变化量；ΔV_S为相邻地层的横波速度的变化量；Δρ为相邻地层的密度的变化量。

②与杨氏模量、泊松比和地层密度有关的Zong-YPD反射系数近似方程为：

其中，杨氏模量反射系数项为

泊松比反射系数项为

地层密度反射系数项为

ΔE为相邻地层的杨氏模量的变化量；Δυ为相邻地层的泊松比的变化量。

③将泊松比υ分解为两个新的泊松比变量：υ₁＝1-υ，

④纵波速度V_P由杨氏模量E和泊松比υ表示为

代入分解后的υ₁和υ₂则表示为：

⑤同样地，横波速度V_S由杨氏模量E和泊松比υ表示为

代入分解后的泊松比υ₁和υ₂则表示为：

⑥将公式(6)和(7)代入公式(4)，得到基于泊松比分解的YPD反射系数近似方程为：

其中，杨氏模量反射系数项为

泊松比υ₁的反射系数项为

泊松比υ₂的反射系数项为

地层密度反射系数项为

⑦根据公式(4)、(5)和(8)，在前三步骤的数据基础上，即可进行AVO模拟，给定不同K值进行AVA扰动分析，作为后续YPD参数反演稳定性质控。

五、基于泊松比分解的K值鲁棒YPD叠前同时反演的步骤，该反演过程在于求解以下目标函数：

其中，d为实际叠前道集数据；待反演参数向量为

R_i为依次对应E、υ₁、υ₂和ρ的反射系数；ξ_i为依次对应E、υ₁、υ₂和ρ的背景模型；λ_i为依次对应E、υ₁、υ₂和ρ的背景模型约束权值；积分算子L和正演映射算子G分别为：

其中，C_i(θ_j)为待反演参数的角度系数；θ_j为实际叠前道集的角度，j＝1,2,3,…,M；W(θ_j)为各个角度地震子波矩阵；n为采样点个数；

对于实际叠前道集的角度θ_j：

求取目标函数(9)的最优反演解为：

由公式(12)得到E、υ₁、υ₂和ρ的反射系数，通过积分即可得到各个参数的绝对值：

由于υ₁＝1-υ，最后仅需要做一步代换回得到泊松比υ＝1-υ₁，此时即完成杨氏模量、泊松比和地层密度的反演。

由于只引入了一步与常值差值运算，并未引入额外数据，也避免了乘、除运算对噪声的增益，这也是用未选用υ₂的原因。

本发明由于采取以上技术方案，其具有以下优点：本发明提供了一种适用于纵、横波速度变化显著的复杂地质条件的地下岩石杨氏模量、泊松比和密度的反演方法。相比Zong-YPD反射系数近似反演，本发明通过泊松比分解降低了反射近似方程对K值的依赖度，构建了一套K值鲁棒的YPD直接反演方法及流程，有效提高了杨氏模量、泊松比和密度参数的反演精度。

附图说明

图1是本发明的流程示意图；

图2是本发明实施例中测井曲线地震尺度和井旁叠前道集处理结果对比图；

图3是本发明实施例中目标层段近似反射系数对K值扰动的鲁棒性分析图；

图4是本发明实施例中基于Zong-YPD近似方程的反演结果图；

图5是本发明实施例中提出基于泊松比分解YPD近似方程的反演结果图；

图6是本发明实施例中叠前地震道集处理前和处理后结果图；

图7是本发明实施例中全叠加地震数据图；

图8是本发明实施例中杨氏模量数据背景模型图；

图9是本发明实施例中泊松比数据背景模型图；

图10是本发明实施例中基于泊松比分解YPD近似方程的杨氏模量反演结果图；

图11是本发明实施例中基于泊松比分解YPD近似方程的泊松比反演结果图。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面对本发明中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

如图1所示，本实施例提供的基于泊松比分解的K值鲁棒性YPD叠前同时反演方法，包括以下步骤：

步骤1)测井曲线地震尺度匹配处理：

①针对井上目的层段的测井数据进行预处理，具体为对测井曲线进行环境校正，剔除泥岩或含煤地层井眼垮塌引起的测量异常值，其中测井曲线包括地层密度曲线、纵波时差曲线和横波时差曲线。

②当缺少横波测井资料，需使用环境校正后的测井曲线进行预测横波速度。

③通过使用Backus平均及滑动窗函数法地层块化平均，进行测井曲线尺度向地震尺度的匹配处理，匹配处理前、后结果如图2所示，可以看到处理前测井曲线(灰色实线)中小尺度信息远远超过地震可分辨尺度范围，且K值变化剧烈。这种K值突变会影响后续敏感性分析和反演，匹配处理后的测井曲线保留了有效中、小尺度信息，且处理后测井合成道集与井旁道集吻合较好。

④在尺度匹配处理后测井曲线的数据基础上，计算杨氏模量、泊松比和K值，并生成相应的测井曲线，用于后续分析和反演。

步骤2)实际叠前道集时深标定：

结合校验炮或零偏VSP测井的时深关系，进行叠后井震标定，以此时深关系对实际叠前道集进行井震标定，以避免二类AVO引起叠后剖面中“暗点”的漏标；图2中三个道集依次为尺度匹配处理后井曲线合成的叠前角度域道集、原始井旁叠前道集和处理后的叠前道集，通过与合成道集对比，原始叠前道集需要做随机噪声压制、剩余时差校正和残余多次波压制。

步骤3)实际叠前道集保AVA处理：

在前两步数据基础上，需要对实际叠前角道集进行剩余时差校正、随机噪声压制、残余多次波压制，得到图2中处理后道集，可以看到剩余时差和噪声得到了有效压制。

步骤4)K值扰动AVO鲁棒性分析：

选取图2中井上目的层段储层和非储层的弹性参数构建双层介质模型，求解单界面反射近似方程：

和

得到

等总反射系数和

R_E(θ)等参数的反射系数项数值；再扰动K值重新计算，结果如图3所示；图3中黑色点线为Zoeppritz方程计算反射系数，虚线为基于泊松比分解近似方程中的密度反射系数项，曲线颜色表示K值大小；对比看出在不同K值条件下，基于泊松比分解的YPD近似方程和Aki-Richards近似方程总反射系精度相当；Zong-YPD近似方程得到的总反射系数随K值变化较大，这主要是由于泊松比反射系数对K值的依赖性引起的，而本发明给出的基于泊松比分解的YPD近似方程中υ₁不受K值影响，即使K值有误的情况下也不受影响，如同Aki-Richards近似方程中的P波速度反射系数项；本发明基于泊松比分解的YPD近似方程和Zong-YPD近似方程的杨氏模量反射系数项受K值的影响相当。本步骤作为本发明中质控环节，分析目的层反射系数方程对K值扰动的鲁棒性，判断是否适用于叠前AVA反演杨氏模量、泊松比和密度。

步骤5)基于泊松比分解的K值鲁棒YPD叠前同时反演，该反演过程在于求解以下目标函数：

在实际实施过程中，首先提取地震子波，并构建低频背景模型；然后对井旁道集进行反演质控分析，图4和图5分别给出了图2中叠前道集的Zong-YPD反演和本发明基于泊松比分解的K值鲁棒YPD叠前同时反演结果，黑色实线为反演结果，灰色细线为测井曲线，灰色粗线为反演所用背景模型。在相同模型和反演参数的条件下，两种反演方法都能得到较为合理的杨氏模量、泊松比和密度结果。然而，对于箭头所指的K值变化剧烈的位置，本发明反演结果与井更吻合，特别是泊松比的反演结果精度和井吻合度改善显著。结合步骤4)和井旁道的反演质控，说明本发明的适用性，确定反演参数。最后，将地震子波和反演参数应用到整个工区。图6为工区部分叠前道集(处理前、后结果)，与图2井旁道集处理一致，处理后的道集中随机噪声和参与多次波得到压制，剩余时差校正后同相轴得到拉平。图7为工区过井线叠后地震剖面，黑色曲线为纵横波速度比(K值的倒数)，虚线为井轨迹。可以看到，纵横波速度比变化较为剧烈，在叠后地震剖面上难以识别弹性参数的改变。图8和图9分别为步骤1)匹配处理后井曲线和解释层位构建的杨氏模量和泊松比低频背景模型，黑色曲线分别为杨氏模量和泊松比，低频背景模型只给出了大套地层参数变化背景趋势，并不反映测井曲线中的高频变化。图10和图11分别为本发明反演得到的杨氏模量和泊松比结果，可以看到反演结果有效刻画出了地层杨氏模量和泊松比的细节变化，取得了良好的效果。

最后应说明的是：以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和范围。

Claims (5)

1.一种基于泊松比分解的K值鲁棒YPD叠前同时反演方法，其特征在于，包括：

一、测井曲线地震尺度匹配处理的步骤；

二、实际叠前道集时深标定的步骤；

三、实际叠前道集保AVO处理的步骤：

四、K值扰动AVA鲁棒性分析的步骤，具体包括：

首先，选取井上目的层段的弹性参数构建双层介质模型，模拟对比Aki-Richards近似方程、Zong-YPD近似方程和基于泊松比分解的YPD近似方程的AVA曲线；

然后，扰动K值模拟AVA曲线的变化，以此分析各个近似方程对K值的确定度和鲁棒性，具体包括以下步骤：

①与纵、横波速度和地层密度变化相关的Aki-Richard反射系数近似方程为：

其中，θ为纵波入射角；纵波速度反射系数项为

横波速度反射系数项为

密度反射系数项为

为相邻地层的纵波速度的变化量；△V_S为相邻地层的横波速度的变化量；△ρ为相邻地层的密度的变化量；

为相邻地层的纵波速度的平均值；

为相邻地层的横波速度的平均值；

为相邻地层的地层密度的平均值；K为地层横纵波速度比的平方；

②与井上杨氏模量、泊松比和地层密度有关的Zong-YPD反射系数近似方程为：

其中，井上杨氏模量反射系数项为

泊松比反射系数项为

地层密度反射系数项为

△E为相邻地层的杨氏模量的变化量；△υ为相邻地层的泊松比的变化量；

为相邻地层的杨氏模量的平均值；

为相邻地层的泊松比的平均值；

③将泊松比υ分解为两个新的泊松比变量：υ₁＝1-υ，

④纵波速度V_P由井上杨氏模量E和泊松比υ表示为

代入分解后的υ₁和υ₂则表示为：

其中，ρ为地层密度；

⑤同样地，横波速度V_S由井上杨氏模量E和泊松比υ表示为

代入分解后的泊松比υ₁和υ₂则表示为：

⑥将公式(6)和(7)代入公式(4)，得到基于泊松比分解的YPD反射系数近似方程为：

其中，杨氏模量反射系数项为

泊松比υ₁的反射系数项为

泊松比υ₂的反射系数项为

地层密度反射系数项为

为相邻地层的分解后的泊松比υ₁的平均值；

相邻地层的分解后的泊松比υ₂的平均值；

⑦根据公式(4)、(5)和(8)，在前三步骤的数据基础上，即可进行AVO模拟，给定不同K值进行AVA扰动分析，作为后续YPD参数反演稳定性质控；

五、基于泊松比分解的K值鲁棒YPD叠前同时反演的步骤。

2.根据权利要求1所述的基于泊松比分解的K值鲁棒YPD叠前同时反演方法，其特征在于，所述测井曲线地震尺度匹配处理的步骤具体包括：

①针对井上目的层段的测井数据进行预处理，具体为对测井曲线进行环境校正，剔除泥岩或含煤地层井眼垮塌引起的测量异常值，其中测井曲线包括地层密度曲线、纵波时差曲线和横波时差曲线；

②当缺少横波测井资料时，使用环境校正后的测井曲线进行预测横波速度；

③使用Backus平均及滑动窗函数法地层块化平均，进行测井曲线尺度向地震尺度的匹配处理，确保Backus处理后的地层密度、声波资料能够有效反映地震尺度响应，且消除K值中的高频异常对后续计算分析的影响；

④在尺度匹配处理后测井曲线的数据基础上，通过公式(1)-(3)分别求取出井上杨氏模量、泊松比和K值曲线，并生成相应的测井曲线，用于后续分析和反演：

式中，E为井上杨氏模量；υ为泊松比；ρ为地层密度；V_S为横波速度；V_P为纵波速度。

3.根据权利要求1所述的基于泊松比分解的K值鲁棒YPD叠前同时反演方法，其特征在于，所述实际叠前道集时深标定的步骤为：

结合校验炮或零偏VSP测井的时深关系，进行叠后井震标定，以此时深关系对实际叠前道集进行井震标定，以避免二类AVO引起叠后剖面中“暗点”的漏标，并确定是否需要进行剩余时差校正：当叠前道集存在显著未拉平或抖动时，则需要进行剩余时差校正。

4.根据权利要求1所述的基于泊松比分解的K值鲁棒YPD叠前同时反演方法，其特征在于，所述实际叠前道集保AVO处理的步骤为：

以井合成叠前道集及井上目的层段的AVO信息为标准，对实际叠前道集进行剩余时差校正、随机噪声压制以及残余多次波压制。

5.根据权利要求1所述的基于泊松比分解的K值鲁棒YPD叠前同时反演方法，其特征在于，所述基于泊松比分解的K值鲁棒YPD叠前同时反演的步骤在于求解以下目标函数：

其中，d为实际叠前道集数据；待反演参数向量为

R_i为依次对应E、υ₁、υ₂和ρ的反射系数；ξ_i为依次对应E、υ₁、υ₂和ρ的背景模型；λ_i为依次对应E、υ₁、υ₂和ρ的背景模型约束权值；积分算子L和正演映射算子G分别为：

其中，C_i(θ_j)为待反演参数的角度系数；θ_j为实际叠前道集的角度，j＝1,2,3,…,M；W(θ_j)为各个角度地震子波矩阵；n为采样点个数；

对于实际叠前道集的角度θ_j：

求取目标函数(9)的最优反演解为：

由公式(12)得到E、υ₁、υ₂和ρ的反射系数，通过积分得到各个参数的绝对值：

其中，r为反演的反射系数向量；由于υ₁＝1-υ，最后仅需要做一步代换回得到泊松比υ＝1-υ₁，此时即完成杨氏模量、泊松比和地层密度的反演。

Images