冷泉羽状流地震波场频谱特征与气含量关系研究

刘一林，李灿苹，勾丽敏，汪洪涛，曾宪军，陈凤英，郭子豪，田鑫裕

(1.广东海洋大学 电子与信息工程学院，广东 湛江 524088；2.中国地质大学(北京)海洋学院，北京 100083；3.中油测井天津分公司解释评价中心，天津 300457；4.广州海洋地质调查局，广东 广州 510075)

0 引 言

天然气水合物，俗称可燃冰，存在于海底沉积层或陆地永久冻土带，具有能量高、分布广、规模大、清洁等特点，被认为是21世纪的重要新能源[1-3]。因此，天然气水合物具有重要能源战略意义，受到世界各国的密切关注[4]。近年来中国在天然气水合物的勘探和开采方面开展了大量的工作，并且取得了显著的成绩，不仅在海域和陆地上都成功钻获实物样品，还在试采上连续取得重大突破[5-6]。2017年5月我国在南海神狐海域首次成功试采天然气水合物[7]，2020年2月第二轮海上作业试采点火成功[8]。两次成功试采标志着我国天然气水合物实现了从“探索性试采”向“试验性试采”的重大跨越，迈出天然气水合物产业化进程中极其关键的一步[9]。

在水合物赋存区域经常发现一种海洋地质现象，即大量海底活动冷泉[10-11]，它是继现代海底热泉活动之后又一新的海洋地质研究领域[12-13]。据估计，全球范围内发育的冷泉活动区有九百多处[14]，墨西哥湾北部布什山(Bush Hill)是冷泉渗漏活动发育的典型地区[15]，卡斯卡迪亚水合物脊(Hydrate Ridge)[16]、新西兰北岛西库朗伊边缘海[17]、尼日尔三角洲边缘[18]等均发现有大量的活动冷泉。我国南海海域、东海海域也有大量活动冷泉发育，其中南海海域的冷泉研究较为成熟，如海马冷泉、台西南冷泉。当冷泉中大量甲烷气体喷涌到海水中，便形成了海底冷泉羽状流。所以，冷泉羽状流是海底气体渗漏的直接证据，同时也对水合物勘探识别起到间接指示作用。研究发现，世界范围内发育的活动冷泉均有高浓度甲烷溢出，其浓度值高出海水背景值数百倍甚至数千倍。如新西兰希库朗伊陆架边缘的Tui冷泉渗漏出的甲烷浓度经常超过500 nmol/L[19],黑海近海底海水甲烷浓度为 50～5.5 μmol/L[20],美国圣塔芭芭拉海峡的甲烷羽状流中甲烷浓度高达1 600 nmol/L[21]。我国南海北部陆坡冷泉出露区的近海底甲烷浓度达3.8～4.2 μmol/L，是上覆水体浓度和海水背景值的300～1 000倍[22]。因而，冷泉羽状流对海洋生态环境以及大气环境也会产生重要影响。

冷泉羽状流的识别可以通过海底可视技术、声呐探测和地震探测技术实现，海底可视技术可以捕捉到羽状流气泡上升的真实照片[23]；声呐技术由于探测频率高从而可以获得羽状流清晰图像[12]。然而随着地震海洋学的发展，地震探测方法不仅可以识别海水中羽状流，还能够将海水与下伏地层连接起来，从而分析羽状流、冷泉和水合物的内在联系。近几年，国内众多学者在羽状流数值模拟方面进行了深入研究，如李灿苹等[2,24-25]通过建立羽状流水体模型，研究了冷泉活动区气泡羽状流的地震响应特征；段沛然等[26]使用交错网格有限差分法进行了羽状流数值模拟并观测其地震响应特征，其正演结果表明地震响应能够准确描述海底冷泉羽流；张闪闪等[27]采用含气泡液体声波方程进行了海底冷泉高频地震波数值模拟，实现了冷泉羽状流地震响应的高精度数值模拟。在羽状流地震探测及与声呐探测相结合方面也取得了一定进展，如刘斌和刘胜旋[28]利用多波束声呐发现了琼东南海域的羽状流，并结合浅层剖面分析了气体渗漏与水合物系统之间的相互作用；韩同刚等[29]总结了羽状流在海底可视技术、声呐系统、地震方法的表现形态特征，分析了三种探测方法的适用性和局限性；杨力等[30]利用多波束声呐数据、多道地震数据以及底质取样结果研究了琼东南海域活动冷泉系统，分析了活动冷泉的羽状流特征、海底地貌、底质特征以及流体活动构造特征。

