目前我国的能源消费结构中煤炭约占67%，预计到2050年仍将占50%以上(闫江伟等, 2013)，因此煤炭在很长一段时间内仍为国民经济发展依赖的最主要能源.60%以上的矿井开采的煤层形成于石炭-二叠纪，形成地质时期早，演化历史与地质条件复杂，煤与瓦斯突出灾害较为严重.随着煤矿开采深度的增加，这种情况将愈加明显.因此对煤矿进行地下结构探测和瓦斯含量预测成为减少煤矿安全事故的关键问题(孙广珍, 2007).

瓦斯是煤在历史时期变质过程中形成的，生于煤层、储存于煤层或者围岩中的气体.导致瓦斯聚集及突出的地质机理主要有三种形式(王德斌和王建涛, 2010)：一是在煤层及其顶底板围岩的裂隙带中瓦斯贮存且比较集中，形成富集区；二是在软煤层中，伴随地层压力的作用，以及煤层采动引起瓦斯动态后移，遇到一定封闭性质地质构造或高应力带阻挡聚集；三是在岩石巷道进入煤层巷道时，原完整岩石对煤层瓦斯具有自然封闭性，但由于采掘活动破坏了自然封闭富集区的完整性，造成瓦斯的涌出或动力突出.

当煤层中的瓦斯压力达到或超过0.74 MPa时，该煤层可被定义为突出煤层(范喜生等, 2013).煤层中瓦斯的存在会降低煤的受力强度，在受到采掘活动影响时，富含瓦斯的煤层易产生破裂，并使瓦斯从高应力区向低应力区转移，进而可能引发煤和瓦斯的突出事故.因此，在煤矿采区规划及工作面布置前，如何有效的探测瓦斯富集区成为高瓦斯矿井安全生产的必要保障(贾兴旺, 2015).瓦斯含量探测主要包括物理采样测试及通过地球物理属性间接探测方法.物理采样测量主要采用原位吸附实验方法(Fu et al., 2009)、钻屑法等.直接采样方法测量的瓦斯含量相对较准确，但该类方法存在两个方面的局限性，包括：(1) 必须在有限钻孔位置才能测量，不能给出宏观的瓦斯分布特征；(2) 只能在已开采区域，或有地面钻孔的位置，此时对于还未形成巷道开采区域的瓦斯含量信息则无法探测.地球物理类方法探测瓦斯是根据煤层在具有不同瓦斯含量时所表现出来的岩石物理属性不同进行探测.相比较于正常煤层，富含瓦斯气体的煤层具有因富含气体、孔隙导致煤层密度减小、对应的地震波速降低、电阻率降低等属性(Afonso, 2015; Lundberg and Sundqvist, 1986)，该特征成为采用地球物理方法探测煤层瓦斯分布的地质地球物理前提(汪志军等, 2011).在构造相对简单，煤层起伏平缓时，顺煤层地震波速变化较小，若煤层出现一定趋势和范围的速度变化，且对应区域无明显的构造异常，则该异常区域可能为裂隙及瓦斯富集区域.因此，可以通过获得顺煤层的地震波速变化来预测瓦斯含量分布.

获得顺煤层速度的地球物理类方法主要包括矿井内工作面地震槽波探测(任亚平, 2015)、工作面地震CT等井下物探方法.但该类方法只能在已布置好运输巷、回风巷和切眼的情况下才能采集数据，且只能测量一个单独的工作面，无法从井田及矿井尺度进行区域探测.而瓦斯赋存、运移特征需要从区域宏观的尺度进行研究.针对目前瓦斯物理采样探测采样点少、矿井探测区域小、无法从区域尺度进行探测的问题，本文提出采用密集地震台阵的背景噪声面波成像方法，对采区尺度的顺煤层速度进行成像，并结合岩石物理实验建立瓦斯含量与地震波速之间的关系，进而研究预测采区内煤层瓦斯含量的分布.

