孫 苗,吳 靜,吳 立,楊鈞凱,覃亞男

(1.湖北國(guó)土資源職業(yè)學(xué)院 環(huán)境與工程學(xué)院,武漢 430090;2.中國(guó)地質(zhì)大學(xué)(武漢) a.巖土鉆掘與防護(hù)教育部工程研究中心;b.工程學(xué)院,武漢 430074;3.湖北工程學(xué)院 土木工程學(xué)院,孝感 432000)

希爾伯特-黃變換(Hilbert-Huang Transform,HHT)由經(jīng)驗(yàn)?zāi)B(tài)分解(Ensemble Empirical Mode,EMD)[1]和希爾伯特變換(Hilbert Transform,HT)[2]組成,其中EMD以其能識(shí)別數(shù)據(jù)內(nèi)在屬性并根據(jù)數(shù)據(jù)本身特征進(jìn)行分解在信號(hào)分析領(lǐng)域得到了廣泛應(yīng)用[3]。對(duì)EMD分解結(jié)果進(jìn)行Hilbert變換,更是實(shí)現(xiàn)了將時(shí)域信號(hào)轉(zhuǎn)變?yōu)轭l域信號(hào),幫助爆破工程技術(shù)人員對(duì)爆破地震波信號(hào)傳播特征、內(nèi)在屬性進(jìn)行識(shí)別,對(duì)爆破地震波危害控制起到了一定的指導(dǎo)作用[4-6]。

但EMD對(duì)含有噪聲的爆破地震波信號(hào)相對(duì)敏感,而實(shí)際爆破地震波監(jiān)測(cè)又無(wú)法避免噪聲信號(hào)的混入[7],這便造成了實(shí)際爆破地震波信號(hào)EMD得到的固有模態(tài)函數(shù)(Intrinsic Mode Function,IMF)[8,9]存在嚴(yán)重模態(tài)混淆[10-12],而Hilbert變換處理這類IMF分量會(huì)產(chǎn)生錯(cuò)誤的時(shí)頻分析結(jié)果,這樣的結(jié)果對(duì)爆破地震波危害效應(yīng)分析意義不大,有時(shí)候甚至?xí)蔀樽璧K爆破工程技術(shù)人員作出判斷的干擾項(xiàng)。因此如何對(duì)HHT處理含噪爆破地震波信號(hào)遇到導(dǎo)致分析精度受損的問(wèn)題進(jìn)行處理,是目前爆破地震波危害效應(yīng)分析亟待解決的問(wèn)題。

鑒于此,對(duì)影響EMD-Hilbert時(shí)頻分析精度的因素進(jìn)行逐一改進(jìn)。對(duì)EMD進(jìn)行改進(jìn),使其可降低對(duì)噪聲的敏感性。具體實(shí)現(xiàn)通過(guò)兩步,首先在EMD中添加自適應(yīng)白噪聲,得到自適應(yīng)補(bǔ)充集合經(jīng)驗(yàn)?zāi)B(tài)分解(Complete Ensemble Empirical Mode Decomposition with Adaptive Noise,CEEMDAN)[13-15],用于抵抗監(jiān)測(cè)中混入噪聲的低頻部分;再者引入多尺度排列熵(Multiscale Permutation Entropy,MPE)[16-18]到CEEMDAN中,充分發(fā)揮MPE隨機(jī)性檢測(cè)能力,抑制高頻噪聲對(duì)分解精度的影響。以上兩步得到的IMF可認(rèn)為是降噪處理后的IMF,對(duì)CEEMDAN·MPE得到的IMF進(jìn)行歸一化Hilbert變換(Normalized Hilbert transform,NHT)[19,20]使得傳統(tǒng)的Hilbert變換不受Bedrosian定理的約束,降低負(fù)值瞬時(shí)頻率出現(xiàn)的概率。通過(guò)上述三步可實(shí)現(xiàn)HHT時(shí)頻分析精度提升。

