一种基于地震反射振幅方位各向异性差异的裂缝预测方法

1.本发明属于油气地球物理勘探技术领域，具体涉及一种基于地震反射振幅方位各向异性差异进行致密砂岩气储层裂缝预测方法。


背景技术：

2.随着石油天然气资源的勘探和开发及各向异性技术的发展与完善，裂缝型油气储层成为当前油气勘探的重要方向。裂缝不仅能作为油气的储层空间，还可为油气的运移提供重要的通道。因此，对裂缝型储层的勘探研究具有十分重要的意义。
3.thomsen提出了弱各向异性理论及弱各向异性参数并给出其表达式。ruger推导出hti介质的纵波反射系数近似表达式。vavrycuk和psencik推导出由各向异性参数表示的两个弱各向异性介质弱对比界面的纵波反射和透射系数公式。mallick研究了振幅随方位角的变化(avaz)，并将该方法应用于实际的方位数据中预测裂缝方位。朱培民详细的讨论了hti介质中弹性参数、各向异性参数和裂缝参数之间的关系，并且提出了利用地震数据反演裂缝密度的方法。martins将各向异性效应加入到弹性阻抗公式中，推导出弱各向异性介质弹性阻抗方程。bachrach推导了hti介质弹性参数、各向异性参数和反射系数的关系，并使用贝叶斯反演方法估计裂缝密度和孔隙度。尹志恒介绍了各向异性技术的发展，总结了应用于裂缝检测的纵波各向异性新的方法与技术。张广智和陈怀震研究了hti介质振幅随偏移距和方位角的变化，并对反射系数公式进行简化，利用叠前数据反演弹性参数和各向异性参数。王康宁(2015)针对任意对称面各向异性介质，对pp波反射系数近似公式进行整理，并结合线性滑动理论对其进行傅里叶级数展开，利用纵波方位地震数据进行裂缝预测。downton提出描述avaz响应的系数可以用方位傅里叶系数表示，并通过线性近似解释p波反射情况下这些系数对于各向异性和裂缝参数的意义。cui基于非接触边界条件推导了新的pp波反射和透射系数公式，得到裂缝介质avo反演的近似公式。barone明确针对垂直裂缝各向同性介质，提出不存在90
°
模糊准确估计裂缝走向方位的方法，并反演与裂缝密度有线性相关的新参数。赵才顺通过岩石物理正演模拟及成像测井资料建立裂缝方位和纵波反射系数的关系，基于叠前纵波方位地震资料利用裂缝相对密度半定量化研究技术预测致密砂岩储层裂缝密度。宋维琪将hti介质pp波反射系数近似公式展开为傅里叶级数形式，并针对不同方位角范围进一步展开，利用最小二乘拟合提取傅里叶级数系数，使各向异性参数和裂缝预测的估计更加准确。
4.目前，各向异性反演仍存在固有问题，即在常规同时反演各向同性速度与各向异性参数的反演机制中，各向同性项加权系数高于各向异性项加权系数几个数量级，使得各向异性参数反演的稳定性与精确性受到极大影响。基于此，急需研发一种基于地震反射振幅方位各向异性差异的致密砂岩气储层无井约束裂缝预测方法，以解决上述问题。


技术实现要素：

