基于峰度檢驗和小波包分解的海洋平臺脈沖噪聲處理方法研究

李 偉，黃 焱

(1.天津大學 水利工程仿真與安全國家重點實驗室，天津 300350；2.天津大學 建筑工程學院，天津 300350；3.天津大學 港口與海洋工程天津市重點實驗室，天津 300350)

海洋平臺作為海上生產和生活的基地[1]，為了確保平臺安全作業，對平臺結構進行實時健康監測，及時辨識損傷并對結構狀態進行評估與預警，具有重要的工程意義。這一方法需要在結構服役期間以傳感系統采集結構的動態響應，通過動態響應提取損傷敏感性特征進行統計分析，進而判定結構的健康狀況。考慮到腐蝕性鹽水環境會對傳感器造成不良影響，極易造成傳感器的損壞，因此，振動傳感器通常安裝在平臺甲板上[2]。

雖然將傳感器布置于平臺甲板對提高傳感器的壽命以及控制結構健康監測的成本起到了積極的作用，然而平臺上的很多設備都布置于甲板位置，包括鉆井與井口設備、油氣處理設備等，這就致使傳感器處于較為復雜的噪聲環境中，振動響應信號將不可避免的受到噪聲的影響。圖1展示了采集得到的現場加速度時程曲線，從圖中可以觀察到曲線存在很多脈沖噪聲，這將影響損傷敏感性特征的提取。由此可見，如何對脈沖噪聲進行降噪，是數據預處理的主要內容之一。

圖1 現場監測加速度時程曲線Fig.1 Acceleration time history curve of field monitoring

脈沖噪聲時域上的特征主要表現為振動響應值通常高于正常信號的幾倍甚至十倍以上，同時噪聲的頻譜具有很寬的頻帶,往往和信號頻率重疊在一起,常用的數字濾波器無法將兩者區分開來[3]。童進等[4]提出了一種時域降噪方法，該方法首先確定脈沖噪聲的位置,并用同相位點的平均值代替該點的值的方法消除脈沖噪聲，這一方法取得了一定效果，但降噪后的曲線仍存在峰值很高的脈沖噪聲。

從上述分析可以發現，僅僅依靠頻域或時域的降噪方法對于脈沖噪聲的處理并不十分理想，而小波包變換可以用于對信號進行時-頻分析[5-6]，目前在信號處理領域取得了廣泛的應用。嚴家斌等[7]提出了利用小波變換模極大值特性的改正方法，并利用實測數據驗證了該降噪算法。何志勇等[8]在小波包域內估計脈沖噪聲的能量分布并依據估計結果計算每段脈沖噪聲的濾波閥值，利用這一算法對軸承故障檢測的數據進行處理。但這一方法對脈沖噪聲的特性有一些假設條件，不適用于較為復雜的脈沖噪聲。盡管目前開展了一些基于小波包的脈沖噪聲降噪方法，但對于海洋平臺脈沖噪聲的研究較少，因此，本文將發展一套基于峰度檢測和小波包分解相結合的海洋平臺脈沖噪聲處理方法，并利用此方法對現場監測數據進行處理。

1 峰度檢驗與小波包聯合降噪方法的提出

小波包變換是小波變換的進一步改進，克服了小波變換對信號的高頻部分沒有細分的缺點，提供了一種更為精細的分析方法。這一方法利用遞歸濾波操作將輸出信號在由高到低的較寬頻帶上同步連續分解，根據頻帶的變化自適應的確定信號分辨率的取值，形成既有低頻又有高頻的能量組成[9]。小波包變換應用于信號降噪的基本思想是未受噪聲污染的信號和噪聲信號的小波包系數具有不同的小波包變換特性。幅值較大的小波包系數代表未受噪聲污染信號的信息，幅值較小的小波包系數代表噪聲信號的信息。通過合理的閾值設置可以調整較大的系數，而將幅值較小的系數設置為零，隨后再進行重構得到降噪后的信號[10]。

然而，傳統的小波包降噪方法主要針對白噪聲和寬帶噪聲，并且具有較好的降噪效果。對于脈沖噪聲，降噪效果并不理想[11]，具體表現為噪聲能量在較短的時域區間內依然比較顯著，因此，如何快速尋找這些時域區間并在區間內對小波包系數進行合理的閾值量化是小波包進行脈沖降噪的技術難點。由此目標出發，本文提出了兩點改進：峰度檢驗，即利用峰度確定脈沖噪聲出現的時域區間和基于K-S檢驗和小波包系數概率統計的閾值設置方法。

1.1 峰度檢測

由于噪聲和正常信號是相對的，本文主要針對脈沖噪聲的消除，因此將非脈沖成分暫視為正常信號。在工程上，通常利用峰度作為特征指標[12]，以表征振動信號中的脈沖信號，假設采集得到的數據為ai(i=1,2,3…，n)，峰度的計算式為

