井间地震时频域联合速度和衰减系数同步层析反演

刘浩杰， 王延光， 魏国华， 尹文笋， 张建中， 郑静静， 陈雨茂

(1. 中国石化 胜利油田分公司物探研究院, 东营 257022;2. 中国海洋大学 海洋地球科学学院，青岛 266100)

0 引言

井间地震是将震源和检波器都置入井中进行探测的一种地球物理方法。由于避开了地表低速带对地震信号高频成分的吸收，可以获得较高分辨率的地震信号，能够进行井间构造、储层物性、储层连通性、注气效果等复杂油气藏问题的评价。

对井间地震资料的处理涉及初至拾取、波场分离、速度层析反演、反射波成像等关键环节[1]。其中井间地震速度层析成像也被称为井间地震CT技术[2]，它能够利用井间地震直达波的时域旅行时信息反演两井间较为精细的速度剖面，为井间地震反射波成像提供了必要的速度模型。理论上，如果井间地震速度层析成像足够精确，则可以直接作为井间岩性、构造和储层等综合解释的成果资料来解决井间构造、储层物性、储层连通性等方面的地质问题。由于井间地震观测角度的局限性、射线分布的不均匀性、波场特征的复杂性和反演问题的多解性等原因[3]，高精度的井间地震速度层析成像是较为困难的。同时，通过对井间地震直达波时域振幅信息和频域频率信息的处理，可以实现两井间的吸收衰减层析成像[4]。由于地下介质的非完全弹性性质，地震波能量在地下介质中传播会受到介质吸收作用的影响。随着地震波的传播，其能量不断衰减，振幅不断减小，高频信息逐渐减少。相对于地震波传播速度，地震波的衰减更能反映岩石物性及流体特征，如饱和度、孔隙度、渗透性、粘性等[4-5]。通过井间地震衰减层析成像，可以反演两井间储层衰减特征的分布情况。

国内、外学者在衰减层析成像方面已进行了较为广泛地研究， Quan等[6]提出了一种计算地震波初至频率变化的算法(质心频移法 Centroid Frequency Shift, CFS) , 通过分析信号的频率变化达到研究岩石的地震波衰减特性的目的，但仅用一维零偏VSP资料进行衰减系数估算，而且求取的衰减量易受噪声干扰等非地层固有衰减因素影响；严又生等[4]、Liao等[5]提出了利用井间地震走时和振幅同时反演速度和衰减的层析成像思路，较好改善了衰减层析成像效果，但在具体方法的系统实现以及在储层流体检测实践应用方面略显不足。为了进一步提高井间衰减层析成像的精度和分辨率，提升衰减成像在含流体性识别方面的应用效果，并降低其对井间地震资料噪声干扰和非地层固有衰减的因素影响，笔者基于胜利油田实际井间地震资料的数据，在井间地震直达波走时精确计算基础上，基于射线路径的一致性，研发了时域速度层析和频域衰减层析反演方法，从井间地震资料中同步反演速度和衰减信息。该方法既充分利用了速度层析对地层结构的反演精度和分辨率高的特点，又联合了衰减层析对储层物性、流体性质响应更为敏感的优势，实现了构造与储层流体性质的同步反演，提高了井间地震成果剖面解决地质和油藏开发问题的能力，扩展了井间地震的应用范围和效益。

1 地震射线路径追踪

井间地震直达波层析反演一般需要四个关键的步骤[1-2]：①旅行时拾取；②反演初始模型建立；③地震射线追踪；④速度层析反演。其中前三步是速度层析反演的基础，旅行时拾取可获取直达波的走时信息，反演初始模型是为了建立网格化的离散速度模型作为反演迭代的初始条件，射线追踪则是为了获得直达波射线路径和理论最小走时。

在对反演初始模型离散化时，由于地层速度线性变化，一个离散单元内上、下速度界面需要用斜面近似。因此，相比于传统的矩形规则网格模型离散方法，根据地层先验信息进行不规则网格剖分的离散化方法则更接近实际地层模型[7-8]。首先，在水平面上进行等间距地质模型剖分，模型将被离散成规则的矩形单元；其次，根据测井和地质综合解释得到的层位信息，在垂向上划分地层，把模型离散成规则分布的六面体单元，即每个单元左右前后四个侧面均平行于垂向轴，顶面和底面为斜面。

