基于地震动信号分析的滑坡堵江过程研究——以白格滑坡堰塞湖为例

林 子 钰,范 刚,2,陈 骎,周 家 文,2

(1.四川大学 水利水电学院,四川 成都 610065; 2.四川大学 水力学与山区河流开发保护国家重点实验室,四川 成都 610065)

0 引 言

滑坡是最常见的地质灾害之一,滑坡体堵截河床或河谷后会形成堰塞湖,堰塞湖一旦发生溃决,将会引发更为严重的次生灾害,中国西南山区是滑坡堰塞湖灾害链发生的重灾区[1-2]。系统揭示滑坡运动过程和堰塞坝溃决过程对灾害监测预警和灾后救援工作都具有十分重要的现实意义,但由于滑坡堰塞湖灾害往往发生在人烟稀少的山区,具有明显的突发性和隐蔽性,导致应急救援部门在灾害发生的时间和地点识别上具有滞后性。例如2018年发生在西藏的加拉滑坡和白格滑坡,均发生在深夜且地处无人区,现有监测手段均无法给出其灾害发生的具体时间[3-4]。滑坡堰塞湖灾害过程的动态信息难以掌握是制约人们研究其成生过程和及时制定应急抢险方案的关键因素。

近几十年来,随着环境地震学的发展,从地震动信号本身出发来识别滑坡堰塞湖的产生过程逐渐成为一种新的方法。灾害发生过程中产生的地震动信号是目前少有的能直观反映其在运动过程中相关信息的数据,目前中国已经建成了覆盖全国的地震监测台网,为地震动信号分析提供了充足的数据来源[5-6]。近些年,学者们在地震动信号频谱分析方面已开展了诸多研究,如陆菜平等[7]利用短时Fourier变换分析了一处矿区的顶板和煤层卸压爆破微震动信号的频谱特性,但该方法无法准确描述此类信号的局部特征;赵国彦等[8]尝试利用频率切片小波变换技术对典型的矿山岩体微震信号和爆破振动信号进行时频分析,但该方法仍然受制于频率切片函数的选择。随着新算法的不断涌现,Huang 等提出了Hilbert-Huang变换方法(HHT),该方法在分析和处理非线性非平稳的瞬态信号方面较传统方法展现出了极大的优越性[9-10]。Fan等[11-12]借助HHT方法在处理地震波这一随机信号方面的优势,对地震作用下含软弱夹层岩质边坡的震损过程进行了识别;Song等[13]基于HHT方法,分析了含不连续面岩质边坡地震破坏机理;徐强等[14]则利用该方法对汶川地震及其余震的地磁数据进行了分析。HHT作为一种新兴的时频分析方法,为研究滑坡堰塞湖成生过程提供了有效的技术手段。

在基于地震动信号的滑坡灾害识别方面,目前国内外学者也已进行了诸多尝试。Vilajosana等[15]的研究证实了地震方法监测滑坡的可行性和有效性,并准确地确定了滑坡的位置和规模。Chen等[16]利用实时记录的宽频带地震信号快速确定了发生在台湾的12次岩崩的时间和位置,并通过对高频数据进行频谱分析和阶段划分来描述岩崩的发生过程。赵娟等[17]基于地震动信号分析揭示了千将坪高速滑坡启动后的滑体运动特征。此外,众多学者通过震相分析,揭示了2017年发生的茂县新磨滑坡的运动过程[18-20]。虽然现有研究在利用地震动信号识别滑坡成灾过程方面进行了一些有益的探索,取得了一些实践经验和科学认识,但是在利用地震动信号识别滑坡-堰塞湖灾害链形成过程方面的研究相对缺乏。

为此,本文选取了具有代表性的2018年白格滑坡堰塞湖作为研究对象,基于台站实测地震动信号,采用HHT信号分析方法获得白格滑坡运动和堰塞湖溃决过程所激发的地震动信号的频谱特征,进而精准还原了白格滑坡堵江和堰塞坝溃决过程,以为今后此类滑坡堰塞湖灾害的精准识别提供技术参考。

1 研究场地和数据

1.1 白格滑坡堰塞湖概况

2018年10月10日和11月3日,在西藏自治区江达县波罗乡白格村金沙江右岸先后发生两次大型山体滑坡,堵塞金沙江上游干流河段,形成堰塞湖。此处河谷呈“V”形谷,滑坡区边坡基岩为元古界雄松群片麻岩组,岩性为浅灰、深灰色黑云斜长片麻岩、角闪片岩、灰绿色蛇纹岩以及夹白云石英片岩,片理产状变化较大。区内构造活动强烈,地震多发,岩体完整性差[4,21]。

