基于HHT方法分析西安地震台VP垂直摆观测背景特征＊

王建昌，颜文华，朱宣明，严珊

基于HHT方法分析西安地震台VP垂直摆观测背景特征＊

王建昌，颜文华，朱宣明，严珊

（陕西省地震局，陕西 西安 710068）

以西安地震台VP垂直摆倾斜仪2018年观测资料为研究对象，用HHT（Hilbert-Huang Transform）方法对VP垂直摆正常时段和干扰时段秒采样观测数据进行处理，获得HHT时频图和Hilbert边际谱图，总结VP垂直摆观测背景的频谱特征，有助于认识形变资料的正常观测背景和干扰特征，为识别地震异常信息提供客观的背景依据。

VP垂直摆；HHT方法；边际谱；背景特征

1 引言

布设在地表或近地表的定点形变观测仪器能监测地壳内部活动构造的微小运动，用于获得前兆异常信息，然而高灵敏度的数字化地形变观测仪器除了记录与地震前兆异常相关的信息外还包含各种干扰信息，主要存在如下几类干 扰[1]：①同震效应，数字化形变仪器可清晰记录到一定震级以上的地方震、近震及远震的地震波，即同震形变波。同震形变波能反映观测仪器的精密度，但对于地震前兆异常信息而言则为干扰。②自然环境干扰，如风扰、气温、气压、降雨等引起的干扰，通常会造成观测曲线加粗或短时畸变等。③观测系统、场地环境以及人为因素造成观测资料的突跳、台阶、畸变等。为了避免对前兆异常信息提取的干扰，准确识别地形变仪器中的干扰信息和正常观测背景，这是从观测数据中提取出可靠前兆异常信息的前提，也是进行地震预测必不可少的一项工作。

目前有众多学者对地形变观测资料的干扰信息和正常观测背景特征进行了分析，包括从原始观测资料曲线形态特征进行分析识别各种干扰因素和利用数字信号处理方法[2-3]，如冯红武等[4]将HHT方法引入到地电场数据处理中，周莉娟等[5]通过HHT方法对连续重力台站背景噪声进行了分析，李杰等[6]采用小波变换方法分析了形变资料的正常观测背景，赵莹[7]利用功率谱密度对VP垂直摆观测背景进行分析，进一步认识了测数据的正常背景和干扰特征。

本文以西安地震台VP宽频带垂直摆2018年观测资料为研究对象，通过HHT方法对VP垂直摆正常时段和干扰时段秒采样观测数据进行分析，依据西安地震台VP垂直摆观测背景中各种干扰因素与正常观测背景时段时频特征的差异，总结不同类型时段观测资料的频谱特征以进行区分，为识别地震异常信息提供客观的背景依据。

2 观测资料与处理方法

2.1 观测资料

西安地震台（简称西安台，经纬度：108.92°E，34.02°N）位于汾渭地震带秦岭山前大断裂南缘，海拔约630 m，VP垂直摆倾斜仪位于西安台形变山洞，形变山洞呈“十”字形，岩性为花岗岩，岩石风化得较为破碎，植被稀少，山体洞长总计72 m，NS向和EW向洞深分别为42 m和30 m，山洞年温差大于1°C。洞口与仪器室之间有四道隔门，保证稳定的观测环境，VP宽频带垂直摆倾斜仪用泡沫板密封，其墩基为混凝土墩。

2.2 HHT方法

美国国家宇航局HUANG等[8]最先提出Hilbert-Huang Transform（HHT）的概念，HHT方法对信号处理具有自适应性，其由经验模态分解（Empirical Mode Decomposition，简称EMD）和希尔伯特谱分析（Hilbert Spectrum Analysis）两部分组成，核心是经验模态分解。EMD分解可以将信号数据分离为有限个本征模态函数（Intrinsic Mode Function，简称IMF），然后利用希尔伯特变换构造解析信号，得到数据的瞬时频率和振幅，进而得到希尔伯特谱及边际谱。

2.2.1 经验模态分解

首先找出信号的所有极大值点和极小值点，用三次样条函数插值，分别得到信号的上、下包络线，对每一时刻的上包络和下包络求均值即得到一条均值线1（），用原始信号（）减去1（）得到时间序列1（），即：

然后检查1（）是否符合IMF条件，如果符合，则1（）为一个IMF分量；否则将1（）当作原始信号重复以上的操作步骤进行不同频率信号的筛选。重复筛选次后，直到k（）符合IMF条件，即得到个IMF分量和一个余量信号k（），公式表示如下：

式（1）（2）（3）中：i（）为第个IMF分量；k（）为余量信号。

2.2.2 希尔伯特谱分析

