基于CEEMD和小波包閾值的組合降噪及泄流結(jié)構(gòu)的模態(tài)識別方法

胡劍超， 練繼建， 馬 斌， 董霄峰

(天津大學(xué) 水利工程仿真與安全國家重點實驗室，天津 300072)

基于CEEMD和小波包閾值的組合降噪及泄流結(jié)構(gòu)的模態(tài)識別方法

胡劍超， 練繼建， 馬 斌， 董霄峰

(天津大學(xué) 水利工程仿真與安全國家重點實驗室，天津 300072)

水流荷載激勵下，泄流結(jié)構(gòu)模態(tài)識別的實測振動信號常受噪聲干擾。以往EEMD算法有添加白噪聲造成非標(biāo)準(zhǔn)IMF導(dǎo)致的模態(tài)分裂問題及IMF有噪聲殘余不能完整重構(gòu)信號的完備性問題。而完備總體經(jīng)驗?zāi)B(tài)分解(CEEMD)法，通過在信號分解的每一層面添加特定高斯白噪聲，利用分解后第一階分量加總平均得到唯一余量計算固有模態(tài)函數(shù)，克服了EEMD的缺點；同時提出CEEMD與小波包閾值結(jié)合的組合降噪方法，運用到向家壩水彈性模型實測振動信號降噪中，驗證了該組合方法降噪的有效性；為了提高帶噪振動響應(yīng)模態(tài)識別的精度，基于數(shù)據(jù)相關(guān)技術(shù)，利用Markov參數(shù)構(gòu)造相關(guān)矩陣R，用該相關(guān)矩陣重構(gòu)Hankel矩陣后SVD分解得到系統(tǒng)最小實現(xiàn)，即數(shù)據(jù)相關(guān)特征系統(tǒng)實現(xiàn)法(ERA/DC)。將濾波降噪結(jié)合模態(tài)識別的整套方法，應(yīng)用到錦屏一級拱壩的泄流實測振動響應(yīng)中，得到了較好的應(yīng)用效果。

泄流結(jié)構(gòu)；濾波降噪；完備總體經(jīng)驗?zāi)B(tài)分解；小波包；模態(tài)識別；相關(guān)特征系統(tǒng)實現(xiàn)法

水工結(jié)構(gòu)的模態(tài)參數(shù)識別是確認(rèn)結(jié)構(gòu)運行健康狀態(tài)的重要技術(shù)手段。由于數(shù)據(jù)采集環(huán)境復(fù)雜以及大型泄流結(jié)構(gòu)激勵難等狀況，必須消除或降低振動響應(yīng)信號中的噪聲影響并在未知輸入激勵條件下進(jìn)行水工結(jié)構(gòu)的模態(tài)識別。

李成業(yè)等[1-3]研究的EMD/EEMD(Empirical Mode Decomposition With Ensemble EMD)+小波的濾波降噪方法具有諸多優(yōu)點：EEMD具有克服模態(tài)混疊、完全自適應(yīng)及IMFs都是單分量函數(shù)的特點，能彌補小波在處理非平穩(wěn)信號上的不足；同時小波變換具有多分辨率分析的優(yōu)點，利用小波閾值函數(shù)濾波，克服了EMD/EEMD濾除噪聲時，容易損失有用信號的缺點。

然而該組合方法仍有不足之處。EEMD通過原始信號中添加高斯白噪聲改進(jìn)EMD算法[4]。鄭近德等[5]指出這種方法使得分解后的信號包含殘余噪聲；每步構(gòu)造添加的白噪聲不同造成IMF分解數(shù)量不同，使得總體平均后的IMF不符合標(biāo)準(zhǔn)定義產(chǎn)生模態(tài)分裂；EEMD不能完整的重構(gòu)原始信號，不具完備性?；诙喑叨确治龅男〔ㄗ儞Q，通過尺度由大到小、時域-頻域變換的過程，由粗到精的展現(xiàn)信號的特征。然而由文獻(xiàn)[6]中小波算法可知，其只能對低頻小波細(xì)節(jié)子空間進(jìn)一步分解，而信號的高頻部分則沒有予以考慮。

由于方法自身的局限性，使得EEMD+小波閾值組合去噪方法在泄流結(jié)構(gòu)振動響應(yīng)分析中受到一定的限制。為了克服EEMD的缺點，Torres等[7]提出完備總體模態(tài)經(jīng)驗分解法(CEEMD)，而小波包變換[8](Wavelet Dacket Transform, WPT)既對信號低頻尺度又對高頻部分進(jìn)行精細(xì)分解。本文提出的基于CEEMD和小波包閾值的組合濾波降噪方法，不僅繼承了EEMD加小波方法的優(yōu)點，還改進(jìn)了該方法不足之處。文中通過向家壩水彈性模型實測振動信號降噪，證明改進(jìn)組合降噪方法更加有效。

