地震分辨率分为纵向分辨率和横向分辨率两个方面。纵向分辨率的概念最早由Rayleigh^[16]给出, 即相邻两个反射界面的分辨率极限为1/4波长, 厚度小于1/4波长的地质体即可以定义为薄层, 利用地震数据无法分辨。为解决薄储集层预测问题, 许多学者开展了大量研究, Widess^[17]提出在理想情况下, 根据振幅横向变化能够识别任意厚度的薄层; 曾洪流^[18]提出了“ 横向检测率” 的概念, 利用地层切片可以检测厚度小于1/4波长的地质体横向变化。

通过一个正演模型来讨论薄层的纵向分辨率和横向分辨率。地质模型为：在泥岩背景下, 连续发育一套厚度为3 m的砂岩, 砂岩之上发育一套厚度为0~3 m的透镜状薄砂岩（见图1a）, 泥岩夹层厚度为3 m, 砂泥岩地震波速度分别为3 500 m/s和2 800 m/s, 密度分别为2.65 g/cm³和2.26 g/cm³。利用35 Hz零相位雷克子波进行褶积。从得到的正演剖面（见图1b）可以看出, 纵向上无法直接分辨出上覆薄砂岩。但由于薄砂岩的发育, 地震波形横向上发生了变化, 提取图1b中最大波峰处的振幅和频率（见图1c）, 可以看出地震振幅和频率横向上发生了很大变化, 即地震波形包含了薄层信息。因此可以得出结论：纵向上地震振幅无法“ 分辨” 薄层, 但可以利用横向波形变化“ 识别” 薄层。

图1

	Figure Option ViewDownloadNew Window
	图1 地质模型与正演剖面

地震波形包含了地震运动学和动力学多种信息, 是地质沉积作用、岩性岩相组合, 储集层物性以及流体等多种地质信息的综合响应。曾洪流^[19]提出了“ 地震沉积学” 的概念, 充分利用地震资料的横向高分辨率, 提高了薄层的预测精度, 在薄储集层预测中得到了广泛应用。

本文在地震沉积学技术的基础上, 引入测井资料, 将地震资料的横向高分辨率和测井曲线的纵向高分辨率有机结合起来, 首先要解决的问题是如何建立地震波形和测井高频信息之间的内在联系。

在实际资料中, 相似的沉积特征往往具有相似的岩性组合, 相似的岩性组合往往具有相似的地震波形特征。选取一个陆相薄储集层实际资料进行分析, 图2a为A、B两口井的地震剖面, 该地震资料信噪比高, 但频率较低, 主频约为20 Hz, 频宽约为10~35 Hz。提取过两口井同一层段井旁道地震波形, 由于两口井存在深度差异, 因此首先对两口井的地震波形经过顶面对齐和厚度一致性校正, 然后进行叠合对比（见图2b）, 从图中可以看出, 两者相关系数达94%。因此, 两口井的井旁地震道具有相似的地震波形特征。

图2

	Figure Option ViewDownloadNew Window
	图2 A井和B井地震波形特征

再进一步分析测井曲线的相似性特征：首先分析与地震反射直接相关的纵波波阻抗曲线, 两口井原始纵波波阻抗曲线相关系数为75%（见图3a）, 远低于地震波形的相似性, 分析其原因, 测井曲线具有非常高的频率, 高频成分的差异降低了其相似性。依次对测井曲线进行滤波, 分别保留0~500 Hz、0~400 Hz、0~300 Hz、0~200 Hz以及0~100 Hz（见图3b— 图3f）, 相关系数逐步提高到88%、92%、93%、94%和95%, 达到了非常高的相似性。由此可见, 由于测井曲线高频信息的差异性导致原始测井曲线相似性较低, 对测井曲线逐步降低频率进行滤波, 测井曲线相似性逐步提高, 当测井曲线频率达200~300 Hz时, 相关系数就达到了地震波形的相似性, 因此可以建立起地震波形与测井高频信息之间的联系, 为高分辨率地震反演提供了依据。开展地震波形指示反演时, 可以通过分析不同频率下测井曲线相似性和目标储集层的厚度共同确定反演的最佳截止频率。