背景噪声面波成像是近十几年来发展起来的有效的地层速度结构成像方法(Shapiro et al., 2005; Lin et al., 2008)，与常用的三维地震勘探方法相比，是一种被动源成像技术(Socco and Strobbia, 2004).它的基本思想是通过计算台站对之间的噪声互相关函数来近似获得台站间的经验格林函数(Lobkis and Weaver, 2001)，每一个台站均可以被认为是面波传播到其他台站的虚拟震源(Huang et al., 2010).对于布设在地表的台站而言，背景噪声的互相关主要由面波组成，面波传播的走时面可以在密集台阵覆盖的区域内构建.面波能量主要集中在地表以下1个波长内(马振宁等, 2019)，随着深度增加，高频面波逐渐衰减，低频面波占主导成分.因此，对应的高频面波反映浅部地层速度结构，低频面波反映深部地层结构(Yao et al., 2008, 2009).早期的背景噪声成像研究主要集中在5~40 s的周期段，用于反演地壳和上地幔的速度结构(Young et al., 2011; 陈凯琪, 2019).近年来，运用背景噪声提取高频面波成分的研究取得了显著进展，使得这种方法越来越多的运用于城市地下空间、矿产资源勘查等浅地表的研究中(Mordret et al., 2013; Liu et al., 2018; Du et al., 2020; 王娟娟等, 2018; 曾求等, 2020)，而目前尚未见基于密集地震台阵的背景噪声成像方法在煤田探测中的应用.

基于此，本文利用布设在阳泉寺家庄煤矿的密集台阵，从背景噪声的互相关函数中获得短周期的基阶瑞利面波的群速度与相速度频散曲线(Sabra et al., 2005).通过面波直接成像的方法(Fang et al., 2015)获得寺家庄井田区域精细的浅部三维横波速度结构.并基于岩石物理实验建立的地震波速度和瓦斯含量关系，进行瓦斯含量分布预测.该研究为煤矿瓦斯探测及后续地面瓦斯抽放和利用提供了新的技术手段，同时为煤矿区域的安全生产提供了技术支持.

寺家庄井田位于山西省晋中市昔阳县境内(图 1a)，井田发育地层主要有太古界、元古界震旦系以及寒武系、奥陶系、石炭系、二叠系、三叠系.第三系、第四系地层在煤田内有不同程度分布.寺家庄矿区位于华北断块吕梁—太行断块沁水、武乡—昔阳北北东向褶带.武乡—昔阳北北东向褶带主要出露二叠系、三叠系地层，是由一系列不同级别褶皱组成的复式向斜.次级褶曲的轴向为北北东向，向斜宽阔，背斜相对较窄.矿区位于沁水煤田东北部，主要含煤地层为石炭系的太原组和二叠系的山西组.区内含煤地层埋藏较深，保存完整.可采煤层为8₁、8₄、9、15号煤层.

图 1

Fig. 1

图 1 (a) 研究区区域地质背景；(b) 寺家庄井田位置；(c) 地震台站分布，其中蓝色三角形表示使用的96个台站，红色三角为图 2中频谱分析台站.AA′和BB′是东西方向的剖面，CC′和DD′是南北方向的剖面 Fig. 1 (a) Regional geological background of the study area; (b) The location of the Sijiazhuang Coal Mine; (c) The distribution of seismic stations. The blue triangles represent the 96 stations of the seismic array and the red triangles denote the stations used for the spectrum analysis in Fig. 2. AA′ and BB′ are east-west sections, and CC′ and DD′ are north-south sections

阳泉煤业主产矿区是瓦斯灾害严重矿区之一，瓦斯大、易自然、难抽采、瓦斯突出问题较为严重(刘宝军等, 2017).区域内大中型断裂构造较少，矿区内高瓦斯区段面积较大，富存高含量、高压力瓦斯.瓦斯灾害突出显现，对矿井的安全生产构成了严重的威胁.随着工作面的回采，地质条件日益复杂，煤层瓦斯含量和压力大幅上升，迫切需要利用地球物理方法手段，对现有瓦斯富藏条件进行探测，对瓦斯进行合理抽放或利用.

2019年3月16日至4月17日，我们在阳泉煤矿寺家庄井田布设了96个短周期地震仪，频率范围为0.2~200 Hz.台站间隔约为0.5 km，对整个研究区进行了一个月的连续背景噪声记录，采样率为100 Hz.背景噪声的数据处理流程参考Bensen等(2007)的背景噪声数据提取面波频散曲线的方法，包括原始数据的分析与预处理，计算噪声互相关函数和提取面波的群速度与相速度频散曲线.