通过对井间模型不规则网格的离散化，就可以较为准确地追踪地震射线路径，得到任一从震源点出发到井间模型所有网格单元的波前走时。我们采用基于Fermat原理的三维射线追踪方法，进行井间地震任一震源点到接收点的射线路径追踪[8]。本文方法能够适应斜井状况，更符合地震波在地层中的传播特征，能够较为高效、准确地追踪到地震波到达接收点的全局最小旅行时路径。井间地震射线路径地确定，也就得到了不同接收点射线走时变化和衰减变化在两井间所经历的路径，这就为速度和衰减的同时层析反演提供了借鉴。

2 井间地震速度层析反演方法

井间地震速度层析成像是井间地震资料处理的关键技术，不同专家基于不同的假设条件提出了一些特色的层析成像方法(如直达波最小旅行时反演、初至时间和振幅信息联合反演)、直达波和反射波联合反演等[2-3,7]，其中最常用和最主要的是直达波最小旅行时层析反演。

根据地震射线理论，旅行时t与慢度s有如下关系

(1)

其中:t为直达波初至时间；L(s)为直达波射线路径；S(x,y,z)为三维空间离散地质模型的慢度，dl为直达波射线在每个网格单元的路径。式(1)中，只有直达波初至时间已知，三维空间的慢度信息和直达波射线路径均需要求解。对该非线性关系表达式进行网格离散化，可以得到如下表达式：

(2)

式中：sj为离散模型第个网格单元的慢度；lij表示第条射线在第个网格单元的长度；m为射线条数；n为离散模型网格单元个数。通过对条射线路径方程的联立求解，根据特定的误差准则，即可利用最优化反演方法反演得到两井间速度剖面。

3 井间地震衰减层析反演方法

地震信号的衰减特征越来越多应用于储层的描述和油气预测，出现了各具特色的地震衰减属性计算方法(如振幅衰减法、上升时间法、频谱比法和质心频移法等[4,9])，但不同方法有各自的适应条件和局限性，需要根据实际资料的情况进行方法的选择和优化。在时域中，振幅衰减法主要考虑了地层对地震波的吸收衰减作用导致地震振幅的变化，但时间域的振幅信息容易受几何扩散、仪器响应特征、反射和透射作用等因素的影响，难以得到准确的地震信号振幅信息。在频率域中，仪器响应特征、几何扩散等因素对地震信号频率变化特征影响有限，能够得到较为准确的地层对地震波吸收衰减特征。质心频移法是利用震源子波和接收地震信号频谱的质心频率变化量，能够相对方便可靠的计算衰减特征[4]。

定义震源地震子波的振幅谱为s(f)，其质心频率f(s)和方差σs分别为

(3)

同理，定义接收地震波振幅谱为R(f)，则其质心频率fR和方差σR分别为

(4)

假定地层吸收衰减系数为α，根据质心频率法原理，由于地层对地震波的吸收衰减作用，地震波的质心频率从f(s)减小到f(R)。因此，根据震源子波和接收地震波频谱质心频率的变化就可以估计沿该射线路径L的地层平均吸收衰减系数。因此，可以建立地层吸收衰减系数、地震波射线路径与地震波质心频率变化量之间如下关系:

(5)

由于井间地震的直达波射线有几千条以上，需要得到质心频率法的二维离散化公式。按照前述的模型离散方法，把井间介质离散成小单元网格，可得如下离散化公式：

(6)

同速度层析反演公式一样，式中αj为离散模型第j个网格单元的吸收衰减系数；lij表示第i条射线在第j个网格单元中的长度；m为直达波旅行时个数，即射线条数；n为离散模型网格单元个数。

式(6)中，接收地震波质心频率fR已知，而震源地震质心频率fs和方差σs均未知。因此，取震源点对应的所有接收点地震频谱方差的平均值为震源点频率方差σs，震源子波质心频率为该炮点对应的所有接收点质心频率的最大值，则震源子波质心频率有如下表达式：

fs= max(fR) + Δf

(7)

其中：max(fR)接收点地震信号质心频率最大值；Δf为需要确定的量。将式(7)代入式(6)，则可以得式(8)。

(8)

同样，通过对m条射线路径方程的联立求解，根据特定的误差准则，即可反演得到两井间吸收衰减剖面。

4 速度和衰减系数同步层析方法原理

