1. 前言

在我国许多地区存在煤层、生物灰岩或油页岩等地层产生的强反射屏蔽层 [1] [2]，此类地层下的薄储层预测难度较大。在地震资料解释中，传统描述强屏蔽层之下薄储层的方法是开时窗提取属性，直接进行储层描述和预测。但由于强屏蔽层的覆盖作用，传统方法的解释结果不能突出下伏弱反射信息，储层预测效果受到较大影响，因此这种强反射的有效剥离对后续解释的顺利进行至关重要。Lilly等 [3] 提出基于频率域特征值分解的多子波分解与重构方法的强屏蔽处理方法。2013年，佘刚等 [4] 将多子波分解重构技术应用到提高储层预测精度上来，并对部分频率子波剥离煤层技术进行了分析。相较于传统的信号分解方法，匹配追踪算法 [5] 由于它灵活的自适应性、稀疏性和简洁性，更具优势 [5] [6] [7]。王仰华(2007) [8] 系统阐述了基于morlet子波的匹配追踪分解算法，设计了快速高精度匹配追踪算法的三步流程，并将该方法应用到煤层剥离上，取得了预期的效果。在此基础上，王仰华(2010) [9] 结合地震波横向约束线性算子，提出了多道匹配追踪的算法，有效的克服了单道匹配追踪的分解非唯一性问题，并将该方法成功地用于去除煤层强反射和含气层下方地层的低频探测。张军华 [10]、李海山 [11] 及朱博华 [12] 基于匹配追踪去强屏蔽算法，分别对不整合面下隐蔽油藏预测、煤层下含气层识别以及油页岩层下河道砂展布特征刻画进行了应用，提高了储层预测的精度。目前，针对深层–超深层系礁滩相储层类型的去强屏蔽处理研究较少。本文利用基于匹配追踪算法针对元坝地区茅三段顶进行强反射分离处理，并利用处理后的地震数据进行礁滩相储层地震沉积微相识别，有效提高了该区储层定性预测的精度和可靠性。

2. 基于Morlet子波族的匹配追踪原理与算法

匹配追踪算法是一个迭代过程，每次迭代都自适应地从地震信号中提出一个最佳匹配波原子 $g_{{\gamma }_n}$，这里n是迭代次数。在第N此迭代后，地震信号 $f\left(t\right)$可以被表示成如下形式 [13]：

$f\left(t\right)={\displaystyle \underset{n=0}{\overset{N-1}{\sum }}{a}_n{g}_{{\gamma }_n}\left(t\right)}+R^{\left(N\right)}f$(1)

其中， $a_n$为第n此迭代提取的波原子 $g_{{\gamma }_n}$的振幅系数， $R^{\left(N\right)}f$是第N − 1次迭代后地震信号的残差，且 $R^{\left(0\right)}f=f$。本文选择Morlet子波进行原子库构建，Morlet子波函数如下：

$g_{{\gamma }_n}\left(t\right)=m_{{\gamma }_n}\left(t\right)=\exp \left[-\frac{\ln \left(2\right)}{{\pi }^2}\frac{{\omega }_m^2{\left(t-u\right)}^2}{{\sigma }^2}\right]\exp \left[i\left({\omega }_m\left(t-u\right)\right)+\varphi \right]$(2)

其中u为时移时间， ${\omega }_m$为主频， $\varphi$为相位， $\sigma$为尺度因子。

具体计算步骤可分为以下三步 [10]：

1) 利用复数道分析技术，可以估计出匹配波原子的三个参数 $\left(u_n,{\omega }_n,{\phi }_n\right)$，并选择最大瞬时振幅处的时间位置作为 $u_n$，该时间位置对应的瞬时频率为平均中心频率 ${\omega }_n$，瞬时相位作为相移 ${\phi }_n$。对于时间轴尺度因子 ${\sigma }_n$，我们通过扫描搜索得到，其依据的评价标准满足公式：

$g_{{\gamma }_n}=\arg \begin{array}{c}\underset{g_{{\gamma }_n}\in D}{\max }\frac{\left|\langle R^{\left(n\right)}f,g_{{\gamma }_n}\rangle \right|}{\Vert g_{{\gamma }_n}\Vert }\end{array}$(3)

其中， $D={\left\{g_{\gamma }\left(t\right)\right\}}_{r\in \Gamma }$是Morlet波原子族， $\Gamma$为其结构参数空间，符号 $\langle •\rangle$表示向量内积， $\Vert •\Vert$表示向量长度。在这里需要强调的是时间尺度因子满足关系式 $\sigma =\frac{{\omega }_m\Delta t}{\pi }$。对于地震信号来说， ${\sigma }_n$一般取0.5~4之间。特别地，当 $\sigma \to \infty$，波原子为正(余)弦函数；当 $\sigma \to 0$，波原子为一脉冲( $\delta$函数)。

