2016, Vol. 59
PP波和PS波联合反演是目前地震储层预测的有利方式，能够提供有用的储层特征信息.自诞生以来，纵横波联合反演在方程参数化、反演算法优化、纵横波资料匹配、实际应用等方面获得长足的进展.Aki和Richards(1980)对Zoeppritz方程简化得到基于纵、横波速度和密度的PP波和PS波AVO近似公式，该简化公式的出现为纵横波联合反演理论的发展奠定了基础.在方程参数化方面，Fatti等(1994)在Aki近似的基础上将纵波和转换波AVO近似方程表示为纵波阻抗、横波阻抗和密度的表达式，利用加权叠加的方法精确地反演纵、横波阻抗；侯栋甲等(2014)提出了基于模量和密度表示的AVO近似公式，采用最小二乘准则构建目标函数以求取相应的参数；张广智等(2014)将纵横波联合反演应用到页岩气储层预测，认为杨氏模量与密度乘积Eρ和泊松比能够增强页岩气储层的异常特征和含油气性，对纵波和转换波公式变换直接反演上述参数来定性指导页岩气储层的勘探.在反演算法优化方面，利用奇异值分解方法从多分量地震资料同时反演不同的弹性参数，反演算法相对稳定且反演 精度高于单纯纵波反演(Mahmoudian and Margrave，2004; Veire and Landrø，2006)；Agullo等(2004)加入测井和地质等先验信息的约束，运用贝叶斯理论从多波地震数据体中反演获取精确的构造特征，在实际应用中取得了很好的效果；常规叠前反演的密度项的反演通常是不稳定的，正演发现转换波反射系数的对密度的敏感程度比纵波反射系数更加强烈，选择全局寻优的算法实现纵横波同步联合反演以提高密度估计的准确性(Leiceaga et al.，2010).纵横波速度比的精确求取是纵波和转换波资料的匹配的基础，Knapp等(2002)利用全角度的PP波偏移道集用反演纵波阻抗，从PS波叠前时间偏移道集求取横波阻抗，用单纯波数据反演的VP/VS提高PP波和PS波数据体同相轴匹配精度；Garotta等(2002)从纵波和转换波截距、梯度中提取最优的纵横波速度比，利用纵波和转换波精确的层位解释结果进行纵波和转换波资料匹配，为纵横波联合反演提供可靠的基础资料.实际多波多分量地震资料在国内外探区的取得了成功的应用，纵波阻抗、横波阻抗、λρ、μρ等弹性参数能够反映储层岩性的变化规律，利用联合算法得到纵、横波速度和密度后求取上述属性体，有效地识别储层特征(Wang et al.，2009，2012; 陈志德等，2011)；Chaveste等(2013)运用PP-PS联合反演预测岩石的脆性特征及裂缝密度，定量地评价TOC含量，将其应用到Marcellus页岩中，降低了甜点预测的风险.尽管纵横波联合反演在上述领域取得了相应的发展，但是，常规的纵横波联合反演只能获取反映储层特征的弹性参数，而缺乏直接有效的流体指示剂进行流体检测.
地震反演技术的应用从勘探阶段深入到开发阶段，储层流体识别逐渐成为研究的热点，近年来取得了快速的发展，Russel等(2003，2011)在Boit-Gassmann 方程的基础上提出Gassmann流体项f，直接反演流体因子、剪切模量和密度等参数，判别储层孔隙流体类型.在深部储层或者缺失大角度信息的情况下，采用两参数叠前反演方法可以直接得到反映储层类型的流体因子，为流体识别提供有力的地球物理支撑(印兴耀等，2014; 李超等，2014).但是，目前的地震流体识别都是基于单纯纵波叠前反演方法，存在精度不高的问题、反演结果不够稳定的特点.
结合纵横波联合反演和储层流体识别优势及存在的问题与不足，本文提出了基于f-μ-ρ的纵横波联合反演方法.首先，对纵波和转换波资料进行同相轴匹配处理，然后从测井、地质等先验信息出发，运用PP和PS波地震资料，建立叠前联合反演流程，直接反演流体因子f、剪切模量和密度等参数. 模型测试和实际资料结果表明方法具有一定的可行性和适用性.
2 方法原理 2.1 转换波反射系数近似公式推导 2.1.1 转换波反射系数公式推导PP和PS波AVO联合反演以Zoeppritz方程为理论依据，但是该方程求解比较复杂，不便应 用，不同的学者对其做线性化近似.Aki和Richards(1980)对其做一定的假设，假定反射界面两侧的弹性参数变化不大，且入射角较小时，近似地将PP波和PS波反射系数表示为纵、横波速度以及密度变化率和入射角的函数关系式，如公式(1)和(2)所示:
|
(1) |
|
(2) |
其中：θ、φ分别是纵波入射角和转换波反射角，VP、VS、ρ分别是纵波速度、横波速度和密度，γsat2=[VP/VS]2是饱和岩石纵横波速度比的平方.
为了直接估计流体特征、骨架信息和密度等参数，Russel等(2011)在Aki近似的基础上，利用孔隙弹性介质理论推导了基于Gassmann流体项f的纵波线性近似方程，表示形式为
|
(3) |
其中：γdry2=[VP/VS]2是干燥岩石的纵横波速度比的平方.
为了研究基于Russell近似的纵横波联合反演方法，首先根据孔隙介质理论推导基于f-μ-ρ的新PS波AVO近似方程.在公式(2)两边同时乘以ρVS2，得到
|
(4) |
Russell等(2003)将Gassmann流体项表示为f=VP2ρ－γVdry2VS2ρ.同时，剪切模量与横波速度和密度的关系为μ=ρVVS2，通过整理可以得到
|
(5) |
|
(6) |
将公式(5)和(6)带入式(4)中可得
|
(7) |
饱和状态下的纵横波速度比的平方表示为γsat2=V2P/V2S，所以，可将剪切模量表示为密度与纵波速度的乘积形式，即μ=ρV2S=ρVV2P/γsat2，化简得到，γsat2Δμ=2ΔVPVPρ+ΔρV2P，将其带入式(7)中，并在方程两边同时除以ρV2S，得到基于Russell近似的转换波AVO反射系数公式
|
(8) |
需要注意，式中θ表示纵波入射角，是弹性界面之上纵波入射角和弹性界面以下纵波透射角的平均值；φ表示转换波反射角，是弹性界面之上转换波反射角和弹性界面以下转换波透射角的平均值.在已知纵横波速度比和纵波入射角θ的情况下利用Snell定律计算转换波反射角φ.
上述公式建立了转换波地震响应特征与储层流体因子、剪切模量和密度的直接关系，将转换波地震资料与流体因子有效的联系起来，通过PS波地震记录进行AVO反演估算流体因子，因此，该方程具有直接反演流体因子的优点；同时该方程相对常规的转换波AVO公式具有较高近似精确的优点(见2.1.2分析)，有利于增强地震反演的稳定性，进而提高储层预测的精度.
2.1.2 近似精度分析建立层状地质模型验证新转换波AVO公式的近似精度，该模型第一层岩性为页岩，第二层岩性为含气砂岩，模型相关参数如表 1所示.对新的AVO近似公式、Zoeppritz精确解、Aki近似进行对比，纵波和转换波AVO曲线如图 1所示，黑色曲线是精确Zoeppritz方程曲线、红色虚线是Aki近似曲线、蓝色圈线是Rusesll近似曲线.曲线对比显示，基于Russell近似的转换波反射系数更接近精确解，入射角小于30°时，新转换波方程的近似精度与精确解的相对误差为2.9%，而Aki近似的转换波方程与精确解的相对误差为3.4%.可见，新近似公式精度在联合反演的要求范围之内.
|
表 1 砂岩模型参数 Table 1 Parameters of sandstone models |
|
图 1 不同近似反射系数曲线对比，(a)纵波反射系数对比，(b)转换波反射系数对比(黑色实线是由Zoeppritz方程所得，红色虚线是Aki近似公式所得，蓝色圈线是Russell近似所得) Fig. 1 Comparison of different approximate reflection coefficients(a)PP wave(b)PS wave(Black solid curves represent exact Zoeppritz equation，Black circle curves represent Aki approximation，Black star curves represent Russell approximation) |
由于贝叶斯反演能结合了地震信息的概率密度特征和待求参数的先验信息，反映所估计参数的统计特性，能够很好地改善反演的“不适定”性问题.所以，本研究根据纵横波联合反演理论，结合PP波和PS波地震数据，在贝叶斯框架下，加入测井解释、先验地质信息的约束，按照图 2所示的流程进行纵横波联合反演方法，获取流体因子f、剪切模量μ和密度ρ等弹性参数的估计值.
|
图 2 基于Russell近似的纵横波联合反演流程图 Fig. 2 Workflow of joint pre-stack inversion based on Russell approximation |
联立方程(3)和方程(8)，不同入射角的纵波和转换波反射系数表示成线性矩阵的关系
|
(9) |
将上式简写为
|
(10) |
给定观测数据 d(包含不同的入射角)和先验地质信息 I 的条件下，参数向量 m 的后验概率分布为
|
(11) |
其中， P(d|m，I)是似然函数， P(m，I)是先验信息.当后验概率分布函数的形状不变的情况下， P(d|I)是一个常数，可以忽略.
如果给定参数向量 m，假设观测数据的第k次测量的噪音是服从方差为σk的高斯分布，那么，第k个测量数据dk的随机分布为
|
(12) |
由于每次测量的结果是相对独立的，且假定高 斯分布的方差是均匀分布的且σm=σ，根据乘积准则可以得到似然函数为
|
(13) |
式中：σ2是噪音方差，M是参数个数，N是每个参数的采样点.对于叠前AVO多参数同时反演，似然函数是一个多变量的高斯分布函数，它的等概率解是一个椭球，所以，不同的参数具有不同的反演不确定性(据陈建江，2007).
假定参数的扰动量服从高斯分布，有
|
(14) |
式中： C m 是待反演参数之间的协方差矩阵，
|
最后，得到参数的估计值
|
(15) |
以实际测井数据建立纵波和转换波合成角度道集.用式(3)和式(8)计算出1°到30°纵波和转换波反射系数，与40 Hz的雷克子波褶积，得到纵波和转换波合成角度道集，并在地震道集中加入随机噪音，合成角度道集如图 3所示，图 4是按照图 2所示的流程进行纵横波联合反演得到的结果，其流体因子f、剪切模量μ和密度ρ的估测值(红色曲线)与 真实值(黑色曲线)非常接近，尤其在目的层段(2200 ms左右)反演结果与测井测量值完全吻合，并将其与单纯的纵波反演(蓝色曲线)比较，结果显示，联合反演得精度更高，更加接近真实值，尤其是密度的反演结果，从而验证了方法取得较好的效果.
|
图 3 SNR=3，正演得到的角度道集(a)纵波角度道集，(b)转换波角度道集 Fig. 3 Pre-stack synthetic seismogram(a)PP wave(b)PS wave with SNR=3 |
|
图 4 SNR=3时，叠前纵横波联合反演所得的流体因子f(N·m-2)、剪切模量(N·m-2)和密度(kg·m-3)的反演结果(红色曲线)与单纯纵波的反演结果(蓝色曲线)以及原始值(黑色曲线)的对作对比比 Fig. 4 Comparison of joint inversion and single PP wave inversion of fluid factor f，shear modulus and density，SNR=3 |
为进一步验证方法的实用性，笔者选取中国西部某实际工区的多波多分量地震资料进行方法试验，首先，对纵波和转换波资料进行匹配处理.由于纵波速度和横波速度不同，所以，地下同一地质体的PP波和PS波的双程旅行时不同，为此要对PP波和PS波地震数据进行同相轴匹配处理.
在零偏移距的情况下，纵波入射，纵波反射，有PP波的双程旅行时表示为
|
(16) |
式中：H是单个地质反射体的厚度.
同样，在零偏移距情况下，假设P波入射，S波反射，可以得到PS波传播时间表示为
|
(17) |
所以，纵波和转换横波的传播时间比值为
|
(18) |
从上式可以看出，纵波和转换波资料匹配的关键在于求取准确的纵横波速度比γ，在匹配前根据井震标定的结果，在纵波和转换波资料上分别解释同一反射体的PP波和PS波旅行时，同时求取精确的纵横波速度比剖面，利用上式将PS波地震数据转换到PP波时间中，可以实现纵波和转换波资料的匹配，图 3是纵波和转换波资料的匹配结果，结果表明纵波和转换波资料匹配程度较高，极性相反.
3.2.2 新的纵横波联合反演匹配前后纵波和转换波资料如图 5、图 6所示，图 6a是PP时间的纵波叠加剖面，该工区的目的层段在3.50~3.54 s(PP时间)，在测井解释的结果是油层，图 6b是匹配前的转换波叠加剖面(PS时间)，图 6c是匹配后的转换波叠加剖面(PP时间)，对比发现，纵波资料与转换波资料匹配程度较高，能够反映同一反射体的地质信息，这为后来的纵横波联合反演提供了数据基础.图 7是通过联合反演得到的参数估计结果，图 7a是流体因子的反演结果，图 7b是剪切模μ的反演结果，图 7c是密度ρ的反演结果.图 7(a，b，c)显示的测井曲线分别是流体因子f、剪切模量μ和密度ρ曲线，井点处测井解释的油层在3.52~3.53 s之间，反演获取的三个弹性参数反演的结果(红色值)与测井解释结果基本符合，能很好地指示该工区的流体特征，进一步验证了新型的联合反演算法的有效性与实用性.
|
图 5 纵波和转换波资料匹配结果 Fig. 5 Result of PP wave and PS wave event matching |
|
图 6 纵波和转换波地震剖面(a)匹配前纵波地震剖面，(b)匹配前转换波地震剖面，(b)匹配转换波地震剖面 Fig. 6 PP and PS wave section，(a)PP wave section before event matching，(b)PS wave section before event matching(b)PS wave section after event matching |
|
图 7 叠前纵横波联合反演所得的流体因子(N·m-2)、剪切模量(N·m-2)和密度(kg·m-3)的反演结果 Fig. 7 Inverted fluid factor，Lame parameter and density |
纵横波联合反演作为地震储层预测和流体识别的有效手段，为油气田的勘探开发提供有力的地球物理依据，本文推导了以流体因子、剪切模量和密度为参量的PS波反射系数公式，页岩和含气砂岩模型之间弹性界面的反射证明了新推导的转换波反射系数的精度能够满足叠前纵横波联合反演的要求.在此基础上，建立能够直接面向流体指示因子的联合反演流程，选择实际井模型资料进行方法测试，测试结果显示，新反演方法获得参数估测值与真实值匹配程度较高，且较单纯纵波反演的精度更高.最后，选取实际的多波多分量叠前地震资料进行试验，通过纵波和转换波资料的匹配进行联合反演，估测结果表明流体因子、剪切模量和密度的反演结果与测井解释一致，验证了基于Russell近似的纵横波联合反演方法可行性，可为叠前纵横波联合反演提供了一种新的方法思路.
利用PP波和PS波联合反演估算孔隙度、饱和度、裂缝密度等储层特征参数，可以提高储层预测的可靠性，多波联合反演方法在致密储层和裂缝储 层等非常规油气藏勘探中应用是其下一步发展方向.
致谢感谢审稿专家和编辑部的指导和帮助.
| Agullo Y, Macè D, Labat K, et al. 2004. Joint PP and PS stratigraphic inversion for prestack time migrated multicomponent data.//2004 SEG Annual Meeting. Denver, Colorado, 889-892. | |
| Aki K, Richards P G. 1980. Quantitative Seismology:Theory and Methods. San Francisco:W. H. Freeman and Co. | |
| Chaveste A, Zhao Z Y, Altan S, et al. 2013. Robust rock properties through PP-PS processing and interpretation-Marcellus Shale. The Leading Edge , 32(1): 86–92. | |
| Chen J J. 2007. Study of Three-term AVO Inversion Method[Ph. D. thesis] (in Chinese). Qingdao:China University of Petroleum (East China). | |
| Chen Z D, Chen S M, Wang J M, et al. 2011. Research to the application of 3D3C seismic data in Lamadian oilfield. Scientia Sinica Terrae (in Chinese) , 41(7): 927–935. | |
| Fatti J L, Smith G C, Vail P V, et al. 1994. Detection of gas in sandstone reservoirs using AVO analysis:A 3-D seismic case history using the Geostack technique. Geophysics , 59(9): 1362–1376. | |
| Garotta R, Granger P Y, Dariu H. 2002. Combined interpretation of PP and PS data provides direct access to elastic rock properties. The Leading Edge , 21(6): 532–535. | |
| Hou D J, Liu Y, Hu G Q, et al. 2014. Prestack multiwave joint inversion for elastic moduli based on Bayesian theory. Chinese Journal of Geophysics , 57(4): 1251–1264. doi: 10.6038/cjg20140422. | |
| Knapp S, Data S, Medley T, et al. 2002. Inversion of PP and PS multicomponent seismic volumes:a case history from the Gulf of Mexico.//2002 SEG Annual Meeting. Salt Lake City, Utah. | |
| Leiceaga G G, Silva J, Artola F, et al. 2010. Enhanced density estimation from prestack inversion of multicomponent seismic data. The Leading Edge , 29(10): 1220–1226. | |
| Li C, Yin X Y, Zhang G Z, et al. 2014. Two-term elastic impedance inversion based on the incident-angle approximation. Chinese Journal of Geophysics , 57(10): 3442–3452. doi: 10.6038/cjg20141030. | |
| Mahmoudian F, Margrave G F. 2004. Three parameter AVO inversion with PP and PS data using offset binning.//2004 SEG Annual Meeting. Denver, Colorado, 240-243. | |
| Russell B H, Hedlin K, Hilterman F J, et al. 2003. Fluid-property discrimination with AVO:A Biot-Gassmann perspective. Geophysics , 68(1): 29–39. | |
| Russell B H, Gray D, Hampson D P. 2011. Linearized AVO and poroelasticity. Geophysics , 76(3): C19–C29. | |
| Veire H H, Landrø M. 2006. Simultaneous inversion of PP and PS seismic data. Geophysics , 71(3): R1–R10. | |
| Wang D X, Zhao Y H, Zhang M B, et al. 2012. Characterizing tight sand gas reservoirs using 3D3C seismic data in Sulige gas field.//2012 SEG Annual Meeting. Las Vegas, Nevada, 1-5. | |
| Wang J S, Huo L, Tian C C, et al. 2009. Lithology and fluid prediction using multiwave seismic data.//2009 SEG Annual Meeting. Houston, Texas, 1237-1241. | |
| Yin X Y, Cao D P, Wang B L, et al. 2014. Research progress of fluid discrimination with prestack seismic inversion. Oil Geophysics Prospecting (in Chinese) , 49(1): 22–34. | |
| Zhang G Z, Du B Y, Li H S, et al. 2014. The method of joint pre-stack inversion of PP and P-SV waves in shale gas reservoirs. Chinese Journal of Geophysics , 57(12): 4141–4149. doi: 10.6038/cjg20141225. | |
| 陈建江. 2007. AVO三参数反演方法研究[博士论文]. 青岛:中国石油大学(华东). | |
| 陈志德, 陈树民, 王建民, 等. 2011. 喇嘛甸油田3D3C地震资料的应用探索. 中国科学:地球科学 , 41(7): 927–935. | |
| 侯栋甲, 刘洋, 胡国庆, 等. 2014. 基于贝叶斯理论的叠前多波联合反演弹性模量方法. 地球物理学报 , 57(4): 1251–1264. | |
| 李超, 印兴耀, 张广智, 等. 2014. 基于入射角的两项流体阻抗反演方法. 地球物理学报 , 57(10): 3442–3452. | |
| 印兴耀, 曹丹平, 王保丽, 等. 2014. 基于叠前地震反演的流体识别方法研究进展. 石油地球物理勘探 , 49(1): 22–34. | |
| 张广智, 杜炳毅, 李海山, 等. 2014. 页岩气储层纵横波叠前联合反演方法. 地球物理学报 , 57(12): 4141–4149. | |