石油物探  2020, Vol. 59 Issue (1): 108-113  DOI: 10.3969/j.issn.1000-1441.2020.01.012
0
文章快速检索     高级检索

引用本文 

宋维琪, 徐月森, 张云银, 等. 基于傅里叶级数分析的各向异性参数估计及裂缝预测[J]. 石油物探, 2020, 59(1): 108-113. DOI: 10.3969/j.issn.1000-1441.2020.01.012.
SONG Weiqi, XU Yuesen, ZHANG Yunyin, et al. Anisotropy parameter estimation and fracture prediction based on Fourier series analysis[J]. Geophysical Prospecting for Petroleum, 2020, 59(1): 108-113. DOI: 10.3969/j.issn.1000-1441.2020.01.012.

基金项目

渤海湾盆地济阳坳陷致密油开发示范工程(2017ZX05072)任务一“致密油藏储层地震预测方法及地应力研究”资助

作者简介

宋维琪(1964—), 男, 教授, 博士生导师, 主要从事地球物理方法教学和科研工作。Email:swq1123@126.com

文章历史

收稿日期:2019-01-14
改回日期:2019-08-12
基于傅里叶级数分析的各向异性参数估计及裂缝预测
宋维琪1 , 徐月森1 , 张云银2 , 张营革2 , 高秋菊2 , 魏欣伟2 , 张明秀2     
1. 中国石油大学(华东)地球科学与技术学院, 山东 青岛 266555;
2. 中国石油化工股份有限公司胜利油田分公司物探研究院, 山东 东营 257055
摘要:利用Ruger反射系数近似公式进行地震方位各向异性椭圆拟合, 研究裂缝密度与走向, 在实际应用时因为待拟合点中存在奇异点, 且椭圆拟合存在90°旋转问题, 因而计算的方位角无法准确指示椭圆的长短轴方向。为此, 将HTI模型纵波反射系数近似公式展开成傅里叶级数形式, 利用最小二乘拟合提取傅里叶级数系数, 从而提高各向异性参数估计和裂缝预测的准确性。针对不同方位角范围的数据, 利用二倍角公式将纵波反射系数傅里叶级数展开式进一步展开, 避免拟合提取傅里叶级数系数时因为同一个cosnφ项对应两个各向异性反射系数而导致拟合结果出错。将该方法应用于HTI介质模型和实际地震资料, 结果表明:对于HTI介质模型, 展开后提取的拟合傅里叶级数系数与理论公式计算结果一致; 对于实际地震资料, 展开后提取傅里叶级数系数与工区裂缝发育及井资料指示的结果吻合, 效果优于椭圆拟合结果和未展开时的傅里叶级数系数拟合结果。
关键词裂缝预测    各向异性参数    傅里叶级数展开    最小二乘拟合    方位各向异性    椭圆拟合    HTI模型    
Anisotropy parameter estimation and fracture prediction based on Fourier series analysis
SONG Weiqi1, XU Yuesen1, ZHANG Yunyin2, ZHANG Geying2, GAO Qiuju2, WEI Xinwei2, ZHANG Mingxiu2     
1. School of Geosciences, China University of Petroleum, Qingdao 266555, China;
2. Geophysical Research Institute of Shengli Oilfield Branch of Sinopec, Dongying 257055, China
Foundation item: This research is financially supported by the Task 1 of Demonstration Project for Tight Oil Development in Jiyang Depression, Bohai Bay Basin (Grant No.2017ZX05072)
Abstract: In practical applications, an approximated formula of Ruger reflection coefficient is used for the elliptic fitting of seismic azimuth anisotropy; this is aimed at studying the density and direction of the fracture in rock masses.As a result, because there are singularities in the fitting points and the elliptic fitting has a 90° rotation issue, the calculated angle of azimuth does not accurately indicate the direction of the elliptical axis.In this study, the approximate formula for the P-wave reflection coefficient of a horizontal transverse isotropic (HTI) medium was expanded through Fourier series; the coefficient of which was extracted by least square fitting.This procedure improved the accuracy in the estimation of the anisotropic parameter, and thus provided a better crack prediction.For data with different azimuth ranges, the Fourier series expansion of the P-wave reflection coefficient was further expanded with the double angle formula; this was to avoid incorrect fitting results caused by the same items corresponding to two values.Further, the method was applied to a HTI model and actual seismic data.The results for the HTI model showed that the fitting Fourier series coefficients obtained after the expansion were consistent with the predicted results of the theoretical formula.For the actual seismic data, the coefficients obtained after the expansion were in agreement with the observations of fracture development, and wells in the applied working area.Furthermore, the results were more accurate than those obtained through elliptic fitting or from the Fourier series coefficients without expansion; thereby it was demonstrated that the proposed approach exhibited a superior capability.
Keywords: fracture prediction    anisotropic parameter    Fourier series expansion    least squares fitting    azimuthal anisotropy    ellipse fitting    HTI model    

