劉 昉， 張魯豐， 龐博慧， 梁 超， 姚 燁

(1.天津大學 水利工程仿真與安全國家重點實驗室，天津 300354；2.華能瀾滄江水電股份有限公司，云南 昆明 650214)

0 引言

導墻作為輕型泄流結構，長期受水流、風等環(huán)境因素的干擾，會由于疲勞和腐蝕發(fā)生開裂損傷，而且導墻水下部分結構損傷不易被直接發(fā)現(xiàn)[1]。Taxakana、Navajo等水利樞紐的導墻都因水流誘發(fā)振動而破壞[2]。對導墻實施振動監(jiān)測是保證水電站安全運行的必要手段，而監(jiān)測數(shù)據(jù)受現(xiàn)場工作環(huán)境影響，會存在噪聲干擾，使得結構振動特征信息被淹沒，因此需要對信號進行降噪以獲得真實信號。

國內(nèi)外學者對信號降噪進行了很多研究，徐國賓等[3]利用小波變換有效地去除了水墊塘底板振動信號的噪聲，同時保留原始信號的細節(jié)部分。李成業(yè)等[4]利用經(jīng)驗模態(tài)分解(empirical mode decomposition，EMD)[5]與小波變換相結合的濾波方法，對拱壩振動信號進行降噪，提取振動特征。EMD無須像小波變換一樣設置基函數(shù)，克服了依賴主觀經(jīng)驗的缺點，但其存在模態(tài)混疊和端點效應問題[6]。胡劍超等[7]對實測振動信號進行處理，采用完備總體經(jīng)驗模態(tài)分解(complete ensemble empirical mode decomposition，CEEMD)與小波包閾值降噪結合的方法對信號進行濾波降噪，CEEMD[8]通過向原信號中添加正負成對的高斯白噪聲以解決EMD存在的模態(tài)混疊現(xiàn)象，但加入的白噪聲難以徹底消除，重構信號中含有殘余噪聲。自適應噪聲的完備總體經(jīng)驗模態(tài)分解(complete ensemble empirical mode decomposition with adaptive noise，CEEMDAN)[9]做出進一步改進，在EMD的每個階段自適應地添加白噪聲，計算唯一余量信號獲得各IMF分量，分解過程完備。重構誤差很小，消除了添加噪聲對信號的影響，提高了計算效率及準確性。Bai等[10]使用CEEMDAN方法、排列熵(permutation entropy，PE)和峰值濾波聯(lián)合去噪算法，能有效去除齒輪箱振動信號噪聲。

基于上述情況，本文使用CEEMDAN和多尺度排列熵(multi-scale permutation entropy，MPE)即CEEMDAN-MPE聯(lián)合降噪方法，對振動數(shù)據(jù)進行CEEMDAN處理；利用MPE檢驗各IMF分量的隨機性，篩分為含噪聲的IMF和純凈IMF；然后采用小波閾值降噪處理含噪聲的IMF分量；最后將處理后的含噪聲的IMF分量與純凈IMF分量進行重構，達到降噪目的，并將其應用于導墻的實測振動位移信號，驗證了該方法的有效性。

1 基本原理

1.1 CEEMDAN算法

CEEMDAN通過在EMD的每個階段添加自適應白噪聲，在較少平均次數(shù)下其重構誤差仍然很小，可以解決EEMD算法中的噪聲殘留問題，運算效率也有相應提高，具體步驟總結如下。

(1)在原始信號x(t)中添加白噪聲ε0ni(t)，其中ε0為第一次添加的白噪聲幅值系數(shù)，則第i次信號可以表示為xi(t)=x(t)+ε0ni(t)，對各xi(t)通過EMD，得到各個信號第一階分量imfi1和殘差ri1(t)。

(1)

(2)定義算子Ei(·)為使用EMD后的第j個分量，對r1(t)中添加白噪聲ε1E1[ni(t)]，然后進行分解，得到第二階分量：

(2)

以此類推，第k個殘余信號為

rk(t)=rk-1(t)-imfk。

(3)

(3)第(k+1)階分量為

(4)

(4)重復步驟(2)和(3)，直到信號無法繼續(xù)分解，運算終止，共得到K個分量，最終殘余量為

(5)

此時，原始信號可以表示為

(6)

1.2 MPE算法

MPE能夠檢測信號的突變性和隨機性，是一種判斷信號復雜度的指標[11]。與排列熵PE只能檢測單一尺度的隨機性相比，MPE能夠衡量信號不同尺度下的隨機性[12]，而泄流振動信號在不同尺度均包含重要信息，因此MPE更適用于該信號分析。計算步驟如下。

對時間序列Z={z1,z2,…,zL}進行多尺度粗?；幚恚?/p>

(7)

N=(L/s-(m-1)t)；

(8)

(9)