白格“10·10”滑坡是一个高位、高剪出口、高速非完全楔形体基岩滑坡,滑坡物源区产生约2.1×107m3的岩土体,在重力作用下失稳掉落并刮铲裹挟约3×106m3的岩土体,形成的堰塞体顶宽约270～320 m,顺河向长约1 330 m,堰塞坝总体左侧高,右侧低,平均坡度约38°。在堵塞金沙江近2 d后,堰塞湖蓄水量达到峰值,约2.9亿m3,湖水漫顶溢流后自然泄洪,逐渐冲刷形成泄流槽,至10月13日18:00左右,险情得到有效控制。

白格“11·3”滑坡是牵引区的部分岩土体在起阻滑作用的碎裂岩体渐进解体后下滑的结果,滑源物质体积约3.7×106m3,沿第一次滑坡的运动路径向下刮铲带走约4.8×106m3的松散堆积体,堵塞原堰塞坝泄流槽并形成规模更大的堰塞体。形成的堆积体横河向最大宽度约300 m,顶部顺河向长约195 m,底部顺河向长度约300 m,颗粒总体较细,与第一次堆积物接近。在堵塞江水近9 d后,堰塞湖蓄水量达到峰值(5.24 亿m3),经人工干预开挖导流槽泄流后,至11月15日08:00,险情缓解[22-24],如图1所示。两次白格滑坡和堰塞坝溃决的具体过程如表1所列。

表1 白格滑坡堰塞湖灾害过程

图1 2019年应急处置后的白格滑坡堰塞湖现场

1.2 地震动数据

本次研究使用的地震动数据均从国家测震台网数据备份中心(http:∥www.seisdmc.ac.cn)获取,共计选取了灾害发生点周围10个测震台站记录到的地震动信号进行分析,各台站具体信息详见图2。在进行地震动信号分析前,首先需要利用SAC对测震台站记录的原始地震动信号进行必要的预处理:为降低信号干扰,首先需要去除台站仪器响应;同时为便于波形识别,本研究将地震仪记录到的东西分量(E)和南北分量(N)的地震动信号旋转到灾害发生点到各测震台站点的径向(R)和切向(T)上。

2 原始地震动信号频谱分析

2.1 HHT基本理论和方法

信号分析技术的基础是Fourier变换(FT),其本质是将信号分解成不同频率的谐信号叠加,但该方法受分析的信号必须是严格周期或平稳的所限。事实上,滑坡堰塞湖灾害产生的地震动信号是一种典型的非线性、非平稳信号,FT方法并不适用,而在其基础上发展出的一系列新的信号分析方法如短时Fourier变换、Wigner-Ville分布、小波变换等,都没有从根本上摆脱FT理论的束缚,仍然需要假定信号在一个很短的时间窗内是伪平稳的。Huang等[9]提出的Hilbert-Huang变换很好地规避了这一问题,可以直接处理非平稳和瞬态问题。

由于非平稳信号的频率是随时间变化的,因此HHT方法中引入了“瞬时频率”的概念。瞬时频率揭示了某一时刻信号能量在频率集中程度的一个度量,基于此而对信号所施加的限制条件就构成了HHT方法的核心步骤——经验模态分解(EMD)。

经验模态分解的本质是对一个信号进行平稳化处理,其结果是将信号分解成一系列具有不同特征尺度的数据序列,每个序列称为一个固有模态函数(IMF),即原始信号X(t)可由n阶IMF分量和一个残余值rn(t) 构成:

再对单个IMF分量ci(t)进行Hilbert变换并叠加后即可得到原始信号的Hilbert时频谱,这是一种联合幅值-频率-时间的三维分布图像,其物理意义是坐标中的某个点代表着在整个时间历程中的该时刻,局部出现了此频率的波动。再将时频谱对时间轴进行积分可得到边际谱,Hilbert边际谱表征了信号能量幅值在频率轴上的分布,其物理意义是:在某一频率上存在着能量意味着具有该频率的波在信号的整个持续时间内某一时刻出现的可能性较高,即边际谱具有一定的概率意义[9]。

2.2 第一次白格滑坡地震信号频谱分析

以“10·10”白格滑坡产生的垂向地震动信号为例,各测震台站记录的原始地震动信号如图3所示,其中信噪比效果最好的是GZI台站。本文选取GZI台站在滑坡发生前后共400 s的数据进行EMD分解,可以得到7阶IMF分量以及1个残差,如图4所示。

图4 EMD分解结果

