长堤油田位于济阳坳陷东北部^[1]，紧邻黄河口凹陷，具有良好的油源条件，中生界与上覆古近系呈角度不整合接触，形成不整合圈闭，发育不整合油藏，且具有埋深浅、自然产能高的特点，是现实的效益储量阵地。

地层油藏勘探已经有较为成熟的理论和行之有效的技术^[2-7]，对于长堤油田不整合油藏的描述难点在于：中生界广泛发育大套致密高速砂岩，但有效储层仅为不整合面之下20m以内的低速疏松砂岩。高速砂岩与上覆古近系低速泥岩形成连续低频强反射T_r同相轴，有效储层的反射波形被T_r同相轴淹没，常规地震资料预测缺乏有效技术手段^[8-9]。

针对强反射屏蔽有效信息的问题，多位学者展开了研究，已经明确了强反射背景下地震反射现象的形成机理^[10-13]，即不整合点附近，由于地层厚度减薄发生调谐效应，且在高速层屏蔽、低地震分辨率的共同影响下，使得地震反射能量减弱，反射波形提前终止，形成空白反射。张军华^[9-10]、刘保国^[14]等分别利用岩性反演、提高分辨率及多属性融合等常规技术描述了强反射屏蔽的砂体展布；秦雪霏^[15]等利用多子波分解与重构的方法去除了煤层强反射，提高储层描述精度；Wang^[16]、李海山^[11]、李传辉^[17]等利用匹配追踪算法分离煤层形成的强反射。以上方法应用于长堤地区发现，常规技术虽然有一定效果，但是描述精度不高；与多子波分解与重构相比，匹配追踪算法时频分辨率更高，但是匹配追踪算法由于搜索原子库寻找匹配原子方法的不同，计算效率和适用性有待提高。

本文针对以上问题，提出了一种强反射背景下描述不整合圈闭有效储层的方法。为改进匹配追踪算法，加入局部频率约束算子以提高运算稳定性和速度；结合实际钻测资料构建一维和二维模型分别验证了该方法的可行性，并明确了振幅属性可以作为预测有效储层的敏感参数。在模型试验的基础上，分离长堤地区中生界顶面强反射，得到能够凸显有效储层的地震数据体，利用振幅属性描述了不整合面之下的有效储层平面展布规律。

1 基于局部频率约束的动态快速匹配追踪方法

匹配追踪方法是一种有效的信号稀疏分解方法，目前最常用的是动态匹配追踪算法，其参数搜索方式可以表示为^[18]：

$ \begin{array}{l} {\gamma _{\rm{n}}} = \{ {u_0}, \omega \in U[\omega ({u_0}), {\delta _{\rm{ \mathsf{ ω} }}}], \varphi \in U[\varphi ({u_0}), {\delta _{\rm{ \mathsf{ φ} }}}]\} \\ {u_0} = \{ {t_0}:A{\rm{ }}({t_0}) = {\rm{max}}[A(t)]\} \end{array} $

式中  γ_n——小波字典的控制参数集合；

A(t)——瞬时振幅；

A(t₀)——t₀时刻瞬时振幅；

φ——瞬时相位；

ω——瞬时频率；

u₀——最大瞬时振幅包络A(t)对应的时间中心；

U[ω(u₀), δ_ω]——以瞬时属性为中心的频率搜索邻域；

U[φ(u₀), δ_φ]——以瞬时属性为中心的相位搜索邻域；

δ_ω、δ_φ——根据信号时频特征定义的参数搜索半径。

这种算法优势在于搜索频率范围时使用动态局部搜索比全局搜索方式增加了计算速度。采用对瞬时相位求导得到瞬时频率，然而此方法求得的瞬时频率不稳定，当信号信噪比较低时会出现不切实际的数值。

本文提出局部频率替代瞬时频率，它能快速且系统地处理噪声强的数据甚至是部分缺失的数据。用瞬时振幅A（t）和瞬时相位φ（t）来表示复数道，如下式：

$ \begin{array}{l} c\left( t \right) = A(t){{\rm{e}}^{{\rm{i}}\varphi \left( t \right)}}\\ A(t) = \sqrt {{x^2}\left( t \right) + {h^2}\left( t \right)} \\ \varphi \left( t \right) = {\rm{arctan}}\left[ {h{\rm{ }}\left( t \right){\rm{ }}{x^{ - 1}}\left( t \right)} \right] \end{array} $

式中  c（t）——地震道的复数道；

x（t）——地震道；

h（t）——实际地震道x（t）的希尔伯特变换。

由此得到对于瞬时频率属性的地震道的数学表达式为：

$ \omega \left( t \right) = \varphi '\left( t \right) = \frac{{x\left( t \right)h'\left( t \right) - x'\left( t \right)h\left( t \right)}}{{{x^2}\left( t \right) + {h^2}\left( t \right)}} $

从上式可以看出瞬时频率是瞬时相位的导数，且瞬时频率是两个信号进行相除，写出瞬时频率的矩阵表达式，即：

