CN104459777A - 基于流体体积模量avo反演的流体识别的方法及系统 - Google Patents

Abstract

本发明提供一种基于流体体积模量AVO反演的流体识别的方法及系统，所述方法包括：采集目的储层的叠前地震资料、测井资料以及地质成果资料；构建表征流体体积模量的AVO近似公式；根据所述的叠前地震资料、测井资料、地质成果资料以及所述的表征流体体积模量的AVO近似公式方程提取子波；根据所述的叠前地震资料、测井资料、地质成果资料以及子波确定流体体积模量数据体；根据所述的流体体积模量数据体识别当前储层孔隙中流体的类型。利用表征流体体积模量的AVO近似公式进行AVO反演可以估算出流体体积模量参数，避免流体识别假象，提高储层孔隙流体类型预测的可靠性。

Description

基于流体体积模量AVO反演的流体识别的方法及系统

技术领域

本发明关于石油地球物理勘探领域，特别是关于叠前地震的反演技术，具体的讲是一种基于流体体积模量AVO反演的流体识别的方法及系统。

背景技术

从20世纪70年代开始，由于地震技术数字化和计算机技术的飞速发展，利用地震波动力学问题进行储层预测和流体识别成为石油工业界的研究热点。其中，AVO(AmplitudeVersus Offset，振幅随偏移距的变化)反演技术作为一项根据振幅随偏移距的变化规律预测地下岩性及含流体性质的主流技术，其核心思想是通过研究反射系数随入射角度的变化规律提取需求的弹性参数，理论基础是描述平面波在水平分界面反射和透射的Zoeppritz方程。由于基于Zoeppritz方程的参数提取属非线性反问题，计算效率低和受噪音影响大等问题制约了其在实际生产中的应用。

为了克服由Zoeppritz方程导出的反射系数形式复杂及不易进行数值计算的困难，许多学者对Zoeppritz方程进行了简化。Koefoed将原来7个独立变量简化为5个独立变量；Bortfeld详细论述了垂直入射的平面纵波反射系数近似计算方法，并给出了区分流体和固体的简化方程；Aki和Richards在假设相邻地层介质弹性参数变化较小的情况下，给出了包含纵横波速度和密度相对变化量的纵波反射系数线性近似公式。在此基础上，许多学者对Aki-Richards方程进行了重新推导、归纳，分别提出了以不同弹性参数表示的纵波反射系数。其中，Shuey给出了突出泊松比的相对反射系数近似表达形式；Fatti等人给出了以相对波阻抗变化表示的反射系数近似公式；Mallick给出了用射线参数表示的反射系数近似形式；Goodway推导了以拉梅参数相对变化量表示的反射系数近似公式。

随着石油勘探开发的深入，地下储层的含流体识别问题成为研究的热点，迫切需要从地震资料中提取更多可以反映地下介质中的流体特征的参数。Russell等人结合饱和流体多孔弹性介质理论对Aki-Richards近似进行重组，提出了包含Gassmann流体项f的Russell反射系数近似公式，突出体现了岩石孔隙流体效应对反射系数的影响，推动了基于AVO反演的流体识别技术的发展。但是由于地下岩石的固液耦合特性，现阶段基于AVO反演方法估算的弹性参数(如波阻抗，拉梅常数，Gassmann流体项等)均会受岩石固体效应(岩石基质，孔隙度等)影响产生流体识别假象，降低了地质情况复杂地区储层含油气预测的可靠性。

因此，为了有效的提高地质情况复杂地区的储层流体表征敏感性和流体识别可靠性，现有技术中急需一种包含对孔隙流体信息更为敏感弹性参数的反射系数近似公式以及相应的AVO反演方法来识别储层流体。

发明内容

为了克服现有技术存在的基于常规AVO近似估算的岩石弹性参数受固液耦合特性影响易产生储层含流体判识假象的问题，本发明提供了一种基于流体体积模量AVO反演的流体识别的方法及系统，在反射系数近似中实现流体体积模量参数分离，流体体积模量对储层流体信息更为敏感，利用表征流体体积模量的AVO近似公式进行AVO反演可以估算出流体体积模量参数，避免流体识别假象，提高储层孔隙流体类型预测的可靠性。

本发明的目的之一是，提供一种基于流体体积模量AVO反演的流体识别的方法，包括：采集目的储层的叠前地震资料、测井资料以及地质成果资料；构建表征流体体积模量的AVO近似公式；根据所述的叠前地震资料、测井资料、地质成果资料以及所述的表征流体体积模量的AVO近似公式方程提取子波；根据所述的叠前地震资料、测井资料、地质成果资料以及子波确定流体体积模量数据体；根据所述的流体体积模量数据体识别当前储层孔隙中流体的类型。

本发明的目的之一是，提供了一种基于流体体积模量AVO反演的流体识别的系统，包括：资料采集装置，用于采集目的储层的叠前地震资料、测井资料以及地质成果资料；近似公式构建装置，用于构建表征流体体积模量的AVO近似公式；子波提取装置，用于根据所述的叠前地震资料、测井资料、地质成果资料以及所述的表征流体体积模量的AVO近似公式方程提取子波；流体体积模量数据体确定装置，用于根据所述的叠前地震资料、测井资料、地质成果资料以及子波确定流体体积模量数据体；流体类型识别装置，用于根据所述的流体体积模量数据体识别当前储层孔隙中流体的类型。

本发明的有益效果在于，提供了一种基于流体体积模量AVO反演的流体识别的方法及系统，基于常规AVO近似估算的岩石弹性参数受固液耦合特性影响易产生储层含流体判识假象，本发明基于孔隙弹性介质理论以及临界孔隙度模型，在反射系数近似中实现流体体积模量参数分离，流体体积模量对储层流体信息更为敏感，且最大限度的消除了储层岩石孔隙度等固体效应造成的流体识别假象，利用表征流体体积模量的AVO近似公式进行AVO反演可以估算出流体体积模量参数，避免流体识别假象，提高储层孔隙流体类型预测的可靠性。

为让本发明的上述和其他目的、特征和优点能更明显易懂，下文特举较佳实施例，并配合所附图式，作详细说明如下。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为本发明实施例提供的一种基于流体体积模量AVO反演的流体识别的方法的实施方式一的流程图；

图2为本发明实施例提供的一种基于流体体积模量AVO反演的流体识别的方法的实施方式二的流程图；

图3为本发明实施例提供的一种基于流体体积模量AVO反演的流体识别的方法的实施方式三的流程图；

图4为本发明实施例提供的一种基于流体体积模量AVO反演的流体识别的方法的实施方式四的流程图；

图5为图1中的步骤S103的具体流程图；

图6为图1中的步骤S104的具体流程图；

图7为本发明实施例提供的一种基于流体体积模量AVO反演的流体识别的系统的实施方式一的结构框图；

图8为本发明实施例提供的一种基于流体体积模量AVO反演的流体识别的系统的实施方式二的结构框图；

图9为本发明实施例提供的一种基于流体体积模量AVO反演的流体识别的系统的实施方式三的结构框图；

图10为本发明实施例提供的一种基于流体体积模量AVO反演的流体识别的系统的实施方式四的结构框图；

图11为本发明实施例提供的一种基于流体体积模量AVO反演的流体识别的系统中的子波提取装置300的具体结构框图；

图12为本发明实施例提供的一种基于流体体积模量AVO反演的流体识别的系统中的流体体积模量数据体确定装置400的具体结构框图；

图13为孔隙流体类型为水和气情况下，纵波速度随孔隙度与含水饱和度的变化趋势图；

图14为孔隙流体类型为水和气情况下，Gassmann流体项随孔隙度与含水饱和度的变化趋势图；

图15为孔隙流体类型为水和气情况下，流体体积模量随孔隙度与含水饱和度的变化趋势图；