式中：m為嵌入維數(shù)；t為時間延遲。

(10)

式中：Pj為第j次符號出現(xiàn)的概率。

對式(10)計算排列熵進行歸一化處理：

(11)

熵值大小代表時間序列的隨機和復雜程度，值越接近于1，隨機性越強，非平穩(wěn)程度越高。參考文獻[13]，將排列熵閾值設置為0.6，大于或等于閾值的IMF分量視為含噪聲的分量，需要進一步處理。

1.3 小波閾值去噪

小波閾值去噪的原理為將信號經(jīng)過小波變換產(chǎn)生的小波系數(shù)中會含有信號的重要信息，真實信號的系數(shù)較大，噪聲信號的系數(shù)較小，通過對小波系數(shù)設置合適的閾值，小于閾值的小波系數(shù)認為是由噪聲引起，將其置0從而達到去噪目的[14]。具體步驟主要包括：

(1)選擇合適的小波基及分解層數(shù)，對含噪信號進行小波分解；

(2)選擇合適的閾值方法對各分解尺度下的小波系數(shù)進行閾值處理；

(3)對小波系數(shù)進行重構，完成降噪。

小波閾值處理主要有軟閾值去噪和硬閾值去噪。軟閾值去噪與硬閾值去噪相比，克服了不連續(xù)的缺陷，去噪后結果更平滑，因此本文采用軟閾值去噪方法，其數(shù)學表達式為

(12)

2 CEEMDAN-MPE聯(lián)合降噪

CEEMDAN-MPE聯(lián)合降噪方法的步驟如下。

首先，利用CEEMDAN將導墻振動信號進行分解，獲得若干IMF分量，計算各IMF的MPE值ZMPE，檢驗各IMF分量的隨機性，若MPE值大于等于設定的閾值0.6，為含噪聲的IMF分量；MPE值小于0.6則為純凈IMF。其次，對于含噪聲的IMF分量使用小波軟閾值方法進行去噪。最后，將小波軟閾值處理后的結果與純凈IMF分量進行重構，得到去噪后的信號，流程見圖1。

圖1 CEEMDAN-MPE降噪方法流程Figure 1 Flow chart of CEEMDAN-MPE denoising method

利用仿真信號驗證聯(lián)合降噪方法的效果。設置采樣時間為10 s，采樣頻率為100 Hz的仿真信號，其表達式為

y1=0.3sin πt+0.5sin 2πt+0.4sin 8πt。

(13)

加噪后信號表達式為

y2=y1+w。

(14)

式中:w=0.5randn(1,1 000)為白噪聲分量,由MATLAB軟件生成，白噪聲隨時間起伏變化快，其功率譜密度在頻帶內(nèi)是均勻的或者變換很小。

采用CEEMDAN-MPE方法對信號進行處理，得到降噪后結果見圖2。由圖2可知，添加的白噪聲已被去除，信號波形與原始信號的相似程度較高。

圖2 降噪后信號波形對比Figure 2 Comparison of signal waveform after denoising

為驗證降噪效果，分別使用CEEMDAN-相關系數(shù)方法及EEMD-MPE方法對信號進行處理，并與本文方法進行對比。含噪聲IMF分量與原信號相關性很差，因此可根據(jù)各IMF分量與原信號的相關系數(shù)來判定虛假分量[16]，相關系數(shù)計算公式為

(15)

式中：X為IMF分量；Y為原始信號；N為信號長度。

選擇信噪比(SNR)、均方根誤差(RMSE)及降噪信號與原始信號的相關系數(shù)作為評價降噪效果的指標。信噪比及相關系數(shù)越大，均方根誤差越小，則有較好的降噪效果。由表1可知本文采用方法降噪效果最優(yōu)。

表1 不同降噪方法的效果評價Table 1 Effect evaluation of different denoising methods

3 工程應用

3.1 工程概況

某水電站泄水建筑物由主河床的4孔泄洪閘和右岸導流明渠內(nèi)3孔泄洪閘組成。以泄流導墻為研究對象，在其頂面布置10個測點，具體布置見圖3，各測點布置三向低頻振動位移傳感器，導墻橫河向振動量大，順河向及垂向振動相對較小，所以主要對導墻橫河向振動進行研究。通過分析發(fā)現(xiàn),距離閘門較近的測點振動量相對較大，且易受到噪聲干擾，淹沒有用信號，影響監(jiān)測結果。因此，選取位于導墻斜坡段的測點2進行降噪處理。

圖3 導墻及測點布置Figure 3 Guide wall and measuring point layout

3.2 降噪處理與分析

測點2的振動位移時程線及功率譜密度如圖4所示。由圖4可知，振動信號有較多“毛刺”，存在噪聲，采用CEEMDAN-MPE聯(lián)合降噪方法對其進行處理。首先進行CEEMDAN分解，共得到13個IMF分量，IMF1～IMF13頻率依次降低。