将GZI台站垂直分量的单个IMF分量进行Hilbert变换并叠加后得到三维Hilbert时频谱,如图5所示。图5表明,信号频率以低频为主,集中出现在110～220 s之间。

图5 第一次滑坡激发的垂向地震动信号的Hilbert时频谱

将时频谱对时间轴积分后得到Hilbert边际谱,如图6所示。滑坡产生的地震动信号频率主要出现在0.001 2～0.376 2 Hz范围内,其中0.008～0.06 Hz为主要组成成分;垂向上的信号峰值出现在0.017 Hz,径向和切向上的峰值均出现在0.035 Hz,说明该频率的信号在滑坡的整个持续时间内出现的概率最大。

图6 第一次滑坡的Hilbert边际谱

2.3 第二次白格滑坡地震信号频谱分析

以同样方法对“11·3”白格滑坡地震动信号进行分析,以GZI台站垂直分量为例,其Hilbert时频谱如图7所示。两次滑坡产生的信号的频谱特征基本一致,信号集中出现在140～240 s之间,但由于第二次滑坡的规模较小,因此信号的信噪比更低。

图8为GZI台站的Hilbert边际谱,从图中可以看出信号频率主要出现在0.001 2～0.379 9 Hz范围内,其中0.008～0.06 Hz为主要组成成分;垂向上的信号峰值出现在0.009 Hz,径向上的信号峰值出现在0.016 Hz,说明该频率的信号在滑坡的整个持续时间内出现的概率最大;切向上的信号峰值不明显,在0.003～0.025 Hz范围内的信号出现概率大致相等。

3 滑坡堰塞湖运动特征分析

3.1 “10·10”白格滑坡运动特征分析

地质灾害在发生过程中会引起地表震动,产生明显的地震动信号,但不同的地质灾害产生的振动信号的频率范围不同。根据频谱分析结果,本次研究将地震动信号划分为低频地震动信号(f≤0.03 Hz)、中频地震动信号(0.03 Hz0.1 Hz)。滑坡体作为一个整体在斜坡上做加速和减速运动过程中会产生低频和中频地震动信号,而高频地震动信号的激发源一般认为与滑坡内部颗粒之间、滑体与滑动面之间的摩擦和碰撞以及弹性介质中局部小尺度的速度扰动有关[25]。

选取信噪比最高的GZI台站作为研究对象,截取滑坡发生前后共200 s的地震动信号进行分析。对原始地震动信号进行两组无限脉冲带通滤波,选择4阶双通道Butterworth滤波器,以垂直分量为例,如图9所示。第一组的滤波范围为0.008～0.030 Hz,此处未选择更低频的地震动信号,这是因为其超出了所有地震台站响应曲线的衰减范围,且噪声在这个周期段被放大从而降低了滑坡地震动信号的信噪比[25],第二组的滤波范围为0.1～1 Hz。为校正滑坡发生时间,取地震波传播速度为3.24 km/s。

图9 “10.10”GZI台站垂向地震动信号

白格第一次滑坡共持续约112 s,从图9(b)可以看出,滑坡产生的地震动信号先以低频信号的形式出现并逐渐增大。将低频信号开始出现的时刻定义为0点,随后伴随高频信号的出现迅速增大到峰值后再缓慢降低到噪声水平。结合高、低频信号的对比,可将滑坡运动分为3个阶段:高位滑坡连续滑落阶段、滑体破碎解体阶段、散落堆积阶段。

第①阶段(第0～28 s):高位滑坡连续滑落阶段,持续时间约28 s。不稳定岩土体在重力作用下加速下滑,尽管滑坡体从启动阶段就开始逐渐分解,但由于此时滑动的土石体没有在空间上被拉长,仍然作为一个整体在斜坡上运动,因此产生的地震动信号以低频为主。

第②阶段(第28～78 s):滑体破碎解体阶段,持续时间约50 s。该阶段滑坡体在下滑过程中受到侧翼山体阻挡,在空间上被拉长,滑体开始大量破碎,整体性被破坏。随着滑坡体对侧翼的不断刮铲裹挟,侵蚀区域山体逐渐失稳塌陷,滑坡体内部小尺度块体之间的碰撞随侧翼被侵蚀而显著增加,开始激发明显的高频信号并迅速增大到峰值。随后,滑体呈碎屑流状态迅速滑入金沙江,前端物质在后方物质的冲击作用下继续向对岸坡爬高[23]。

第③阶段(第78～112 s):散落堆积阶段,共持续约34 s。此时滑坡主运动基本停止,因此高频信号的振幅显著降低。部分爬高块体在重力作用下回落到金沙江,未冲上对岸的物质沿金沙江上下游滑动并最终停止运动,地震动信号逐渐恢复到噪声水平。