图16为一维地层模型曲线示意图；

图17为未加噪声时的叠前CMP道集示意图；

图18为利用未加噪声时的叠前CMP道集进行反演得到的反演结果示意图；

图19为信噪比为1:1的叠前CMP道集示意图；

图20为利用信噪比为1:1的叠前CMP道集进行反演得到的反演结果示意图；

图21为研究测线的叠加地震剖面示意图；

图22为估算的流体体积模量K_f的剖面示意图。

具体实施方式

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

由于地下岩石的固液耦合特性，现阶段基于AVO反演方法估算的弹性参数(如波阻抗，Gassmann流体项等)均会受岩石固体效应(岩石基质，孔隙度等)影响产生流体识别假象，降低了地质情况复杂地区储层含油气预测的可靠性。本发明基于孔隙弹性介质理论以及临界孔隙度模型推导新的AVO近似公式，在反射系数近似中实现流体体积模量参数分离，基于新的近似公式进行AVO反演可以得到对储层流体信息更为敏感的弹性参数——流体体积模量。

图1为本发明提出的一种基于流体体积模量AVO反演的流体识别的方法的实施方式一的具体流程图，为了提高流体体积模量估算的稳定性，本发明基于贝叶斯理论建立反演目标函数，利用重加权最小二乘迭代算法循环迭代求解目标函数，实现地震尺度流体体积模量的定量估算。由图1可知，在实施方式一中，所述的方法包括：

S101：采集目的储层的叠前地震资料、测井资料以及地质成果资料。

叠前地震资料的采集在具体的实施方式中可通过如下方式：在勘探工区中，根据目的层的深度和地质特点，考虑动校正拉伸、干扰波、多次波等因素来设计宽方位角的三维观测系统，保证有足够的偏移距、方位角。经过激发、接收，得到可以满足AVO分析需求的宽方位大角度的叠前地震资料。

测井资料的采集在具体的实施方式中可通过如下方式：在勘探工区中进行全波列测井，得到测井资料，主要包含纵波速度、横波速度、密度等全波列测井曲线，孔隙度等解释成果曲线，测井层位，录井资料等。

S102：构建表征流体体积模量的AVO近似公式。

AVO(Amplitude Versus Offset，振幅随偏移距的变化)近似公式是进行AVO反演的基础，本发明基于孔隙弹性介质理论以及临界孔隙度模型，推导得到的表征流体体积模量的AVO近似公式，其为基于AVO反演的流体识别提供了理论基础。

Russell等人从Aki-Richards近似公式出发，推导了突出储层流体特征的反射特征近似公式：