圖4 測點2振動信號時程圖與功率譜密度圖Figure 4 Time history and power spectral density of vibration signal at measuring point 2

通過MPE方法計算各IMF分量的MPE值(ZMPE)，需要設置合適的參數(shù)，本文使用互信息法和偽近鄰法[17]計算時間延遲τ，選擇嵌入維數(shù)m。信號的互信息與時間延遲τ的關系見圖5。τ=1時，互信息值取得的第一個極小值，為最優(yōu)延遲時間。

圖5 互信息隨時間延遲變化曲線Figure 5 Curve of mutual information with time delay

在時間延遲τ基礎上，嵌入維數(shù)通過分析偽近鄰率隨嵌入維數(shù)的變化規(guī)律來確定。判據(jù)1閾值設置為15，判據(jù)2閾值設置為2，最大嵌入維數(shù)設置為10。偽近鄰率計算結果如圖6所示。

圖6 偽近鄰率隨嵌入維數(shù)變化曲線Figure 6 Curve of FNNP with embedding dimension

圖6中，判據(jù)1的偽鄰近率逐漸降低，判據(jù)2的偽鄰近率在嵌入維數(shù)為6之后不再隨嵌入維數(shù)的增加而減少，因此，最佳嵌入維數(shù)m=6。選定嵌入維數(shù)及時間延遲后，經(jīng)過試算，選擇尺度因子s=12。各IMF分量的ZMPE如圖7所示。IMF1～IMF13的ZMPE逐漸降低，其中IMF1～IMF6的ZMPE大于0.6，視為含噪聲的IMF分量；IMF7～IMF13的ZMPE值小于0.6，為純凈IMF分量。

圖7 IMF分量的ZMPEFigure 7 ZMPE of the IMF component

采用小波軟閾值去噪方法對IMF1～IMF6分量進行處理，選擇去噪效果較好的db4小波基，進行4層分解，固定閾值為9.673 0。將降噪后的IMF1～IMF6與IMF7～IMF13進行重構，得到降噪后的導墻振動信號如圖8所示。

圖8 降噪后信號時程圖與功率譜密度圖Figure 8 History diagram and power spectral density diagram of signal after denoising

對比圖4和圖8可知，降噪后信號噪聲減少，更好地反映了導墻振動的波形特征。從功率譜密度圖中可知，大于6 Hz的高頻白噪聲成分已被去除，低頻成分信息保存完好，保留了信號的關鍵信息。

3.3 降噪效果評價

為進一步驗證本文算法的有效性，分別用EEMD-MPE算法及CEEMDAN-相關系數(shù)算法對信號進行處理，降噪結果見圖9。

由圖9可知，CEEMDAN-MPE聯(lián)合降噪的效果優(yōu)于另外2種方法，EEMD-MPE降噪效果一般，仍存在有部分高頻噪聲，是因為EEMD分解時加入的白噪聲無法徹底消除；CEEMDAN-相關系數(shù)降噪后信號峰值點處仍存在“毛刺”，而且部分低頻成分被濾除，是由于CEEMDAN分解能夠消除加入白噪聲的影響，但相關系數(shù)判斷各IMF分量屬性不準確，使得效果不佳。在本文采用的方法中，通過計算各IMF分量的MPE值，能更好地判斷出含噪聲的IMF分量及純凈IMF分量。

圖9 降噪效果對比圖Figure 9 Comparison diagram of denoising effect

由于在實際工程中，實測信號中有效信號和噪聲信號的功率均是未知的，無法采用信噪比判斷降噪效果。因此本文引入降噪誤差比(dnSNR)判斷降噪質(zhì)量[18]，其計算公式如下：

dnSNR=10lg(Ps/Pg)。

(16)

式中：Ps為含噪信號的功率；Pg為去除噪聲信號的功率，dnSNR值越小，表明降噪效果越顯著。各方法的降噪誤差比計算結果見表2。

由表2可知，CEEMDAN-MPE算法的降噪誤差比最小，表明該算法降噪效果最好。

表2 降噪誤差比計算結果Table 2 Calculation results of denoising error ratio

4 結論

使用CEEMDAN和MPE聯(lián)合降噪的方法，對導墻振動信號進行CEEMDAN分解，利用MPE方法計算各IMF分量隨機程度，結合小波變換對含噪聲的IMF進行降噪，并與余下純凈IMF分量進行重構，完成降噪。通過該方法對實測信號進行分析，結果表明該方法在有效降噪的同時，更好地保留了結構特征信息。比較于EEMD-MPE方法及CEEMDAN-相關系數(shù)方法，該方法更適用于導墻振動信號降噪，提高了導墻振動的監(jiān)測精度，為實現(xiàn)泄流導墻健康運行及故障診斷提供依據(jù)。