图3

	Figure Option ViewDownloadNew Window
	图3 A、B井纵波波阻抗曲线不同频带范围滤波

当储集层特征复杂、纵波波阻抗无法区分岩性时, 需要借助自然伽马等曲线区分岩性, 但是基于褶积模型的地震反演无法预测自然伽马曲线。借助分析纵波波阻抗曲线相似性的思路来分析自然伽马曲线, 如图4a, 两口井原始自然伽马曲线相关系数只有71%, 依次对自然伽马曲线进行滤波, 保留0~500 Hz、0~400 Hz、0~300 Hz、0~200 Hz以及0~100 Hz频率成分（见图4b— 图4f）, 相关系数逐步提高到88%、90%、92%、93%和95%, 达到了和纵波波阻抗基本一致的相关性。

图4

	Figure Option ViewDownloadNew Window
	图4 A、B井自然伽马曲线不同频带范围滤波

从沉积学的角度分析, 在等时地层格架内, 地震波形的变化反映了岩性组合的变化, 岩性组合是沉积相或者地震相的表现形式。从测井曲线的角度分析, 波阻抗和自然伽马等其他曲线均可以反映岩性组合或沉积相的变化。因此, 低频地震信息与高频测井信息的对应关系不仅适用于波阻抗曲线, 同样适用于自然伽马等非波阻抗曲线, 这就奠定了利用低频地震资料和高频测井资料信息进行高分辨率反演（或模拟）的基础。需要说明的是, 对于波阻抗曲线为地震波形指示反演, 而对于自然伽马等非波阻抗曲线为地震波形指示模拟。

地震波形参数和井点属性的对应关系可以定义为一个$n\times m$阶矩阵A, 则有：

$\mathbf{A}\text{= }\mathbf{U}\Lambda \ {{\mathbf{V}}^{T}}$ （1）

式中, U和V为地震波形数据和井点属性数据, 分别为$n\times n$阶和$m\times m$阶正交矩阵, V^T为V的共轭转置, Λ 为$n\times m$阶非负实数对角矩阵, 表示地震波形数据和井点数据的相关性。

奇异值分解即对A进行正交分解, 当矩阵的秩为r时, 则矩阵A可分解为r个本征向量的代数和, 则矩阵A的总能量可表示为：

${{\left\| \mathbf{A} \right\|}^{2}}=\sum\limits_{i=1}^{r}{\delta _{i}^{2}}$ （2）

矩阵A经过（2）式进行奇异值分解之后, 可以用前r个非零奇异值对应的奇异向量表示矩阵A的主要特征, 即实现地震波形动态聚类分析, 得到不同储集层类型地震波形与测井曲线特征的对应关系, 建立初始的样本集。

小波变换将数字信号从时间域变换到频率域, 可以表征信号在时间域和频率域的特征, 能够同时定量预测信号的低频稳定部分和高频突变部分。

对前文建立的样本集中的测井曲线分别用不同的截止频率l利用（3）式开展离散小波变换, 将测井曲线分解为低— 中频宏观特征信息、高频细节信息和超高频噪声信息3部分, 需要说明的是, 由于测井曲线具有非常高的频率, 以声波时差曲线为例, 其频率高达20 kHz^[21], 因此, 这里定义的低— 中频大致为100~200 Hz, 甚至300~400 Hz的频率范围, 即本文中相对于地震频率范围定义的“ 测井高频信息” , 提取的低— 中频宏观特征信息即该样本集中所有测井曲线的共性结构, 该共性结构可以作为波形指示反演的初始模型。

