引言

震动波在地层传播时，由于地层介质的热传导、黏滞性等因素而产生吸收衰减作用，导致震动波能量减小、振幅减弱的同时，相位产生畸变。对地震资料的后续处理会产生严重的影响。因此，对地震资料的高频部分进行合理有效地补偿是提高地震资料分辨率的主要有效手段之一（丁进杰等，2013）。

反Q滤波方法在地震资料处理和解释中有着重要的应用，该方法主要针对地层吸收衰减特性，对吸收衰减能量进行补偿和相位畸变的校正，可以提高地震资料的分辨率。反Q滤波有助于对薄互层识别、预测优良储层和合理解释振幅（Zhao et al，2019）。同时，震动波衰减是储层预测和流体识别的重要属性。作为表示介质衰减特性的重要参数，品质因子Q可以用来衡量黏弹性介质的黏滞性对震动波能量造成的衰减程度。因此，更高精确度地提取品质因子Q和对信号进行反Q滤波处理对于高分辨率地震勘探、深部信号增强和油气藏预测具有重要意义。

一般来说对实际地震资料衰减补偿时，反Q滤波一般用于补偿深部地震波的反射能量（Hale，1981）在Futterman模型的基础上，提出的大地滤波的逆运算方法，即早期的反Q滤波方法，用于补偿地震波的能量，但该方法的运行时间较长。裴江云和何樵登（1994）基于Kjartansson速度散射关系提出一种新的反Q滤波算法，对衰减的震源响应进行泰勒级数展开。白桦和李鲲鹏（1999）基于短时傅里叶变换（short-time Fourier transform）提出一种地震波能量衰减补偿方法，该方法增加了一个时窗函数，然而，短时傅里叶变换缺点是时窗宽度的频率是固定的，它的时频分辨率不能有效地适应地震信号的变化特征，最终也不能有效地提高地震资料的分辨率。为了克服短时傅里叶变换时窗的缺陷，Wang （2002）提出了一种基于波场延拓理论的稳定反Q滤波方法，并对该方法进行了改进（Wang，2006），扩展到地层Q值随着地层深度的不断变化的情况。姚振兴等（2003）结合在非弹性介质中地震波传播的理论提出了一种新的反Q滤波方法，有效地补偿了深部地震波的能量衰减。随着时频分析技术的发展，刘喜武等（2006）采用广义S变换时频分析技术来补偿地震波能量的衰减。Braga和Moraes （2013）基于对尺度与频率关系的研究，通过对一维连续小波的变换来补偿地震波的吸收和衰减。但它仍然存在一些缺陷，即其尺度宽度不能根据地震信号的特性灵活改变。通过研究截止频率法和稳定性因子法，张固澜等（2015）基于地震资料动态范围对地震资料分辨率的影响，提出了自适应增益限反Q滤波方法。王小杰和栾锡武（2017）基于小波分频技术提出了一种Q值补偿方法，在该方法中利用小波分频技术得到Q值，然后用反Q滤波法补偿衰减地震波的能量。赵岩等（2021）基于时频域信噪比提出了一种自适应增益限的反Q滤波方法，根据时频域内的信噪比阈值确定地震资料的有效频带。反Q滤波的幅值增益极限适应有效频带的截止频率，对不同时间的地震记录进行有限的可变增益幅值随时间的补偿。

品质因子Q是进行反Q滤波处理的基础，Q的计算准确性直接影响到反Q滤波处理的补偿效果，目前反Q滤波法中估算Q值方法的精度和稳定性欠佳。此外，目前的反Q滤波法大多应用于对叠后地震资料的衰减补偿，对叠前地震资料的补偿较少，而对叠前地震资料的衰减补偿可以更好地进行储层预测和高分辨率地震勘探。为此，本文拟采用基于震动波传播的褶积原理求出品质因子Q值，用于补偿衰减的震动波，并在改进广义S变换的基础上，对叠前地震资料进行反Q滤波处理，以期论证基于震动波的褶积原理的Q值提取与改进广义S变换结合的反Q滤波方法的优越性与实用性。

1.   方法原理

1.1   改进的广义S变换

由于高频部分对地层的吸收衰减作用反馈更为敏感，而改进广义S变换在高频部分具有较高的频率分辨率，可以满足吸收衰减补偿的要求。张固澜等（2010）改进广义S变换选用了时窗宽度与频率f成正比的可变宽度的高斯窗函数。改进广义S变换的时窗函数为：