3.2 “11·3”白格滑坡运动特征分析

图10为GZI台站垂直分量上的地震动信号,滤波范围同上。白格第二次滑坡持续时间为102 s,同样以低频信号开始出现的时刻定义为0点,比较图10(b)、(c)可以看出第二次滑坡产生的低频信号的振幅明显比高频信号小,但两次滑坡产生的地震动信号的特征基本相同。将滑坡运动分为3个阶段:滑坡启动阶段、滑体破碎解体阶段、散落堆积阶段。

第①阶段(0～24 s):滑坡启动阶段,持续时间约24 s。原滑坡点后缘的不稳定岩土体失稳后沿第一次滑坡的滑床向下崩落,对滑床产生猛烈的冲击碰撞作用,产生明显的低频信号,信号振幅持续增大。

第②阶段(24～68 s):滑体破碎解体阶段,持续时间约44 s。滑体裹挟第一次滑坡残留的物质后呈碎屑流状滑入金沙江,但由于第二次滑坡的动能较小,未出现向对岸陡坡爬高的现象[23],此阶段的特征是产生强烈的高频地震动信号。

第③阶段(68～102 s):散落堆积阶段,持续时间约34 s。此时滑坡堆积形态已初步形成,以小规模块体局部缓慢减速运动为主,最终停止运动并再次阻断金沙江,信号逐渐恢复到噪声水平。

3.3 第一次堰塞坝溃决运动特征分析

由于GZI台站在这段时间内有部分信号缺失,因此选择效果更好的BTA台站进行分析,以垂直地震动信号为例,如图11所示,选定的滤波范围为0.008～0.03 Hz。

根据地震动信号分析结果,结合调查资料,可将白格堰塞坝溃决过程分为4个阶段:过流孕育阶段、侵蚀阶段、溃口加速扩展阶段以及水沙再平衡阶段[26]。

第①阶段:过流孕育阶段。此阶段以10月10日22:07滑坡阻断金沙江为起始,至10月12日17:30溃口缓慢发展至堰塞体自然溢流为止,共持续约43.5 h。期间堰塞坝迎水坡水位逐渐上升,由于水流对坝体的侵蚀作用较小,水流挟沙能力弱,因此无明显的地震动信号。

第②阶段:侵蚀阶段。此阶段以12日17:30溃口贯通为起始,至13日00:45堰塞湖水位回落为止,共持续约7 h,如图11(b)所示。在溃口刚刚贯通的一段时间内,坝体顶部的上游来水贴于溃口底部徐徐流出,下泄时会携带背水坡坝顶位置的颗粒物质到背水坡下游,但尚不足以携带大颗粒物质,与该阶段前期无明显的振动信号的情况相符;由于此时天然流量仍然大于下泄流量,堰塞湖水位持续上升,且随着迎水坡水位不断上升,下泄流量也会不断增大,在这段时间内溃口共发生两次较大规模的坍塌和两次小规模的局部坍塌,产生4 段明显的振动信号。

第③阶段:溃口加速扩展阶段。此阶段以13日00:45堰塞湖水位回落为起始,到13日09:30溃口形态基本不再变化为止,共持续约8.5 h,如图11(c)所示。随着下泄流量不断增大,溃决水流的挟沙能力变强,水流不断冲刷坡脚,使堰塞坝稳定性降低。从06:00 左右开始,溃口不断发生明显的坍塌,持续产生较为明显的地震动信号,在07:00左右达到下泄流量峰值,堰前水位快速降低,在09:00左右,溃口共发生了2次大规模的坍塌,产生2段非常强烈的振动信号,与现场右岸拢口被完全冲开、坝体高度骤降的实际情况一致,此时的信号振幅在堰塞坝溃决的整个持续时间内达到最大。

第④阶段:水沙再平衡阶段。此阶段以13日09:30 堰前水位迅速下降为起始,至当日18:00坝下水位平稳为止,共持续约8.5 h,如图11(d)所示。溃决水流对坝体物质持续冲刷,期间产生了1段明显的振动信号,但由于水流的侵蚀能力逐渐减弱,信号振幅也随之减小,最后逐渐稳定在噪声水平。

3.4 第二次堰塞坝溃决运动特征分析

采用同样的滤波方法对GZI台站记录的地震动信号进行处理,以垂直分量为例,如图12所示。尽管受到人工干预,但两次堰塞坝溃决过程大致相似,仍然分为4个阶段:过流孕育阶段、侵蚀阶段、溃口加速扩展阶段、水沙再平衡阶段。