(1)

式中：k為峰度系數；n為信號的總采樣點數；ai為第i個采樣點；μ和σ分別為采集得到的數據的均值和標準差。

峰度是用于描述分布形態的陡緩程度的統計參量，正態分布的峰度為3。Kline[13]指出，當樣本數據峰度的絕對值小于10時，可以認為觀測變量基本服從正態分布，因此，本文在判斷信號是否存在脈沖噪聲時，將峰度的閾值設置為10。為了獲取信號中的脈沖噪聲出現的時間區域，需要對信號劃分時間段以便求得各個時間段的峰度。具體方法如下：首先需要將原始信號進行Q等分，每段長度為L，用bm表示各段，記為{bm}，m=1,2,…，Q。需要注意的是，若L取值太小，原始信號的變化趨勢不能較好的反映出來，會導致脈沖噪聲和其它信號的混疊；若L取值太大，則對于脈沖噪聲的時間分辨率較低。因此，L的取值應在兩者之間進行折中。由于脈沖噪聲持續時間通常較短，為了能有比較好的時間分辨率，本文將L的最大取值設置為15 s。當L取某一值時，可根據式(1)得到一系列峰度值{km}，m=1,2,…，Q。對這一系列峰度值進行統計，求得在此L值條件下的變異度Cvm，其定義為

(2)

1.2 閾值選取

小波包降噪中閾值選取十分重要，依據峰度檢驗的結果可以確定脈沖噪聲以及非脈沖信號的時間區域，因此這兩類信號在每個小波包節點上映射的時域區間的小波包系數也可以相應的確定。假設第j層第k個小波包節點處的非脈沖信號對應的小波包系數G(j,k)可以表示為G(j,k)={α1,α2，α3，α4,…,αn}，對于閾值的估計，采用K-S檢驗來估計非脈沖信號對應的小波包系數的概率分布函數，該方法可以檢驗一個經驗分布是否符合某種理論分布，本文考慮3種理論分布：①正態分布；②極值分布；③瑞利分布。確定了非脈沖信號段的小波包系數的分布規律后，閾值Th可根據99.7%的置信水平估計。文獻[14]指出，如果脈沖噪聲的小波包系數W(j,k)大于閾值Th，則將其視為脈沖系數；為了抑制脈沖，將該系數衰減為新系數W′(j,k)，由

(3)

根據式(3)的閾值修正方法，可以將某個節點脈沖噪聲段的小波包系數進行修正，其它節點的小波包系數與上述處理方式類似，這里就不再贅述。對每個節點的小波包系數進行閾值量化后，再對信號進行重構，即可得到降噪后的曲線。

由此可見，基于第一點改進，可以確定脈沖噪聲在每個小波包節點上映射的時域區間，基于第二點改進，可利用非脈沖信號時段的小波包系數的分布規律確定系數調整閾值，再對脈沖噪聲時段的小波包系數進行閾值處理，因此，本文方法可以很好的適應脈沖噪聲的時域分布，圖2為本文算法的流程圖，對這一方法的主要步驟進行了總結。下面將結合數值模擬的計算結果，對本文提出的方法進行應用。

圖2 信號中含有脈沖噪聲的降噪步驟Fig.2 De-noising procedure for time series polluted by impulse noise

2 數值模擬信號驗證

2.1 數值模擬信號的獲取

本文主要研究信號中的脈沖噪聲消除方法，目標海洋平臺是位于渤海中部海域的NB35-2平臺(見圖3)。該平臺由導管架結構作為支撐，上部組塊分為3層，平臺作業水深12.6 m，平臺總質量為9 610.4 t。

圖3 NB35-2海洋石油平臺Fig.3 The NB35-2 offshore oil platform

數值模擬信號的獲取是通過在ANSYS中對導管架平臺數值模型(見圖4)開展相應的動力分析計算得到的，在計算中主要考慮了波流載荷的影響，為模擬各種海況下平臺所受的波浪沖擊，利用隨機波浪理論模擬現實海況，對模型平臺施加隨機波浪激勵。圖5(a)展示了計算得到的1.0 m波高的加速度時程曲線，時間為10 min。在獲得原始加速度信號的基礎上，隨后利用MATLAB軟件生成在時域中起始位置和寬度隨機的脈沖噪聲，用來污染圖5(a)中純凈的加速度信號，利用這一方法得到的含有脈沖噪聲的信號，能夠較好的檢驗本方法的有效性，含脈沖噪聲的信號如圖5(b)所示。

圖4 平臺有限元模型Fig.4 Finite element model of the platform

2.2 信號降噪與分析

