Hilbert-Huang 变换在爆破震动信号分析中的应用

 张义平, 李夕兵

(中南大学 资源与安全工程学院, 湖南 长沙, 410083)

摘要: 介绍了Hilbert-Huang变换(HHT)法的原理、 内容和优越性, 并用仿真信号进行实例分析, 以验证其关键技术经验模态分解(EMD)的高效性、 自适应性, 以及其时频图能定量地描述时间与瞬时频率的关系。 用HHT法对爆破震动信号进行分析与处理。 研究结果表明: EMD能很好地按不同的时间尺度对信号进行分解, 分解后的固有模态函数能反映信号本身所固有的特性; 能将Hilbert能量谱中的信号能量清晰地表示在时间-频率-能量的分布图上; HHT法能有效地提取爆破震动信号的时频特征; HHT法比小波分析更具适应性, 为爆破震动信号的分析与处理提供了新的研究思路与方向。

关键词: 爆破震动; 时频分析; Hibert-Huang变换; 经验模态分解; 瞬时频率; Hilbert谱

中图分类号:TB41; TD235 文献标识码:A 文章编号: 1672-7207(2005)05-0882-06

Application of Hilbert-Huang transform in blasting vibration signal analysis

ZHANG Yi-ping, LI Xi-bing

(School of Resources and Safety Engineering, Central South University, Changsha 410083, China)

Abstract: The principle, contents and superiorities of Hilbert-Huang transform (HHT) were briefly introduced and example analysis was carried out with simulative signal on computer, and blasting vibration signal was analyzed and processed as an example as well. The results show that the empirical mode decomposition method which is the key of HHT method with different time scales can efficiently decompose the blasting vibration signal into intrinsic mode functions reflecting the intrinsic characteristics of the signal, and the signal energy can be clearly shown on the Hilbert energy spectrum with an energy-frequency-time distribution. Besides, the HHT method, which is more adaptive than other methods, can draw time-frequency characteristics from blasting vibration signal and provides a new research approach to blasting vibration analyzing and processing.

Key words: blasting vibration; Hilbert-Huang transform; empirical mode decomposition; instantaneous frequency; Hilbert spectrum

   爆破作用引起的地震波是一种极为复杂的宽带非平稳随机波, 所产生的地震信号具有短时、 突变且衰减快等特点^[1,2]。 从爆破地震信号中提取峰值和主震频率等有限的特征信息的处理爆破地震信号手段已不能满足需要。 随着时频分析技术的发展, 越来越多的处理手段已运用到爆破信号的处理中, 其中小波分析已成为研究热点, 它具有诸多优点, 正逐渐取代传统的爆破震动信号分析方法^[3-6]。 但是小波本质上是一种窗口可调的傅里叶变换,其小波窗内的信号必须是平稳的, 因而没有根本摆脱Fourier分析的局限, 同时还存在选择小波基和频谱扩散的问题^[7-9]。

  Hilbert-Huang变换(HHT)是N. E. Huang等于1998年提出的一种新的信号处理方法, 是近年来对以傅里叶变换为基础的线性和稳态谱分析的一个重大突破^[10], 它主要由经验模态分解(Empirical Mode Decomposition, EMD)方法和Hilbert变换(Hilbert Transform)两部分组成。 固有模态函数(Intrinsic Mode Function, IMF)和EMD法比小波分析更具良好的适应性^[8,10]。 尽管这种新方法提出的时间不长, 其理论和方法有待研究和完善, 但已在不同领域得到了研究和应用^[11-15]。 本文引入HHT方法, 通过计算机仿真处理, 对现场爆破震动信号进行实例分析, 以期为爆破震动信号的分析与处理提供一种新的方法。

1 Hilbert-Huang变换

1.1 EMD方法

EMD方法, 即Huang变换, 是HHT方法的关键, 它依据信号本身的时间尺度(信号相邻峰值点之间的时差)特征来进行分解, 将信号分解为含有不同时间尺度且满足以下2个定义条件的IMF:

a. 整个数据序列中, 极值点的数量与过零点的数量相等或至多相差1;

b. 信号上任意一点,由局部极大值点确定的包络线和由局部极小值点确定的包络线的均值均为0,即信号关于时间轴局部对称。

