孫 苗, 楊鈞凱, 吳 立

(1. 湖北國(guó)土資源職業(yè)學(xué)院 環(huán)境與工程學(xué)院,武漢 430090; 2. 武漢華中科大建筑規(guī)劃設(shè)計(jì)研究院有限公司,武漢 430070; 3. 中國(guó)地質(zhì)大學(xué)(武漢) 巖土鉆掘與防護(hù)教育部工程研究中心,武漢 430074; 4. 中國(guó)地質(zhì)大學(xué)(武漢) 工程學(xué)院, 武漢 430074)

自然界中絕大部分信號(hào)都是非線性、非穩(wěn)態(tài)信號(hào),振動(dòng)信號(hào)屬于此類(lèi),是典型的復(fù)雜多分量信號(hào)[1-3],希爾伯特-黃變換(Hilbert-Huang Transform, HHT)[4]具有的自適應(yīng)性使其很適合振動(dòng)信號(hào)時(shí)頻分析[5-6].通過(guò)時(shí)頻分析,可實(shí)現(xiàn)對(duì)以時(shí)域形式呈現(xiàn)的信號(hào)進(jìn)行“翻譯”,使其展現(xiàn)更多的隱藏信息.以此幫助相關(guān)工作人員研究振動(dòng)信號(hào)產(chǎn)生的時(shí)頻能量,是否會(huì)引起受迫振動(dòng)結(jié)構(gòu)產(chǎn)生不良反應(yīng)[7-9].

HHT由經(jīng)驗(yàn)?zāi)B(tài)分解(Empirical Mode Decomposition, EMD)和Hilbert變換組成.振動(dòng)信號(hào)監(jiān)測(cè)過(guò)程中不可避免會(huì)因測(cè)試儀器、監(jiān)測(cè)環(huán)境等因素而導(dǎo)致實(shí)際測(cè)得的信號(hào)中混有噪聲.噪聲的混入使得EMD得到的固有模態(tài)函數(shù)(Intrinsic Mode Function, IMF)具有嚴(yán)重的模態(tài)混淆[10-11],而傳統(tǒng)Hilbert變換受Bedrosian定理的約束在處理此類(lèi)IMF分量會(huì)得到負(fù)值瞬時(shí)頻率,導(dǎo)致實(shí)測(cè)振動(dòng)信號(hào)的時(shí)頻特征難以識(shí)別,或者識(shí)別結(jié)果缺乏準(zhǔn)確性[12-13].

基于此,提出CEEMDAN·MPE-INHT時(shí)頻分析算法,該算法實(shí)現(xiàn)分兩步.第1步抑制噪聲干擾,對(duì)EMD進(jìn)行雙重改進(jìn),為控制低頻噪聲改進(jìn)EMD得到自適應(yīng)補(bǔ)充集合經(jīng)驗(yàn)?zāi)B(tài)分解(Complete Ensemble Empirical Mode Decomposition with Adaptive Noise,CEEMDAN),為控制高頻噪聲引入多尺度排列熵(Multiscale Permutation Entropy, MPE),即通過(guò)CEEMDAN·MPE控制噪聲干擾.第2步排除Bedrosian定理的限制,優(yōu)化Hilbert變換,得到改進(jìn)歸一化Hilbert變換(Improved Normalized Hilbert Transform, INHT).

首先,通過(guò)混有噪聲的仿真振動(dòng)信號(hào)驗(yàn)證該算法不僅可有效抑制EMD固有的模態(tài)混淆現(xiàn)象,同時(shí)得到具有實(shí)際物理意義的瞬時(shí)頻率.再將該算法用于實(shí)際工程獲得的含噪振動(dòng)信號(hào)時(shí)頻分析中,發(fā)現(xiàn)該算法得到的時(shí)頻譜在時(shí)域和頻域均具有較高的分辨率,可有效提取振動(dòng)信號(hào)時(shí)頻能特征參數(shù),提高了時(shí)頻分析參數(shù)提取精度,對(duì)振動(dòng)信號(hào)危害控制具有一定的指導(dǎo)作用.

1 CEEMDAN·MPE-INHT時(shí)頻分析算法

1.1 CEEMDAN