$O(l)=\arg \left( \min \left\| \mathbf{W}-\overline{\mathbf{W}} \right\| \right)=\arg \left( \min \left\| \int_{\text{0}}^{l}{\psi }(\omega , t)\text{d}\omega -\overline{\mathbf{W}} \right\| \right)$ （3）

地震反演的基础是褶积模型为：

$\mathbf{d}=\mathbf{Gm}+\mathbf{n}$ （4）

假设噪声n满足高斯分布：

$P(\mathbf{n})=\frac{1}{\sqrt{2\text{ }\!\!\pi\!\!\text{ }{{\sigma }^{2}}}}\exp \left( -\frac{1}{2{{\sigma }^{2}}}{{\mathbf{n}}^{\text{T}}}\mathbf{n} \right)$ （5）

将（5）式代入（4）式, 建立地震数据似然函数：

$P\left( \mathbf{d}|\mathbf{m}, \mathbf{I} \right)=\frac{1}{{{\left( \sigma \sqrt{2\text{ }\!\!\pi\!\!\text{ }} \right)}^{N}}}\exp \left[ -\frac{\sum\limits_{i=1}^{N}{{{\left( \Delta {{\mathbf{d}}_{i}}-\mathbf{G}\centerdot \Delta {{\mathbf{m}}_{i}} \right)}^{2}}}}{2{{\sigma }^{2}}} \right]$ （6）

在贝叶斯反演中, 假设弹性参数模型m也符合高斯分布, 可以得到模型的先验分布为：

$P\left( \mathbf{m}|\mathbf{I} \right)=\frac{1}{\sqrt{2\text{ }\!\!\pi\!\!\text{ }\left| {{\sigma }_{m}} \right|}}\exp \left( -\frac{{{\mathbf{m}}^{\text{T}}}\mathbf{m}}{2{{\sigma }_{m}}} \right)$ （7）

将数据条件概率分布与模型先验概率分布的乘积作为模型的后验概率分布函数, 其具体表达式为：

$\begin{align} & P\left( \mathbf{d}|\mathbf{m}, \mathbf{I} \right)=\frac{1}{{{\left( \sigma \sqrt{2\text{ }\!\!\pi\!\!\text{ }} \right)}^{N}}}\exp \left[ -\frac{\sum\limits_{i=1}^{N}{{{\left( \Delta {{\mathbf{d}}_{i}}-\mathbf{G}\centerdot \Delta {{\mathbf{m}}_{i}} \right)}^{2}}}}{2{{\sigma }^{2}}} \right]\times \\ & \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \frac{1}{2{{\text{ }\!\!\pi\!\!\text{ }}^{\frac{3}{2}}}\sqrt{{{\left| {{\sigma }_{\Delta m}} \right|}^{3}}}}\exp \left( -\frac{\Delta {{\mathbf{m}}^{\text{T}}}\Delta \mathbf{m}}{2{{\sigma }_{\Delta m}}} \right) \\ \end{align}$ （8）

这样对于给定的地震波形d, 可按照后验概率应用Gibbs抽样法计算得到模型m的期望值。

（8）式中, 概率最大时得到的解, 即为反演的最终解。对（8）式两端求对数, 得到目标函数：

$O\left( \mathbf{m}|\mathbf{d}, \mathbf{I} \right)=-\frac{1}{2{{\sigma }^{2}}}\sum\limits_{i=1}^{N}{{{\left( \Delta {{\mathbf{d}}_{i}}-\mathbf{G}\centerdot \Delta {{\mathbf{m}}_{i}} \right)}^{2}}}-\frac{\Delta {{\mathbf{m}}^{\text{T}}}\Delta \mathbf{m}}{2{{\sigma }_{\Delta m}}}$ （9）

为了使后验概率最大, 对式（8）关于模型参数Δ m求导得到：

