1. 引言
随着我国现代化建设的不断加快,城市地下资源逐渐成为我国经济发展的重要战略资源。因此,目前我国对地下空间的利用逐渐往深部地区发展,然而很多城市地区对深部地层的现有地质资料缺失严重,导致开发地下空间资源时容易出现安全隐患。
在城市地区的地质勘探中,重力勘探、高密度电阻率法以及钻孔探测等多种工程物探方法被广泛应用。例如,罗国平[1]于2000年在华北平原的城市区域使用重力勘探技术成功地识别了地下热水储存的有利位置——基岩坳陷中的凸起构造。这种方法不仅在地球深部结构研究、石油与煤炭资源普查、固体矿产开发及水文调查等多个领域发挥着日益重要的作用,还展示了其在寻找特定地质构造方面的有效性。苏永军等人[2]在2023年采用高密度电阻率法对雄安新区已填埋坑塘的空间特征进行了详细探测。通过该方法获得的数据经过钻孔验证后显示了高度的一致性,表明高密度电阻率法在横向和纵向分辨率上均表现出色,能够提供浅层丰富的地质信息,满足精细化探测的需求。然而,这些传统工程物探方法在近地表浅地表结构勘探方面也存在一定的局限性。重力勘探虽然有效,但其分辨率相对较低,且由于不同地质结构可能产生相似的重力异常,可能导致结果的不确定性增加。高密度电阻率法等电磁类物探方法在城市环境中可能会受到来自各种电磁信号的干扰,从而影响探测结果的精度。钻孔探测作为一种直接获取地下介质信息的方法,尽管准确性高,但施工成本高昂,并且难以揭示钻孔之间的地质情况,导致钻孔之间地下结构的不确定性较大。相比之下,面波勘探因其施工简便、成本低廉而成为近几十年来广泛使用的近地表勘探方法。
面波勘探能够有效地探测浅层地质结构[3],为城市地区的地质调查提供了另一种可靠的选择。根据震源类型不同,面波勘探可分为主动源面波勘探方法和被动源面波勘探方法。Park等人[4]在1998年提出的多道面波分析法(Multichannel Analysis of Surface Waves, MAPS)因其具有非侵入性、非破坏性和精度相对较高的特点,已经成为主动源面波勘探中常用方法。但这种依靠人工震源的主动源面波方法易受其他震源的干扰、场地的限制且探测深度相对较浅。随着面波勘探理论的不断成熟,近十年来被动源面波勘探方法在地球结构的多个尺度上发挥了作用。Aki [5]于1957提出,通过使用空间自相关方法(Spatial Auto-Correlation method, SPAC)在瑞雷波中提取频散曲线进行反演,得到横波速度结构。原始的SPAC对排列要求比较严格,要求检波器呈圆形排列,这在实际应用中较为困难。Lin和Okada在1993年提出了扩展空间自相关(Extended Spatial Auto-Correlation method, ESPAC)方法,这使得使用任意形状的阵列进行空间自相关面波频散分析成为可能。Cheng等人[6] [7]在2016年提出了基于互相关函数的被动源面波多道分析方法(Multi-channel Analysis of Passive Surface Waves, MAPS)。针对线性排列的多道背景噪声记录,MAPS方法利用地震干涉法计算得到大量互相关函数对,按照台站间距排列,形成类似主动源面波信号的虚拟炮集记录。汪利民等人[8]在2022年介绍了一种背景噪声多阶面波频散成像方法,称为频率贝塞尔变换法(Frequency-Bessel transform method, F-J)。席超强在2021年提出了改进的频率贝塞尔变换方法[9],使用单方向速度扫描算子汉克尔函数替换双方向速度扫描算子贝塞尔函数,能够有效避免“交叉”假频影响,获取更高质量的多阶瑞雷波频散能量图像。与MAPS方法类似,该方法也需要首先计算得到互相关函数。随后利用频率贝塞尔变换,将数据从时间–空间域转换到频率–速度域。上述的这类干涉类频散分析方法被广泛用于通过互相关背景噪声来提取两个检波器之间的近似格林函数,而格林函数表示在一个位置上观测到的信号是由另一个位置脉冲源的响应。在检波器数量较多的情况下,可以根据距离对格林函数进行排序,形成虚拟炮集。近年来,随着背景噪声成像技术的成熟,越来越多的研究者关注于利用城市中人为背景噪声探测城市地区的近地表结构。Nakata在2016年基于交通噪声使用地震干涉法成功提取出了体波和面波[10]。Quiros等人在2016年从铁路交通噪声信号中恢复了体波和面波[11]。Zheng等人在2021年利用背景噪声层析成像技术,检测城市地区的隐藏断层[12]。Mi等人在2022年利用城市道路沿线的交通噪声对近地表结构进行了刻画[13]。刘红蕾[14]利用采集自杭州市的三分量实测交通地震背景噪声数据,分别提取出垂向分量和径向分量中的瑞利波频散曲线,联合两个分量的面波频散曲线反演出测区地下横波速度结构,以此验证多分量交通噪声在实际应用中的适用性、有效性。
以往,由于受到单分量检波器只能采集垂向数据的限制,主要的研究集中在利用地震干涉法在从垂向分量的地震数据中,获取瑞雷波信号,反演其频散曲线得到地下横波速度。现在随着多分量检波器的出现,可以同时获取垂向、径向和切向三个分量的地震数据。其中在垂向分量和径向分量地震数据中可以提取到瑞雷波信号,在切向分量中可以提取出勒夫波信号[15]。利用瑞雷波频散曲线和勒夫波频散曲线进行联合反演,能够有效避免单一瑞雷波频散曲线反演中常见的“模式接吻”现象。Chmiel等人[16]在2019年通过使用瑞雷波与勒夫波联合反演格洛宁根气田区域的横波速度证明了与单独反演瑞雷波频散曲线相比,使用瑞雷波频散曲线和勒夫波频散曲线联合反演可以有效地提高精度并减少非唯一性。
本文以当涂某地区地下结构为研究对象,采用三分量检波器布置密集线性阵列记录背景噪声。首先,分析该地区的背景噪声的频谱特征,判断噪声的主要来源。然后利用地震干涉法提取高质量的瑞雷波和勒夫波频散曲线。最后结合联合反演得到了横波速度剖面。与现有地质资料进行了对比,揭示了土壤和基岩的地下结构。
2. 数据采集与处理
2.1. 调查区域和数据采集
2.1.1. 调查区域概述
乌溪镇位于安徽省当涂县大公埒东南角的水阳江北岸,处于两省(安徽、江苏)四县(当涂、湾沚区、宣州、高淳)交汇处。东邻江苏省高淳县,南与安徽省宣城和芜湖隔河相望,西接黄池镇,北连塘南镇。乌溪镇自古以来便是商业繁华之地,边贸发达,是当涂县著名的八大古镇之一。
本次调查的重点区域位于乌溪镇乌溪大巷附近。该区域计划建设住宅区,因此,在施工前进行准确细致的岩土结构和基岩深度调查显得尤为重要。根据之前的地质调查和钻孔测量,该地区主要为平原地貌,包括冲积平原、冲积湖积平原和湖积平原等地貌类型。平均海拔在10至20米之间,地势平坦,由粘质砂土和砂质粘土构成。
2.1.2. 数据采集
调查区域被道路和城镇环绕,如图1所示,周边有多条公路。道路上频繁行驶的汽车和公交车等交通工具提供了充足的背景噪声源。采用了主频为5 Hz的地球脉动三分量检波器进行数据采集。这些检波器沿道路方向布置,测线包含了103个站点。从2024年3月20日中午12点00分开始,连续记录了一个半小时的噪声数据。节点之间的部署间隔为5米,测线总长度约为500米。
为了便于后续的数据处理和提高工作效率,所有节点均沿道路方向平行部署。这种部署方式有利于提高数据的一致性和质量,从而更好地进行地质结构分析。
Figure 1. Satellite map of the line area
图1. 测线区域卫星图
2.2. 三分量背景噪声特征
为了更好地了解观测到的背景噪声记录的频谱特征,本文对原始的三分量数据进行了频谱分析。本文随机选择该条测线上的48号站点作为示例,并在图中显示出了时频信号图。从得到的图上看能量主要集中在2~15 Hz这个频段内,有两个峰值,分别为2 Hz和8 Hz。这些频率范围可能与周边噪声源类型有关。
另外,如图2所示,在时频信号图中本文发现Z分量功率一般比N和E分量的强。根据之前的研究,环境噪声中包含瑞雷波和勒夫波。然而,环境噪声中许多噪声源是垂直激发的,比如人类行走。车辆的运动虽然既会产生垂直信号,也会产生水平信号,但垂直方向的能量会更加集中。这可能会使观测到的垂直分量能量相对于水平分量的能量更强。而且,从时频信号图上可以看出,12点到12点30以内噪声能量更强,频率更高,与人类的活动变化相一致。结合调查区域周边存在多条公路可以判断该地区背景噪声主要是由交通噪声主导的。
Figure 2. Three-component time-frequency signal diagrams with Z (a), E (b) and N (c) components, respectively
图2. 三分量时频信号图,分别为Z (a)、E (b)、N (c)分量
2.3. 数据预处理和噪声互相关
在提取背景噪声中的地震信号进行互相关之前,本文按照工作流程对环境噪声进行预处理。首先,将数据分成30 s的片段,再对每个台站的三个分量应用相同的预处理参数,包括去除片段的均值和趋势,进行时间归一化和频率谱白化以平衡所有频率的能量。在本研究中,本文使用1个半小时的噪声记录来提取高质量的面波信号。
2.3.1. 噪声互相关
为了得到虚拟炮集,本文使用地震干涉法对提取的三分量地震信号进行互相关。
在均匀入射条件下,ZZ分量瑞雷波互相关函数表达式表示为:
(1)
RR分量互相关函数为:
(2)
TT分量互相关函数为:
(3)
其中,
为瑞雷波的垂直相对水平方向的振幅比,
为第一类零阶贝塞尔函数,
为第一类第二阶贝塞尔函数,
为瑞雷波功率谱,
为勒夫波功率谱,k为波数,r为垂直分量站台对之间的距离,ω为角频率。
在互相关之后,本文使用线性叠加来平均时间噪声变化的影响,并改善相关信号。最后,本文得到了三个分量各自台站对之间的3个互相关函数。结果如图3所示。在ZZ分量中可以观察到清晰的瑞雷波信号(图3(a))。EE分量(图3(b))也显示出明显的勒夫波信号。对比图3(c),本文注意到勒夫波比瑞雷波到达得更早。
Figure 3. Three-component cross-correlation functions with ZZ (a), EE (b) and NN (c) components, respectively
图3. 三分量互相关函数,分别为ZZ (a)、EE (b)、NN (c)分量
2.3.2. 三分量频散曲线的提取
根据采样定理,检波器的布设长度应该足够长,以提取瑞雷波和勒夫波的可靠速度。水平分辨率随着检波器的布设长度的增加而降低。根据最佳检波器展布长度(约为最大调查深度的两倍)的选择经验,每32个检波器为一个子阵列,以获得80 m深度以内的横波速度剖面。对于每个子阵列,可以获得以第一个台站作为虚拟震源,其余台站作为接收机的虚拟炮集。将虚拟炮集的因果支和非因果支叠加后平均,并应用相移法从ZZ、EE和NN分量虚拟炮集中获得频散光谱。
2.3.3. 相移法成像原理
相移法是Park等人于1998年提出的一种频散曲线提取方法。
原理是设时间–空间域中的一个炮集记录为
,当对该记录在时间域上做傅里叶变换时,其频谱
可以表达为振幅谱
和相位谱
的乘积,即:
(4)
相位谱携带着面波频散的所有信息,而振幅谱则反映了除频散外的其他地震信号信息。对于二维单频平面波,其表达形式为:
(5)
其中,t为时间;f为频率;x为偏移距;k为波数;
为初相位。
经傅里叶变换后,相位谱记录了不同频率简谐波在
时的相位
,因此相位谱仅依赖于偏移距x,并随着位置线性变化,与初至时间无关。所以频谱表达式可写为:
(6)
其中,
,
为圆频率;
为对应的相速度。
将式(6)中特定频率的波场经相位平移和振幅标准化后沿炮检距方向叠加:
(7)
若面波沿地表传播,则存在下面关系:
(8)
对于目标频率范围,通过使用不同的相位值
值进行搜索,并在
域进行积分可以得到能量分布图。当积分结果接近最大值(即相当于所有道的总和)时,表明已经成功拟合了相位线,意味着找到了与数据最匹配的速度值。利用下面公式:
(9)
可以进一步转换至
域频散能量图。进而提取出所需的频散能量曲线。在图4(a)中,在2~15 Hz频段可以观察到清晰连续的瑞雷波频散能量。与ZZ分量相比,RR分量中基阶模态的频谱带更短,在低频部分消失,而第一个高阶模态出现,能量较强,如图4(c)所示。此外,在图4(b)中,从EE分量获得了勒夫波,包括基阶模态和高阶模态,考虑到提取的高质量频散曲线,本文将ZZ分量和EE分量互相关函数频散分析中得到的频散能量分别表示为瑞雷波和勒夫波的频散能量。各子阵列道间距均为5 m,由此得到72条瑞雷波和72条勒夫波频散曲线。
Figure 4. Dispersion diagrams of background noise signal with ZZ (a), EE (b) and NN (c) components, respectively
图4. 背景噪声信号频散图,分别为ZZ (a)、EE (b)、NN (c)分量
3. 反演及结果
本文利用嵌入库恩蒙克雷斯算法的模式搜索(The Pattern Search with Embedded Kuhn-Munkres Algorithm, PSEKM)算法对频散曲线进行联合反演,该算法是由Yan等人于2022年提出的一种多阶模式频散曲线反演方法。该方法将反演过程分为两个阶段:第一阶段采用模式搜索(Pattern Search, PS)来反演基阶频散曲线,以获取一维横波速度;在第二阶段,将第一阶段得到的反演结果设置为初始模型,并将高阶频散曲线添加到反演系统中,采用PSEKM算法来最小化拟合误差。
为了获得每个子阵中部地下结构的一维横波速度剖面,本文联合反演了从每个子阵提取的瑞雷波和勒夫波频散曲线。由于频散谱质量的限制,只选择2 Hz以上的频散曲线。
与瑞雷波单独反演相比,联合反演的优点是可以减少反演的非唯一性,提高反演模型的可靠性。同时为了方便反演结果与以往地质资料的对比,笔者将一维深度模型参数设为4层,其中底部为半空间。为了评估反演的准确性,本文展示了三个一维反演的例子(图5)。基于瑞雷波和勒夫波反演结果的理论频散曲线同时收敛于实测频散曲线(图5)。联合反演得到的瑞雷波和勒夫波的平均均方根误差分别为7.92 m/s和10.87 m/s,表明频散曲线拟合良好。
Figure 5. Inversion of Rayleigh wave and Love wave dispersion curves: with (a) Fitting of Love wave dispersion curve at point 24, (d) Fitting of Rayleigh wave dispersion curve at point 24, (g) Inversion velocity profile at point 24, (b) Fitting of Love wave dispersion curve at point 45, (e) Fitting of Rayleigh wave dispersion curve at point 45, (h) Inversion velocity profile at point 45 and (c) Fitting of Love wave dispersion curve at point 60, (f) Fitting of Rayleigh wave dispersion curve at point 60, (i) Inversion velocity profile at point 60
图5. 瑞雷波与勒夫波频散曲线反演情况:分别为(a) 24号点勒夫波频散曲线拟合情况、(d) 24号点瑞雷波频散曲线拟合情况、(g) 24号点反演速度剖面、(b) 45号点勒夫波频散曲线拟合情况、(e) 45号点瑞雷波频散曲线拟合情况、(h) 45号点反演速度剖面和(c) 60号点勒夫波频散曲线拟合情况、(f) 60号点瑞雷波频散曲线拟合情况、(i) 60号点反演速度剖面
分别将瑞雷波单独反演与联合反演得到的一维横波速度模型结合起来,形成2维横波剖面(见图6),单一使用瑞雷波进行反演时,所获得的地下结构信息显得较为模糊,尤其是在40米处,模型揭示出一个明显的低速异常区。相比之下,联合反演方法提供了更加精细和清晰的地下结构图像,不仅消除了低速异常区,而且与地质资料的高度一致性,有力地验证了瑞雷波与勒夫波联合反演的有效性和准确性。
联合反演结果成功地揭示了地表至80米深度范围内连续且可靠的2维横波速度分布特征。根据地质资料显示,该研究区域的地表覆盖层主要由粉砂质粘土和细粒粉砂构成,这是河海交互作用下的典型沉积产物。在剖面中,相对横波较高速度的区域对应于砂砾。横波速度超过500 m/s的区域对应基岩部分,主要由不同风化程度的泥质粉砂岩构成。随着深度的增加,岩石的风化程度逐渐减弱,导致S波速度呈现递增趋势。在水平方向上,由于研究区域位于坡度极缓的地形上,基岩界面基本平坦,这一发现与联合反演的结果相吻合。
Figure 6. 2D shear wave velocity profile: (a) Rayleigh wave inversion shear wave velocity profile; (b) Joint inversion shear wave velocity profile
图6. 2维横波速度剖面:(a) 瑞雷波反演横波速度剖面;(b) 联合反演横波速度剖面
4. 结论
本文介绍了一个三分量背景噪声成像技术在浅层地表勘探中的应用案例。本研究选取了中国马鞍山当涂某地区作为研究对象,旨在揭示该区域的近地表结构特征。通过地震干涉法,成功从连续1.5小时的交通噪声记录中提取出了瑞雷波和勒夫波。从切向分量中提取的勒夫波能够为地下介质的特性分析提供额外的信息,有助于减少反演过程中的不确定性,并为背景噪声成像提供更为严格的约束条件。通过对背景噪声产生的瑞雷波和勒夫波进行联合反演,获得了更加精确和可靠的结果,探测深度达到了80米。与仅采用瑞雷波进行反演的方法相比,联合反演方法所得到的二维S波速度剖面与地质资料的吻合度更高,同时也能更好地匹配频散曲线。这些发现表明,在近地表介质的刻画过程中,利用多分量背景噪声不仅能够提高精度和效率,而且由于在高度城市化区域所需的观测周期较短,这种方法还具备地下监测的应用潜力。
基金项目
安徽省教育厅项目(2022AH050836):基于多分量背景噪声的废弃矿井地下空间精细探测研究。
NOTES
*通讯作者。