如上所述，利用井间地震直达波初至信息，通过层析反演可以得到两井之间的速度剖面。通过对井间地震直达波频域信息的处理，可以实现两井之间的吸收衰减层析成像。速度层析是在时域利用直达波的初至，建立矩阵最优化求解方程组。质心频移衰减层析法则是通过计算接收地震波在频率上的变化量来计算衰减系数，计算精度较高，但求取的衰减量易受噪声干扰影响而产生假的衰减异常。因此可以把时域的速度层析与频域的衰减层析两种方法联合起来，同时反演速度和衰减系数，既能提高效率，同时又相互增加了约束条件，提高了反演结果的可靠性，得到更加准确的速度和地层吸收衰减变化。

井间地震时域速度层析公式(2)，如果有m条射线，则用矩阵表示为式(9)。

AvS=T

(9)

其中，

S是由离散模型网格单元慢度所组成的n维矢量，sj为第j个网格单元上的慢度，n为单元总数；T为由接收直达波初至走时组成的m维矢量，ti为第个i直达波走时数据，m为观测数据的总道数；而矩阵Av则是m×n系数矩阵，其元素为直达波射线在该网格单元中的长度，lij表示第j条射线在第j个单元中的长度。

井间地震速度层析反演问题是，已知直达波走时求取射线路径的和每个网格节点的慢度，即已知T求S和Av。式(9)中慢度矢量s和射线路径矩阵Av均未知，而慢度又是直达波射线路径的非线性函数，因此，式(9)实际上非线性方程组。这种非线性反演问题实际求解时，需要采用逐次迭代的方法最优化求解。即先给定一个初始模型s0，用该初始模型计算射线路径和理论旅行时，根据式(9)求出慢度的扰动量δs，修正初始模型s=s0+δs，得到新的慢度模型。如此反复进行，直到计算的理论旅行时与拾取的初至旅行时之差满足给定的条件为止，这时得到的模型便作为最终反演结果。因此在求解过程中，如果有更多的确定性约束条件，则层析反演得到的速度就更加准确。

对于井间地震衰减层析公式(8)，若共有m条射线，则可形成一个线性方程矩阵：

AqΧ=B

(10)

其中，

利用最小二乘原理，通过将式(9)、式(10)联立求解，可得如下方程：

(11)

联立的反演方程组反演方程是超定的，数据量一般都很大，所以反演算法需要能节省内存，和较快的计算效率。这里同样采用阻尼最小二乘法进行求解[10-11]。

具体求解过程如下：①建立初始速度模型和吸收衰减模型；②计算射线路径，并计算理论旅行时和理论质心频率衰减值；③计算理论和实际旅行时之差和质心频率衰减差；④通过修改初始模型，通过一系列的优化迭代；⑤同时反演得到井间地震的速度层析剖面和吸收衰减层析剖面。通过两者的同时反演，减小噪音因素对反演结果的影响，提高了反演的准确度。

5 方法效果分析

5.1 正演模型的速度与衰减联合层析试验

利用典型薄互储层地质模型的井间地震粘弹正演模拟结果，对时频域联合速度和衰减同步层析方法的效果进行了分析。图1(a)为该模型对应的油藏地质剖面，从图中可以看出薄油层尖灭、油水过渡层、油水同层等，储层厚度在5 m-15 m左右，包含丰富的地质信息。根据实际地层速度填充并网格化，可得到纵波速度模型(图1(b))，依据岩石物理分析结果，对油层、水层、含气层分别填充不同的衰减值，衰减系数分别为0.02、0.01、0.03，可得到纵波衰减系数模型(图1(c))。采用粘弹介质波动方程正演方法进行模拟，正演主要观测系统参数为：两井间距离400 m，左井为炮点阵列，右井为接收点阵列；炮点检波点均为168个，炮点间距和检波器间距均为3 m，均匀分布在左右井中；正演子波为雷克子波，激发主频200 Hz，采样间隔0.5 ms，记录长度1 s。图2为第1炮各接收点的初至波及其频谱，随着传播距离增大，初至波的振幅和频率皆逐渐减小，表明了同时利用直达波初至和频率变化信息进行联合反演的可行性。按照时频域联合速度和衰减同步层析方法流程对正演数据进行了处理，得到了层析反演速度模型和衰减系数模型。由图3可以看出，速度和衰减反演的结果与模型剖面均基本一致，而且仅用质心频率反演的衰减系数(图3(c))，其模型上部出现了局部较强的衰减异常值，与真实模型有偏差。其原因是上部地层衰减系数值较小，质心频率移动量较小，导致质心频率插值计算的较大误差，而用联合反演方法则可消除顶部假衰减异常，而且还能将层间砂体尖灭引起的衰减异常精确反演出来，具有更高的分辨率。

