A Wavelet Comprehensive Threshold Seismic Signal Denoising Method
-
摘要: 小波阈值方法中硬、软阈值方法是地震信号降噪常用方法,但容易造成信号中高频信息丢失导致地震误判和漏判情况发生。小波综合阈值方法继承和发展了硬、软阈值降噪方法的优点,对信号高频部分用硬阈值方法,以提高高频信号能量,对信号低频部分用软阈值方法,提高信号降噪能力的同时保证信号连续性和光滑性。利用噪声信号小波系数小和地震信号小波系数大的特征,进行雷克子波降噪仿真实验和实际地震信号降噪实验。仿真实验表明,小波综合阈值方法降噪后波形MSE值最小,且降噪后与原信号波形最近似,降噪后波形高频部分能量增强且抑制低频部分能量。最后,对实际采集的地震信号进行降噪处理,处理后信号中能量增强被压制,利用处理后的信号可得到地震的初至时间。Abstract: The wavelet threshold method with hard or soft threshold is common for noise reduction of seismic signal, but it is often to cause the loss of high frequency information in the signal to lead in earthquake misjudgment or missing judgment. Taking the advantages of hard threshold and soft threshold denoising method, the wavelet comprehensive threshold improves the high frequency signal energy by the hard threshold treatment, and keeps the signal smoothness and reduce the signal degradation of noise capacity while ensuring signal continuity by lowness of the signal. Based on the small wavelet coefficients of the noise signal and the large wavelet coefficients of the seismic signal, the simulation experiment of the Ricker signal and the noise reduction experiment is carried out. The simulation results show that the MSE of the wavelet comprehensive threshold is the closest to the original signal waveform after noise reduction, and the energy of the high frequency part of the waveform is reduced and the low frequency part is suppressed. Finally, taking the noise reduction of the natural waveform as an example, we found that noise signal energy is suppressed and the first time of the waveform can be obtained.
-
Key words:
- Comprehensive threshold /
- Seismic signal /
- Wavelet method
-
引言
地震预警是利用电磁波传播速度远大于地震动传播速度以及破坏性地震波(S波或面波)的传播速度小于破坏性较小的P波速度的原理,在地震发生后尽可能短的时间内,确定地震基本参数或直接估计预警目标区的地震动和破坏程度,在破坏性地震波到来前,对可能破坏区发出地震警报。地震预警系统可以快速地进行地震定位并根据震中距不同,提供几秒、十几秒、甚至几十秒的时间,将预警信息在破坏性地震波到达之前发送至用户,以帮助用户采取应急措施,达到减灾的目的(马强,2008;殷海涛等,2012;张红才等,2013)。目前,日本、墨西哥、罗马和中国台湾地区已经建立了地震预警系统并提供预警服务,美国、瑞士、意大利和中国大陆则已建立了试验性质的地震预警系统(Allen等,2009;赵兵,2011)。
天津是中国北方最大的港口城市和环渤海经济中心,是国家防震减灾的重点设防城市,具有发生强烈地震的背景和条件,历史上曾多次遭遇特别重大的地震灾害。1679年9月2日三河-平谷8级地震和1976年7月28日唐山7.8级地震都使天津市遭受了严重的破坏。其中,唐山7.8级地震使宁河和汉沽的地震烈度达到Ⅸ度,塘沽及市区大部分地区地震烈度达到Ⅷ度,从而使天津成为中国唯一遭受Ⅷ度地震灾害的特大型城市。因此,建设天津简易烈度计地震预警试验区,能够建立完善有效的地震预警机制,进一步提高天津地区的地震灾害防御能力,尽可能地避免未来强震造成的人员伤亡和经济损失。
1. 试验区建设概况
天津简易烈度计地震预警试验区的建设目标是在天津地区新建80个简易烈度计观测点,与天津行政区内具备实时传输能力的测震台站和强震动台站共同组建天津地震预警观测网络,其网络规模如图 1所示。观测网络东西孔径约127km,南北孔径约160km,平均台站间距约10km。结合天津地震观测台网现有布局特征,选择观测条件较好的中小学校内平房或二层楼房的一层作为简易烈度仪观测场址。仪器固定在观测室不高于30cm的承重墙上,采用锚固方式安装,如图 2所示。数据传输采用移动4G无线VPDN方式,VPDN为虚拟专用拨号网络,是建立在VPN基础上且隶属与于VPN的1种业务。根据天津市远程网络传输需求,烈度计端路由设备首先通过4G无线网络接入附近移动运营商的无线网,再通过移动运营商VPDN专用网接入地震行业网,实现节点端与地震系统内网的互联互通。移动4G无线VPDN组网流程如图 3所示。
2. 简易烈度计性能
综合考虑仪器的安装方式、性能指标以及京津冀地震烈度速报与预警示范区的整体规划,天津简易烈度计地震预警试验区的台站观测仪器选用GL-P2B型和MI3000型三分向型简易烈度计,均具备16位以上分辨率,水平向测量范围-2—2g,垂直向测量范围-1—3g,动态范围大于80dB。传感器采用微机械加速度传感器(MEMS),MEMS加速度传感器为微电子技术与微机械工程结合发展的1种惯性器件,与传统力加速度传感器相比,具有体积小、成本低等特点。将MEMS传感器应用于地震仪器中,可以较低的成本提高台网密度,并可实现烈度速报与地震预警功能(Chung等,2011;Wu,2013;D'Alessandro等,2014)。2种仪器的性能指标见表 1,仪器的幅频特性曲线见图 4。
表 1 简易烈度计性能指标Table 1. Performance index of simple intensity meter型号 数量 采集延迟 分辨率 动态范围 加速度平坦频率 线性误差 GL-P2B 40 <1s 28位 >80dB 0.1—20Hz <1% MI3000 40 <1s 16位 >80dB 0.1—20Hz <1% 3. 地震预警数据处理系统
首先,简易烈度计台站、强震台站和测震台站的数据实时汇集到地震数据汇集系统,再通过行业内网将数据实时传输至地震预警数据处理系统,数据传输流程如图 5所示。当地震发生时,地震预警数据处理系统利用多个台站收到的地震初至波信号,使用STA/LTA算法检查震相触发,使用AIC算法精确拾取震相到时,使用“着未着”算法快速分析多个台站震相到时数据,确定地震发生的时刻和位置,使用P波最大振幅和台站震中距估算地震震级,从而获得地震基本参数,并向地震波尚未到达的地区发布准确的地震预警信号。
2018年8月16日13时50分20.97秒,天津市蓟州区发生2.2级地震,距离震中较近的简易烈度计台站、强震实时台站和测震台站均监测到了P波信号。地震预警数据处理系统采用256字节的数据打包模式进行数据传输,地震发生后8s发出第一报,发震时刻为2018年8月16日13时50分20.7秒,震中地点为天津蓟州区,震级为2.0级(转换震级),震中地点距天津市中心约94km,横波预计到达时间为12s,触发台站39个(图 6)。由于此次地震的震级较小,尚需经历更大的地震对天津简易烈度计地震预警试验区的效能进行检验。
4. 简易烈度计数据记录
简易烈度计观测的数据存储在地震数据处理中心服务器上,采用定长512字节MiniSeed格式,对于任一时段的数据,可利用JOPENS6.0数据处理系统转换为SEED、ASCII等多种数据文件。在数据处理过程中,采用python语言编写相关的计算程序,主要使用了python中的Numpy、Matplotlib、Basemap、Obspy等扩展包,有效提高了系统的开发效率和数据处理的工作效率。
在2018年8月16日天津蓟州2.2级地震中,位于天津蓟州附近的简易烈度计台站A0055、A0056、A0058、A0063均获取了观测数据。根据《仪器地震烈度计算(征求意见稿)》1的要求,对0.1—10Hz频率范围内的观测数据进行带通滤波处理,然后三分向合成计算最大峰值加速度。计算得到距离震中4km的A0058台站的最大峰值加速度为7.4gal;距离震中10km的A0055台站的最大峰值加速度为6.46gal;距离震中10km的A0063台站的最大峰值加速度为6.78gal;距离震中15km的A0056台站的最大峰值加速度为6.27gal。
1 中国地震局,2015.仪器地震烈度计算暂行规程.
由于简易烈度计台站没有同址的测震或强震台站,故选择台站距离相邻、方位相差较小的强震台站BAJ与简易烈度计台站A0056进行波形数据对比。BAJ台距离震中12km,三分向合成的最大峰值加速度为6.89gal,使用的地震计为SLJ-100型加速度计,测量的频带范围为0.1—30Hz,采样率为100Hz。简易烈度计台站A0056使用GL-P2B型MEMES加速度计,测量频带范围为0.1—30Hz,采样率为100Hz。观测台站分布如图 7所示。
由于观测频带不同,采用4阶巴特沃斯(Buttorworth)带通滤波器对强震台BAJ与简易烈度计台A0056的0.1—20Hz观测数据进行滤波,结果如图 8所示。强震台BAJ三分向记录的波形峰值加速度(PGA)分别为3.99gal、3.88gal和1.43gal,简易烈度计台A0056三分向记录的波形峰值加速度(PGA)分别为1.23gal、1.28gal和0.78gal。2个台站的峰值加速度和波形数据时程均不相同,这是由于台站位置不同造成的,其中强震台BAJ距离震中12km,简易烈度计台A0056距离震中15km,二者直线距离相差3km,且相对震中的方位有差异,但二者NS向与EW向波形存在一定的相似性。
观测仪器的记录噪声水平在一定程度上能够直接反应仪器的记录能力,尤其是记录中小地震事件的能力。为此,截取了简易烈度计台站A0056无明显干扰的10分钟的数据记录,利用Welch算法计算加速度噪声功率谱(Peterson,1993),结果如图 9所示。由图可以看出,简易烈度计台A0056的加速度噪声功率谱处于高噪声模型(New High Noise Model)和低噪声模型(New Low Noise Model)之上,且远高于高噪声模型。由此可知,相对于测震台站对地脉动噪声的完整记录(在NHNM和NLNM之间),简易烈度计台站记录的噪声数据完全是仪器的自噪声,天然噪声则淹没在仪器的自噪声中无法分辨。由于简易烈度计的仪器自噪声已经超出地脉动水平,无法观测记录到有效地脉动,且仪器通常布设在人员比较密集的建筑物内,因此,简易烈度计台站的记录数据信噪比较低,地震数据初至震相会淹没在背景噪声中。
各类地震观测仪器的观测特点不同,测震台站主要服务于地球物理学的研究,旨在了解地球内部的构造和地震活动性,为此提供各种地震参数,包括地震基本参数、地震活动性、地震波传播、地震趋势等,地震观测的场地条件要求为基岩、深井且背景噪声应尽量小。强震观测主要为地震工程学研究和结构抗震提供资料,旨在测定地面和建筑物在地震作用下的运动过程,包括地震影响场、地面及结构的地震动特征,主要针对的是近场强震,仪器频带范围为-2—2g,通常要求背景噪声小于1cm/s2,城区或乡村均可布设观测仪器。而简易烈度计是以MEMS传感器为主的非专业强震仪,主要针对的是近场强震,要求背景噪声小于2cm/s2,记录质量较低,一般布设在自由地表或2层以下低矮结构建筑物内1。简易烈度计具有成本低、安装要求不高、维护简单等特点,适于高密度布设,以提高区域地震预警和烈度速报能力。
1 马强,2017.地震预警与烈度速报系统的试验示范.浙江省地震局内部报告.
5. 讨论与结论
本文对天津简易烈度计地震预警试验区的建设、台站分布情况、仪器性能指标、数据处理及数据记录情况进行了介绍,探讨了简易烈度计在预警台网中的应用,获得了以下初步认识:
(1)通过建设天津简易烈度计地震预警试验区,加密了天津强震台网的密度,初步实现了天津及邻区的地震预警功能。
(2)由测震台站、强震动实时传输台站及简易烈度计台站组成的天津地震预警观测网络,拓宽了各台网的数据应用范围,不仅可应用于地震应急快速响应,还可应用于地震工程和科学研究。
(3)简易烈度计台站造价低、建设周期短,易于大面积布设或规模化应用,可用较少的投资实现更好的应用(姚会琴等,2012),是现有测震台站和强震台站的有益补充。
(4)简易烈度计记录的噪声数据完全是仪器的自噪声,高于高噪声模型,但对于震中附近满足一定信噪比条件的台站,可获得较准确的初至震相到时信息,对于地震预警信息的准确产出具有意义。
(5)今后还需进一步积累地震数据,结合天津及邻区的地震分布特点优化地震预警参数模型,以期提升该地区的防震减灾能力以及公共服务水平。
-
表 1 噪声分类
Table 1. Noise classification
噪声类型 频率范围/Hz 降噪难度 环境噪声 3—20 难 次生噪声 5—30 难 仪器噪声 1—2 易 表 2 仿真降噪后所得结果的SNR和MSE值
Table 2. SNR and MSE from the de-noised signal in simulation
降噪方法 SNR MSE 硬阈值 18.5736 0.2907 软阈值 19.3285 0.2897 小波综合阈值 19.4136 0.2784 -
耿冠世, 俞言祥, 2015.中国西部地区震源破裂尺度与震级的经验关系.震灾防御技术, 10(1):68-76. doi: 10.11899/zzfy20150107 侯跃伟, 赵兵, 田韬, 2015.基于Daubechies小波分析的南京数字化钻孔形变震前变化特征研究.震灾防御技术, 10(2):388-396. doi: 10.11899/zzfy20150219 孔祥茜, 吴继伟, 岳继光, 2005.地震信号小波变换的去噪方法.计算机辅助工程, 14(3):52-56. http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=jsjfzgc200503012 范涛, 赵兆, 吴海, 2014.矿井瞬变电磁多匝回线电感影响消除及曲线偏移研究.煤炭学报, 39(05):932-940. http://www.oalib.com/paper/4232079 刘霞, 潘洪屏, 高晓春, 2010.基于小波阈值的地震信号去噪处理.科学技术与工程, 10(29):7251-7254. doi: 10.3969/j.issn.1671-1815.2010.29.030 李英, 张淑贞, 许康生, 2006.小波降噪方法在地震信号处理中的应用.西北地震学报, 28(2):159-162. http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=xbdzxb200602013 曲中党, 吴蔚, 贺日政等, 2015.基于S变换的软阈值滤波在深地震反射数据处理中的应用.地球物理学报, 58(9):3157-3168. doi: 10.6038/cjg20150912 唐守锋, 童敏明, 潘玉祥等, 2011.煤岩破裂微震信号的小波特征能谱系数分析法.仪器仪表学报, 32(7):1521-1527. http://www.cnki.com.cn/Article/CJFDTOTAL-YQXB201107012.htm 魏学强, 袁洪克, 秦晶晶等, 2016.广义S变换地震信号时频分析.震灾防御技术, 11(4):808-813. doi: 10.11899/zzfy20160411 曾宪伟, 赵卫明, 师海阔等, 2010.利用小波包变换对地震信号进行时频分析时小波基函数的选取.地震研究, 33(4):323-328. http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=dzyj201004012 郑作亚, 卢秀山, 李克行, 2007.一类小波基函数的构造及其在测量数据处理中的应用探讨.测绘科学, 32(3):9-11. http://industry.wanfangdata.com.cn/dl/Detail/Periodical?id=Periodical_chkx200703003 Bruni V., Vitulano D., 2006. Wavelet-based signal de-noising via simple singularities approximation. Signal Processing, 86(4):859-876. doi: 10.1016/j.sigpro.2005.06.017 Mousavi S. M., Langston C. A., 2016. Hybrid seismic denoising using higher-order statistics and improved wavelet block thresholding. Bulletin of the Seismological Society of America, 106(4):1380-1393. doi: 10.1785/0120150345 Xia L. M., Zheng W., Liang M. D., 2017. Application of a comprehensive one dimensional wavelet threshold denoising algorithm in blasting signal. In: Balas V. E., Jain L. C., Zhao X. M., Information Technology and Intelligent Transportation Systems. Cham: Springer, 517-527. -