$ \mathit{\boldsymbol{W}} = {D^{ - 1}}n $

式中  W——瞬时频率向量；

n——匹配原子数目；

D——对角矩阵算子。

为了避免出现“除以零”导致计算错误的现象，在分母上加一阻尼项ε，成为：

$ \mathit{\boldsymbol{W}}{\rm{inst}} = {\left( {D + \varepsilon \mathit{\boldsymbol{I}}} \right)^{ - 1}}n $

其中

$ \mathit{\boldsymbol{I = }}\left[ {\begin{array}{*{20}{c}} 1&0&0\\ 0&1&0\\ 1&0&{ \cdots 1} \end{array}} \right] $

式中  W_inst——加阻尼项的瞬时频率向量；

ε——阻尼项；

I——单位矩阵。

阻尼项ε可以防止瞬时频率受噪声和不稳定性因素的影响，却不能使其变得稳定。故选取形状正则化算子，结合整形光滑算子S，进行进一步约束，用迭代反演法进行求解，如下式：

$ {\mathit{\boldsymbol{W}}_{{\rm{loc}}}} = [{\lambda ^2}\mathit{\boldsymbol{I}} + S{(D - {\lambda ^2}\mathit{\boldsymbol{I}}]^{ - 1}}Sn $

式中  W_loc——加入尺度因子的瞬时频率向量；

S——整形光滑算子；

λ——尺度因子。

利用迭代反演方法求解时，通过尺度因子λ进行缩放，可以保留原始物理维度，从而有效加快收敛速度。

2 强反射层分离 2.1 强反射分离流程

在地层结构复杂、沉积类型多变的地区，地震同相轴通常是多个地层或岩性界面的叠加反映^[19]。匹配追踪算法可实现对每个界面子波的匹配识别，将干扰信息进行分离，从而得到有效信息。具体流程如下：

(1）利用测井资料进行合成记录标定，分析目的层段地震资料时频特征，确定强反射界面的主频。

(2）根据时频分析特征，提取地震数据体子波，构建超完备子波字典。

(3）利用局部频率约束的匹配追踪算法搜索子波字典，实现强反射特征原子的匹配拾取，并获得最终的强反射干扰记录。

(4）结合实际井区地质背景，建立模型，分析反射系数、频率、相位、尺度因子λ，实现强反射分离。

(5）结合实钻资料，检验分离方法的合理性。

2.2 模型试算

为验证强反射分离方法的可行性，依据长堤地区中生界的实际钻测资料构建模型进行正演。统计工区内38口井的速度特征，古近系底部为沙三段泥岩，速度为3300m/s；中生界砂岩速度在4500m/s左右，泥岩速度为3800m/s。以此速度特征为基础数据，首先建立一个不含不整合强反射系数的一维地质模型（图 1a列），图中红框内上部0.1s处的反射系数代表了一套储层，正演得到的原始地震记录（图 1b列）能够反映该套储层。在该套储层处加入强反射系数，使得该套储层的反射系数包含于不整合强反射系数之中（如图 1c列），由于波阻抗差异较大，在0.1s处出现一个较强的反射系数（图 1d列），在地震记录(图 1e列)中出现强反射同相轴，将储层的反射信息屏蔽。利用局部频率约束的匹配追踪算法，匹配出强反射（图 1f列），该强反射分离后得到地震记录（图 1g列)，其中黑色波形显示为分离强反射后的地震记录，红色为不含强反射系数模型得到的地震记录，可以看出两种地震记录波形基本吻合。正演结果验证了该方法对分离不整合地层强反射的有效性与准确性，可以利用该方法得到不整合面附近被屏蔽的有效反射信息。

2.3 尺度因子λ确定

从模型试算中可以看出，在强反射分离过程中，有两个关键点：第一，强反射特征原子的最优匹配和快速拾取。运用局部频率约束的匹配追踪算法得以解决，该方法快速搜索子波字典，实现强反射特征原子的最优匹配拾取，并获得最终的强反射干扰记录。第二，为了使分离后的地震数据，既不被强反射所屏蔽又能最大限度凸显有效信息，尺度因子λ的选择至关重要^[20]。尺度因子越大，则分离效果越清楚，但如果完全去除，势必会同时去除掉部分有效信息。为解决这一问题，继续设计模型对尺度因子进行实验。如图 2所示，a列为构建的模型，有效储层和强反射位于0.25~0.35s之间，b列为原始地震数据，c列为匹配拾取的强反射记录，d列、e列、f列、g列分别为λ取值为0.2、0.4、0.6、0.8时分离后的地震数据。可以看出d、e两列强反射分离不够，有效储层没有显现；而g列分离过度，同时去除了有效储层的信息；只有f列红色波形（无强反射时波形）与黑色波形（分离后的波形）基本重合，做到了强反射合理分离，故而长堤地区尺度因子λ选为0.6。

3 储层描述

为准确描述储层的平面展布规律，进行二维模型的试算。依据工区实际地质情况建立二维地质模型（图 3a），砂岩厚度小于20m，且砂体之间横向有变化点、纵向叠置，其上稳定发育高速层。图 3b为理论模型的原始正演地震记录的波形特征，由于受到高速层的影响，底部砂体的有效反射信息被强反射屏蔽，砂体之间不能分辨。图 3c为强反射匹配识别分离后的结果，在砂体尖灭处波形特征有变化，可分辨砂体期次，但追踪描述单砂体时砂体边界依然有多解性，不能完全满足识别需求。

进一步对多种属性进行分析，发现振幅属性能够体现砂体之间的变化点。从图 3d可以看出，振幅属性能够清楚地表示砂体之间的尖灭点。由图 4可以看出，原始地震资料振幅与砂体厚度没有相关性（图 4a）；分离强反射之后的振幅与砂体厚度呈良好的正相关关系（图 4b）。因此可以利用振幅类属性预测被强反射屏蔽的储层的平面展布规律。

综上所述，强反射背景下的不整合圈闭储层描述流程如图 5所示。

4 应用效果

长堤地区中生界储层发育程度受不同地层岩性和不整合面控制较为明显，中生界三台组和蒙阴组储层较为发育，岩性多为粗砂岩到细砂岩，并且不整合面附近储层由于长期遭受风化剥蚀，储层改造作用明显。通过已钻井分析，有效储层段多集中于紧邻不整合面之下0~20m之内，孔隙度为10%~15%，渗透率为80~150mD，大于20m相同岩性储层孔渗性变差，孔隙度只有5%~8%，渗透率仅为20~50mD。经过钻井统计，不整合面附近疏松储层速度为3900m/s，而远离不整合面相同层位中生界速度为4500m/s，有效储集层段速度明显降低了600m/s，但是由于该段厚度较薄，经过精细标定可知，在地震资料上不能形成单独反射，淹没于中生界高速砂岩之中。实际地震资料频谱分析表明，目的层段主频仅为20Hz。若平均速度为4500m/s，则可分辨储层的厚度为60m，而长堤地区有效储层厚度小于20m，故用常规地震属性无法进行准确描述^[21-23]。

首先利用局部频率约束的匹配追踪算法分离强反射，效果明显。T_r同向轴位于2100ms之上，原始剖面表现为连续稳定、低频强振幅，桩213-斜10井在3821~3828m发育一套有效储层（图 6c红框内），标定后位于T_r轴之上（图 6a红框内），被T_r同相轴所屏蔽，不能识别该套储层；强反射分离后，从剖面(图 6b)可以看出，T_r反射同相轴成为中弱振幅、时断时续，尤其在井点处由于储层发育波形有明显变化，可对该套储层进行识别描述。

以分离强反射的地震数据体为基础，沿层包含地震同相轴提取均方根振幅属性^[24]。为对比局部频率约束算子的算法优势，选取同样的属性、开取相同的时窗、应用相同的成图参数进行储层预测（图 7），图 7a是依据常规匹配追踪算法分离强反射后的地震数据所做的储层预测图，图 7b为加入局部频率约束算子的匹配追踪算法分离强反射后的地震数据所做的储层预测图。总体看来，两张图都反映了长堤地区三台组Ⅳ砂组有效储层的平面展布规律，整体呈北西向条带状展布，主要发育在临近尖灭点处。但是图 7b中有效储层的边界刻画更为清晰，同时也更符合实钻井况。Ⅳ砂组共有5口直井钻遇，其中4口井（桩205、桩12、桩205-8、桩205-5）有效储层厚度在6~21m范围，桩11井有效储层不发育，而图 7a显示为储层发育区，图 7b更为符合。桩205块为该层的探明区块，而图 7a预测范围比实钻范围明显减小，图 7b与实际钻探效果吻合。据图 7b部署的桩205-平10井在水平段钻遇有效储层2层102.2m，即为桩205-5井南部红色区域的沿桩205-平10井轨迹的水平宽度，钻探结果与预测范围完全吻合，为长堤地区中生界不整合油藏勘探提供了部署依据。

5 结论

本文以长堤地区中生界不整合油藏为例，针对被强反射屏蔽的不整合圈闭有效储层，形成局部频率约束的动态快速匹配追踪算法分离强反射。构建符合长堤地区地质特征的模型分析和实际资料处理应用表明：

(1）这一算法能够快速、稳定地对数据平滑处理，有效避免常规匹配追踪算法可能会出现的瞬时频率“异常值”。同时，应用该算法处理后的数据对于储层的边界刻画更为精确。

(2）强反射分离时，尺度因子λ的选择需要反复实验，保证最大限度地凸显储层有效信息，达到强反射波组的合理分离。

(3）强反射分离后的数据体能够识别砂体的变化点，明确砂体的发育期次。利用振幅类属性可以预测储层的平面展布规律，与实钻井吻合较好，具有良好的应用前景。但是，要得到准确度更高的有效储层厚度预测结果，还有待于进一步研究。