Lian等[9]運用特征系統(tǒng)實現(xiàn)法(Eigensystem Realization Algolithm, ERA)，根據(jù)水工結(jié)構(gòu)工作狀態(tài)下只受環(huán)境激勵(如水流脈動荷載)的特點，僅用拱壩泄流響應(yīng)成功進(jìn)行了模態(tài)識別。ERA實質(zhì)是利用脈沖或自由衰減響應(yīng)數(shù)據(jù)，構(gòu)造廣義Hankel矩陣進(jìn)行奇異值SVD分解，尋找系統(tǒng)的最小實現(xiàn)。Hankel矩陣的維數(shù)與模型階次相關(guān)聯(lián)，影響著系統(tǒng)矩陣A是否為系統(tǒng)最小實現(xiàn)，最終決定模態(tài)識別精度。若系統(tǒng)輸入響應(yīng)中混有噪聲時，對Hankel矩陣影響較大[10]。Juang[11-12]用相關(guān)響應(yīng)數(shù)據(jù)R代替脈沖響應(yīng)數(shù)據(jù)，構(gòu)造數(shù)據(jù)相關(guān)Hankel矩陣進(jìn)行SVD分解，即數(shù)據(jù)相關(guān)特征系統(tǒng)實現(xiàn)法(ERA/DC)。本文結(jié)合泄流激勵下水工結(jié)構(gòu)振動特點，將ERA/DC引入泄流結(jié)構(gòu)模態(tài)參數(shù)識別中。并將基于組合降噪的ERA/DC應(yīng)用于錦屏一級拱壩實測振動響應(yīng)分析中，得到了精確的識別結(jié)果。

1 完備總體模態(tài)經(jīng)驗分解法CEEMD

以EEMD為基礎(chǔ)，CEEMD算法原理如下：

步驟1：x[n]為離散數(shù)n的原始信號；定義算子Ej(·)表示信號EMD分解后第j階IMF；wi[n](i=1,…,I)是單位方差零均值高斯白噪聲；εk系數(shù)是各階段信噪比；

(1)

步驟4：對于一階余量r1[n]中添加白噪聲信號，構(gòu)造信號r1[n]+ε1E1(wi[n])并進(jìn)行EMD分解，加總平均后得到第二階固有模態(tài)函數(shù)分量：

(2)

步驟6：對k階余量rk[n]中添加白噪聲信號，有構(gòu)造信號rk[n]+εkEk(wi[n])，進(jìn)行EMD分解，得到第k+1個固有模態(tài)函數(shù)分量：

(3)

步驟7：回到步驟(5)，計算下一個k。算法構(gòu)造循環(huán)取第一階IMFs分量加總平均，避免了EEMD同時分解的IMFs數(shù)量不同造成的模態(tài)分裂問題；

步驟8：重復(fù)步驟5～步驟7，直到余量不能再分解為止，最終得到：

(4)

式中：R是最終余量；K是固有模態(tài)分量總數(shù)。目標(biāo)信號x[n]可以表示為

(5)

式(5)為原始信號x[n]得到精確重構(gòu)，即該方法是具有完備性的分解方法。

下面結(jié)合一簡單算例驗證CEEMD方法的優(yōu)越性。設(shè)采樣頻率為200 Hz、采樣時間5 s的原始信號x，混入白噪聲k后得到含噪信號z，表達(dá)式如下：

z=x+k=3sin(πt)+5sin(2πt)+4.5sin(8πt)+5randn(1,1 000)

信號x，z如圖1，染噪信號z分別用EMD、EEMD、CEEMD方法濾波降噪(EEMD與CEEMD添加相同標(biāo)準(zhǔn)差與數(shù)量的白噪聲，一般取0.01倍～0.5倍原信號標(biāo)準(zhǔn)差，添加次數(shù)100～300次)，濾波后與x對比，如圖2所示。

(a)原始信號x(b)混入白噪聲后信號2

圖1 原始信號x與混入白噪聲后信號z

Fig.1 Original signalxand signalzmixed white noise

圖2 三種方法濾波降噪信號圖(局部)

從圖2(a)、圖2(b)可知，EMD或EEMD 方法濾波降噪都有不足之處，而圖(c)中CEEMD濾波降噪的重構(gòu)信號與原信號x波形最為接近。CEEMD方法改進(jìn)了EEMD不足之處，得到了顯著優(yōu)于EMD以及EEMD的濾波結(jié)果。

2 小波包(WPT)閾值濾波方法