式中：t为信号的传播时间；s和r为两个调节窗函数的因子，用于调控窗函数的衰减趋势和时间延续长度。

设时间序列为

式中，L²（R）表示实数域上的平方可积函数空间，则输入信号h（t）的NGST为：

式中，τ为窗函数的时间位置，i为虚数单位。

反变换表达式

1.2   褶积原理求取品质因子

以震动波传播过程是从震源的击发位置到传感器的接收位置的褶积过程的理论为基础（贾晓东等，2021），假设震源击发的位置在A位置，传感器位置为B，代表震动波受黏弹性介质影响产生衰减程度的品质因子Q根据震动波传播的原理，有

式中，Q⁻¹为震动波受地层吸收衰减作用的参数，*为褶积符号，A_i（t）和B_j（t）为传感器采集到的时域信号。等式两边同时进行傅里叶变换，得到

其中式（6）和式（7）均为傅里叶域上的表达式。

Q是经过Q⁻¹经过傅里叶反变换得到，取Q⁻¹的反算子，则

计算出品质因子Q可用于补偿震动波能量，消除由介质吸收引起的震动波能量衰减。

1.3   基于改进广义S变换的反Q滤波法

本文提出的基于褶积原理求取品质因子Q方法和改进广义S变换相结合的反Q滤波方法的基本思想是：假设震动波在地层介质中传播时，不会因为地下介质的衰减特性而产生能量损耗，深部反射波和地表反射波的振幅谱具有相同的带宽，它们的相位谱只差一个线性相位。假设将地震信号划分为不同的频率，它们之间对应于时间的振幅谱是相似的。然而，由于地层介质的黏弹性和不均匀性，不同时间段的震动波能量分布不同，因此可以通过地层乘以时变品质因子来补偿地震信号的衰减（高军等，1996；高静怀等，1996，2003）。

通过大地滤波后，当震源子波通过地下介质时，其能量随着传播时间的增加而逐渐减弱，这个衰减过程被称为Q滤波作用（熊晓军等，2006；何书梅等，2020）。假设地层反射系数为r_k （k＝1，2，···，k），由于子波时变，每一个反射系数都会经过不同的Q滤波，r_k对应时间t_k经Q滤波后的振幅响应为:

式中，p（0，ω）为初始时刻的等效振幅，Q_eq（t_k）为t_k处的等效Q值，ω为信号的数字频率，j为虚数单位。

由于介质吸收特性，需要对震动波能量的衰减进行补偿，信号表示如下：

式中：r_k为等效反射系数（r_k对应时间t_k，r₀对应初始时刻t₀），表示在没有能量衰减的情况下，振幅之间的相对关系；x（t_k，ω）为t_k时刻对应的地震记录。式（9）和式（10）均为时频域上的表达式。

将t_k时刻各频率的反射波振幅与初始时刻振幅的比称为衰减比率，表示如下：

用计算得到的品质因子Q对t_k时刻的信号进行加权处理，补偿后的信号表示为：

式中，$r'_k $为与常数合并后的等效反射系数，|x（t_k，ω）|′为补偿后的地震记录。

将t₀时刻的时频谱乘以等效反射系数即得到t_k时刻的时频谱，即恢复了衰减的震动波能量。因此我们只需对各时刻的时频谱进行加权处理，然后进行广义S反变换，最终得到补偿后的重构信号。

进行反Q滤波法补偿震动波的步骤如下：① 对传感器采集到的信号通过改进广义S变换进行时频分析；② 对变换后的信号基于褶积原理求取品质因子Q；③ 用品质因子Q加权相应时刻的信号；④ 进行广义S反变换得到重构信号。

2.   震动波衰减补偿模型试验

本文从理论研究的角度出发，试验模型模拟了实际工程中常见的地质结构，通过室内搭建理想模型的方法研究了基于褶积原理求取品质因子Q的准确性以及在改进广义S变换基础上对叠前地震资料衰减补偿的效果。

试验模型根据研究内容包含的地质条件进行搭建，利用多通道高频数字信号采集系统对试验数据进行采集。本次试验搭建了一个复杂模型（图1）。模型按介质不同分为三层，主要材料为石英砂、石灰、石膏和水，按预先计算好的材料用料进行配料，搅拌均匀后倒入模具并进行夯实和抹平。