CEEMDAN[14]是在分解的每一階段添加有限次自適應(yīng)白噪聲,實(shí)現(xiàn)在減少分解次數(shù)的同時(shí)減小重構(gòu)誤差.通過(guò)算法原理分析,發(fā)現(xiàn)CEEMDAN對(duì)監(jiān)測(cè)中混入低頻噪聲的控制源于所添加自適應(yīng)白噪聲的抵消能力.

1.2 MPE

MPE[15]相比排列熵(Permutation Entropy, PE)[16]的優(yōu)勢(shì)在于多尺度粗?；幚?這一處理是將時(shí)間序列(本文時(shí)間序列為CEEMDAN得到的IMF)進(jìn)行分段并在各段進(jìn)行平均化,進(jìn)一步提高處理精度.再通過(guò)設(shè)置排列熵閾值,控制處理后時(shí)間序列的熵值,進(jìn)而實(shí)現(xiàn)抑制振動(dòng)信號(hào)在監(jiān)測(cè)中混入的高頻噪聲.

1.3 改進(jìn)歸一化Hilbert變換

傳統(tǒng)的Hilbert變換受Bedrosian定理的約束,在處理非平穩(wěn)、多分量信號(hào)時(shí)會(huì)得到負(fù)值瞬時(shí)頻率.針對(duì)此現(xiàn)象,Huang等[17]提出歸一化希爾伯特-黃變換(Normalized Hilbert Transform, NHT),以此來(lái)提高瞬時(shí)頻率的計(jì)算精度.本文在NHT的基礎(chǔ)上做出了改進(jìn),Huang等[17]用于NHT的IMF來(lái)自于EMD或集合經(jīng)驗(yàn)?zāi)B(tài)分解(Ensemble Empirical Mode Decomposition, EEMD)[18-19],但EMD和EEMD在進(jìn)行含噪信號(hào)分解時(shí)不可避免會(huì)產(chǎn)生模態(tài)混淆的問(wèn)題.而本文IMF由CEEMDAN·MPE得到,對(duì)IMF模態(tài)混淆進(jìn)行了抑制.對(duì)此類(lèi)IMF進(jìn)行NHT,便可實(shí)現(xiàn)改進(jìn)歸一化Hilbert變換(INHT).INHT的操作和NHT十分相似,區(qū)別在于,NHT所用的IMF來(lái)源是EMD或EEMD,而INHT使用的IMF來(lái)自于CEEMDAN·MPE,因此這里不加累述,具體可參考文獻(xiàn)[17],將其中的IMF替換成CEEMDAN·MPE得到的IMF即可實(shí)現(xiàn)INHT.

1.4 CEEMDAN·MPE-INHT時(shí)頻分析算法流程圖解

CEEMDAN·MPE-INHT算法的建立思路是對(duì)含噪振動(dòng)信號(hào)時(shí)頻分析時(shí)遇到影響其精度的問(wèn)題進(jìn)行逐一改進(jìn),前述1.1節(jié)～1.3節(jié)簡(jiǎn)單說(shuō)明了組成CEEMDAN·MPE-INHT的每一部分單獨(dú)算法的運(yùn)行原理,并未詳細(xì)展開(kāi),為增加論文的可讀性,將CEEMDAN·MPE-INHT時(shí)頻分析算法的運(yùn)算流程圖進(jìn)行展示,如圖1所示.

圖1 CEEMDAN·MPE-INHT時(shí)頻分析算法流程圖Fig.1 Flow chart of CEEMDAN·MPE-INHT time-frequency analysis algorithm

2 仿真信號(hào)CEEMDAN·MPE-INHT時(shí)頻分析

2.1 仿真信號(hào)的建立

構(gòu)建含有噪聲的仿真振動(dòng)信號(hào),噪聲的添加是為了模擬在實(shí)際振動(dòng)信號(hào)監(jiān)測(cè)時(shí)混入的噪聲.仿真信號(hào)波形圖如圖2所示.仿真信號(hào)S(t)=x1(t)+x2(t),x1(t)=wgn(1,N, 0.1),即功率為0.1的噪聲信號(hào),N為采樣點(diǎn)數(shù),如圖2(a)所示;x2(t)=sin(2π×50t),即頻率f=50 Hz的正弦信號(hào),如圖2(b)所示;仿真信號(hào)如圖2(c)所示.采樣點(diǎn)數(shù)N=512,采樣時(shí)間t=1.0 s.