最后進(jìn)行含噪仿真爆破振動(dòng)信號(hào)和實(shí)測(cè)含噪水下鉆孔爆破地震波信號(hào)CEEMDAN·MPE-NHT時(shí)頻分析算法和HHT算法對(duì)比研究,以驗(yàn)證CEEMDAN·MPE-NHT算法能有效提高HHT時(shí)頻精度,得到反映真實(shí)爆破振動(dòng)屬性的時(shí)頻特征參數(shù),對(duì)爆破振動(dòng)特征識(shí)別和制定科學(xué)的防護(hù)措施具有一定的指導(dǎo)作用。

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

CEEMDAN·MPE-NHT時(shí)頻分析算法可分兩步構(gòu)建,第一步建立CEEMDAN·MPE算法,其后對(duì)CEEMDAN·MPE算法得到的IMF進(jìn)行NHT。

CEEMDAN·MPE算法的核心是通過(guò)剔除噪聲來(lái)抑制EMD模態(tài)混淆。它將MPE的隨機(jī)性檢測(cè)功能與CEEMDAN的自適應(yīng)性相結(jié)合,可有效地減少噪聲的干擾。有關(guān)CEEMDAN和MPE的詳細(xì)介紹,請(qǐng)參閱文獻(xiàn)[13-15]和[16-18]。在這里,僅通過(guò)流程圖1分析CEEMDAN·MPE的操作過(guò)程。

圖1 CEEMDAN·MPE算法流程圖Fig. 1 CEEMDAN·MPE algorithm flow chart

Hilbert變換得到具有實(shí)際物理意義的瞬時(shí)頻率通常需滿足十分苛刻的要求,模態(tài)混淆的IMF分量往往不能滿足要求。鑒于此Huang等提出NHT[20],以此來(lái)提高瞬時(shí)頻率的計(jì)算精度。

NHT的原始操作是對(duì)IMF的調(diào)頻分量進(jìn)行Hilbert變換,IMF來(lái)源是EMD或集合經(jīng)驗(yàn)?zāi)B(tài)分解(Ensemble Empirical Mode Decomposition,EEMD)。不難發(fā)現(xiàn),本文得到的IMF物理意義較EMD和EEMD得到的IMF更加清晰,真實(shí)性更高。將文獻(xiàn)[20]中的IMF替換成CEEMDAN·MPE得到的IMF,便是本文的實(shí)現(xiàn)途徑。進(jìn)一步分析,得到CEEMDAN·MPE-NHT時(shí)頻分析算法運(yùn)算流程圖,見(jiàn)圖2。

圖2 CEEMDAN·MPE-NHT算法流程圖Fig. 2 CEEMDAN·MPE-NHT algorithm flow chart

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

為突出CEEMDAN·MPE-NHT時(shí)頻分析算法相比HHT可有效提高含噪爆破地震波信號(hào)時(shí)頻特征參數(shù)的提取精度,特進(jìn)行CEEMDAN·MPE-NHT和HHT的含噪仿真振動(dòng)信號(hào)時(shí)頻分析對(duì)比研究。

含噪仿真振動(dòng)信號(hào)S(t)=x1(t)+x2(t),其中x1(t)=wgn(1,N,0.1),即功率為0.1的噪聲信號(hào),如圖3所示;x2(t)=sin(2×pi×100×t),即頻率為100 Hz的正弦信號(hào),如圖4所示;仿真信號(hào)如圖5所示。采樣點(diǎn)數(shù)N=512,采樣時(shí)間t=1/N∶1/N∶1。

圖3 功率為0.1的噪聲信號(hào)圖Fig. 3 Noise signal with power of 0.1

圖4 頻率為100 Hz的正弦信號(hào)圖Fig. 4 Sine signal with frequency of 100 Hz