希尔伯特谱分析就是对IMF进行希尔伯特变换，得到时间序列的瞬时频率及振幅，也就是Hilbert谱，表示为：

式（4）中：i（）为瞬时振幅；j（）为瞬时频率；为取实部，式（4）中省略了余量信号k（）。

Hilbert谱（，）准确地表示了时间序列的振幅随频率和时间的变化规律。将（，）对时间进行积分得到Hilbert边际谱（），公式为：

Hilbert边际谱表示在某时间段内振幅随频率的变化情况，从边际谱上可以得到所研究的信号的优势频率范围。

3 观测资料背景特征

本文以西安台秒采样VP宽频带垂直摆倾斜仪EW分量观测资料为研究对象，根据奈奎斯特采样定理，VP垂直摆倾斜仪能可靠观测到周期大于2 s的信号，即计算结果频带范围为0～0.5 Hz。

3.1 正常背景特征

选取西安台2018年任一时段VP垂直摆EW分量秒采样无干扰的正常波形资料，对选取的2 880个数据点进行EMD分解得到5个IMF及一个余项，如图1所示，经过希尔伯特谱分析后得到希尔伯特二维时频图，如图2所示，再分别对5个IMF进行Hilbert边际谱分析。

从图2的时频图中可以看出，能量在频率域内的分布比较零散，整体上看主要集中于0.05～0.3 Hz频段内。

如图3所示，IMF1至IMF5分量边际谱表征的优势频率范围从高到低展布，其中IMF1分量的边际谱能量分布较发散，优势频率范围为0.05～0.35 Hz，IMF1分量的边际谱能量主要集中于0.05～0.15 Hz，IMF3分量的边际谱能量主要集中于0.05～0.1 Hz，IMF4和IMF5分量的优势频带宽度分别为0.003～0.04 Hz、0.003～0.02 Hz。整体从原始数据边际谱图上可以看出，在正常背景下，VP垂直摆EW分量优势能量的频率范围为0.05～0.15 Hz、0.05～0.1 Hz及低频范围0.003～0.02 Hz。其中VP垂直摆记录到的高频范围0.05～0.1 Hz为第1类地脉动信号（10～20 s），而0.05～0.15 Hz除了包含一部分第2类（5～10 s）地脉动信号外，还包含了其他成分的信号。

图1 西安台VP垂直摆EW分量正常时段原始曲线、各个IMF和余项

图2 西安台VP垂直摆EW分量正常时段Hilbert时频图

图3 西安台VP垂直摆EW分量正常时段边际谱示意图

3.2 干扰因素特征分析

高灵敏度的数字化地形变观测仪器除了记录与地震前兆异常相关的信息外还包含各种干扰信息，常见的干扰一般可分为自然环境干扰、人为干扰（检修仪器、调零等）、同震效应等，以上干扰因素会使仪器记录的观测曲线变化形态主要呈现为加粗或固体潮畸变（自然环境干扰）、突跳或阶变（人为干扰、同震响应）。本文根据对观测曲线变化形态的不同而分别对各种干扰因素进行分析，因同震形变波与人为干扰造成的突跳或台阶在振幅和持续时间上有差异而单独分析。

3.2.1 曲线加粗

2018-08-31T23：00—09-01T00：15，因刮风使得气压短时变化较大而造成西安台VP垂直摆倾斜仪原始观测曲线表现为加粗形态，呈现为高频干扰，与气压变化同步发生。图4为风扰时段原始曲线、各个IMF和余项示意图，风扰对观测曲线造成的干扰在各个IMF分量上都有表现。通过对此时段的观测曲线进行EMD分解并进行边际谱分析，如图5所示。对比图3正常时段边际谱图，可以看出在风扰时段内，在0.025～0.1 Hz频段，能量显著升高。风扰主要是大风以负载形式作用于安装仪器的整个山体的结果，风力作用于山体，导致山体承受一定的荷载而发生微弱倾斜[2]。而赵莹[7]利用功率谱密度方法对常熟台VP摆EW分量风扰时段观测资料分析显示，在0.2～0.3 Hz频段内，风扰时段内背景噪声能量显著升高。两种方法得到的风扰时段能量优势频率的差异可能与VP垂直摆仪器响应差异以及倾斜仪器所处山体山洞的深度、山体的岩性等有关。

图5 西安台VP垂直摆EW分量风扰时段边际谱示意图

3.2.2 突跳、阶变

2018-05-01因VP垂直摆倾斜仪EW分量电压值靠近量程范围而进洞对此调零，干扰时间较短。如图6所示，调零期间造成VP垂直摆倾斜仪产生突跳和小幅阶变，调零结束后数据恢复正常。通过HHT方法对调零时段的原始观测数据处理，从得到的图6时频图及结合图7边际谱图可以看到，调零期间的数据能量位于0.05 Hz频段以下，主要是长周期成份的干扰，与正常背景下的频谱特征有显著差异。