圖2 仿真信號(hào)波形圖Fig.2 Waveform of simulation signal

2.2 仿真信號(hào)模態(tài)分解

為突出CEEMDAN·MPE對(duì)EMD模態(tài)混淆的抑制作用,分別采用EMD、EEMD、CEEMDAN和CEEMDAN·MPE對(duì)S(t)進(jìn)行模態(tài)分解,分解結(jié)果如圖3所示.圖中:IMF為不同模態(tài)分解方法得到的固有模態(tài)函數(shù);R為不同模態(tài)分解方法得到的余項(xiàng).圖3(a)為EMD結(jié)果,不難發(fā)現(xiàn)IMF1是難以除去的噪聲信號(hào);IMF高頻模態(tài)混淆嚴(yán)重,如IMF2在t=0.24,0.42,0.54,0.80 s附近都出現(xiàn)高頻向中頻發(fā)展的現(xiàn)象;中頻模態(tài)混淆有所緩解,如IMF3在t=0.42 s及右端點(diǎn)附近存在向低頻發(fā)展的趨勢(shì);低頻相對(duì)穩(wěn)定,如IMF6.圖3(b)為EEMD結(jié)果,高頻IMF模態(tài)混淆相比EMD有所緩解,但依舊存在,如IMF1在t=0.42 s處向中頻發(fā)展;中頻相對(duì)穩(wěn)定,如IMF4和IMF5;低頻出現(xiàn)模態(tài)分裂[17],該現(xiàn)象是模態(tài)混淆的第2種現(xiàn)象,即相同的模態(tài)分量出現(xiàn)在不同的固有模態(tài)函數(shù)中,如IMF7和IMF8.圖3(c)為CEEMDAN結(jié)果,低頻分布穩(wěn)定,高頻IMF1出現(xiàn)了模態(tài)混淆,出現(xiàn)在t=0.42,0.64 s附近;圖3(d)為CEEMDAN·MPE結(jié)果,分解得到的IMF從高頻向低頻依次排列,未見(jiàn)明顯模態(tài)混淆現(xiàn)象.

圖3 仿真信號(hào)經(jīng)不同方法的分解結(jié)果Fig.3 Results of decomposition of simulation signals by different methods

2.3 仿真信號(hào)Hilbert變換和INHT對(duì)比分析

為驗(yàn)證CEEMDAN·MPE-INHT含噪振動(dòng)信號(hào)時(shí)頻分析的準(zhǔn)確性,依據(jù)圖3不同模態(tài)分解方法所得IMF結(jié)果,分別進(jìn)行EMD-HT、EEMDNHT、CEEMDAN-INHT和CEEMDAN·MPEINHT,得到的變換結(jié)果即時(shí)頻譜如圖4所示.圖中:E為能量譜密度.

圖4 仿真信號(hào)時(shí)頻譜圖Fig.4 Time-frequency spectrum of simulated signals

圖4(a)為EMD-HT時(shí)頻譜,可發(fā)現(xiàn)時(shí)頻譜出現(xiàn)了50 Hz以上的虛假分量,該分量較發(fā)散,難以識(shí)別,時(shí)頻譜在時(shí)間和頻率這兩個(gè)維度的分辨率都不高.圖4(b)為EEMD-HT時(shí)頻譜,可發(fā)現(xiàn)50 Hz以上的虛假分量相比EMD-HT少,頻率分辨率有所提高.圖4(c)為CEEMDAN-INHT得到的時(shí)頻譜圖,基本未見(jiàn)低頻模態(tài)混淆,50 Hz以下時(shí)頻譜頻率分辨率高,但是在t=0.12,0.36,0.80 s附近出現(xiàn)了高頻噪聲干擾導(dǎo)致頻率分辨率降低.對(duì)比圖4(a)、圖4(b)和圖4(c)時(shí)頻譜在高頻和低頻處的分辨率,可以發(fā)現(xiàn)CEEMDAN相比EMD和EEMD,對(duì)低頻噪聲具有很好的抑制作用,但對(duì)高頻噪聲缺乏控制力度.圖4(d)所示為CEEMDAN·MPE-INHT得到的時(shí)頻譜圖,該時(shí)頻譜在時(shí)間和頻率維度均具有較高分辨率,未見(jiàn)高頻和低頻模態(tài)混淆.再次表明CEEMDAN對(duì)噪聲信號(hào)引起的低頻模態(tài)混淆具有良好的抑制作用,MPE也能有效控制噪聲信號(hào)引起的高頻模態(tài)混淆,同時(shí)CEEMDAN·MPE得到的IMF經(jīng)過(guò)INHT能夠得到具有實(shí)際物理意義的時(shí)頻信息,即CEEMDAN·MPE-INHT時(shí)頻分析算法不僅可有效抑制噪聲信號(hào)引起的EMD模態(tài)混淆,同時(shí)得到時(shí)頻分辨率雙高的信號(hào)頻譜圖.