图1 正演地质及其地球物理参数模型Fig.1 Geological model and its parameters model(a)薄互层油藏地质模型；(b)纵波速度模型；(c)纵波衰减系数模型

图2 第一炮各接收点的初至波及其频谱Fig.2 Direct wave and its amplitude spectrum for the first gather(a)初至波；(b)频谱

图3 不同方法衰减层析反演比较Fig.3 Comparison between different attenuation coefficient tomography(a) 走时层析反演速度；(b)质心频率法反演的衰减系数剖面；(c)联合反演法反演的衰减系数剖面

5.2 实际数据联合层析

对胜利油田某井间地震实际采集数据进行了速度和衰减的同步层析反演。该实际资料共有43 568道采样数据，道采样点数为2 000，道采样间隔为0.5 ms，炮检距为3 m，道间距为3 m，共有炮点329个，检波点300个。

图4 第49炮的地震记录与提取的初至波波形Fig.4 Seismic record and its direct wave and amplitude spectrum for the 49 gather(a)地震记录；(b)初至波波形；(c)初至波对应的频谱

图5 速度和衰减联合层析反演和油藏剖面比较Fig.5 Simultaneous inversion of velocity and attenuation coefficient and the associate reservoir secion(a)井间油藏地质剖面；(b)速度层析反演剖面；(c)衰减层析剖面

图6 两种方法计算的衰减系数与声波测井速度、油藏描述的对比Fig.6 Comparison of attenuation tomography for different method and its reservois section

图4为第49炮的地震记录、初至波波形及其对应的频谱，从图4(c)可以看出，随着偏移距的增大，初至波振幅与主频逐渐降低。图5为联合层析反演得到的速度和衰减系数剖面与过井实际油藏地质剖面的对比。从图5中可以看出，反演得到的速度剖面和衰减系数剖面与油藏地质剖面一致，特别是衰减特征能够反映较好地区分气层及油层特征。图6为利用VSP计算的衰减系数、时频域联合层析得到的衰减系数、声波测井速度和油藏描述的对比图。由图6可见，井间地震资料反演的衰减系数剖面具有较高的分辨率，与VSP计算的井旁衰减系数曲线趋势较为一致，但分辨率明显更高，可识别出厚度为5 m～8 m的薄气层。图6中通过与油藏剖面的综合对比可以看出，联合反演的衰减特征与储层的含流体性质有较好的对应关系，表明了方法的可行性和可靠性。特别是，1 670 m～1 730 m的多个含气层衰减幅度明显，与水层差异较大，其衰减系数达到0.025，为水层的2倍～3倍。总体而言，不同含流体性质的储层衰减特征表现为：气层强衰减，水层弱衰减，油层介于两者之间，与水层相比衰减特征更明显。将衰减系数反演结果结合油藏岩性剖面进行联合解释，可进行储层的流体识别和含油气性检测，获得了较好的地质效果。

6 结论

基于井间地震直达波射线路径的一致性，本研究提出了井间地震时域速度层析和频域衰减层析联合反演的方法，在模型和实际资料处理中取得了较好效果，可以得到以下认识。

1)井间地震资料中蕴含了丰富的地下储层岩性和物性的信息，是高精度油气描述和预测重要的基础资料。

2)井间地震时域速度层析和频域衰减层析的联合反演，兼顾了时域直达波初至信息和频域振幅衰减信息，同时利用井间地震资料的运动学和动力学信息，增加了约束条件，增强了抗干扰能力，能够提取到更为准确和丰富的储层岩性和物性信息。

3)相比于仅利用直达波和反射波的走时信息来获取速度反演结果的方法，将地震资料的动力学特征中包含的储层和流体信息通过联合反演的方式获取新的井间地震成果剖面，提高了井间地震成果解决地质和油藏开发问题的能力。

4)多尺度资料、多域信息的联合是今后油气精细描述和预测重要方向。