首先对原始数据进行分析与预处理，将每个台站连续波形数据分割成以天为单位的文件；再对各个台站每天的数据进行频谱分析，如图 2所示.通过频谱分析，我们发现连续采集的数据在某些固定频率能量过强，例如5 Hz、15 Hz、25 Hz.通过对比各个台站数据不同时间的频谱图，发现这些5 Hz、15 Hz、25 Hz的单频噪声只出现在某些台站的部分时间段，说明不是因为仪器本身问题造成的.在数据采集期间，位于工区西南部有部分掘进巷道在施工，推断该单频信号可能来源于地下施工机械的震动.为了防止单一频谱对背景噪声信号的干扰，在做互相关之前，对这些频率做陷波处理.完成频谱分析后，进而对数据做谱白化处理，并对数据做时间域去均值、去趋势等预处理.

图 2

Fig. 2

图 2 所选择的四个台站(图 1c中红色三角形)的频谱.台站分别为YQ001、YQ004、YQ053、YQ079 Fig. 2 The spectrum of selected stations (red triangles in Fig. 1c) of YQ001, YQ004, YQ053, and YQ079

原始时间域数据预处理后，我们对每个台站对数据的Z分量数据做互相关运算(Roux et al., 2005).在利用互相关函数进行叠加时，我们分别使用了线性叠加和相位加权叠加方法(Stockwell et al., 1996; Schimmel and Gallart, 2007; Schimmel et al., 2011)对互相关函数进行叠加.图 3是两种方法叠加后得到的YQ001号台站和其他台站之间的互相关函数(CF).从图 3可以看出，采用基于S变换的相位加权叠加方法对干扰信号有更明显的压制效果，其互相关函数具有更好的信噪比.这是因为在密集台阵的小区域噪声研究中，短周期的噪声信号较复杂，观测时间相对较短，致使目标信号常与其他高频信号混合在一起.使用相位加权叠加可以通过减弱非相干性噪声信号的方式达到增强相干信号的目的.该方法能够获得更高信噪比的互相关函数，同时可以减弱观测时长不足的影响，相对短的时间数据叠加就可以达到线性叠加的效果(李选涛, 2019).

图 3

Fig. 3

图 3 台站YQ001与其他台站之间不同叠加方式的互相关函数 (a) 线性叠加；(b) 相位加权叠加. Fig. 3 Cross-correlation functions between YQ001 and other stations using different stacking methods (a) The linear stacking; (b) The phase weighted stacking.

对于叠加的互相关函数，我们使用基于多重滤波的图像分析法(姚华建等, 2004; Yao et al., 2006)提取台站对的群速度与相速度频散曲线.同时，为了满足面波传播远场近似，拾取准则为台间距(L)大于1.5倍的波长(λ)(Yao et al., 2011).首先对中心周期为0.2~2 s、间隔为0.05 s、通带宽度为0.025 s的互相关函数进行带通滤波，然后构建面波群速度的速度-周期(v-T)图像.图 4为一个台站对群速度与相速度频散曲线图.利用速度-周期的图像提取相速度时，为避免拾取相速度受到2π相位影响，利用同一台站对群速度与相速度频散曲线通常具有相似趋势的性质(Luo et al., 2011)，将该台站对对应的群速度绘制在相速度图上(图 4b中粉色曲线)，并以此为参照拾取相速度频散曲线.最终获得基阶瑞利面波在0.2~1.4 s的群速度和相速度频散曲线(图 5)，其中相速度频散曲线2091条、群速度频散曲线2360条.从频散曲线可以明显看出，群速度与相速度主要集中在约1.0~1.5 km·s^-1.

图 4

Fig. 4

图 4 瑞利面波群速度和相速度频散曲线图 红色线代表台间距=1.5倍波长，红点代表最后拾取的频散曲线点.(a) 群速度频散曲线分析；(b) 相速度频散曲线分析，粉色线为相速度频散曲线. Fig. 4 Rayleigh wave group and phase velocity dispersion curves The red line represents the two stations spacing =1.5 times the wavelength, and the red dots represent the picked dispersion curves. (a) Group velocity dispersion curve; (b) Phase velocity dispersion curve, and pink line is the phase velocity dispersion curve.

图 5

Fig. 5

图 5 周期0.2~1.4 s间的(a)群速度和(b)相速度的频散曲线 Fig. 5 Dispersion curve of (a) group velocity and (b) phase velocity between periods of 0.2~1.4 s