小波包變換閾值濾波效果主要受閾值的影響。閾值函數(shù)可分為硬閾值函數(shù)和軟閾值函數(shù)。軟閾值降噪可實現(xiàn)最優(yōu)估計，保證降噪后的重構(gòu)信號與原信號有同樣的平滑，且不會產(chǎn)生附加震蕩，其表達(dá)式為

(6)

3 CEEMD與WPT閾值的組合濾波方法

壩體振動主要受表中孔水流脈動荷載激發(fā)，該脈動荷載是具有一定頻帶寬度的平穩(wěn)隨機(jī)過程，表現(xiàn)為有色白噪聲譜的性質(zhì)。壩體動力響應(yīng)通常是低頻或平穩(wěn)的振動信號，而環(huán)境噪聲常為非平穩(wěn)的高頻信號。

CEEMD將實測振動響應(yīng)信號的不同頻率模量由高到低自適應(yīng)清晰分離，分離出的帶有噪聲污染的信號固有模態(tài)分量，可利用小波包閾值的濾波方法有效去除噪聲。該方法適用于水工結(jié)構(gòu)的振動信號濾波降噪處理，有利于下一步水工結(jié)構(gòu)模態(tài)精確識別。

基于CEEMD和WPT的組合濾波降噪方法步驟如下：

下面結(jié)合向家壩消力池內(nèi)導(dǎo)墻水彈性模型驗證文中組合濾波方法的有效性，如圖3所示。

圖3 向家壩水彈模型隔振系統(tǒng)與右消力池內(nèi)導(dǎo)墻模型

Fig.3 Vibration isolation system of Xiang JiaBa hydroelastic model and Guide wall hydroelastic model in the right stilling basin

為了研究向家壩水電站泄洪誘發(fā)場地振動規(guī)律，天津大學(xué)水利工程系建造了1∶80比尺的向家壩水彈性模型。為了減小外界振源對壩體振動信號的干擾，采取在水彈模型混凝土基礎(chǔ)四周設(shè)置隔振溝、模型上下游連接段由止水橡膠軟連接等隔振措施，見圖3。試驗證明，該減振隔振系統(tǒng)有效減小了隔振溝外振動源的干擾[13]。

向家壩采用消力池消能，見圖3，右消力池內(nèi)設(shè)置一導(dǎo)墻水彈模型。導(dǎo)墻長寬高為99 cm×4 cm×38 cm。材料模擬碾壓混凝土，彈模為105 MPa，泊松比0.167。布設(shè)4個動位移測點。選擇測點4在上游水位370 m、表孔開度6 m工況實測導(dǎo)墻動位移響應(yīng)數(shù)據(jù)進(jìn)行分析，采樣頻率200 Hz，取10 s數(shù)據(jù)分析。

圖4 實測導(dǎo)墻振動信號與添加噪聲后振動信號

由于設(shè)置減振隔振系統(tǒng)，該動位移實測響應(yīng)數(shù)據(jù)x(n)的“毛刺”很少，信號較為光滑，可認(rèn)為其基本不含噪聲信號，如圖4(a)。在實測信號x(n)中加入白噪聲k，受噪聲污染信號z=x(n)+k(SNR=3.39 dB)，如圖4(b)所示，針對染噪信號z，分別用CEEMD、小波包(WPT)閾值、EEMD+小波、CEEMD+WPT四種方法進(jìn)行濾波降噪。

圖5為各方法濾波重構(gòu)波形與原信號x對比，表1為信噪比SNR、濾波信號與原始信號x互相關(guān)系數(shù)。

圖5 四種方法濾波降噪對比圖(局部)

信號類型降噪方法信噪比SNR/dB互相關(guān)系數(shù)構(gòu)造染噪信號z(SNR=3.39)CEEMD9.9070.9475小波包(WPT)11.2420.9626EEMD+小波11.8310.9666CEEMD+WPT13.0060.9747

由四種方法濾波降噪后的振動信號與原始信號x(n)對比圖5看出：CEEMD分解加帶通濾波的做法，得到的SNR較低，重構(gòu)信號仍有較多毛刺，說明加入的噪聲有殘留且部分有用信號也被一并去除；小波包閾值濾波的信噪比大于CEEMD，但重構(gòu)波形仍存在毛刺；EEMD+小波的濾波方法，SNR較高，信號波形重構(gòu)較好、毛刺基本消除，說明噪聲基本被濾除，原始信號保留較好；而CEEMD+WPT的組合濾波降噪振動信號，白噪聲“毛刺”不僅被很好的濾除干凈，從信噪比SNR=13.006 dB、相關(guān)系數(shù)=0.974 7來看，原始振動信號x(n)保留最好、降噪效果最佳。

4 基于ERA/DC的泄流結(jié)構(gòu)模態(tài)參數(shù)識別方法