与声呐相比，地震方法除了可以识别冷泉羽状流，还可以通过地震属性分析羽状流地震响应特征，由地震属性反演羽状流气含量，从而进一步探讨气含量与海底地层中水合物赋存状态之间存在的内在机理。本文将利用我国南海某测区羽状流地震偏移剖面数据进行地震波场频谱特征分析，通过已有气含量与振幅属性关系模型反演气含量，最终获得羽状流地震信号频谱参量与气含量的关系。本研究为羽状流的地震探测以及气含量反演的实际应用奠定基础。

1 地震信号频谱分析及频谱特征参量

信号频谱分析，就是计算信号的傅立叶变换，在频域中对信号进行处理，达到预期要求[31]。由于地震波的频谱特征包含了地下地层的岩性及构造特点信息，所以研究地震波的频谱特征有众多意义，比如，分析有效波和干扰波在频谱上的差异，可有效地指导地震仪器设计、地震资料处理参数的选择；利用频谱特征可以对地震资料进行地层岩性解释和特殊岩体解释等[32]。羽状流地震信号是地震波在海水介质中传播时遇到大量气泡产生的散射波，因此，研究羽状流地震信号频谱特征有助于进一步分析海水中气泡的分布特点以及气泡含量的分布状态。

设离散地震信号x(nΔt)={x(0·Δt),x(Δt),x(2Δt), ……,x((N-1)·Δt)}，通过离散傅立叶变换，其频谱和频率的计算公式[33]如下：

(1)

(2)

式中：x(nΔt)为时域地震信号；Δt为信号的时间采样间隔，ms；m为离散采样点；N为采样长度；fm为第m个样点的频率；X(fm)为频域频谱。

根据以上地震信号频谱和频率计算公式，本文提取以下频谱特征参量[34]：

(1)主频，频谱曲线的最大峰值所对应的频率，表示地震脉冲中能量最大的简谐分量。

(2)振幅谱最大值，主频对应的振幅值。

(3)优势频宽，谱值超过一定门坎值T的频带宽度，对归一化振幅谱而言，T=0.5。

(4)有限带宽能量，给定频率范围(优势频宽)内能量，计算公式为：

(3)

式中：X(m)为频谱；A1为有限带宽能量；m1为起始样点；m2为终止样点。

(5)总能量，频域振幅的平方和，它与地层的反射系数等因素有关，在含油气范围内，波形总能量一般响应较强。在所选时窗内其计算公式为：

(4)

式中：X(m)为频谱；A2为有限带宽能量；m为样点位置；N为样点长度。

2 羽状流地震数据来源

本文所用地震数据采集于南海北部琼东南海域[35]，其位置如图1所示。该海域位于海南岛东南向的新生代被动大陆边缘型盆地[30]。在琼东南海域，不仅在地层中发现了天然气水合物，还通过ROV(遥控无人潜水器)在海底发现了大面积冷泉系统(海马冷泉)，同时结合多波束影像数据观测到多个气泡羽状流[28]。

地震方法不仅是探测天然气水合物的有效手段，还可以探测到水体中的气泡羽状流[35-36]。对研究区(图1黑色矩形框区域)的水体地震数据进行处理，其主要流程包括：滤波，去异常值，振幅补偿，切除海底，去除次一级噪声，速度分析，炮域去除水平轴，切除直达波，去除类似海底反射，第二次振幅补偿以及叠前时间偏移。经过以上步骤的处理，最终获得水体中气泡羽状流的地震偏移剖面，如图2所示。