本研究使用基于射线追踪的面波频散直接反演方法(Fang et al., 2015)确定台阵下方的三维速度结构.相比传统的两步法面波反演，直接反演方法可以利用不同周期下的面波走时与频散确定三维速度结构，同时拟合面波走时数据和频散数据.面波直接反演方法使用了快速行进算法(Rawlinson and Sambridge, 2004)计算每个周期台站间的面波走时和射线路径进行层析反演，并同时进行频散反演.在反演中考虑了浅层速度结构复杂时的弯曲射线传播路径(Yin et al., 2016)，可以更可靠的获取精细的三维速度结构，具体反演流程如下：

通过线性化的反演方法，沿路径i不同频率的面波传播时间扰动可以表示为：

(1)

式(1)中，t_i^obs(ω)为台站间某一周期面波的观测走时，t_i(ω)是模型正演计算的走时，U_k(ω)和δU_k(ω)分别为二维平面上第k个网格点的群速度(或相速度)及其扰动, 其中ω为角频率.群速度(或相速度)的扰动可以表示为：

(2)

式(2)中，α_k(z_j)，β_k(z_j)和ρ_k(z_j)分别为深度上第j个网格的纵波速度、横波速度和密度.由于面波频散对横波速度最为敏感，我们通过经验公式，将纵波速度α和密度ρ的扰动通过横波速度β的扰动来表示.由此，将公式(2)代入公式(1)可以得到：

(3)

式(3)中，R_α(z_j)与R_ρ(z_j)是基于经验公式的比例因子(Fang et al., 2015).公式(3)写成矩阵的形式为：

(4)

式(4)中的d表示的是所有路径的面波走时残差向量，G表示数据敏感度矩阵，m为模型横波速度扰动向量.通过求解式(4)可以得到模型的更新量，然后在新的模型下重新进行正演计算射线路径与敏感核矩阵，反复迭代直到反演收敛.为了使反演具有更好的稳定性，在反演中加入了阻尼正则化和一阶平滑正则化(Aster et al., 2018).

为了获得更好的三维横波速度反演的初始速度模型，我们首先使用整个区域的平均群速度与相速度的频散曲线(图 6a、6b)进行横波速度反演得到一维速度模型，并将该一维模型拓展成三维模型作为面波一步法反演的初始速度模型.一维频散采用CPS程序包(Herrmann, 2013)完成线性最小二乘迭代反演.

图 6

Fig. 6

图 6 (a) 平均群速度频散曲线; (b) 平均相速度频散曲线，误差棒表示标准偏差; (c) 黑线为随机生成的一系列模型作为初始模型对应的反演结果，红线为平均速度模型; (d) 基于平均速度模型正演得到的群速度和相速度频散曲线(虚线)与平均频散曲线(实线)对比.蓝色代表相速度，黑色代表群速度, 实线为观测值，虚线为理论值 Fig. 6 (a) Average group velocity dispersion curve; (b) Average phase velocity dispersion curve, error bars represent standard deviations; (c) The black lines are the inverted results from a series of randomly generated initial models, the red line is the average velocity model; (d) The comparison between the average group velocity and phase velocity dispersion curves (solid line) and the predicted dispersion curves (dashed line) derived from the average model. Blue lines represent phase velocity, and black lines represent group velocity, the solid line is the observed value, the dashed line is the theoretical value

在对平均相速度和平均群速度进行反演时，初始一维速度模型采用类蒙特卡罗方法随机生成.首先在0.8~2.5 km·s^-1范围内随机生成3000个最大深度为3 km的随机速度模型，然后对所有的速度模型进行一维反演，并且设置残差小于0.02为接受阈值.由此我们获得了1600个一维横波速度模型(图 6c)，最后将所有的速度模型取平均作为最后的一维横波速度模型，如图 7中的红线所示.为了验证获得的平均一维速度模型的合理性(Gu et al., 2019)，基于此模型我们正演了频散曲线.与平均频散曲线对比(图 6d)，可以看出基于平均一维模型正演得到的频散曲线可以较好的拟合观测数据的平均频散曲线，说明了使用该一维模型作为三维反演初始模型的合理性.

图 7

Fig. 7

图 7 不同周期(0.2 s，0.4 s，0.7 s，0.9 s，1.2 s，1.4 s)的基阶瑞利面波群速度与相速度对横波速度在不同深度上的敏感核 (a) 相速度敏感核；(b) 群速度敏感核. Fig. 7 The depth-varying sensitivities of fundamental-mode Rayleigh wave group velocity and phase velocity at different periods (0.2 s, 0.4 s, 0.7 s, 0.9 s, 1.2 s, 1.4 s) to shear wave velocity (a) Phase velocity sensitivity kernels; (b) Group velocity sensitivity kernels.