4.1數(shù)據(jù)相關(guān)特征系統(tǒng)實現(xiàn)法(ERA/DC)算法原理

根據(jù)振動理論，線性系統(tǒng)的脈沖響應(yīng)函數(shù)矩陣，即馬可夫參數(shù)(Markov)為

Y0=D,Y1=CB,Y2=CAB,…,

Y(k)=CAk-1B

(7)

式中：Y(k)為l×m階矩陣，又是脈沖響應(yīng)函數(shù)矩陣。A為系統(tǒng)矩陣、B為輸入矩陣、C為輸出矩陣。由有限維線性時不變系統(tǒng)狀態(tài)方程知，系統(tǒng)模態(tài)識別問題即為系統(tǒng)的最小實現(xiàn)問題，即求得矩陣A，B，C并使其階次最小。由Markov參數(shù)Y(k)，構(gòu)造分塊Hankel矩陣：

(8)

式中：H(k-1)為lξ×mη階矩陣，ξ和η為任意整數(shù)。當(dāng)令k取1時，有Hankel矩陣H(0)。在ERA方法中，會針對H(0)進(jìn)行SVD分解。而在ERA/DC方法中，則構(gòu)建γ=lξ維方陣R(q)：

R(q)=H(q)HT(0)

(9)

式中：數(shù)據(jù)相關(guān)矩陣R(q)為lξ×lξ維，Hankel矩陣H(0)、H(q)為lξ×mη階矩陣。當(dāng)Hankel矩陣的列數(shù)很大時，R(q)的維數(shù)將?H(q)。若當(dāng)q=0時，有：

R(0)=H(0)HT(0)=

(10)

構(gòu)建(α+1)×(β+1)階塊相關(guān)矩陣U(q)：

U(q)=

(11)

式中：q為避免有噪聲存在時相關(guān)矩陣產(chǎn)生偏移的整數(shù)；r為防止相鄰R塊產(chǎn)生重疊的整數(shù)；α和β為分析中包含多少相關(guān)延遲。針對U(0)矩陣SVD分解：

U(0)=PDQT

(12)

式中：左奇異矩陣P和右奇異矩陣Q的維數(shù)為(α+1)γ×n和(β+1)γ×n，奇異值對角矩陣D為n×n階方陣。構(gòu)造U(0)偽逆矩陣U?，代入式(16)～式(18)推得：

A=D-1/2PTU(1)QD1/2

B=[EγTPD1/2]?H(0)EmC=EpTPD1/2

(13)

式中：EγT=[Iγ0γ]，EmT=[Im0m]，EpT=[Ip0p]，Iγ表示γ階單位矩陣，0γ表示γ階零矩陣。

針對系統(tǒng)矩陣A進(jìn)行奇異值分解，求得特征值矩陣Z=diag(Z1,Z2,…,Zn)，系統(tǒng)特征值與矩陣A特征值有關(guān)系：λi=ln(Zi)/T，T為采樣間隔。求得模態(tài)參數(shù)為

有阻尼頻率：ωdi=imag(λi)

阻尼比：ξi=-Re(λi)/|λi|

(14)

式中：Re(·)為取實部，imag(·)為取虛部，|·|表示取模。ERA/DC與ERA算法的主要區(qū)別在于首先進(jìn)行數(shù)據(jù)相關(guān)運算后重構(gòu)Hankel矩陣U(0)，這樣可壓縮矩陣規(guī)模，減少奇異值分解的計算時間，并減少噪聲影響。相較于ERA，ERA/DC具有更好的識別精度。

4.2基于組合降噪的ERA/DC模態(tài)識別方法

ERA/DC通過Markov參數(shù)(即為脈沖響應(yīng)函數(shù))進(jìn)行模態(tài)參數(shù)識別。水工結(jié)構(gòu)在泄流激勵下振動時，可假定為線性時不變系統(tǒng)，利用自然激勵技術(shù)(NExT)將參考點和響應(yīng)點數(shù)據(jù)變換為測點間的互相關(guān)函數(shù)，得到結(jié)構(gòu)的脈沖響應(yīng)函數(shù)。

在未知輸入激勵下，結(jié)構(gòu)的振動模態(tài)階次通常是未知的。原始信號由延時嵌陷技術(shù)重構(gòu)的吸引子軌道矩陣D，與ERA方法中Hankel矩陣H(0)的構(gòu)造形式一致，基于奇異熵定階技術(shù)，H(0)的定階(即系統(tǒng)的定階)可根據(jù)奇異熵增量的變化進(jìn)行。觀察式11可知，ERA/DC用相關(guān)R矩陣構(gòu)造的Hankel矩陣，當(dāng)q=0時，U(0)的矩陣形式與H(0)相同，因此由ERA/DC進(jìn)行結(jié)構(gòu)模態(tài)識別，同樣可根據(jù)奇異熵增量譜確定系統(tǒng)階次。

