近年来, 随着中国油气勘探开发进入中—高勘探阶段, 油气勘探技术逐渐成熟, 油气勘探逐步进入隐蔽油气藏勘探阶段, 滩坝砂岩储层逐渐引起人们的关注^[1]。国内外许多学者对滩坝的地质特征和识别预测技术进行了大量研究, 滩坝砂具有单层砂体厚度薄、横向变化快、地震资料纵向分辨率低等特点。朱筱敏等^[2]、陈世悦等^[3]根据滩坝的主要物质组成将其分为砂质滩坝和生物碎屑滩坝, 杨勇强等^[4]建立了陆相断陷湖盆的滩坝沉积模式, 刘书会等^[5]利用地震属性对滩坝薄砂体进行了预测, 才巨宏^[6]、韩宏伟^[7]等利用波形分析技术和地震特征反演技术对滩坝储层进行了预测, 黄捍东等^[8]通过非线性随机反演方法对陆相薄砂岩储层进行了预测。这些研究成果均较好地指导了滩坝勘探, 但仍然存在较多问题:地震属性多解性较强; 波阻抗反演分辨率太低, 无法满足薄层需求; 测井约束反演高频分量主要来源于初始阻抗模型, 多解性强; 地质统计学反演随机性较强, 对井位分布要求较高。笔者基于相控模型的高精度反演, 在层序地层格架约束下, 运用沉积学基本原理, 充分利用地震波形的横向变化反映储层空间的相变规律, 实现滩坝储层和沉积相带的预测。

在地震勘探中, 由声波计算得到的波阻抗曲线通常表示为

式中, Y为由声波计算得到的波阻抗值; X为待求解的地下地层实际波阻抗值; N为随机噪音, 一般情况下, 假设N为服从高斯分布的平衡白噪音, 其数学期望为0, 协方差为σ。

由于${\left\| {\mathit{\boldsymbol{Y}}-\mathit{\boldsymbol{X}}} \right\|^2} $与N一样均符合高斯分布, 因此可以确定其目标函数为

波阻抗反演过程就是通过公式(2)研究X为何值时使J₁达到最小。由于地震资料带限宽度的限制, 只能提供波阻抗反演的中频成分, 而高频和低频成分则主要通过测井数据、地质信息获得^[9]。波阻抗反演存在不稳定性和多解性, 仅使用一种方法很难获得较理想解, 因此须在目标函数中引入先验信息对其进行约束, 从而获得较稳定的反演解。在统计方法中, 通常采用最大后验估计, 则目标函数为

即可表示为测量数据估算可信度项与先验项的和。其中Z为波阻抗, 先验项J₂(Z)可以定义为势函数的和, 是与测井资料、地质信息等先验信息有关的函数。J₁(Z)为与某些后验信息有关的函数。λ为平滑参数, 用于协调J₁(Z)和J₂(Z)之间相互影响, 为常数。公式(3)进一步改写为

式中, φ为势函数; λ和δ为正则参数, λ用于平衡数据项和先验项之间的相互影响, δ是一种刻度参数, 调节不连续处梯度值^[10]。公式(4)由两部分组成, 第一部分是关于测量模型和数据之间的某种一致性, 当其值达到最小时, 一致性会变得很好, 但解稳定性变差; 第二部分是关于解“光滑性”的测量标准, 当其值取最小时, 获得的近似解较稳定。

马尔科夫链蒙特卡罗(MCMC)方法是在贝叶斯框架下, 用已有资料进行约束, 既可使最优解满足参数的统计特性, 同时融入了先验信息, 解的精度得到提高^[11]。其基本思想是:通过构建一个平稳分布为π(x)的马尔科夫链来获得π(x)的样本, 被估计的模型参数值是链的状态空间, 被估计的后验分布为链的贝叶斯后验分布。

利用MCMC方法对贝叶斯理论下的后验概率分布进行随机抽样, 反演问题希望得到的是所估参数的后验概率分布, 因此各个参数的马尔科夫链应收敛于所估计参数的后验概率分布^[11]。