另外，为获得研究区域面波频散数据在深度方向上对横波速度的约束能力，图 7为基于图 6c获得的该区域一维速度模型计算得到的基阶瑞利面波群速度与相速度对横波速度的深度敏感核.从图 7中可以看出，在选取的面波频带内(0.2~1.4 s)，提取的频散曲线可以对800 m以上的地层速度结构具有较好的约束.

研究区域内台站间面波射线路径的覆盖密度很大程度决定了反演结果的可靠性与分辨率.图 8a和图 8b分别给出了研究区域不同周期(0.2 s，0.6 s，1 s，1.2 s，1.4 s)的群速度与相速度的射线覆盖情况.从图中可以看到，整体射线的覆盖密度是中间密，四周相对稀疏，且随着周期增加，射线数量减少，射线分布密度变小.但总体来看，周期0.2~1.0 s范围内，由于台站的均匀分布，各周期群速度与相速度对于研究区域都有良好的覆盖，能满足走时层析成像要求.本研究区最大可靠成像深度位于1.0 s周期瑞利面波的敏感深度处.

图 8

Fig. 8

图 8 不同周期(0.2 s, 0.6 s, 1.0 s, 1.2 s, 1.4 s)下射线路径分布图 第一行为群速度射线路径分布; 第二行为相速度射线路径分布. Fig. 8 Ray path distributions at different periods (0.2 s, 0.6 s, 1.0 s, 1.2 s, 1.4 s) First row, group velocity ray path distribution; Second row, phase velocity ray path distributions.

此外，为了测试当前数据分布(台站间射线分布)和反演网格下的模型分辨率，我们进行了棋盘分辨率测试(图 9a)，首先将研究区域在水平面上参数化为14×18个网格，南北方向与东西方向的网格间距均为0.005°(~500 m).在深度方向上的网格点数为20个.通过前面深度敏感核的分析，设置模型的深度范围为0~2 km.合成测试的棋盘模型以图 6c所示一维平均横波速度模型为基础，设置异常大小为水平和垂直方向分别为2×3×4个网格.我们使用该棋盘模型计算与实际数据分布一致的台站间不同周期的瑞利波走时，在合成数据中添加2%的随机噪声，并使用与实际数据相同的反演策略进行反演.图 9b显示四个不同深度的棋盘模型恢复的水平速度切片.结果表明，因为台站分布及射线数据的覆盖程度较好，反演区域不同深度的异常均能得到较好的恢复，其中300~500 m深度的恢复程度最好(台间距~500 m提取的频散最高频为5 Hz，面波分辨率随深度减小).

图 9

Fig. 9

图 9 棋盘分辨率测试(切片深度：300 m，400 m，500 m，600 m) (a) 合成测试棋盘模型; (b) 棋盘恢复模型. Fig. 9 Checkerboard resolution test (slice depths: 300 m, 400 m, 500 m, 600 m) (a) Synthetic checkerboard model; (b) Recovered checkerboard model.

在实际观测数据反演中，我们将前面获得的一维平均速度模型扩充成三维初始速度模型，利用群速度和相速度频散曲线联合反演三维横波速度结构.反演结束后，大部分射线路径的走时残差都接近于零(图 10a).通过7次迭代走时残差的均方根由反演前的1.64 s下降到0.21 s(图 10b)，所有残差的均值为0.013 s，表明反演得到的三维模型可以很好地拟合观测数据.

图 10

Fig. 10

图 10 (a) 反演前后的走时残差分布图和(b)走时残差随迭代次数变化图 Fig. 10 (a) Comparison of traveltime residuals before and after inversion; (b) The variation of traveltime residuals with iterations

图 11和图 12分别为三维横波速度模型沿不同测线(AA′、BB′、CC′、DD′)的垂直速度剖面图和在不同深度(300 m、400 m、500 m、600 m)的水平速度切片图，其中剖面位置如图 1c所示.从图 11剖面图中可以看出，研究区域浅部的横波速度相对较低，约为1.0 km·s^-1.在400~800 m深度范围内速度变化整体上呈水平层状分布.研究区15号煤层埋深在400~500 m深，沿煤层的地震横波速度分布在1.0~1.5 km·s^-1左右，整体背景速度在2.0 km·s^-1.500 m以深为奥陶系灰岩，整体呈高速特征.速度在垂向上从浅至深符合地层的岩性特征，这也进一步说明密集台阵的背景噪声在浅地表矿产尺度可以获得较为可靠的速度模型.