5.本发明的目的就在于提供一种基于地震反射振幅方位各向异性差异的致密砂岩
气储层无井约束裂缝预测方法，以解决提高各向异性参数反演的稳定性和精度，直观表征致密砂岩气储层裂缝发育程度的问题。在地震avaz反射系数方法的基础上通过与参考方位数据做差消除各向同性部分的影响，提出方位反射系数差的反演方法对各向异性参数进行反演，该方法充分利用各向异性部分的响应并且提高了反演的稳定性和精度；并且，通过岩石物理理论，进一步计算得到更具有物理意义的裂缝密度参数，直观描述裂缝的发育程度；最后，在裂缝预测与弹性反演压缩率预测的基础上，构建用于致密砂岩气储层裂缝及流体综合评价新的指示因子。
6.本发明的目的是通过以下技术方案实现的：
7.一种基于地震反射振幅方位各向异性差异的裂缝预测方法，包括以下步骤：
8.步骤a、输入方位角道集地震数据；
9.步骤b、对hti介质pp波方位反射系数公式进行线性化，由纵波速度、横波速度、密度和各向异性参数等线性表达，选取裂缝面法向为参考方位，构建以裂缝法向为参考单位的方位反射系数差公式；
[0010][0011]
其中，θ表示入射角，表示观测方位，表示垂直于裂缝走向的方位，v
p
和vs分别ε
(v)
,γ和δ
(v)
表示thomsen各向异性参数，符号
ˉ
和
△
分别表示上层和下层介质性质的平均值和差值；
[0012]
步骤c、基于方位反射系数差公式，通过理论模型计算并分析线性化的pp波方位反射系数，以及方位反射系数对储层裂缝的敏感性；
[0013]
步骤d、基于方位反射系数差公式构建方位地震反演方法，计算各向异性参数：
[0014][0015]
其中，ai＝sin2θi,bi＝tan2θi,i＝1,
…
,m和j＝1,
…
,n，m和n分别为地震数据中入射角、方位角道集个数；
[0016]
步骤e、应用傅里叶级数展开方法对步骤a得到的方位角道集数据预测裂缝方位，并以裂缝面法向为参考方位构建地震方位差道集，作为步骤d反演方法的输入数据；
[0017]
步骤f、基于步骤d构建的方位地震反演方法、输入步骤e中构建的地震方位差道集为输入数据进行反演，计算含裂缝的致密砂岩气储层hti各向异性参数；
[0018]
步骤g、基于schoenberg线性滑动和hudson币状裂缝岩石物理模型，将步骤f反演
得到的各向异性参数转换为裂缝切向弱度参数δ
t
，将切向裂缝弱度δ
t
进一步转换为裂缝密度参数e，直观描述致密砂岩气储层裂缝发育程度；
[0019]
δ
t
＝2γ
ꢀꢀꢀ
(23)
[0020][0021]
其中，g为s波和p波速度比的平方；
[0022]
步骤h、将步骤g得到的裂缝密度e与地震弹性反演得到的压缩率1/k相乘并归一化，构建用于致密砂岩气储层裂缝及流体综合评价的指示因子e
×
1/k，其中，裂缝密度e反映致密砂岩气储层裂缝发育程度，压缩率1/k反映致密砂岩气储层含气性；
[0023]
步骤i、由测井数据及储层开发数据对步骤h预测结果进行标定。
[0024]
进一步地，步骤a，对宽方位地震数据进行方位角道集处理作为基于地震反射振幅方位各向异性差异的裂缝预测方法的输入数据。
[0025]
与现有技术相比，本发明的有益效果是：
[0026]
本发明提出基于地震反射振幅方位各向异性差异的方法，在裂缝方位预测的基础上，以裂缝面法向为参考构建新的方位各向异性反演机制，通过消除各向同性项、保留各向异性项的方式进行反演；与常规方法相比，本发明方法能够显著提高各向异性参数反演的稳定性和精度。此外，通过岩石物理理论，将反演得到的各向异性参数转换为更具物理含义的裂缝密度参数，能够直观表征致密砂岩气储层裂缝发育程度；并且，提出综合预测因子，同时考虑储层裂缝与流体，能够为致密砂岩气储层评价提供重要依据。
附图说明
[0027]
为了更清楚地说明本发明实施例的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，应当理解，以下附图仅示出了本发明的某些实施例，因此不应被看作是对范围的限定，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他相关的附图。
[0028]
图1是本发明实施方法流程图；
[0029]
图2是hti介质地震方位反射示意图；
[0030]
图3是入射角为25
°
时，反射系数权重随方位角的变化；
[0031]
图4a和图4b是反射系数和方位振幅差随方位角和入射角的变化图；
[0032]
图4c和图4d是入射角为15
°
、20
°
和25
°
时反射系数和方位振幅差随方位角变化的极坐标图；
[0033]
图5是入射角为25
°
时，过c井的六个方位叠加地震剖面；
[0034]
图6a为裂缝密度e连井剖面；图6b为压缩率1/k连井剖面；
[0035]
图6c为裂缝密度e与压缩率1/k乘积归一化的连井剖面；
[0036]
图7a为裂缝密度e平面图；
[0037]
图7b为压缩率1/k平面图；
[0038]
图7c为裂缝密度e与压缩率1/k乘积归一化的平面图。
具体实施方式
[0039]
下面结合具体实施案例对本发明作进一步的说明，本实施案例在以本发明技术为前提下进行实施，给出了详细的实施方式，但本发明的保护范围不限于以下的实施例。
[0040]
本发明基于地震反射振幅方位各向异性差异的裂缝预测方法，即在地震avaz反射系数方法的基础上通过与参考方位数据做差消除各向同性部分的影响，得到方位反射系数差的反演公式并对各向异性参数进行反演，再将各向异性参数转换成能直接表征裂缝的裂缝密度参数，最后，在裂缝预测与弹性反演压缩率预测的基础上，构建用于致密砂岩气储层裂缝及流体综合评价新的指示因，本发明提出的方法能充分利用各向异性部分的响应提高了反演过程的稳定性和精度，并对致密砂岩气储层裂缝发育及流体有很好的估计，为储层评价提供重要依据。
[0041]
具体原理如下：
[0042]rü
ger(1997,1998)提出弱各向异性hti介质pp波反射系数近似公式：
[0043][0044]
其中，θ表示入射角，为方位角,v
p
和vs分别表示p波和s波速度。z＝ρv
p
表示p波阻抗，g＝ρv
s2
表示横波模量。ε
(v)
,γ和δ
(v)
表示thomsen各向异性参数(thomsen,1986)。符号
ˉ
和
△
分别表示上层和下层介质性质的平均值和差值。
[0045]
thomsen各向异性参数的表达式为：
[0046][0047][0048][0049]
其中，c
ij
表示hti介质的弹性刚度系数
[0050]
将公式(1)进行线性化，由纵波速度、横波速度、密度和各向异性参数表示：
[0051][0052]
令：
[0053][0054]
公式(5)简化为：
[0055][0056]
hti介质pp波反射系数公式可以表示为各向同性项和各向异性项之和：
[0057][0058]
其中
[0059][0060][0061]
由于公式(6)中各向同性项和各向异性项加权系数的大小不同，直接同时反演可能会导致各向异性参数估计的不稳定性。为了提高各向异性参数估计的鲁棒性，通过减去参考方位反射系数消除公式(7)中各向同性项的影响，利用方位反射系数差反演各向异性参数的公式表示为：
[0062][0063]
公式(11)利用垂直于裂缝走向方位上的反射振幅与不同方位角的反射振幅之差来消除公式(7)中的各向同性项。其中表示观测方位，表示垂直于裂缝走向的方位。
[0064]
对于m个入射角θ、n个方位角的叠前方位地震集，根据公式(11)表示的各向异性参数反演矩阵的形式为：
[0065][0066]
其中，ai＝sin2θi,bi＝tan2θi,i＝1,
…
,m和j＝1,
…
,n，m和n分别为地震数据中入射角、方位角道集个数。
[0067]
公式(12)的简化形式为：
[0068]
d＝gx
ꢀꢀꢀ
(13)
[0069]
其中
[0070][0071]
利用阻尼最小二乘法对公式(13)进行求解，反演得到各向异性参数：
[0072]
x＝(g
t
·
g+σ2i)-1
·gtdꢀꢀꢀ
(14)
[0073]
其中，σ表示阻尼因子，i是单位矩阵。
[0074]
由式公式(11)-公式(14)可知，利用方位反射系数差反演各向异性参数的初始步骤是估算垂直于裂缝走向的方位ikelle(1996,1997)针对弱各向异性情况，将反射系数公式按方位角展开成傅里叶级数形式：
[0075][0076]
其中，系数d、en和fn(n＝2,4,
…
)是入射角的函数。
[0077]
根据三角函数关系：
[0078][0079]
公式(15)可整理成反射系数的级数形式：
[0080][0081]
公式(17)中的相位项可以表示为：
[0082][0083]
根据downton和roure(2015)，垂直于裂缝走向的方位可以由相位项估计
[0084][0085]
schoenberg和douma(1998)在线性滑动模型和hudson模型的基础上，提出了在描述裂缝时具有更明确物理性质的裂缝弱度参数，其表达式为：
[0086]
[0087][0088]
其中，δn和δ
t
分别表示法向和切向裂缝弱度，a和c是便士状裂隙或裂缝的长半轴和短半轴，e为裂缝密度，k’和μ’为裂缝充填物的体积模量和剪切模量，μ为母岩的剪切模量，g为s波和p波速度比的平方。
[0089]
bakulin(2000)给出裂缝弱度和各向异性参数之间的关系为：
[0090][0091]
δ
t
＝2γ
ꢀꢀꢀ
(23)
[0092]
切向裂缝弱度δ
t
可以进一步与裂缝密度相关：
[0093][0094]
本发明基于地震反射振幅方位各向异性差异的裂缝预测方法，包括以下步骤：
[0095]
步骤a、对宽方位地震数据进行方位角道集处理将其转换为分角度、分方位叠加数据，作为基于地震反射振幅方位各向异性差异的裂缝预测方法的输入数据。
[0096]
步骤b、将ruger的hti介质pp波方位反射系数公式进行线性化，由纵波速度、横波速度、密度和各向异性参数等线性表达，选取裂缝面法向为参考方位，构建方位反射系数差公式，在消除各向同性项影响的同时保留待反演的各向异性项，增强各向异性参数反演稳定性和精度。
[0097][0098]
其中，θ表示入射角，表示观测方位，表示垂直于裂缝走向的方位，v
p
和vs分别ε
(v)
,γ和δ
(v)
表示thomsen各向异性参数，符号
ˉ
和
△
分别表示上层和下层介质性质的平均值和差值。
[0099]
步骤b中，方位反射系数差方法与传统的反射系数方法不同，该方法构建以裂缝法向为参考方位的方位反射系数差表达式，通过消除地震反射中与裂缝响应无关的各向同性项，充分提取利用振幅方位各向异性差异，提高了各向异性参数及裂缝预测的稳定性和精度。
[0100]
步骤c、基于步骤b所述公式，通过理论模型计算并分析线性化的pp波方位反射系数，以及方位反射系数对储层裂缝的敏感性。
[0101]
步骤c，基于线性化pp波反射系数公式以及方位反射系数差公式的理论模型测试，分析了各向同性、各向异性项权重系数随方位角的变化规律，明确了在裂缝面方向附近，各向异性项权重系数比各向同性项权重系数的数量级小，该差异会造成的反演结果的不稳定
性，同时验证了本发明提出的方位反射系数差方法通过消除各向同性项避免同时反演各向同性和各向异性参数的合理性与有效性。
[0102]
步骤d、基于步骤b推导的的方位反射系数差公式构建方位地震反演方法，计算各向异性参数。
[0103][0104]
其中，ai＝sin2θi,bi＝tan2θi,i＝1,
…
,m和j＝1,
…
,n，m和n分别为地震数据中入射角、方位角道集个数。
[0105]
步骤d，基于本发明的方位反射系数差公式构建的方位地震反演方法，在裂缝方位预测的基础上，为hti各向异性参数的稳定与高精度反演提供了理论基础。
[0106]
步骤e、基于傅里叶级数展开的方法对步骤a得到的方位角道集数据预测裂缝方位，并以裂缝面法向为参考方位构建地震方位差道集，作为步骤d反演方法的输入数据。
[0107]
步骤f：基于步骤d构建的方位地震反演方法、输入步骤e中构建的地震方位差道集为输入数据进行反演，计算含裂缝的致密砂岩气储层hti各向异性参数。
[0108]
步骤f，基于本发明构建的方位地震反演方法，对于m个入射角、n个方位角的方位地震数据进行反演，为致密砂岩气储层hti各向异性参数反演提供了方法。
[0109]
步骤g、基于schoenberg线性滑动和hudson币状裂缝岩石物理模型，将步骤f反演得到的各向异性参数转换为裂缝切向弱度参数δ
t
，并根据bakulin(2000)的理论，将切向裂缝弱度δ
t
进一步转换为裂缝密度参数e，直观描述致密砂岩气储层裂缝发育程度。
[0110]
δ
t
＝2γ
ꢀꢀꢀ
(23)
[0111][0112]
其中，g为s波和p波速度比的平方。
[0113]
步骤g，基于schoenberg线性滑动和hudson币状岩石物理模型，根据各向异性参数与裂缝弱度、以及裂缝弱度与裂缝密度的关系，将反演得到的各向异性参数转换为裂缝密度，描述裂缝的发育程度的直观描述。
[0114]
步骤h、将步骤g得到的裂缝密度e与地震弹性反演得到的压缩率1/k相乘并归一化，构建用于致密砂岩气储层裂缝及流体综合评价的指示因子e
×
1/k。其中，裂缝密度e反映致密砂岩气储层裂缝发育程度，压缩率1/k反映致密砂岩气储层含气性。
[0115]
步骤h，在裂缝预测与弹性反演压缩率预测的基础上，构建用于致密砂岩气储层裂缝及流体综合评价新的指示因子e
×
1/k，其中裂缝密度e反映致密砂岩气储层裂缝发育程度，压缩率1/k反映致密砂岩气储层含气性。通过测井及储层开发数据标定，新的综合评价因子能够为致密砂岩气储层评价提供重要依据。
[0116]
步骤i、由测井数据及储层开发数据对步骤h预测结果进行标定。
[0117]
实施例1
[0118]
理论模型：设计了一个发育垂直定向排列裂缝的致密砂岩层被各向同性泥岩层覆盖的模型，致密砂岩层模型等效为hti介质，模型示意图如图2所示，表1为模型参数。在入射角固定为25
°
，方位角从0
°
到180
°
变化的情况下，计算公式(7)中的a-e的权重系数，如图3所示。各向同性项中的权重系数a、b和c是入射角θ的函数不随方位角变化，而各向异性项的权重系数d、e和f与入射角θ和方位角都有关，权重系数随方位角的变化而变化，并且在方位角接近90
°
(裂缝方向)时，其幅度要小得多。因此，直接用公式(7)同时反演各向同性性质和各向异性参数反演结果可能会不稳定，特别是涉及到靠近裂缝方向的地震数据时，由于各向同性项和各向异性项的加权系数数量级不同，造成反演结果的不稳定性。
[0119]
为了提高各向异性参数反演的鲁棒性，本发明提出了一种新的反演方法，消除各向同性项，从而减小公式(7)中各权重系数大小不同的影响。分别利用公式(5)和公式(11)计算pp波反射系数和相应的方位反射系数差(图4a)随入射角和方位角的变化(图4b)。极坐标图显示了入射角分别为15
°
、20
°
和25
°
时反射系数随方位角的变化裂缝，裂缝走向和垂直于裂缝走向的方向也在图中标出(图4c和图4d)。从图中可以看出入射角越大，反射系数和反射系数差异随方位角的变化就越明显。此外，反射系数差的方位变化(图4d)比反射系数的方位变化(图4c)更明显，这表明着通过消除各向异性背景的影响，可以提高方位反射系数差对各向异性程度或裂缝密度的敏感性。因此，利用方位反射系数差能够更可靠地估计各向异性参数或裂缝密度。
[0120]
实施例2
[0121]
实际应用：图1为本发明实施方法流程图。实际数据资料为自鄂尔多斯盆地致密砂岩气储层，根据地质背景资料和fmi成像测井资料显示，该地区发育近垂直的高角度裂缝，因此可以将该地区致密砂岩等效为hti介质。对地震数据进行整理得到三个入射角和六个方位角的方位地震叠加数据，三个入射角分别为7
°
、15
°
和25
°
，六个方位角分别为40
°
、60
°
、80
°
、100
°
、120
°
、140
°
。图5展示了过产气井c的方位叠加地震剖面，其中入射角固定为25
°
，不同方位角叠加地震剖面不同。
[0122]
在反演之前根据傅里叶级数展开的方法对裂缝方位进行估计，并对地震数据进行校正，以保证反演结果的准确性。将18个方位数据作为输入值进行反演，得到各向异性参数。根据公式(22)和公式(23)由各向异性参数计算得到法向裂缝弱度和切向裂缝弱度。利用公式(24)进一步得到裂缝密度，来描述裂缝的发育程度。
[0123]
图6是过井a、b、c和d的连井剖面。图6a是裂缝密度连井剖面图，层位和测井渗透率曲线在图中标出，计算得到的裂缝密度与测井渗透率存在一定程度的相关性，井b和井c在致密砂岩含气层he4段具有较高的渗透率和裂缝密度，同时4口井的测井渗透率与预测裂缝密度在部分位置也表现出明显的相关性。这说明垂向裂缝的发育可以提高研究区致密砂岩的渗透率。然而，由于垂向裂缝并不是决定储层渗透率的唯一因素，测井渗透率与裂缝密度不一定完全吻合。
[0124]
由于致密砂岩的压缩率1/k是识别天然气的良好指标，其中k代表体积模量。高1/k值(低k值)一般对应高致密砂岩含气饱和度sg值。因此，利用1/k连井剖面(图6b)来估算致密砂岩气中sg的分布，其中黑色曲线为含气饱和度测井曲线。井b和井c均为产气井对应的
1/k值较高，并且与sg的测井数据一致。然而，天然气饱和度和天然裂缝的发育程度对致密砂岩气藏的有效开发均具有重要意义。因此，如图6c连井剖面所示，本发明提出了一种同时考虑裂缝密度和含气饱和度影响的归一化参数e
×
1/k，用于测量致密砂岩气有效开发潜力。新的归一化参数e
×
1/k可以更好地预测潜在的天然气面积分布。
[0125]
在致密砂岩气he4段预测的裂缝密度e平面图(图7a)和渗透率1/k平面图(图7b)上，标出了直井a-h的平台位置及由平台b、c和e延伸的斜井位置，并绘制了直井a、b、c和d的井位置连线，其中红点表示产气井，白点表示干井。图7a与图7b不完全吻合说明裂缝发育的分布(图7a)与致密砂岩含气的分布(图7b)不一定一致。同时，与含气饱和度相比，he4段致密砂岩裂缝密度的空间变化表现出更强的非均质性。图7a和图7b表明，与干井相比，产气井的e值和1/k值普遍较高(sg值较高)。本发明综合考虑裂缝发育程度和含气性的影响，得到e
×
1/k的归一化平面图(图7c)，综合因子e
×
1/k的归一化与开采井产量有良好的吻合，能识别出he4段致密砂岩的产气区，预测结构为致密砂岩气藏的开发部署提供重要信息。
[0126]
本发明首先对各向异性hti介质pp波反射系数近似公式进行线性化，推导由纵波速度、横波速度、密度和各向异性参数表征的线性化公式。在此基础上，构建以裂缝法向为参考方位的方位反射系数差表达式，通过消除地震反射中与裂缝响应无关的各向同性项，充分提取利用振幅方位各向异性差异，提高了各向异性参数及裂缝预测的稳定性和精度。通过理论模型测试明确了在裂缝面方向附近，各向异性项权重系数比各向同性项权重系数的数量级小，该差异会造成的反演结果的不稳定性，同时验证了本发明提出的方位反射系数差方法的合理性与有效性。并且，应用schoenberg线性滑动和hudson币状岩石物理理论，将本发明反演的hti各向异性参数转换为裂缝密度参数，直观描述裂缝的发育程度。最后，在裂缝预测与弹性反演压缩率预测的基础上，构建用于致密砂岩气储层裂缝及流体综合评价新的指示因子e
×
1/k，其中裂缝密度e反映致密砂岩气储层裂缝发育程度，压缩率1/k反映致密砂岩气储层含气性。通过测井及储层开发数据标定，新的综合评价因子能够为致密砂岩气储层评价提供重要依据。
[0127]
注意，上述仅为本发明的较佳实施例及所运用技术原理。本领域技术人员会理解，本发明不限于这里所述的特定实施例，对本领域技术人员来说能够进行各种明显的变化、重新调整和替代而不会脱离本发明的保护范围。因此，虽然通过以上实施例对本发明进行了较为详细的说明，但是本发明不仅仅限于以上实施例，在不脱离本发明构思的情况下，还可以包括更多其他等效实施例，而本发明的范围由所附的权利要求范围决定。