该分解算法称为筛选过程(Sifting process), 其步骤如下: 首先找出原始信号x(t)上所有的极值点, 用三次样条函数曲线对所有的极大值点进行插值, 从而连接各极大值点拟合出x(t)的上包络线, 同理得到下包络线, 并确保2条包络线包含了所有的信号数据。 接着将两条包络线的均值定义为m₁, x(t)与m₁的差定义为h₁₁=x(t)-m₁。 若h₁₁满足IMF的定义条件, 则为第一个IMF, 否则继续对h₁₁进行重复筛选。 一般h₁₁不会满足IMF的条件, 假定经过k次筛选后的结果h_1k满足IMF的定义, 则x(t)第一个IMF分量记为c₁=h_1k。 然后令x(t)与c₁的差r(t)=x(t)-c₁(t)为新的信号数据, 重复以上筛选过程, 再依次得到IMF分量c₂, c₃, …, c_n, 当分量c_n或余量r_n小于预先设定之值, 或余量r_n成为单调函数时, 就结束筛选过程。 经上述分解步骤后, x(t)就可分解为n个IMF分量及r_n的和:

1.2 Hilbert变换

信号经分解后得到多个IMF的组合, 对每个IMF分量进行Hilbert变换, 即可得到每个IMF分量的瞬时频率, 综合所有IMF分量的瞬时频谱就可获得Hilbert谱。 先对信号x(t)的IMF分量c(t)作Hilbert变换:

其中, V为柯西主值。 据此构造解析信号z(t):

其中: a(t)为幅值函数, Φ(t)为相位函数。

再在此基础上求出的角瞬时频率ω(t)和瞬时频率f(t)分别为

将每一个IMF分量进行Hilbert变换之后, 则可把x(t)表示成:

式(8)中除去了余差r_n, Re表示取实部。 上式把信号幅度在三维空间中可表达成时间与瞬时频率的函数, 即Hilbert谱。 若振幅的平方对时间积分, 则可以得到Hilbert能量谱:

希尔伯特能量谱提供了每个频率的能量计算式,表达了每个频率在整个时间长度内所累积的能量。

1.3 HHT方法的优越性

HHT变换具有诸多优越性: 其关键方法EMD分解没有固定的先验基底,是自适应的, 比依赖于先验函数基的傅里叶及小波等分析方法更适合处理非平稳信号; 首次定义了IMF, 指出其幅值允许改变, 突破了传统的将幅值不变的简谐信号定义为基底的局限,使信号分析更加灵活多变; 每一个IMF可以看作是信号中一个固有的振动模态,通过Hilbert变换得到的瞬时频率具有清晰的物理意义,能够表达信号的局部特征; 它能精确地做出能量—时间—频率图, 这是小波等其他信号分析难以做到的^[10]。

2 计算机仿真实例

为验证HHT方法的有效性, 利用3个正弦信号与1个指数趋势向相加作为仿真信号(如图1所示)进行仿真, 该信号由以下方程表示:

S(t)=sin(200πt)+sin(100πt)+sin(40πt)+e^t。(10)

对该信号(采样率1 kHz, 数据长度1024点)进行EMD分解, 得到3个分量c₁, c₂, c₃和一个趋势量c₄(如图2所示)以及相应的各IMF分量的时频图(如图3所示)。 从图2可以看出, EMD方法是按不同的时间尺度分解信号, 先分解出高频, 再依次分解出低频, 次低频, 最后得到趋势项。 分解出的3个IMF分量和趋势向正是仿真信号的4个原始信号, 表明了分解的高效性。 EMD分解过程是信号本身所决定的, 是一个自适应分解过程, 能很快地提取信号特征并分解出信号的分量, 其过程类似小波分解但又不同。 小波变换中需要选择具体的小波基, 小波基不同其结果也不一样, 有时甚至相差很大, 这就决定其分解过程会出现较多的分量及导致频谱扩散。 时频图定量地描述了时间与瞬时频率的关系。 从图3可以看出, 瞬时频率分别很好地集中在100 Hz, 50 Hz, 20 Hz和0 Hz周围, 只是端点由于污染有点摆动, 但其结果与预期结果有很好的一致性, 表明HHT法是有效的。

图 1   仿真信号

Fig. 1   Simulative signal

图 2   信号的EMD分量

Fig. 2   EMD components of signal

图 3   IMF分量的时频图

Fig. 3   Time-frequency distribution of IMF components

3 应 用

在某开山填海工程中采用中深孔爆破, 由于周边环境复杂, 对一些特殊结构物进行了爆破震动监测, 以确保正常生产和人员设备的安全。 图4所示为其中的一个典型的多段爆破地震动波原始信号, 其采样间隔时间为0.4 ms。 现用HHT法对信号进行分析。

图 4   爆破震动信号

Fig. 4   Blasting vibration signal

EMD分解结果的分量及其频谱分别如图5~7所示, 原始信号重构和重构误差结果如图8所示。 从图5~8可以看出:

a. 原始信号已被分解为8个IMF分量c₁~c₈和1个余量, 并依次按时间尺度(相邻两显著波峰间的距离)从大至小的顺序分解出, 即先分解高频再分解低频。 c₁频率最高, 波长最短, 随着分解的进行, 所得IMF分量频率逐渐变低, 波长越来越长,直到分解出频率已经很低的最后一个余量。 各IMF分量包含了不同的时间特征尺度, 可以不同的分辨率显示信号特征, 说明EMD分解中分辨率是自适应的, 相对小波的多分辨率分析, EMD方法在分解信号时更加简单。

图 5   EMD分解结果

Fig. 5   Results of EMD decomposition

b. EMD分解中没有固定的基函数, 每个IMF分量的提取都是由信号本性所决定的, 具有自适应性, 使信号分析更加灵活多变。 而其他分析方法如小波和Fourier法需预先设定基函数, 基函数本身有一定的局限, 并且选择的不同其结果也不一样, 有时甚至相差较大, 另外也不可避免地引入了频率扩散。

c. EMD分解出的IMF分量大都有清晰的物理意义。 分量c₁频率最高, 所占能量非常小, 并均匀分布于整个信号中, 表明它是监测中引入的高频噪声, 需在进一步的分析中去噪。 c₂为震动信号中的高频部分, 其能量很小, 说明在地震波的传播过程中高频已大幅度衰减。 分量c₃, c₄, c₅和c₆频率较低, 振幅增大, 包含了信号的大部分能量, 为原始信号的优势频段, 对结构物的破坏主要由这些分量所造成, 应该加以重点考虑。 另外, c₇和c₈是分解后频率更小的分量, 可能是信号固有的, 也可能是其他情况引起的。 最后的余量表明监测仪器的漂零或信号微弱的变化趋势。 这表明EMD方法是一种新的主成分分析方法, 它按频率从高至低提取信号本身所固有的全部模态函数, 并且分解出的几个IMF分量集中了原信号中最显著的信息。

图 6   IMF分量c₁~c₄的频谱

Fig. 6   Spectrum of IMF component c₁-c₄

图 7   IMF分量c₅~c₈的频谱

Fig. 7   Spectrum of IMF component c₅-c₈

d. 信号的频谱很丰富, 分布在0~1250 Hz之间, 但其分布又是有规律的, 信号频谱大部分分布在150 Hz以下。 可见多段爆破震动信号的优势频率主要集中在10~50 Hz, 并可分成多个子震频率段。 这和图9所示原始信号的频谱基本一致。

图 8   重构信号(a)和重构误差(b)

Fig. 8   Signal after reconstruction (a) and error of reconstruction (b)

图 9   爆破震动信号的频谱

Fig. 9   Spectrum of blasting vibration signal

e. 通过各IMF分量重构后的信号和原始信号相比是非常接近的, 其误差很小。 从图6可看出, 重构后的信号和原始信号的误差量级在10^-16以上, 可完全满足工程计算和分析要求, 说明EMD分解的可靠性。

将EMD分解得到的IMF分量通过Hilbert变换,再利用公式(9)可以得出Hilbert能量谱, 结果如图10所示。 从Hilbert能量谱图中可以看出, 由EMD得到的各IMF分量很清晰以频率—时间—振幅图的形式表现在颜色编码图上(能量以对数形式表示), 其分辨率很高。 各IMF分量的瞬时频率是一个围绕中心频率上下小幅波动的量,相互间很少交叉重叠。 波动能量基本上都很清晰地分布在频率为100 Hz左右, 并主要集中在50Hz的低频区。 而高于100 Hz的频段其能量分布很少, 这和信号经EMD分解后所分析的结果是完全吻合的, 进一步表明HHT法在爆破震动信号分析中的有效性。

图 10   信号的Hilbert能量谱

Fig. 10   Hilbert energy spectrum of signal

4 结 论

a. 通过Matlab语言编程, 用仿真信号对其进行计算机算例分析, 显示了EMD分解的高效性、 自适应性, 时频图能定量地描述时间与瞬时频率的关系, 为Hilbert谱提供了变换的基础。

b. 利用HHT法对爆破震动信号进行分析与处理, 结果表明HHT法能更好地适应非平稳信号的特点, 能依据爆破震动信号本身固有的特征按不同的时间尺度快速地分解信号, 并将信号能量定量地表示在时间-频率-能量分布的Hilbert能量谱图上, 能有效地提取爆破震动信号的时频特征, 较小波分析方法更具适应性和优越性。

c. HHT 变换是分析爆破地震信号的有效的新方法, 为进一步认识地震波的传播机理、 破坏原因和危害的预测提供了新的研究思路与方向。

参考文献:

[1]张雪亮, 黄树棠. 爆破地震效应[M]. 北京: 地震出版社, 1981.

ZHANG Xue-liang, HUANG Shu-tang. Blasting Seismic Effect[M]. Beijing: Earthquake Press, 1981.