图 11

Fig. 11

图 11 沿图 1c中AA′、BB′、CC′、DD′的垂直剖面横波速度图 Fig. 11 The results of shear wave velocity along the vertical sections of AA′, BB′, CC′ and DD′ in Fig. 1c

图 12

Fig. 12

图 12 三维横波速度结构水平切面图.(a)、(b)、(c)和(d)分别为300 m，400 m，500 m和600 m深度的横波速度结构. 黑色三角形表示台站位置. Fig. 12 Map view of the three-dimensional shear wave velocity model at different depths of (a) 300 m, (b) 400 m, (c) 500 m and (d) 600 m The black triangles indicate the locations of seismic stations.

从图 12速度水平切面可以看出，低速异常呈北西方向展布，但地表三维地震反射勘探显示，研究区煤层相对比较完整且起伏较小，断层不发育，地层结构相对简单.当煤层横向连续、地层结构相对简单时，地震波速度应在较小范围内变化.从速度平面图可以看出，低速异常有一定连续性与空间展布范围.

陈信平等(2013)研究发现煤层的瓦斯含量与其密度、波速之间存在负相关关系，即煤层瓦斯含量高，对应煤系地层密度小、地震波速低；煤层瓦斯含量低，对应煤系地层密度大，地震波速高.并认为煤层气储层的含气量与其弹性参数(例如速度，密度)之间的负相关关系可能是煤层这类储层内在的固有的规律性的关系.同时也有大量的岩石物理实验的结果支持该观点(申振华，2011；陈信平等，2013；刘盛东等，2015).煤系地层的地震波速度受到不同因素影响，比如温度、流体、压力以及地层完整性等影响.依据三维地震结果及巷道实际揭露，研究区煤层横向连续性较好、煤层厚度变化及起伏均较小，因此，影响煤系地层速度的主要因素为瓦斯气体含量及地层裂隙.地层裂隙同时又成为瓦斯富集及运移的通道，在煤系地层地质结构相对简单的区域，低速异常可能即为瓦斯富集区对应的地球物理特征.因此，我们在获得顺煤层速度分布的基础上，结合岩石物理实验建立的煤层速度与瓦斯含量的关系，进而预测研究区瓦斯含量的分布具有可行性.

阳泉寺家庄矿井的含煤地层为石炭系上统太原组(C3t)和二叠系下统山西组(P1x)，含煤地层总厚度约160 m，分别为8号、9号、11号、12号、13号、14号、15号煤.其中15号煤层在本区域内稳定沉积，煤层结构简单，在全区内均有分布，厚度稳定，属于全区稳定可采煤.这也为根据地震波速度预测瓦斯富集区提供了地质基础.

实测数据表明，以甲烷(瓦斯)为代表的烃类气体的密度为0.7168 mg·cm^-3，声波时差为2260 μs·m^-1；而烟煤的密度为1.25~1.35 g·cm^-3，声波时差为400~560 μs·m^-1.通过阳泉煤矿区域的超声波岩石物理实验(Wang et al., 2019)，得出煤样吸附能力(Q)与气压值的关系如图 13a所示.随着气压值的增大，煤样的瓦斯吸附能力呈现上升趋势.超声波煤样弹性测试的结果如图 13b所示，随着气压的升高，煤样的超声波P波速度呈现下降趋势，这为根据地震波速度预测瓦斯含量提供定量分析的基础.

图 13

Fig. 13

图 13 煤样气体吸附能力及地震波速随压力变化关系.数据来源于Wang等(2019) (a) 煤样吸附能力(Q)与气压值(P)关系; (b) 煤样地震波(VP)速度随压力(P)变化关系. Fig. 13 The variations of the gas absorption capacity and seismic wave velocity of coal samples with pressures. Data from Wang et al.(2019) (a) The relationship between coal sample absorption capacity (Q) with the pressure (P); (b) The relationship between seismic velocity and pressure of the coal samples.