為突出CEEMDAN·MPE對(duì)EMD模態(tài)混淆的抑制作用,分別采用EMD和CEEMDAN·MPE對(duì)S(t)進(jìn)行模態(tài)分解。圖6為EMD分解結(jié)果,不難發(fā)現(xiàn)IMF1是難以除去的噪聲信號(hào);高頻模態(tài)混淆嚴(yán)重,如IMF2在0.53 s,0.64 s及0.85 s附近均出現(xiàn)了高頻向中頻發(fā)散的趨勢(shì);中頻模態(tài)混淆有所緩解,如IMF3在0.41 s及右端點(diǎn)附近存在向低頻發(fā)展的趨勢(shì);低頻相對(duì)穩(wěn)定,如IMF5和IMF6。圖7為CEEMDAN·MPE分解結(jié)果,分解得到的IMF從高頻向低頻依次排列,未見(jiàn)明顯模態(tài)混淆現(xiàn)象,模態(tài)顯示相比EMD結(jié)果清晰且穩(wěn)定。

圖7 仿真信號(hào)經(jīng)CEEMDAN·MPE得到的IMFFig. 7 IMFs of simulated signal by CEEMDAN·MPE

為突出NHT對(duì)瞬時(shí)頻率具備清晰物理意義的貢獻(xiàn),進(jìn)行Hilbert變換和NHT對(duì)比分析。對(duì)EMD得到的IMF進(jìn)行Hilbert變換,得到的時(shí)頻譜見(jiàn)圖8。對(duì)CEEMDAN·MPE得到的IMF進(jìn)行NHT,得到的時(shí)頻譜見(jiàn)圖9。

圖8 EMD-HT得到的仿真信號(hào)時(shí)頻譜圖Fig. 8 Time-frequency spectrum of simulated signal obtained by EMD-HT

圖9 CEEMDAN·MPE-NHT得到的仿真信號(hào)時(shí)頻譜圖Fig. 9 Time-frequency spectrum of simulated signal obtained by CEEMDAN-INHT

圖8為EMD-HT時(shí)頻譜,模態(tài)混淆分量經(jīng)Hilbert變換得到的時(shí)頻譜分辨率不高,出現(xiàn)了100 Hz以上難以識(shí)別的虛假分量。圖9為CEEMDAN·MPE-NHT時(shí)頻譜,該算法通過(guò)CEEMDAN·MPE來(lái)抑制S(t)中的x1(t);CEEMDAN·MPE得到的IMF經(jīng)NHT處理后得到的頻率分布在0～100 Hz,和x2(t)具有對(duì)應(yīng)性,時(shí)頻參數(shù)可清晰識(shí)別。側(cè)面說(shuō)明了CEEMDAN·MPE得到的IMF經(jīng)過(guò)NHT能夠得到具有實(shí)際物理意義的時(shí)頻信息,即CEEMDAN·MPE-NHT時(shí)頻分析算法不僅可有效抑制噪聲引起的EMD模態(tài)混淆,同時(shí)得到時(shí)頻分辨率雙高的信號(hào)頻譜圖。

3 炸礁爆破地震波信號(hào)CEEMDAN·MPE-NHT時(shí)頻分析

以三峽-葛洲壩兩壩間蓮沱河段航道整治炸礁工程為研究對(duì)象,選擇下岸溪～丁頭鎮(zhèn)區(qū)間LT5炸礁區(qū)為研究對(duì)象。該區(qū)總的炸礁工程量2561.0 m3,清碴工程量2484.2 m3。LT5區(qū)周邊環(huán)境復(fù)雜,炸礁區(qū)距離最近的民房70 m、距離公路(陡紙線)81 m,還有滑坡、崩塌和泥石流等不良地質(zhì)情況,對(duì)工程影響較大。