图  1  室内试验模型

（a） 模型搭建；（b） 模型A面

Figure  1.  The model of indoor experiment

（a） Model building；（b）A-side of the model

下载: 全尺寸图片 幻灯片

根据团队已有试验资料成果和试配试验结果，最终确定各岩层的配比参数。其中第一层模拟细砂岩层，厚度为0.35 m，砂胶比为3 ∶ 1；第二层模拟煤层进行搭建，厚度为0.60 m，砂胶比为5 ∶ 1；第三层模拟泥岩层搭建，厚度为0.35 m，砂胶比为4 ∶ 1。各层原料的用量如表1所示。

表  1  相似材料用量表

Table  1.  The list of similar material consumption

岩层名称	模型厚度/cm	各相似材料用量/kg
		沙子	石灰	石膏	水
泥岩层	35	963.2	160.6	160.6	183.4
煤层	60	1694.0	101.7	237.3	290.4
细砂岩层	35	1133.3	149.9	133.6	202.4

下载: 导出CSV 

| 显示表格

模型的尺寸为2.30 m×1.00 m×1.30 m，模型长边和高对应的两个侧面即为A面和B面。其中A面左侧面记做C面，B面右侧的面记为D面，模型顶面记为E面。

本次试验采用的多通道高频数字信号采集系统进行监测和数据采集工作，传感器的命名方法为序号排列而成，传感器记为“1”，“2”，“3”···。传感器布置位置如图2所示。试验采用人工击发小球模拟震源，在预设位置（图2）手动放置钢球垂直撞击试验模型表面，详细记录震源位置及各传感器接收信号的情况。

图  2  传感器布置图

Figure  2.  Sensor layout diagram

下载: 全尺寸图片 幻灯片

3.   震动波衰减补偿

3.1   时频特性分析

本试验采用的多通道高频数字信号采集系统的采样频率为50 kHz，选取1，4号传感器采集到Z₁震源的信号，利用改进广义S变换对传感器采集的信号进行时频特性分析。

根据图3中的原始信号，可以看出：信号持续时间较短，振动速度衰减较快。经过改进广义S变换后的原始信号时频谱如图4所示，从频率分布范围来看，人工击发小球模拟震源的信号包含多个频率成分，频率主要集中在0.5—1 kHz。1号传感器信号在0.01 s附近时频率较大、具有的能量较高；4号传感器信号在0.1 s附近时频率较大、具有的能量较高。同时，原始信号的能量损耗随着传播时间而增加，导致信号的主频减弱，有效频段变窄，其中夹杂着由于实验过程中的随机噪音干扰造成的个别不规则抖动，地震信号的改进广义S变换时频分辨率与原始信号的对应关系很准确，可以标定地震信号在不同时刻频率的变化状况。

图  3  1号（a）和4号（b）传感器接收的原始信号

Figure  3.  Original microseism signals of sensor 1 （a） and 4 （b）

下载: 全尺寸图片 幻灯片

图  4  1号（a）和4号（b）传感器接收到的原始信号经过改进广义S变换后的时频谱

Figure  4.  Time-frequency spectrum of the original signals received by sensor 1 （a） and 4 （b） after the improved generalized S transform

下载: 全尺寸图片 幻灯片

通过分析原始信号经过改进广义S变换前后的时频谱如图4所示，由图4a可以看出，频率范围分布广，高频成分明显，并且图中清晰地展示出在400 Hz附近的能量集中区，同时400 Hz的信号存在时间较长，约为0.025 s，信号在0.01 s附近时频率较大，可以达到650 Hz；图4b中能量集中区位于700—750 Hz范围内，震动波在0.09 s附近的频率较大，可以达到1250 Hz，频率为650 Hz的信号存在时间较长，约为0.13 s。通过分析改进广义S变换后的时频图（图4）可以明显的看出地震信号模型的时频分布，时频谱中能量最大的区域为颜色较深的区域，对应的时间与时域信号波峰中心出现的时间也是一致的，可以从时频谱中读取到有效信息，同时更好地保留原始信号的特征。

改进广义S变换其时窗函数可以随着频率变化而灵活改变，能区分某些低频信号，具有较高的时频聚集性，信号时间和频率对应关系准确，较好地展示了各频段的地震信号。同时，改进广义S变换不仅可以描述信号随时间变化的频率特征，而且可以分离同一时间不同的频率分量，对信号频率特征进行定量研究，丰富了地震信号波形分析方法。这些优势极大地提高了地震资料的分辨率，对于准确提取品质因子Q和有效补偿震动波能量衰减有很大帮助。