[2]汪旭光,于亚伦. 关于爆破震动安全判据的几个问题[J]. 工程爆破, 2001, 7(2): 88-92.

WANG Xu-guang, YU Ya-lun. On several problems of safety criterion for blasting vibration engineering blasting[J]. Engineering Blasting, 2001, 7(2): 88-92.

[3]张贤达, 保铮 . 非平稳信号分析与处理[M]. 北京: 国防工业出版社, 1998.

ZHANG Xan-da, BAO Zhen. Analysis and processing of nonstationary signal[M]. Beijing: National Defence Industry Press, 1998.

[4]林大超, 施惠基, 白春华, 等.爆炸地震效应的时频分析[J].爆炸与冲击, 2003, 23(1): 31-35.

LIN Da-chao, SHI Hui-ji, BAI Chun-hua, et al. Time-frequency analysis of explosion seismic effects[J]. Explosion and Shock Waves, 2003, 23(1): 31-35.

[5]何军, 于亚伦, 梁文基. 爆破地震信号的小波分析[J]. 岩土工程学报, 1998, 20(1): 47-50.

HE Jun, YU Ya-lun, LIANG Wen-ji. Wavelet analysis for blasting seismic signals[J]. Chinese Journal of Geotechnical Engineering, 1998, 20(1): 47-50.

[6]宋光明. 爆破振动小波包分析理论与应用研究[D]. 长沙: 中南大学资源环境与建筑工程学院, 2001.

SONG Guan-ming. Wavelet Packet Analysis of Blasting Vibration: Theory and Applications[D]. Changsha: School of Resources and Safety Engineering, Central South University, 2001.

[7]凌同华. 爆破震动效应及其灾害的主动控制[D]. 长沙: 中南大学资源与安全工程学院, 2004.

LING Tong-hua. Blast Vibration Effect and Initiative Control of Vibrational Damage[D]. Changsha: School of Resources and Safety Engineering, Central South University, 2004.

[8]崔锦泰. 小波分析导论[M]. 西安: 西安交通大学出版社, 1995.

CUI Jin-tai. Theoretical System of Wavelet Packet[M]. Xi′an: Xi′an Jiaotong University Press, 1995.

[9]石春香, 罗奇峰. 时程信号的Hilbert-Huang变换与小波分析[J]. 地震学报, 2003, 25(4): 398-405.

SHI Chun-xiang, LUO Qi-feng. Hilbert-huang transform and wavelet analysis of time history signal[J]. Journal of Earthquake, 2003, 25(4): 398-405.

[10]仲佑明, 秦树人, 汤宝平. 一种振动信号新变换法的研究[J]. 振动工程学报, 2002, 15(2): 233-238.

ZHONG You-ming, QIN Shu-ren, TANG Bao-ping. Study on a new transform method for vibration signal[J]. Journal of Vibration Engineering, 2002, 15(2): 233-238.

[11]Huang N E, Shen Z, long S R, et al. The empirical mode decomposition and the hilbert spectrum for nonlinear and nonstationary time series analysis[J] . Pro Roy Soc London A, 1998, 454: 903-995.

[12]ZHAO Jin-ping, HUANG Da-ji. Mirror extending and circular spline function for empirical mode decomposition method[J]. Journal of Zhejiang University(Science), 2001, 2(3): 247-252.

[13]ZHANG Rui-chong, MA Shuo, et al. Interpretation and application of hilbert-huang transformation for seismic performance analyses[J]. Technical Council on Lifeline Earthquake Engineering Monograph, 2003, 25: 657-666.

[14]丁康, 陈健林, 苏向荣, 平稳和非平稳振动信号的若干处理方法及发展[J]. 振动工程学报, 2003, 16(1): 1-10.

DING Kang, CHEN Jian-lin, SU Xiang-rong. Development in vibration signal analysis and processing methods[J]. Journal of Vibration Engineering, 2003, 16(1): 1-10.

[15]胡昌华, 张军波, 夏军. 基于MATLAB的系统分析与设计——小波分析[M]. 西安: 西安电子科技大学出版社, 1999.

HU Chang-hua, ZHANG Jun-bo, XIA Jun. System Analysis and Design Based on MATLAB6.X-Wavelet Analysis[M]. Xi′an: Xidian University Press, 1999.

收稿日期:2005-01-14

基金项目:国家自然科学基金资助项目(50490274, 50490272, 10472134)

作者简介:张义平(1970-), 男, 湖南邵东人, 博士研究生, 从事岩土灾害控制、 震动信号分析与处理等研究

论文联系人: 张义平, 男, 博士研究生; 电话: 0731-8836628(O); E-mail: zyp_stone@163.com