本文选取图2中一部分地震数据进行地震信号频谱分析，具体数据如下：CDP2201—3100，共900道，时间上从600～1 800 ms的地震数据，采样间隔为2 ms，即600个样点，如图2中蓝色框区域。从图2中可以看出，该区域上部为羽状流气泡稀疏区，下部为羽状流气泡密集区，说明整个区域从上到下(对应深度增加)羽状流气含量越来越高。为了分析羽状流气含量的变化是否会引起地震信号频谱特征参量的改变，将所选蓝色框区域分为18个模块(图2中红色虚线区域)，每一模块都是100×300的二维数据体，18个模块所提取的地震数据依次命名为plume_mn，其中m取值1、2、3、4、5、6，n取值1、2、3。为便于研究，首先抽取整个剖面中CDP2550和CDP2750两道(图2中黄色虚线)数据分析该剖面地震信号时域波形特征以及频谱特征；然后再抽取每个模块中第100道和第200道数据提取频谱参量和反演气含量，尝试分析频谱特征参量与气含量是否相关。最后对每个模块的所有道地震数据进行频谱分析，从而获得羽状流气含量与地震信号频谱特征参量之间的关系。

3 羽状流气含量反演

根据第1节地震信号频谱参量提取方法，利用图2中的实际羽状流地震数据，即可获得羽状流地震信号频谱参量。为了探讨羽状流气含量是否与地震信号频谱特征参量相关，还需要利用实际羽状流地震数据反演出气含量。根据李灿苹等人研究[36]，可知羽状流气含量与均方根振幅存在如下关系模型：

y=0.048x+7.99×10-4

(5)

式中：x是气含量；y是均方根振幅。均方根振幅可以通过如下公式计算获得：

(6)

即，在所选时窗内计算起始样点n1和终止样点n2内地震信号x(nΔt)的均方根A3，其中Δt是采样间隔，n为样点，n2-n1为样点总数。

因此，利用实际羽状流地震数据先由式(6)计算出均方根振幅，然后再由式(5)即可反演出羽状流气含量。当利用单道地震数据时，将反演出的气含量一维数组进行平均后作为该模块该地震道的气含量；当利用多道地震数据(每个模块的所有数据)时，将反演出的气含量二维数组进行平均后作为该模块的气含量。

4 羽状流地震信号频谱参量与气含量的关系

4.1 羽状流单道地震数据

首先抽取图2中整个剖面的CDP2550和CDP2750两道数据分析该剖面地震信号时域波形特征以及频谱特征，然后再抽取每个模块中第100道和第200道数据提取频谱参量和反演气含量，尝试分析频谱特征参量与气含量是否相关。

4.1.1 单道羽状流时域波形及频谱特征

选取图2中整个剖面的CDP2550道和CDP2750道地震数据进行时域波形特征分析及其频谱特征分析，其时域波形及频谱如图3所示，两道地震信号所在深度相同，只是横向位置不同。

从图3中可以看出，两道信号时域上都呈现出散射波动特征，特别是0.8 s之后，散射波动愈加强烈，这可能是由于在0.8 s之后所对应深度区域存在更多的羽状流气泡。含气泡的羽状流水体介质属于气-液双相介质，也是随机介质，作为散射体的气泡，其半径相对于地震波长来说很小，所以会产生散射波动现象。相比之下，CDP2750道信号在深部振幅稍强，说明此处气泡分布较密集或者气泡半径稍大，故波动幅度较大。两道信号的频谱带宽基本相同，都在250 Hz以内；由于是散射波频谱，信号频率成分较多，所以开叉较多，但在主频两边频谱分布较对称(尤其是CDP2550道)。CDP2550道频谱主要能量集中在60～200 Hz以内，突出的尖峰较少；CDP2750道频域振幅也稍强，与时域相对应，频谱主要能量集中在40～210 Hz以内，在160 Hz附近具有较强单峰特征，说明该处羽状流气含量较大，产生振动较强。两道信号频谱在50～200 Hz范围内跳动剧烈，说明处于该位置的气泡遇到声波振动较强烈，而且是不同频率的振动。

4.1.2 单道信号各频谱参量与气含量的关系