選取炸礁爆破施工期間處于較危險(xiǎn)地帶的民房作為研究對(duì)象,對(duì)民房進(jìn)行現(xiàn)場(chǎng)監(jiān)測(cè)得到一系列監(jiān)測(cè)數(shù)據(jù)。觀察監(jiān)測(cè)數(shù)據(jù)發(fā)現(xiàn)水下鉆孔爆破地震波信號(hào)三個(gè)方向振動(dòng)速度呈現(xiàn)出有規(guī)律的一致性,水平徑向峰值振動(dòng)速度最大,水平切向次之,而垂直方向振動(dòng)速度最小?？紤]文章篇幅,僅對(duì)水下鉆孔爆破地震波監(jiān)測(cè)信號(hào)的水平徑向峰值振動(dòng)速度進(jìn)行CEEMDAN·MPE-NHT時(shí)頻分析。有代表性水平徑向地震波監(jiān)測(cè)信號(hào)如圖10所示。

圖10 水平徑向地震波監(jiān)測(cè)信號(hào)波形圖Fig. 10 Waveform of monitored radial seismic signal

對(duì)圖10信號(hào)進(jìn)行CEEMDAN·MPE-NHT時(shí)頻分析。先對(duì)圖10信號(hào)進(jìn)行CEEMDAN·MPE,得到如圖11所示的IMF分量?？砂l(fā)現(xiàn),IMF從高頻到低頻排列,每個(gè)IMF都攜帶了爆破地震波信號(hào)的特定時(shí)頻信息。其中IMF3攜帶圖10所示信號(hào)的主要能量,其次是IMF4和IMF5。

圖11 水平徑向地震波信號(hào)CEEMDAN·MPE結(jié)果Fig. 11 The results of CEEMDAN·MPE of the radial seismic signal

進(jìn)一步分析,對(duì)圖11中的IMF1～I(xiàn)MF6執(zhí)行NHT,獲得每個(gè)IMF的時(shí)頻能量特征圖,其中圖12～圖17為對(duì)應(yīng)IMF的邊際譜;圖18～圖23為對(duì)應(yīng)IMF的時(shí)頻譜;這里未對(duì)IMF7和R進(jìn)行變換,是因?yàn)榻Y(jié)合圖11不難發(fā)現(xiàn)IMF7和R占有的能量可忽略不計(jì)。觀察圖12～圖23可發(fā)現(xiàn),CEEMDAN·MPE-NHT獲得的時(shí)頻譜在時(shí)域和頻域都具有高分辨率,這點(diǎn)和仿真信號(hào)時(shí)頻分析得到的結(jié)果一致。IMF1的頻率范圍為100～250 Hz,持續(xù)時(shí)間為0.24～0.38 s;IMF2的頻率范圍為50～100 Hz,持續(xù)時(shí)間為0.23～0.42 s;IMF3的頻率范圍為20～50 Hz,持續(xù)時(shí)間為0.22～0.47 s;IMF4的頻率范圍為10～30 Hz,持續(xù)時(shí)間為0.28～1.17 s;IMF5的頻率范圍為5～10 Hz,持續(xù)時(shí)間為0.00～1.20 s;IMF6頻率最低,約為0～5 Hz,持續(xù)時(shí)間為0.00～1.20 s。

圖12 IMF1邊際譜Fig. 12 Marginal spectrum of IMF1

圖13 IMF2邊際譜Fig. 13 Marginal spectrum of IMF2

圖14 IMF3邊際譜Fig. 14 Marginal spectrum of IMF3

圖15 IMF4邊際譜Fig. 15 Marginal spectrum of IMF4

圖16 IMF5邊際譜Fig. 16 Marginal spectrum of IMF5

圖17 IMF6邊際譜Fig. 17 Marginal spectrum of IMF6

圖18 IMF1時(shí)頻譜Fig. 18 Time-frequency spectrum of IMF1

圖19 IMF2時(shí)頻譜Fig. 19 Time-frequency spectrum of IMF2

圖20 IMF3時(shí)頻譜Fig. 20 Time-frequency spectrum of IMF3

圖21 IMF4時(shí)頻譜Fig. 21 Time-frequency spectrum of IMF4

圖22 IMF5時(shí)頻譜Fig. 22 Time-frequency spectrum of IMF5

圖23 IMF6時(shí)頻譜Fig. 23 Time-frequency spectrum of IMF6