3 基于 CEEMDAN·MPE-INHT 的含噪振動(dòng)信號(hào)時(shí)頻分析

振動(dòng)信號(hào)在自然界分布極廣,本文實(shí)際應(yīng)用中,以含噪爆破地震波振動(dòng)信號(hào)作為研究對(duì)象,重點(diǎn)研究CEEMDAN·MPE-INHT算法的實(shí)際含噪振動(dòng)信號(hào)時(shí)頻分析.

為滿足交通運(yùn)輸需求,將福建某原有雙向四車(chē)道隧道原位擴(kuò)建為雙向八車(chē)道隧道,相關(guān)設(shè)計(jì)參數(shù)如表1所示.施工要求右側(cè)隧道在封閉施工過(guò)程中左側(cè)隧道依然保持正常通車(chē)狀態(tài),因此,在右側(cè)隧道進(jìn)行爆破施工時(shí),有必要對(duì)左側(cè)通車(chē)隧道進(jìn)行動(dòng)態(tài)跟蹤監(jiān)測(cè).

表1 隧道擴(kuò)建前后設(shè)計(jì)參數(shù)Tab.1 Design parameters before and after tunnel expansion

采用UBOX-5016爆破振動(dòng)智能監(jiān)測(cè)儀對(duì)左側(cè)通車(chē)隧道進(jìn)行監(jiān)測(cè),選取監(jiān)測(cè)點(diǎn)中一條典型爆破振動(dòng)信號(hào)作為研究對(duì)象,如圖5所示.圖中:v為速度;采樣時(shí)間t=1.0 s;采樣點(diǎn)為 1 500.

圖5 實(shí)際爆破地震波監(jiān)測(cè)信號(hào)波形圖Fig.5 The waveform of the actual blasting seismic wave monitoring signal

對(duì)圖5所示的爆破地震波監(jiān)測(cè)信號(hào)進(jìn)行基于CEEMDAN·MPE-INHT的時(shí)頻分析,分析結(jié)果如圖6～圖13所示.圖6～圖13中的(a)圖是圖5信號(hào)經(jīng)CEEMDAN·MPE得到的IMF;(b)圖是每個(gè)IMF對(duì)應(yīng)的邊際譜;(c)圖是每個(gè)IMF的時(shí)頻譜圖,這里忽略了余項(xiàng)R.

圖6 IMF1時(shí)頻能量特征圖Fig.6 Time-frequency energy characteristic of IMF1

圖7 IMF2時(shí)頻能量特征圖Fig.7 Time-frequency energy characteristic of IMF2

圖8 IMF3時(shí)頻能量特征圖Fig.8 Time-frequency energy characteristic of IMF3

圖9 IMF4時(shí)頻能量特征圖Fig.9 Time-frequency energy characteristic of IMF4

圖10 IMF5時(shí)頻能量特征圖Fig.10 Time-frequency energy characteristic of IMF5

圖11 IMF6時(shí)頻能量特征圖Fig.11 Time-frequency energy characteristic of IMF6

圖12 IMF7時(shí)頻能量特征圖Fig.12 Time-frequency energy characteristic of IMF7

圖13 IMF8時(shí)頻能量特征圖Fig.13 Time-frequency energy characteristic of IMF8

觀察圖6～圖13的(a)圖可發(fā)現(xiàn),圖5信號(hào)經(jīng)CEEMDAN·MPE得到了8個(gè)IMF分量,IMF從高頻到低頻依次排列,每個(gè)IMF分量都比較穩(wěn)定,未見(jiàn)明顯的模態(tài)混淆.