图12 GZI台站垂向地震动信号加速度时程

第①阶段:过流孕育阶段。此阶段以11月3日17:22滑坡阻断金沙江为起始,到11月12日10:50堰塞湖通过人工开挖导流槽开始过流为止。由于此次堰塞坝的规模较大,因此从8日起,通过人工开挖导流槽进行抢险工作,在8～10日期间,GZI台站记录的地震动信号会突然出现非常强烈的振动,振幅范围远远超过在自然条件下产生的振动信号,因此为便于观察,仅展示从12日2:00～10:50的地震动信号,如图12(b)所示。由于受到人为因素的干扰,堰塞坝已经不断出现局部小规模的坍塌,能够监测到明显的振动信号,期间堰塞湖水位持续上涨,以溯源侵蚀为主。在12日04:45之后,水位已经高于导流槽底坎,直至10:50白格堰塞湖正式通过人工开挖的导流槽开始过流,由于水流的冲击作用,此时溃口发生了一次较大规模的坍塌。

第②阶段:侵蚀阶段。此阶段以12日10:50人工导流槽开始过流为起始,到13日08:55导流槽被大幅冲开为止,共持续约22 h,如图12(c)所示。从图中可以观察到4段明显的振动信号,说明溃口在此期间共发生了4次较大规模的坍塌,其中振幅最大的一次大致出现在13日07:00～09:00,与现场导流槽被大幅冲开的情况吻合。现场专家曾估算,湖水从过流到冲刷导流槽形成较大缺口大约需要20 h,将形成较大规模泄流,信号振动结果与专家预测结果相符。

第③阶段:溃口加速扩展阶段。此阶段以13日08:55导流槽被大幅冲开为起始,到14日20:30导流槽被完全冲开为止,共持续约35.5 h,如图12(d)所示。这一阶段内堰塞坝共发生了8次规模较大的坍塌,其中规模最大的一次发生时间接近14日20:30,此时人工导流槽被完全冲开,堰前水位迅速下降,因此产生了强烈的振动信号。

第④阶段:水沙再平衡阶段。此阶段以14日20:30 堰前水位急速下降为起始,到15日08:00坝下水位平稳为止,共持续约11.5 h。期间坝前水位持续降低,水流的侵蚀能力逐渐减弱,不再产生强烈的振动信号,最终金沙江水位基本稳定,溃坝险情解除,信号恢复到噪声水平。

4 结 论

本文对2018年两次白格滑坡堰塞湖灾害激发的地震动信号进行了频谱特征和灾害发展过程分析,得到如下结论。

(1) 地震动信号时程分析结果表明:第一次白格滑坡的发生时间是10月10日22:05:41,结束时间是22:07:33,共持续112 s;第二次滑坡的发生时间是11月3日17:21:15,结束时间是17:22:57,共持续102 s。

(2) 两次白格滑坡激发的地震动信号的频率主要分布在0.008～0.06 Hz范围内,以低频信号为主。其中第一次滑坡垂向信号峰值出现在0.017 Hz,径向和切向峰值均出现在0.035 Hz;第二次滑坡垂向信号峰值出现在0.009 Hz,径向信号峰值出现在0.016 Hz,切向信号峰值不明显,在0.003～0.025 Hz范围内的信号出现概率大致相等。

(3) 两次白格滑坡产生的地震动信号均先以低频信号的形式出现并逐渐增大,随后出现强烈高频信号,并迅速增大到峰值后再缓慢降低到噪声水平。以强烈高频地震动信号的出现和消失为节点,可将滑坡的运动过程分为3个阶段:第一次滑坡包括高位滑坡连续滑落阶段(28 s)、滑体破碎解体阶段(50 s)、散落堆积阶段(34 s);第二次滑坡包括滑坡启动阶段(24 s)、滑体破碎解体阶段(44 s)、散落堆积阶段(34 s)。

(4) 两次堰塞坝溃决的过程均可分为4个阶段:过流孕育阶段、侵蚀阶段、溃口加速扩展阶段、水沙再平衡阶段,第一次堰塞湖漫顶后自然泄流与第二次人工引导泄流在过流孕育阶段产生的地震动信号存在较大差异。

(5) 基于测震台站实时记录的连续波形数据可以对此类地质灾害进行远程监测,直观反映灾害发生过程中的一些动态信息,可为滑坡堰塞湖灾害的精准识别和灾后救援工作提供参考,也为后续深入揭示滑坡堰塞湖灾害激发低频地震动信号的物理源机制及灾害过程动力特征参数反演奠定了基础。