希尔伯特_黄变换谱及其在地震信号分析中的应用

第34卷第2期福州大学学报(自然科学版)Vol.34No.2 2006年4月Journal of Fuzhou University(Natural Science)Apr.2006

文章编号:1000-2243(2006)02-0260-05希尔伯特-黄变换谱及其在地震信号分析中的应用

陈子雄,吴琛,周瑞忠

(福州大学土木建筑工程学院,福建福州350002)

摘要:介绍了希尔伯特-黄变换(HHT)这一非线性、非平稳信号处理方法,并利用HHT处理了地震工程中

常用的El Centro地震波,得到了该信号的Hilbert谱、边际谱和能量谱,提取了该信号的主要动力特性,并与

该信号的Fourier分析结果进行了对比,显示出HHT这一方法的优越性.

关键词:希尔伯特-黄变换;经验模态分解;固有模态函数;地震信号

中图分类号:TU311.3文献标识码:A

Hilbert-Huang transform spectru m and its application in seismic signal analysis

CHEN Zi-xiong,W U Chen,ZHOU Rui-zhong

(College of Civil Engineering and Architecture,Fuzhou University,Fuzhou,Fujian350002,China) Abstract:HHT is a ne w method to deal with non-linear and non-stationary data.El Centro earth-

quake wave is analyzed by HHT.Through the way,Hilbert spectrum,marginal spectrum and energy

spec trum are got and dynamic property is extrac ted.The comparison between HHT spectrum and Fourier

spec trum is made and the superiority of HHT is demonstrated.

Keyw ords:Hilbert-Huang transform;empirical mode decomposition;intrinsic mode function;seismic signal

地震信号具有短时、突变等特点,是一种典型的非平稳随机信号,必须对其进行分析与处理,才可以提取信号的主要特征.传统的Fourie r变换能够表述信号的频率特性,但不提供任何时域信息[1],而小波分析虽然在时域和频域都具有很好的局部化性质,但本质上仍是一种窗口可调的Fourier变换,在小波窗内的信号必须是平稳的,因而没有根本摆脱Fourier分析的局限[2].小波基的选择也是信号分析中的一个重要问题,另外,小波基的有限长会造成信号能量的泄漏,使信号的能量-频率-时间分布很难定量表述.

Hilbert-Huang变换(HH T)的信号处理方法被认为是近年来对以Fourier变换为基础对线性和稳态谱分析的一个重大突破[2].它由经验模态分解(E mpirical Mode Decomposition,E MD)方法和Hilbert变换(H T)两部分组成,其核心是E MD分解.该方法采用了固有模态函数(Intrinsic Mode Function,I MF)概念以及将任意信号分解为I MF组成的思想,即E MD法,使得瞬时频率具有实际的物理意义[3].它不受Fourier分析的局限,可依据数据本身的时间尺度特征进行模态分解,分解过程中保留了数据本身的特性,再对各I MF分量进行Hilbert变换,得到信号能量在时间尺度上的分布规律,实现地震动力特性的提取.

1Hilbert-Huang变换

1.1经验模态分解和固有模态函数

经验模态分解(EMD)的目的是通过对非线性非平稳信号的分解获得一系列表征信号特征时间尺度的固有模态函数(I MF),使得各个I MF是窄带信号,可以进行Hilbert分析.首先设定两个条件:¹整个时间序列的极大极小值数目与过零点数目相等或最多相差一个;º时间序列的任意点上,由极大值确

收稿日期:2005-07-27

作者简介:陈子雄(1981-),男,硕士研究生;通讯联系人:周瑞忠,教授.

基金项目:教育部博士点专项科研基金资助项目(20040386004)

定的包络与由极小值确定的包络的均值始终为零.能满足以上两个条件的信号称为I MF 信号,用Hilbert 变换求I MF 信号的瞬时参数,其结果是准确的.

求取I MF 的具体步骤如下:

1)找出原始序列x (t)的各个局部极大值.这里,为更好地保留原序列的特性,局部极大值定义为时间序列中的某个时刻的值,其前一时刻的值不比它大,后一时刻的值也不比它大.然后用三阶样条函数进行插值,得到原序列的上包络序列值x max (t).同理,可以得到下包络序列值x min (t ).

2)对每个时刻的x ma x (t)和x min (t)取平均,得到瞬时平均值m(t):

m (t)=[x max (t)+x min (t)]/2

(1) 3)用原序列x (t)减去瞬时平均值m(t),得到一个去掉低频的新序列h(t ):

h(t)=x (t)-m(t)

(2)

此时,如果h(t)满足I MF 的两个条件,那它就是固有模态函数.否则,把h (t)当作原序列,重复以上步骤,直至满足I MF 的定义,得到第一个固有模态函数,并以c 1(t)记之.一般来说,c 1(t)代表了原始序列中的高频部分.然后,用原序列减去c 1(t),得到剩余值序列r 1(t):

r 1(t)=x (t)-c 1(t)

(3)

至此,提取第1个固有模态函数的过程全部完成.然后,把r 1(t)作为一个新的原序列,按照以上步骤,依次提取第2,第3,,,,直至第n 个固有模态函数c n (t).之后,r n (t)变成一个单调序列,再也没有内在模函数能被提取出来.

如果把分解后的各分量合并起来,就得到原序列x (t):

x (t )=

E n

j =1

c j (t)+

r n (t )(4)

1.2 Hilbert 变换和Hilbert 谱

通过E MD 分解得到的I MF 分量非常适合作Hilbert 变换,进而求出瞬时频率,得到HHT 谱.Hilbert 变换是一种线性变换,如果输入信号是平稳的,那么输出信号也应该是平稳的;Hilbert 变换强调局部属性,这避免了Fourier 变换时为拟合原序列而产生的许多多余的、事实上并不存在的高、低频成分.

对于任一固有模态函数c(t),其Hilbert 变换^c (t)定义为[4]:

^c (t )=

1P P V Q

c(S )

t -S

d S (5)

式中:P V 代表柯西主值,则对于c(t)的解析信号z (t)为:

z (t)=c(t)+i^c (t )=a(t )e

i H (t )

(6)式中:a(t)和H (t)分别为信号x (t)的瞬时振幅和瞬时相位,按下式计算:

a (t)=

c 2(t)+^c 2(t)

(7)H (t)=arctan (^c (t)/c(t))

(8)

由瞬时相位可得到信号的瞬时频率:

X (t)=d H (t)/d t (9)

可见,由Hilbert 变换得到的振幅和频率都是时间的函数,如果把振幅显示在频率-时间平面上,就可以得到Hilbert 幅值谱H (X ,t),简称Hilbert 谱,记作:

H (X ,t)=Re E n

j=1a j (t)e i Q X j

(t)d t

(10)

H (X ,t)精确地描述了信号的幅值随时间和频率的变化规律.将H (X ,t)对时间积分,就得到Hilbert

边际谱h (X ):

h(X )=

Q T

H (X ,

t)d t (11)

边际谱提供了对每个频率的振幅量测,表达了在整个时间长度内振幅的累积.将振幅的平方对频率积

#

261#第2期陈子雄,等:希尔伯特-黄变换谱及其在地震信号分析中的应用