地应变、地倾斜等定点形变观测主要为了监测固体潮汐的变化,其测量指标能直接反应地壳介质的微动态变化,进而捕捉到地壳介质破裂前的力学变化信息,适用于监测中短和短临阶段的地震前兆(张雁滨等,2001、2002;孙伶俐等,2013)。但在形变观测过程中,不可避免地受到各种因素的干扰,因此如何准确识别干扰或前兆异常显得尤为重要(刘建明等,2016;狄樑等,2017;赵莹,2018)。
时频分析方法是研究非平稳随机信号的有效工具,通过时频分析,可以在时频空间得到信号的频率、强度随时间的分布(Stockwell et al,1996;张燕等,2004;周挚等,2005;陈学华,2006)。目前,时频分析方法已逐渐应用于形变数据分析中,例如,吕品姬等(2011)研究了小波分解STFT方法在地形变观测数据中的应用;王宁等(2014、2015)研究了时频分析方法在形变数据中的应用和形变数据中天然扰动在时频域的响应特征;戴勇等(2013、2016)将时频分析方法应用于包头台形变和地震前兆数据分析;方燕勋等(2019)研究了S变换在定点形变观测中的应用。
本文用S变换,对内蒙古中部地区受自然环境干扰和区内、区外地震波影响(近震、远震)的形变数据时频响应特征进行分析,以便进一步认识形变数据的干扰特征,为识别内蒙古中部地区地震前兆的异常信息提供更加客观的参考。
1 方法原理S变换以Morlet小波为基本小波,是连续小波变换的延伸。信号s(t)的S变换定义为
$ s(\tau, f)=\int_{-\infty}^{\infty} s(t) \frac{|f|}{\sqrt{2 \pi}} \mathrm{e}^{-\frac{(\tau-t)^{2 / 2}}{2}} \mathrm{e}^{-i 2 \pi f t} \mathrm{d} t $ | (1) |
式中,τ和f分别表示时间和频率,均为实数。基本小波由简谐波和高斯函数的乘积构成,定义为
$ w(t, f)=\frac{|f|}{\sqrt{2 \pi}} \exp \left(\frac{-t^{2} f^{2}}{2}-2 \pi i f t\right) $ | (2) |
S变换采用宽度可变的高斯函数,其小波基函数可随频率变化而自动调节分析时宽,在低频段的时窗较宽,可获得较高的频率分辨率,而高频段的时窗较窄,可获得较高的时间分辨率。它既解决了短时傅立叶变换中时窗大小不能调节的问题,又具备了小波变换的多分辨率特性,同时还避免了小波变换无法与傅立叶变换保持联系的问题,这些特点使S变换广泛应用于非平稳信号分析(姚家骏等,2011;王宁等,2014)。
2 不同干扰因素的时频响应特征在内蒙古中部地区定点形变观测中,常见的干扰有自然环境(降雨、大风、气压)干扰和地震波影响。本文用S变换,选取该地区受干扰突出且具代表性的呼和浩特台、乌加河台和包头台的形变数据进行分析,分别总结受不同干扰因素影响的形变观测数据的时频域特征。
文中时频图的频率为归一化频率fn,即实际频率f与采样频率fs的比值
$ f_{n}=f / f_{s} $ | (3) |
归一化频率fn的范围为0~0. 5(戴勇等,2016)。时频图中的色标由相应频率和时间对应的S变换结果转换得到,与原始观测数据中不同频率成分振幅值的变化特征相一致,其反映了观测数据中各频率信号能量随时间的变化和相互之间能量的对比(万永革,2012;刘学谦等,2015)。
2.1 自然环境干扰 2.1.1 降雨干扰在实际观测中,呼和浩特台的水管倾斜仪、洞体应变仪存在明显的降雨干扰。选取该台站水管倾斜仪NS分量2018年整时值数据以及洞体应变仪NS分量2018年8月29日~9月6日分钟值数据进行分析,其中水管倾斜仪NS分量5~6月、7~8月和9~10月数据受雨季降雨影响,固体潮出现3次加速下降变化,洞体应变仪NS分量2018年9月1日12:00~2日18:00受强降雨干扰,固体潮曲线出现转折性快速上升变化,如图 1(a)、(b)、(d)、(e)所示。
![]() |
图 1 降雨干扰时频图 (a)水管仪原始曲线;(b)降雨量;(c)水管仪S变换图;(d)洞体应变仪原始曲线(8月29日~9月6日);(e)降雨量(8月29日~9月6日);(f)洞体应变仪S变换图 |
由时频图 1(c)、(f)可看出:①水管倾斜仪全年整时值数据在归一化频率0.083和0.042处,即周期约为12h和24h处,有两处明显的周期性波块,分别是半日波、日波;②水管倾斜仪整时值数据5~6月、7~8月和9~10月在归一化频率0~0.1频段内、洞体应变仪分钟值数据2018年9月1~2日在归一化频率0~0.025频段内均存在干扰信号,与原始信号受降雨干扰时间段一致,且干扰信号的能量强度高于正常背景值;③降雨影响表现为三角状分布,表明该干扰信号的频率先增大后减小。
降雨影响常会导致形变观测曲线整体线性趋势发生转折或上升(下降)速率发生变化(赵莹,2018)。由图 1(a)、(b)、(d)、(e)可看出,观测曲线上升(下降)速率随着持续降雨的影响逐渐增大,随着降雨的结束,其速率逐渐减小,观测曲线恢复稳定,该过程与时频域频率变化过程一致。
王宁等(2015)认为对于整时值形变数据,降水量较多时会产生周期大于日波的干扰信号。对于本文整时值水管倾斜仪数据,其受降雨的影响既有周期大于日波的干扰信号,如5~6月,也有周期小于日波的干扰信号,如7~8月和9~10月。降雨对形变观测的影响与降雨时间、降雨量及降雨过程等有密切关系(汪翠枝等,2010),对比3次降雨量(图 1(b))可知,3次降雨的持续时间、累积降雨量、降雨过程均不同,因此干扰信号周期的不同可能与降雨量的多少、降雨持续时间的长短等引起观测数据变化趋势、速率的不同有关。呼和浩特台洞体应变仪分钟值数据受降雨影响,其干扰信号的归一化频率小于0.025,即周期大于40min。
2.1.2 大风干扰包头台、乌加河台水管倾斜仪受季节性大风干扰明显,选用包头台水管倾斜仪NS分量2018年5月25~26日和乌加河台水管倾斜仪NS分量2018年3月13~15日的分钟值数据,其中包头台水管倾斜仪NS分量2018年5月25日21:21~26日7:05、乌加河台水管倾斜仪NS分量2018年3月14日07:13~15日01:13因大风干扰导致观测曲线加粗,出现毛刺现象,其它时段观测正常,如图 2(a)、(c)所示。
![]() |
图 2 大风干扰时频图 (a)包头台水管仪原始曲线;(b)包头台水管仪S变换图;(c)乌加河台水管仪原始曲线;(d)乌加河台水管仪S变换图 |
由时频分析结果(图 2(b)、(d))可看出:①在与时域对应的干扰时段内出现高频扰动信号,最大归一化频率约为0.4;②干扰时段的能量强度明显高于正常时段,干扰结束后能量强度逐渐恢复至正常水平,且能量强度与观测曲线受风扰影响变化幅度成正比,幅度变化越大,能量强度越高,例如包头台水管倾斜仪NS分量2018年5月26日0时、6时的观测曲线变化幅度明显大于其它干扰时段,其时频域对应能量集中在色标值0.1左右,高于其它时段。
对比大风干扰的时频域和时域特征,时频域的频率大小和能量强度随横轴时间变化过程,可很好地反映出时域观测曲线受风扰影响的过程,且大风干扰对形变数据的影响集中在高频段,能量强度与观测曲线受风扰影响变化幅度成正比。
2.1.3 气压干扰呼和浩特台、乌加河台洞体应变仪易受气压干扰,选用呼和浩特台洞体应变仪NS分量2018年6月12日和乌加河台洞体应变仪NS分量2018年6月24日的分钟值数据进行分析,其中呼和浩特台洞体应变仪9~13时、乌加河台洞体应变仪16~19时观测曲线受气压短周期骤变影响,出现明显的固体潮畸变,其它时段观测正常,如图 3(c)、(d)、(e)、(f)所示,其时频结果见图 3((a)、(b))。
![]() |
图 3 气压干扰时频图 (a)呼和浩特台洞体应变仪S变换图;(b)乌加河台洞体应变仪S变换图;(c)呼和浩特台洞体应变仪原始曲线;(d)呼和浩特台气压值;(e)乌加河台洞体应变仪原始曲线;(f)乌加河台气压值 |
由时频分析结果(图 3(a)、(b))可以看出:①在与时域对应的干扰时段内有扰动信号,其中图 3(a)扰动信号的最大归一化频率约为0. 2,图 3(b)扰动信号的最大归一化频率约为0.15;②干扰时段的能量强度明显高于正常时段,且最高能量区域集中在归一化频率0~0.05频段内;③干扰信号的时频域能量变化形态与时域变化形态具有一致性。
综合本区域气压干扰的时域、时频域特征,对于分钟值形变数据,在气压短周期骤变的地方,观测数据在时域会同步出现形态类似的变化;在时频域上,干扰信号的最高能量区域集中在归一化频率0~0.05内,即干扰信号的归一化优势频率分布在0~0.05之内,该结果与王宁等(2015)获得的结果(对于分钟值形变数据,在气压变化处,观测数据有相同的变化趋势,在时频域上,干扰信号的归一化优势频率低于0.05)具有一致性,且本分析区域内干扰信号的最大归一化频率可达到0.2左右。
2.2 地震波影响形变仪器受地震影响时,大部分仪器能记录到地震波形,也有部分仪器会出现阶变等同震响应(赵莹,2018)。而地震对形变数据的影响相较其正常背景下的观测数据而言,也可看做是一种干扰信号(李继业等,2015)。
选取乌加河台秒采样VP垂直摆NS分量记录的2017年6月3日内蒙古阿拉善左旗M5.0地震和2018年9月8日云南普洱墨江县M5.9地震的数据进行时频分析,并与同台JCZ地震计记录的同时段地震波数据的时频响应特征进行对比,具体地震参数见表 1。
![]() |
表 1 地震参数 |
VP垂直摆和JCZ地震计记录的M5.0、M5.9地震波的原始曲线见图 4(a)、(b)和图 5(a)、(b),其中图 4(b)VP垂直摆记录的M5.0地震出现阶变同震响应。为便于比较VP垂直摆和JCZ地震计记录地震波信息的频段,将时频图中的归一化频率转化为真实频率值。
![]() |
图 4 内蒙古阿拉善左旗M5.0地震的震扰时频图 (a)JCZ原始曲线;(b)VP原始曲线;(c)JCZ S变换图;(d)VP S变换图 |
由时频分析(图 4(c)、(d)和图 5(c)、(d))可看出:①VP垂直摆和JCZ地震计在地震波到达初期,高频信号瞬间增多,且频带较宽,之后随着大幅值地震波到达,能量强度明显增强。随着地震波的衰减,影响频段逐渐变窄,能量强度也逐渐减弱。VP垂直摆和JCZ地震计反映的地震波频率、强度随时间变化过程类似,二者具有较好的一致性。②对于M5.0内蒙古区内近震,JCZ地震计记录地震波的最大影响频段为0~4.5Hz、高能量强度集中在0.5~2.5Hz,VP垂直摆记录的最大频段为0~0.5Hz、高能量强度集中在0~0.3Hz。③对于M5.9内蒙古区外远震,JCZ地震计记录地震波的最大影响频段为0.08~0.64Hz、高能量值集中在0.08~0.32Hz,VP垂直摆记录M5.9影响频段为0.03~0.22Hz,高能量值集中在0.06~0.16Hz。④无论近震还是远震,VP垂直摆记录地震波的频段均小于JCZ地震计记录的频段,且VP垂直摆主要记录的是低频地震波信号。
![]() |
图 5 云南普洱墨江县M5.9地震的震扰时频图 (a)JCZ原始曲线;(b)VP原始曲线;(c)JCZ S变换图;(d)VP S变换图; |
本文基于S波变换,对内蒙古中部地区定点形变观测中受降雨、大风、气压干扰以及地震波影响的典型干扰数据进行时频分析,获得如下结果和认识:
(1) 对呼和浩特台水管倾斜仪、洞体应变仪受降雨干扰的时频特征分析认为,降雨干扰信号主要集中在低频区域,呈三角状分布,其频率先增大后减小,与时域观测曲线受降雨影响时速率变化过程一致。对于整时值水管倾斜仪数据,干扰信号的周期可能与降雨量的多少、持续时间的长短等引起观测数据变化趋势、速率不同有关,其周期既可能大于日波,也可能会小于日波;呼和浩特台洞体应变仪分钟值数据受降雨影响时,干扰信号的归一化频率在0~0.025频段,即干扰信号的周期大于40min。
(2) 对包头台、乌加河台水管倾斜仪分钟值数据受大风干扰的时频特征分析认为,大风干扰的时频响应特征表现为高频干扰,最大归一化频率约为0.4,即最小周期约为2.5min,时频域能量强度与观测曲线受风扰影响变化幅度成正比,且干扰信号在时频域的频率、能量强度随时间变化过程与时域观测曲线受风扰影响过程一致。
(3) 对呼和浩特台、乌加河台洞体应变仪分钟值数据受气压干扰的时频特征分析认为,气压干扰信号的归一化优势频率分布在0~0.05频段内,即优势周期大于20min,干扰信号的最大归一化频率可达到0.2左右,即最小周期为5min,且干扰信号在时频域变化形态与其时域变化形态具有一致性。
(4) 在时频域中,受降雨、大风、气压影响的频段存在交叉部分,但不同干扰信号在时频域中能量强度、最高频率存在一定差别,因此由时频特征来判断具体干扰因素时,可综合参考这两点进行区分。并且,时频图中降雨、大风、气压干扰的色标范围不一致,总体上,降雨较大,气压次之,大风较小,表明降雨干扰信号优势频率的振幅与正常背景信号频率振幅的比值较大,因此相比于正常背景信号,降雨干扰信号的能量可能要大于气压、大风干扰。
(5) 对比分析乌加河台VP垂直摆和JCZ地震计记录地震波的时频特征,二者的频率、能量强度随时间的变化过程均可直观反映出地震波的衰减过程,具有较好的一致性。但VP垂直摆记录的地震波频段要小于JCZ地震计,且VP垂直摆主要记录的是低频地震波。
由于时频特征分析所选取的干扰样本量有限,对形变观测中不同干扰因素的时频响应特征总结不完全。随着样本量的增加以及对各类干扰因素影响机理的深入认识,今后将进一步研究总结,以期得出更加全面的内蒙古中部地区形变观测中不同干扰因素的时频响应特征。
陈学华, 2006, 时频分布与地震信号谱分析研究, 硕士学位论文, 成都: 成都理工大学.
|
戴勇、高立新、陈立峰等, 2016, 地震前兆数据时频分析, 地震地磁观测与研究, 37(6): 95-101. |
戴勇、高立新、杨彦明等, 2013, 基于小波变换方法的包头台形变分析, 华南地震, 33(4): 39-46. |
狄樑、陆德明、丁建国等, 2017, 常熟台倾斜观测资料的异常特征分析, 震灾防御技术, 12(2): 415-422. |
方燕勋、卞根发, 2019, S变换在定点形变观测中的应用, 华南地震, 39(1): 25-30. |
李继业、武晓军、刘长生等, 2015, 前兆观测干扰信号频谱特征分析, 中国地震, 31(3): 574-583. |
刘建明、李志海、孙甲宁等, 2016, 基于小波分析提取地倾斜异常特征, 地震, 36(1): 38-48. |
刘学谦、薛秀秀、杨龙翔等, 2015, 基于S变换的河南地区形变及重力观测数据时频特征研究, 地震, 35(4): 127-135. |
吕品姬、赵兵、陈志遥等, 2011, 小波分解-STFT方法在地形变观测数据中的应用, 大地测量与地球动力学, 31(5): 136-140. |
孙伶俐、李明、蒋玲霞等, 2013, 湖北省潮汐形变观测异常及干扰识别, 大地测量与地球动力学, 33(增刊Ⅰ): 36-40. |
万永革, 2012, 数字信号处理的MATLAB实现, 2版, 北京: 科学出版社.
|
汪翠枝、张磊、刘双庆等, 2010, 定点形变观测的降雨干扰及排除方法研究, 华北地震科学, 28(1): 42-47. |
王宁、吴云、张燕, 2014, 时频分析方法在形变数据中的应用研究, 地震工程学报, 36(2): 413-420. |
王宁、吴云、张燕等, 2015, 形变数据中天然扰动在时频域的响应研究, 内陆地震, 29(1): 65-70. |
杨澄雨, 2018, 时频分析方法在地震资料处理中的研究与应用, 硕士学位论文, 成都: 成都理工大学.
|
姚家骏、杨立明、冯建刚, 2011, 常用时频分析方法在数字地震波特征量分析中的应用, 西北地震学报, 33(2): 105-110. |
张雁滨、蒋俊、陈绍绪等, 2001, 连续形变的前兆参量、判别方法及实用化研究, 内陆地震, 15(1): 1-10. |
张雁滨、蒋骏、钱家栋等, 2002, 地壳介质微形变异常与强震短临前兆, 地震学报, 24(1): 103-108. |
张燕、吴云、刘永启等, 2004, 小波分析在地壳形变资料处理中的应用, 地震学报, 26(增刊): 103-109. |
赵莹, 2018, VP垂直摆观测背景的功率谱密度特征分析, 大地测量与地球动力学, 38(10): 1080-1085, 1095. |
周挚、山秀明、梁红等, 2005, 固体潮时频分析新方法, 地震研究, 28(4): 334-339. |
Stockwell R G, Mansinha L, Lowe R P, 1966, Localization of the complex spectrum:the S transform, IEEE Trans Signal Precess, 44(4): 998-1001. |