台站测震资料干扰分析处理软件研究

何思源，李贵元，刘华姣，钟李彬，徐文海

(四川省地震局成都地震基准台,四川 成都 611730)

地震观测是研究地震和地球内部结构的重要手段，也是地震学理论体系建立和发展的重要基础，数字地震记录已成为现代地震研究的基石(李兴泉等，2018)。目前，MSDP软件广泛应用于台站,该软件提供了波形分析功能，缺少波形干扰处理功能,无法满足实际工作需求。虽然通过滤波程序可以达到滤掉部分干扰的效果，但整体缺乏综合分析处理能力，难以推广应用。地震波形中包含了众多的信息，其中不乏各种干扰，这些干扰对数字地震记录的质量造成了影响，便会对于震相的识别与分析造成干扰，为测震分析工作增加难度和工作量，因此如何快速有效地完成对干扰的去除就显得十分必要。本研究基于MATLAB开发了一款测震资料干扰分析处理软件，其内部不仅集成了快速傅里叶变换、希尔伯特-黄变换、小波变换、速度功率谱估计这些分析算法，还包括了不同类型的去噪方法，可以满足装配速度型地震计台站观测资料的干扰分析与处理需求。

1 总体设计思路

如何去除地震波形中的干扰信息，保留其有效信息，首先需要了解有效信息和干扰的特点，才能正确地将二者加以区分。傅里叶变换、希尔伯特-黄变换、小波变换和速度功率谱估计则是整个“波形分析”功能的基础，整体流程如图1所示。“波形处理”采用的是FIR滤波器和小波去噪两种方法，FIR滤波器由于具有线性相位的特性而被广泛使用，小波变换去噪实质是软阈值去噪法。阈值法去噪的优点是噪声几乎完全得到有效抑制，且反映原始信号的特征尖峰点也得到很好的保留，能够使估计信号实现最大，均方误差最小化，即去噪后的估计信号是原始信号的近似最优估计。小波变换具有时频局部化特性、多分辨特性、去相关特性和选基灵活性的特点(周宇峰等，2008)。前两个性质决定了小波变换去噪相比传统方法具有独特的优势，是信号处理的得力工具。在利用小波变换分析和处理信号时，不同领域的学者对如何选取合适的小波基有不尽相同的看法，而且对同一信号选取不同的小波基进行处理时会出现完全不同的结果(张华等，2011)。因此，需要根据信号的特点，采取折衷的办法进行小波基的选取(高玉宝等，2008)。软件中给出了一些常用类型的小波基供用户进行选择，如图2所示。

图1 波形分析与处理总流程

图2 软件所提供的小波基类型

2 功能演示

2.1 干扰的分析

选取成都地震台JCZ-1T地震计记录的2016年7月27日19时59分四川平武2.5级地震垂直分向波形数据，见如图3a。选择“波形分析”功能中的“频谱分析”功能，可以分别得到FFT频谱图和希尔伯特边际谱，如图3b所示。“Welch速度功率谱”可以根据设置不同的参数(FFT长度、窗函数长度、重采样精度)得到相应的速度功率谱估计图，如图3c。根据FFT频谱图及希尔伯特边际谱可知，在2～4 Hz的频率区间内存在较高的幅值，而速度功率谱估计也在此区间出现了高幅值。

图3 波形分析功能示例

为了研究该2～4 Hz的频率区间所表示的内容，文中首先采用FIR滤波器对地震波形中的2～4 Hz频率范围进行带阻滤波，以此研究2～4 Hz频率范围是否与地震波形中的干扰有关。选择软件中的“带阻滤波器”功能，通过设置通、阻带频率，便可以实现带阻滤波，如图4所示。滤波前后的对比结果显示：滤除2～4 Hz频率范围的干扰后，干扰被压制的效果十分明显，原始波形中的“毛刺”干扰也得到了有效的改善，表明滤除的2～4 Hz频段波正是干扰所在的频率范围。