${O}'\left( \Delta \mathbf{m} \right)=\frac{1}{{{\sigma }^{2}}}\left( {{\mathbf{G}}^{\text{T}}}\mathbf{G}\Delta \mathbf{m}-{{\mathbf{G}}^{\text{T}}}\Delta \mathbf{d} \right)-\frac{\Delta \mathbf{m}}{{{\sigma }_{\Delta m}}}$ （10）

令${O}'\left( \Delta \mathbf{m} \right)$=0, 则模型扰动量为：

$\Delta \mathbf{m}={{\left( {{\mathbf{G}}^{\text{T}}}\mathbf{G}+\frac{{{\sigma }^{2}}}{{{\sigma }_{\Delta m}}}\mathbf{I} \right)}^{-1}}{{\mathbf{G}}^{\text{T}}}\Delta \mathbf{d}$ （11）

使用迭代模型扰动量的方法逼近样本数据, 得到最终的高分辨率反演结果。

地震波形指示反演方法具有以下特点：①地震波形指示反演利用地震波形横向变化驱动高频测井信息特征实现高分辨率反演。反演结果纵向上与测井高频信息吻合, 具有纵向高分辨率（当测井曲线高频共性结构为200~300 Hz时, 薄层预测的分辨率为2~3 m）。同时反演结果横向上遵循地震波形的变化, 也具有横向高分辨率, 因此, 地震波形指示反演可以同时提高反演结果的纵横向分辨率, 是一种高精度的反演方法。②井间唯一存在的数据为地震数据, 地震波形的横向变化体现了岩相组合的变化。地震波形指示反演通过地震波形横向变化代替变差函数空间域插值模拟, 实现了地震相自动控制下的反演, 克服了传统相控反演需要人为事先给定沉积相而导致的主观性, 是真正意义上的相控反演。③地震波形指示反演方法不仅可以实现高分辨率波阻抗反演, 还可以实现自然伽马、电阻率和孔隙度等非波阻抗参数的相控高分辨率模拟, 突破了地震反演只能得到波阻抗结果的局限性, 对于利用地震信息进行储集层参数模拟是一次巨大的进步。

研究区位于松辽盆地大庆长垣北部白垩系姚家组葡萄花油层。前人研究认为, 该地层沉积时期大庆长垣北部地区受北部水系控制形成大型河流— 三角洲复合体, 广泛发育近南北向水下分流河道沉积^[22]。砂体埋藏深度为900~1 200 m, 单砂体厚度为0.2~15.0 m, 47%的砂体厚度为1~3 m（见图8a）, 呈现出典型的薄互层特征。油田进入开发后期阶段, 挖潜的主要对象从厚层砂体转变为厚度小、横向变化快的砂体, 因此需要精细刻画薄互层砂体。

图8

	Figure Option ViewDownloadNew Window
	图8 砂体厚度特征和岩性解释图版

利用测、录井资料开展岩性测井解释, 建立砂泥岩解释图版（见图8b）。砂岩纵波波阻抗小于7 600m· g/(s· cm³), 因此, 纵波波阻抗可以识别砂泥岩, 可以开展地震波形指示反演得到纵波波阻抗结果预测砂岩分布。

本次研究地震三维工区面积6.25 km², 工区内有钻井39口, 且分布较为均匀。为验证波形指示反演算法的合理性, 开展不同数量井参与反演实验：第1组实验39口井全部参与反演（见图9a、图9b）; 第2组实验模拟勘探阶段少井的情况, 从39口井中随机选取10口井参与反演, 其余29口井为验证井（见图9c、图9d）; 第3组实验模拟评价阶段不均匀井网的情况, 北边25口井参与反演, 南部14口井为验证井（见图9e、图9f）。从上述3组实验的反演剖面（见图9a、图9c、图9e）上可以看出, 无论是参与井还是验证井, 反演结果与钻井吻合程度都比较高, 并且3组实验反演结果高度相似。按照岩性直方图确定的砂岩门槛值, 提取了上油层组砂岩（反演剖面中黑色虚线框, 见图9a、图9c、图9e）厚度图（见图9b、图9d、图9f）, 从图中可见, 3个厚度图形态基本一致, 都表现出了近南北向发育的两个河道砂体分支, 主河道边界清晰, 河道展布形态符合地质规律。

