金森博雄1972年首次发现震前低频波(Kanamori et al,1974),也有人称之为“形变波”“静地震”等(高金哲等,2005);随后越来越多的研究者证实了这种低频短临异常信息的存在(李世愚等,2000;王庆良等,2006;张淑亮等, 2008, 2009, 2011)。近年来,国内外一些学者对这种短临异常进行了大量探索研究(冯德益等,1984;Dieterich,1992;Rice,1993;Monastersky,1994;郭增建等,1996;刘万琴等,1999),认为该异常现象是一种较为可靠的地震短临前兆信息,有可能成为地震预报的一个突破口。模拟记录时代,定点形变学科的震前短临异常形态主要以脉冲式扰动、脉动加粗等为主(陈德福等,1993;谢觉民等,1993;王贤仁等, 1994, 2002;李瑞浩等,1995;许昭永等,2001;范桂英等,2002);数字化形变观测投入运行以来,由于仪器的分辨率以及采样率的提高,记录到的震前短临异常形态表现更加复杂。高振强等(2009)针对临汾水平摆记录的远震前异常现象进行分析,认为远场效应的异常形态以“阶跃型和脉动加粗”为主;李艳等(2013)总结了近场短临异常的波动特性,发现近场短临异常首先以连续脉动变粗以及突跳异常信号开始,逐渐发展为间断性脉冲信号;李继业等(2014)统计分析了3次中强地震前近场短临异常图像,并归纳了与地震三要素相关的经验公式;杨立明等(2015)对汶川、玉树以及青川等3次震例进行研究,发现均存在临震波动现象,具有重现性。
目前,针对这种短临异常信息的识别多数停留在原始数据直观对比上,有些微弱的前兆异常信号不能被很好地识别,而Hilbert-Huang变换可以将任意信号依据信号自身的时间尺度特征,进行多尺度时频特征的刻画,得到的时频谱有较高的时间分辨率与频率分辨率。张立等(2007)和周挚等(2008)基于HHT提取震前重力固体潮异常信息,发现震前瞬时频率特征参数的确有变化。李东生等(2013)应用HHT方法分析地电场数据,认为可以很好地提取主要周期成分;王宁等(2014)对比了小波变换、S变换以及HHT变换3种时频分析方法,结果认为,采用HHT方法分析信号组成成分特征时效果较好。基于此,本文利用HHT变换方法,对河津ML 5.2地震前临汾盆地数字化洞体形变观测资料进行分析,以期探索震前洞体形变观测资料的短临异常时频特征演化。
1 研究方法Hilbert—Huang变换(HHT)最初由科学家Huang N E(1998)等提出,主要包含2个步骤:经验模态分解(Empirical mode decomposition,EMD)与Hilbert谱分析(Hilbert spectrum analysis,HSA)。EMD可以自适应地将原始信号分解为不同频率的固有模态函数(Intrinsic Mode Function,IMF),对IMF分量进行希尔伯特变换,得到瞬时振幅以及二维时频谱。
1.1 经验模态分解找出原始信号S(t)所有的极大值点与极小值点,用3次样条函数插值拟合得到原始数据序列的上、下包络线;计算上、下包络线均值得到平均包络线m(t),用S(t)减去m(t),得到一个去掉低频的新序列ci(t)。如果ci(t)满足IMF的数学条件,定义ci(t)为IMF的一个分量,否则以ci(t)为分析信号继续重复以上步骤,直到满足要求为止。信号经EMD被分解成n个IMF分量和一个信号余量,即
$s\left( t \right) = \sum\limits_{i = 1}^n {{c_i} + {r_n}} $ | (1) |
其中,ci表示第i个IMF分量,rn为残余项。
1.2 Hilbert谱分析对原始信号进行EMD分解,得到IMF序列,利用希尔伯特变换,对每个IMF进行计算,得到瞬时振幅、瞬时频率,综合IMF所有分量的瞬时频率,得到二维时频谱及边际谱。希尔伯特变换公式为
$Y\left( t \right) = \frac{1}{{\rm{ \mathit{ π} }}}\int_{ - \infty }^{ + \infty } {\frac{{{c_i}\left( T \right)}}{{t - T}}{\rm{d}}T} $ | (2) |
自2007年数字化洞体形变观测以来,临汾盆地仅发生一次较为显著的地震事件,即2010年1月24日河津ML 5.2地震。距离此次地震较近的形变观测站有临汾和侯马地震台,震中距分别为89 km和65 km,且均位于临汾盆地。其中河津ML 5.2地震与临汾地震台同处临汾盆地西缘罗云山山前断裂带附近,侯马地震台地处临汾盆地南缘峨嵋台地北缘断裂带(图 1)。因此,研究此次地震前临汾与侯马地震台形变观测资料的短临异常信息,对地震短临预测可能有一定参考价值。
![]() |
图 1 测点及河津ML 5.2地震空间分布 Fig.1 The distribution of stations and Hejin ML 5.2 earthquake |
根据已有研究成果(王庆良等,2005),该短临异常现象出现的有利时段一般为震前7天内。为了有效提取异常信息,选取2009年12月26日至2010年1月24日临汾地震台水平摆和垂直摆倾斜仪NS测向、EW测向的分钟值数据以及侯马地震台伸缩仪和水管仪的NS测向、EW测向的分钟值数据为研究对象,分别对其进行HHT分析。
2.2 震前短临异常现象2010年1月10日临汾地震台水平摆、垂直摆同步出现低频波动变化,且在水平摆NS向叠加高频脉冲现象,1月15日后脉冲突跳现象消失,于1月24日发生河津ML 5.2地震(图 2)。以上现象揭示了河津ML 5.2地震短临异常特征,表现为不但有低频信号扰动,且伴随高频脉冲信号。
![]() |
图 2 河津ML 5.2地震水平摆与垂直摆原始观测曲线 Fig.2 The original observation curve of the horizontal pendulum and vertical pendulum for Hejin ML 5.2 earthquake |
如前人所述,具有低频波动或高频脉冲信号的短临异常现象出现时间一般为震前7天左右(王庆良等,2005),为了有效提取异常信号以及考虑计算量,本文选取河津ML 5.2地震前30天的分钟值观测数据进行分析。首先进行数据预处理,通过查看观测日志去除人为干扰数据,对少量干扰数据采用3次样条插值方法进行插值补全,对连续缺数或连续干扰较长时间段的观测资料直接舍弃;为了凸显高频短临异常信息,采用高通滤波方法将固体潮及以下频段滤除。利用HHT方法,对上述预处理数据进行分析,得到瞬时振幅和时频谱参数。此外,针对临汾水平摆和垂直摆在河津地震前同步出现固体潮以下频段的低频扰动信息,利用HHT变换,计算得到震前30天边际谱参数。在实际计算中,考虑到计算机硬件限制,将观测数据按5天一组进行分组处理。
3.1 瞬时振幅和时频谱参数分析利用HHT变换,对临汾、侯马地震台4种观测手段8个测项数据进行计算,得到瞬时振幅,分析发现,河津地震前侯马伸缩仪NS向、临汾水平摆NS向以及侯马水管EW向在瞬时振幅参数上存在高值异常,见图 3—图 5的(a)图(限于篇幅只展示前3阶瞬时振幅结果)。具体演化特征为:震前21天左右首先在侯马伸缩仪NS向上出现间断性高值异常,持续约10天,之后消失;震前14天左右在临汾水平摆NS向上出现间断性脉冲型信号,持续约5天,之后消失;震前2天在侯马水管EW向上出现微弱的持续性高值异常现象,直至发震,震后逐渐减弱。由图 1可知,临汾台地处罗云山山前断裂附近,其走向为NNE,侯马台与震中连线近东西向,因此从测点、震中以及断裂走向的空间分布特征推测,临汾台监测仪器对来自震源的信息响应能力应表现为NS向优于EW向,而侯马台监测仪器对来自震源的信息响应能力应表现为EW向优于NS向。实际计算结果表明,临汾台水平摆和侯马台水管倾斜仪在测向的响应能力上与理论分析较为一致,而侯马台伸缩仪在测向响应能力上与理论分析结果有所矛盾。对此,对侯马台伸缩仪NS向瞬时振幅异常信度进行回溯分析,将数字化观测以来的分钟值数据进行高通滤波,得到高频信息时序曲线,见图 3(b),发现2011年之前侯马伸缩仪NS向背景噪声比较大,2011年更换传感器之后背景噪声明显减小,认为侯马伸缩仪NS向在2010年1月24日河津ML 5.2地震前瞬时振幅参数上的高值异常与仪器传感器性能有一定关系,因此侯马伸缩仪NS向瞬时振幅高值异常信度较低。
![]() |
图 3 伸缩仪NS瞬时振幅Am与高通滤波时序曲线 Fig.3 The instantaneous amplitude and sequence diagram after high-pass filtering of NS extensometer |
![]() |
图 4 水平摆NS瞬时振幅Am与时频谱 Fig.4 The instantaneous amplitude and time-frequency spectrum of NS horizontal pendulum |
![]() |
图 5 侯马水管EW瞬时振幅Am与时频谱 Fig.5 The instantaneous amplitude and time-frequency spectrum of EW water tube tiltmeter at Houma station |
对临汾水平摆NS向与侯马水管EW向瞬时振幅高值异常时段进行时频谱计算,结果显示:水平摆NS向存在频带较宽的不连续高频成分,且高频成分能量相对较高,见图 4(b);侯马水管EW向存在频带较宽的连续高频成分,但高频成分能量相对较低,见图 5(b)。由于水平摆的固有频率比水管倾斜仪高,因此水平摆记录高频信息的能力强于水管倾斜仪,上述实际监测结果与理论分析较为一致。
以上研究结果表明,震前高频信息的传播和接收与监测仪器的频率响应、测点附近断裂走向以及测点与震中连线的方位有关。
3.2 边际谱参数分析由临汾水平摆和垂直摆边际谱参数可知,边际谱图像在日波和半日波频率范围存在明显的优势频率,低于日波频段也存在明显的低频扰动信息,且在NS分量和EW分量上的优势频率分布较为一致(图 6),垂直摆与水平摆有类似结果,说明可能记录的是同一低频异常信号。但这种低频扰动信息只在摆式倾斜仪上有所记录,且二分量均有记录,该现象与高频信息的传播和接收特征有明显区别,可能揭示了低频异常信息的传播路径较为复杂,且接收主要受限于测点仪器的频率响应特征。
![]() |
图 6 水平摆边际谱 Fig.6 The marginal spectrum of the horizontal pendulum |
综合上述分析结果得出:地震前兆短临异常的频率范围较宽,有高于固体潮频段的高频信息,也有低于固体潮频段的低频扰动信息,地震前低频扰动信息主要出现在水平摆和垂直摆观测曲线上,而高频异常信号主要体现在水平摆NS向和水管倾斜仪EW向上,随着发震时刻的逼近,震源区激发的短临异常信息以高频成分为主。
4 讨论与结论 4.1 讨论上述结果表明,不同观测仪器、同一观测仪器不同测向所记录的短临异常起始时间、频率成分以及携带能量存在的差异性,可能与不同监测仪器的频率响应、地震波传播路径以及测点位置的基岩完整性等有关。根据断裂力学中波的传播理论(冯德益,1984),波在传播过程中具有明显的频散现象,因此,随着震中距的增大,短临异常信息中携带的高频信息逐渐减弱,若时频谱包含相对较高的频率成分,则间接表明可能是近场信息。
国内外许多实验和震例研究结果均证实粘滑前预滑的存在,而这种震前短临异常信息可能是地震发生前的预滑造成的。非均匀断层滑动失稳成核过程的实验结果(马胜利,2003)也表明:粘滑失稳过程中产生2种不同性质的弹性波,即低频、低能量的弹性波从断裂上部向下部传播及高频、高能量的弹性波从断裂向两侧传播。
4.2 结论本文在已有研究基础上,采用HHT变换方法,对河津ML 5.2地震前30天的临汾盆地洞体形变观测数据进行处理,从瞬时振幅、时频谱以及边际谱参数的角度,对震前记录的短临异常信息进行归纳总结,初步得出以下结论。
(1) 形态特征。地震前兆短临异常的传播过程具有一定频散现象,且传播路径和仪器响应频率的不同可能导致不同监测点记录的异常形态具有一定差异性。
(2) 时间特征。震前14天左右在水平摆和垂直摆同步出现低频扰动信息,且水平摆NS向伴随高频异常信息;震前2天水管倾斜仪EW向出现脉动加粗现象,反映了震源体在不同阶段所激发的异常信息随时间的演化特征。
(3) 频率特征。摆式倾斜仪首先同步出现低频扰动信号,然后水平摆NS向出现间断性高频信号,水管倾斜仪临震阶段出现连续性高频成分,反映了震源体随时间的演化,震源处所激发的频率信息由低频向高频过渡。因此,临震阶段的高频成分可能指示的是近场前兆信息,对判定是否为近震的前兆短临异常提供重要参考,倘若在同一监测点布设多种不同响应频段的观测仪器,对捕捉短临预报信息具有一定意义。
感谢山西省地震局张淑亮研究员的热心帮助和指导。陈德福, 罗荣祥, 刘国培. 地壳形变动力学观测与研究[M]. 北京: 海洋出版社, 1993, 174-187. | |
陈德福. 潮汐形变前驱波的时空特征[J]. 大地测量与地球动力学, 2006, 26(2): 24-30. | |
范桂英, 拾礼勤, 李兰生. 2001年11月14日昆仑山口西8.1级大震前徐州台数字化形变观测仪器记录的异常[J]. 地震学刊, 2002, 22(2): 7-11. | |
冯德益, 潘琴龙, 郑斯华, 等. 长周期形变波及其所反应的短期和临震地震前兆[J]. 地震学报, 1984, 6(1): 41-56. | |
高金哲, 吕政, 张洪艳, 等. 地震前驱波观测与进展[J]. 华南地震, 2005, 25(1): 53-58. | |
高振强, 冯建琴, 张仰辉, 等. 山西临汾中心地震台石英摆倾斜仪前驱波特征分析[J]. 山西地震, 2009, 25(1): 1-6. | |
郭增建, 秦保燕. 震源成核过程研究综述[J]. 国际地震动态, 1996, 2: 1-4. | |
国家地震局科技监测司编. 大同-阳高地震研究[M]. 北京: 地震出版社, 1993, 77-87. | |
李继业, 武晓军, 任建辉, 等. 低频前驱波频谱特征分析[J]. 地震地磁观测与研究, 2014, 35(Z2): 38-46. | |
李瑞浩, 刘贵安, 陈冬生, 等. 滇西试验场重力潮汐观测及其地震响应特征的研究[J]. 中国地震, 1995, 11(2): 161-167. | |
李世愚, 腾春凯, 卢振业, 等. 典型构造微破裂集结的实验研究[J]. 地震学报, 2000, 22(3): 278-286. | |
李艳, 高振强, 程冬焱, 等. 河津4.8级地震前山西临汾台倾斜前驱波特征分析[J]. 地震研究, 2013, 36(2): 160-164. | |
刘琦, 张晶. S变换在汶川地震前后应变变换分析中的应用[J]. 大地测量与地球动力学, 2011, 31(4): 6-9. | |
刘琦, 张晶, 池顺良, 等. 2013年芦山MS 7.0地震前后姑咱台四分量钻孔应变时频特征分析[J]. 地震学报, 2014, 36(5): 770-779. | |
刘万琴, 李世愚, 郑志真, 等. 破坏性矿震震前短临阶段震源过程研究[J]. 地震学报, 1999, 21(1): 57-64. | |
马胜利, 刘力强, 马瑾, 等. 均匀和非均匀断层滑动失稳成核过程的实验研究[J]. 中国科学, 2003, 33(Z1): 45-52. | |
王庆良, 张晓东, 崔笃信, 等. 理解前兆异常变化机理和地震短临前兆[J]. 国际地震动态, 2006(5): 131-144. | |
王宁, 吴云, 张燕. 时频分析方法在形变数据中的应用研究[J]. 地震工程学报, 2014, 36(2): 413-420. | |
王贤仁, 卢艳清. 澜沧-耿马地震前FSQ仪记录到的粘滑和慢地震及对短临预报探讨[J]. 地震研究, 1994, 17(2): 157-162. | |
王贤仁, 卢艳清. 倾斜扰动与短临地震预报[J]. 地震研究, 2002, 25(1): 20-24. | |
许昭永, 王贤仁, 黄雍, 等. 倾斜扰动——一种可能的有效地震前兆[J]. 地震研究, 2001, 24(2): 131-135. | |
杨立明, 郝臻, 王建军, 等. 汶川、玉树地震临震波动现象的震例研究[J]. 国际地震动态, 2015(9): 62 | |
张淑亮, 程紫燕, 唐垒黎, 等. 井水位前驱波现象与震源成核之间关系的研究[J]. 中国地震, 2008, 24(2): 159-171. | |
张淑亮, 刘瑞春, 宁亚灵, 等. 汶川MS 8.0地震前山西前兆低频前驱波特征分析[J]. 大地测量与地球动力学, 2009, 29(6): 35-39. | |
张淑亮, 吕芳, 李斌. 井水位前驱波普遍性问题的探讨[J]. 中国地震, 2011, 27(4): 348-362. | |
张立, 傅荣珊, 周挚, 等. 基于HHT提取重力固体潮的地震前兆信息[J]. 地震学报, 2007, 29(2): 222-226. | |
周挚, 山秀明, 张立, 等. 基于HHT提取昆明、下关重力固体潮的地震前兆信息[J]. 地球物理学报, 2008, 51(3): 836-844. | |
Dieterich J H. Earthquake nucleation on faults with rate and state dependent strength[J]. Tectonophysics, 1992, 211: 115-134. DOI:10.1016/0040-1951(92)90055-B | |
Kanamori J, Cipar J J. Focal process of the great Chilean earthquake May 22, 1960[J]. Phys Earth Planet Inter, 1974, 9: 128-136. DOI:10.1016/0031-9201(74)90029-6 | |
Huang N E, Zheng S, Long S R. The empirical mode decomposition and the Hilbert spectrum for non-1iner and non-stationary time series analyses[J]. Proe R Sco Lond, 1998, 545: 903-995. | |
Monastersky R. Before the earthquake:detecting the slow groan[J]. Science News, 1994, 146(23): 374 | |
Rice J R. Spatio-temporal complexity of slip on a fault[J]. J Geophys Res, 1993, 98: 9 885-9 908. DOI:10.1029/93JB00191 |