圖5(a)中的信號被脈沖噪聲所干擾，為了對圖5(b)中的信號進行降噪，首先應求取信號的峰度，進而得到脈沖噪聲所在的時域區間。表1展示了不同L值條件下的變異度，隨著L值增加變異度呈現出逐漸增加的變化趨勢，在L取15時，變異程度最高，因此在進行峰度檢驗時，可將每段統計時間長度設置為15 s。圖6展示了峰度值隨時間的變化，可以觀察到峰度的變化和脈沖噪聲出現的時間具有很好的對應性，由此可見，利用這一指標可以準確的檢測到脈沖噪聲出現的時域區間。

表1 變異度隨L值的變化Tab.1 Plotted results of Cvm against L

圖5 加速度時程信號Fig.5 Time history of acceleration in the numerical simulation

圖6 峰度檢驗結果Fig.6 Plotted results of kurtosis against time

峰度檢驗完成后，下一步將對信號進行小波包分解，采用具有多種小波基優良特性的db5小波函數進行，圖7展示了經3層小波包分解后的完全二叉小波包節點樹，樹中節點下的數字為節點編號。經過3層小波包分解后，即可得到各個頻段的小波包系數，如圖8所示，展示了分解后第3層各個節點的小波包系數。從分解的小波包系數可以看到，在脈沖噪聲發生的時間區域，小波包系數的絕對值較大，顯著超過了非脈沖噪聲信號的小波包系數。除此之外，脈沖噪聲不僅影響了信號低頻部分的小波包系數，同時對高頻部分的小波包系數也有顯著影響，這正說明脈沖噪聲具有較寬的頻率分布。

圖7 三層分解的完全二叉小波包樹Fig.7 Binary wavelet packet tree using 3 layer decomposition

圖8 第三層各個節點的小波包系數Fig.8 Wavelet packet coefficients of each node in the 3 layer

根據式(3)的閾值調整算法，需要首先對未受脈沖噪聲影響的信號段的小波包系數進行K-S檢驗，以(3,1)節點的小波包系數為例，小波包系數較好的服從正態分布，如圖9所示，因此可根據正態分布的統計規律設置這一節點的閾值為小波包系數的均值加減三倍標準差。其它節點的小波包系數的閾值設置方法也按這一節點的處理方法進行統計，從而得到各個節點的閾值。

圖9 (3,1)節點小波包系數正態分布概率圖Fig.9 Probability plots of wavelet packet coefficients of node (3,1) for normal distribution

對各節點脈沖噪聲段的小波包系數進行閾值處理后即可通過小波包逆變換重構回原始數據，結果如圖10所示。從圖中可以觀察到，盡管脈沖導致的顯著能量集中有了改善，但仍然可以從曲線中觀察到脈沖噪聲的跡象，這主要是由于式(3)的閾值調整方法是將幅值較大的小波包系數調整至閾值水平，這對于本文的脈沖噪聲的消除不理想，因此需要進一步改進小波包系數的調整算法。

圖10 數值模擬加速度時程降噪結果Fig.10 The results of de-noising for the time history of acceleration in the numerical simulation

2.3 改進方法與效果評價

當小波包系數超過閾值時，乘以一個加權函數q(j,k)[15]可使得其值迅速減小。因此，通過對脈沖噪聲時段內的小波包系數進行加權，是十分重要的系數調整方法。考慮到某一節點的閾值Th和這一節點脈沖噪聲時段的小波包系數W(j,k)的比值可以用來描述脈沖噪聲時段的小波包系數偏離正常系數的程度，因此加權函數q(j,k)可以表示為

(4)

新的小波包系數設定方法可以在式(3)的基礎上調整為

(5)

從式(5)的表達式可以看出，當脈沖噪聲時段的小波包系數顯著超過閾值時，經加權處理會對原始系數賦予較大的權重，而當脈沖噪聲時段的小波包系數與閾值相差不顯著時，加權處理后的小波包系數與閾值較為接近。經上述修正后，對調整后的小波包系數進行重構，得到的降噪曲線如圖11所示，此時脈沖噪聲作用的跡象已經消失，與前述的處理方法相比，降噪效果較好。

圖11 新閾值條件加速度時程降噪結果Fig.11 The results of de-noising method using the new threshold

為了定量評價本文提出的降噪算法，首先定義輸入信噪比(RSN,IN)和輸出信噪比(RSN,OUT)

(6)

(7)

除了上述參數外，原始信號與降噪信號之間的均方根誤差(root mean square error，RMSE)也是評價降噪效果的重要指標，其定義為

(8)

(9)

式中，T為信號時長。根據式(9)，選取R(τ)的最大值作為兩組信號的相關性，即

(10)