图9

	Figure Option ViewDownloadNew Window
	图9 不同数量井参与反演结果对比

为进一步定量描述反演预测的精度, 提取反演预测的上油组砂体厚度, 与测井解释砂岩厚度进行对比分析。图10a为39口井全部参与反演对比结果, 反演预测砂岩厚度与测井解释厚度相关系数达90.8%; 图10b为模拟勘探阶段随机选取10口井参与反演的结果, 红色点为参与反演井, 相关系数达91.2%, 蓝色点为验证井结果, 相关系数达80%。图10c模拟评价阶段北边25口井参与反演结果, 红色点为参与反演井, 相关系数达90.2%, 蓝色点为验证井结果, 相关系数达80.8%。

图10

	Figure Option ViewDownloadNew Window
	图10 反演预测砂体厚度与测井解释砂岩厚度交会图

通过上述不同方法反演结果对比以及不同数量井参与反演结果表明, 波形指示反演结果具有两方面的优势：①分辨率高, 能识别厚度为2~3 m的薄互层, 同时反演精度高, 参与井吻合率达90%, 验证井吻合率达80%; ②地震波形驱动测井曲线实现高分辨率反演, 具有相控思想, 不受井数量和井位分布影响, 反演结果符合地质规律。

地震波形指示反演对大庆长垣薄互层的识别精度达到了2 m, 通过地震波形指示反演与主流的确定性反演方法— — 稀疏脉冲反演结果的对比分析, 进一步证明地震波形指示反演具有更高的精度。

目的层三维地震资料主频为40 Hz, 纵波速度为3 000 m/s左右, 按照地震分辨率1/4波长理论, 地震资料可识别的最大厚度约为18.75 m, 显然不满足薄互层砂体识别的要求。从过井的地震剖面上（见图11a）可以看到, 地震资料无法识别葡萄花油层的薄互层砂体。首先开展常规稀疏脉冲反演, 从得到的纵波波阻抗剖面上（见图11b）可以清楚地看到, 稀疏脉冲反演分辨率与地震资料大致相当, 只能大致识别大套的砂层组, 无法预测薄互层砂体。然后开展波形指示反演, 从得到的纵波波阻抗反演剖面上（见图11c）, 可以直观地看到, 波形指示反演结果的纵向分辨率明显高于稀疏脉冲反演, 能够识别薄互层砂岩, 同时, 砂岩横向变化特征符合地震特征, 能够较好地刻画砂体横向边界, 即波形指示反演也具有较高的横向分辨率。

图11

	Figure Option ViewDownloadNew Window
	图11 地震剖面与不同方法反演效果对比

为进一步定量论证两种反演方法的精度, 分别从两个反演数据体上提取井点处曲线（以图11中W3井为例）, 如图12中第1、2道红色曲线为实测纵波波阻抗曲线（为便于和反演结果对比, 实际纵波波阻抗曲线以2 m采样率进行了滤波）, 可以看出, 波形指示反演结果与实测曲线吻合度远高于稀疏脉冲反演。分别对两种反演结果按照砂岩门槛值进行岩性解释, 与实际测井岩性解释结果对比（见图12）可以看出, 波形指示反演解释岩性与测井解释岩性吻合度高, 可识别厚度为1.8~3.0 m的砂岩; 而稀疏脉冲反演结果与测井解释岩性吻合度低, 可识别的砂体厚度为10 m。由此可见, 波形指示反演无论是吻合率还是反演精度都远高于稀疏脉冲反演。

图12

	Figure Option ViewDownloadNew Window
	图12 稀疏脉冲反演与波形指示反演结果对比