图4 带阻滤波前后波形

2.2 干扰的处理

选取成都地震台JCZ-1T地震计记录的2017年12月10日08点14分九寨沟2.9级地震记录作为处理实例。原始波形自身存在明显的波形“毛刺”，同时波形中还存在高频干扰信息，可以采用小波去噪和2～4 Hz带阻滤波法综合去除噪声的办法进行处理，其中2～4 Hz带阻滤波法是为了消除测震资料中的日常普遍性干扰(何思源等，2019)；再将处理的结果进行仿真，如图5所示。其中，小波去噪法选择的小波基为bior2.4小波，因为bior2.4小波基重构精度高，去噪效果良好，能够更加完整地保留了有效信号(何思源等，2019)，其分解层数为4。原始波形经过处理后， “毛刺”及高频干扰得到了有效的压制，仿真后波形的垂直分向上原本被干扰所掩盖的初至震相变得清晰可辨，信噪比得到了明显的提升。

图5 去噪前后波形及仿真结果

3 实例应用

由于经济发展、城市规模扩大以及交通建设等对部分台站的观测条件产生了影响。泸州地震台由于地处市区且长期受到周边城建施工的影响，台站观测数据整体呈“昼差夜良”的特点，如图6所示。泸州地震台的观测数据质量保持“昼夜一致”对提高该地震台测震数据质量有着积极的意义。分别选取一段2019年3月12日泸州地震台白天和夜晚记录的波形数据进行波形分析。通过“波形分析”功能中的“频谱分析”可以分别得到白天和夜晚观测数据的FFT频谱图，白天相比夜晚主要在高频段具有更高的幅值水平(见图7)，再采用小波去噪与滤波器相结合的方式进行处理。结果表明，白天的观测数据去噪处理后，波形质量得到明显改善，能够与夜晚的观测数据质量保持一致，如图8所示。

图6 泸州地震台波形数据

图7 泸州地震台白天及夜晚波形数据FFT频谱图

图8 泸州地震台去噪前后波形质量对比

4 疑似普遍性干扰提示功能

此项功能的作用是通过记录用户使用带阻滤波器滤除的频率区间，当某一频率区间累计出现了20次，那么软件会自动将其识别为“疑似普遍性干扰”，并及时提醒用户注意此频率区间信息的干扰。该功能的原理是用户在第一次使用“带阻滤波器”功能时，会在D盘自动建立一个文件夹且命名为“滤波器结果文件”，里面会自动生成一个名为“带阻滤波”的Excel表格并记录下本次所选择滤除的频率区间。之后用户每使用一次“带阻滤波器”功能，软件都会对选择处理的频率区间进行记录并统计出现的次数。当有一组频率区间第20次出现在Excel表格中时，软件就会弹出消息提醒框，提示用户注意该频率区间。例如，软件在进行带阻滤波处理时，用户有20次选择滤除2～4 Hz的频率区间，那么软件便会弹出提示信息，如图9所示。

图9 “疑似普遍性干扰频率”提示功能示例

5 认识与讨论

在测震资料干扰分析处理软件的开发和使用过程中，有以下几点认识：

(1)软件的核心是“波形分析”与“波形处理”功能。“波形分析”将傅里叶变换、希尔伯特-黄变换、小波变换和速度功率谱估计作为基础，对于测震资料中的干扰能够实现时域、频率的分析；“波形处理”功能以FIR滤波器及小波去噪为手段，能够去一般的干扰信号。(2)“疑似普遍干扰频率”提示功能作为附加功能，可以及时地为用户反馈测震资料中可能存在的普遍性干扰，提示用户对干扰产生的原因进行核实。(3)在实际应用中分析地震波形中存在的干扰，并能进行有效的处理，改善波形数据的质量。(4)经过分析发现，成都地震台测震资料中普遍存在频率区间为2～4 Hz的干扰，需要对该频段可疑干扰源进行核实。