2. 中国石油化工股份有限公司西北油田分公司研究院, 新疆乌鲁木齐 830011;
3. 中国石油集团东方地球物理勘探有限责任公司大庆物探一公司, 黑龙江大庆 163300
2. Research Institute of Petroleum Exploration & Production, Sinopec Northwest Company, Urumqi 830011, China;
3. No.1 Geophysical Exploration Company of BGP, CNPC, Daqing 163300, China
随着勘探开发的不断深入, 勘探目标逐渐转向“深、隐、薄、碎、小”的地质体, 构造解释和储层描述的难度越来越大, 对地震成像精度的要求越来越高, 地震资料的各向异性处理已经成为常规处理流程中不可或缺的重要环节。对具有垂向对称轴的横向各向同性(vertical transverse isotropy, VTI)介质进行各向异性处理主要解决了层状介质引起的垂向各向异性, 对具有水平对称轴的横向各向同性(horizontal transverse isotropy, HTI)介质进行各向异性处理主要解决了垂直裂缝引起的方位各向异性, 这两种处理方法均可改善成像聚焦性和减小井震深度误差。但当层状介质与垂直裂缝正交时, 地震波表现为正交各向异性特征, 造成速度描述困难、断裂成像精度低、绕射波不收敛和方位成像道集呈螺旋状弯曲等现象, 低品质的成像剖面对后续构造精细解释和储层精细描述影响很大。因此, 研究正交各向异性参数建模及地震成像非常必要[1-2]。
现阶段, 各向异性处理技术主要针对薄互层TTI和垂直裂缝HTI各向异性介质, 其理论发展较快, 相关技术较成熟, 并取得了良好的应用效果。正交各向异性处理技术的实用化发展相对滞后, 尽管国外大型石油公司的商业软件在近几年推出了相关处理模块, 但技术不够完善, 实际应用受限。TSVANKIN等[3-4]基于Thomsen参数的思路, 给出了HTI介质各向异性参数, 并建立了HTI介质各向异性参数估计思路; 随后TSVANKIN[5-7]提出了表征正交介质的各向异性参数, 介绍了各个参数的物理意义, 推动了正交各向异性处理技术的发展; 郭恺等[8]通过分析Thomsen各向异性参数的角度敏感性, 提出了多参数反演策略, 取得了良好的效果; BAKULIN等[9-11]对HTI各向异性介质、正交各向异性介质和单斜各向异性介质等进行了研究, 并利用地震反射波资料反演裂缝参数; GRECHKA等[12]在多分量宽方位采集数据的基础上, 利用地震反射波资料对正交介质各向异性参数进行了估计。
本文在前人研究的基础上, 首先研究了正交介质各向异性参数初始建模方法, 将数据进行分方位处理从而预测了裂缝的主要发育方向, 采用分方位速度分析和参数融合方法建立了空间分布合理、数值可靠的正交介质各向异性参数初始模型; 其次, 研究了正交介质各向异性参数精细建模方法, 引入方位正交各向异性介质层析公式提高了层析反演的稳定性, 制定了多参数反演策略, 在方位成像道集的基础上进行参数更新, 建立了精细的正交介质各向异性参数模型, 有效提高了断裂和缝洞储层的成像精度; 最后利用多层模型数据和某探区实际资料验证了本文方法的适用性和有效性。
1 方法原理正交各向异性介质可用9个各向异性参数描述:vP0, ε1, ε2, δ1, δ2, δ3, α,θ,φ, 其中, vP0是地震波垂向速度; ε1, δ1是平行于垂直裂缝方向的各向异性参数; ε2, δ2是垂直于垂直裂缝方向的各向异性参数; δ3是水平方向(即岩层方向)的各向异性参数; α,θ,φ分别是裂缝方位角、对称轴倾角、对称轴方位角。为了获取上述正交各向异性参数, 首先通过方位速度分析和参数提取建立了参数vP0, ε1, ε2, δ1, δ2, δ3, α的模型, 然后通过地层扫描得到参数θ,φ, 利用上述参数成像并提取方位道集后, 再进行层析反演进一步更新各向异性参数, 最后建立正交各向异性参数精细模型。
1.1 正交各向异性参数初始建模 1.1.1 方位速度分析及裂缝方位提取在实际应用中, 可根据工区前期的勘探情况, 将工区分为多个扇区(图 1), 在各个扇区内进行速度分析(图 2), 提取螺旋状成像道集所在地层的平均速度, 按方位排列, 拟合成椭圆形态(图 3)。椭圆的长轴平行于垂直裂缝, 由参数ε1, δ1描述; 短轴垂直于垂直裂缝, 由参数ε2, δ2描述。裂缝方向与零方位的夹角为α, 当α≠0时, 需要对坐标系统进行旋转, 公式如下:
$ \begin{array}{l} s = x\cos \alpha + y\sin \alpha \\ t = y\cos \alpha - x\sin \alpha \end{array} $ | (1) |
式中:(x,y)为直角坐标系中的某一点; (s,t)为旋转后坐标系中的某一点。
1.1.2 方位参数提取及初始建模采用上述方位速度分析方法, 按方位提取各向异性参数, 具体公式如下:
$ \left\{ \begin{array}{*{35}{l}} v_{\text{NMO}}^{-2}(\varphi )=\frac{{{\cos }^{2}}(\varphi -\alpha )}{v_{\text{P}, \text{NMO}}^{2}}+\frac{{{\sin }^{2}}(\varphi -\alpha )}{v_{\text{Po}}^{2}} \\ {{v}_{\text{P}, \text{NMO}}}={{v}_{\text{PO}}}\sqrt{1+2{{\delta }^{v}}} \\ \end{array} \right. $ | (2) |
式中:vNMO是表征各向异性的方位动校正速度; vP, NMO是表征各向异性的动校正速度垂向分量; vP0是纵波速度垂向分量; φ是方位角; δv是方位各向异性参数。
此时, 提取的方位各向异性参数描述的是扇区内各向异性, 不是正交各向异性, 故需要利用以下公式建立正交各向异性参数模型:
$ \left\{ {\begin{array}{*{20}{l}} {{\varepsilon _1} = {\varepsilon _{{\rm{vT1}}}} + {\varepsilon _{{\rm{HTI}}}}, {\varepsilon _2} = {\varepsilon _{{\rm{VTI}}}}}\\ {{\delta _1} = {\delta _{{\rm{VTI}}}} + {\delta _{{\rm{HTI}}}}, {\delta _2} = {\delta _{{\rm{VTI}}}}}\\ {{\delta _1} = {\delta _{{\rm{VTI}}}} + {\delta _{{\rm{HTI}}}}, {\delta _2} = {\delta _{{\rm{VTI}}}}}\\ {{\varepsilon _{{\rm{HTI}}}} = - \frac{{{\varepsilon _{{\rm{HTI}}}}}}{{\left( {1 + 2{\varepsilon _{{\rm{HTI}}}}} \right)}}}\\ {{\delta _{{\rm{HT}}1}} = - \frac{{{\delta _{{\rm{HT}}1}} - 2{\varepsilon _{{\rm{HT}}1}}\left( {1 + {\varepsilon _{{\rm{HTI}}}}} \right)}}{{{{\left( {1 + 2{\varepsilon _{{\rm{HTI}}}}} \right)}^2}}}}\\ {{\delta _3} = \frac{1}{2}\left[ {{{\left( {\sqrt {1 + 2{\delta _1}} \sqrt {1 + 2{\delta _2}} + 2\sqrt {{\varepsilon _1} - {\delta _1}} \cdot \frac{1}{{\sqrt {{\varepsilon _2} - {\delta _2}} }}} \right)}^2}/{{\left( {1 + 2{\varepsilon _2}} \right)}^2} - 1} \right]} \end{array}} \right. $ | (3) |
式中:εHTI, δHTI是平行于垂直裂缝方向的各向异性参数; εVTI, δVTI是垂直于垂直裂缝方向的各向异性参数。
通过公式(2)提取方位各向异性参数, 再通过公式(3)融合并建立正交各向异性参数模型, 可得到方位信息准确、空间分布合理的正交各向异性参数初始模型。
1.2 正交各向异性参数精细建模 1.2.1 正交各向异性层析矩阵的建立正交各向异性参数与方位相关, 因此不包含方位信息的常规层析反演矩阵不适于正交各向异性介质。我们记录不同方位的核函数和旅行时误差, 代入以下公式, 再根据不同方位的旅行时反演对应的参数:
$ \left\{ \begin{array}{*{35}{l}} \mathit{\boldsymbol{A}}(\varphi )\Delta \mathit{\boldsymbol{m}}=\Delta \tau (\varphi ),\left( \varphi =0,\cdots ,{{45}^{{}^\circ }},\cdots ,{{90}^{{}^\circ }} \right) \\ \mathit{\boldsymbol{A}}(\varphi )=\mathit{\boldsymbol{L}}(\varphi )\cdot \\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left[ \frac{{{\partial }_{{{S}_{\text{g}}}}}}{\partial {{s}_{\text{Po}}}}\frac{\partial {{s}_{\text{g}}}}{\partial {{\text{ }\!\!\varepsilon\!\!\text{ }}_{1}}}\frac{\partial {{s}_{\text{g}}}}{\partial {{\text{ }\!\!\varepsilon\!\!\text{ }}_{2}}}\frac{\partial {{s}_{\text{g}}}}{\partial {{\delta }_{1}}}\frac{\partial {{s}_{\text{g}}}}{\partial {{\delta }_{2}}}\frac{\partial {{s}_{\text{g}}}}{\partial {{\delta }_{3}}} \right] \\ \Delta \mathit{\boldsymbol{m}}=\left[ \Delta {{s}_{\text{Po}}}\quad \Delta {{\text{ }\!\!\varepsilon\!\!\text{ }}_{1}}\quad \Delta {{\text{ }\!\!\varepsilon\!\!\text{ }}_{2}}\quad \Delta {{\delta }_{1}}\quad \Delta {{\delta }_{2}}\quad \Delta {{\delta }_{3}} \right] \\ \end{array} \right. $ | (4) |
式中:A(φ)是方位角φ的核函数矩阵; Δm是待反演参数矩阵; Δτ(φ)是方位角φ的旅行时残差矩阵; sg是群速度; L(φ)为方位角φ的射线长度矩阵; Δ代表更新量。
1.2.2 正交各向异性多参数反演策略正交各向异性参数多, 数量级相差大, 因此同时反演所有参数非常不稳定。此外参数之间相互干扰, 反演精度也会受到影响。因此, 针对正交各向异性介质的多参数反演, 我们制定了如下的反演策略:
1) 通过正交各向异性偏移提取全方位道集, 根据分方位数据偏移结果提取分方位角道集;
2) 因为vP0参数对所有方位角道集的旅行时影响一致, 所以采用全方位道集反演vP0;
3) ε1, δ1只影响裂缝方向旅行时, ε2, δ2只影响垂直裂缝方向旅行时, 因此, 在准确反演vP0的基础上, 利用两个方向的道集分别反演ε1, δ1和ε2, δ2;
4) 裂缝方向和垂直裂缝方向的数据量小, 反演精度有限, 可采用除裂缝方向和垂直裂缝方向外的方位道集进行反演, 将采用参数扫描及全局寻优方法更新后的ε1, δ1和ε2, δ2代入如下的Thomsen参数ε(φ)和δ(φ)的计算公式:
$ \varepsilon (\varphi ) = {\varepsilon _1}{\sin ^4}\varphi + {\varepsilon _2}{\cos ^4}\varphi + \left( {2{\varepsilon _2} + {\delta _3}} \right) \\ \cdot {\sin ^2}\varphi {\cos ^2}\varphi \\ \delta (\varphi ) = {\delta _1}{\sin ^2}\varphi + {\delta _2}{\cos ^2}\varphi $ | (5) |
5) 将上述步骤得到的ε1, δ1和ε2, δ2, 代入公式(6)估算δ3:
$ {\delta _3} = \frac{1}{2}\left[ {{{(\sqrt {1 + 2{\delta _1}} \sqrt {1 + 2{\delta _2}} + 2\sqrt {{\varepsilon _1} - {\delta _1}} \cdot \sqrt {{\varepsilon _2} - {\delta _2}} )}^2}/{{\left( {1 + 2{\varepsilon _2}} \right)}^2} - 1{\rm{ }}} \right] $ | (6) |
我们采用正交各向异性多参数反演方法, 建立了正交各向异性参数精细模型, 为后续偏移成像提供了精确的参数模型数据。
2 模型试算设计的正交各向异性参数多层模型如图 4所示, 模型横向长度为3000m, 包括300个CDP, 深度为4000m, 采样间隔为10m。模型由浅到深分为4层, 第1层为各向同性地层, vP0=3000m/s, ε1=ε2=δ1=δ2=δ3=0;第2层为TTI各向异性地层, vP0=3250m/s, ε1=ε2=0.1, δ1=δ2=0.05, δ3=0;第3层为正交各向异性地层, vP0=3500m/s, ε1=0.3, ε2=0.035, δ1=0.266, δ2=0.001, δ3=0.22;第4层为各向同性地层, vP0=4000m/s, ε1=ε2=δ1=δ2=δ3=0, 该地层作用是提供反射界面, 不参与反演参数的更新。该模型包含各向同性、TTI各向异性和正交各向异性介质, 可模拟实际地质情况并验证本文方法的适用性。
如图 5所示, 从正交各向异性参数初始模型(vP0=3000m/s, ε1=ε2=δ1=δ2=δ3=0)开始反演迭代更新, 采用本文方法得到的正交各向异性参数精细模型和数值曲线分别如图 6和图 7所示, 可以看出, 背景模型不包含任何梯度信息, 更新难度大。我们采用了递进式的更新思路和针对性的反演策略, 建立的精细模型数值准确、构造界面明显、空间分布合理, 能精确地分辨和描述各向同性、TTI各向异性、正交各向异性介质, 刻画各个介质的各向异性特征, 并验证了本文方法的有效性和适用性。
对某探区实际资料采用本文方法进行处理以验证方法的有效性。该探区目的层深, 地层倾角小, 复杂断裂和缝洞储层发育, 构造复杂, 存在较强的正交各向异性特性, 断裂和缝洞是油气有利聚集区。对探区采用TTI各向异性处理方法得到的成像结果质量较高, 但断裂和缝洞成像精度较低, 不利于储层精细描述和断裂准确预测。
图 8为采用本文方法建立的正交各向异性参数精确模型。模型参数分布合理, 符合断裂和缝洞发育规律, 参数数值准确, 能够描述地震波速度的方位特性。图 9a和图 9b分别采用TTI各向异性处理方法和本文方法得到的偏移成像结果。TTI介质各向异性建模及成像整体质量高, 但是碳酸盐岩目的储层的断裂和缝洞成像精度低, 存在断裂不清晰、缝洞绕射波聚焦性差等现象。经本文方法处理后, 断裂清晰、绕射波收敛, 断裂和缝洞的成像精度进一步得到了提升。图 10a和图 10b分别为对探区采用TTI各向异性处理方法和本文方法提取的OVT域成像道集。可以看出, 对TTI介质采用各向异性处理方法时, 存在方位信息的OVT域成像道集适应性差, 无法将不同方位的道集完全拉平的问题, 说明TTI各向异性处理方法不能准确描述地震波速度的方位差异。经本文方法处理后, OVT域成像道集被完全拉平, 这是因为本文方法可通过多参数准确刻画地震波速度, 再利用速度将弯曲的成像道集拉平。
本文通过正交各向异性介质多参数建模方法的研究及其在实际地震资料中的应用, 得出以下结论:
对正交各向异性介质, 采用基于分方位速度分析的裂缝方位预测和各向异性参数提取方法, 建立了正交各向异性参数初始模型, 再通过建立方位正交各向异性层析矩阵, 并制定合理的多参数反演策略, 最终建立精细的正交各向异性参数模型, 为后续偏移成像和道集提取提供了可靠的模型数据。
理论模型数据和实际资料测试结果表明, 正交各向异性介质多参数建模方法能够提高断裂的成像精度、增强绕射波收敛、改善缝洞体的聚焦性, 同时能够拉平OVT域成像道集, 提升成像效果。虽然很多通过消除道间时差拉平OVT域成像道集的方法可以提高叠加效果, 改善剖面分辨率, 但是对于断裂和缝洞的成像精度改善有限。因为断裂和缝洞对空间速度非常敏感, 仅拉平道集无法从根本上提高成像效果, 所以必须依赖正交各向异性参数准确刻画速度的方向差异, 消除旅行时误差, 提升整体成像效果。
[1] |
刘瑞合, 赵金玉, 印兴耀, 等. VTI介质各向异性参数层析反演策略与应用[J]. 石油地球物理勘探, 2017, 52(3): 484-490. LIU R H, ZHAO J Y, YIN X Y, et al. Strategy of anisotropic parameter tomography inversion in VTI medium[J]. Oil Geophysical Prospecting, 2017, 52(3): 484-490. |
[2] |
傅旦丹, 何樵登. 正交各向异性介质地震弹性波场的伪谱法正演模拟[J]. 石油物探, 2001, 40(3): 8-14. FU D D, HE Q D. Elastic wavefield forward modeling with pseudo-spectral method in orthorhombic anisotropy media[J]. Geophysical Prospecting for Petroleum, 2001, 40(3): 8-14. DOI:10.3969/j.issn.1000-1441.2001.03.002 |
[3] |
TSVANKIN I, THOMSEN L. Non-hyperbolic reflection moveout in anisotropic media[J]. Geophysics, 1994, 59(8): 1290-1304. DOI:10.1190/1.1443686 |
[4] |
TSVANKIN I, THOMSEN L. Inversion of reflection traveltimes for transverse isotropy[J]. Geophysics, 1995, 60(4): 1096-1108. |
[5] |
TSVANKIN I. Seismic signatures and analysis of reflection data in anisotropic media[M]. Pergamon: Elsevier Science Press, 2001: 1-57.
|
[6] |
TSVANKIN I. Reflection moveout and parameter estimation for horizontal transverse isotropy[J]. Geophysics, 1997, 62(2): 614-629. DOI:10.1190/1.1444170 |
[7] |
TSVANKIN I. Anisotropic parameters and P-wave velocity for orthorhombic media[J]. Geophysics, 1997, 62(4): 1292-1309. DOI:10.1190/1.1444231 |
[8] |
郭恺, 杨林. 一种新的TTI介质多参数联合层析反演方法[J]. 石油物探, 2019, 58(3): 412-418. GUO K, YANG L. A new multi-parameter joint tomography inversion method for TTI medium[J]. Geophysical Prospecting for Petroleum, 2019, 58(3): 412-418. DOI:10.3969/j.issn.1000-1441.2019.03.010 |
[9] |
BAKULIN A, GRECHKA A, TSVANKIN I. Estimation of fracture parameters from reflection seismic data-Part Ⅰ:HTI model due to a single fracture set[J]. Geophysics, 2000, 65(6): 1788-1802. DOI:10.1190/1.1444863 |
[10] |
BAKULIN A, GRECHKA A, TSVANKIN I. Estimation of fracture parameters from reflection seismic data-Part Ⅱ:Fracture models with orthorhombic symmetry[J]. Geophysics, 2000, 65(6): 1803-1817. DOI:10.1190/1.1444864 |
[11] |
BAKULIN A, GRECHKA A, TSVANKIN I. Estimation of fracture parameters from reflection seismic data-Part Ⅲ:Fracture models with monoclinic symmetry[J]. Geophysics, 2000, 65(6): 1818-1830. DOI:10.1190/1.1444865 |
[12] |
GRECHKA V, PECH A, TSVANKIN I. Parameter estimation in orthorhombic media using multicomponent wide-azimuth reflection data[J]. Geophysics, 2005, 70(2): 441-452. |