为分析羽状流单道地震信号各频谱参量是否与气含量相关，首先提取plume_mn各模块中第100道和第200道地震数据；然后按照本文第3节中气含量的反演方法，先反演出每个模块中第100道和第200道的气含量；再利用本文第1节中提取频谱参量的方法提取每个模块中第100道和第200道的频谱参量，即总能量、有限带宽能量、振幅谱最大值、主频和优势频宽；最后将每个模块中单道平均气含量作为横坐标，每个模块中单道各频谱参量值作为纵坐标，绘制出散点图，并添加趋势线，最终结果如图4和图5所示。

4.1.2.1 第100道羽状流地震信号各频谱参量与气含量相关性分析

从图4可以看出，气含量改变时会引起各频谱参量的变化，各振幅类频谱参量(总能量、有限带宽能量和振幅谱最大值)与气含量都呈现出较明显的正相关关系，即随着气含量的增大，各振幅类频谱参量也随之明显增大，其中总能量和有限带宽能量与气含量的相关程度较高，散点基本分布在拟合线附近。主频和气含量也呈现出正相关的趋势，只是散点分布较零散，主频在40～160 Hz之间，主要集中在80～150 Hz之间。优势频宽分布零较散，总体呈下降趋势，即随着气含量的增大，优势频宽减小。这说明当气含量增大时，气泡的振动愈发强烈，但气泡振动的频率更加集中，频带变窄。图4中优势频宽出现一个零值，这是由于该模块中第100道归一化后的振幅谱在门槛值0.5之上的只有一个点，所以优势频宽结果为0。

4.1.2.2 第200道羽状流地震信号各频谱参量与气含量相关性分析

从图5可以看出，气含量改变时同样会引起各频谱参量的变化，但气含量的变化范围较第100道的略大。与第100道类似，各振幅类频谱参量(总能量、有限带宽能量和振幅谱最大值)与气含量也都呈现出较明显的正相关关系，即随着气含量的增大，各振幅类频谱参量也随之明显增大，其中总能量和有限带宽能量与气含量的相关程度也较高，散点基本分布在拟合线附近。主频和气含量的关系也与第100道类似，呈现出正相关的趋势，只是散点分布较零散，但主频范围缩小，在60～160 Hz之间，也主要集中在80～150 Hz之间。优势频宽分布仍然较零散，但与气含量的关系总体略呈上升趋势，这说明当气含量增大时，气泡的振动愈发强烈，但气泡振动的频率范围略微加宽。优势频宽与气含量的关系之所以与第100道有差异，这是由于每道信号强弱不同，不同位置气泡的密集程度以及气泡大小均有所不同。图5中优势频宽也出现一个0值，同样也是由于该模块中第200道归一化后的振幅谱在门槛值0.5之上的只有一个点，故优势频宽结果为0。

4.2 羽状流多道地震数据

从单道地震数据分析可知，气含量的改变会引起各频谱参量的变化，而且气含量与振幅类频谱参量呈现出较好的正相关关系，频率类参量与气含量的相关性不是很明显，同时两道的优势频宽也有差异。为进一步验证以上结果，下面增加统计样本数量，利用18个模块中所有道地震数据进行频谱参量特征分析。

4.2.1 多道羽状流时域波形及频谱特征

选取模块plume_42中所有道地震数据，其时域波形用三维图显示，如图6(a)所示，x、y坐标分别为CDP号和双程旅行时，z坐标为时域振幅。对plume_42中二维地震数据进行二维傅立叶变换，幅度频谱用三维图显示，如图6(b)所示，x、y坐标分别为CDP号和频率，z坐标为频域振幅。为方便显示，此处CDP号为plume_42模块中的重新标号。图6(a)的情况与单道相似，时域上仍然都呈现出散射波动特征，振幅在±0.10范围变化，说明该模块所对应的区域内含有大量气泡，当地震波入射到含有气泡的海水介质时产生了散射波动现象。图6(b)的三维频谱图显示出plume_42中所有道信号的幅度频谱，其中CDP1—CDP150的幅度较小，CDP151—CDP280的幅度较大；同单道类似，各道频谱带宽基本相同，都在250 Hz以内，呈现开叉较多的散射波频谱特征，说明不同区域气泡遇到声波发生不同频率的振动。

4.2.2 多道信号各频谱参量与气含量的关系