3.2   反Q滤波处理

利用本文的基于褶积原理求取品质因子Q的方法与改进广义S变换结合的反Q滤波法对传感器采集到的信号进行衰减补偿。对原始信号进行改进广义S变换处理，然后对处理后的信号根据本文中基于褶积原理求出品质因子Q值，将该值加权到经改进广义S变换后的信号，经过广义S反变换得到重构信号（图5）。

图  5  1号（a）和4号（b）传感器的重构信号

Figure  5.  Reconstructed signal of sensor 1 （a） and sensor 4 （b）

下载: 全尺寸图片 幻灯片

图6为经过反Q滤波处理后的幅值图，图中红色曲线为反Q滤波处理之前的信号，蓝色曲线为经过反Q滤波处理后的信号。图6a补偿前的有效频带范围为200—500 Hz，经过补偿后的有效频带范围拓宽为200—600 Hz，补偿后的主频为420 Hz；图6b补偿前的有效频带范围为410—700 Hz，经过补偿后的有效频带范围扩宽为250—900 Hz，补偿后的主频为610 Hz。分析补偿后幅值谱可知，在频率范围内的波形幅值都有了一定提高，说明低频部分和高频部分的能量都得到了补偿，波形的波峰也变得更加尖锐说明主频能量得到恢复。通过对比图7与图4补偿前后的时频谱，发现经过补偿原始信号的频带明显拓宽，分辨率得到了提高。证明了该方法可以较好地补偿震动波吸收导致的衰减，同时，该方法无需手动设置调节效果的参数，对原始信号的干预程度低，适合应用于信噪比较高的地震资料，能够有效地拓宽地震资料的频率，提高分辨率，为更准确地解释地震资料奠定了良好的基础。

图  6  1号（a）和4号（b）传感器接收到原始信号补偿前后的幅值谱

Figure  6.  Amplitude spectrum of original signal received by sensor 1 （a） and sensor 4 （b） before and after compensation

下载: 全尺寸图片 幻灯片

图  7  1号（a）和4号（b）传感器接收的原始信号补偿后时频谱

Figure  7.  Time-frequency spectrum of original signal received by the sensor 1 （a） and sensor 4 （b） after compensation

下载: 全尺寸图片 幻灯片

4.   讨论与结论

本文旨在还原震源信息，补偿震动波衰减信号。将基于褶积原理求取品质因子Q方法与基于改进广义S变换相结合，建立了震动波衰减补偿方法。通过开展震动波衰减补偿模型试验，采用改进广义S变换进行时频特性分析，并通过褶积原理计算品质因子Q值，最后对信号进行衰减补偿得到重构信号，建立了新的反Q滤波方法，可以使衰减能量得到有效的恢复，为获得复杂构造的地质剖面信息打下基础。主要结论如下：

1） 信号经过改进广义S变换进行时频分析后，发现试验得到的信号频率集中在0.5—1 kHz，同时可以清楚的看出信号的时频分布情况，对应关系准确，能量主要集中在350—700 Hz的高频段，清晰地展示了能量集中区，保留了信号的原始特征。

2） 采用基于褶积原理求取品质因子Q方法与基于改进广义S变换结合的反Q滤波法对震动波衰减补偿模型试验的数据进行反Q滤波处理，在拓宽信号频带和补偿衰减部分的效果显著，信号的有效频带拓宽和振速幅值的衰减得到补偿，进而有效地恢复了衰减的能量。

当前的反Q滤波法对于大数据量叠前的地震资料处理，计算效率较低并且Q值的计算模型过于理想化，在实际应用中，造成工作量的增加。本文将基于褶积原理求取品质因子Q方法与基于改进广义S变换的反Q滤波法相结合，高效地提取了时变Q值，近似地按照波的传播路径进行反Q滤波处理。但该方法对提取品质因子Q值的精度有较高要求，如何保证该方法在不同地质情况的计算精确度以及进一步增强对噪声的压制效果是下一步研究的重点。此外，随着多波和多分量探测技术的发展，如何实现结合P波和S波的震动波衰减能量补偿将是未来的一个热点。Q值计算和多源信息融合技术的研究有望成为未来Q值准确估计和应用的发展方向之一。