0 引言

二叠系茅口组是四川盆地最早的油气勘探层系之一，其勘探历程表现出多类型、多领域的特征，是四川盆地油气勘探的热门层系^[1-3]。20世纪50年代至今，以高陡构造裂缝储集体、岩溶缝洞储集体、热液白云岩为目标开展了多次勘探^[4-5]。胡东风^[6]提出四川盆地元坝地区茅口组三段为缓坡型台地边缘相，发育沉积能量较高的浅滩，在东吴构造运动时期抬升并受岩溶作用改造，形成相控岩溶储层^[7]。

四川盆地元坝地区茅口组三段相控岩溶储层成因机理复杂，受沉积环境和流体作用叠加控制，储层表现为非均质性强、结构复杂、溶蚀孔洞充填严重等特征。该地区地震响应特征复杂，岩溶地震异常差异大，导致岩溶异常带刻画多解性强，有效储层定量预测难度大。以往学者在相控岩溶发育模式及储层特征的基础上，对四川盆地茅口组复杂岩溶目标描述及定量预测等难题进行了研究。李源等^[8]通过研究子波分解重构技术增强了岩溶地震异常；徐敏等^[9]、赵虎等^[10]利用分频相干、曲率等地震结构类属性刻画出岩溶异常带分布^[9-10]；熊晓军等^[11]通过研究相控随机反演获得了岩溶低速异常体，并定量预测了岩溶储层的分布。元坝地区茅口组岩溶纵向上表现为结构特征复杂、隐蔽性强，应用上述岩溶刻画及预测技术描述的岩溶储层匹配性差，精度低。基于此，以岩溶发育模式为基础，通过改进基于倾角约束的增强相干技术来降低元坝地区茅口组相控岩溶地震异常刻画的多解性，结合地震结构属性研究岩溶纵、横向结构变化的相控岩溶有效储层的反演方法，以期为该区茅口组岩溶目标描述及勘探部署提供借鉴。

1 地质概况

元坝地区位于四川省苍溪县，地处四川盆地北部，构造上处于川北低缓构造带与川中平缓构造带结合部，其西北部与九龙山背斜构造带相接，东北部与通南巴背斜构造带相邻，南部与元坝低缓褶皱带相连（图 1a）。元坝地区茅口组早中期（茅一段、茅二段）以缓坡沉积为主，水体由深变浅，茅一段地层厚度为90~100 m，岩性为泥质灰岩、泥晶灰岩；茅二段地层厚度为70~90 m，水体变浅，生屑灰岩含量增加，岩性为泥晶灰岩、生屑灰岩；目的层段茅三段沉积时期，受峨眉地裂运动控制，开江—梁平陆棚初具雏形，元坝地区逐渐演化为缓坡型台地沉积，地层厚度为80~160 m，岩性为泥晶灰岩、生屑灰岩夹薄层灰质云岩（图 1b）。

2 岩溶储层发育模式及特征

元坝地区茅口组三段处于台缘浅滩沉积相，为相对高能的沉积环境，受东吴运动时期构造抬升大气淡水淋滤和溶蚀作用，形成相控岩溶储层。早期浅滩沉积是储层发育的基础，后期岩溶改造是储层形成的关键。

淡永等^[12]依据现代岩溶理论分析钻井岩溶特征并建立岩溶发育模式，纵向上把茅口组划分为表层风化壳、垂直渗流带、水平潜流带3个岩溶作用带。其中，表层风化壳位于茅三段顶部，距茅口组顶20 m处，以片流和散流为主，溶蚀作用强，发育孔隙-孔洞型储层；垂直渗流带位于茅三段中上部，距茅口组顶20~60 m处，受重力及流体势能控制，以垂向渗流为主，溶蚀作用较强，发育孔隙-孔洞型和半充填溶洞型储层；水平潜流带位于茅三段中部，距茅口组顶60~100 m处，受底部致密隔挡层影响，以径向管流为主，溶蚀作用变弱，局部发育半充填溶洞型储层（图 2a）。综合钻井、地震及岩溶地貌等资料，依据水动力强弱、流体势、饱和度等的差异，把茅三段上亚段划分为岩溶残丘和溶丘平原2个三级地貌单元^[13]。其中，岩溶残丘包括垄脊、峰丛、坡地等3个微地貌单元；溶丘平原包括槽谷、丘丛2个微地貌单元（表 1）。分析研究区过YB7井—YB701井—YB8井的岩溶发育模式，认为垄脊和峰丛区储层发育较好。垄脊区（YB7井，YB701井）顶部发育孔隙-孔洞型储层，中下部发育半充填溶洞型储层；峰丛区主要发育半充填溶洞型储层；坡地区（YB8井）顶部发育薄层孔隙-孔洞型储层，中下部半充填溶洞储层较发育（图 2b）。