设M为某一空间, n为产生的总样本数, m为马尔科夫链达到平稳时的样本数, 则MCMC的思路为^[12]:

(1) 构建一条Markov链, 使其收敛至平稳分布π(x), 常用的构造转移方法是Gibbs抽样和Metropolis-Hastings抽样。

(2) 产生样本。由M中的某一点x⁽⁰⁾出发, 用(1)中的Markov链进行抽样模拟, 产生点序列:x₁, …, x_n。

(3) 蒙特卡罗积分。任一函数f(x)的期望估计为

Metropolis-Hastings抽样的基本思路^[12]是:选择一转移函数q(x; x^(i-1))和初始值x⁽⁰⁾, 若第i次迭代开始时的参数值为x^(i-1), 则第i次迭代过程为

(1) 从q(x; x^(i-1))中抽取一个备选值x′;

(2) 计算接受概率:$\alpha \left( {{x^{\left( {i-1} \right)}}, x'} \right) = \min \left\{ {\frac{{\pi \left( {{x^{\left( {i-1} \right)}};x'} \right)}}{{\pi \left( {{x^{\left( {i-1} \right)}}q\left( {x';{x^{\left( {i - 1} \right)}}} \right)} \right)}}} \right\} $;

(3) 以概率α(x^(i-1), x′), 置x⁽ⁱ⁾=x′, 以概率1-α(x^(i-1), x′), 置x⁽ⁱ⁾=x^(i-1);

(4) 重复(1)~(3)步n次, 则可以得到后验样本x⁽¹⁾, x⁽²⁾, x⁽³⁾…x⁽ⁿ⁾, 从而可计算后验分布的各阶矩以进行统计推理。

波形指示反演技术流程如图 1所示。首先对测井资料处理, 包括测井曲线标准化、重构和敏感曲线分析, 建立层序地层格架, 沉积特征及沉积相特征分析; 然后分析测井数据和地震波形特征, 参照样本空间分布距离和曲线分布特征, 优选与待判别道波形特征相似程度较高的样本井建立初始模型, 统计其纵波阻抗建立先验概率函数; 将初始模型与地震波阻抗体进行匹配滤波, 求得似然函数, 对样本进行多尺度分解, 逐步滤除高频成分; 最后基于贝叶斯理论, 联合似然分布与先验分布得到后验概率分布, 并将其作为目标函数; 基于波形指示优选的样本在空间上具有较好的相关性, 采用Metropolis-Hastings抽样算法对后验概率分布抽样, 选取目标函数最大值时的解作为可行随机实现, 求取多次可行实现的均值作为期望输出。

图 1

图 1 波形指示反演流程 Fig.1 Flow chart of waveform indication inversion

东营凹陷王家岗地区沙四上亚段主要发育滩坝相沉积, 滩坝砂岩油藏是其中一种重要的岩性油藏类型^[13], 分布面积广、认识程度较低, 但勘探潜力巨大。实际勘探开发中往往会遇到一系列困难:岩性组合为砂泥岩薄互层, 单层砂体多小于2 m, 储层物性差, 横向连续性差, 分布规律难以掌握; 同时地震资料主频低, 目的层段约为30 Hz, 有效频带范围10~70 Hz, 目的层速度约为2 000 m/s, 地震资料对储层的分辨能力偏低, 对滩坝砂岩储层地震相特征认识不清, 导致滩坝砂体识别和预测难度大^[14-16]。

渤海湾盆地是中国东部最重要的含油气盆地, 济阳坳陷位于渤海湾盆地的东南部, 其中东营凹陷位于济阳坳陷的东南部, 是渤海湾盆地中一个重要的三级构造单元。王家岗地区位于东营凹陷南斜坡的东段, 其西北部为牛庄洼陷, 东北部为广利洼陷, 东部为八面河断裂鼻状构造带, 南部为广饶凸起。构造位置处于近东西向陈官庄-王家港断裂带与东北向八面河断裂鼻状构造带交汇处, 西南部为近北西向纯化草桥断裂鼻状构造带(图 2)。

图 2

图 2 东营凹陷区域构造 Fig.2 Structure units in Wangjiagang area, Dongying Depression

建立层序地层格架为开展储层预测和沉积相研究提供了地震地质的等时格架, 在地震沉积相研究中不可缺少。王家岗地区沙四上亚段为一个完整的三级层序^[17], 应用初次湖泛面和最大湖泛面将其划分为低水位体系域、湖侵体系域和高水位体系域三个体系域^[18]。沙四上亚段顶界面在地震剖面上对应于T6反射层, 反射特征为强振幅、强连续反射, 全区内可以较大范围追踪, 为沙四段与沙三段的分界面; 沙四上亚段层序内部以沉积作用转换面T7为界, 为低水位体系域向湖侵体系域过渡面, 地震剖面上中强振幅、中连续反射特征; T7′为沙四上亚段和沙四下亚段的分界面。湖泊低水位体系域对应准层序组1、准层序组2和准层序组3, 准层序组1湖泊水体较浅, 准层序组2时期水体开始逐渐变浅至准层序组3时最浅, 主要发育砂泥间互的滩坝相沉积; 湖侵体系域对应准层序组4, 湖泊水体开始快速上升, 发育一套分布稳定高导泥岩; 高水位体系域对应于准层序5、准层序组6和准层序组7, 早期湖泊水体深且湖平面稳定, 沉积厚层泥岩, 晚期水位逐渐下降, 有一定量的沉积物供应, 形成一定分布范围的滩坝相沉积(图 3)。

图 3

图 3 东营凹陷南北向层序地层剖面 Fig.3 North-south profile of sequence stratigraphic in Wangjiagang area, Dongying Depression

由于测井仪器误差等因素导致测井数据出现偏差, 采用频率直方图法进行测井资料的标准化处理; 统计分析寻找对岩性变化比较敏感的测井曲线, 通过交会分析认为波阻抗曲线对砂泥岩区分较差, 砂岩和泥岩部分阻抗值重叠(图 4(a)), 无法用现有波阻抗曲线进行反演, 但自然电位曲线对砂泥岩区分较好。通过对声波曲线重构, 提高地震反演的分辨率和精度, 将声波中体现地层背景速度的低频信息与岩性敏感曲线自然电位曲线的高频信息通过信息融合技术进行融合, 形成重构声波曲线。通过重构声波曲线得到的波阻抗对砂岩和泥岩的区分较好(图 4(b))。

图 4

图 4 王家岗地区滩坝砂岩和泥岩波阻抗概率直方图 Fig.4 Histogram of impedance of beach-bar sandstone and mudstone in Wangjiagang area, Dongying Depression

子波的准确提取是地震反演的基础和关键。本文中根据地震资料的统计特性和测井资料的特点, 采用循环迭代法子波提取。其一般提取流程为:首先使用地震道的自相关性计算初始子波的振幅谱, 结合测井曲线确定初始子波的相位谱(假定为线性相位); 然后利用振幅谱和相位谱的信息合成一个理论雷克子波, 用理论子波合成的记录初步进行标定。在初步标定的基础上, 提取井旁地震道实际地震子波, 用该子波制作合成记录并重新调整时深关系。如此反复, 直至得到相位振幅谱变化稳定的精细子波。

频率参数影响反演的分辨率, 控制反演结果的有效频带范围。地震资料有效信息主要集中在中频; 测井曲线属于宽频带信息, 能够为反演提供低频信息和高频信息; 在地震波形指示反演中, 低频成分主要从测井资料中获取, 中频成分通过地震资料求取, 高频成分是在地震波形指示下进行的随机模拟结果, 整个过程是由确定到随机的逐步过渡, 随机成分较少。本次反演的低频参数和高频参数优选如下:低频参数中高通频率取8 Hz, 高截频率取15 Hz; 高频参数中低截频率取70 Hz, 低通频率取110 Hz, 高通频率取300 Hz, 高截频率取350 Hz。

根据东营凹陷沙四段上亚段岩心、录井、测井等资料将碎屑岩滩坝划分为滩砂和坝砂, 其具有不同的沉积特征(图 5)。

图 5

图 5 东营凹陷沙四上亚段滩坝砂岩沉积及地震波形特征 Fig.5 Sedimentary and seismic wave characteristics of beach-bar sand of Es4s in Wangjiagang area, Dongying Depression

滩砂亚相发育于滨浅湖地区, 平行岸线分布^[1], 砂层多而厚度薄, 平面上呈席状或带状广泛分布, 物性较差, 沉积构造主要有波状层理、浪成砂纹层理、冲洗交错层理及压扁层理等, 可见植物根化石。根据沉积砂体厚度、粒度、分布特征等进一步分为滩席和滩脊两种微相^[19]。滩席微相砂体粒度细, 岩性为灰色泥质粉砂岩、粉砂岩, 砂泥岩频繁薄互层, 砂体厚度小于2 m, 自然电位曲线为中低幅锯齿状; 滩脊微相岩性以粉砂岩、泥质粉砂岩、粉砂质泥岩为主, 自然电位曲线为中高幅锯齿状。

坝砂亚相分布在滩砂中, 粒度相对较粗。砂体单层厚, 呈条带状或透镜状, 分布面积较小, 物性较好, 沉积构造主要有平行层理、波状层理、浪成砂纹层理和块状层理, 见植物根化石、介壳化石等。坝亚相包括坝主体和坝边缘两种微相类型。坝主体微相是主体部分, 粒度粗, 岩性以厚层相对均质的灰褐色粉砂岩、粉细砂岩为主, 粒序呈反韵律或复合韵律, 自然电位曲线为中高幅漏斗形、箱形组合型。坝边缘微相岩性主要为灰色粉砂岩、泥质粉砂岩、粉砂质泥岩, 粒序呈反韵律或复合韵律, 自然电位曲线为中幅漏斗形或指形。

滨浅湖泥以灰色、深灰色泥岩、页岩为主, 砂体不发育, 自然电位曲线为低幅或微幅齿形或平滑线形。

碳酸盐滩坝多分布于靠近物源区且无大量陆源碎屑物质供给的湖湾地区^[20]。岩性主要为黄色页状灰岩、深灰色碎屑石灰岩和灰质页岩, 滩单层厚度多小于2 m, 坝单层厚度多大于2 m, 感应电导率、声波时差为低值尖峰, 电阻为高值, 自然电位曲线值较高, 低幅齿形。

根据正演模拟及样本井处沉积特征, 分析不同岩性组合特征的波形特征:①坝砂亚相单砂层厚度较大, 地震反射波形特征为半极值宽度大、幅度小、反射能量中等的中幅单峰波形, 且坝砂反射振幅与坝砂厚度呈正相关, 坝砂发育处地震反射振幅明显大于坝砂不发育处^[21](图 5(a)); ②滩砂亚相砂泥岩薄互层, 地震反射波形特征为中—强单峰波形和复合波波形。当厚层泥岩夹薄层砂岩时, 反射特征为半极值宽度小、幅度大、反射能量较强的中—强单峰波形(图 5(b)); 当双层砂岩夹薄层泥岩时, 反射特征为复波波形, 中等偏下振幅, 随上、下两侧砂岩的多少而呈不对状, 且振幅强处砂岩总厚度小而振幅弱处砂岩厚度反而大^[7](图 5(c)、(d))。

采用约束稀疏脉冲反演和波形指示反演2种反演方法对目的层滩坝砂岩进行预测, 反演结果分别见图 6和图 7。低阻值代表泥岩, 高阻值代表砂岩, 约束稀疏脉冲反演和地震波形指示反演的结果整体趋势相同, 但约束稀疏脉冲反演分辨率较低, 对薄层和横向变化较快的储层预测难度较大; 而地震波形指示反演分辨率高, 反演剖面中—高阻发育规律与连井剖面上砂体发育规律基本一致, 反演阻抗尖灭点可以代表砂岩尖灭点, 砂岩尖灭自然, 同时可以体现出砂岩的内部结构和沉积演化规律。地震波形指示反演具有宽频阻抗信息, 其中低频成分通过井插值模型获得, 中频成分来自地震资料的相对阻抗, 而高频成分是地震波形指示下的随机模拟结果。地震波形是地下地质体地震响应参数的综合反映, 不同地层岩性、物性的变化都会引起反射特征的变化。该反演方法在提高反演分辨率的同时, 充分利用地震波形的横向变化来反映储层空间的相变特征, 体现了相控思想, 是一种真正的井震结合反演方法。反演结果垂向分辨率和横向分辨率明显得到提高, 符合地质规律, 反演效果好。该反演方法在沉积相稳定且断层发育较少, 测井资料丰富, 地震资料信噪比高, 主频越高越好, 地震波形与沉积相变有明显对应关系的地区应用效果较好, 特别适用于成熟探区和开发区块的薄层精细预测。

图 6

图 6 稀疏脉冲反演波阻抗剖面 Fig.6 Parse spike inversion section of impedance

图 7

图 7 波形指示反演剖面 Fig.7 Waveform indication inversion section of impedance

沙四上纯下沉积时期东营湖盆处于断陷初期, 气候由干旱向潮湿转变, 凹陷南斜坡相对平缓的古地貌背景, 决定了广阔的湖区以滨浅湖滩坝沉积为主^[22]。基于层序地层格架下对反演数据体提取地层切片, 从平面分布上预测砂体的展布范围和演化规律(图 8)。低位域时期湖盆面积小, 沉积水体整体较浅。准层序组1沉积期刚进入滨浅湖沉积初期, 物源主要来自东南部, 物源供应不足, 以滩砂沉积为主; 准层序组2沉积期水体稍有加深, 之后水体维持相对稳定, 沉积物由东南向西北推进, 在水下高地处沉积滩坝砂体, 坝砂非常发育, 且砂体分布稳定, 由北部向南砂体厚度逐渐减薄; 准层序组3沉积期基准面持续上升, 基准面上升缓慢, 沉积物以加积为主, 滩砂和坝砂间互发育, 滩坝砂体相对发育, 砂体分布面积增大而泥岩面积缩小, 西北部远离物源处的水下高地沉积碳酸盐岩滩坝。湖侵域时期基准面迅速上升并相对稳定, 只零星发育碳酸盐岩滩坝和滩砂, 以碳酸盐岩滩坝为主。高位域时期湖盆面积较湖侵域变化不大, 水体相对变浅。准层序组5和准层序组6沉积期, 滩坝砂呈孤立椭圆状分布, 远离物源处的水下高地沉积碳酸盐岩滩坝; 准层序组7沉积期基准面略有下降, 沉积物整体具进积特征, 水下高地处沉积滩坝砂体, 西北部远离物源处沉积碳酸盐岩滩坝, 碳酸盐岩滩坝分布面积较大(图 9)。

图 8

图 8 东营凹陷沙四上亚段滩坝反演体切片 Fig.8 Inversion-slice of beach bar sandstone of Es4s in Wangjiagang area, Dongying Depression

图 9

图 9 东营凹陷沙四上亚段滩坝沉积微相 Fig.9 Sedimentary microfacies of beach-bar sandstone of Es4s in Wangjiagang area, Dongying Depression

(1) 本文中反演方法是在层序地层格架下利用沉积学原理, 基于贝叶斯理论下的反演过程, 充分利用地震波形的横向变化反映储层空间的相变规律, 更好地体现了相控的思想, 分辨率明显提高, 预测结果更符合地质规律, 是一种新的井震结合反演方法。

(2) 东营凹陷沙四上亚段在低位域时期准层序组2时期以发育坝砂为主, 砂体厚且横向分布稳定, 准层序组1沉积期和准层序组3沉积期以滩砂和坝砂间互发育, 砂体厚度相对较小; 湖侵域和高位域时期主要发育碳酸盐岩滩坝, 砂质滩坝呈零星状分布。