页岩气储层脆性地震预测方法及装置
阅读：10发布：2021-12-19
专利汇可以提供页岩气储层脆性地震预测方法及装置专利检索，专利查询，专利分析的服务。并且本 发明 提供了一种 页岩 气储层脆性 地震 预测方法及装置，该方法包括：利用一工区的实际地震数据和多个弹性参数的先验分布函数及Zoeppritz方程，基于Bayesian原理，借助广义线性反演思想求解得到所述弹性参数的扰动量的求解表达式；根据所述求解表达式和所述弹性参数的初始模型计算得到所述弹性参数的扰动量值；根据所述扰动量值和所述初始模型计算得到所述弹性参数的当前模型；基于所述实际地震数据、所述当前模型及所述多个弹性参数的Zoeppritz方程计算得到地震反演残差；若所述地震反演残差小于一设定 阈值 ，将所述当前模型作为所述弹性参数的值，根据所述弹性参数的值计算得到所述工区的页岩气储层脆性指数。本发明能够提高页岩气储层脆性指数预测的准确度。，下面是页岩气储层脆性地震预测方法及装置专利的具体信息内容。
1.一种页岩气储层脆性地震预测方法，其特征在于，包括：
利用一工区的实际地震数据和多个弹性参数的先验分布函数及Zoeppritz方程，基于Bayesian原理，借助广义线性反演思想求解得到所述弹性参数的扰动量的求解表达式，所述多个弹性参数包括杨氏模量、泊松比及密度；
根据所述求解表达式和所述弹性参数的初始模型计算得到所述弹性参数的扰动量值；
根据所述扰动量值和所述初始模型计算得到所述弹性参数的当前模型；
基于所述实际地震数据、所述当前模型及所述多个弹性参数的Zoeppritz方程计算得到地震反演残差；
若所述地震反演残差小于一设定阈值，将所述当前模型作为所述弹性参数的值，根据所述弹性参数的值计算得到所述工区的页岩气储层脆性指数。
2.如权利要求1所述的页岩气储层脆性地震预测方法，其特征在于，利用一工区的实际地震数据和多个弹性参数的先验分布函数及Zoeppritz方程，基于Bayesian原理，借助广义线性反演思想求解得到所述弹性参数的扰动量的求解表达式之前，还包括：
利用所述工区的测井数据统计分析得到多个弹性参数的先验分布函数。
3.如权利要求1所述的页岩气储层脆性地震预测方法，其特征在于，利用一工区的实际地震数据和多个弹性参数的先验分布函数及Zoeppritz方程，基于Bayesian原理，借助广义线性反演思想求解得到所述弹性参数的扰动量的求解表达式之前，还包括：
对初始Zoeppritz方程进行等价变换，得到所述多个弹性参数的Zoeppritz方程。
4.如权利要求1所述的页岩气储层脆性地震预测方法，其特征在于，根据所述扰动量值和所述初始模型计算得到所述弹性参数的当前模型之前，还包括：
利用地震构造解释数据和所述工区的测井数据建立所述弹性参数的初始模型。
5.如权利要求1所述的页岩气储层脆性地震预测方法，其特征在于，基于所述实际地震数据、所述当前模型及所述多个弹性参数的Zoeppritz方程计算得到地震反演残差，包括：
基于所述实际地震数据提取得到地震子波；
利用所述地震子波、所述初始模型、所述多个弹性参数的Zoeppritz方程及所述实际地震数据，计算得到地震反演残差。
6.如权利要求1所述的页岩气储层脆性地震预测方法，其特征在于，还包括：
若所述地震反演残差大于或等于所述设定阈值，则根据所述求解表达式和所述当前模型重新计算所述弹性参数的扰动量值，并根据重新计算的扰动量值和所述当前模型重新计算所述弹性参数的当前模型，以及基于所述实际地震数据、重新计算的当前模型及所述多个弹性参数的Zoeppritz方程重新计算地震反演残差，如此重复进行直到地震反演残差小于所述设定阈值，并将此时重新计算的当前模型作为所述弹性参数的值。
7.如权利要求5所述的页岩气储层脆性地震预测方法，其特征在于，基于所述实际地震数据提取得到地震子波，包括：
通过统计方法从所述实际地震数据的地震叠前道集提取波形依赖于地震波入射角度的多个所述地震子波。
8.如权利要求5所述的页岩气储层脆性地震预测方法，其特征在于，利用所述地震子波、所述初始模型、所述多个弹性参数的Zoeppritz方程及所述实际地震数据，计算得到地震反演残差，包括：
以所述初始模型作为输入，利用所述多个弹性参数的Zoeppritz方程计算各地震入射角度的地震反射系数；
将所述地震子波和所述地震反射系数进行褶积，得到所述工区的模拟的地震角道集；
将所述工区的模拟的地震角道集与所述实际地震数据中的地震道角集作差，得到所述地震反演残差。
9.如权利要求5所述的页岩气储层脆性地震预测方法，其特征在于，利用所述地震子波、所述初始模型、所述多个弹性参数的Zoeppritz方程及所述实际地震数据，计算得到地震反演残差之前，还包括：
利用所述工区内生产井的测井数据、初始Zoeppritz方程及所述实际地震数据三者中的一个或多个对所述地震子波进行修正，修正后的地震子波用于计算所述地震反演残差。
10.如权利要求9所述的页岩气储层脆性地震预测方法，其特征在于，利用所述工区内生产井的测井数据、初始Zoeppritz方程及所述实际地震数据三者中的一个或多个对所述地震子波进行修正，包括：
以所述测井数据作为输入模型，利用初始Zoeppritz方程正演得到模拟的地震角道集；
根据所述模拟的地震角道集和所述实际地震数据中的所述生产井周围的地震角道集，计算得到地震振幅缩放因子；
利用所述地震振幅缩放因子对所述地震子波进行振幅缩放处理。
11.如权利要求2所述的页岩气储层脆性地震预测方法，其特征在于，利用所述工区的测井数据统计分析得到多个弹性参数的先验分布函数，包括：
基于所述弹性参数服从柯西分布的假设和所述工区内所有生产井的测井数据，计算得到所述弹性参数的均值和所述多个弹性参数的扰动量的协方差矩阵，作为所述先验分布函数。
12.如权利要求3所述的页岩气储层脆性地震预测方法，其特征在于，对初始Zoeppritz方程进行等价变换，得到所述多个弹性参数的Zoeppritz方程，包括：
基于岩石物理理论和地层介质的各向同性假设，对初始Zoeppritz方程进行等价变换，得到所述多个弹性参数的Zoeppritz方程。
13.如权利要求4所述的页岩气储层脆性地震预测方法，其特征在于，利用地震构造解释数据和所述工区的测井数据建立所述弹性参数的初始模型，包括：
依据所述地震构造解释数据进行层位散点插值，得到地质模型；
根据所述地质模型和所述测井数据进行弹性参数横向插值，得到所述初始模型。
14.如权利要求1所述的页岩气储层脆性地震预测方法，其特征在于，利用一工区的实际地震数据和多个弹性参数的先验分布函数及Zoeppritz方程，基于Bayesian原理，借助广义线性反演思想求解得到所述弹性参数的扰动量的求解表达式，包括：
利用所述先验分布函数、所述实际地震数据、所述初始模型、所述多个弹性参数的Zoeppritz方程所对应的正演算子，基于Bayesian原理，构建得到最大后验概率意义下的关于所述弹性参数的目标函数；
利用广义线性反演思想对所述目标函数进行求解，得到所述求解表达式。
15.如权利要求1所述的页岩气储层脆性地震预测方法，其特征在于，所述实际地震数据包括PP波和PS波，所述地震反演残差为包含PP波和PS波的多波的地震联合反演残差。
16.一种页岩气储层脆性地震预测装置，其特征在于，包括：
扰动量表达式生成单元，用于利用一工区的实际地震数据和多个弹性参数的先验分布函数及Zoeppritz方程，基于Bayesian原理，借助广义线性反演思想求解得到所述弹性参数的扰动量的求解表达式，所述多个弹性参数包括杨氏模量、泊松比及密度；
扰动量值生成单元，用于根据所述求解表达式和所述弹性参数的初始模型计算得到所述弹性参数的扰动量值；
当前模型生成单元，用于根据所述扰动量值和所述初始模型计算得到所述弹性参数的当前模型；
反演残差生成单元，用于基于所述实际地震数据、所述当前模型及所述多个弹性参数的Zoeppritz方程计算得到地震反演残差；
脆性指数生成单元，用于若所述地震反演残差小于一设定阈值，将所述当前模型作为所述弹性参数的值，根据所述弹性参数的值计算得到所述工区的页岩气储层脆性指数。
17.如权利要求16所述的页岩气储层脆性地震预测装置，其特征在于，还包括：
先验分布生成单元，用于利用所述工区的测井数据统计分析得到多个弹性参数的先验分布函数。
18.如权利要求16所述的页岩气储层脆性地震预测装置，其特征在于，还包括：
Zoeppritz方程变换单元，用于对初始Zoeppritz方程进行等价变换，得到所述多个弹性参数的Zoeppritz方程。
19.如权利要求16所述的页岩气储层脆性地震预测装置，其特征在于，还包括：
初始模型生成单元，用于利用地震构造解释数据和所述工区的测井数据建立所述弹性参数的初始模型。
20.如权利要求16所述的页岩气储层脆性地震预测装置，其特征在于，所述反演残差生成单元，包括：
地震子波提取模块，用于基于所述实际地震数据提取得到地震子波；
反演残差生成模块，用于利用所述地震子波、所述初始模型、所述多个弹性参数的Zoeppritz方程及所述实际地震数据，计算得到地震反演残差。
21.如权利要求16所述的页岩气储层脆性地震预测装置，其特征在于，还包括：
弹性参数值更新单元，用于若所述地震反演残差大于或等于所述设定阈值，则根据所述求解表达式和所述当前模型重新计算所述弹性参数的扰动量值，并根据重新计算的扰动量值和所述当前模型重新计算所述弹性参数的当前模型，以及基于所述实际地震数据、重新计算的当前模型及所述多个弹性参数的Zoeppritz方程重新计算地震反演残差，如此重复进行直到地震反演残差小于所述设定阈值，并将此时重新计算的当前模型作为所述弹性参数的值。
22.如权利要求20所述的页岩气储层脆性地震预测装置，其特征在于，所述地震子波提取模块，包括：
地震子波筛选模块，用于通过统计方法从所述实际地震数据的地震叠前道集提取波形依赖于地震波入射角度的多个所述地震子波。
23.如权利要求20所述的页岩气储层脆性地震预测装置，其特征在于，所述反演残差生成模块，包括：
地震反射系数生成模块，用于以所述初始模型作为输入，利用所述多个弹性参数的Zoeppritz方程计算各地震入射角度的地震反射系数；
模拟地震角道集生成模块，用于将所述地震子波和所述地震反射系数进行褶积，得到所述工区的模拟的地震角道集；
地震反演残差确定模块，用于将所述工区的模拟的地震角道集与所述实际地震数据中的地震道角集作差，得到所述地震反演残差。
24.如权利要求20所述的页岩气储层脆性地震预测装置，其特征在于，所述反演残差生成单元，还包括：
地震子波修正模块，用于利用所述工区内生产井的测井数据、初始Zoeppritz方程及所述实际地震数据三者中的一个或多个对所述地震子波进行修正，修正后的地震子波用于计算所述地震反演残差。
25.如权利要求24所述的页岩气储层脆性地震预测装置，其特征在于，所述地震子波修正模块，包括：
模拟地震角道集正演模块，用于以所述测井数据作为输入模型，利用初始Zoeppritz方程正演得到模拟的地震角道集；
地震振幅缩放因子生成模块，用于根据所述模拟的地震角道集和所述实际地震数据中的所述生产井周围的地震角道集，计算得到地震振幅缩放因子；
地震振幅缩放模块，用于利用所述地震振幅缩放因子对所述地震子波进行振幅缩放处理。
26.如权利要求17所述的页岩气储层脆性地震预测装置，其特征在于，所述先验分布生成单元，包括：
先验分布生成模块，用于基于所述弹性参数服从柯西分布的假设和所述工区内所有生产井的测井数据，计算得到所述弹性参数的均值和所述多个弹性参数的扰动量的协方差矩阵，作为所述先验分布函数。
27.如权利要求18所述的页岩气储层脆性地震预测装置，其特征在于，所述Zoeppritz方程变换单元，包括：
Zoeppritz方程变换模块，用于基于岩石物理理论和地层介质的各向同性假设，对初始Zoeppritz方程进行等价变换，得到所述多个弹性参数的Zoeppritz方程。
28.如权利要求19所述的页岩气储层脆性地震预测装置，其特征在于，所述初始模型生成单元，包括：
地质模型生成模块，用于依据所述地震构造解释数据进行层位散点插值，得到地质模型；
初始模型生成模块，用于根据所述地质模型和所述测井数据进行弹性参数横向插值，得到所述初始模型。
29.如权利要求16所述的页岩气储层脆性地震预测装置，其特征在于，所述扰动量表达式生成单元，包括：
目标函数生成模块，用于利用所述先验分布函数、所述实际地震数据、所述初始模型、所述多个弹性参数的Zoeppritz方程所对应的正演算子，基于Bayesian原理，构建得到最大后验概率意义下的关于所述弹性参数的目标函数；
扰动量表达式生成模块，用于利用广义线性反演思想对所述目标函数进行求解，得到所述求解表达式。
30.如权利要求16所述的页岩气储层脆性地震预测装置，其特征在于，所述扰动量表达式生成单元，还用于执行：所述实际地震数据包括PP波和PS波；反演残差生成单元，还用于执行：所述地震反演残差为包含PP波和PS波的多波的地震联合反演残差。
页岩气储层脆性地震预测方法及装置
技术领域
[0001] 本发明涉及油气勘探开发技术领域，尤其涉及一种页岩气储层脆性地震预测方法及装置。
背景技术
[0002] 页岩气作为全球重要接替资源，其勘探开发技术成为目前人们的研究热点和难点。脆性是岩石重要的力学特征之一，是页岩气储层可压裂性的关键参数，根据页岩气储层的可压裂性，可以确定有利井位或“甜点”位置。
[0003] 通常，利用脆性矿物含量估算岩石脆性最为直观而且准确，但必须以岩样测试数据和测井分析为基础，并且局部的岩样测试结果和单口井的分析结果不能客观反映周围岩石的脆性分布特征。当储层岩样测试数据和测井数据有限时，基于岩样测试和测井分析的脆性矿物含量估算结果横向分辨率低，难于实现整个工区储层可压裂性评价。
[0004] 岩石力学参数杨氏模量和泊松比是重要的岩石脆性指示因子，在非常规油气藏地层脆性预测中发挥着重要作用，而获取岩石力学参数的有效途径之一就是叠前地震反演。地震数据具有较好的空间连续性，在储层表征中起着非常重要的作用。高精度反演方法使反演结果具有重要价值，可以提高地震数据解释的可信度。利用反演结果可以更好地实现储层空间描述，很大程度地提高储层表征的价值。
[0005] 通常，杨氏模量和泊松比是通过由叠前反演得到的纵波速度、横波速度以及密度间接计算或者直接基于Zoeppritz方程的近似公式反演而得到。然而，间接计算存在累计误差，大大降低了杨氏模量和泊松比的计算精度，并且基于Zoeppritz方程近似公式的反演方法在利用大偏移距地震振幅(大角度)信息和纵横波联合反演等方面均存在较大问题，严重制约着叠前反演精度的提高。因此，上述获取杨氏模量和泊松比的途径在计算精度上不能满足实际油藏储层岩石脆性精细表征要求。
发明内容
[0006] 本发明提供一种页岩气储层脆性地震预测方法及装置，以提高页岩气储层脆性地震预测精度。
[0007] 本发明提供一种页岩气储层脆性地震预测方法，包括：利用一工区的实际地震数据和多个弹性参数的先验分布函数及Zoeppritz方程，基于Bayesian原理，借助广义线性反演思想求解得到所述弹性参数的扰动量的求解表达式；根据所述求解表达式和所述弹性参数的初始模型计算得到所述弹性参数的扰动量值；根据所述扰动量值和所述初始模型计算得到所述弹性参数的当前模型；基于所述实际地震数据、所述当前模型及所述多个弹性参数的Zoeppritz方程计算得到地震反演残差；若所述地震反演残差小于一设定阈值，将所述当前模型作为所述弹性参数的值，根据所述弹性参数的值计算得到所述工区的页岩气储层脆性指数。
[0008] 一个实施例中，利用一工区的实际地震数据和多个弹性参数的先验分布函数及Zoeppritz方程，基于Bayesian原理，借助广义线性反演思想求解得到所述弹性参数的扰动量的求解表达式之前，还包括：利用所述工区的测井数据统计分析得到多个弹性参数的先验分布函数。
[0009] 一个实施例中，利用一工区的实际地震数据和多个弹性参数的先验分布函数及Zoeppritz方程，基于Bayesian原理，借助广义线性反演思想求解得到所述弹性参数的扰动量的求解表达式之前，还包括：对初始Zoeppritz方程进行等价变换，得到所述多个弹性参数的Zoeppritz方程。
[0010] 一个实施例中，根据所述扰动量值和所述初始模型计算得到所述弹性参数的当前模型之前，还包括：利用地震构造解释数据和所述工区的测井数据建立所述弹性参数的初始模型。
[0011] 一个实施例中，基于所述实际地震数据、所述当前模型及所述多个弹性参数的Zoeppritz方程计算得到地震反演残差，包括：基于所述实际地震数据提取得到地震子波；利用所述地震子波、所述初始模型、所述多个弹性参数的Zoeppritz方程及所述实际地震数据，计算得到地震反演残差。
[0012] 一个实施例中，还包括：若所述地震反演残差大于或等于所述设定阈值，则根据所述求解表达式和所述当前模型重新计算所述弹性参数的扰动量值，并根据重新计算的扰动量值和所述当前模型重新计算所述弹性参数的当前模型，以及基于所述实际地震数据、重新计算的当前模型及所述多个弹性参数的Zoeppritz方程重新计算地震反演残差，如此重复进行直到地震反演残差小于所述设定阈值，并将此时重新计算的当前模型作为所述弹性参数的值。
[0014] 一个实施例中，利用所述地震子波、所述初始模型、所述多个弹性参数的Zoeppritz方程及所述实际地震数据，计算得到地震反演残差，包括：以所述初始模型作为输入，利用所述多个弹性参数的Zoeppritz方程计算各地震入射角度的地震反射系数；将所述地震子波和所述地震反射系数进行褶积，得到所述工区的模拟的地震角道集；将所述工区的模拟的地震角道集与所述实际地震数据中的地震道角集作差，得到所述地震反演残差。
[0015] 一个实施例中，利用所述地震子波、所述初始模型、所述多个弹性参数的Zoeppritz方程及所述实际地震数据，计算得到地震反演残差之前，还包括：利用所述工区内生产井的测井数据、初始Zoeppritz方程及所述实际地震数据三者中的一个或多个对所述地震子波进行修正，修正后的地震子波用于计算所述地震反演残差。
[0016] 一个实施例中，利用所述工区内生产井的测井数据、初始Zoeppritz方程及所述实际地震数据三者中的一个或多个对所述地震子波进行修正，包括：以所述测井数据作为输入模型，利用初始Zoeppritz方程正演得到模拟的地震角道集；根据所述模拟的地震角道集和所述实际地震数据中的所述生产井周围的地震角道集，计算得到地震振幅缩放因子；利用所述地震振幅缩放因子对所述地震子波进行振幅缩放处理。
[0017] 一个实施例中，利用所述工区的测井数据统计分析得到多个弹性参数的先验分布函数，包括：基于所述弹性参数服从柯西分布的假设和所述工区内所有生产井的测井数据，计算得到所述弹性参数的均值和所述多个弹性参数的扰动量的协方差矩阵，作为所述先验分布函数。
[0018] 一个实施例中，对初始Zoeppritz方程进行等价变换，得到所述多个弹性参数的Zoeppritz方程，包括：基于岩石物理理论和地层介质的各向同性假设，对初始Zoeppritz方程进行等价变换，得到所述多个弹性参数的Zoeppritz方程。
[0019] 一个实施例中，利用地震构造解释数据和所述工区的测井数据建立所述弹性参数的初始模型，包括：依据所述地震构造解释数据进行层位散点插值，得到地质模型；根据所述地质模型和所述测井数据进行弹性参数横向插值，得到所述初始模型。
[0020] 一个实施例中，利用一工区的实际地震数据和多个弹性参数的先验分布函数及Zoeppritz方程，基于Bayesian原理，借助广义线性反演思想求解得到所述弹性参数的扰动量的求解表达式，包括：利用所述先验分布函数、所述实际地震数据、所述初始模型、所述多个弹性参数的Zoeppritz方程所对应的正演算子，基于Bayesian原理，构建得到最大后验概率意义下的关于所述弹性参数的目标函数；利用广义线性反演思想对所述目标函数进行求解，得到所述求解表达式。
[0021] 一个实施例中，所述实际地震数据包括PP波和PS波，所述地震反演残差为包含PP波和PS波的多波的地震联合反演残差。
[0022] 本发明还提供一种页岩气储层脆性地震预测装置，包括：扰动量表达式生成单元，用于利用一工区的实际地震数据和多个弹性参数的先验分布函数及Zoeppritz方程，基于Bayesian原理，借助广义线性反演思想求解得到所述弹性参数的扰动量的求解表达式；扰动量值生成单元，用于根据所述求解表达式和所述弹性参数的初始模型计算得到所述弹性参数的扰动量值；当前模型生成单元，用于根据所述扰动量值和所述初始模型计算得到所述弹性参数的当前模型；反演残差生成单元，用于基于所述实际地震数据、所述当前模型及所述多个弹性参数的Zoeppritz方程计算得到地震反演残差；脆性指数生成单元，用于若所述地震反演残差小于一设定阈值，将所述当前模型作为所述弹性参数的值，根据所述弹性参数的值计算得到所述工区的页岩气储层脆性指数。
[0023] 一个实施例中，还包括：先验分布生成单元，用于利用所述工区的测井数据统计分析得到多个弹性参数的先验分布函数。
[0024] 一个实施例中，还包括：Zoeppritz方程变换单元，用于对初始Zoeppritz方程进行等价变换，得到所述多个弹性参数的Zoeppritz方程。
[0025] 一个实施例中，还包括：初始模型生成单元，用于利用地震构造解释数据和所述工区的测井数据建立所述弹性参数的初始模型。
[0026] 一个实施例中，所述反演残差生成单元，包括：地震子波提取模块，用于基于所述实际地震数据提取得到地震子波；反演残差生成模块，用于利用所述地震子波、所述初始模型、所述多个弹性参数的Zoeppritz方程及所述实际地震数据，计算得到地震反演残差。
[0027] 一个实施例中，还包括：弹性参数值更新单元，用于若所述地震反演残差大于或等于所述设定阈值，则根据所述求解表达式和所述当前模型重新计算所述弹性参数的扰动量值，并根据重新计算的扰动量值和所述当前模型重新计算所述弹性参数的当前模型，以及基于所述实际地震数据、重新计算的当前模型及所述多个弹性参数的Zoeppritz方程重新计算地震反演残差，如此重复进行直到地震反演残差小于所述设定阈值，并将此时重新计算的当前模型作为所述弹性参数的值。
[0028] 一个实施例中，所述地震子波提取模块，包括：地震子波筛选模块，用于通过统计方法从所述实际地震数据的地震叠前道集提取波形依赖于地震波入射角度的多个所述地震子波。
[0029] 一个实施例中，所述反演残差生成模块，包括：地震反射系数生成模块，用于以所述初始模型作为输入，利用所述多个弹性参数的Zoeppritz方程计算各地震入射角度的地震反射系数；模拟地震角道集生成模块，用于将所述地震子波和所述地震反射系数进行褶积，得到所述工区的模拟的地震角道集；地震反演残差确定模块，用于将所述工区的模拟的地震角道集与所述实际地震数据中的地震道角集作差，得到所述地震反演残差。
[0030] 一个实施例中，所述反演残差生成单元，还包括：地震子波修正模块，用于利用所述工区内生产井的测井数据、初始Zoeppritz方程及所述实际地震数据三者中的一个或多个对所述地震子波进行修正，修正后的地震子波用于计算所述地震反演残差。
[0031] 一个实施例中，所述地震子波修正模块，包括：模拟地震角道集正演模块，用于以所述测井数据作为输入模型，利用初始Zoeppritz方程正演得到模拟的地震角道集；地震振幅缩放因子生成模块，用于根据所述模拟的地震角道集和所述实际地震数据中的所述生产井周围的地震角道集，计算得到地震振幅缩放因子；地震振幅缩放模块，用于利用所述地震振幅缩放因子对所述地震子波进行振幅缩放处理。
[0032] 一个实施例中，所述先验分布生成单元，包括：先验分布生成模块，用于基于所述弹性参数服从柯西分布的假设和所述工区内所有生产井的测井数据，计算得到所述弹性参数的均值和所述多个弹性参数的扰动量的协方差矩阵，作为所述先验分布函数。
[0033] 一个实施例中，所述Zoeppritz方程变换单元，包括：Zoeppritz方程变换模块，用于基于岩石物理理论和地层介质的各向同性假设，对初始Zoeppritz方程进行等价变换，得到所述多个弹性参数的Zoeppritz方程。
[0034] 一个实施例中，所述初始模型生成单元，包括：地质模型生成模块，用于依据所述地震构造解释数据进行层位散点插值，得到地质模型；初始模型生成模块，用于根据所述地质模型和所述测井数据进行弹性参数横向插值，得到所述初始模型。
[0035] 一个实施例中，所述扰动量表达式生成单元，包括：目标函数生成模块，用于利用所述先验分布函数、所述实际地震数据、所述初始模型、所述多个弹性参数的Zoeppritz方程所对应的正演算子，基于Bayesian原理，构建得到最大后验概率意义下的关于所述弹性参数的目标函数；扰动量表达式生成模块，用于利用广义线性反演思想对所述目标函数进行求解，得到所述求解表达式。
[0036] 一个实施例中，所述扰动量表达式生成单元，还用于执行：所述实际地震数据包括PP波和PS波；反演残差生成单元，还用于执行：所述地震反演残差为包含PP波和PS波的多波的地震联合反演残差。
[0037] 本发明的页岩气储层脆性地震预测装置，基于关于多个弹性参数的Zoeppritz方程可直接通过叠前地震数据反演得到能够直接表征岩石脆性的岩石力学参数，能够有效避免由现有弹性参数间接计算岩石力学参数带来的累积误差。相对于现有的Zoeppritz近似方程而言，本发明地震叠前AVO反演基于关于多个弹性参数的Zoeppritz方程能够克服因入射角(或偏移距)过大而带来计算误差的问题，从而使得我们能够利用大偏移距信息进行储层参数预测。本发明关于多个弹性参数的精确Zoeppritz方程的AVO反演方法，具有更好的适应性，能够满足泥砂岩薄互层结构，纵向弹性参数差异大的储层条件，如页岩气储层，有效提高叠前反演精度。本发明可以获得弹性参数的较佳模型，以此计算岩石脆性指数，能够提供高可靠性的页岩气“甜点”地震识别方法，相比于通过岩样测试和测井分析估算工区储层的岩石脆性分布特征，预测“甜点”位置的方法，能提供横向连续性更好，且更加真实的反映地下页岩气储层脆性分布特征脆性指数数据体。附图说明
[0038] 为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。在附图中：
[0039] 图1是本发明一实施例的页岩气储层脆性地震预测方法的流程示意图；
[0040] 图2是本发明另一实施例的页岩气储层脆性地震预测方法的流程示意图；
[0041] 图3是本发明又一实施例的页岩气储层脆性地震预测方法的流程示意图；
[0042] 图4是本发明再一实施例的页岩气储层脆性地震预测方法的流程示意图；
[0043] 图5是本发明一实施例中计算地震反演残差方法的流程示意图；
[0044] 图6是本发明再一实施例的页岩气储层脆性地震预测方法的流程示意图；
[0045] 图7是本发明一实施例中计算地震反演残差方法的流程示意图；
[0046] 图8是本发明另一实施例中计算地震反演残差方法的流程示意图；
[0047] 图9是本发明一实施例中对地震子波进行修正方法的流程示意图；
[0048] 图10是本发明一实施例中建立弹性参数的初始模型的流程示意图；
[0049] 图11是本发明一实施例中计算得到弹性参数的扰动量求解表达式方法的流程示意图；
[0050] 图12是本发明又一实施例的页岩气储层脆性地震预测方法的流程示意图；
[0051] 图13是本发明一实施例中输入的不含噪地震叠前AVO PP波角道集地震数据；
[0052] 图14是本发明一实施例中输入的不含噪叠前AVO PS波角道集地震数据；
[0053] 图15是本发明一实施例中输入的含噪叠前AVO PP波角道集地震数据；
[0054] 图16是本发明一实施例中输入的含噪叠前AVO PS波角道集地震数据；
[0055] 图17是本发明一实施例中无噪情况下多波AVO反演得到的弹性参数；
[0056] 图18是本发明一实施例中含噪情况下多波AVO反演得到的弹性参数；
[0057] 图19和图20分别是本发明一实施例中不含噪情况下和含噪情况下计算得到的岩石脆性指数；
[0058] 图21是本发明一实施例的页岩气储层脆性地震预测装置的结构示意图；
[0059] 图22是本发明另一实施例的页岩气储层脆性地震预测装置的结构示意图；
[0060] 图23是本发明又一实施例的页岩气储层脆性地震预测装置的结构示意图；
[0061] 图24是本发明再一实施例的页岩气储层脆性地震预测装置的结构示意图；
[0062] 图25是本发明一实施例中反演残差生成单元的结构示意图；
[0063] 图26是本发明再一实施例的页岩气储层脆性地震预测装置的结构示意图；
[0064] 图27是本发明一实施例中反演残差生成模块的结构示意图；
[0065] 图28是本发明另一实施例中反演残差生成单元的结构示意图；
[0066] 图29是本发明一实施例中地震子波修正模块的结构示意图；
[0067] 图30是本发明一实施例中初始模型生成单元的结构示意图；
[0068] 图31是本发明一实施例中扰动量表达式生成单元的结构示意图；
[0069] 图32是本发明另一实施例的页岩气储层脆性地震预测装置的结构示意图。
具体实施方式
[0070] 为使本发明实施例的目的、技术方案和优点更加清楚明白，下面结合附图对本发明实施例做进一步详细说明。在此，本发明的示意性实施例及其说明用于解释本发明，但并不作为对本发明的限定。
[0071] 发明人经过创造性劳动发现，目前基于叠前地震数据的页岩气储层脆性地震预测方法研究存在以下问题：
[0072] (1)传统弹性参数，如纵波速度、横波速度等不能很好的表征岩石脆性。
[0073] (2)通过纵波速度、横波速度和密度间接计算得到的杨氏模量和泊松比存在累积误差。
[0074] (3)油气藏储层特别是页岩气储层最常见的结构是泥砂岩薄互层结构，纵向弹性参数差异较大，极大地降低了Zoeppritz近似公式的求解精度，基于Zoeppritz方程近似表达式的AVO反演方法不再适用于这类储层。
[0075] (4)目前常规的AVO反演方法主要是基于近似公式对杨氏模量、泊松比、密度等岩石力学参数的反射系数进行反演，然后通过反演得到的对应参数的反射系数间接地求解得到，降低了反演的精度。
[0076] (5)传统的叠前AVO反演方法一般只考虑PP波数据，而PS波数包含大量的岩石力学参数信息，单独利用PP地震AVO数据反演杨氏模量、泊松比具有较大的不确定性。
[0077] 针对现有的杨氏模量和泊松比是由叠前反演得到的纵波速度、横波速度以及密度间接计算获得，或者直接基于Zoeppritz方程的近似公式反演获得，这种间接计算导致累计误差，降低了杨氏模量和泊松比的计算精度的问题，且基于近似公式的反演存在一系列的精度问题，发明人提出了一种页岩气储层脆性地震预测方法，创造性地由现有Zoeppritz方程得到关于多个弹性参数的Zoeppritz方程，并基于新的Zoeppritz方程进行页岩气储层脆性地震预测，能够提高预测精度。并基于现有技术中存在的诸多其他问题提出各实施例的方法以进一步提高预测精度。
[0078] 图1是本发明一实施例的页岩气储层脆性地震预测方法的流程示意图。如图1所示，本发明实施例的页岩气储层脆性地震预测方法，可包括步骤：
[0079] S110：利用一工区的实际地震数据和多个弹性参数的先验分布函数及Zoeppritz方程，基于Bayesian原理，借助广义线性反演思想求解得到上述弹性参数的扰动量的求解表达式；
[0080] S120：根据上述求解表达式和上述弹性参数的初始模型计算得到上述弹性参数的扰动量值；
[0081] S130：根据上述扰动量值和上述初始模型计算得到上述弹性参数的当前模型；
[0082] S140：基于上述实际地震数据、上述当前模型及上述多个弹性参数的Zoeppritz方程计算得到地震反演残差；
[0083] S150：若上述地震反演残差小于一设定阈值，将上述当前模型作为上述弹性参数的值，根据上述弹性参数的值计算得到上述工区的页岩气储层脆性指数。
[0084] 在上述步骤S110中，上述弹性参数可以指介质的岩石力学参数，可以包括杨氏模量、泊松比及密度等参数。
[0085] 该实际地震数据可以是通过地震测量的实际记录，可以是实际角度域地震道集。该实际地震数据可以包括多种地震波，例如PP波和PS波，PS波包含大量的岩石力学信息，基于包含PP波和PS波的实际地震数据计算弹性参数的扰动量，有利于后续得到确定性更高的弹性参数。
[0086] 多个弹性参数的先验分布函数及Zoeppritz方程具体是指多个弹性参数的先验分布函数和关于多个弹性参数的Zoeppritz方程，关于多个弹性参数的Zoeppritz方程可以是以弹性力学参数为变量的新形式的Zoeppritz方程，可以通过对现有的Zoeppritz方程进行改写或变换得到。
[0087] 现有的Zoeppritz方程不是关于上述多个弹性参数的方程，可以是关于横波波速、纵波波速及密度的方程。现有的Zoeppritz方程可以具有多种形式，例如Zoeppritz方程近似表达式、精确(exact)的Zoeppritz方程，可如后述公式(3)所示。
[0088] 在上述步骤S120中，上述弹性参数的初始模型可以是初始弹性参数模型，可通过多种不同方法得到，例如，利用地震构造解释资料和测井数据经过建模得到。
[0089] 在上述步骤S130中，该当前模型可以是上述弹性参数的最新模型，可通过在初始模型上增加一定的扰动量得到。
[0090] 在上述步骤S140中，可以利用上述当前模型和上述多个弹性参数的Zoeppritz方程得到模拟的地震数据，再利用模拟的地震数据和实际地震数据计算得到地震反演残差。其中，地震反演残差是基于弹性参数的Zoeppritz方程计算得到，可以克服现有的Zoeppritz方程因入射角或偏移距过大带来的计算误差的问题，从而能够利用大偏移距信息进行储层参数预测，可以提高岩石脆性指数预测结果的准确度。
[0091] 一个实施例中，上述实际地震数据可包括PP波和PS波，利用上述当前模型和上述多个弹性参数的Zoeppritz方程所得到的模拟的地震数据可包括PP波和PS波，则该地震反演残差可以为包含PP波和PS波的多波的地震联合反演残差。
[0092] 在上述步骤S150中，该设定阈值可以根据需要进行设置。可以采用多种方法根据上述弹性参数的值计算上述工区的页岩气储层脆性指数，例如，可以是Rickman于2008年提出的基于弹性参数的脆性指数(BrittlenessIndex，BI)计算公式，或者是根据工区实际情况对Rickman于2008年提出的基于弹性参数的脆性指数BI计算公式进行合理的修正，然后利用修正后的计算公式根据上述弹性参数的值计算整个工区的岩石脆性指数，以此可以预测有利压裂区，可为“甜点”位置的确定提供重要的依据。
[0093] 一个实施例中，若上述地震反演残差不小于该设定阈值，可基于当前模型重复步骤S120～S140重新计算，得到新的地震反演残差，并再次判断新的地震反演残差是否小于该设定阈值，如此迭代进行，直到地震反演残差小于该设定阈值。通过地震反演残差控制迭代计算的次数，有利于获得弹性参数的最优模型，从而能够获得最优的岩石脆性指数预测结果。
[0094] 本发明实施例中，基于关于多个弹性参数的Zoeppritz方程计算得到上述弹性参数的扰动量的求解表达式，并基于上述求解表达式计算得到上述弹性参数的扰动量值，再计算得到上述弹性参数的当前模型，并基于上述当前模型及上述多个弹性参数的Zoeppritz方程计算得到地震反演残差，无需利用由叠前反演得到的纵波速度、横波速度以及密度间接获得或者基于现有Zoeppritz方程的近似公式反演获得弹性参数的模型，避免产生累计误差，提高了弹性参数的计算精度。相对于现有的Zoeppritz近似方程而言，关于多个弹性参数的Zoeppritz方程，能够克服了因入射角(或偏移距)过大而带来计算误差的问题，从而使得能够利用大偏移距信息进行储层参数预测。通过地震反演残差控制是否将计算得到的当前模型作为弹性参数的最终模型，有利于得到弹性参数的最优模型，利用该最优模型预测岩石脆性指数，可以提高预测精度。
[0095] 进一步，本发明实施例中，上述实际地震数据可包括PP波和PS波，该地震反演残差可以为包含PP波和PS波，联合PP波与PS波进行多个弹性参数反演，与现有通常仅考虑PP波，本实施例中不仅考虑了PP波，还考虑了PS波，由于PS波地震数据包含更多岩石物理属性和密度信息，能有效提高反演的稳定性。
[0096] 图2是本发明另一实施例的页岩气储层脆性地震预测方法的流程示意图。如图2所示，图1所示的页岩气储层脆性地震预测方法，利用一工区的实际地震数据和多个弹性参数的先验分布函数及Zoeppritz方程，基于Bayesian原理，借助广义线性反演思想求解得到上述弹性参数的扰动量的求解表达式(步骤S110)之前，还可包括步骤：
[0097] S160：利用上述工区的测井数据统计分析得到多个弹性参数的先验分布函数。
[0098] 在上述步骤S160中，较佳地，可以基于工区内所有测井数据统计弹性参数的先验分布函数，以此利用充足的测井数据可以得到更准确的先验分布函数。该先验分布函数可以包括各弹性参数的均值及所有弹性参数的扰动量相关性的协方差矩阵。上述多个弹性参数可包括杨氏模量、泊松比和密度。
[0099] 本实施例中，除了根据上述扰动量值和上述初始模型计算得到上述弹性参数的当前模型，还利用工区的测井数据获得地下介质地震弹性参数模型的信息作为先验信息，可以降低步骤S140中地震反演的不确定性，提高反演过程的稳定性，提高地震反演残差的精度，进而提高页岩气储层脆性指数预测的准确度。
[0100] 一个实施例中，在上述步骤S160中，利用上述工区的测井数据统计分析得到多个弹性参数的先验分布函数的方法，具体实施方式可为：基于上述弹性参数服从柯西分布的假设和上述工区内所有生产井的测井数据，计算得到上述弹性参数的均值和上述多个弹性参数的扰动量的协方差矩阵，作为上述先验分布函数。
[0101] 本实施例中，可通过测井数据分析获得需要的弹性参数，在利用多变量柯西分布函数求取各弹性参数的自相关系数和互相关系数，并构建多个弹性参数扰动量的协方差矩阵，形成符合该工区的弹性参数的模型参数的先验分布函数。本实施例中的先验分布函数包括多个弹性参数的扰动量的协方差矩阵，且基于上述弹性参数服从柯西分布的假设获得，利用该先验分布函数计算弹性参数的扰动量的求解表达式，可以很好地消除待反演弹性参数扰动量之间的相关性，从而降低反演的不确定性，提高反演的精度，进而提高脆性指数预测准确度。
[0102] 一个具体实施例中，上述多个弹性参数可包括杨氏模量、泊松比及密度三个参数，以此为例，上述多变量柯西分布函数可以是三变量柯西分布函数。可采用该三变量柯西分布函数作为上述三个弹性参数的先验分布关系式，并可基于工区内所有测井数据统计得到的各弹性参数的模型信息，求取各弹性参数的均值、各弹性参数的自相关系数及三个弹性参数的互相关系数，进而构建三弹性参数的扰动量的协方差矩阵，利用求取的均值和协方差矩阵形成该工区的弹性参数的先验分布函数。
[0103] 上述三变量柯西分布函数所对应的规则化表达式可为：
[0104]
[0105] 其中，R(m)表示柯西分布函数，m＝[E,σ,den]T为弹性参数的模型参数，E为杨氏模量，σ为泊松比，den为密度，μ为由各弹性参数的均值构成的向量，N是弹性参数的模型参数的长度。矩阵Φi＝(Di)TΨ-1Di，其中，Ψ可以是3×3的协方差矩阵，可包含三个待反演的弹性参数之间统计相关性，可用于在计算过程中消除三个弹性参数的统计相关性。由于研究算法最终是转化为求解三个弹性参数的扰动量，因此在实际计算过程中协方差矩阵Ψ为扰动量的协方差矩阵，用于去除不同参数扰动量之间的相关性，使反演结果变得更加稳定。矩阵D可以是3×3N的矩阵，具体形式可如下所示：
[0106]
[0107] 图3是本发明又一实施例的页岩气储层脆性地震预测方法的流程示意图。如图3所示，图1所示的页岩气储层脆性地震预测方法，利用一工区的实际地震数据和多个弹性参数的先验分布函数及Zoeppritz方程，基于Bayesian原理，借助广义线性反演思想求解得到上述弹性参数的扰动量的求解表达式(步骤S110)之前，还可包括步骤：
[0108] S170：对初始Zoeppritz方程进行等价变换，得到上述多个弹性参数的Zoeppritz方程。
[0109] 在上述步骤S170中，该初始Zoeppritz方程可以是包括近似方程在内的多种现有形式，较佳地是精确(exact)的Zoeppritz方程。得到的上述多个弹性参数的Zoeppritz方程是新形式的Zoeppritz方程。
[0110] 一个实施例中，在上述步骤S170中，对初始Zoeppritz方程进行等价变换，得到上述多个弹性参数的Zoeppritz方程的方法，具体实施方式可为：基于岩石物理理论和地层介质的各向同性假设，对初始Zoeppritz方程进行等价变换，得到上述多个弹性参数的Zoeppritz方程。
[0111] 本实施例中，根据岩石物理理论，对现有形式的Zoeppritz方程进行等价改写，可得到关于多个弹性参数的新形式的Zoeppritz方程。相对于Zoeppritz近似方程而言，能够克服了因入射角(或偏移距)过大而带来计算误差的问题，基于该种关于多个弹性参数的Zoeppritz方程在上述步骤S140中进行的AVO反演，具有更好的岩石适应性，能够满足泥砂岩薄互层结构，纵向弹性参数差异大的储层条件，如页岩气储层，能够有效提高叠前反演精度，进而提高脆性指数预测准确度。
[0112] 具体地，例如，首先基于各向同性假设，利用杨氏模量、泊松比、纵横波速度以及密度之间的关系，对传统形式的Zoeppritz方程进行等价转换，最终得到关于杨氏模量、泊松比和密度这三个弹性参数的新形式的Zoeppritz方程。
[0113] 一个具体实施例中，上述多个弹性参数可以包括杨氏模量、泊松比和密度三个参数，以此为例，根据岩石物理理论将初始Zoeppritz方程转化为基于杨氏模量、泊松比和密度的新形式的Zoeppritz方程的过程可如下所述：
[0114] 首先，初始Zoeppritz方程描述平面波入射到两个半无限均匀各向同性弹性介质的分界面时的反射和透射情况，其表达式可如下所示：
[0115]
[0116] 其中，下标1和2表示分界面两侧的介质，θ1为入射角， 为反射角，θ2、 为透射角，VP1、VP2为纵波波速，VS1、VS2为横波波速，RPP为PP波反射系数，RPS为PS波反射系数，TPP为PP波透射系数，TPS为PS波透射系数，ρ1、ρ2表示介质密度。
[0117] 然后，根据Snell定律有：
[0118]
[0119] 利用公式(4)将公式(3)中的反射角 和透射角θ2、 用入射角θ1表示，公式(3)可变换为如下形式：
[0120]
[0121] 根据岩石物理理论，在均匀各向同性完全弹性无孔介质中，纵波速度Vp1、Vp2和横波速度Vs1、Vs2可以分别表示为：
[0122]
[0123] 其中，M1和M2分别表示分界面两侧介质中的纵波模量，μ1和μ2分别表示分界面两侧介质中的横波模量。
[0124] 将公式(6)代入公式(5)可得到基于纵波模量、横波模量和密度的Zoeppritz方程表达式：
[0125]
[0126] 在各向同性介质中，纵波模量和横波模量与杨氏模量和泊松比之间的关系可如下所示：
[0127]
[0128] 其中，E1和E2分别表示分界面两侧介质的杨氏模量，σ1和σ2分别表示分界面两侧介质的泊松比，M1和M2分别表示分界面两侧介质中的纵波模量，μ1和μ2分别表示分界面两侧介质中的横波模量。
[0129] 假设用参数η1、η2、τ1、τ2代替泊松比σ1和σ2表示如下：
[0130]
[0131] 则公式(8)可变换为：
[0132]
[0133] 将公式(10)代入公式(7)便可得到关于杨氏模量、泊松比和密度三个弹性参数的Zoeppritz方程，如下所示：
[0134]
[0135] 基于关于杨氏模量、泊松比和密度三个弹性参数的Zoeppritz方程的正演过程可表示如下：
[0136]
[0137] 其中，下标pp表示PP波，下标ps表示PS波，G表示关于杨氏模量、泊松比和密度三个弹性参数的Zoeppritz方程所对应的正演算子，W表示地震子波的矩阵。
[0138] 图4是本发明再一实施例的页岩气储层脆性地震预测方法的流程示意图。如图4所示，图1所示的页岩气储层脆性地震预测方法中，根据上述扰动量值和上述初始模型计算得到上述弹性参数的当前模型(步骤S130)之前，还可包括步骤：
[0139] S180：利用地震构造解释数据和上述工区的测井数据建立上述弹性参数的初始模型。
[0140] 本实施例中，利用地震构造解释数据和上述工区的测井数据可以通过多种方法建立上述弹性参数的初始模型。例如，可利用三维空间插值方法建立弹性参数的初始模型。具体地，可利用地震构造解释资料，基于沉积模式建立地质模型，并将测井资料，按构造模式进行插值和外推，得到每条测线的初始弹性参数模型，初始弹性参数模型可包括杨氏模量模型、泊松比模型和密度模型，以此可以快速建立弹性参数的初始模。
[0141] 图5是本发明一实施例中计算地震反演残差方法的流程示意图。如图5所示，在上述步骤S140中，基于上述实际地震数据、上述当前模型及上述多个弹性参数的Zoeppritz方程计算得到地震反演残差的方法，可包括步骤：
[0142] S141：基于上述实际地震数据提取得到地震子波；
[0143] S142：利用上述地震子波、上述初始模型、上述多个弹性参数的Zoeppritz方程及上述实际地震数据，计算得到地震反演残差。
[0144] 在上述步骤S141中，该实际地震数据可以是实际的地震叠前道集，可以包括角度域的PP波道集和PS波道集。一个实施例中，可基于实际的地震叠前道集和测井数据采取统计方法提取地震子波，该测井数据可以对基于上述实际地震数据提取得到地震子波进行修正，以提高提取地震子波的稳定性。
[0145] 在上述步骤S142中，上述初始模型可以是时间域的初始弹性参数模型。可利用该初始模型、基于多个弹性参数的Zoeppritz方程及上述实际地震数据计算得到工区的模拟的地震角道集，例如角度域模拟的PP波道集和PS波道集，再由模拟的地震角道集和角度域的实际地震道集，计算得到例如PP波与PS波的地震反演残差。
[0146] 一个实施例中，在上述步骤S141中，基于上述实际地震数据提取得到地震子波的方法，具体实施方式可为：通过统计方法从上述实际地震数据的地震叠前道集提取波形依赖于地震波入射角度的多个上述地震子波。
[0147] 本实施例中，地震子波在传播过程中受地层的影响会发生波形或频率变化，提取依赖于入射角度的地震子波，能有效提高振幅匹配程度。
[0148] 图6是本发明再一实施例的页岩气储层脆性地震预测方法的流程示意图。如图6所示，图1所示的页岩气储层脆性地震预测方法，还可包括步骤：
[0149] S190：若上述地震反演残差大于或等于上述设定阈值，则根据上述求解表达式和上述当前模型重新计算上述弹性参数的扰动量值，并根据重新计算的扰动量值和上述当前模型重新计算上述弹性参数的当前模型，以及基于上述实际地震数据、重新计算的当前模型及上述多个弹性参数的Zoeppritz方程重新计算地震反演残差，如此重复进行直到地震反演残差小于上述设定阈值，并将此时重新计算的当前模型作为上述弹性参数的值。
[0150] 在上述步骤S190中，根据弹性参数扰动量的求解表达式公式和反演残差，可采用迭代重加权最小二乘算法计算弹性参数的扰动量，重复迭代上述步骤S120～S140，并通过反演残差控制最大迭代次数，可以获得最优模型参数，进而提高岩石脆性指数的预测精度。
[0151] 图7是本发明一实施例中计算地震反演残差方法的流程示意图。如图7所示，在上述步骤S142中，利用上述地震子波、上述初始模型、上述多个弹性参数的Zoeppritz方程及上述实际地震数据，计算得到地震反演残差的方法，具体地，可包括步骤：
[0152] S1421：以上述初始模型作为输入，利用上述多个弹性参数的Zoeppritz方程计算各地震入射角度的地震反射系数；
[0153] S1422：将上述地震子波和上述地震反射系数进行褶积，得到上述工区的模拟的地震角道集；
[0154] S1423：将上述工区的模拟的地震角道集与上述实际地震数据中的地震道角集作差，得到上述地震反演残差。
[0155] 在上述步骤S1421中，可以时间域的初始模型作为输入，利用上述多个弹性参数的Zoeppritz方程计算每一个角度地震波例如PP波与PS波的反射系数向量。
[0156] 在上述步骤S1422中，将上述地震子波和上述地震反射系数进行褶积，可得到上述工区的模拟的地震角道集，例如角度域的PP波角道集与PS波角道集。
[0157] 在上述步骤S1423中，将模拟的地震角道集(例如模拟的PP波角道集与PS波角道集)与实际的地震角道集(例如实际的PP波角道集与PS波角道集)作差，可得到上述反演残差(例如PP波与PS波的多波联合反演残差)。本实施例的反演残差弹性参数的初始模型计算得到，可以有效控制上述迭代次数，获取弹性参数的最优模型。
[0158] 图8是本发明另一实施例中计算地震反演残差方法的流程示意图。如图8所示，图5所示的计算地震反演残差方法，利用上述地震子波、上述初始模型、上述多个弹性参数的Zoeppritz方程及上述实际地震数据，计算得到地震反演残差(步骤S142)之前，还可包括步骤：
[0159] S143：利用上述工区内生产井的测井数据、初始Zoeppritz方程及上述实际地震数据三者中的一个或多个对上述地震子波进行修正，修正后的地震子波用于计算上述地震反演残差。
[0160] 本实施例中，可以通过多种方式单独或一同对上述地震子波进行修正。例如，可以仅利用上述工区内生产井的测井数据对上述地震子波进行修正，或者，可以上述工区内生产井的测井数据、初始Zoeppritz方程及上述实际地震数据三者求得地震振幅缩放因子，在利用地震振幅缩放因子对上述地震子波进行修正。修正后的地震子波有利于使上述模拟的地震数据和上述实际地震数据具有更好的匹配度。
[0161] 图9是本发明一实施例中对地震子波进行修正方法的流程示意图。如图9所示，在上述步骤S143中，利用上述工区内生产井的测井数据、初始Zoeppritz方程及上述实际地震数据三者中的一个或多个对上述地震子波进行修正的方法，具体地，可包括步骤：
[0162] S1431：以上述测井数据作为输入模型，利用初始Zoeppritz方程正演得到模拟的地震角道集；
[0163] S1432：根据上述模拟的地震角道集和上述实际地震数据中的上述生产井周围的地震角道集，计算得到地震振幅缩放因子；
[0164] S1433：利用上述地震振幅缩放因子对上述地震子波进行振幅缩放处理。
[0165] 在上述步骤S1431中，该初始Zoeppritz方程可以是现有精确的Zoeppritz方程。正演得到模拟的地震角道集可以包括角度域的PP波角道集和PS波角道集。
[0166] 在上述步骤S1432～S1433中，一个实施例中，实际地震数据的信噪比较高，可为角道集的每一道使用统一的振幅缩放因子，以保证振幅随偏移距的变化关系。另一实施例中，实际地震数据的信噪比低，可分别利用近偏移距、中偏移距、远偏移距计算地震振幅缩放因子，以保证模拟记录(模拟的地震数据)与实际记录(实际的地震数据)的最佳振幅匹配，减少噪声对反演过程的影响。
[0167] 本实施例中，实际的地震振幅往往是相对值，采用初始Zoeppritz方程正演模拟的地震数据振幅与实际地震振幅存在一定数值差异。利用地震振幅缩放因子对提取的地震子波进行振幅缩放处理，能够使得模拟的地震数据的振幅与实际的地震数据具有较好的振幅匹配度。
[0168] 图10是本发明一实施例中建立弹性参数的初始模型的流程示意图。如图10所示，在上述步骤S180中，利用地震构造解释数据和上述工区的测井数据建立上述弹性参数的初始模型的方法，具体地，可包括步骤：
[0169] S181：依据上述地震构造解释数据进行层位散点插值，得到地质模型；
[0170] S182：根据上述地质模型和上述测井数据进行弹性参数横向插值，得到上述初始模型。
[0171] 本实施例中，可利用三维空间插值方法建立弹性参数模型。在上述步骤S181中，首先利用散点插值的方法对各个层位的数据进行插值，完成地质层位建模，在上述步骤S181中，根据地质层位进行弹性参数横向插值，即将测井信息进行横向插值，计算得到地下每个点上的弹性参数值，完成初始弹性参数建模的任务。本实施例中建立初始模型的方法较简便，有助于提高初始模型的创建速度，进而可提高页岩气储层脆性地震预测的效率。
[0172] 图11是本发明一实施例中计算得到弹性参数的扰动量求解表达式方法的流程示意图。如图11所示，在上述步骤S110中，利用一工区的实际地震数据和多个弹性参数的先验分布函数及Zoeppritz方程，基于Bayesian原理，借助广义线性反演思想求解得到上述弹性参数的扰动量的求解表达式方法，可包括步骤：
[0173] S111：利用上述先验分布函数、上述实际地震数据、上述初始模型、上述多个弹性参数的Zoeppritz方程所对应的正演算子，基于Bayesian原理，构建得到最大后验概率概率意义下的关于上述弹性参数的目标函数；
[0174] S112：利用广义线性反演思想对上述目标函数进行求解，得到上述求解表达式。
[0175] 在上述步骤S111中，可基于Bayesian原理、先验分布函数和多个弹性参数的精确Zoeppritz方程正演算子构建最大后验概率意义下的反演目标函数。在上述步骤S112中，可利用广义线性反演思想对反演目标函数进行求解，得到模型参数扰动量的求解表达式。
[0176] 一个实施例中，可利用上述初始模型对反演目标函数进行约束，并假设地震数据噪声服从高斯分布，先验分布函数服从柯西分布，那么，反演似然函数和先验概率分布可分别满足高斯分布和柯西分布。根据Bayesian原理，综合反演似然函数和先验分布函数可得到后验概率分布函数。根据后验概率确定反演的目标函数。
[0177] 一个实施例中，可考虑借助于泰勒展开对上述目标函数进行简化，变成关于模型参数扰动量的函数，然后将目标函数对扰动量进行求导并令该导数等于零，可得到扰动量的迭代求解公式，以此可以克服弹性参数求解困难的问题。其中，泰勒展开可求取关于弹性参数的Zoeppritz方程正演算子关于模型参数的一阶偏导数，该导数可以通过解析或数值方法求得。
[0178] 一个具体实施例中，设弹性参数的模型参数为mT＝(m1,m2,…,mn)T，观测地震数据即实际地震数据为dT＝(d1,d2,…,dn)T，由贝叶斯理论可以得到，在已知实际叠前地震数据的情况下，反演地下介质弹性参数的问题可以归结为求解一个后验概率函数：
[0179]
[0180] 其中，P(d)＝∫p(d|m)p(m)dm为归一化因子，可以看作是常数，P(d|m)为似然函数，P(m)为先验分布函数。假设似然函数和先验分布函数分别满足高斯分布和柯西分布，则似然函数P(d|m)和先验分布函数P(m)可分别为：
[0181]
[0182]
[0183] 其中，G为正演算子，CD是实际地震数据中噪声的协方差矩阵，Nd是观测数据即地震数据的长度。多波地震数据(包括PP波和PS波)的后验概率分布可以表示为：
[0184]
[0185] 一般可通过求解公式(16)的最大值作为最优解，该求解过程可等价于求解下式(公式(17))所示目标函数的最小值所对应的解，该目标函数可表示为：
[0186]
[0187] 一个实施例中，可假设噪声相对于观测数据而言是不相关的，那么噪声的协方差矩阵就可以简化为一个对角矩阵，在这种情况下，噪声的协方差矩阵可表示为，CD＝σn2I，其中，σn2是噪声的方差，I是Nd×Nd的单位阵，Nd是观测数据的长度。在上述目标函数中，可假设PP波数据的噪声方差为σpp，PS波数据的噪声方差为σps，则目标函数即公式(17)可以简化为：
[0188]
[0189] 其中，α＝σpp/σps为控制纵横波数据的使用比重，β＝σpp为控制先验信息的权重，下标pp表示PP波，下标ps表示PS波。
[0190] 借助广义线性反演思想对该目标函数进行求解，首先将关于杨氏模量、泊松比及0 1
密度三个弹性参数的Zoeppritz方程正演算子在初始模型参数即弹性参数的初始模型m (m＝m0+Δm)处进行泰勒展开，如下所示：
密度三个弹性参数的Zoeppritz方程正演算子在初始模型参数即弹性参数的初始模型m (m＝m0+Δm)处进行泰勒展开，如下所示：
[0191]
[0192] 可对公式(19)取一阶近似，即利用线性关系近似表达非线性关系，将公式(19)带入公式(18)，并对约束项进行改写可得到：
[0193]
[0194] 将上面得到的目标函数即公式(20)对弹性参数扰动量Δm进行求导并令导数等于零，可以得到扰动量为：
[0195]
[0196] 其中：
[0197]
[0198]
[0199] 考虑到等式即公式(21)两边都含有Δm，可采用迭代重加权最小平方算法求解扰动量，初始扰动量Δm0可取为零。由于该方法收敛速度快，并且该方法属于多层迭代，一般在此计算扰动量时可只迭代一次，以此可减少计算量。最终扰动量的求解表达式可以简化为：
[0200] Δmk＝(Hk)-1γk (22)
[0201] 其中：
[0202]
[0203]
[0204] 其中，正演算子关于模型参数的一阶导数 与 可通过解析的方法求得：
[0205]
[0206] 其中，W是地震子波矩阵，rpp(m)和rps(m)分别为纵波反射系数序列和转换横波反射系数序列。对于其中某一个角度θi有，
[0207]
[0208] 将公式(24)推广到n个输入角度的情况，正演算子G(m)关于模型参数m的一阶偏导数可变成一个(n×N)×3N矩阵，然后再推广到纵横波联合反演的情况，一阶偏导数最终就变成了一个(2n×N)×3N的矩阵。
[0209] 将公式(11)所示的Zoeppritz方程的新形式写成矩阵形式，可如下所示：
[0210] Ar＝b (25)
[0211] 其中，A表示公式(11)等号左侧的第一个矩阵，r表示公式(11)等号左侧的第二个矩阵，b表示公式(11)等号右侧的矩阵。
[0212] 对公式(25)两边同时对界面两侧的杨氏模量、泊松比和密度求偏导，可以得到：
[0213]
[0214] 由公式(26)进一步整理，便可得到某一界面反射系数和透射系数关于界面两侧的杨氏模量、泊松比和密度的偏导数：
[0215]
[0216] 最终，根据上述推导，可以得到弹性参数的模型参数m的更新迭代公式：
[0217] mk+1＝mk+ηkΔmk,k＝0,1,2... (28)
[0218] 其中，ηk是第k次迭代的步长因子，具体可根据需要设定。
[0219] 一个实施例中，在上述步骤S150中，可根据工区实际情况对Rickman于2008年提出的基于弹性参数的脆性指数(Brittleness Index，简称BI)计算公式进行合理的修正，然后上述反演得到的杨氏模量和泊松比数据体定量计算整个工区的岩石脆性指数，预测有利压裂区，为“甜点”位置的确定提供重要的依据。
[0220] 页岩脆性指数的计算公式可为：
[0221]
[0222] 其中，E为杨氏模量，单位为GPa；σ为泊松比，无量纲；BI为脆性指数，用％表示。
[0223] 一个实施例中，本发明实施例的方法具体地可采取以下过程来实现：基于地震数据提取角度依赖的子波，通过测井数据正演模拟和井旁道地震数据确定振幅缩放因子→基于工区内所有测井数据统计模型参数的先验信息→利用地震构造解释资料和测井数据，基于沉积模式建立初始弹性参数模型→基于岩石物理理论推导新形式的Zoeppritz方程，利用新形式的精确Zoeppritz方程正演模拟角度域的多波地震叠前道集并计算反演残差→利用广义线性反演思想将目标函数改写成关于弹性参数扰动量的函数，然后对扰动量求导并令导数等于零得到扰动量的迭代求解公式→利用迭代重加权最小二乘算法对扰动量进行求解，并利用求解得到的扰动量对模型参数进行更新→重复以上工作步骤直至反演残差达到要求或达到最大迭代次数→输出最终模型参数，包括杨氏模量、泊松比、密度模型→利用反演得到的岩石力学参数计算岩石的脆性指数，指导页岩气储层的开发。
[0224] 图12是本发明又一实施例的页岩气储层脆性地震预测方法的流程示意图。如图12所示，本发明实施例的页岩气储层脆性地震预测方法，具体地，可包括步骤：
[0225] S301：基于地震数据提取角度依赖的子波，基于测井数据和Zoeppritz正演模拟地震角道集并结合实际井旁地震数据确定振幅缩放因子；
[0226] S302：基于工区内所有测井数据统计模型参数的先验信息及其均值，包括杨氏模量、泊松比和密度，并计算包含三参数扰动量相关性的协方差矩阵；
[0227] S303：利用地震构造解释资料和测井数据，基于沉积模式建立时间域的初始弹性参数模型；
[0228] S304：根据岩石物理理论，对传统形式的Zoeppritz方程进行等价改写，得到关于杨氏模量、泊松比和密度的新形式的Zoeppritz方程，基于时间域的初始弹性参数模型和新形式的精确Zoeppritz方程正演模拟角度域的多波地震叠前道集，由正演模拟记录和实际记录直接计算PP波与PS波反演残差；
[0229] S305：基于Bayesian原理、先验信息和新形式的精确Zoeppritz方程正演算子构建最大后验概率意义下的反演目标函数，利用广义线性反演思想对反演目标函数进行求解，得到弹性参数扰动量的求解表达式；
[0230] S306：根据模型参数扰动量的求解表达式公式和反演残差，采用迭代重加权最小二乘算法计算模型参数的扰动量，重复迭代以上步骤S304～S306，并通过反演残差控制最大迭代次数，获得最优岩石力学参数，例如包括杨氏模量、泊松比和密度；
[0231] S307：根据工区实际情况对Rickman于2008年提出的基于弹性参数的脆性指数BI计算公式进行合理的修正，然后上述反演得到的杨氏模量和泊松比数据体定量计算整个工区的岩石脆性指数，预测有利压裂区，为“甜点”位置的确定提供重要的依据。
[0232] 图13是本发明一实施例中输入的不含噪地震叠前AVO PP波角道集地震数据，图14是本发明一实施例中输入的不含噪叠前AVO PS波角道集地震数据，图15是本发明一实施例中输入的含噪叠前AVO PP波角道集地震数据，图16是本发明一实施例中输入的含噪叠前AVO PS波角道集地震数据。如图13所示，地震数据为不含噪的PP波。如图14所示，地震数据为不含噪的PS波。如图15所示，地震数据为含噪的PP波。如图16所示，地震数据为含噪的PS波。在图13至图16中，纵轴表示时间，单位为秒，横轴表示地震入射角度。
[0233] 图17是本发明一实施例中无噪情况下多波AVO反演得到的弹性参数，图18是本发明一实施例中含噪情况下多波AVO反演得到的弹性参数。如图17和图18所示，采用本发明实施例的关于弹性参数的Zoeppritz方程的进行多波AVO反演得到多个弹性参数，包括杨氏模量、泊松比及密度。纵轴表示时间，单位为秒，横轴自左至右表示杨氏模量(单位为牛/平方米)、泊松比(无量纲)和密度(单位为千克/立方米)。由相同的时间坐标(纵坐标)可以对比各弹性参数的变化情况。如图18所示，在加入信噪比为2的随机噪声时的反演结果中，三变量柯西先验约束对保持反演过程稳定起到了关键作用。
[0234] 图19和图20分别是本发明一实施例中不含噪情况下和含噪情况下计算得到的岩石脆性指数。如图19和图20所示，根据Rickman于2008年提出的基于弹性参数的脆性指数计算公式利用模型数据计算得到的岩石脆性指数，图19为不含噪声时的计算结果，图20为加入信噪比为2的随机噪声时的计算结果，纵轴表示时间，单位为秒，横轴表示脆性指数，单位为％。由图13至图20对比可知，本发明的方法可以取得较佳的效果。
[0235] 本实施例中，基于精确Zoeppritz方程的多波AVO反演方法不仅对杨氏模量、泊松比有精确的预测结果，由于引入了包含密度信息的先验分布，并采用PS横波进行联合反演，其对密度模型也预测准确。
[0236] 本发明的页岩气储层脆性地震预测方法，基于关于多个弹性参数的Zoeppritz方程可直接通过叠前地震数据反演得到能够直接表征岩石脆性的岩石力学参数，能够有效避免由现有弹性参数间接计算岩石力学参数带来的累积误差。相对于现有的Zoeppritz近似方程而言，本发明地震叠前AVO反演基于关于多个弹性参数的Zoeppritz方程能够克服因入射角(或偏移距)过大而带来计算误差的问题，从而使得我们能够利用大偏移距信息进行储层参数预测。本发明关于多个弹性参数的精确Zoeppritz方程的AVO反演方法，具有更好的适应性，能够满足泥砂岩薄互层结构，纵向弹性参数差异大的储层条件，如页岩气储层，有效提高叠前反演精度。本发明可以获得弹性参数的较佳模型，以此计算岩石脆性指数，能够提供高可靠性的页岩气“甜点”地震识别方法，相比于通过岩样测试和测井分析估算工区储层的岩石脆性分布特征，预测“甜点”位置的方法，能提供横向连续性更好，且更加真实的反映地下页岩气储层脆性分布特征脆性指数数据体。进一步，本发明联合PP波与PS波进行多个弹性参数的反演，PS波地震数据包含更多岩石物理属性和密度信息，能有效提高反演的稳定性，克服现有技术中仅考虑PP波而只能丢弃包含大量岩石力学信息的PS波的缺点。再进一步，本发明通过多变量柯西先验分布引入模型参数扰动量的协方差矩阵，能够很好地消除待反演的多个弹性参数扰动量之间的相关性，降低了反演的不确定性，提高了反演的精度。
[0237] 基于与图1所示的页岩气储层脆性地震预测方法相同的发明构思，本申请实施例还提供了一种页岩气储层脆性地震预测装置，如下面实施例所述。由于该页岩气储层脆性地震预测装置解决问题的原理与页岩气储层脆性地震预测方法相似，因此该页岩气储层脆性地震预测装置的实施可以参见页岩气储层脆性地震预测方法的实施，重复之处不再赘述。
[0238] 图21是本发明一实施例的页岩气储层脆性地震预测装置的结构示意图。如图21所示，本发明实施例的页岩气储层脆性地震预测装置，可包括：扰动量表达式生成单元210、扰动量值生成单元220、当前模型生成单元230、反演残差生成单元240及脆性指数生成单元250，上述各单元顺序连接。
[0239] 扰动量表达式生成单元210用于利用一工区的实际地震数据和多个弹性参数的先验分布函数及Zoeppritz方程，基于Bayesian原理，借助广义线性反演思想求解得到上述弹性参数的扰动量的求解表达式。
[0240] 扰动量值生成单元220用于根据上述求解表达式和上述弹性参数的初始模型计算得到上述弹性参数的扰动量值。
[0241] 当前模型生成单元230用于根据上述扰动量值和上述初始模型计算得到上述弹性参数的当前模型。
[0242] 反演残差生成单元240用于基于上述实际地震数据、上述当前模型及上述多个弹性参数的Zoeppritz方程计算得到地震反演残差。
[0243] 脆性指数生成单元250用于若上述地震反演残差小于一设定阈值，将上述当前模型作为上述弹性参数的值，根据上述弹性参数的值计算得到上述工区的页岩气储层脆性指数。
[0244] 在扰动量表达式生成单元210中，上述弹性参数可以指介质的岩石力学参数，可以包括杨氏模量、泊松比及密度等参数。
[0245] 该实际地震数据可以是通过地震测量的实际记录，可以是实际角度域地震道集。该实际地震数据可以包括多种地震波，例如PP波和PS波，PS波包含大量的岩石力学信息，基于包含PP波和PS波的实际地震数据计算弹性参数的扰动量，有利于后续得到确定性更高的弹性参数。
[0246] 多个弹性参数的先验分布函数及Zoeppritz方程具体是指多个弹性参数的先验分布函数和关于多个弹性参数的Zoeppritz方程，关于多个弹性参数的Zoeppritz方程可以是以弹性力学参数为变量的新形式的Zoeppritz方程，可以通过对现有的Zoeppritz方程进行改写或变换得到。
[0247] 现有的Zoeppritz方程不是关于上述多个弹性参数的方程，可以是关于横波波速、纵波波速及密度的方程。现有的Zoeppritz方程可以具有多种形式，例如Zoeppritz方程近似表达式、精确(exact)的Zoeppritz方程，可如上述公式(3)所示。
[0248] 在扰动量值生成单元220中，上述弹性参数的初始模型可以是初始弹性参数模型，可通过多种不同方法得到，例如，利用地震构造解释资料和测井数据经过建模得到。
[0249] 在当前模型生成单元230中，该当前模型可以是上述弹性参数的最新模型，可通过在初始模型上增加一定的扰动量得到。
[0250] 在反演残差生成单元240中，可以利用上述当前模型和上述多个弹性参数的Zoeppritz方程得到模拟的地震数据，再利用模拟的地震数据和实际地震数据计算得到地震反演残差。其中，地震反演残差是基于弹性参数的Zoeppritz方程计算得到，可以克服现有的Zoeppritz方程因入射角或偏移距过大带来的计算误差的问题，从而能够利用大偏移距信息进行储层参数预测，可以提高岩石脆性指数预测结果的准确度。
[0251] 在脆性指数生成单元250中，该设定阈值可以根据需要进行设置。可以采用多种方法根据上述弹性参数的值计算上述工区的页岩气储层脆性指数，例如，可以是Rickman于2008年提出的基于弹性参数的脆性指数(BrittlenessIndex，BI)计算公式，或者是根据工区实际情况对Rickman于2008年提出的基于弹性参数的脆性指数BI计算公式进行合理的修正，然后利用修正后的计算公式根据上述弹性参数的值计算整个工区的岩石脆性指数，以此可以预测有利压裂区，可为“甜点”位置的确定提供重要的依据。
[0252] 一个实施例中，若上述地震反演残差不小于该设定阈值，可基于当前模型重复利用扰动量值生成单元220、当前模型生成单元230及反演残差生成单元240重新计算，得到新的地震反演残差，并再次判断新的地震反演残差是否小于该设定阈值，如此迭代进行，直到地震反演残差小于该设定阈值。通过地震反演残差控制迭代计算的次数，有利于获得弹性参数的最优模型，从而能够获得最优的岩石脆性指数预测结果。
[0253] 本发明实施例中，基于关于多个弹性参数的Zoeppritz方程计算得到上述弹性参数的扰动量的求解表达式，并基于上述求解表达式计算得到上述弹性参数的扰动量值，再计算得到上述弹性参数的当前模型，并基于上述当前模型及上述多个弹性参数的Zoeppritz方程计算得到地震反演残差，无需利用由叠前反演得到的纵波速度、横波速度以及密度间接获得或者基于现有Zoeppritz方程的近似公式反演获得弹性参数的模型，避免产生累计误差，提高了弹性参数的计算精度。相对于现有的Zoeppritz近似方程而言，关于多个弹性参数的Zoeppritz方程，能够克服了因入射角(或偏移距)过大而带来计算误差的问题，从而使得能够利用大偏移距信息进行储层参数预测。通过地震反演残差控制是否将计算得到的当前模型作为弹性参数的最终模型，有利于得到弹性参数的最优模型，利用该最优模型预测岩石脆性指数，可以提高预测精度。
[0254] 一个实施例中，上述扰动量表达式生成单元，还用于执行：上述实际地震数据包括PP波和PS波；反演残差生成单元，还用于执行：上述地震反演残差为包含PP波和PS波的多波的地震联合反演残差。上述实际地震数据可包括PP波和PS波，利用上述当前模型和上述多个弹性参数的Zoeppritz方程所得到的模拟的地震数据可包括PP波和PS波，则该地震反演残差可以为包含PP波和PS波的多波的地震联合反演残差。
[0255] 本发明实施例中，上述实际地震数据可包括PP波和PS波，该地震反演残差可以为包含PP波和PS波，联合PP波与PS波进行多个弹性参数反演，与现有通常仅考虑PP波，本实施例中不仅考虑了PP波，还考虑了PS波，由于PS波地震数据包含更多岩石物理属性和密度信息，能有效提高反演的稳定性。
[0256] 图22是本发明另一实施例的页岩气储层脆性地震预测装置的结构示意图。如图22所示，图21所示的页岩气储层脆性地震预测装置，还可包括：先验分布生成单元260，与上述扰动量表达式生成单元210连接。
[0257] 先验分布生成单元260用于利用上述工区的测井数据统计分析得到多个弹性参数的先验分布函数。
[0258] 在先验分布生成单元260中，较佳地，可以基于工区内所有测井数据统计弹性参数的先验分布函数，以此利用充足的测井数据可以得到更准确的先验分布函数。该先验分布函数可以包括各弹性参数的均值及所有弹性参数的扰动量相关性的协方差矩阵。上述多个弹性参数可包括杨氏模量、泊松比和密度。
[0259] 本实施例中，除了根据上述扰动量值和上述初始模型计算得到上述弹性参数的当前模型，还利用工区的测井数据获得地下介质地震弹性参数模型的信息作为先验信息，可以降低地震反演的不确定性，提高反演过程的稳定性，提高地震反演残差的精度，进而提高页岩气储层脆性指数预测的准确度。
[0260] 一个实施例中，上述先验分布生成单元260，可包括：先验分布生成模块261。
[0261] 先验分布生成模块261用于基于上述弹性参数服从柯西分布的假设和上述工区内所有生产井的测井数据，计算得到上述弹性参数的均值和上述多个弹性参数的扰动量的协方差矩阵，作为上述先验分布函数。
[0262] 本实施例中，可通过测井数据分析获得需要的弹性参数，在利用多变量柯西分布函数求取各弹性参数的自相关系数和互相关系数，并构建多个弹性参数扰动量的协方差矩阵，形成符合该工区的弹性参数的模型参数的先验分布函数。本实施例中的先验分布函数包括多个弹性参数的扰动量的协方差矩阵，且基于上述弹性参数服从柯西分布的假设获得，利用该先验分布函数计算弹性参数的扰动量的求解表达式，可以很好地消除待反演弹性参数扰动量之间的相关性，从而降低反演的不确定性，提高反演的精度，进而提高脆性指数预测准确度。
[0263] 一个具体实施例中，上述多个弹性参数可包括杨氏模量、泊松比及密度三个参数，以此为例，上述多变量柯西分布函数可以是三变量柯西分布函数。可采用该三变量柯西分布函数作为上述三个弹性参数的先验分布关系式，并可基于工区内所有测井数据统计得到的各弹性参数的模型信息，求取各弹性参数的均值、各弹性参数的自相关系数及三个弹性参数的互相关系数，进而构建三弹性参数的扰动量的协方差矩阵，利用求取的均值和协方差矩阵形成该工区的弹性参数的先验分布函数。
[0264] 图23是本发明又一实施例的页岩气储层脆性地震预测装置的结构示意图。如图23所示，图21所示的页岩气储层脆性地震预测装置，还可包括：Zoeppritz方程变换单元270，与上述扰动量表达式生成单元210连接。
[0265] Zoeppritz方程变换单元270用于对初始Zoeppritz方程进行等价变换，得到上述多个弹性参数的Zoeppritz方程。
[0266] 在上述Zoeppritz方程变换单元270中，该初始Zoeppritz方程可以是包括近似方程在内的多种现有形式，较佳地是精确(exact)的Zoeppritz方程。得到的上述多个弹性参数的Zoeppritz方程是新形式的Zoeppritz方程。
[0267] 一个实施例中，上述Zoeppritz方程变换单元270，可包括：Zoeppritz方程变换模块271。
[0268] Zoeppritz方程变换模块271用于基于岩石物理理论和地层介质的各向同性假设，对初始Zoeppritz方程进行等价变换，得到上述多个弹性参数的Zoeppritz方程。
[0269] 本实施例中，根据岩石物理理论，对现有形式的Zoeppritz方程进行等价改写，可得到关于多个弹性参数的新形式的Zoeppritz方程。相对于Zoeppritz近似方程而言，能够克服了因入射角(或偏移距)过大而带来计算误差的问题，基于该种关于多个弹性参数的Zoeppritz方程进行的AVO反演，具有更好的岩石适应性，能够满足泥砂岩薄互层结构，纵向弹性参数差异大的储层条件，如页岩气储层，能够有效提高叠前反演精度，进而提高脆性指数预测准确度。
[0270] 具体地，例如，首先基于各向同性假设，利用杨氏模量、泊松比、纵横波速度以及密度之间的关系，对传统形式的Zoeppritz方程进行等价转换，最终得到关于杨氏模量、泊松比和密度这三个弹性参数的新形式的Zoeppritz方程。
[0271] 图24是本发明再一实施例的页岩气储层脆性地震预测装置的结构示意图。如图24所示，图21所示的页岩气储层脆性地震预测装置，还可包括：初始模型生成单元280，可连接于扰动量值生成单元220和当前模型生成单元230之间。
[0272] 初始模型生成单元280用于利用地震构造解释数据和上述工区的测井数据建立上述弹性参数的初始模型。
[0273] 本实施例中，利用地震构造解释数据和上述工区的测井数据可以通过多种方法建立上述弹性参数的初始模型。例如，可利用三维空间插值方法建立弹性参数的初始模型。具体地，可利用地震构造解释资料，基于沉积模式建立地质模型，并将测井资料，按构造模式进行插值和外推，得到每条测线的初始弹性参数模型，初始弹性参数模型可包括杨氏模量模型、泊松比模型和密度模型，以此可以快速建立弹性参数的初始模。
[0274] 图25是本发明一实施例中反演残差生成单元的结构示意图。如图25所示，上述反演残差生成单元240可包括：地震子波提取模块241和反演残差生成模块242，二者相互连接。
[0275] 地震子波提取模块241用于基于上述实际地震数据提取得到地震子波。
[0276] 反演残差生成模块242用于利用上述地震子波、上述初始模型、上述多个弹性参数的Zoeppritz方程及上述实际地震数据，计算得到地震反演残差。
[0277] 在上述地震子波提取模块241中，该实际地震数据可以是实际的地震叠前道集，可以包括角度域的PP波道集和PS波道集。一个实施例中，可基于实际的地震叠前道集和测井数据采取统计方法提取地震子波，该测井数据可以对基于上述实际地震数据提取得到地震子波进行修正，以提高提取地震子波的稳定性。
[0278] 在上述反演残差生成模块242中，上述初始模型可以是时间域的初始弹性参数模型。可利用该初始模型、基于多个弹性参数的Zoeppritz方程及上述实际地震数据计算得到工区的模拟的地震角道集，例如角度域模拟的PP波道集和PS波道集，再由模拟的地震角道集和角度域的实际地震道集，计算得到例如PP波与PS波的地震反演残差。
[0279] 一个实施例中，上述地震子波提取模块241可包括：地震子波筛选模块2411。
[0280] 地震子波筛选模块2411用于通过统计方法从上述实际地震数据的地震叠前道集提取波形依赖于地震波入射角度的多个上述地震子波。
[0281] 本实施例中，地震子波在传播过程中受地层的影响会发生波形或频率变化，提取依赖于入射角度的地震子波，能有效提高振幅匹配程度。
[0282] 图26是本发明再一实施例的页岩气储层脆性地震预测装置的结构示意图。如图26所示，图21所示的页岩气储层脆性地震预测装置，还可包括：弹性参数值更新单元290，与上述反演残差生成单元240连接。
[0283] 弹性参数值更新单元290用于若上述地震反演残差大于或等于上述设定阈值，则根据上述求解表达式和上述当前模型重新计算上述弹性参数的扰动量值，并根据重新计算的扰动量值和上述当前模型重新计算上述弹性参数的当前模型，以及基于上述实际地震数据、重新计算的当前模型及上述多个弹性参数的Zoeppritz方程重新计算地震反演残差，如此重复进行直到地震反演残差小于上述设定阈值，并将此时重新计算的当前模型作为上述弹性参数的值。
[0284] 在上述弹性参数值更新单元290中，根据弹性参数扰动量的求解表达式公式和反演残差，可采用迭代重加权最小二乘算法计算弹性参数的扰动量，重复迭代计算反演残差，并通过反演残差控制最大迭代次数，可以获得最优模型参数，进而提高岩石脆性指数的预测精度。
[0285] 图27是本发明一实施例中反演残差生成模块的结构示意图。如图27所示，上述反演残差生成模块242可包括：地震反射系数生成模块2421、模拟地震角道集生成模块2422及地震反演残差确定模块2423，上述各模块顺序连接。
[0286] 地震反射系数生成模块2421用于以上述初始模型作为输入，利用上述多个弹性参数的Zoeppritz方程计算各地震入射角度的地震反射系数。
[0287] 模拟地震角道集生成模块2422用于将上述地震子波和上述地震反射系数进行褶积，得到上述工区的模拟的地震角道集。
[0288] 地震反演残差确定模块2423用于将上述工区的模拟的地震角道集与上述实际地震数据中的地震道角集作差，得到上述地震反演残差。
[0289] 在地震反射系数生成模块2421中，可以时间域的初始模型作为输入，利用上述多个弹性参数的Zoeppritz方程计算每一个角度地震波例如PP波与PS波的反射系数向量。
[0290] 在模拟地震角道集生成模块2422中，将上述地震子波和上述地震反射系数进行褶积，可得到上述工区的模拟的地震角道集，例如角度域的PP波角道集与PS波角道集。
[0291] 在地震反演残差确定模块2423中，将模拟的地震角道集(例如模拟的PP波角道集与PS波角道集)与实际的地震角道集(例如实际的PP波角道集与PS波角道集)作差，可得到上述反演残差(例如PP波与PS波的多波联合反演残差)。本实施例的反演残差弹性参数的初始模型计算得到，可以有效控制上述迭代次数，获取弹性参数的最优模型。
[0292] 图28是本发明另一实施例中反演残差生成单元的结构示意图。如图28所示，上述反演残差生成单元240还可包括：地震子波修正模块243，连接于地震子波提取模块241和反演残差生成模块242之间。
[0293] 地震子波修正模块243用于利用上述工区内生产井的测井数据、初始Zoeppritz方程及上述实际地震数据三者中的一个或多个对上述地震子波进行修正，修正后的地震子波用于计算上述地震反演残差。
[0294] 本实施例中，可以通过多种方式单独或一同对上述地震子波进行修正。例如，可以仅利用上述工区内生产井的测井数据对上述地震子波进行修正，或者，可以上述工区内生产井的测井数据、初始Zoeppritz方程及上述实际地震数据三者求得地震振幅缩放因子，在利用地震振幅缩放因子对上述地震子波进行修正。修正后的地震子波有利于使上述模拟的地震数据和上述实际地震数据具有更好的匹配度。
[0295] 图29是本发明一实施例中地震子波修正模块的结构示意图。如图29所示，上述地震子波修正模块243可包括：模拟地震角道集正演模块2431、地震振幅缩放因子生成模块2432及地震振幅缩放模块2433，上述各模块顺序连接。
[0296] 模拟地震角道集正演模块2431用于以上述测井数据作为输入模型，利用初始Zoeppritz方程正演得到模拟的地震角道集。
[0297] 地震振幅缩放因子生成模块2432用于根据上述模拟的地震角道集和上述实际地震数据中的上述生产井周围的地震角道集，计算得到地震振幅缩放因子。
[0298] 地震振幅缩放模块2433用于利用上述地震振幅缩放因子对上述地震子波进行振幅缩放处理。
[0299] 在上述模拟地震角道集正演模块2431中，该初始Zoeppritz方程可以是现有精确的Zoeppritz方程。正演得到模拟的地震角道集可以包括角度域的PP波角道集和PS波角道集。
[0300] 在上述步骤地震振幅缩放因子生成模块2432和地震振幅缩放模块2433中，一个实施例中，实际地震数据的信噪比较高，可为角道集的每一道使用统一的振幅缩放因子，以保证振幅随偏移距的变化关系。另一实施例中，实际地震数据的信噪比低，可分别利用近偏移距、中偏移距、远偏移距计算地震振幅缩放因子，以保证模拟记录(模拟的地震数据)与实际记录(实际的地震数据)的最佳振幅匹配，减少噪声对反演过程的影响。
[0301] 本实施例中，实际的地震振幅往往是相对值，采用初始Zoeppritz方程正演模拟的地震数据振幅与实际地震振幅存在一定数值差异。利用地震振幅缩放因子对提取的地震子波进行振幅缩放处理，能够使得模拟的地震数据的振幅与实际的地震数据具有较好的振幅匹配度。
[0302] 图30是本发明一实施例中初始模型生成单元的结构示意图。如图30所示，上述初始模型生成单元280，可包括：地质模型生成模块281和初始模型生成模块282，二者相互连接。
[0303] 地质模型生成模块281用于依据上述地震构造解释数据进行层位散点插值，得到地质模型。
[0304] 初始模型生成模块282用于根据上述地质模型和上述测井数据进行弹性参数横向插值，得到上述初始模型。
[0305] 本实施例中，可利用三维空间插值方法建立弹性参数模型。上述地质模型生成模块281，首先利用散点插值的方法对各个层位的数据进行插值，完成地质层位建模，初始模型生成模块282，根据地质层位进行弹性参数横向插值，即将测井信息进行横向插值，计算得到地下每个点上的弹性参数值，完成初始弹性参数建模的任务。本实施例中建立初始模型的方法较简便，有助于提高初始模型的创建速度，进而可提高页岩气储层脆性地震预测的效率。
[0306] 图31是本发明一实施例中扰动量表达式生成单元的结构示意图。如图31所示，上述扰动量表达式生成单元210，可包括：目标函数生成模块211和扰动量表达式生成模块212，二者相互连接。
[0307] 目标函数生成模块211用于利用上述先验分布函数、上述实际地震数据、上述初始模型、上述多个弹性参数的Zoeppritz方程所对应的正演算子，基于Bayesian原理，构建得到最大后验概率意义下的关于上述弹性参数的目标函数。
[0308] 扰动量表达式生成模块212用于利用广义线性反演思想对上述目标函数进行求解，得到上述求解表达式。
[0309] 在上述目标函数生成模块211中，可基于Bayesian原理、先验分布函数和多个弹性参数的精确Zoeppritz方程正演算子构建最大后验概率意义下的反演目标函数。在上述扰动量表达式生成模块212中，可利用广义线性反演思想对反演目标函数进行求解，得到模型参数扰动量的求解表达式。
[0310] 一个实施例中，可利用上述初始模型对反演目标函数进行约束，并假设地震数据噪声服从高斯分布，先验分布函数服从柯西分布，那么，反演似然函数和先验概率分布可分别满足高斯分布和柯西分布。根据Bayesian原理，综合反演似然函数和先验分布函数可得到后验概率分布函数。根据后验概率确定反演的目标函数。
[0311] 一个实施例中，可考虑借助于泰勒展开对上述目标函数进行简化，变成关于模型参数扰动量的函数，然后将目标函数对扰动量进行求导并令该导数等于零，可得到扰动量的迭代求解公式，以此可以克服弹性参数求解困难的问题。其中，泰勒展开可求取关于弹性参数的Zoeppritz方程正演算子关于模型参数的一阶偏导数，该导数可以通过解析或数值方法求得。
[0312] 图32是本发明另一实施例的页岩气储层脆性地震预测装置的结构示意图。如图32所示，该装置可包括：地震子波获取及修正单元401、弹性参数均值及协方差矩阵获取单元402、初始弹性参数模型获取单元403、PP波与PS波反演残差获取单元404、弹性参数扰动量求解表达式获取单元405、最优弹性参数获取单元406及岩石脆性指数预测单元407。
[0313] 地震子波获取及修正单元401用于：基于地震数据提取角度依赖的子波，基于测井数据和Zoeppritz正演模拟地震角道集并结合实际井旁地震数据确定振幅缩放因子。
[0314] 弹性参数均值及协方差矩阵获取单元402用于：基于工区内所有测井数据统计模型参数的先验信息及其均值，包括杨氏模量、泊松比和密度，并计算包含三参数扰动量相关性的协方差矩阵。
[0315] 初始弹性参数模型获取单元403用于：利用地震构造解释资料和测井数据，基于沉积模式建立时间域的初始弹性参数模型。
[0316] PP波与PS波反演残差获取单元404，用于：根据岩石物理理论，对传统形式的Zoeppritz方程进行等价改写，得到关于杨氏模量、泊松比和密度的新形式的Zoeppritz方程，基于时间域的初始弹性参数模型和新形式的精确Zoeppritz方程正演模拟角度域的多波地震叠前道集，由正演模拟记录和实际记录直接计算PP波与PS波反演残差。
[0317] 弹性参数扰动量求解表达式获取单元405用于：基于Bayesian原理、先验信息和新形式的精确Zoeppritz方程正演算子构建最大后验概率意义下的反演目标函数，利用广义线性反演思想对反演目标函数进行求解，得到弹性参数扰动量的求解表达式。
[0318] 最优弹性参数获取单元406用于：根据模型参数扰动量的求解表达式公式和反演残差，采用迭代重加权最小二乘算法计算模型参数的扰动量，重复迭代以上步骤S304～S306，并通过反演残差控制最大迭代次数，获得最优岩石力学参数，例如包括杨氏模量、泊松比和密度。
[0319] 岩石脆性指数预测单元407用于：根据工区实际情况对Rickman于2008年提出的基于弹性参数的脆性指数BI计算公式进行合理的修正，然后上述反演得到的杨氏模量和泊松比数据体定量计算整个工区的岩石脆性指数，预测有利压裂区，为“甜点”位置的确定提供重要的依据。
[0320] 本发明的页岩气储层脆性地震预测装置，基于关于多个弹性参数的Zoeppritz方程可直接通过叠前地震数据反演得到能够直接表征岩石脆性的岩石力学参数，能够有效避免由现有弹性参数间接计算岩石力学参数带来的累积误差。相对于现有的Zoeppritz近似方程而言，本发明地震叠前AVO反演基于关于多个弹性参数的Zoeppritz方程能够克服因入射角(或偏移距)过大而带来计算误差的问题，从而使得我们能够利用大偏移距信息进行储层参数预测。本发明关于多个弹性参数的精确Zoeppritz方程的AVO反演方法，具有更好的适应性，能够满足泥砂岩薄互层结构，纵向弹性参数差异大的储层条件，如页岩气储层，有效提高叠前反演精度。本发明可以获得弹性参数的较佳模型，以此计算岩石脆性指数，能够提供高可靠性的页岩气“甜点”地震识别方法，相比于通过岩样测试和测井分析估算工区储层的岩石脆性分布特征，预测“甜点”位置的方法，能提供横向连续性更好，且更加真实的反映地下页岩气储层脆性分布特征脆性指数数据体。进一步，本发明联合PP波与PS波进行多个弹性参数的反演，PS波地震数据包含更多岩石物理属性和密度信息，能有效提高反演的稳定性，克服现有技术中仅考虑PP波而只能丢弃包含大量岩石力学信息的PS波的缺点。再进一步，本发明通过多变量柯西先验分布引入模型参数扰动量的协方差矩阵，能够很好地消除待反演的多个弹性参数扰动量之间的相关性，降低了反演的不确定性，提高了反演的精度。
[0321] 在本说明书的描述中，参考术语“一个实施例”、“一个具体实施例”、“一些实施例”、“例如”、“示例”、“具体示例”、或“一些示例”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施例或示例中。在本说明书中，对上述术语的示意性表述不一定指的是相同的实施例或示例。而且，描述的具体特征、结构、材料或者特点可以在任何的一个或多个实施例或示例中以合适的方式结合。
[0322] 本领域内的技术人员应明白，本发明的实施例可提供为方法、系统、或计算机程序产品。因此，本发明可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且，本发明可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。
[0323] 本发明是参照根据本发明实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。
[0324] 这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。
[0325] 这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。
[0326] 以上所述的具体实施例，对本发明的目的、技术方案和有益效果进行了进一步详细说明，所应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限定本发明的保护范围，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。
| 标题 | 发布/更新时间 | 阅读量 |
|---|---|---|
| 地震波传播速度获取方法、测深方法及地震波观测系统 | 2020-05-12 | 327 |
| 一种对地震波屏蔽与波导的抗震建筑 | 2020-05-13 | 1000 |
| 基于地震波的地下入侵监测系统 | 2020-05-13 | 467 |
| 地震波正演方法及系统 | 2020-05-11 | 820 |
| 一种利用斜井直达地震波走时计算地层地震波速度的方法 | 2020-05-12 | 481 |
| 一种地震波分解方法 | 2020-05-11 | 670 |
| 地震波立体演示台 | 2020-05-11 | 111 |
| 勘测三维地震波的方法及设备 | 2020-05-12 | 120 |
| 一种地震波形匹配方法及装置 | 2020-05-13 | 226 |
| 用于地震波传感器的耦合装置 | 2020-05-13 | 230 |
高效检索全球专利专利汇是专利免费检索,专利查询,专利分析-国家发明专利查询检索分析平台，是提供专利分析，专利查询，专利检索等数据服务功能的知识产权数据服务商。
我们的产品包含105个国家的1.26亿组数据，免费查、免费专利分析。
分析报告
专利汇分析报告产品可以对行业情报数据进行梳理分析，涉及维度包括行业专利基本状况分析、地域分析、技术分析、发明人分析、申请人分析、专利权人分析、失效分析、核心专利分析、法律分析、研发重点分析、企业专利处境分析、技术处境分析、专利寿命分析、企业定位分析、引证分析等超过60个分析角度，系统通过AI智能系统对图表进行解读，只需1分钟，一键生成行业专利分析报告。
QQ群二维码
意见反馈
意见反馈