元坝地区茅三段相控岩溶储层受浅滩相带和岩溶作用共同控制，储层类型主要为孔隙-孔洞型和半充填溶洞型。其中，孔隙-孔洞型储层发育于茅三段顶部，岩性主要为云质生屑灰岩、灰质白云岩、细—中晶白云岩；半充填溶洞型储层发育于茅三段中上部，岩性主要为生屑灰岩，具泥晶生屑结构，生屑主要包括棘屑、有孔虫和腕足类。茅口组储层储集空间类型多样，包括孔隙、溶洞和裂缝三大类。其中，孔隙包括晶间孔、晶间溶孔，生屑溶孔、粒间溶孔；溶洞主要为次生溶洞，形态极不规则，洞径为0.2~5.0 cm，溶洞内充填物为方解石、硅质、黄铁矿。

3 岩溶有利带精细刻画

四川盆地二叠系茅口组晚期岩溶作用时间短，岩溶的规模和尺度较塔里木盆地奥陶系小。元坝地区茅口组储层受沉积作用与岩溶作用的共同控制，岩溶非均质性强、结构复杂，地震异常识别及刻画难度较大。

3.1 岩溶古地貌研究 3.1.1 岩溶古地貌恢复方法

元坝地区茅三段晚期岩溶作用主要受东吴运动影响。由于东吴运动以整体抬升为主，岩溶古地貌主要受茅口组晚期沉积地貌及剥蚀厚度的共同影响，利用常规的古地貌恢复方法（如印模法、残厚法）难以准确恢复岩溶古地貌^[14-15]。综合考虑沉积地貌后形成的地表趋势面和岩溶作用后的微古地貌，提出基于趋势面和微地貌的岩溶古地貌恢复方法。该方法恢复茅口组晚期岩溶古地貌可分为3步：①通过低通滤波计算茅三段残余地层厚度的宏观地貌趋势面；②对印模地层厚度进行高通滤波，计算出岩溶作用后的微地貌信息；③融合宏观地貌趋势面和岩溶微地貌信息，恢复出较准确的岩溶古地貌（图 3）。

3.1.2 岩溶古地貌单元划分

元坝地区茅口组岩溶作用时期主要发育地表水系^[16]，地表水由YB7—YB701—YB8井区高地向南部和北部散流，地下水经潜流带底界汇聚运移至岩溶盆地。YB7—YB701—YB8井区为台地边缘沉积相，整体地貌形态相对较高，微地貌形态及信息丰富，能表征岩溶作用后的地貌差异特征，为该区茅口组岩溶发育模式及岩溶异常精细刻画奠定了基础（图 4）。

受岩溶古地貌高差的影响，各岩溶地貌单元发育的岩溶储层类型不同，表现的地震响应特征也具有明显差异。岩溶残丘地貌单元相对较高，为岩溶作用供水区域；溶丘平原地貌单元相对较低，为岩溶作用汇水区域（图 5）。

3.2 地震精细刻画

元坝地区茅口组岩溶地貌单元类型多样，岩溶结构变化快，地震异常响应复杂，因此利用常规的相干、曲率及结构张量等属性刻画岩溶发育带的精度较低。研究区位于川北九龙山构造倾末端，区内发育低缓褶皱构造，影响了地震结构类属性的表征能力。

3.2.1 增强相干岩溶异常刻画技术

相干是求取相邻地震道之间的相似系数，再组成相干数据体。计算过程中主要受分析窗口的空间组合方式、时间长度及相干时移等因素影响。目前，最新的相干技术为第三代相干，是利用地震道组成协方差矩阵的特征值来计算相干值。该技术在分辨率及噪声压制方面的应用效果较好^[17-18]，但其并未考虑地层倾角的影响。

在相干分析窗口中考虑地层倾角的变化，沿倾角提取相邻地震道的数据，再利用增强相干算法求解相干体，形成基于倾角约束的增强相干算法。首先，沿三维地震数据主测线、联络线方向求取目标点位的视倾角（p，q），其中，p代表主测线方向视倾角，q代表联络线方向视倾角；其次，通过相干分析窗口计算各地震道的协方差矩阵：

$ M=\left[\begin{array}{cccc} \sum\nolimits_{n=1}^N u_{n 1}^2 & \sum\nolimits_{n=1}^N u_{n 1} u_{n 2} & \cdots & \sum\nolimits_{n=1}^N u_{n 1} u_{n J} \\ \sum\nolimits_{n=1}^N u_{n 2} u_{n 1} & \sum\nolimits_{n=1}^N u_{n 2}^2 & \cdots & \sum\nolimits_{n=1}^N u_{n 2} u_{n J} \\ \cdots & \cdots & \cdots & \cdots \\ \sum\nolimits_{n=1}^N u_{n J} u_{n 1} & \sum\nolimits_{n=1}^N u_{n J} u_{n 2} & \cdots & \sum\nolimits_{n=1}^N u_{n J}^2 \end{array}\right] $