$R_{\mathrm{PP}}\left(\mathrm{\θ}\right)=\left[(1-\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{{\mathrm{\γ}}_{\mathrm{sat}}^2})\frac{{\sec }^2\mathrm{\θ}}{4}\right]\frac{\mathrm{\Δf}}{f}+[\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{4{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\θ}-\frac{2}{{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\θ}]\frac{\mathrm{\Δ\μ}}{\mathrm{\μ}}+[\frac{1}{2}-\frac{{\sec }^2\mathrm{\θ}}{4}]\frac{\mathrm{\Δ\ρ}}{\mathrm{\ρ}}---\left(1\right)$

Dehua Han等人提出了Gassmann流体项的经验公式：

f＝G(φ)K_f      (2)

其中， $G\left(\mathrm{\φ}\right)=\frac{{(1-K_n)}^2}{\mathrm{\φ}}{K}_n=\frac{K_{\mathrm{dry}}}{K_m}$

将公式(2)代入公式(1)，并考虑到剪切模量不受孔隙流体的影响，利用干岩石剪切模量μ_dry替换μ，进行相应变换，可以得到：

$R_{\mathrm{PP}}\left(\mathrm{\θ}\right)=\left[(1-\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{{\mathrm{\γ}}_{\mathrm{sat}}^2})\frac{{\sec }^2\mathrm{\θ}}{4}\right]\frac{\mathrm{\Δ}\left(G\left(\mathrm{\φ}\right)K_f\right)}{G\left(\mathrm{\φ}\right)K_f}+[\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{4{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\θ}-\frac{2}{{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\θ}]\frac{\mathrm{\Δ\μ}}{\mathrm{\μ}}+[\frac{1}{2}-\frac{{\sec }^2\mathrm{\θ}}{4}]\frac{\mathrm{\Δ\ρ}}{\mathrm{\ρ}}---\left(3\right)$

Nur提出的临界孔隙度模型，表达式如下式所示：

$\left\{\begin{array}{c}{K}_{\mathrm{dry}}=K_m(1-\frac{\mathrm{\φ}}{{\mathrm{\φ}}_c})\\ {\mathrm{\μ}}_{\mathrm{dry}}={\mathrm{\μ}}_m(1-\frac{\mathrm{\φ}}{{\mathrm{\φ}}_c})\end{array}\right.---\left(4\right)$

其中，φ_c表示临界孔隙度，K_dry表示干岩石的体积模量，μ_dry表示干岩石的剪切模量，K_m表示固体矿物基质的体积模量，μ_m表示矿物基质的剪切模量。

结合Nur公式(4)对公式(3)进一步展开可以得到：

$\begin{array}{c}{R}_{\mathrm{PP}}\left(\mathrm{\θ}\right)=\left[(1-\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{{\mathrm{\γ}}_{\mathrm{sat}}^2})\frac{{\sec }^2\mathrm{\θ}}{4}\right](\frac{\mathrm{\ΔG}\left(\mathrm{\φ}\right)}{G\left(\mathrm{\φ}\right)}+\frac{\mathrm{\Δ}{K}_f}{K_f})+[\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{4{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\θ}-\frac{2}{{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\θ}](\frac{\mathrm{\Δ}{\mathrm{\μ}}_m}{{\mathrm{\μ}}_m}+\frac{\mathrm{\Δ}({\mathrm{\φ}}_c-\mathrm{\φ})}{{\mathrm{\φ}}_c-\mathrm{\φ}})\\ +[\frac{1}{2}-\frac{{\sec }^2\mathrm{\θ}}{4}]\frac{\mathrm{\Δ\ρ}}{\mathrm{\ρ}}\end{array}---\left(5\right)$

将代入G(φ)，进一步展开可以得到：

$\left.\begin{array}{c}{R}_{\mathrm{PP}}\left(\mathrm{\θ}\right)=\left[(1-\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{{\mathrm{\γ}}_{\mathrm{sat}}^2})\frac{{\sec }^2\mathrm{\θ}}{4}\right]\frac{\mathrm{\Δ}{K}_f}{K_f}[+\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{4{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\θ}-\frac{2}{{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\θ}]\frac{\mathrm{\Δ}{\mathrm{\μ}}_m}{{\mathrm{\μ}}_m}+[\frac{{\sec }^2\mathrm{\θ}}{4}-\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{2{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\θ}+\frac{2}{{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\θ}]\frac{\mathrm{\Δ\φ}}{\mathrm{\φ}}\\ +[\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{4{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\θ}-\frac{2}{{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\θ}]\frac{\mathrm{\Δ\φ}({\mathrm{\φ}}_c-\mathrm{\φ})}{\mathrm{\φ}({\mathrm{\φ}}_c-\mathrm{\φ})}+[\frac{1}{2}-\frac{{\sec }^2\mathrm{\θ}}{4}]\frac{\mathrm{\Δ\ρ}}{\mathrm{\ρ}}\end{array}\right]---\left(6\right)$

令f_m＝φμ，并利用公式(4)，最终构建的表征流体体积模量的AVO近似公式如下所示：

$\begin{array}{c}{R}_{\mathrm{PP}}\left(\mathrm{\θ}\right)=\left[(1-\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{{\mathrm{\γ}}_{\mathrm{sat}}^2})\frac{{\sec }^2\mathrm{\θ}}{4}\right]\frac{\mathrm{\Δ}{K}_f}{K_f}+[\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{4{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\θ}-\frac{2}{{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\θ}]\frac{\mathrm{\Δ}\left(f_m\right)}{f_m}+[\frac{1}{2}-\frac{{\sec }^2\mathrm{\θ}}{4}]\frac{\mathrm{\Δ\ρ}}{\mathrm{\ρ}}\\ +(\frac{{\sec }^2\mathrm{\θ}}{4}-\frac{{\mathrm{\γ}}_{\mathrm{dry}}^2}{2{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\θ}+\frac{2}{{\mathrm{\γ}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\θ})\frac{\mathrm{\Δ\φ}}{\mathrm{\φ}}\end{array}---\left(7\right)$

其中，f_m＝φμ为地下岩石的固体项，θ为入射角度，K_f为地下岩石的流体体积模量，f_m为地下岩石的固体项，φ为地下岩石的孔隙度，ρ为地下岩石的密度，μ为地下岩石的剪切模量，△K_f为界面两侧的流体体积模量的差值，△f_m为界面两侧固体项的差值，△ρ为界面两侧的密度的差值，△φ为界面两侧的孔隙度的差值，为干岩石纵横波速度比的平方，为饱和岩石纵横波速度比的平方。

S103：根据所述的叠前地震资料、测井资料、地质成果资料以及所述的表征流体体积模量的AVO近似公式提取子波。图5为步骤S103的具体流程图。

S104：根据所述的叠前地震资料、测井资料、地质成果资料以及子波确定流体体积模量数据体，图6为步骤S104的具体流程图。

S105：根据所述的流体体积模量数据体识别当前储层孔隙中流体的类型。

图2为本发明实施例提供的一种基于流体体积模量AVO反演的流体识别的方法的实施方式二的流程图，由图2可知，在实施方式二中，该方法具体包括：

S201：采集目的储层的叠前地震资料、测井资料以及地质成果资料。

叠前地震资料的采集在具体的实施方式中可通过如下方式：在勘探工区中，根据目的层的深度和地质特点，考虑动校正拉伸、干扰波、多次波等因素来设计宽方位角的三维观测系统，保证有足够的偏移距、方位角。经过激发、接收，得到可以满足AVO分析需求的宽方位大角度的叠前地震资料。

测井资料的采集在具体的实施方式中可通过如下方式：在勘探工区中进行全波列测井，得到测井资料，主要包含纵波速度、横波速度、密度等全波列测井曲线，孔隙度等解释成果曲线，测井层位，录井资料等。

S202：对所述的叠前地震资料进行保幅处理，得到CMP道集(即共中心点道集)资料。

对所述的叠前地震资料进行保幅处理，以保证形成可用于进行AVO反演的高品质CMP道集资料。具体的保幅处理包括：波前扩散补偿、震源组合与检波器组合效应的校正、反Q滤波、地表一致性处理、多次波消除、叠前去噪等。

S203：构建表征流体体积模量的AVO近似公式方程。

S204：根据所述的叠前地震资料、测井资料、地质成果资料以及所述的表征流体体积模量的AVO近似公式提取子波。图5为步骤S204的具体流程图。

S205：根据所述的叠前地震资料、测井资料、地质成果资料以及子波确定流体体积模量数据体，图6为步骤S205的具体流程图。

S206：根据所述的流体体积模量数据体识别当前储层孔隙中流体的类型。

图3为本发明实施例提供的一种基于流体体积模量AVO反演的流体识别的方法的实施方式三的流程图，由图3可知，在实施方式三中，该方法具体包括：

S301：采集目的储层的叠前地震资料、测井资料以及地质成果资料。

叠前地震资料的采集在具体的实施方式中可通过如下方式：在勘探工区中，根据目的层的深度和地质特点，考虑动校正拉伸、干扰波、多次波等因素来设计宽方位角的三维观测系统，保证有足够的偏移距、方位角。经过激发、接收，得到可以满足AVO分析需求的宽方位大角度的叠前地震资料。

测井资料的采集在具体的实施方式中可通过如下方式：在勘探工区中进行全波列测井，得到测井资料，主要包含纵波速度、横波速度、密度等全波列测井曲线，孔隙度等解释成果曲线，测井层位，录井资料等。

S302：对所述的叠前地震资料进行保幅处理，得到CMP道集资料。

对所述的叠前地震资料进行保幅处理，以保证形成可用于进行AVO反演的高品质CMP道集资料。具体的保幅处理包括：波前扩散补偿、震源组合与检波器组合效应的校正、反Q滤波、地表一致性处理、多次波消除、叠前去噪等。

S303：对所述的测井资料进行环境校正；

S304：对经过环境校正后的测井资料进行奇异值消除；

S305：根据经过奇异值消除后的测井资料确定所述目的储层的干岩石纵横波速度比、饱和岩石纵横波速度比。具体的实施方式中，在勘探工区中进行全波列测井，得到的测井资料作为先验信息，需要经过环境校正、奇异值消除等预处理后，结合岩石物理经验公式或者实验统计结果计算干岩石纵横波速度比γ_dry、饱和岩石纵横波速度比γ_sat的数值。

S306：根据所述的目的储层的干岩样纵横波速度比以及所述的测井资料确定与所述测井资料对应的流体体积模量曲线以及剪切模量曲线。在具体的实施方式中，流体体积模量为K_f，剪切模量为μ曲线。

S307：构建表征流体体积模量的AVO近似公式方程。

S308：根据所述的叠前地震资料、测井资料、地质成果资料以及所述的表征流体体积模量的AVO近似公式提取子波。图5为步骤S308的具体流程图。

S309：根据所述的叠前地震资料、测井资料、地质成果资料以及子波确定流体体积模量数据体，图6为步骤S309的具体流程图。

S310：根据所述的流体体积模量数据体识别当前储层孔隙中流体的类型。

图4为本发明实施例提供的一种基于流体体积模量AVO反演的流体识别的方法的实施方式四的流程图，由图4可知，在实施方式四中，该方法具体包括：

S401：采集目的储层的叠前地震资料、测井资料以及地质成果资料。

叠前地震资料的采集在具体的实施方式中可通过如下方式：在勘探工区中，根据目的层的深度和地质特点，考虑动校正拉伸、干扰波、多次波等因素来设计宽方位角的三维观测系统，保证有足够的偏移距、方位角。经过激发、接收，得到可以满足AVO分析需求的宽方位大角度的叠前地震资料。

测井资料的采集在具体的实施方式中可通过如下方式：在勘探工区中进行全波列测井，得到测井资料，主要包含纵波速度、横波速度、密度等全波列测井曲线，孔隙度等解释成果曲线，测井层位，录井资料等。

S402：对所述的叠前地震资料进行保幅处理，得到CMP道集资料。

对所述的叠前地震资料进行保幅处理，以保证形成可用于进行AVO反演的高品质CMP道集资料。具体的保幅处理包括：波前扩散补偿、震源组合与检波器组合效应的校正、反Q滤波、地表一致性处理、多次波消除、叠前去噪等。

S403：对所述的测井资料进行环境校正；

S404：对经过环境校正后的测井资料进行奇异值消除；

S405：根据经过奇异值消除后的测井资料确定所述目的储层的干岩石纵横波速度比、饱和岩石纵横波速度比。具体的实施方式中，在勘探工区中进行全波列测井，得到的测井资料作为先验信息，需要经过环境校正、奇异值消除等预处理后，结合岩石物理经验公式或者实验统计结果计算干岩石纵横波速度比γ_dry、饱和岩石纵横波速度比γ_sat的数值。

S406：根据所述的目的储层的干岩样纵横波速度比以及所述的测井资料确定与所述测井资料对应的流体体积模量曲线以及剪切模量曲线。在具体的实施方式中，流体体积模量为K_f，剪切模量为μ曲线。

S407：根据所述的叠前地震资料、测井资料、地质成果资料对目的储层进行精细层位追踪，得到地震层位数据。

在具体的实施方式中，结合工区的地质背景，在保证目的层位符合工区地质沉积模式与层序格架的基础上，对研究的目的层进行精细地质解释，解释得到的地震层位资料不仅用来标识时深关系，还用于建立反演过程的模型约束。

S408：构建表征流体体积模量的AVO近似公式方程。

S409：根据所述的叠前地震资料、测井资料、地质成果资料以及所述的表征流体体积模量的AVO近似公式提取子波。图5为步骤S409的具体流程图。由图5可知，该步骤具体包括：

S501：从所述的测井资料中提取出地下岩石的密度以及孔隙度。在具体的实施方式中，地下岩石的密度通过ρ来表示，地下岩石的孔隙度通过φ来表示。

S502：根据所述的剪切模量曲线以及所述的孔隙度确定地下岩石的固体项，在具体的实施方式中，地下岩石的固体项通过f_m来表示，则f_m＝φμ。

S503：根据所述地下岩石的密度、固体项、孔隙度、流体体积模量曲线、表征流体体积模量的AVO近似公式、叠前CMP道集以及地震层位数据提取子波。

由图4可知，在实施方式四中，该方法还包括：

S410：根据所述的叠前地震资料、测井资料、地质成果资料以及子波确定流体体积模量数据体。图6为步骤S410的具体流程图。由图6可知，该步骤具体包括：

S601：利用所述的叠前CMP道集、测井资料、地质成果资料以及子波，基于贝叶斯理论，使用高斯分布函数为似然函数、柯西分布为先验约束建立AVO反演目标函数。建立的AVO反演目标函数为：

(G^ΤG+ηQ+αC^ΤC)m＝G^Τd+αC^Τξ      (8)

$\mathrm{\η}=\mathrm{\λ}\frac{{\mathrm{\σ}}_n^2}{{\mathrm{\σ}}_m^2}$

$\mathrm{\ξ}=\frac{1}{2}\ln \left(\frac{M}{M_0}\right)$

其中，为CMP道集资料中噪音的方差，为弹性参数相对变化量的方差，λ、α为加权系数，m为待反演的弹性参数相对变化量矩阵，G为正演算子矩阵，C为积分矩阵，d为CMP道集矩阵，Q为斜对角加权矩阵，M为模型约束矩阵，M₀为M的均值。

S602：利用重加权最小二乘迭代算法求解所述的AVO反演目标函数，得到流体体积模量的相对变化量。在具体的实施方式中，该步骤具体包括：

(1)、建立初始反射系数序列(初始模型)m₀；

(2)、根据地震资料的品质，选择合适的加权系数λ、α；

(3)、根据先验信息、模型数据及测井数据计算G^TG、G^Τd、C^ΤC、C^Τξ，利用第k-1次迭代的反演结果m_k-1(第一次迭代利用初始模型m₀)计算η和Q；

(4)、对公式(8)进行求逆计算，得到第k次反演结果：

m_k＝(G^ΤG+ηQ+αC^ΤC)^-1(G^Τd+αC^Τξ)

在本发明的其他实施方式中，除了上述步骤之外，还可以设计如下收敛条件：

$\frac{|J\left(m_k\right)-J\left(m_{k-1}\right)|}{\left(\right|J\left(m_k\right)|+|J\left(m_{k-1}\right)\left|\right)/2}<e$

_e为收敛误差。利用m_k计算得到目标函数值J(m_k)，当满足上述收敛条件时终止迭代，当不满足收敛条件时，跳回第(2)步，迭代继续进行。

S603：利用道积分思想将所述的相对变化量转化为流体体积模量数据体。

由图4可知，在实施方式四中，该方法还包括：

S411：根据所述的流体体积模量数据体识别当前储层孔隙中流体的类型。

估算得到流体体积模量数据体以后，与储层流体体积模量的理论值进行对比，即可以识别储层孔隙中流体的类型。

如上所述，即为本发明提供的一种基于流体体积模量AVO反演的流体识别的方法，实现了基于流体模量分离AVO近似公式的流体体积模量参数估算，流体体积模量参数对储层含流体类型敏感度高，且最大限度的消除了储层岩石孔隙度等固体效应造成的流体识别假象，提高了储层孔隙流体类型预测的可靠性。

图7为本发明提出的一种基于流体体积模量AVO反演的流体识别的系统的实施方式一的具体结构框图，为了提高流体体积模量估算的稳定性，本发明基于贝叶斯理论建立反演目标函数，利用重加权最小二乘迭代算法循环迭代求解目标函数，实现地震尺度流体体积模量的定量估算。由图7可知，在实施方式一中，所述的系统包括：

资料采集装置100，用于采集目的储层的叠前地震资料、测井资料以及地质成果资料。

叠前地震资料的采集在具体的实施方式中可通过如下方式：在勘探工区中，根据目的层的深度和地质特点，考虑动校正拉伸、干扰波、多次波等因素来设计宽方位角的三维观测系统，保证有足够的偏移距、方位角。经过激发、接收，得到可以满足AVO分析需求的宽方位大角度的叠前地震资料。

测井资料的采集在具体的实施方式中可通过如下方式：在勘探工区中进行全波列测井，得到测井资料，主要包含纵波速度、横波速度、密度等全波列测井曲线，孔隙度等解释成果曲线，测井层位，录井资料等。

近似公式构建装置200，用于构建表征流体体积模量的AVO近似公式。

AVO(Amplitude Versus Offset，振幅随偏移距的变化)近似公式是进行AVO反演的基础，本发明基于孔隙弹性介质理论以及临界孔隙度模型，推导得到的表征流体体积模量的AVO近似公式，其为基于AVO反演的流体识别提供了理论基础。

Russell等人从Aki-Richards近似公式出发，推导了突出储层流体特征的反射特征近似公式：

$R_{\mathrm{PP}}\left(\mathrm{\&theta;}\right)=\left[(1-\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2})\frac{{\sec }^2\mathrm{\&theta;}}{4}\right]\frac{\mathrm{\&Delta;f}}{f}+[\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{4{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\&theta;}-\frac{2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\&theta;}]\frac{\mathrm{\&Delta;\&mu;}}{\mathrm{\&mu;}}+[\frac{1}{2}-\frac{{\sec }^2\mathrm{\&theta;}}{4}]\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}}---\left(1\right)$

Dehua Han等人提出了Gassmann流体项的经验公式：

f＝G(φ)K_f        (2)

其中， $G\left(\mathrm{\&phi;}\right)=\frac{{(1-K_n)}^2}{\mathrm{\&phi;}}{K}_n=\frac{K_{\mathrm{dry}}}{K_m}$

将公式(2)代入公式(1)，并考虑到剪切模量不受孔隙流体的影响，利用干岩石剪切模量μ_dry替换μ，进行相应变换，可以得到：

$R_{\mathrm{PP}}\left(\mathrm{\&theta;}\right)=\left[(1-\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2})\frac{{\sec }^2\mathrm{\&theta;}}{4}\right]\frac{\mathrm{\&Delta;}\left(G\left(\mathrm{\&phi;}\right)K_f\right)}{G\left(\mathrm{\&phi;}\right)K_f}+[\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{4{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\&theta;}-\frac{2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\&theta;}]\frac{\mathrm{\&Delta;\&mu;}}{\mathrm{\&mu;}}+[\frac{1}{2}-\frac{{\sec }^2\mathrm{\&theta;}}{4}]\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}}---\left(3\right)$

Nur提出的临界孔隙度模型，表达式如下式所示：

$\left\{\begin{array}{c}{K}_{\mathrm{dry}}=K_m(1-\frac{\mathrm{\&phi;}}{{\mathrm{\&phi;}}_c})\\ {\mathrm{\&mu;}}_{\mathrm{dry}}={\mathrm{\&mu;}}_m(1-\frac{\mathrm{\&phi;}}{{\mathrm{\&phi;}}_c})\end{array}\right.---\left(4\right)$

其中，φ_c表示临界孔隙度，K_dry表示干岩石的体积模量，μ_dry表示干岩石的剪切模量，K_m表示固体矿物基质的体积模量，μ_m表示矿物基质的剪切模量。

结合Nur公式(4)对公式(3)进一步展开可以得到：

$\begin{array}{c}{R}_{\mathrm{PP}}\left(\mathrm{\&theta;}\right)=\left[(1-\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2})\frac{{\sec }^2\mathrm{\&theta;}}{4}\right](\frac{\mathrm{\&Delta;G}\left(\mathrm{\&phi;}\right)}{G\left(\mathrm{\&phi;}\right)}+\frac{\mathrm{\&Delta;}{K}_f}{K_f})+[\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{4{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\&theta;}-\frac{2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\&theta;}](\frac{\mathrm{\&Delta;}{\mathrm{\&mu;}}_m}{{\mathrm{\&mu;}}_m}+\frac{\mathrm{\&Delta;}({\mathrm{\&phi;}}_c-\mathrm{\&phi;})}{{\mathrm{\&phi;}}_c-\mathrm{\&phi;}})\\ +[\frac{1}{2}-\frac{{\sec }^2\mathrm{\&theta;}}{4}]\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}}\end{array}---\left(5\right)$

将代入G(φ)，进一步展开可以得到：

$\left.\begin{array}{c}{R}_{\mathrm{PP}}\left(\mathrm{\&theta;}\right)=\left[(1-\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2})\frac{{\sec }^2\mathrm{\&theta;}}{4}\right]\frac{\mathrm{\&Delta;}{K}_f}{K_f}[+\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{4{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\&theta;}-\frac{2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\&theta;}]\frac{\mathrm{\&Delta;}{\mathrm{\&mu;}}_m}{{\mathrm{\&mu;}}_m}+[\frac{{\sec }^2\mathrm{\&theta;}}{4}-\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{2{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\&theta;}+\frac{2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\&theta;}]\frac{\mathrm{\&Delta;\&phi;}}{\mathrm{\&phi;}}\\ +[\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{4{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\&theta;}-\frac{2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\&theta;}]\frac{\mathrm{\&Delta;\&phi;}({\mathrm{\&phi;}}_c-\mathrm{\&phi;})}{\mathrm{\&phi;}({\mathrm{\&phi;}}_c-\mathrm{\&phi;})}+[\frac{1}{2}-\frac{{\sec }^2\mathrm{\&theta;}}{4}]\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}}\end{array}\right]---\left(6\right)$

令f_m＝φμ，并利用公式(4)，最终构建的表征流体体积模量的AVO近似公式如下所示：

$\begin{array}{c}{R}_{\mathrm{PP}}\left(\mathrm{\&theta;}\right)=\left[(1-\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2})\frac{{\sec }^2\mathrm{\&theta;}}{4}\right]\frac{\mathrm{\&Delta;}{K}_f}{K_f}+[\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{4{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\&theta;}-\frac{2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\&theta;}]\frac{\mathrm{\&Delta;}\left(f_m\right)}{f_m}+[\frac{1}{2}-\frac{{\sec }^2\mathrm{\&theta;}}{4}]\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}}\\ +(\frac{{\sec }^2\mathrm{\&theta;}}{4}-\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{2{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\&theta;}+\frac{2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\&theta;})\frac{\mathrm{\&Delta;\&phi;}}{\mathrm{\&phi;}}\end{array}---\left(7\right)$

其中，f_m＝φμ为地下岩石的固体项，θ为入射角度，K_f为地下岩石的流体体积模量，f_m为地下岩石的固体项，φ为地下岩石的孔隙度，ρ为地下岩石的密度，μ为地下岩石的剪切模量，△K_f为界面两侧的流体体积模量的差值，△f_m为界面两侧固体项的差值，△ρ为界面两侧的密度的差值，△φ为界面两侧的孔隙度的差值，为干岩石纵横波速度比的平方，为饱和岩石纵横波速度比的平方。

子波提取装置300，用于根据所述的叠前地震资料、测井资料、地质成果资料以及所述的表征流体体积模量的AVO近似公式提取子波。图11为子波提取装置的具体结构框图。

流体体积模量数据体确定装置400，用于根据所述的叠前地震资料、测井资料、地质成果资料以及子波确定流体体积模量数据体，图12为流体体积模量数据体确定装置的具体结构框图。

流体类型识别装置500，用于根据所述的流体体积模量数据体识别当前储层孔隙中流体的类型。

图8为本发明实施例提供的一种基于流体体积模量AVO反演的流体识别的系统的实施方式二的结构框图，由图8可知，在实施方式二中，该系统具体包括：

资料采集装置100，用于采集目的储层的叠前地震资料、测井资料以及地质成果资料。

叠前地震资料的采集在具体的实施方式中可通过如下方式：在勘探工区中，根据目的层的深度和地质特点，考虑动校正拉伸、干扰波、多次波等因素来设计宽方位角的三维观测系统，保证有足够的偏移距、方位角。经过激发、接收，得到可以满足AVO分析需求的宽方位大角度的叠前地震资料。

测井资料的采集在具体的实施方式中可通过如下方式：在勘探工区中进行全波列测井，得到测井资料，主要包含纵波速度、横波速度、密度等全波列测井曲线，孔隙度等解释成果曲线，测井层位，录井资料等。

保幅处理装置600，用于对所述的叠前地震资料进行保幅处理，得到CMP道集资料。

对所述的叠前地震资料进行保幅处理，以保证形成可用于进行AVO反演的高品质CMP道集资料。具体的保幅处理包括：波前扩散补偿、震源组合与检波器组合效应的校正、反Q滤波、地表一致性处理、多次波消除、叠前去噪等。

近似公式构建装置200，用于构建表征流体体积模量的AVO近似公式方程。

子波提取装置300，用于根据所述的叠前地震资料、测井资料、地质成果资料以及所述的表征流体体积模量的AVO近似公式提取子波。图11为子波提取装置的具体结构框图。

流体体积模量数据体确定装置400，用于根据所述的叠前地震资料、测井资料、地质成果资料以及子波确定流体体积模量数据体，图12为流体体积模量数据体确定装置的具体结构框图。

流体类型识别装置500，用于根据所述的流体体积模量数据体识别当前储层孔隙中流体的类型。

图9为本发明实施例提供的一种基于流体体积模量AVO反演的流体识别的系统的实施方式三的结构框图，由图9可知，在实施方式三中，该系统具体包括：

资料采集装置100，用于采集目的储层的叠前地震资料、测井资料以及地质成果资料。

叠前地震资料的采集在具体的实施方式中可通过如下方式：在勘探工区中，根据目的层的深度和地质特点，考虑动校正拉伸、干扰波、多次波等因素来设计宽方位角的三维观测系统，保证有足够的偏移距、方位角。经过激发、接收，得到可以满足AVO分析需求的宽方位大角度的叠前地震资料。

测井资料的采集在具体的实施方式中可通过如下方式：在勘探工区中进行全波列测井，得到测井资料，主要包含纵波速度、横波速度、密度等全波列测井曲线，孔隙度等解释成果曲线，测井层位，录井资料等。

保幅处理装置600，用于对所述的叠前地震资料进行保幅处理，得到CMP道集资料。

对所述的叠前地震资料进行保幅处理，以保证形成可用于进行AVO反演的高品质CMP道集资料。具体的保幅处理包括：波前扩散补偿、震源组合与检波器组合效应的校正、反Q滤波、地表一致性处理、多次波消除、叠前去噪等。

环境校正装置700，用于对所述的测井资料进行环境校正；

奇异值消除装置800，用于对经过环境校正后的测井资料进行奇异值消除；

纵横速度比确定装置900，用于根据经过奇异值消除后的测井资料确定所述目的储层的干岩石纵横波速度比、饱和岩石纵横波速度比。具体的实施方式中，在勘探工区中进行全波列测井，得到的测井资料作为先验信息，需要经过环境校正、奇异值消除等预处理后，结合岩石物理经验公式或者实验统计结果计算干岩石纵横波速度比γ_dry、饱和岩石纵横波速度比γ_sat的数值。

剪切模量曲线确定装置1000，用于根据所述的目的储层的干岩样纵横波速度比以及所述的测井资料确定与所述测井资料对应的流体体积模量曲线以及剪切模量曲线。在具体的实施方式中，流体体积模量为K_f，剪切模量为μ曲线。

近似公式构建装置200，用于构建表征流体体积模量的AVO近似公式方程。

子波提取装置300，用于根据所述的叠前地震资料、测井资料、地质成果资料以及所述的表征流体体积模量的AVO近似公式提取子波。图11为子波提取装置的具体结构框图。

流体体积模量数据体确定装置400，用于根据所述的叠前地震资料、测井资料、地质成果资料以及子波确定流体体积模量数据体，图12为流体体积模量数据体确定装置的具体结构框图。

流体类型识别装置500，用于根据所述的流体体积模量数据体识别当前储层孔隙中流体的类型。

图10为本发明实施例提供的一种基于流体体积模量AVO反演的流体识别的系统的实施方式四的结构框图，由图10可知，在实施方式四中，该系统具体包括：

资料采集装置100，用于采集目的储层的叠前地震资料、测井资料以及地质成果资料。

叠前地震资料的采集在具体的实施方式中可通过如下方式：在勘探工区中，根据目的层的深度和地质特点，考虑动校正拉伸、干扰波、多次波等因素来设计宽方位角的三维观测系统，保证有足够的偏移距、方位角。经过激发、接收，得到可以满足AVO分析需求的宽方位大角度的叠前地震资料。

测井资料的采集在具体的实施方式中可通过如下方式：在勘探工区中进行全波列测井，得到测井资料，主要包含纵波速度、横波速度、密度等全波列测井曲线，孔隙度等解释成果曲线，测井层位，录井资料等。

保幅处理装置600，用于对所述的叠前地震资料进行保幅处理，得到CMP道集资料。

对所述的叠前地震资料进行保幅处理，以保证形成可用于进行AVO反演的高品质CMP道集资料。具体的保幅处理包括：波前扩散补偿、震源组合与检波器组合效应的校正、反Q滤波、地表一致性处理、多次波消除、叠前去噪等。

环境校正装置700，用于对所述的测井资料进行环境校正；

奇异值消除装置800，用于对经过环境校正后的测井资料进行奇异值消除；

纵横速度比确定装置900，用于根据经过奇异值消除后的测井资料确定所述目的储层的干岩石纵横波速度比、饱和岩石纵横波速度比。具体的实施方式中，在勘探工区中进行全波列测井，得到的测井资料作为先验信息，需要经过环境校正、奇异值消除等预处理后，结合岩石物理经验公式或者实验统计结果计算干岩石纵横波速度比γ_dry、饱和岩石纵横波速度比γ_sat的数值。

剪切模量曲线确定装置1000，用于根据所述的目的储层的干岩样纵横波速度比以及所述的测井资料确定与所述测井资料对应的流体体积模量曲线以及剪切模量曲线。在具体的实施方式中，流体体积模量为K_f，剪切模量为μ曲线。

层位追踪装置1100，用于根据所述的叠前地震资料、测井资料、地质成果资料对目的储层进行精细层位追踪，得到地震层位数据。

在具体的实施方式中，结合工区的地质背景，在保证目的层位符合工区地质沉积模式与层序格架的基础上，对研究的目的层进行精细地质解释，解释得到的地震层位资料不仅用来标识时深关系，还用于建立反演过程的模型约束。

近似公式构建装置200，用于构建表征流体体积模量的AVO近似公式方程。

子波提取装置300，用于根据所述的叠前地震资料、测井资料、地质成果资料以及所述的表征流体体积模量的AVO近似公式提取子波。图11为子波提取装置300的具体结构框图。由图11可知，该子波提取装置具体包括：

孔隙度提取模块301，用于从所述的测井资料中提取出地下岩石的密度以及孔隙度。在具体的实施方式中，地下岩石的密度通过ρ来表示，地下岩石的孔隙度通过φ来表示。

固体项确定模块302，用于根据所述的剪切模量曲线以及所述的孔隙度确定地下岩石的固体项，在具体的实施方式中，地下岩石的固体项通过f_m来表示，则f_m＝φμ。

子波提取模块303，用于根据所述地下岩石的密度、固体项、孔隙度、流体体积模量曲线、表征流体体积模量的AVO近似公式、叠前CMP道集以及地震层位数据提取子波。

由图4可知，在实施方式四中，该系统还包括：

流体体积模量数据体确定装置400，用于根据所述的叠前地震资料、测井资料、地质成果资料以及子波确定流体体积模量数据体。图12为流体体积模量数据体确定装置的具体结构框图。由图12可知，流体体积模量数据体确定装置具体包括：

反演目标函数建立模块401，用于利用所述的叠前CMP道集、测井资料、地质成果资料以及子波，基于贝叶斯理论，使用高斯分布函数为似然函数、柯西分布为先验约束建立AVO反演目标函数。建立的AVO反演目标函数为：

(G^ΤG+ηQ+αC^ΤC)m＝G^Τd+αC^Τξ       (8)

$\mathrm{\&eta;}=\mathrm{\&lambda;}\frac{{\mathrm{\&sigma;}}_n^2}{{\mathrm{\&sigma;}}_m^2}$

$\mathrm{\&xi;}=\frac{1}{2}\ln \left(\frac{M}{M_0}\right)$

其中，为CMP道集资料中噪音的方差，为弹性参数相对变化量的方差，λ、α为加权系数，m为待反演的弹性参数相对变化量矩阵，G为正演算子矩阵，C为积分矩阵，d为CMP道集矩阵，Q为斜对角加权矩阵，M为模型约束矩阵，M₀为M的均值。

相对变化量确定模块402，用于利用重加权最小二乘迭代算法求解所述的AVO反演目标函数，得到流体体积模量的相对变化量。在具体的实施方式中，该步骤具体包括：

(1)、建立初始反射系数序列(初始模型)m₀；

(2)、根据地震资料的品质，选择合适的加权系数λ、α；

(3)、根据先验信息、模型数据及测井数据计算G^TG、G^Τd、C^ΤC、C^Τξ，利用第k-1次迭代的反演结果m_k-1(第一次迭代利用初始模型m₀)计算η和Q；

(4)、对公式(8)进行求逆计算，得到第k次反演结果：

m_k＝(G^ΤG+ηQ+αC^ΤC)^-1(G^Τd+αC^Τξ)

在本发明的其他实施方式中，除了上述步骤之外，还可以设计如下收敛条件：

$\frac{|J\left(m_k\right)-J\left(m_{k-1}\right)|}{\left(\right|J\left(m_k\right)|+|J\left(m_{k-1}\right)\left|\right)/2}<e$

_e为收敛误差。利用m_k计算得到目标函数值J(m_k)，当满足上述收敛条件时终止迭代，当不满足收敛条件时，跳回第(2)步，迭代继续进行。

转化模块403，用于利用道积分思想将所述的相对变化量转化为流体体积模量数据体。

由图10可知，在实施方式四中，该系统还包括：

流体类型识别装置500，用于根据所述的流体体积模量数据体识别当前储层孔隙中流体的类型。

估算得到流体体积模量数据体以后，与储层流体体积模量的理论值进行对比，即可以识别储层孔隙中流体的类型。

如上所述，即为本发明提供的一种基于流体体积模量AVO反演的流体识别的系统，实现了基于流体模量分离AVO近似公式的流体体积模量参数估算，流体体积模量参数对储层含流体类型敏感度高，且最大限度的消除了储层岩石孔隙度等固体效应造成的流体识别假象，提高了储层孔隙流体类型预测的可靠性。

下面结合具体的实施例，详细介绍本发明的技术方案。

图13为孔隙流体类型为水和气情况下，纵波速度随孔隙度与含水饱和度的变化趋势图，图14为孔隙流体类型为水和气情况下，Gassmann流体项随孔隙度与含水饱和度的变化趋势图，图15为孔隙流体类型为水和气情况下，流体体积模量随孔隙度与含水饱和度的变化趋势图。为了研究流体指示敏感性，采用流体替代方法设计了三种模型，即常规用来指示流体类型的岩石弹性参数纵波速度V_P，Gassmann流体项f以及流体体积模量K_f随孔隙流体类型以及孔隙度的变化模型。从图中可以看到，相比较于孔隙流体类型，孔隙度对常规弹性参数V_P影响大，f次之，K_f的取值完全依赖于储层流体类型，即在实际应用中利用常规弹性参数进行流体识别时会存在储层岩石的固体效应(特别是孔隙度)造成的假象，而K_f与含水饱和度成完全线性变化趋势，可以作为一种更为敏感的弹性参数参与流体识别。

图16为一维地层模型曲线示意图，用来验证反演方法的可行性与抗噪性。褶积40Hz零相位Ricker子波生成叠前CMP道集(入射角度范围是0—30度，时间采样间隔是2ms)，以此为基础对不同信噪比条件下的反演结果进行分析。

图17为未加噪声时的叠前CMP道集示意图，图18为利用未加噪声时的叠前CMP道集进行反演得到的反演结果示意图。其中实线表示实际模型，虚线表示反演结果。通过比较发现，当无噪音时，在先验信息的约束下，流体体积模量K_f的反演结果与模型数据吻合较好。

图19为信噪比为1:1的叠前CMP道集示意图，图20为利用信噪比为1:1的叠前CMP道集进行反演得到的反演结果示意图。其中实线表示实际模型，虚线表示反演结果。当信噪比为1:1时，反演的流体体积模量K_f虽然受到一定的影响，但总体趋势仍然接近实际模型，误差在可接受范围内。因此通过模型测试知道，采用本研究提出的流体模量提取方法在测井先验信息的约束下对噪音具有较好的压制，能够较好的实现流体体积模量参数的定量提取，为流体识别提供可靠的弹性参数数据支撑。

图21为研究测线的叠加地震剖面示意图。图中所投曲线为测井解释柱状图，三角形表示气层，菱形表示水层。其2.58s与2.64s附近发育两套砂岩，测井资料显示2.58s处砂岩含气，而2.64s处砂岩含水，但由于受古地形影响，下覆含水砂岩孔隙度高于上覆含气砂岩的孔隙度，在流体识别中易造成含流体异常假象。

图22为估算的流体体积模量K_f的剖面示意图。可以看到其与测井解释结果吻合较好，较为准确的反映出了井周围储层含流体特征及分布情况，验证了流体体积模量参数在流体识别中的有效性。

综上所述，本发明提出的一种基于流体体积模量AVO反演的流体识别的方法及系统，为了提高流体体积模量估算的稳定性，本方案基于贝叶斯理论建立反演目标函数，利用重加权最小二乘迭代算法循环迭代求解目标函数，实现地震尺度流体体积模量的定量估算。在具体实现过程中首先需要进行叠前CMP道集资料和测井资料的预处理，然后进行子波提取，在实现流体体积模量相对变化量提取的基础上最终得到流体体积模量参数。

基于常规AVO近似估算的岩石弹性参数受固液耦合特性影响易产生储层含流体判识假象，本发明基于孔隙弹性介质理论以及临界孔隙度模型，在反射系数近似中实现流体体积模量参数分离，流体体积模量对储层流体信息更为敏感，利用表征流体体积模量的AVO近似公式进行AVO反演可以估算出流体体积模量参数，避免流体识别假象，提高储层孔隙流体类型预测的可靠性。

本领域普通技术人员可以理解实现上述实施例方法中的全部或部分流程，可以通过计算机程序来指令相关的硬件来完成，所述的程序可存储于一般计算机可读取存储介质中，该程序在执行时，可包括如上述各方法的实施例的流程。其中，所述的存储介质可为磁碟、光盘、只读存储记忆体(Read-Only Memory，ROM)或随机存储记忆体(Random AccessMemory，RAM)等。

本领域技术人员还可以了解到本发明实施例列出的各种功能是通过硬件还是软件来实现取决于特定的应用和整个系统的设计要求。本领域技术人员可以对于每种特定的应用，可以使用各种方法实现所述的功能，但这种实现不应被理解为超出本发明实施例保护的范围。

本发明中应用了具体实施例对本发明的原理及实施方式进行了阐述，以上实施例的说明只是用于帮助理解本发明的方法及其核心思想；同时，对于本领域的一般技术人员，依据本发明的思想，在具体实施方式及应用范围上均会有改变之处，综上所述，本说明书内容不应理解为对本发明的限制。

Claims (16)

1.一种基于流体体积模量AVO反演的流体识别的方法，其特征是，所述的方法具体包括：

采集目的储层的叠前地震资料、测井资料以及地质成果资料；

构建表征流体体积模量的振幅随偏移距的变化AVO近似公式；

根据所述的叠前地震资料、测井资料、地质成果资料以及所述的表征流体体积模量的AVO近似公式提取子波；

根据所述的叠前地震资料、测井资料、地质成果资料以及子波确定流体体积模量数据体；

根据所述的流体体积模量数据体识别当前储层孔隙中流体的类型。

2.根据权利要求1所述的方法，其特征是，所述的方法在构建表征流体体积模量的AVO近似公式之前还包括：

对所述的叠前地震资料进行保幅处理，得到共中心点CMP道集资料。

3.根据权利要求2所述的方法，其特征是，所述的方法在构建表征流体体积模量的AVO近似公式之前还包括：

对所述的测井资料进行环境校正；

对经过环境校正后的测井资料进行奇异值消除；

根据经过奇异值消除后的测井资料确定所述目的储层的干岩石纵横波速度比、饱和岩石纵横波速度比；

根据所述的目的储层的干岩样纵横波速度比以及所述的测井资料确定与所述测井资料对应的流体体积模量曲线以及剪切模量曲线。

4.根据权利要求3所述的方法，其特征是，所述的方法在构建表征流体体积模量的AVO近似公式之前还包括：

根据所述的叠前地震资料、测井资料、地质成果资料对目的储层进行精细层位追踪，得到地震层位数据。

5.根据权利要求4所述的方法，其特征是，所述的表征流体体积模量的AVO近似公式为：

$\begin{array}{c}{R}_{\mathrm{PP}}\left(\mathrm{\&theta;}\right)=\left[(1-\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2})\frac{{\sec }^2\mathrm{\&theta;}}{4}\right]\frac{{\mathrm{\&Delta;K}}_f}{K_f}+[\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{{4\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\&theta;}-\frac{2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\&theta;}]\frac{\mathrm{\&Delta;}\left(f_m\right)}{f_m}+[\frac{1}{2}-\frac{{\sec }^2\mathrm{\&theta;}}{4}]\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}}\\ +(\frac{{\sec }^2\mathrm{\&theta;}}{4}-\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{{2\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\&theta;}+\frac{2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\&theta;})\frac{\mathrm{\&Delta;\&phi;}}{\mathrm{\&phi;}}\end{array}$

其中，f_m＝φμ为地下岩石的固体项，θ为入射角度，K_f为地下岩石的流体体积模量，f_m为地下岩石的固体项，φ为地下岩石的孔隙度，ρ为地下岩石的密度，μ为地下岩石的剪切模量，△K_f为界面两侧的流体体积模量的差值，△f_m为界面两侧固体项的差值，△ρ为界面两侧的密度的差值，△φ为界面两侧的孔隙度的差值，为干岩石纵横波速度比的平方，为饱和岩石纵横波速度比的平方。

6.根据权利要求5所述的方法，其特征是，根据所述的叠前地震资料、测井资料、地质成果资料以及表征流体体积模量的AVO近似公式提取子波具体包括：

从所述的测井资料中提取出地下岩石的密度以及孔隙度；

根据所述的剪切模量曲线以及所述的孔隙度确定地下岩石的固体项；

根据所述地下岩石的密度、固体项、孔隙度、流体体积模量曲线、表征流体体积模量的AVO近似公式、叠前CMP道集以及地震层位数据提取子波。

7.根据权利要求6所述的方法，其特征是，根据所述的叠前地震资料、测井资料、地质成果资料以及子波确定流体体积模量数据体具体包括：

利用所述的叠前CMP道集、测井资料、地质成果资料以及子波，基于贝叶斯理论，使用高斯分布函数为似然函数、柯西分布为先验约束建立AVO反演目标函数；

利用重加权最小二乘迭代算法求解所述的AVO反演目标函数，得到流体体积模量的相对变化量；

利用道积分思想将所述的相对变化量转化为流体体积模量数据体。

8.根据权利要求7所述的方法，其特征是，建立的AVO反演目标函数为：

(G^ΤG+ηQ+αC^ΤC)m＝G^Τd+αC^Τξ

$\mathrm{\&eta;}=\mathrm{\&lambda;}\frac{{\mathrm{\&sigma;}}_n^2}{{\mathrm{\&sigma;}}_m^2}$

$\mathrm{\&xi;}=\frac{1}{2}\ln \left(\frac{M}{M_0}\right)$

其中，为CMP道集资料中噪音的方差，为弹性参数相对变化量的方差，λ、α为加权系数，m为待反演的弹性参数相对变化量矩阵，G为正演算子矩阵，C为积分矩阵，d为CMP道集矩阵，Q为斜对角加权矩阵，M为模型约束矩阵，M₀为M的均值。

9.一种基于流体体积模量AVO反演的流体识别的系统，其特征是，所述的系统具体包括：

资料采集装置，用于采集目的储层的叠前地震资料、测井资料以及地质成果资料；

近似公式构建装置，用于构建表征流体体积模量的振幅随偏移距的变化AVO近似公式；

子波提取装置，用于根据所述的叠前地震资料、测井资料、地质成果资料以及所述的表征流体体积模量的AVO近似公式提取子波；

流体体积模量数据体确定装置，用于根据所述的叠前地震资料、测井资料、地质成果资料以及子波确定流体体积模量数据体；

流体类型识别装置，用于根据所述的流体体积模量数据体识别当前储层孔隙中流体的类型。

10.根据权利要求9所述的系统，其特征是，所述的系统还包括：

保幅处理装置，用于对所述的叠前地震资料进行保幅处理，得到共中心点CMP道集资料。

11.根据权利要求10所述的系统，其特征是，所述的系统还包括：

环境校正装置，用于对所述的测井资料进行环境校正；

奇异值消除装置，用于对经过环境校正后的测井资料进行奇异值消除；

纵横速度比确定装置，用于根据经过奇异值消除后的测井资料确定所述目的储层的干岩石纵横波速度比、饱和岩石纵横波速度比；

剪切模量曲线确定装置，用于根据所述的目的储层的干岩样纵横波速度比以及所述的测井资料确定与所述测井资料对应的流体体积模量曲线以及剪切模量曲线。

12.根据权利要求11所述的系统，其特征是，所述的系统还包括：

层位追踪装置，用于根据所述的叠前地震资料、测井资料、地质成果资料对目的储层进行精细层位追踪，得到地震层位数据。

13.根据权利要求12所述的系统，其特征是，所述的表征流体体积模量的AVO近似公式为：

$\begin{array}{c}{R}_{\mathrm{PP}}\left(\mathrm{\&theta;}\right)=\left[(1-\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2})\frac{{\sec }^2\mathrm{\&theta;}}{4}\right]\frac{{\mathrm{\&Delta;K}}_f}{K_f}+[\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{{4\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\&theta;}-\frac{2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\&theta;}]\frac{\mathrm{\&Delta;}\left(f_m\right)}{f_m}+[\frac{1}{2}-\frac{{\sec }^2\mathrm{\&theta;}}{4}]\frac{\mathrm{\&Delta;\&rho;}}{\mathrm{\&rho;}}\\ +(\frac{{\sec }^2\mathrm{\&theta;}}{4}-\frac{{\mathrm{\&gamma;}}_{\mathrm{dry}}^2}{{2\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sec }^2\mathrm{\&theta;}+\frac{2}{{\mathrm{\&gamma;}}_{\mathrm{sat}}^2}{\sin }^2\mathrm{\&theta;})\frac{\mathrm{\&Delta;\&phi;}}{\mathrm{\&phi;}}\end{array}$

其中，f_m＝φμ为地下岩石的固体项，θ为入射角度，K_f为地下岩石的流体体积模量，f_m为地下岩石的固体项，φ为地下岩石的孔隙度，ρ为地下岩石的密度，μ为地下岩石的剪切模量，△K_f为界面两侧的流体体积模量的差值，△f_m为界面两侧固体项的差值，△ρ为界面两侧的密度的差值，△φ为界面两侧的孔隙度的差值，为干岩石纵横波速度比的平方，为饱和岩石纵横波速度比的平方。

14.根据权利要求13所述的系统，其特征是，所述的子波提取装置具体包括：

孔隙度提取模块，用于从所述的测井资料中提取出地下岩石的密度以及孔隙度；

固体项确定模块，用于根据所述的剪切模量曲线以及所述的孔隙度确定地下岩石的固体项；

子波提取模块，用于根据所述地下岩石的密度、固体项、孔隙度、流体体积模量曲线、表征流体体积模量的AVO近似公式、叠前CMP道集以及地震层位数据提取子波。

15.根据权利要求14所述的系统，其特征是，所述的流体体积模量数据体确定装置具体包括：

反演目标函数建立模块，用于利用所述的叠前CMP道集、测井资料、地质成果资料以及子波，基于贝叶斯理论，使用高斯分布函数为似然函数、柯西分布为先验约束建立AVO反演目标函数；

相对变化量确定模块，用于利用重加权最小二乘迭代算法求解所述的AVO反演目标函数，得到流体体积模量的相对变化量；

转化模块，用于利用道积分思想将所述的相对变化量转化为流体体积模量数据体。

16.根据权利要求15所述的系统，其特征是，建立的AVO反演目标函数为：

(G^ΤG+ηQ+αC^ΤC)m＝G^Τd+αC^Τξ

$\mathrm{\&eta;}=\mathrm{\&lambda;}\frac{{\mathrm{\&sigma;}}_n^2}{{\mathrm{\&sigma;}}_m^2}$

$\mathrm{\&xi;}=\frac{1}{2}\ln \left(\frac{M}{M_0}\right)$

其中，为CMP道集资料中噪音的方差，为弹性参数相对变化量的方差，λ、α为加权系数，m为待反演的弹性参数相对变化量矩阵，G为正演算子矩阵，C为积分矩阵，d为CMP道集矩阵，Q为斜对角加权矩阵，M为模型约束矩阵，M₀为M的均值。