随着各向异性介质理论研究的发展与完善, 非常规致密油气储层作为当前油气勘探的重点领域之一, 相关的地震预测技术也在不断深入发展, 并逐渐应用到生产实际领域[1-2]。实际地球物理介质都不同程度地存在着各向异性的特性。作为非常规致密油气藏, 裂缝储层是典型的各向异性介质[3-4]。地下地层中发育的近垂直裂缝, 不仅是良好的油气运输通道, 也是良好的油气存储空间, 因此研究裂缝储层具有重要的实用价值[5]

各向异性介质理论在19世纪末被引入地球物理勘探领域, 并在近二、三十年得到了长足的发展与应用。CRAMPIN[6-7]提出方位各向异性和横波分裂等概念, 极大地推动了各向异性理论在应用领域的发展; THOMSEN[8], HUDSON[9]先后提出不同的各向异性理论模型, 对该理论的实际应用做出了重要贡献; RUGER[10]提出具有水平对称轴的横向各向同性介质(HTI)纵波反射系数近似方程, 推动了利用纵波方位各向异性信息预测裂缝技术的发展; 王顺昌等[11]研究了倾斜裂缝地层qP波方位反射系数椭圆特征, 将qP波方位反射系数拟合成椭圆, 得到了椭圆参数与裂缝密度、裂缝倾向和裂缝倾角的关系; 杨敏等[12]研究了裂缝参数对地震方位各向异性特征的影响, 发现HTI介质、裂缝介质的各向异性特征随裂缝密度增加而增强, P波速度、反射系数随方位角呈现椭圆、周期性变化, 各向异性椭圆方位指示裂缝发育方向; 王洪求等[13]分析了不同地震属性的方位各向异性并展开裂缝预测; 杨勤勇等[14]研究了纵波方位各向异性及其在裂缝检测中的应用, 证实了垂直裂缝存在方位各向异性特征, 提出了基于二维多方位预测的精确解法和基于三维多方位预测的最小二乘法; 王康宁等[15]研究了基于HTI模型的纵波方位各向异性反射系数近似公式的傅里叶级数展开式; 多位学者研究了具有两组正交对称轴的各向同性的纵波反射系数近似公式的傅里叶级数展开[16-20]。本文将纵波反射系数近似公式展开成傅里叶级数形式, 并针对不同方位角范围的数据, 利用二倍角公式进一步展开, 避免了在最小二乘拟合条件下提取傅里叶级数系数, 造成同一个cosnφ项对应两个各向异性反射系数而导致的拟合结果非唯一, 将该傅里叶级数系数拟合方法应用于HTI模型和实际地震资料, 可实现各向异性参数估计和裂缝预测。

1 基于傅里叶级数分析的各向异性参数预测 1.1 基于HTI模型的纵波反射系数公式傅里叶级数展开式

随地震观测方位发生变化的纵波反射系数序列可以表示为周期为2π的函数, 傅里叶级数展开形式为[15]:

$ R(\varphi)=a_{0}+\sum\limits_{n=1}^{\infty}\left[a_{n} \cos n \varphi+b_{n} \sin n \varphi\right] $ (1)

式中:φ为方位角; a0, an, bn(n=1, 2)为常数。假设裂缝介质为HTI介质, 结合SCHOENBERG等[21]裂缝介质理论模型, 得到傅里叶级数形式的纵波方位各向异性反射系数RPP(i, φ)的近似公式[15]为:

$ R_{\mathrm{PP}}(i, \varphi)=a_{0}+a_{2} \cos 2 \varphi+a_{4} \cos 4 \varphi $ (2)

式中:a2=1/2{g[ΔδT-(1-2gδN]+g(g-1)ΔδNtan2i}sin2i; a4=1/8g2δT-gΔδN)sin2itan2i; δT, δN分别表示切向弱度和法向弱度[22]; $g=\bar{\beta}^{2} / \bar{\alpha}^{2}$; ,分别表示反射界面附近平均纵波速度、横波速度; i表示入射角。