為了進(jìn)一步減少實測信號中噪聲影響，模態(tài)識別之前須對響應(yīng)數(shù)據(jù)運用CEEMD+WPT的組合濾波方法進(jìn)行降噪，這就是基于組合降噪方法的ERA/DC泄流結(jié)構(gòu)模態(tài)識別方法，其流程如圖6所示。

圖6 CEEMD加小波包組合濾波后的信號

為了說明ERA/DC方法針對帶噪信號識別的優(yōu)越性，建立混凝土懸臂梁有限元模型驗證。梁長20.0 m、截面尺寸1.4 m×0.6 m、動彈模30 GPa。梁的一端全約束固定，如圖7。為了充分激發(fā)梁的各階模態(tài)，沿梁軸線選擇六個點P1～P6施加負(fù)Y向的典型水流脈動荷載，荷載頻率100 Hz，時長30 s，如圖8所示。

圖7 懸臂梁有限元模型與荷載施加方向圖

圖8 典型水流脈動荷載時程線圖(局部)

計算懸臂梁動力響應(yīng)，提取懸臂梁D1～D6測點振動信號。以靠近懸臂梁固端的D6為參考點，其與梁中間D4測點的時程線如圖9，用NExT法計算D4-D6的脈沖響應(yīng)函數(shù)如圖10(a)所示。

得到脈沖響應(yīng)函數(shù)后計算奇異熵增量隨階次的變化如上圖10(b)所示，奇異熵階次為9階時，熵增量開始緩慢增長并趨于平穩(wěn)，說明信號的有效特征信息量已經(jīng)趨于飽和，特征信息已經(jīng)基本完整。剔除特征值中的非模態(tài)項(非共軛根)和共軛項(重復(fù)項)后，利用模態(tài)置信因子MAC準(zhǔn)則判別虛假模態(tài)[14]。得出結(jié)構(gòu)的前3階振動階次進(jìn)行對比，如表2所示。

圖9 測點D4與D6振動時程線圖(局部)

(a)脈沖響應(yīng)函數(shù)(b)脈沖響應(yīng)奇異熵增量

圖10 D4-D6點的脈沖響應(yīng)函數(shù)與奇異熵增量譜圖

由表2可知，ERA法識別頻率誤差<2%，而 ERA/DC方法頻率識別誤差<1%，都準(zhǔn)確識別出了懸臂梁前三階模態(tài)。

圖11 添加噪聲脈沖函數(shù)信號與濾波降噪信號

為了對比振動信號受到噪聲污染情況下，ERA與ERA/DC方法的識別精度，將D4-D6脈沖響應(yīng)函數(shù)添加白噪聲(信噪比SNR=20 dB)，信號如圖11(a)所示。運用兩種方法模態(tài)識別，結(jié)果如表3所示。

表3帶噪脈沖響應(yīng)信號的模態(tài)頻率識別對比(SNR=20)

Tab.3Modalfrequencyidentificationofpulse-responsesignalmixednoise(SNR=20dB)

階數(shù)數(shù)值模擬頻率ERA方法ERA/DC方法頻率/Hz誤差/%頻率/Hz誤差/%一階1.2121.2805.6111.2222.887二階2.8102.9946.5482.8261.281三階7.5418.1027.4397.5521.366

由表3看出，信號受噪聲污染后，ERA識別懸臂梁三階模態(tài)頻率，相較于無干擾情況誤差顯著增大。其一階、三階誤差都>5%，第三階誤差>7%。ERA/DC識別中，三階誤差都<3%，表明其模態(tài)識別抗噪聲干擾能力強(qiáng)于ERA方法。

將該染噪脈沖函數(shù)信號(D4點)用CEEMD+小波包組合降噪，結(jié)果如圖11(b)所示。濾波后信號信噪比SNR=25.538、相關(guān)系數(shù)0.998 6。模態(tài)識別結(jié)果，如表4所示。

表4濾波脈沖函數(shù)響應(yīng)信號模態(tài)頻率識別對比

Tab.4Frequencyidentificationoffilteringpulse-responsesignal

階數(shù)數(shù)值模擬頻率ERA方法ERA/DC方法頻率/Hz誤差/%頻率/Hz誤差/%一階1.2121.2361.9801.2020.825二階2.8102.7472.2422.8441.210三階7.5417.7632.9447.5200.278

染噪信號濾波后識別結(jié)果由表4可知：ERA方法識別中，頻率誤差降至<3%，ERA/DC頻率識別誤差降至<1.5%，證明該組合降噪方法與ERA/DC方法結(jié)合的有效性與必要性。

