2. 中国科学院大学, 北京 100049
2. University of Chinese Academy of Sciences, Beijing 100049, China
滑坡灾害在我国频繁发生,通常会造成巨大的人员伤亡、财产损失和基础设施破坏(殷坤龙等,2001;朱良峰等,2004;黄润秋,2007;季伟峰等,2007;Huang,2009),因此有必要对人类住所附近的滑坡体进行监测。滑坡的表面运动学特征目前可以通过大地测量、摄影测量等遥感手段获取(唐亚明等,2012;Wang,2013),然而对于大型滑坡,表面运动很难反映深部的情况。地震学手段常基于滑坡发生时产生的地震信号分析滑坡的动力学特征(Xie et al,2020;Zhao et al,2020),同时也能通过检测滑坡变形过程中微破裂的活动性来了解滑坡内部的变形情况(Gomberg et al,2011;盛敏汉等,2018)。放置在滑坡上的地震仪不但能记录到滑坡变形辐射出的地震信号,还可能会记录到人类活动信号。
滑坡上的地震仪能记录到大量的震动信号,包括滑坡体外事件的地震信号、与滑坡体变形相关的微震信号、滑坡上的滚石信号;以及非自然信号,包括车辆、行人经过地震仪所引起的震动信号等;甚至暴雨降临时地震仪也有记录(Gomberg et al,2011;盛敏汉等,2018)。与滑坡形变相关的微震称作滑坡震(slidequake),其指滑坡体在运动过程中边界、基底和围岩的相互作用以及局部破裂产生能被地震仪记录到的信号(Helmstetter et al,2010)。Yamada等(2016)在日本的Rausu滑坡附近的地震记录中发现该滑坡在2次坍塌前的数小时内发生了大量的重复事件,认为该事件为滑坡底部边界和基岩之间的粘滑(stick-slip)产生,并认为这一过程对于认识滑坡的触发十分重要。
目前研究滑坡震主要是通过对滑坡震信号的识别、定位,进而分析滑坡震的时空变化来认识滑坡的活动性。滑坡上的滑坡震一般震级较小,Spillmann等(2007)在研究瑞士阿尔卑斯地区滑坡时得到滑坡震震级主要分布在-2~-1级之间,这将对地震事件的识别带来挑战。而中国西南地区部分滑坡附近靠近河岸,土地比较肥沃,滑坡附近或滑坡体上常有居民居住,人类活动的干扰对滑坡震的识别带来进一步挑战。因此,鉴别与检测这些地震事件是利用地震学方法研究滑坡内部形变的重要前提。
本研究利用2017年2月至3月在梭坡滑坡上布设的34台地震流动台站记录到的连续地震波形,识别和检测了地震仪记录到的震动事件,通过进一步的波形分析并利用网格搜索的方法对记录较好的事件进行定位。我们发现地震仪不仅记录到了和梭坡滑坡变形相关的滑坡震事件,还清晰记录到了在滑坡上燃放烟花产生的震动事件。
1 梭坡滑坡概况梭坡滑坡位于四川省丹巴县内大渡河上游,其堆积体总体走向为S30°E,坡面倾角约为30°。滑坡前缘高程约1840m,延伸至大渡河河谷,滑坡体积超过3×108m3(陈菲等,2012)。梭坡滑坡侵蚀地貌较为发育,地势险峻,以高山陡坡为主,地质构造发育,降雨集中,人类工程、经济活动较强烈(阳东,2010)。滑坡东侧为藏族古碉群集中区,近年来,梭坡滑坡灾害频发的趋势尤为明显,造成了巨大经济财产损失(许秦坤,2012)。由于古碉群独特的历史、人文和旅游价值,梭坡滑坡的稳定性对古碉群保护和水库蓄水位的选择均十分重要。
为了监测该滑坡体上的滑坡震及其他震动事件,2017年2月至3月在梭坡滑坡上布设了34台地震仪(图 1),平均台间距约为50m,最大台间距约为600m。所用仪器为5s短周期地震仪,固有频率为4.5Hz,采样率为100Hz。
蓝色三角表示布设在梭坡滑坡上的34台地震流动台 |
在进行微震信号检测前,首先将原始数据进行高通5Hz滤波,以削弱远震信号的干扰。由于滑坡震震级通常比较微弱(Spillmann et al,2007),本文采用拾取微弱信号能力较强的滑动时窗互相关(Sliding-window Cross-correlation,SCC)方法(Yang et al,2009;杨慧等,2018;Ma et al,2020)进行微震信号检测。
由于滑动时窗互相关方法检测地震信号需要依靠模板事件,需要依据到时和波形信息寻找模板事件。首先,确定一个所有台站均能记录到的区域地震,并按照P波到时对齐所有波形,以校正GPS时钟差和台站下方速度结构的影响;同时,由于滑坡土质松软,波速较低,如果事件发生在滑坡上将会导致明显的到时差,我们随后人工查看了大量对齐后的波形数据(一周)观察是否存在明显到时差,并试图找到分布在滑坡上所有类型的事件;最终,确定了两类信噪比较高的事件(图 2):第Ⅰ类事件(A)只能被D04台及研究区域西北部的少数几个台记录到,可能与事件本身能量较小有关;第Ⅱ类事件(B)信号在D12台最先被观测到,推测事件发生在D12台附近。
纵坐标为旋转到滑坡面上的震中距;红色虚线为参考视速度;波形滤波范围5~30Hz |
基于拾取的微震信号到时信息对两类事件进行定位,以更好地获取微震的空间分布。根据各个台站的到时延迟可大致估计这些波形均以一个较低的速度沿滑坡的表面传播,为了确定这一点,我们首先将两类不同事件的ENZ分量旋转到了滑坡面上的RTV分量(其中RTV分别为滑坡面上台站波形记录相对事件位置的径向分量、横向分量以及垂直分量),并将不同事件的R分量进行希尔伯特变换后与对应的V分量进行对比,如果事件为面波信号则变换后R与V分量的相位差会消除。由图 3(b)、3(d)可明显看出各事件的R、V分量信号的相位差经过希尔伯特变换后几乎均被消除,一定程度上说明这些事件为面波信号。此外,我们还分析了两类事件的R-V分量质点运动轨迹(图 3(c)、3(e)),结果显示两类事件的R-V分量质点运动轨迹均近似逆进椭圆,符合Rayleigh面波的质点运动轨迹特征,由此认为这两类事件信号在高通5Hz的频带范围内主要发育面波信号。在确定两类滑坡震事件主要为面波信号后,采用盛敏汉等(2018)一层模型网格搜索技术来确定事件的水平位置,图 3(a)中红色五角星为最终搜索得到的最佳位置。根据图 3中两类事件的定位结果,将波形以沿滑坡表面到各台站的距离进行排列(图 2),可以看出两类事件随着震中距的增大振幅逐渐减小,由此说明定位结果较为可靠。
(a)事件A、B的地震定位图,其中红色五角星表示事件A、B的定位位置,蓝色三角形表示台站;(b)、(c)分别为发生在D04台的事件A的希尔伯特变换图与R-V分量质点运动轨迹图,其中红色波形为V分量波形,黑色虚线与实线分别为希尔变换前后R分量波形,滤波频带范围为5~30Hz;(d)、(e)分别为发生在D12台的事件B的希尔伯特变换图与R-V分量质点运动轨迹图,其中红色波形为V分量波形,黑色虚线与实线分别为希尔变换前后R分量波形,滤波频带范围为5~30Hz |
随后,我们拟以图 2中两类波形为模板,采用滑动时窗互相关(SCC)方法(Yang et al,2009)匹配与上述两类波形相似的滑坡震事件。为了确保所选模板的波形质量以提高模板匹配的可靠性,对于A、B两类事件,分别选取了信噪比最高的D04台以及D12台作为模板进行互相关匹配并对互相关匹配的结果进行人工检查,最终筛选出46个事件A以及48个事件B。
3 滑坡上烟花信号的监测完成事件A、B的拾取与定位工作后,发现除事件A、B外,在滑坡东部还检测到一类长时间尺度位置和波形一致的重复成簇事件C(图 4(a)),该事件在D03台到时最早且波形信噪比最高,推测可能发生在D03附近。为了更好地了解重复事件的发生机理,分别从重复性、波形特征以及时间分布特征等方面对其进行了研究。
(a)事件C地震剖面图,纵坐标为旋转到滑坡面上的震中距,红色虚线为参考视速度,红色框为模板事件窗口,波形滤波范围为5~30Hz;(b)重复事件与模板事件的互相关系数及振幅比 |
首先,为研究事件C的重复性,截取D03台一段波形作为模板,如图 4(a)中红框所示,将其与随后发生的若干重复事件进行互相关计算,得到不同重复信号与模板信号的振幅比及互相关系数(图 4(b))。结果显示,所有重复事件与模板事件的互相关系数、振幅比均在0.8以上,说明该类重复事件之间具有较好的一致性,猜测其可能为由人类活动(如道路打钻、施工或燃放烟花等)产生的一致性较高的重复事件。此外,由图 4(a)可知事件C传播距离超过200m就很难被地震仪记录到,考虑到该滑坡道路主要为土路,无道路施工现象,同时确认在梭坡滑坡内D03附近处仅存在一户人家(图 5中黑色圆圈),而梭坡滑坡外距离滑坡内最近的住户约300m,说明该重复事件主要由D03附近的住户产生。
(a)事件C地震定位图,其中红色五角星表示事件C的定位位置,蓝色三角形表示台站,黑色圈为滑坡内唯一住户;(b)、(c)分别为发生在D03台事件C的希尔伯特变换图与R-V分量质点运动轨迹图,其中红色波形为V分量波形,黑色虚线与实线分别为希尔变换前后R分量波形,滤波频带范围为5~30Hz |
其次,对事件C进行定位,获取其波形特征,结果如图 5 (a)所示。根据得到的定位结果,将D03的三分量波形ENZ分量旋转到滑坡面上的RTV分量,并对R分量进行希尔伯特变换,图 5 (b)表明变换后的R分量与V分量的相位差几乎消除。同时,对事件C的R-V分量质点运动轨迹进行分析(图 5 (c)),结果显示事件C的R-V分量质点运动轨迹也近似逆进椭圆。因此事件C波形在高通滤波5Hz也主要发育面波信号,结合其重复性特征,认为其可能为烟花爆炸引起的空气震动耦合地面运动产生的微震信号(Ishihara et al,2012)。
采用图 4(a)中的模板事件,通过滑动时窗互相关方法对所有波形进行搜索来获取重复事件C的时间分布特征。最终统计得到638个重复事件C(图 6(a)),通过其时间分布图(图 6(a)中红色实线为发生时间段)可知,事件C发生的日期为2017年2月22日、2017年2月27日、2017年3月3日以及2017年3月17日,并且这4天的地震事件较未发生事件C的日期(2017年2月23日)数量明显增多,且主要集中在每天的北京时间中午12点至下午6点间。由于该滑坡区域主要住户为藏民,为此我们首先查找了藏历,发现2017年2月27日为藏历新年,而根据藏族传统,此时为藏民庆祝的时间。另外,我们将得到的互相关函数包络相邻最大振幅的时间差作为事件发生的间隔进行统计,得到图 6(b)所示事件C发生间隔的柱状图,发现大部分事件C的两两时间间隔分布在1.5~2.5s区间内,这一现象与烟花燃放时的爆炸时间间隔吻合。由此推测该事件可能是人工燃放烟花产生的信号。
(a)重复事件C的时间分布图,红色实线表示重复事件发生的时间段;(b)重复事件C的时间间隔图 |
本文识别了梭坡滑坡上的震动事件并利用滑动互相关方法检测了相似事件,通过网格搜索的方法对信噪比高的事件进行了定位。分析滑坡体微震事件的地震学性质能够为理解滑坡的动力学过程、进行滑坡稳定性分析以及灾害防治提供关键信息,因此我们分别从波速以及频谱等地震学性质对拾取到的微震信号进行讨论。
高频面波波速能够反映滑坡浅层沉积层的物质空间分布信息(曾求等,2020)。根据一层模型网格搜索技术获得的微震定位结果,计算了微震位置到不同台站沿滑坡表面的震中距,将各事件按照震中距排列并剔除信噪比较低的信号(图 2(a)、4(a))。通过微震事件的位置和到时信息可以估算出各事件信号的面波视速度,结果显示,不同事件间的面波波速存在一定差异,其中事件A速度约为250m/s,事件B约为375m/s,事件C约为250m/s。这些速度差异体现了滑坡浅层沉积层的物质的横向分布不均匀,事件A和C较事件B速度偏低,这与我们在野外布台时所观察到的地质环境相对应。事件A与C发生在滑坡中段内部,浅层沉积层为较软的风化土壤;事件B发生在滑坡东部边缘处,接近藏族古碉群住宅区,浅层沉积层为稍硬的风化块石夹土(陈菲等,2012)。该速度结果也在一定程度上说明了定位结果的可靠性。
微震事件的频谱特征同样也包含了地震波传播过程中的介质信息(Fan et al,2017)。为了削弱衰减对频谱的影响,我们选取了震中距均约为100m的事件A、B、C进行频谱分析,如图 7所示。其中,事件A与事件C的频谱特征相似,主频均约为20Hz,因此并不能简单地通过波形频谱特征分析将烟花事件与滑坡震事件进行区分。相比之下事件B主频较高(约30Hz)。结合三类事件的面波视速度信息可知,主频较低的事件A与事件C对应的视速度也较小(均约250m/s),而主频较高的事件B对应的面波时速度较大(约为375m/s)。这可能是由于事件A与事件C均发生于位于滑坡内较为松散的风化土壤层,高频面波衰减较快,而事件B发生于滑坡边缘较硬的风化块石夹土层,大于20Hz的高频面波信号得以保留。由此表明,对于微震事件的频谱分析能够获得滑坡体浅层介质的物理性质分布特征。
此外,通过对事件C的重复性、波形特征以及事件分布特征进行分析,确认了事件C为滑坡上烟花爆炸引起的空气震动耦合地面运动产生的微震信号。这一方面有助于与滑坡内部产生的微震信号进行区分,在滑坡活动性分析时能够有效去除人工信号的干扰;另一方面信噪比较高的烟花信号可以作为人工爆炸源提供滑坡浅层介质物性参数,为研究滑坡体的浅层速度结构及其随时间变化等信息提供基础。
本文观测到了滑坡滑动变形过程中辐射出的大量地震信号,对于这些地震事件发生机理、与滑坡变形的关系及与天然断层类比的情况,则需要进一步对这些微震的机制解进行分析。同时,本文在定位过程中所使用的单层模型过于简单,但是对于浅地表区域有利于问题的简化。而对于重复的烟花事件,下一步考虑对所有事件进行定位,并结合速度结构反演研究滑坡体速度随时间的变化情况。
陈菲、邓建辉、高明忠等, 2012, 丹巴县莫洛村滑坡的地质成因与稳定评价, 岩土力学, 33(6): 1781-1786. DOI:10.3969/j.issn.1000-7598.2012.06.028 |
黄润秋, 2007, 20世纪以来中国的大型滑坡及其发生机制, 岩石力学与工程学报, 26(3): 433-454. DOI:10.3321/j.issn:1000-6915.2007.03.001 |
季伟峰、胡时友、宋军, 2007, 中国西南地区主要地质灾害及常用监测方法, 中国地质灾害与防治学报, 18(S1): 38-41. |
盛敏汉、储日升、危自根等, 2018, 四川省理县西山村滑坡运动变形过程中的微震研究, 地球物理学报, 61(1): 171-182. |
唐亚明、张茂省、薛强等, 2012, 滑坡监测预警国内外研究现状及评述, 地质论评, 58(3): 533-541. DOI:10.3969/j.issn.0371-5736.2012.03.014 |
许秦坤, 2012. 狭长通道火灾烟气热分层及运动机制研究. 博士学位论文. 合肥: 中国科学技术大学.
|
阳东, 2010. 狭长受限空间火灾烟气分层与卷吸特性研究. 博士学位论文. 合肥: 中国科学技术大学.
|
杨慧、储日升、盛敏汉, 2018, 2015山东平邑石膏矿塌陷地震震源参数测定, 地球物理学进展, 33(1): 125-132. |
殷坤龙、朱良峰, 2001, 滑坡灾害空间区划及GIS应用研究, 地学前缘, 8(2): 279-284. DOI:10.3321/j.issn:1005-2321.2001.02.010 |
曾求、储日升、曾祥方等, 2020, 利用地震背景噪声提取西山村滑坡高频面波信号, 地球物理学报, 63(5): 1830-1843. |
朱良峰、吴信才、殷坤龙等, 2004, 基于信息量模型的中国滑坡灾害风险区划研究, 地球科学与环境学报, 26(3): 52-56. DOI:10.3969/j.issn.1672-6561.2004.03.011 |
Fan G, Zhang L M, Zhang J J, et al, 2017, Energy-based analysis of mechanisms of earthquake-induced landslide using Hilbert-Huang transform and marginal spectrum, Rock Mech Rock Eng, 50(9): 2425-2441. DOI:10.1007/s00603-017-1245-8 |
Gomberg J, Schulz W, Bodin P, et al, 2011, Seismic and geodetic signatures of fault slip at the Slumgullion landslide natural laboratory, J Geophys Res: Solid Earth, 116(B9): B09404. |
Helmstetter A, Garambois S, 2010, Seismic monitoring of Séchilienne rockslide(French Alps): Analysis of seismic signals and their correlation with rainfalls, J Geophys Res: Earth Surf, 115(F3): F03016. |
Huang R Q, 2009, Some catastrophic landslides since the twentieth century in the southwest of China, Landslides, 6(1): 69-81. DOI:10.1007/s10346-009-0142-y |
Ishihara Y, Hiramatsu Y, Yamamoto M Y, et al, 2012, Infrasound/seismic observation of the Hayabusa reentry: Observations and preliminary results, Earth Planets Space, 64(7): 655-660. DOI:10.5047/eps.2012.01.003 |
Ma H C, Chu R S, Sheng M H, et al, 2020, Template matching for simple waveforms with low signal-to-noise ratio and its application to icequake detection, Earthq Sci, 33: 1-8. DOI:10.29382/eqs-2020-0001-01 |
Spillmann T, Maurer H, Green A G, et al, 2007, Microseismic investigation of an unstable mountain slope in the Swiss Alps, J Geophys Res: Solid Earth, 112(B7): B07301. |
Wang G Q, 2013, Millimeter-accuracy GPS landslide monitoring using precise point positioning with single receiver phase ambiguity(PPP-SRPA)resolution: a case study in Puerto Rico, J Geodetic Sci, 3(1): 22-31. |
Xie J, Chu R S, Ni S D, 2020, Relocation of the 17 June 2017 Nuugaatsiaq(Greenland)Landslide Based on Green's Functions From Ambient Seismic Noises, J Geophys Res: Solid Earth, 125(5): e2019JB018947. |
Yamada M, Mori J, Matsushi Y, 2016, Possible stick-slip behavior before the Rausu landslide inferred from repeating seismic events, Geophys Res Lett, 43(17): 9038-9044. DOI:10.1002/2016GL069288 |
Yang H F, Zhu L P, Chu R S, 2009, Fault-plane determination of the 18 April 2008 Mount Carmel, Illinois, earthquake by detecting and relocating aftershocks, Bull Seismol Soc Am, 99(6): 3413-3420. DOI:10.1785/0120090038 |
Zhao J, Ouyang C J, Ni S D, et al, 2020, Analysis of the June 2017 Maoxian landslide processes with force histories from seismological inversion and terrain features, Geophys J Int, 222: 1956-1976. |