2) 在第一步中估计出了波原子的结构参数 ${\gamma }_n=\left\{u_n,{\omega }_n,{\phi }_n,{\sigma }_n\right\}$后，需要在 $\Gamma$的子空间中优化这些参数，即需要给定一个 ${\gamma }_n$的领域 $\left[{\gamma }_n-\Delta {\gamma }_n,{\gamma }_n+\Delta {\gamma }_n\right]$。

3) 最后，需要估计波原子 $g_{{\gamma }_n}$的振幅系数 $a_n$，其依据是正交投影定理，即

$a_n=\frac{\left|\langle R^{\left(n\right)}f,g_{{\gamma }_n}\rangle \right|}{{\Vert g_{{\gamma }_n}\Vert }^2}$(4)

显然，由(2) (4)式可以导出一个评价关系式：

$\left\{a_n,{\gamma }_n\right\}=\arg \underset{a_n,{\gamma }_n\in \Gamma }{\min }{\Vert R^{\left(n+1\right)}f\left(t\right)\Vert }^2$(5)

式(5)表示最终迭代的残差向量 $R^{N+1}f\left(t\right)$为弱小的噪声数据，不再进行匹配。

如图1所示，地震信号经过数次迭代后，可分解成不同频率、不同相位、不同形状的波原子，当残差达到设定的阀值后停止迭代。对比可知，重构的地震信号与原始信号在形态与结构上基本一致，仅在局部存在较小的误差。

根据这个特点，可选取某一时窗内所需的波原子进行重构，剔除干扰信号对应的波原子，以此弱化干扰信号的影响，达到增强有效信号的目的。

3. 二维模型试算

如图2所示模型，其顶部灰色地层为一套泥岩，纵波阻抗为11.5 kg/m³*m/s；底部绿色地层为灰岩，纵波阻抗为17.5 kg/m³*m/s；黄色部分为储层，纵波阻抗为14.5 kg/m³*m/s；黑色地层为一套相对高速层，纵波阻抗为18 kg/m³*m/s。正演主频为25 Hz，采样间隔为1 ms。当两套储层纵向间隔较小时，靠近泥岩层和高速层反射界面的储层响应消失，仅剩空白反射，储层响应特征和边界无法识别。利用匹配追踪算法识别强反射层波组后进行剔除，可以看到顶部储层响应特征得以突出，呈现亮点反射特征；同时储层尖灭点可以拾取，提高了储层识别能力并改善储层边界识别精度。

4. 实际工区应用

4.1. 工区背景介绍

四川盆地茅口组沉积早期整体为缓坡相沉积，茅三段沉积时期，受峨眉地裂运动活动影响，伴随着强烈的拉张作用，受古地貌控制，梁平–开江陆棚雏形形成，元坝地区演变为镶边台地沉积 [14] [15] [16]。新完钻井A井井深7366米，钻井过程中茅三段见良好油气显示，测试获得105.94万方高产工业气流。

如图3所示，茅三段早期，元坝地区发生相对沉降，以外缓坡泥岩沉积为主；在茅三段晚期，随着拉张作用进一步加强，元坝东拉张下沉形成明显破折，沉积模式由缓坡变为镶边台地，纵向上A井区加积作用形成厚层台地边缘生屑滩沉积 [17] [18]。

元坝地区茅三段开阔台地表现为地层厚度厚，强峰强谷、平行–亚平行反射特征；台地边缘相为碳酸盐岩快速建造，地层厚度最大，向两侧厚度减小，呈中强变振幅丘状反射；陆棚–斜坡相地层厚度薄，强峰强谷、平行–楔状反射，岩性相变为泥灰岩、泥岩。

A井茅三段可分为1、2两个亚段，在2亚段顶部发育一套台缘浅滩储层，厚度为23.7 m，表现为“低伽马、中低速”的特征。由离散合成记录可知，上覆吴二段底部泥岩和吴一段底部沉凝灰岩速度较低，与围岩阻抗差异较大，形成强反射界面，遮挡了下部茅三2亚段储层与围岩地震响应，导致储层段地震响应特征无法有效识别(图4)。

4.2. 针对滩相储层的应用

深层地震资料分辨率较低，信噪比低，且滩相储层分布受相带影响较大，直接去除强反射能量效果较差。结合实际工区情况，对方法进行了两点改进：首先将匹配追踪得到的波原子与原始地震记录进行褶积，然后根据单井分析结果将波原子能量乘以分离系数，如此可有效提高匹配追踪的抗噪性及强反射分离剖面与单井的吻合度，具体表达式如下所示：

$S_{new}=S-\lambda \times m_t\ast S$(6)