5 錦屏一級拱壩實測振動信號應(yīng)用分析

錦屏一級水電站是雅礱江干流下游河段的控制性梯級電站。其最大壩高305 m，是目前世界最高的混凝土雙曲拱壩。2014年天津大學(xué)水利工程系針對錦屏一級樞紐工程進(jìn)行了水力學(xué)原型觀測。沿壩頂1#～24#壩段，水平向順?biāo)鞑荚O(shè)H1～H7共7個壩頂動位移傳感器。水庫汛期泄洪期間，實測了不同泄流工況下，泄流激勵下壩體振動響應(yīng)。傳感器測點布置圖及錦屏一級拱壩有限元模型，如圖12所示。

圖12 傳感器測點布置及錦屏一級拱壩有限元模型圖

5.1濾波降噪

流激振動造成的拱壩動位移均方根通常為幾十個微米。由于振動量級較小，混入的干擾信號和背景噪聲容易淹沒有用信號。分析各工況壩體振動規(guī)律：拱冠梁處振動最大但遠(yuǎn)離干擾源，壩肩振動量小且易受噪聲干擾。因此選取工況7(3#表孔25%開度、泄流量329 m3/s、上游水位1 879.2 m、下游1 645.6 m)時，拱冠梁附近參考測點H3及右壩肩附近測點H6，采樣頻率200 Hz，時長40 s，其時程線及歸一化功率譜，如圖13所示。

圖13 H3、H6測點振動信號(局部)及歸一化功率譜圖

由錦屏有限元模態(tài)計算得拱壩前六階模態(tài)頻率，頻率范圍在1.0～3.5 Hz之間。測點H3位于拱冠梁附近，受到噪聲干擾小，由其功率譜圖13(b)可知其振動能量集中在1.0～3.0 Hz頻帶部分，與有限元計算相符。

在現(xiàn)場泄流觀測期間，發(fā)現(xiàn)右壩肩附近存在少量施工作業(yè)。壩肩測點H6振動信號“毛刺”較多，并出現(xiàn)漂移現(xiàn)象，運用最小二乘擬合法消除趨勢項后如圖13(c)。功率譜如圖13(d)，<0.5 Hz的低頻部分有能量集中，且>4.5 Hz處有寬頻白噪聲存在，直接應(yīng)用H6實測振動數(shù)據(jù)進(jìn)行模態(tài)識別會產(chǎn)生很大誤差。

為了得到拱壩真實振動信號，運用文中提出的組合濾波降噪方法濾除原始信號中的噪聲。將H6測點信號CEEMD分解，添加白噪聲標(biāo)準(zhǔn)差1.494 μm，添加數(shù)300次，最大迭代數(shù)5 000次。分解得IMF1～I(xiàn)MF13分量，IMFs分量從高頻到低頻依次排列如上圖14。

圖14 H6測點實測信號CEEMD分解圖

做出IMF1～I(xiàn)MF12分量歸一化功率譜圖如下圖15所示。前兩階IMF能量集中在高頻且頻帶很寬，為寬頻白噪聲分量，運用低通濾波進(jìn)行濾除。IMF3～I(xiàn)MF6是振動信號分量，但其中混入了大量噪聲信號，需要運用WPT方法進(jìn)一步濾除噪聲。IMF7～I(xiàn)MF9分量振動能量集中在1.0～3.0 Hz附近，判斷其屬于結(jié)構(gòu)真實振動信號。而IMF10～13為低頻趨勢項干擾，予以高通濾除。

IMF3～I(xiàn)MF6混有白噪聲，運用小波包閾值方法進(jìn)一步降噪。小波包分解采用dbN小波基，共N層分解，不同IMF分量采用不同的N值。采用軟閾值方法，默認(rèn)全局閾值。濾波后，由IMF3～I(xiàn)MF9重構(gòu)信號如圖16所示。

從圖16可知，原始信號毛刺被去除。從功率譜曲線可知，低頻干擾噪聲和高頻白噪聲被有效消除了，響應(yīng)中反映結(jié)構(gòu)模態(tài)信息的共振分量所占的比重大大增強(qiáng)。

圖15 H6信號IMF分量歸一化功率譜圖

圖16 H6組合濾波降噪振動信號(局部)及功率譜圖

5.2模態(tài)識別

由圖17可知，H6測點實測振動響應(yīng)濾除噪聲信號后，以拱冠梁附近的測點H3為參照點，運用NExT方法計算H3與H6測點脈沖響應(yīng)函數(shù)；由H6-H3脈沖響應(yīng)函數(shù)，首先構(gòu)造4 000×4 000階Hankel矩陣H(0)以及2 000×2 000階Hankel矩陣U(0)，計算奇異熵增量隨階次變化。確定奇異熵階次為11階，采用文中ERA/DC方法進(jìn)行識別。剔除非模態(tài)項和共軛項后用MAC準(zhǔn)則判別虛假模態(tài)，得到結(jié)構(gòu)前6階振動階次，識別結(jié)果列于表5。表5還列出了拱壩有限元模態(tài)計算結(jié)果，并將該振動信號同時運用ITD、STD、Prony(復(fù)指數(shù)法)、ARMA(時序分析法)、SSI(隨機(jī)子空間法)與ERA(特征系統(tǒng)實現(xiàn)法)進(jìn)行模態(tài)識別，識別結(jié)果一同列于表5中。

圖17 H3-H6測點脈沖響應(yīng)函數(shù)與奇異熵增量譜

表5 錦屏一級拱壩模態(tài)頻率識別對比

拱壩是一種典型的密頻結(jié)構(gòu)，綜合各方法的模態(tài)識別的結(jié)果，可得到錦屏一級拱壩的前六階工作頻率。可以看出，ERA/DC方法對于低頻、密頻具有較強(qiáng)的識別能力。

將各方法識別結(jié)果進(jìn)行對比，發(fā)現(xiàn)識別結(jié)果都較為接近，都能識別出拱壩前六階頻率。而文中所提方法與拱壩有限元計算結(jié)果同樣相吻合。文中方法的計算結(jié)果不僅能夠滿足工程分析的需要，且具有相對更好的計算精度。

6 結(jié) 論

本文針對泄流結(jié)構(gòu)振動響應(yīng)信號，提出一種組合濾波降噪結(jié)合模態(tài)識別的整套方法，應(yīng)用到實際泄流結(jié)構(gòu)的狀態(tài)參數(shù)識別中，得到以下結(jié)論：

(1) CEEMD解決了EEMD算法導(dǎo)致的模態(tài)分裂及完備性問題，小波包能獲得豐富的時頻局部化信息，更適合突變信號、非平穩(wěn)信號的分析與濾波。向家壩水彈性模型實測振動信號降噪說明，CEEMD+小波包閾值的組合降噪方法較EEMD+小波方法具有更好的濾波降噪效果。

(2) 基于數(shù)據(jù)相關(guān)技術(shù)的ERA/DC方法克服了泄流結(jié)構(gòu)輸入激勵條件難獲得的問題，利用懸臂梁振動模型進(jìn)行驗證，表明其對含噪信號具有比ERA更好的模態(tài)識別精度，同時說明了其與組合降噪方法結(jié)合的必要性。

(3) 將CEEMD+小波包閾值組合降噪方法結(jié)合ERA/DC模態(tài)識別方法，應(yīng)用到錦屏一級拱壩泄流實測振動響應(yīng)信號中，得到了很好的實際應(yīng)用效果。

[1] 李成業(yè)，練繼建，劉昉，等. EMD與小波閾值聯(lián)合濾波方法的改進(jìn)及其在泄流結(jié)構(gòu)振動分析中的應(yīng)用[J]. 振動與沖擊，2013，32(9)：63-70.

LI Chengye, LIAN Jijian, LIU Fang, et al. An improved filtering method based on EMD and wavelet-threshold and its application in vibration analysis for a flood discharge structure[J]. Journal of vibration and shock, 2013, 32(9): 63-70.

[2] ZHANG Yan, LIAN Jijian, LIU Fang. An improved filtering method based on EEMD and wavelet-threshold for modal parameter identification of hydraulic structure[J]. Mechanical Systems and Signal Processing, 2016, 68: 316-329.

[3] 邵忍平，曹精明，李永龍. 基于EMD小波閾值去噪和時頻分析的齒輪故障模式識別與診斷[J]. 振動與沖擊，2012，31(8)：96-101.

SHAO Renping, CAO Jingming, LI Yonglong. Gear fault pattern identification and diagnosis using time-frequency analysis and wavelet threshold de-noising based on EMD[J]. Journal of vibration and shock, 2012, 31(8): 96-101.

[4] WU Z H, HUANG N E. Ensemble empirical mode decomposition: A noise assisted data analysis method[J]. Advances in Adaptive Data Analysis, 2009, 1(1): 1-41.

[5] 鄭近德，程軍圣，楊宇. 改進(jìn)的EEMD算法及其應(yīng)用研究[J]. 振動與沖擊，2013，32(21)：21-26.

ZHENG Jinde, CHENG Junsheng, YANG Yu. Modified EEMD algorithm and its applications[J]. Journal of Vibration and Shock, 2013, 32(21): 21-26.

[6] OUADFEUL S A, ALIOUANE L. Random seismic noise attenuation data using the discrete and continuous wavelet transforms[J]. Arabian Journal of Geosciences, 2014, 7(7): 2531-2537.

[7] TORRES M E, COLOMINAS M A, SCHLOTTHAUER G, et al.A complete ensemble mode decomposition with adaptive noise[C] // IEEE International Conference on Acoustics. Speech and signal Processing(ICASSP). USA, NJ: IEEE. 2011. 4144-4147.

[8] HEMMATI F, ORFALI W, GADALA M S. Roller bearing acoustic signature extraction by wavelet packet transform, applications in fault detection and size estimation[J]. Applied Acoustics, 2016, 104(1): 101-118.

[9] LIAN J J,LI H K,ZHANG J W. ERA modal identification method for hydraulic structures based on order determination and noise reduction of singular entropy[J]. Science in China series E-Technological Sciences, 2009, 52(2): 400-412.

[10] 李蕾紅，陸秋海，任革學(xué).特征系統(tǒng)實現(xiàn)法的識別特性研究及算法的推廣[J]. 工程力學(xué)，2002，19(1)：109-115.

LI Leihong, LU Qiuhai, REN Gexue. A study of eigensystem realization algorithm and its generalization[J]. Engineering Mechanics, 2002, 19(1): 109-115.

[11] JUANG J N, COOPER J E, WRIGHT J R. An eigensystem realization algorithm using data correlation (ERA/DC) for modal parameter identification[J]. Control Theory and Advanced Technology, 1988, 4(1): 5-14.

[12] CHIANG D Y, LIN C S. Identification of modal parameters from ambient vibration data using eigensystem realization algorithm with correlation technique[J]. Journal of Mechanical Science and Technology, 2010, 24(12): 2377-2382.

[13] 李淑君，練繼建，歐陽秋平. 高壩泄流誘發(fā)壩區(qū)及場地振動的振源分析研究[J]. 水利水電技術(shù)，2014，45(9)：47-51.

LI Shujun, LIAN Jijian, OUYANG Qiupin. Study on source of dam discharge-induced vibration in dam region and site[J]. Water Resources and Hydropower Engineering, 2014, 45(9): 47-51.

[14] 董霄峰，練繼建，楊敏，等. 海上風(fēng)機(jī)結(jié)構(gòu)工作模態(tài)識別的組合降噪方法[J]. 天津大學(xué)學(xué)報(自然科學(xué)與工程技術(shù)版)，2015，48(3)：203-208.

DONG Xiaofeng, LIAN Jijian, YANG Min, et al. A compound de-noising method for operational modal identification of offshore wind turbine structure[J]. Journal of Tianjin University(Science and Technology),2015,48(3):203-208.

Ade-noisingandmodalidentificationcombinedmethodbasedonCEEMDandwaveletpacketthresholdforflooddischargestructures

HU Jianchao， LIAN Jijian， MA Bin， DONG Xiaofeng

(State Key Laboratory of Hydraulic Engineering Simulation and Safety, Tianjin University, Tianjin, 300072 China)

Vibration signals measured for modal identification of a flood discharge structure under water load excitation are often interfered by noise. The EEMD algorithm has problems of modes splitting caused by nonstandard IMFs produced due to adding white noise and IMFs with residual noise not being able to completely reconstruct original signals. Here, the complete ensemble empirical mode decomposition (CEEMD) was proposed. Using it, the drawbacks of EEMD were overcome through adding a specific white noise to each stage of signal decomposition and using the unique residual obtained by adding the average to the first order component gained after decomposition to calculate IMFs. Furthermore, a de-noising method with CEEMD and the wavelet packet threshold combined was proposed to be applied in de-noising vibration signals measured from Xiangjia dam’s hydro-elastic model, the effectiveness of this method was verified. In order to improve the precision of modal identification with the noised vibration responses, based on the data correlation technology, Markov parameters were used to construct a correlation matrix R. The R matrix was employed to reconstruct a Hankel matrix. After its SVD, the minimum system realization was obtained. This algorithm was called the eigen-system realization algorithm with data correlation (ERA/DC). Finally, the proposed de-noising and modal identification combined method was applied in the measured vibration responses of Jinping arch dam to get a better application effect.

flood discharge structure; noise reduction; CEEMD; wavelet packet; modal identification; ERA/DC

國家重點基礎(chǔ)研究發(fā)展計劃(973計劃)子課題“極端條件下高碾壓混凝土壩動力響應(yīng)機(jī)制”(2013CB035905-1)

2016-06-06 修改稿收到日期：2016-08-08

胡劍超 男，碩士生，1990年生

練繼建 男，教授，博士生導(dǎo)師，1965年生

TV65

： A

10.13465/j.cnki.jvs.2017.17.001