通過(guò)圖12～圖23的分析,可發(fā)現(xiàn)信號(hào)在低頻停留時(shí)間大于高頻,主要能量集中在50 Hz以下。

進(jìn)一步分析,獲得圖10中信號(hào)的三維時(shí)頻能量譜,如圖24所示。從圖24可看出,本次水下炸礁工程監(jiān)測(cè)得到的水下鉆孔爆破地震波信號(hào)的主頻為26.83 Hz,次主頻為19.32 Hz。

圖24 信號(hào)時(shí)間-頻率-能量三維圖Fig. 24 Time-frequency-energy spectrum of the signal

根據(jù)結(jié)構(gòu)抗震知識(shí),當(dāng)爆破地震波的頻率與房屋的固有頻率相同時(shí),結(jié)構(gòu)的振幅將達(dá)到最大,即發(fā)生共振危險(xiǎn)。物體的受迫振動(dòng)是根據(jù)特定規(guī)律進(jìn)行的,即由于形狀和結(jié)構(gòu)的不同,它們具有不同的固有頻率。為了探索建筑結(jié)構(gòu)的固有頻率,有必要分析結(jié)構(gòu)的振動(dòng)特性,如振動(dòng)特征頻率與頻率對(duì)應(yīng)的振型,即模態(tài)分析。通過(guò)PKPM結(jié)構(gòu)軟件的SETWE模塊,對(duì)距離爆區(qū)最近的2層房屋進(jìn)行了模態(tài)分析。房屋的第一個(gè)6階振型圖如圖25所示,每個(gè)振型對(duì)應(yīng)的固有頻率如表1所示。

表1 民房前6階陣型對(duì)應(yīng)的自振頻率(單位:Hz)Table 1 Natural frequency corresponding to the first 6-order array of civil houses(unit:Hz)

圖25 民房前6階陣型圖Fig. 25 The first six formations of civilian houses

從表1中可看出,受保護(hù)房屋的第六階陣型對(duì)應(yīng)的自振頻率為27.08 Hz,本次爆破的主頻為26.83 Hz;被保護(hù)房屋的第五階陣型對(duì)應(yīng)的自振頻率為20.17 Hz,本次爆破的次主頻為19.32 Hz。根據(jù)結(jié)構(gòu)共振的知識(shí),不難發(fā)現(xiàn)此次爆破產(chǎn)生的地震波可能會(huì)引起該房屋的共振。

進(jìn)一步定量分析不同IMF分量對(duì)民房每一階陣型的影響程度,IMF對(duì)不同的建筑結(jié)構(gòu)具有不同的放大倍數(shù),根據(jù)式(1),可計(jì)算IMF分量對(duì)民房每一階陣型的放大倍數(shù),計(jì)算結(jié)果見(jiàn)表2,β為單個(gè)IMF主頻和房屋各階陣型固有頻率之比,λ為阻尼比,一般建筑都λ取0.05?？紤]到本次爆破地震波能量主要在IMF3中,其次在IMF4及IMF5中,所以僅對(duì)IMF3～I(xiàn)MF5對(duì)房屋不同陣型的放大性進(jìn)行討論。

表2 IMF分量對(duì)不同陣型下民房的放大系數(shù)DTable 2 The IMF component corresponds to the magnification factor D under different formations of civil houses

(1)

分析表2可發(fā)現(xiàn),房屋固有頻率和信號(hào)主頻越接近即單個(gè)IMF主頻和房屋各階陣型固有頻率之比β趨近于1時(shí),D越大。D值大小反應(yīng)在一定的外界爆破信號(hào)激勵(lì)下,結(jié)構(gòu)對(duì)應(yīng)質(zhì)點(diǎn)產(chǎn)生的質(zhì)點(diǎn)振動(dòng)速度放大效應(yīng),D值越大對(duì)應(yīng)的爆破作用放大效應(yīng)越強(qiáng),對(duì)結(jié)構(gòu)的損傷越大,結(jié)構(gòu)發(fā)生破壞的概率越大。