其中 $S_{new}$为处理后新的地震信号，S为原始地震信号， $\lambda$为分离系数， $m_t$为最佳匹配波原子。

如图5所示，直接剥离匹配的强反射记录，抗噪性差、横向连续性较差；引入谱分解，将匹配得到的最佳波形原子与地震记录进行褶积运算，提高了匹配的强反射记录与地震记录的相关性，可有效改善强屏蔽处理后地震数据的横向连续性。

如图6所示，分别对A井和B井泥岩段速度、密度曲线进行编辑，去除泥岩低速、低密层，并制作合成记录。对比可知，位于台地边缘的A井，去除泥岩后，储层与下伏致密灰岩之间的反射界面为中–强振幅反射；位于开阔台地的B井，去除泥岩后，储层与下伏致密灰岩之间的反射界面为弱振幅反射。根据参数测试，A井选取分离系数 $\lambda =1.6$，B井选取分离系数 $\lambda =1.2$，强反射分离后的波组与单井正演结果能最佳吻合。

根据已有的测井数据进行单井正演模拟确定去除强反射层影响后的地震波形特征；然后通过调整匹配追踪分离系数将强反射削弱至与正演模拟结果相匹配，从而确定分离系数大小；最后通过相带控制进行插值，获得分离系数 $\lambda$的平面分布。由于该区斜坡–陆棚和台地边缘相带钻井较少，均只有一口，开阔台地内各井分离系数 $\lambda$相近，故而直接根据相带赋予不同的 $\lambda$值。加入分离系数约束后，在A井、C井井点处强屏蔽处理剖面与实际井震标定吻合度更高。位于开阔台地相的C井在加入分离系数约束的强屏蔽处理后，其茅三段顶部为弱振幅反射，而位于台地边缘相的A井在处理后，其茅三段顶部为中–强振幅反射，两口典型井的处理效果都与沉积相带的地震响应特征认识更加匹配(图7)。

4.3. 在沉积微相识别中的应用

结合地质认识和单井替换实验可知，元坝地区茅三段台地边缘相包含三个微相：① 台缘丘滩相，发育于台缘靠斜坡一侧，沉积时期水动力相对较强，以纵向加积为主，地层厚度较大，整体呈丘状，茅口组顶界表现为波峰或复合波峰反射。② 台缘浅滩相，发育于台缘靠台地一侧，水动力相对较弱，发育浅滩和滩间，地层厚度减薄，内部主要由单波峰和上强下弱复波地震反射构成。③ 滩间，发育在滩体之间，表现为茅三2亚内部出现强波峰反射(图8上)。

去除强反射能量后，台缘丘滩相和台缘浅滩相均为靠近茅三段顶部的亮点反射特征，其中台缘丘滩丘状外形明显，台缘浅滩表现为微幅丘状–平行外形特征，台缘丘滩波峰反射能量强于台缘浅滩。滩间为靠近茅三2亚中部的波峰反射，下凹外形明显(图8下)。

对比去除强反射轴前后剖面可知，台地边缘的三种微相在原始剖面上反射特征差异较小，界限不清；去除强反射层能量后，三种微相界限清晰，反射特征和外形特征差异明显，分异性较好。

通过提取相同时窗的地震属性进行对比分析可见，原始地震数据(图9左)在茅三2亚段振幅属性台地边缘与斜坡–陆棚相边界清晰，但台地边缘与开阔台地的边界无法有效识别，同时台地边缘相带内振幅属性的响应受上覆吴家坪泥岩/沉凝灰岩的强反射的影响，内部能量差异小，不能有效表现滩相储层的响应。去强轴后地震数据(图9右)相带边界清晰，靠近斜坡一侧的丘滩与靠近台地一侧的浅滩沿台地边缘成排成带展布，振幅能量差异明显，滩间可能为滩体发育沟通的水道；同时，由振幅能量变化可推测，该区滩体虽大面积分布，但内部横向存在一定的非均质性，整体上靠近台缘外侧的滩体储层可能更好。

5. 结论

基于匹配追踪算法的强反射分离方法可以恢复强反射层下薄储层真实地震响应特征。将匹配波原子与原始地震资料卷积后，可有效改善分离剖面的横向连续性；对于多个薄层反射能量相互干涉的情况，综合地质沉积相带确定分离系数大小，可使强反射分离后剖面地质意义更加明确。将其应用在元坝地区茅三段台缘滩相碳酸盐岩储层综合预测中，相比前期预测结果，沉积微相细节变化更丰富，相带边界展布更清晰，有效提升了该区地震微相识别精度。

基金项目

国家科技重大专项《四川盆地下组合油气富集规律与勘探目标评价(2017ZX05005-003-006)》资助。

NOTES

^*通讯作者。

参考文献