觀察圖6～圖13可發(fā)現(xiàn),基于CEEMDAN·MPE-INHT時(shí)頻分析算法得到的時(shí)頻譜圖在時(shí)域和頻域都擁有很高的分辨率.其中IMF1頻率范圍為100～300 Hz,持續(xù)時(shí)間0.00～1.00 s;IMF2頻率范圍為60～100 Hz,持續(xù)時(shí)間0.10～0.80 s;IMF3頻率范圍為40～60 Hz,持續(xù)時(shí)間0.00～0.80 s;IMF4頻率范圍為20～40 Hz,持續(xù)時(shí)間0.00～1.00 s;IMF5頻率范圍為10～20 Hz,持續(xù)時(shí)間0.08～0.80 s;IMF6頻率范圍為5～10 Hz,持續(xù)時(shí)間0.00～1.00 s;IMF7和IMF8頻率最低,主要在0～5 Hz,持續(xù)時(shí)間0.00～1.00 s.通過(guò)分析圖6～圖13可發(fā)現(xiàn),信號(hào)能量主要集中在IMF1和IMF2中,其次在IMF3和IMF4中,地震波信號(hào)頻率主要分布在100～300 Hz,能量也主要集中在高頻.

進(jìn)一步分析得到圖5中的信號(hào)的三維時(shí)間-頻率-能量譜圖,如圖14所示.觀察圖14,可發(fā)現(xiàn)圖5含噪振動(dòng)信號(hào)的主頻是110.5 Hz,爆破地震波信號(hào)能量也主要集中在100 Hz以上的頻率部分,該結(jié)果和單個(gè)IMF時(shí)頻能量特征圖時(shí)頻能量信息分布一致.現(xiàn)有研究表明[20-21],地下工程中既有隧道的固有頻率在50 Hz左右,結(jié)合本次爆破地震波信號(hào)時(shí)頻圖和三維時(shí)間-頻率-能量分布情況,可以初步判斷本工程分離式隧道一側(cè)隧道擴(kuò)挖爆破產(chǎn)生的地震波信號(hào)不會(huì)引起另一側(cè)隧道的共振.

圖14 爆破地震波監(jiān)測(cè)信號(hào)時(shí)間-頻率-能量三維圖Fig.14 Time-frequency-energy three-dimensional diagram of blasting seismic wave monitoring signal

綜上所述,基于CEEMDAN·MPE-INHT的含噪振動(dòng)信號(hào)時(shí)頻分析算法,不僅有助于抑制含噪監(jiān)測(cè)信號(hào)引起傳統(tǒng)EMD固有的模態(tài)混淆現(xiàn)象,同時(shí)CEEMDAN·MPE得到的IMF經(jīng)INHT得到的時(shí)頻圖和三維時(shí)頻能量圖可清晰展示振動(dòng)信號(hào)內(nèi)在所蘊(yùn)含的頻率-能量之間的對(duì)應(yīng)關(guān)系,即得到的時(shí)頻信息具有清晰的物理意義.對(duì)含噪振動(dòng)信號(hào)特征參數(shù)提取、振動(dòng)危害控制及科學(xué)制定防護(hù)措施具有重要的意義.

4 結(jié)論

(1) 通過(guò)CEEMDAN·MPE-INHT算法在含噪仿真振動(dòng)信號(hào)和實(shí)際含噪振動(dòng)信號(hào)時(shí)頻分析中的應(yīng)用,驗(yàn)證了該算法可抑制EMD固有的模態(tài)混淆現(xiàn)象,并得到時(shí)頻分辨率均高的時(shí)頻譜圖,譜圖對(duì)應(yīng)的時(shí)頻參數(shù)實(shí)際物理意義明確.

(2) INHT得到的瞬時(shí)頻率為IMF的調(diào)頻分量Hilbert變換得到的,使得傳統(tǒng)的Hilbert變換免受Bedrosian定理的限制,因此INHT相比Hilbert變換得到的瞬時(shí)頻率真實(shí)性更高.

(3) 基于CEEMDAN·MPE-INHT含噪振動(dòng)信號(hào)時(shí)頻分析算法,可得出含噪振動(dòng)信號(hào)內(nèi)在所蘊(yùn)含的時(shí)間-頻率-能量之間的對(duì)應(yīng)關(guān)系,該對(duì)應(yīng)關(guān)系有助于含噪振動(dòng)信號(hào)內(nèi)在特征的識(shí)別和實(shí)際工程振動(dòng)危害控制.