(1)

式中：u_nJ为相干分析时窗内第J道的第n个样点元素。

对式（2）进行特征值分解得到J个特征值λ₁，λ₂，…，λ_J，对应的特征向量为V₁，V₂，……，V_J。

最后，以目标点为中心沿视倾角（p，q）计算相干分析窗口所有地震道的增强相干：

$ C_3(p, q)=\frac{\lambda_1}{\sum\limits_{j=1}^J \lambda_j} $

(2)

在增强相干算法中引入地层倾角的约束，其相干分析窗口中更接近实际地层，刻画的岩溶地震异常更加精细。

过元坝地区YB7井茅口组距顶界面70 m范围内，岩溶储层较发育，地震同相轴横向变化快，振幅能量差异明显，表现出岩溶地震异常的差异（图 6a）。由倾角约束增强相干（图 6b）方法计算出YB7井区茅口组的岩溶异常分布，比利用常规增强相干方法计算的结果（图 6c）更清晰，新方法更准确地表征了岩溶异常，该异常强弱可反映岩溶储层的分布。

3.2.2 岩溶有利带地震精细描述

利用倾角约束增强相干岩溶刻画技术计算岩溶异常体，可表征研究区茅口组复杂的岩溶异常展布特征。岩溶残丘地貌单元内的岩溶异常相对发育，主要集中在垄脊和峰丛地貌单元内，且峰丛地貌单元内的岩溶异常强度大，说明岩溶作用改造较强，为岩溶异常有利区；溶丘平原地貌单元内局部地貌高差大的位置有较明显的丘丛带岩溶异常，其中西北部丘丛岩溶异常相对发育，为丘丛带岩溶有利区（图 7）。

4 相控岩溶有效储层地震预测

常规储层预测方法是基于层状介质储层，其无法适用于元坝地区茅口组纵横向结构复杂的岩溶储层预测。为提高岩溶储层预测结果的合理性，基于岩溶结构类属性，对相控岩溶地震反演方法和非线性伽马模拟方法及适应于相控岩溶的有效储层预测技术进行了研究。

4.1 相控岩溶地震反演技术

勘探实践证实，在复杂目标储层反演中，初始模型对反演结果的影响较大。构建常规反演初始模型的方法为确定性插值，横向以反距离加权、自然邻域、克里金等算法对钻井信息进行插值^[19-21]，但数学插值模型化强，不符合强非均质储层发育规律。依据茅口组岩溶储层发育模式，以表征半充填溶洞储层复杂结构的地震结构类属性（倾角约束增强相干、多尺度边缘检测及结构张量）作为约束，建立了岩溶属性构型建模方法，构建符合岩溶储层发育规律的初始模型（图 8）。

$ Z^*\left(u_0\right)=\sum\nolimits_{i=1}^n \alpha_i Z\left(u_i\right)+\beta_1 Y_1\left(u_0\right)+\beta_2 Y_2\left(u_0\right) $

(3)

式中：α_i，β₁和β₂为权重系数，Y₁（u₀），Y₂（u₀）为不同的岩溶结构类属性。以协同变量把Y₁，Y₂融入到克里金插值算法中，构建岩溶相控初始模型。

在相控岩溶建模的基础上，为准确预测岩溶低速异常带，以贝叶斯框架联合似然函数分布和先验分布得到后验概率分布，构建反演目标函数，形成高分辨率岩溶相控地震反演技术。该方法在优选样本点时，利用马尔科夫链蒙特卡洛随机模拟进行最优估计^[22-23]，可使反演结果在横向上较好地表征岩溶异常非均质性。样本点优选模拟过程为

$ z\left(x_0\right)=\sum\nolimits_{i=1}^n \omega_i z\left(x_i\right) $

(4)

式中：z（x₀）为未知点的值；z（x_i）为波形优选已知样本点的值；ω_i为第i个已知样本点对未知样本点的权重，n为优选样本数，个。

由相控岩溶地震反演方法预测出元坝地区茅口组岩溶低速异常带的分布特征，在反演剖面中，低速异常表现出上、中、下3层不同结构类型。其中，茅口组晚期上部为表层风化壳，表现为中低波阻抗值，主要分布在古地貌较高、沉积浅滩能量相对高的位置；中部为垂直渗流带，表现出纵向中—高波阻抗值，包括致密灰岩、方解石填充和部分泥质充填异常等，位于岩溶垂直渗流带内，岩溶异常相对不发育；下部为水平潜流带，以少量中—低波阻抗值的横向延伸储层为主，主要包括半充填溶洞储层（图 9a）。相控岩溶地震反演结果能客观反映岩溶纵横向变化的复杂结构，符合岩溶发育模式。