不難發(fā)現(xiàn)IMF3可使民房第6階陣型房屋質(zhì)點(diǎn)速度產(chǎn)生9.9倍的放大效應(yīng)。當(dāng)IMF主頻和陣型自振頻率相等時(shí),即β=1,此時(shí)D=10,結(jié)構(gòu)振動(dòng)幅度達(dá)到最大,這便解釋了為什么民房處測(cè)得的速度峰值僅0.49 cm/s,但在民房處依舊產(chǎn)生了大量的裂縫,IMF分量對(duì)民房每一階陣型的放大倍數(shù)是不一樣的,實(shí)際結(jié)構(gòu)在共振的作用下振幅呈現(xiàn)出放大的特征,可據(jù)此解釋房屋在監(jiān)測(cè)低振速環(huán)境下發(fā)生開(kāi)裂的原因。因此,實(shí)際工作中安全評(píng)估不能以單一速度峰值作為判別標(biāo)準(zhǔn)。

繼續(xù)觀察表2還可發(fā)現(xiàn),對(duì)于民房結(jié)構(gòu)而言由于結(jié)構(gòu)前六階陣型相對(duì)較小,而高頻信號(hào)對(duì)此種結(jié)構(gòu)的影響相對(duì)50 Hz以下的低頻影響相對(duì)小。針對(duì)此現(xiàn)象建議施工采用降低單段藥量法、微差起爆、優(yōu)化裝藥結(jié)構(gòu)等提高爆破地震波信號(hào)頻率的方法。

綜上分析可發(fā)現(xiàn),基于CEEMDAN·MPE-NHT的爆破地震波信號(hào)時(shí)頻分析方法,不僅有助于抑制含噪監(jiān)測(cè)信號(hào)引起的EMD模態(tài)混淆,同時(shí)CEEMDAN·MPE得到的IMF經(jīng)NHT處理后,得到的時(shí)頻譜在分辨率和精度上都得到了提升,分析結(jié)果有助于爆破振動(dòng)特征的識(shí)別和爆破振動(dòng)危害控制。

4 結(jié)論

(1)為提高HHT含噪爆破地震波信號(hào)時(shí)頻分析精度,得到反映真實(shí)爆破振動(dòng)屬性的時(shí)頻能量參數(shù)。對(duì)HHT進(jìn)行了改進(jìn),得到CEEMDAN·MPE-NHT時(shí)頻分析算法。該算法通過(guò)控制噪聲來(lái)改善EMD模態(tài)混淆,對(duì)IMF調(diào)頻分量進(jìn)行Hilbert變換,降低負(fù)值瞬時(shí)頻率出現(xiàn)的概率,從算法原理上實(shí)現(xiàn)提升信號(hào)時(shí)頻分析精度。

(2) 通過(guò)HHT和CEEMDAN·MPE-NHT算法仿真爆破振動(dòng)信號(hào)時(shí)頻分析對(duì)比研究可發(fā)現(xiàn),CEEMDAN·MPE-NHT算法將CEEMDAN的自適應(yīng)性和MPE隨機(jī)性檢測(cè)能力相結(jié)合,降低EMD對(duì)噪聲的敏感性。不僅可有效抑制噪聲引起的EMD模態(tài)混淆,同時(shí)結(jié)合NHT可得到時(shí)頻分辨率雙高的信號(hào)頻譜圖。

(3) 根據(jù)結(jié)構(gòu)共振放大系數(shù)D的原理,發(fā)現(xiàn)單個(gè)IMF分量主頻和房屋不同陣型下固有頻率越接近時(shí),對(duì)應(yīng)質(zhì)點(diǎn)振動(dòng)速度放大效應(yīng)越明顯。當(dāng)出現(xiàn)極限情況,即當(dāng)β=1,對(duì)應(yīng)質(zhì)點(diǎn)振動(dòng)速度增幅可高達(dá)10倍?？蓳?jù)此解釋房屋在低振速情況下,出現(xiàn)大量的爆破振動(dòng)裂縫的情形。