传统的地震属性方位各向异性椭圆拟合方法通常基于Ruger两项近似方程, 当选定某一个固定的入射角时, 将纵波反射系数公式进行傅里叶级数展开, 提取二阶项系数a2, 可直接得到裂缝发育密度。

1.2 基于二倍角的纵波反射系数公式傅里叶级数展开式

利用纵波反射系数公式的傅里叶级数展开式预测裂缝发育密度和走向, 通常采用反演方法或拟合方法。考虑到实际地震资料信噪比低的特点, 本文采用拟合方法估计各向异性参数。

当地震数据方位角的范围为0~45°时, 针对(2)式, 可采用最小二乘拟合方法计算傅里叶级数系数a2; 但对于方位角范围在0~90°的地震数据, 采用(2)式拟合a2并不合适, 因为此时一个cos4φ项在拟合时会对应两个RPP值, 导致拟合结果非唯一。因此需要利用二倍角公式将(2)式展开成如下形式再进行拟合:

$ \begin{aligned} R_{\mathrm{PP}}(i, \varphi) &=a_{0}+a_{2} \cos 2 \varphi+a_{4} \cos 4 \varphi \\ &=2 a_{4} \cos ^{2} 2 \varphi+a_{2} \cos 2 \varphi+\left(a_{0}-a_{4}\right) \end{aligned} $ (3)

针对实际应用中方位角范围为0~180°的地震数据, 需要将(3)式继续用二倍角公式展开, 否则在拟合时, 一个cos2φ项也会对应两个振幅值, 导致拟合结果不确定。如图 1所示, 拟合傅里叶级数系数a2时, 一个cos2φ项本应对应一个振幅值, 但因为数据方位角范围扩大至0~180°, 使得一个cos(2φ)项对应两个振幅值, 导致拟合结果存在非唯一性。图中蓝色曲线代表 0~180°方位角范围的实际地震数据, 小圆圈代表 9个不同方位的振幅值, 紫色曲线代表拟合数据。

图 1 未用二倍角公式展开时错误的拟合结果

根据二倍角公式将(3)式展开为:

$ \begin{array}{l} {R_{{\rm{PP}}}}(i, \varphi ) = 8{a_4}{\cos ^4}\varphi + \left( {2{a_2} - 8{a_4}} \right){\cos ^2}\varphi + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\left( {{a_0} - {a_2} + {a_4}} \right) \end{array} $ (4)

对于方位角范围在0~180°的数据, 一个cosφ项对应一个RPP值, 拟合结果的唯一性得到保障。

2 HTI介质理论模型试算

HTI介质是具有水平对称轴的横向各向同性介质, 弹性系数矩阵由5个独立的弹性常数表示。本文采用的双层HTI介质理论模型各向异性参数如表 1所示。

表 1 双层HTI介质理论模型各向异性参数

各向同性面平行于纵轴, i为入射角。越靠近各向异性对称轴, 介质的各向异性就会越强。利用最小二乘拟合的方法, 根据双层HTI介质理论模型拟合a2, 结果如图 2a所示, 小圆圈代表不同入射角时的a2值, 紫色曲线为拟合曲线, 可以看出随着入射角的增加, a2的绝对值不断增大, 拟合得到的a2变化趋势与各向异性参数理论公式计算结果(图 2b)一致, 说明本文方法能有效地反映介质的各向异性特征。

图 2 傅里叶级数系数a2随入射角变化情况 a本文方法拟合; b理论公式计算
3 实际应用

将本文方法应用于渤南地区四扣工区致密砂岩的储层研究。研究区内地质构造复杂, 断层、裂缝等发育; 致密砂岩储层基质颗粒杂乱, 分选性差, 渗透率低, 非均质性强, 岩性变化大; 地层的形成和演化过程非常复杂, 各向异性特征明显。利用基于两项Ruger方程的地震属性方位各向异性椭圆拟合方法即传统方法, 对比验证本文方法的有效性和精度。

渤南地区地震资料包含9个方位角(0~180°)数据, 在远、中、近3个固定偏移距上各自形成了9个分方位角道集叠加的地震数据, 每隔20°划分一个方位角道集, 提取每个方位角的平均振幅属性, 将每一个点的振幅属性值A按照其方位角分解为:

$ \left\{\begin{array}{l} {A_{x}=A \cos \varphi} \\ {A_{y}=A \sin \varphi} \end{array}\right. $ (5)

式中:φ为该方位角道集的中心方位角; Ax, Ay分别为沿横轴、纵轴方向的振幅属性值。

因近偏移距地震数据各向异性不明显, 远偏移距地震数据的信噪比低而且数据不完整, 故本文使用中偏移距的地震数据进行计算。首先将某一个坐标点9个方位的每一个地震振幅属性值都分别分解为AxAy, 然后进行振幅属性方位各向异性椭圆拟合, 最后得到的椭圆拟合长短轴之比的结果如图 3所示。采用地震属性方位各向异性椭圆拟合方法预测裂缝发育密度和走向时, 对属性值波动大、信噪比低的数据, 得到的计算结果精度低, 抗噪能力差。

图 3 中偏移距地震数据椭圆拟合长短轴之比的结果

图 4可知该地区发育北西—南东向的断层, 由图 5可知, 本文方法拟合结果与实际一致, 说明本文方法在预测裂缝发育密度和走向时计算精度高, 抗噪能力强, 效果优于传统的椭圆拟合方法。采用本文方法对该地区地震资料进行计算, 拟合傅里叶级数第2项系数a2, 对比未利用二倍角公式展开和利用二倍角公式展开后的傅里叶级数系数a2拟合结果(图 5)发现, 后者的拟合精确明显提高, 更为准确地实现了各向异性参数估计和裂缝预测。

图 4 义170井裂缝玫瑰图
图 5 利用中偏移距地震资料拟合的傅里叶级数系数a2 a未用二倍角公式展开; b利用二倍角公式展开后
4 结论

将HTI理论模型纵波反射系数公式的傅里叶级数展开式, 利用二倍角公式进一步展开, 避免了同一个cosnφ项对应两个RPP值导致的拟合结果非唯一性问题, 提高了基于拟合的各向异性参数估计和裂缝预测的准确性。本文方法在HTI模型和实际地震资料应用结果表明, 相较于传统的椭圆拟合方法, 本文方法抗噪能力更强, 计算精度更高, 更为准确地实现了渤南地区各向异性参数估计和裂缝发育的预测。

参考文献
[1]
刘建英.三维方位各向异性分析及裂缝检测[D].成都: 成都理工大学, 2010
LIU J Y.Three-dimensional azimuthal anisotropy analysis and crack detection[D].Chengdu: Chengdu University of Technology, 2010
[2]
吴萍, 杨长春, 王真理. HTI介质中的反射纵波方位属性[J]. 地球物理学进展, 2009, 24(3): 944-950.
WU P, YANG C C, WANG Z L. Azimuth properties of reflected p-wave in HTI medium[J]. Progress in Geophysics, 2009, 24(3): 944-950. DOI:10.3969/j.issn.1004-2903.2009.03.016
[3]
孙炜, 何治亮, 李玉凤. 改进的方位各向异性裂缝预测方法及其应用[J]. 石油地球物理勘探, 2014, 49(6): 1170-1178.
SUN W, HE Z L, LI Y F. Improved azimuthal anisotropic fracture prediction method and its application[J]. Oil Geophysical Prospecting, 2014, 49(6): 1170-1178.
[4]
李雷涛, 肖秋红, 肖伟. 优化的方位各向异性裂缝预测方法及应用[J]. 断块油气田, 2016, 23(4): 455-459.
LI L T, XIAO Q H, XIAO W. Optimized azimuthal anisotropic fracture prediction method and its application[J]. Fault Block Oil and Gas Field, 2016, 23(4): 455-459.
[5]
程冰洁, 张玉芬. AVO简化方程的物理意义及其在油气识别中的应用[J]. 物探化探计算技术, 2003, 25(1): 26-30.
CHENG B J, ZHANG Y F. Physical significance of simplified AVO equation and its application in oil and gas identification[J]. Geophysical and Geochemical Exploration Calculation Technology, 2003, 25(1): 26-30. DOI:10.3969/j.issn.1001-1749.2003.01.006
[6]
CRAMPIN S. A review of wave motion in an aniso-tropic and cracked elastic media[J]. Wave Motion, 1981, 3(4): 343-391. DOI:10.1016/0165-2125(81)90026-3
[7]
CRAMPIN S. Evaluation of anisotropy by shear-wave splitting[J]. Geophysics, 1985, 50(1): 142-152. DOI:10.1190/1.1441824
[8]
THOMSEN L. Weak elastic anisotropy[J]. Geophysics, 1986, 51(10): 1954-1966. DOI:10.1190/1.1442051
[9]
HUDSON J A. Wave speeds and attenuation of elastic waves in material containing cracks[J]. Geophysical Journal International, 1981, 64(1): 133-150. DOI:10.1111/j.1365-246X.1981.tb02662.x
[10]
RUGER A. Variation of P-wave reflectivity with offset and azimuth in anisotropic media[J]. Geophysics, 1996, 61(3): 935-947.
[11]
王顺昌, 田景春, 贾瀛. 倾斜裂缝地层qP波方位反射系数椭圆特征研究[J]. 物探化探计算技术, 2014, 36(4): 452-457.
WANG S C, TIAN J C, JIA Y. Study on the elliptical characteristics of qP wave azimuth reflection coefficient in inclined fractured strata[J]. Geophysical and Geochemical Exploration Calculation Technology, 2014, 36(4): 452-457. DOI:10.3969/j.issn.1001-1749.2014.04.14
[12]
杨敏, 李明, 兰峰. 裂缝参数对地震方位各向异性特征的影响[J]. 重庆科技学院学报(自然科学版), 2013, 15(2): 31-34.
YANG M, LI M, LAN F. Effects of fracture parameters on seismic azimuth anisotropy[J]. Journal of Chongqing University of Science and Technology (Natural Science Edition), 2013, 15(2): 31-34. DOI:10.3969/j.issn.1673-1980.2013.02.008
[13]
王洪求, 杨午阳, 谢春辉. 不同地震属性的方位各向异性分析及裂缝预测[J]. 石油地球物理勘探, 2014, 49(5): 925-931.
WANG H Q, YANG W Y, XIE C H. Azimuthal anisotropy analysis and fracture prediction of different seismic attributes[J]. Oil Geophysical Prospecting, 2014, 49(5): 925-931.
[14]
杨勤勇, 赵群, 王世星. 纵波方位各向异性及其在裂缝检测中的应用[J]. 石油物探, 2006, 45(2): 177-181.
YANG Q Y, ZHAO Q, WANG S X. Compressional azimuth anisotropy and its application in fracture detection[J]. Geophysical Prospecting for Petroleum, 2006, 45(2): 177-181. DOI:10.3969/j.issn.1000-1441.2006.02.011
[15]
王康宁, 孙赞东, 侯昕晔. 基于傅里叶级数展开的纵波方位各向异性裂缝预测[J]. 石油物探, 2015, 54(6): 755-761.
WANG K N, SUN Z D, HOU X Y. Prediction of longitudinal wave azimuth anisotropic fractures based on Fourier series expansion[J]. Geophysical Prospecting for Petroleum, 2015, 54(6): 755-761. DOI:10.3969/j.issn.1000-1441.2015.06.014
[16]
BARONE A, SEN M K. A new Fourier azimuthal amplitude variation fracture characterization method:Case study in the Haynesville Shale[J]. Geophysics, 2018, 83(1): WA101-WA120. DOI:10.1190/geo2017-0030.1
[17]
ANTHONY B, MRINAL K S. Comparison of HTI and Orthorhombic methods for determining fracture density and fracture azimuth from 3D seismic data[J]. Expanded Abstracts of 85th Annual Internat SEG Mtg, 2015, 2916-2920.
[18]
CHEADLE S P, BROWN R J, LAWTON D C. Orthorhombic anisotropy:A Physical seismic modeling study[J]. Geophysics, 1991, 56(10): 1603-1613. DOI:10.1190/1.1442971
[19]
JONATHAN E D, BENJAMIN R. Interpreting azimuthal fourier coefficients for anisotropic and fracture parameters[J]. Interpretation, 2015, 3(3): ST9-ST27. DOI:10.1190/INT-2014-0235.1
[20]
BENJAMIN R. Statistical moments for azimuthal anisotropy characterization[J]. Expanded Abstracts of 84th Annual Internat SEG Mtg, 2014, 290-294.
[21]
SCHOENBERG M, SAYERS C M. Seismic anisotropy of fractured rock[J]. Geophysics, 1995, 60(1): 204-211. DOI:10.1190/1.1443748
[22]
DOWNTON J, RUSSELL H. Azimuthal Fourier coefficients:A simple method to estimate fracture parameters[J]. Expanded Abstracts of 81st Annual Internat SEG Mtg, 2011, 269-273.