4.2 岩溶有效储层预测

在岩溶低速异常带基础上，为准确识别出岩溶有效储层和泥质充填范围，引入支持向量机（SVM）算法。基于地震振幅与频率的关系（AVF）形成非线性伽马模拟方法，预测出更符合岩溶复杂结构的伽马体，剔除高伽马泥岩充填异常，从而预测出岩溶有效储层的展布。

4.2.1 基于SVM的非线性伽马模拟方法

支持向量机（Support Vector Machine，SVM）是一种解决小样本、非线性函数的机器学习算法，该方法利用有限样本信息提高复杂非线性模型的学习精度，预测及推广能力强^[24]。相比神经网络方法，SVM具有较强的泛化能力，可有效克服过拟合缺陷。

支持向量机运用了回归估计（SVR）的思路：首先把输入数据通过非线性内积函数转换到高维空间，然后把转换后的输入数据在高维空间中进行线性回归并分析构造回归估计函数。假设训练样本集为$\left[\left(x_i, y_i\right), i=1, 2, \cdots, n\right] $，其中$ x_i \in \boldsymbol{R}^N$为N维样本输入值，y_i ∈ R为样本输出值，n为样本数。引入非线性支持向量机回归估计函数表达式：

$ f(x)=\omega \cdot \varphi\left(x_i\right)+b $

(5)

式中：ω·φ（x）表示向量ω和φ（x）的内积；ω的维数即为高维空间维数，b ∈ R为偏置量。

上述待求解问题属于凸二次问题，且特征空间维数较高，为满足Mercer条件，引入任意对称函数K(x_i, x_j)作为核函数代替非线性映射。依据Wolfe对偶定理，采用拉格朗日乘子法求解，将上述问题转换为对偶函数即可稳定求解。

与地球物理弹性参数相比，伽马反演的地球物理意义不够明确。因此，针对茅口组岩溶异常非均质特征，根据振幅随频率的变化信息（AVF）所反映的频散特征^[25]，研究基于SVM的非线性伽马模拟方法。其中，振幅与频率的关系在频率域泰勒展开，可建立非线性伽马模拟目标函数，利用支持向量机算法求取最优解，模拟出更准确的伽马展布规律。

$ r_0=R(t, f)={\mathbf{G}}[m(t, f)] $

(6)

式中：G为非线性算子；r₀为地震数据的矩阵形式；R（t，f）为反射系数的时频谱；t为时间，ms；f频率，Hz；m为模型参数（如波阻抗）。

由非线性伽马模拟结果可看出，研究区茅口组岩溶低速异常段顶部伽马值较高，泥质含量增多，中段泥质逐渐减少，非均质性强，符合岩溶低速异常发育规律（图 9b）。

4.2.2 岩溶有效储层预测及描述

在岩溶相控低速异常反演和非线性伽马模拟的基础上，高伽马岩溶异常体内充填了低速泥岩，以自然伽马大于78 API为门槛对岩溶低速异常带内的泥质充填段进行剔除（图 10），进而进行岩溶有效储层的预测。研究区茅口组岩溶作用的顶部风化壳岩溶低速异常范围广，但泥质充填严重，泥质剔除后残留的中低阻抗为有效储层（图 9c）；中段垂直渗流带岩溶作用以纵向为主，泥质充填减少，整体有效储层较少；水平潜流带岩溶作用减弱，泥质充填较少，有效储层欠发育。

在岩溶有效储层反演剖面中，以小于17 430 g/cm³·m/s的波阻抗值为门槛，预测茅三段岩溶有效储层厚度。研究区茅口组岩溶有效储层主要发育在岩溶残丘的垄脊带和峰丛带，有效储层最大厚度达28 m，平均厚度为13 m，预测结果符合该区茅口组岩溶发育模式及规律（图 11）。

5 结论

（1）元坝地区茅口组岩溶储层受古地貌控制，趋势面及微地貌的岩溶古地貌恢复和倾角约束的增强相干方法，减少了构造对岩溶地震异常信号的干扰，明确了茅口组晚期岩溶残丘地貌单元内垄脊和峰丛是相控岩溶储层发育的有利带。

（2）元坝地区茅口组岩溶储层结构复杂，非均质性强，溶蚀孔洞充填严重。基于岩溶结构类地震属性，建立了相控岩溶有效储层地震预测技术，定量预测了茅口组相控岩溶储层分布范围，认为垄脊带、峰丛带内储层厚度最大。

（3）元坝地区茅口组相控岩溶储层主要发育在古地貌较高位置，处于岩溶残丘单元内的垄脊带和峰丛带，呈北西向展布，储层厚度大，是研究区茅口组相控岩溶领域下一步勘探的有利目标。