背景技术：

2.地震储层预测和流体识别技术主要是通过分析地震波的速度、振幅、相位、频率、波形等参数的变化来预测储集岩层的分布范围以及所含流体的性质等。岩性和储层物性的变化，造成了地震反射波速度、振幅、相位、频率、波形等的相应变化。这些变化是当今地震储层预测的主要依据。近年来，随着勘探程度的不断提高，岩性油气藏和隐蔽油气藏成为了油气田勘探开发的新领域，这类油气田勘探难度大，不仅要能有效的预测储层的空间分布，还须准确的识别储层中所含流体的性质。勘探对象的隐蔽、复杂，要求我们要深度挖掘地震资料，采用先进的现代信号处理技术，进行更精确、更高分辨率的储层预测和流体识别技术研究。随着频谱分析技术、地震各向异性成像、流体因子分析等地震资料处理和解释技术的不断发展，地震资料品质的不断提高，描述储层和流体特征参数的变化在地震剖面上的清晰度将越来越明显，可信度也将越来越高。另一方面，各种地震岩石物理模型和理论的不断发展和完善，也将使我们对地震波在含流体介质的传播机制具有更深入的理解，人们对弹性波速度频散和衰减与岩石物理性质之间关系的认识必将不断深化。3.avo油气预测技术在30年前就受到了业界的高度重视。80年代初，ostrander作为第一人开创性的介绍了具有avo效应的实际应用，该实际应用表明某个页岩底部的含天然气砂岩具有振幅随偏移距变化的性质。其内在原因为，岩石孔隙中的天然气会导致纵波速度会大大降低，但是并不会对横波速度造成影响，于是，纵横波速度的相对比例改变了，这种改变就会反应在不同炮检距的反射的振幅上。但是，伴随着avo技术的不断进化，学者们陆续发现随着炮检距的增大，含气砂岩的振幅并不是一定增加的，也会出现异常情况。因此，研究者开始使用avo交汇图的方法来对含油属性进行判断。ruthorford(1989)将砂岩分为三类，分类的标准是在垂直入射的情况下，反射系数的差别。castagna(1997)把砂岩分为四类，其划分的依据是交汇图的截距和梯度这两个性质。在其划分的四类砂岩中，只有一类砂岩的振幅是伴随炮检距的变大而变大的。为了完善avo理论，以及进一步挖掘avo技术的潜力，还有研究者给出了不同油气的检查因子和参数，其中包括smith-gidlow(1987)和fatti(1994)，他们给出了纵横波速度、阻抗、截距、梯度等参数，这些参数能够较为正确描述含油气层中的流体与其周围岩体性质的不同。以上这些参数均可用于对油气层的岩性以及流体较为准确的刻画。4.然而，常规的avo技术忽视了频率因素。但是研究表明，地震波的反射系数与频率有着不可分割的联系。地震反射系数与频率越相关，那么我们就越有可能观察到地震波通过含油气地层反射后的吸收衰减和速度发散等现象。介质中饱和流体与地震波衰减之间存在的这种关系表明，应该在以流体为主要影响因素的前提下来模拟地震波的衰减。chapman在2003年对thomsen和hudson提出的单尺度模型进行了改进和扩展，其核心思想是将单尺度变成了两种尺度。将双尺度裂隙下的流体互相影响纳入研究范围，在微裂隙孔隙基岩中含有中等尺度裂纹，裂纹长度远小于地震波长，但是同时比粒级尺度孔隙或微裂缝要大得多。chapman在这个模型中增加了一个长度特征变量，这个长度特征量又与中等规模裂缝非常关，并且与地震频率有着千丝万缕的关系。该模型对实际的储集层的介质更能进行准确的描述，根据该模型得到的正演记录，可知，avo与频率有不可分割的关系。该新理论是进行频变avo模拟分析的一个合理选择。5.实际地层介质通常是各向异性的，特别是对于非常规储层中的页岩储层。在vti介质中，各向同性假设条件下推导得到的精确zoeppritz方程计算得到的反射系数与vti介质精确反射系数方程计算得到的反射系数之间存在较大的偏差，并且zoeppritz方程不包含各向异性参数，因此精确zoeppritz方程在vti介质的反演中不再适用。同样，各向同性介质的线性反射系数近似公式也不再适用。6.因此，有必要开发一种vti介质叠前频变avo反演方法、装置、电子设备及介质。7.公开于本发明背景技术部分的信息仅仅旨在加深对本发明的一般背景技术的理解，而不应当被视为承认或以任何形式暗示该信息构成已为本领域技术人员所公知的现有技术。