由式(6)可知原始信號的信噪比(RSN,IN)為1.38，在進行算法比較時，首先根據圖10的降噪結果計算上述比較參數，這一方法是Nongpiur提出的小波包閾值方法。除了這一方法外，中值濾波方法也是一種廣泛采用的濾除脈沖噪聲的方法，其具有計算簡單以及速度快等特點，作為另外一種比較算法與本文方法進行對比，三種方法的比較結果如表2所示。從結果可以看出，本文的方法無論從信噪比還是均方根誤差來說都是最好的，同時，原始信號和經本文建立的方法降噪后信號的相關系數為0.97，具有很高的相關性，降噪后的能量損失僅為3%，表明降噪后的信號是原始信號的高精度逼近。

表2 三種方法降噪后RSN，DRMSE和相關性的對比Tab.2 Comparison RSN,DRMSE and correlation of noise reduction in three kinds of method

3 現場監測信號

3.1 振動響應信號的獲取

現場數據的獲取主要依托“科技部十三五重大專項”的支持，利用低頻加速度傳感器獲取NB35-2平臺的振動響應，加速度傳感器為德國MMF公司的KS823B系列三軸加速度傳感器，每個傳感器安裝在防爆盒內進行保護，共在平臺甲板上布置了18個傳感器，數據采集頻率為50 Hz。圖12展示了不同測點采集得到的加速度傳感器時程曲線，由圖可知，在現實環境激勵下，平臺加速度響應時程曲線存在很多脈沖噪聲。

圖12 不同測點的加速度時程曲線Fig.12 Acceleration time history of different measuring points

以某一天的加速度數據為例，圖13展示了各個加速度傳感器一天內脈沖噪聲的出現次數，圖中的x軸表示時間，y軸表示傳感器編號，z軸表示脈沖噪聲出現次數。可以觀察到，由于不同傳感器處于平臺的不同位置，因此每個傳感器具有不同的脈沖噪聲特點。同時，脈沖噪聲出現頻率較高的時段大致在8:00—12:00、14：00—18:00，這也正對應平臺常規作業時段。由此可見，脈沖噪聲的出現和平臺上的生產作業息息相關。

圖13 一天內各個加速度傳感器脈沖噪聲出現次數Fig.13 Impulse noise for each acceleration sensor in a day

3.2 現場數據降噪結果

將本文建立的降噪算法應用于某一受到脈沖噪聲影響的現場實測數據(見圖14)，以驗證該方法的有效性和合理性。將該數據進行峰度檢驗后確定受到脈沖噪聲干擾的三個時間段分別是100～115 s，120～135 s和180～195 s。隨后對數據進行小波包分解，并利用K-S檢驗非脈沖時段的小波包系數的分布，以(3,1)節點為例，圖15～圖16分別展示了該結點的小波包系數，以及非脈沖時段小波包系數的分布特點。

圖14 受到脈沖噪聲影響的現場實測加速度數據Fig.14 Acceleration data affected by impulse noise in the field monitoring

圖15 現場數據(3,1)節點的小波包系數Fig.15 Wavelet packet coefficients of node (3,1) for field monitoring acceleration data

圖16 現場數據(3,1)節點小波包系數正態分布概率圖Fig.16 Probability plots of wavelet packet coefficients of node (3,1) for field monitoring acceleration data

如圖16所示，非脈沖時段(3,1)節點的小波包系數較好的服從正態分布，其余節點的小波包系數和這一節點的分布特點具有一致性，這里就不再一一展示。根據上文的小波包系數閾值處理方法對脈沖時段的小波包系數進行加權處理，最后再對處理過后的小波包系數進行重構即可得到消除脈沖噪聲后的加速度時程信號，降噪結果如圖17所示。從圖中可以看出，基于峰度檢驗和小波包分解的脈沖噪聲處理方法很好的消除了現場信號中的脈沖噪聲。

圖17 降噪完成后的現場加速度數據Fig.17 Field measured acceleration data after noise reduction

4 結 論

本文針對隨機海洋環境下某導管架平臺監測數據經常出現的脈沖噪聲開展了相關的降噪研究工作，考慮由于脈沖噪聲通常具有較寬的頻帶，單一的頻域或者時域處理方法具有一定的局限性，因此發展了一套基于峰度檢驗和小波包分解相結合的處理方法，經過數值模擬和現場數據的驗證，利用該方法對含有脈沖噪聲的數據處理后，信噪比有較大的提升，并且信號的相關性高，能有效消除脈沖噪聲。

本文所開展的工作主要針對信號中的脈沖噪聲，由于現場采集得到的信號是高度復雜的非線性耦合，其組成十分復雜，如何在脈沖噪聲處理后，進行其余噪聲的識別和降噪工作是下一步研究的重點，相關研究工作正在陸續開展。