图 14为通过岩石物理实验得到的实测瓦斯含量与煤层纵波速度关系图(申振华, 2011)，红色散点为不同瓦斯含量对应的地震纵波速度，黑色实线为通过散点拟合得到的回归表达式.可以看出瓦斯含量与地震纵波速度具有反比例对应关系，即随着瓦斯含量增大，地震纵波速度降低.通过回归计算，以幂指数函数拟合得到瓦斯含量与速度对应关系为y=1+703.3e^-0.0018x(y为瓦斯含量(m³/t), x为地震纵波速度(m·s^-1)).这也进一步说明，在瓦斯含量与地震波速度具有较好的回归关系的基础上，使用地震波速度来预测寺家庄矿区的15号煤层的瓦斯含量是可行的.

图 14

Fig. 14

图 14 瓦斯含量与地震纵波速度的关系图. 数据来源于申振华(2011) Fig. 14 The relationship between methane content and seismic P wave velocity. Data from Shen (2011)

由于背景噪声成像得到的是地震横波速度，而岩石物理实验(图 14)得到的是地震纵波速度与瓦斯含量之间的关系.因此，在利用地震横波速度估算瓦斯含量之前首先要把横波速度转化为纵波速度，转换公式采用Brocher(2005)提出的纵横波速度转换经验公式，V_P=0.9409+2.0947V_S-0.8206V_S²+0.2683V_S³-0.0251V_S⁴，在得到纵波速度的基础上，利用公式y=1+703.3e^-0.0018x，即可估算得到沿15号煤层的瓦斯含量(图 15a).图 15a中的黑色方框为寺家庄矿井某一巷道位置，该巷道均匀分布瓦斯实测点.沿15号煤层的瓦斯含量整体在8 m³/t左右，其中研究区域的西南与东南部较低，最低值约为5 m³/t.研究区域中部瓦斯含量较高，最高约为12.8 m³/t.

图 15

Fig. 15

图 15 瓦斯含量预测与实际巷道测量对比图 (a) 沿15号煤层瓦斯含量分布平面图；(b) 巷道内瓦斯预测含量；(c) 巷道内瓦斯实测含量. Fig. 15 Comparison of methane content between prediction and actual roadway measurement (a) Distribution of methane concentration along No.15 coal seam; (b) Predicted methane concentrations in the roadway and (c) Measured methane concentrations in the roadway.

图 15b与图 15c分别为在巷道位置通过背景噪声成像结果与岩石物理实验预测的瓦斯含量分布与实际测量瓦斯含量分布对比图.图 15b为预测的该巷道位置的瓦斯含量分布图，图 15c为在巷道内实际测量的瓦斯含量分布图.从对比图中可以看出，寺家庄这一巷道的瓦斯含量预测值与实际测量值具有相同的趋势，但是在数值上整体偏低.主要的原因是因为数据中包含噪声，因此背景噪声成像需要使用平滑约束和阻尼等正则化方法使模型稳定，进而影响速度异常的恢复.但相比传统的瓦斯含量估算方法，通过背景噪声速度成像预测瓦斯含量可以一次获得整个研究区的瓦斯分布特征，为后续采区设计工作面布设提供参考依据.同时也进一步说明：通过背景噪声成像得到煤层速度，并结合岩石物理实验预测瓦斯含量具有合理性，也为煤矿瓦斯含量预测研究提供了一种全新的技术方式.

本文利用在原阳泉煤业(集团)有限责任公司寺家庄煤矿布设的96个地震台站记录的一个月的连续波形数据，基于背景噪声成像方法获得了寺家庄煤矿的浅部横波速度结构.根据岩石物理实验建立了地震波速度和煤层瓦斯含量之间的关系.煤层内的瓦斯含量高，会导致煤层的骨架强度降低，从而使地震波速下降.15号煤层低速异常呈北西走向，推测为瓦斯含量较高区域.通过背景噪声获得的速度结构以及岩石物理实验获得的速度-瓦斯含量经验关系对整个研究区的15号煤层的瓦斯含量分布进行了估计，获得了与实测巷道瓦斯含量具有相同趋势的预测结果，为煤矿后续的开采生产提供了技术支持.本次研究表明通过背景噪声成像得到的速度结构，对于煤系地层构造相对简单的高瓦斯矿井，可以获得煤层的瓦斯含量分布，这为研究煤矿瓦斯分布提供了一种全新的技术形式.