为利用多道信号进一步分析羽状流地震信号各频谱参量与气含量的相关性，首先提取plume_mn中每个模块所有道地震数据；然后按照本文第3节中气含量的反演方法，反演出每个模块的平均气含量；再依据本文第1节中提取频谱参量的方法提取每个模块的各频谱参量，即总能量、有限带宽能量、振幅谱最大值、主频和优势频宽；最后将每个模块的平均气含量作为横坐标，每个模块的各频谱参量作为纵坐标，绘制出散点图，并添加趋势线，最终结果如图7所示。

从图7可以看出，平均气含量的变化范围与单道情况大致相似，只是最低值略大；气含量改变时同样会引起各频谱参量的变化，而且各振幅类频谱参量(总能量、有限带宽能量和振幅谱最大值)与气含量都呈现出更加明显的正相关关系，即随着气含量的增大，各振幅类频谱参量也随之明显增大，其中总能量和有限带宽能量与气含量的相关程度很高，散点基本分布在拟合线附近；主频和气含量的关系与单道情况类似，也呈现出正相关的趋势，只是散点分布较零散，主频在60～150 Hz之间，主要集中在90～140 Hz之间；优势频宽分布较零散，总体呈下降趋势，与第100道情况类似，即随着气含量的增大，优势频宽减小，这说明当气含量增大时，气泡的振动愈发强烈，但气泡振动的频率更加集中，频带变窄。

羽状流地震信号频谱参量(除优势频宽)与气含量呈现出如此明显的正相关关系，足以说明甲烷气泡羽状流首先能够引起地震响应，然后当羽状流气含量发生变化时，地震属性也将发生变化。羽状流频谱参量(除优势频宽)与气含量之间为正相关关系，这是由于气含量增大，气泡半径随之增大(气含量为体积含量)，或者气泡个数增多，使得地震波发生散射的散射面增大，散射波的能量也随之增大，从而导致振幅增强，而且气泡分布密度较高部位散射相干加强。因此，通过羽状流地震信号频谱参量与气含量的相关性分析可以推测海水中气泡的分布特点以及气泡含量的分布状态。

5 结 论

本文利用实际羽状流地震偏移剖面中单道和多道地震数据，首先进行了地震波时域波形特征分析和地震波场频谱特征分析，然后通过已有气含量与振幅属性关系模型反演气含量，再提取羽状流地震频谱参量，最后获得羽状流地震信号频谱参量与气含量的关系，具体认识如下：

(1)羽状流时域波形呈现出散射波动特征，不同位置上由于气泡含量的差异导致了振幅强弱变化；羽状流频谱带宽基本在250Hz左右，频谱开叉较多，具有单峰特征，说明不同区域气泡遇到声波发生不同频率的振动。

(2)综合单道和多道地震信号频谱参量与气含量的相关性分析表明：气含量的改变可引起各频谱参量的变化；各频谱参量(除优势频宽)与气含量都呈现出较明显的正相关关系，其中总能量和有限带宽能量与气含量的相关性最好，主频和气含量呈现出正相关的趋势；单道的优势频宽大致呈上升和下降两种趋势，而多道优势频宽则分布较散，总体呈负相关趋势。

(3)羽状流频谱参量(除优势频宽)与气含量呈现出明显的正相关关系，说明甲烷气泡羽状流能够引起地震响应，并且当羽状流气含量发生变化时，地震属性也将发生变化。导致这种结果的原因是，气含量增大时地震波在阻抗界面受到的散射作用更强，则散射波的能量也将随之增大，从而形成能量更强的信号以及更高的频谱参量值。

综上，羽状流地震频谱参量与气含量呈现较好的相关关系，为羽状流的气含量反演做了前期铺垫。羽状流气含量的分布特征及其变化状态与海底气体来源及地质构造直接相关，而水合物分解的甲烷是羽状流气体的重要来源。所以，可以将水体中羽状流气含量地震剖面与海底水合物地震剖面连接起来，分析羽状流赋存状态、海水下伏地层岩性和构造特征、水合物赋存环境特征以及流体运移通道特征等，从而探索羽状流与水合物之间的内在联系，寻找识别水合物的地球物理特征标志，为水合物的勘探识别提供有意义的方法指导。