技术实现要素：

8.本发明提出了一种叠前频变avo反演方法、装置、电子设备及介质，其能够通过对原始地震数据进行时频分析，来完成vti介质频变avo反演，并构建更灵敏的流体因子，可以直接得到地下流体的分布，假像更少，也减少了累计误差。9.第一方面，本公开实施例提供了一种叠前频变avo反演方法，包括：10.根据avo反演，获得纵波速度差异；11.通过频变avo反演，获得频变流体项；12.根据所述纵波速度差异与所述频变流体项，计算绝对纵波速度差频变流体因子。13.优选地，通过频变avo反演，获得频变流体项包括：14.通过频变avo反演，获得频变流体项的表达式；15.根据贝叶斯理论，将子波加入分块矩阵，获得频变流体矩阵；16.通过基于cauchy约束的贝叶斯频变avo反演理论，建立目标函数，进而获得最终反演方程；17.采用迭代重加权最小二乘优化算法求解所述最终反演方程，获得频变流体项。18.优选地，频变流体项的表达式为：[0019][0020]其中，为频变流体项，[0021]优选地，所述频变流体矩阵为：[0022][0023]其中，wj(j＝1,2...,k)表示频率为fj的子波矩阵，的子波矩阵，[0024]优选地，所述目标函数为：[0025][0026]其中，为噪音的方差，为需要求得的模型参数m的方差。[0027]优选地，所述最终反演方程为：[0028]m＝(gtg+λq)-1gtδdꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ(4)[0029]其中，[0030]优选地，通过公式(5)计算绝对纵波速度差频变流体因子：[0031][0032]其中，apdfff为绝对纵波速度差频变流体因子，δvp为纵波速度差异，为频变流体项。[0033]作为本公开实施例的一种具体实现方式，[0034]第二方面，本公开实施例还提供了一种叠前频变avo反演装置，包括：[0035]纵波速度差异计算模块，根据avo反演，获得纵波速度差异；[0036]频变流体项计算模块，通过频变avo反演，获得频变流体项；[0037]绝对纵波速度差频变流体因子计算模块，根据所述纵波速度差异与所述频变流体项，计算绝对纵波速度差频变流体因子。[0038]优选地，通过频变avo反演，获得频变流体项包括：[0039]通过频变avo反演，获得频变流体项的表达式；[0040]根据贝叶斯理论，将子波加入分块矩阵，获得频变流体矩阵；[0041]通过基于cauchy约束的贝叶斯频变avo反演理论，建立目标函数，进而获得最终反演方程；[0042]采用迭代重加权最小二乘优化算法求解所述最终反演方程，获得频变流体项。[0043]优选地，频变流体项的表达式为：[0044][0045]其中，为频变流体项，[0046]优选地，所述频变流体矩阵为：[0047][0048]其中，wj(j＝1,2...,k)表示频率为fj的子波矩阵，的子波矩阵，[0049]优选地，所述目标函数为：[0050][0051]其中，为噪音的方差，为需要求得的模型参数m的方差。[0052]优选地，所述最终反演方程为：[0053]m＝(gtg+λq)-1gtδdꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ(4)[0054]其中，优选地，通过公式(5)计算绝对纵波速度差频变流体因子：[0055][0056]其中，apdfff为绝对纵波速度差频变流体因子，δvp为纵波速度差异，为频变流体项。[0057]第三方面，本公开实施例还提供了一种电子设备，该电子设备包括：[0058]存储器，存储有可执行指令；[0059]处理器，所述处理器运行所述存储器中的所述可执行指令，以实现所述的叠前频变avo反演方法。[0060]第四方面，本公开实施例还提供了一种计算机可读存储介质，该计算机可读存储介质存储有计算机程序，该计算机程序被处理器执行时实现所述的叠前频变avo反演方法。[0061]其有益效果在于：。[0062]本发明的方法和装置具有其它的特性和优点，这些特性和优点从并入本文中的附图和随后的具体实施方式中将是显而易见的，或者将在并入本文中的附图和随后的具体实施方式中进行详细陈述，这些附图和具体实施方式共同用于解释本发明的特定原理。附图说明[0063]通过结合附图对本发明示例性实施例进行更详细的描述，本发明的上述以及其它目的、特征和优势将变得更加明显，其中，在本发明示例性实施例中，相同的参考标号通常代表相同部件。[0064]图1示出了根据本发明的一个实施例的叠前频变avo反演方法的步骤的流程图。[0065]图2示出了根据本发明的一个实施例的利用广义s变换方法的反演结果的示意图。[0066]图3示出了根据本发明的一个实施例的利用稀疏约束反演谱分解方法的反演结果的示意图。[0067]图4示出了根据本发明的一个实施例的一种叠前频变avo反演装置的框图。[0068]附图标记说明：[0069]201、纵波速度差异计算模块；202、频变流体项计算模块；203、绝对纵波速度差频变流体因子计算模块。具体实施方式[0070]下面将更详细地描述本发明的优选实施方式。虽然以下描述了本发明的优选实施方式，然而应该理解，可以以各种形式实现本发明而不应被这里阐述的实施方式所限制。[0071]本发明提供一种叠前频变avo反演方法，包括：[0072]根据avo反演，获得纵波速度差异。[0073]具体地，基于贝叶斯理论的传统avo反演获得纵横波速度，进而得到纵波速度差。[0074]通过频变avo反演，获得频变流体项；在一个示例中，通过频变avo反演，获得频变流体项包括：通过频变avo反演，获得频变流体项的表达式；根据贝叶斯理论，将子波加入分块矩阵，获得频变流体矩阵；通过基于cauchy约束的贝叶斯频变avo反演理论，建立目标函数，进而获得最终反演方程；采用迭代重加权最小二乘优化算法求解最终反演方程，获得频变流体项。[0075]在一个示例中，频变流体项的表达式为：[0076][0077]其中，为频变流体项，[0078]在一个示例中，频变流体矩阵为：[0079][0080]其中，wj(j＝1,2...,k)表示频率为fj的子波矩阵，的子波矩阵，[0081]在一个示例中，目标函数为：[0082][0083]其中，为噪音的方差，为需要求得的模型参数m的方差。[0084]在一个示例中，最终反演方程为：[0085]m＝(gtg+λq)-1gtδdꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ(4)[0086]其中，[0087]其中，apdfff为绝对纵波速度差频变流体因子，δvp为纵波速度差异，为频变流体项。[0088]具体地，假设地层呈现vti弱各项异性；实际地震数据的入射角度范围是中小角度；实际地层无强吸收衰减现象。[0089]rüger(1997)在前人研究的基础上，基于弹性参数为弱不连续性和弱各向异性假设，推导出了vti介质的反射系数近似公式。若所有各向异性参数相同为零，则rppvti(θ)可化简为各向同性反射系数的shuey的近似表达式。式中ε和δ为thomsen各向异性参数中的两个。[0090][0091]其中，z＝ρvp是垂直p波阻抗，g＝ρvs2表示垂直剪切模量。[0092]对rüger(1997)推导的反射系数公式进行整理，将公式整理为纵波速度项，纵波阻抗项，剪切模量项和两个各向异性参数项。之后对公式进行泰勒展开，进而进行频变avo反演。[0093][0094]在参考频率f0处求导可得：[0095][0096]其中，d(θ)δε+e(θ)δδ＝b(f)ani。[0097]整理上式为：[0098][0099]进而得到三个频变流体项：[0100][0101]同上面的频变反演步骤类似，根据贝叶斯理论，将子波加入分块矩阵，得到公式(2)，将公式(2)简写成矩阵形式：[0102]δd＝gmꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ(11)[0103]假设待反演模型参数满足柯西分布，在贝叶斯参数估计框架下，利用基于cauchy约束的贝叶斯频变avo反演理论建立目标函数为公式(3)。[0104]在求解过程中为了使得目标函数达到最小，求取目标函数f(m)对模型参数m的梯度，并让其等于0，就可以得到最终反演方程为公式(4)，选用迭代重加权最小二乘(irls)优化算法求解，即可得到频变流体项[0105]根据纵波速度差异与频变流体项，计算绝对纵波速度差频变流体因子。在一个示例中，通过公式(5)计算绝对纵波速度差频变流体因子：[0106][0107]具体地，基于传统avo反演获得的弹性参数构建的流体因子可以指示流体的存在，但精度有限。当地下地层含流体时，纵波速度会发生变化，造成含流体储层位置与周围岩层的纵波速度差，所以纵波速度差可以作为流体指示因子，但由于反演结果精度的限制，速度差识别流体效果不佳。频变avo反演可以提取地震数据中由于流体的存在而引起的频散信息，因此，频散特性也可以用于识别流体。如果把传统avo反演获得的纵波速度差异与频变avo反演获得的频变属性结合起来，就可以结合两种反演用于识别流体的优势，进一步锁定油气优势区域，提高利用地震数据识别流体的精度。因此，通过公式(5)计算绝对纵波速度差频变流体因子。[0108]本发明还提供一种叠前频变avo反演装置，包括：[0109]纵波速度差异计算模块，根据avo反演，获得纵波速度差异。[0110]具体地，基于贝叶斯理论的传统avo反演获得纵横波速度，进而得到纵波速度差。[0111]频变流体项计算模块，通过频变avo反演，获得频变流体项；在一个示例中，通过频变avo反演，获得频变流体项包括：通过频变avo反演，获得频变流体项的表达式；根据贝叶斯理论，将子波加入分块矩阵，获得频变流体矩阵；通过基于cauchy约束的贝叶斯频变avo反演理论，建立目标函数，进而获得最终反演方程；采用迭代重加权最小二乘优化算法求解最终反演方程，获得频变流体项。[0112]在一个示例中，频变流体项的表达式为：[0113][0114]其中，为频变流体项，[0115]在一个示例中，频变流体矩阵为：[0116][0117]其中，wj(j＝1,2...,k)表示频率为fj的子波矩阵，的子波矩阵，[0118]在一个示例中，目标函数为：[0119][0120]其中，为噪音的方差，为需要求得的模型参数m的方差。[0121]在一个示例中，最终反演方程为：[0122]m＝(gtg+λq)-1gtδdꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ(4)[0123]其中，[0124]其中，apdfff为绝对纵波速度差频变流体因子，δvp为纵波速度差异，为频变流体项。[0125]具体地，假设地层呈现vti弱各项异性；实际地震数据的入射角度范围是中小角度；实际地层无强吸收衰减现象。[0126]rüger(1997)在前人研究的基础上，基于弹性参数为弱不连续性和弱各向异性假设，推导出了vti介质的反射系数近似公式。若所有各向异性参数相同为零，则rppvti(θ)可化简为各向同性反射系数的shuey的近似表达式。式中ε和δ为thomsen各向异性参数中的两个。[0127]对rüger(1997)推导的反射系数公式进行整理，将公式整理为纵波速度项，纵波阻抗项，剪切模量项和两个各向异性参数项。之后对公式进行泰勒展开，进而进行频变avo反演。[0128]在参考频率f0处求导可得公式(8)，整理上式为公式(9)，进而得到三个频变流体项为公式(10)。[0129]同上面的频变反演步骤类似，根据贝叶斯理论，将子波加入分块矩阵，得到公式(2)，将公式(2)简写成矩阵形式为公式(11)。[0130]假设待反演模型参数满足柯西分布，在贝叶斯参数估计框架下，利用基于cauchy约束的贝叶斯频变avo反演理论建立目标函数为公式(3)。[0131]在求解过程中为了使得目标函数达到最小，求取目标函数f(m)对模型参数m的梯度，并让其等于0，就可以得到最终反演方程为公式(4)，选用迭代重加权最小二乘(irls)优化算法求解，即可得到频变流体项[0132]绝对纵波速度差频变流体因子计算模块，根据纵波速度差异与频变流体项，计算绝对纵波速度差频变流体因子。在一个示例中，通过公式(5)计算绝对纵波速度差频变流体因子：[0133][0134]具体地，基于传统avo反演获得的弹性参数构建的流体因子可以指示流体的存在，但精度有限。当地下地层含流体时，纵波速度会发生变化，造成含流体储层位置与周围岩层的纵波速度差，所以纵波速度差可以作为流体指示因子，但由于反演结果精度的限制，速度差识别流体效果不佳。频变avo反演可以提取地震数据中由于流体的存在而引起的频散信息，因此，频散特性也可以用于识别流体。如果把传统avo反演获得的纵波速度差异与频变avo反演获得的频变属性结合起来，就可以结合两种反演用于识别流体的优势，进一步锁定油气优势区域，提高利用地震数据识别流体的精度。因此，通过公式(5)计算绝对纵波速度差频变流体因子。[0135]本发明还提供一种电子设备，电子设备包括：存储器，存储有可执行指令；处理器，处理器运行存储器中的可执行指令，以实现上述的叠前频变avo反演方法。[0136]本发明还提供一种计算机可读存储介质，该计算机可读存储介质存储有计算机程序，该计算机程序被处理器执行时实现上述的叠前频变avo反演方法。[0137]为便于理解本发明实施例的方案及其效果，以下给出四个具体应用示例。本领域技术人员应理解，该示例仅为了便于理解本发明，其任何具体细节并非意在以任何方式限制本发明。[0138]实施例1[0139]图1示出了根据本发明的一个实施例的叠前频变avo反演方法的步骤的流程图。[0140]如图1所示，该叠前频变avo反演方法包括：步骤101，根据avo反演，获得纵波速度差异；步骤102，通过频变avo反演，获得频变流体项；步骤103，根据纵波速度差异与频变流体项，计算绝对纵波速度差频变流体因子。[0141]基于贝叶斯理论的传统avo反演获得纵横波速度，进而得到纵波速度差。[0142]通过频变avo反演，获得频变流体项的表达式为公式(1)；根据贝叶斯理论，将子波加入分块矩阵，获得频变流体矩阵为公式(2)；通过基于cauchy约束的贝叶斯频变avo反演理论，建立目标函数为公式(3)，进而获得最终反演方程为公式(4)；采用迭代重加权最小二乘优化算法求解最终反演方程，获得频变流体项。[0143]根据纵波速度差异与频变流体项，通过公式(5)计算绝对纵波速度差。[0144]图2示出了根据本发明的一个实施例的利用广义s变换方法的反演结果的示意图，其中，stack是叠加剖面，df是纵波频散流体因子，apdfff为新构建的绝对纵波速度差流体因子。从图2可以看出，借助广义s变换对原始地震信号进行时频分析时自动改变窗口大小的自适应性，可以较好地从地震记录提取各个频率的信息，结合基于vti介质的反演算法，计算出的纵波频散流体因子df可以较好的指示地下流体的分布。但是还有一些非速度频散造成的假象，比如动校拉伸和近地表的强吸收衰减作用。新构建的流体因子apdfff结合传统的速度差异大小，将速度差异大小的影响因素一起考虑。这样就减少了部分假象，突显了真正地下速度频散的地方。从图中也可以看出，apdfff的剖面能量显示更集中，更真实。[0145]图3示出了根据本发明的一个实施例的利用稀疏约束反演谱分解方法的反演结果的示意图，其中，stack是叠加剖面，df是纵波频散流体因子，apdfff为新构建的绝对纵波速度差流体因子。图3对比图2，使用了基于稀疏约束反演谱分解的时频分析方法，得到了更为稀疏，能量更聚焦的时频分析结果。使用这样的时频域地震记录进行频变avo反演，得到的结果更为精确。从图中可以看出，与图2的纵波流体因子df相比，图3的df能量更聚焦，位置更精确。再考虑绝对纵波速度差异的贡献，使得结果新流体因子apdfff误差更小，能更准确的指示vti地下介质的流体分布。[0146]实施例2[0147]图4示出了根据本发明的一个实施例的一种叠前频变avo反演装置的框图。[0148]如图4所示，该叠前频变avo反演装置，包括：[0149]纵波速度差异计算模块201，根据avo反演，获得纵波速度差异；[0150]频变流体项计算模块202，通过频变avo反演，获得频变流体项；[0151]绝对纵波速度差频变流体因子计算模块203，根据纵波速度差异与频变流体项，计算绝对纵波速度差频变流体因子。[0152]作为可选方案，通过频变avo反演，获得频变流体项包括：[0153]通过频变avo反演，获得频变流体项的表达式；[0154]根据贝叶斯理论，将子波加入分块矩阵，获得频变流体矩阵；[0155]通过基于cauchy约束的贝叶斯频变avo反演理论，建立目标函数，进而获得最终反演方程；[0156]采用迭代重加权最小二乘优化算法求解最终反演方程，获得频变流体项。[0157]作为可选方案，频变流体项的表达式为：[0158][0159]其中，为频变流体项，作为可选方案，频变流体矩阵为：[0160][0161]其中，wj(j＝1,2...,k)表示频率为fj的子波矩阵，的子波矩阵，[0162]作为可选方案，目标函数为：[0163][0164]其中，为噪音的方差，为需要求得的模型参数m的方差。[0165]作为可选方案，最终反演方程为：[0166]m＝(gtg+λq)-1gtδdꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ(4)[0167]其中，作为可选方案，通过公式(5)计算绝对纵波速度差频变流体因子：[0168][0169]其中，apdfff为绝对纵波速度差频变流体因子，δvp为纵波速度差异，为频变流体项。[0170]实施例3[0171]本公开提供一种电子设备包括，该电子设备包括：存储器，存储有可执行指令；处理器，处理器运行存储器中的可执行指令，以实现上述叠前频变avo反演方法。[0172]根据本公开实施例的电子设备包括存储器和处理器。[0173]该存储器用于存储非暂时性计算机可读指令。具体地，存储器可以包括一个或多个计算机程序产品，该计算机程序产品可以包括各种形式的计算机可读存储介质，例如易失性存储器和/或非易失性存储器。该易失性存储器例如可以包括随机存取存储器(ram)和/或高速缓冲存储器(cache)等。该非易失性存储器例如可以包括只读存储器(rom)、硬盘、闪存等。[0174]该处理器可以是中央处理单元(cpu)或者具有数据处理能力和/或指令执行能力的其它形式的处理单元，并且可以控制电子设备中的其它组件以执行期望的功能。在本公开的一个实施例中，该处理器用于运行该存储器中存储的该计算机可读指令。[0175]本领域技术人员应能理解，为了解决如何获得良好用户体验效果的技术问题，本实施例中也可以包括诸如通信总线、接口等公知的结构，这些公知的结构也应包含在本公开的保护范围之内。[0176]有关本实施例的详细说明可以参考前述各实施例中的相应说明，在此不再赘述。[0177]实施例4[0178]本公开实施例提供一种计算机可读存储介质，该计算机可读存储介质存储有计算机程序，该计算机程序被处理器执行时实现所述的叠前频变avo反演方法。[0179]根据本公开实施例的计算机可读存储介质，其上存储有非暂时性计算机可读指令。当该非暂时性计算机可读指令由处理器运行时，执行前述的本公开各实施例方法的全部或部分步骤。[0180]上述计算机可读存储介质包括但不限于：光存储介质(例如：cd－rom和dvd)、磁光存储介质(例如：mo)、磁存储介质(例如：磁带或移动硬盘)、具有内置的可重写非易失性存储器的媒体(例如：存储卡)和具有内置rom的媒体(例如：rom盒)。[0181]本领域技术人员应理解，上面对本发明的实施例的描述的目的仅为了示例性地说明本发明的实施例的有益效果，并不意在将本发明的实施例限制于所给出的任何示例。[0182]以上已经描述了本发明的各实施例，上述说明是示例性的，并非穷尽性的，并且也不限于所披露的各实施例。在不偏离所说明的各实施例的范围和精神的情况下，对于本技术领域的普通技术人员来说许多修改和变更都是显而易见的。