图6 调零时段原始曲线及Hilbert时频图

图7 调零时边际谱示意图

数字化形变仪器可清晰记录到一定震级以上的地震波，秒采样的地形变观测仪器观测曲线上基本可以识别出类似于地震波的d-P和d-S波[9]。远震同震波在垂直摆倾斜仪观测曲线上主要表现为高频干扰（突跳），近震同震波多表现为阶跃，也有会表现为固体潮畸变等同震效应[10]。本文选取的地震为2018-09-12T19：06发生在陕西省宁强县（32.75°N，105.69°E）s5.3地震，其震源深度11 km，震中距332 km。

西安台VP垂直摆EW分量记录到的地震原始观测曲线及经过EMD分解后的各IMF分量和余项，如图8所示，从图8可以看出从高频率信号的IMF1至低频率信号的IMF5和余项都含有地震同震波信号，经过希尔伯特时频分析得到的图9时频图及结合图10边际谱图，可以看到宁强县s5.3地震同震形变波的频率范围为0.05 Hz以下，宁强地震同震波能量主要集中于低频范围内。

图9 陕西省宁强县Ms5.3地震Hilbert时频图

图10 陕西省宁强县Ms5.3地震边际谱示意图

4 结论

本文通过HHT方法对西安地震台秒采样VP宽频带垂直摆正常时段和干扰时段观测数据进行处理，从非平稳的观测数据中提取Hilbert时频图和边际谱图，分析不同类型时段观测资料的频谱特征，得到以下结果：①在无任何干扰的正常情况下，西安台VP垂直摆倾斜仪能清晰记录到第1类和第2类地脉动信号，正常背景下主要是低频的自然噪声，集中于0.003～0.02 Hz频带内；②VP垂直摆在风扰时段表现为曲线加粗，会导致频率范围0.025～0.1 Hz能量显著升高；③导致原始观测曲线突跳、小幅度阶跃的进洞调零干扰与陕西省宁强县s5.3地震同震波干扰的能量优势频带都集中于0.05 Hz以下。

［1］宋治平，武安绪，王梅，等.小波分析方法在形变数字化资料处理中的应用［J］.大地测量与地球动力学，2003，23（4）：21-27.

［2］王秋宁，古云鹤.西安地震台形变观测自然环境干扰分析［J］.地震地磁观测与研究，2016，37（3）：72-77.

［3］李惠玲，高云峰，程冬焱，等.VP宽频带垂直摆倾斜仪观测干扰识别［J］.地震地磁观测与研究，2018，39（2）：100-107.

［4］冯红武，颜文华，王建昌，等.地电场观测中城市轨道交通干扰剔除［J］.地震地磁观测与研究，2019，40（3）：70-78.

［5］周莉娟，李辉，杨光亮，等.基于HHT的连续重力台站背景噪声分析方法研究［J］.大地测量与地球动力学，2012，32（5）：113-123.

［6］李杰，刘希强，李红，等.利用小波变换方法分析形变观测资料的正常背景变化特征［J］.地震学报，2005，27（1）：33-41.

［7］赵莹.VP垂直摆观测背景的功率谱密度特征分析［J］.大地测量与地球动力学，2018，38（10）：1080-1085.

［8］HUANG N E，ZHENG S，LONG S R，et al. The empirical mode decomposition and Hilbert spectrum for nonliner and non-stationary time series analysis［J］.Proceedings of the royal society a mathematical physical and engineering sciences，1998（454）：903-995.

［9］章鑫，叶青.闽粤地区同震形变波记录序列分析及小波能谱特征［J］.大地测量与地球动力学，2017，37（12）：1302-1307.

［10］胡宝慧，张浩，常金龙，等.鹤岗地震台垂直摆倾斜仪记录同震响应［J］.地震地磁观测与研究，2017，38（4）：137-141.

［11］葛洪魁，陈海潮，欧阳飚，等.流动地震观测背景噪声的台基响应［J］.地球物理学报，2013，56（3）：857-868.

［12］MORLET J，ARENS G，FOURGEAU E，et al. Wave propagation and sampling theory——Part II：Sampling theory and complex waves［J］.Geophysics，1982，47（2）：222-236.

2095－6835（2021）06－0030－05

P315

A

10.15913/j.cnki.kjycx.2021.06.010

陕西省地震局启航与创新基金（编号：QC201815、QC201907）

王建昌（1970—），男，工程师，主要从事地震监测工作。

〔编辑：丁琳〕