馬小剛,王永平,杜 歷
(寧夏水利科學研究院,銀川 750021)
1998年,美國宇航局的華裔科學家黃鍔(N.E.Huang)博士提出在瞬時頻率概念基礎上用希爾伯特-黃變換(Hilbert-Huang Transform,HHT)處理非線性、非平穩信號,包括經驗模式分解(Empirical Mode Decomposition,EMD)和希爾伯特譜分析(Hilbert Spectrum Analysis,HSA)兩個處理過程[1–2]。由于經EMD算法分解信號后得到的固有模態函數(Intrinsic Mode Function,IMF)分量具有非線性,故該方法非常適合用于分析非線性、非平穩信號。然而,該方法采用三次樣條插值來擬合信號數據的包絡曲線,很難實現包絡曲線準確到達兩端邊界處,由此產生分解誤差,這就是所謂的“端點效應”。這種誤差一旦產生,會隨著分解層數的增加而增大,最終使信號的分解結果無意義[3]。Huang提出用特征波法來解決這一難題,已將該方法申請為美國宇航局專利,但他同時指出該方法僅對一些特定的信號波形有效[4]。隨后,其他學者對此問題展開研究,提出一些新方法[5–20],這些新方法對特定信號的處理結果令人滿意。受鏡像延拓法的啟發,本文提出一種能有效抑制端點效應的新方法,從而實現對非線性、非穩定信號的分析和處理。
通常端點效應的產生是由三次樣條擬合信號數據包絡曲線時引起,這是因為擬合包絡曲線時需要用到前后兩個相鄰的極值點序列,如此一來就會引發信號數據的極值序列在兩端邊界處的取值問題,即極值序列在邊界兩端處既有可能是極大值點,也有可能是極小值點,還有可能是非極值點。由于是在信號數據兩端邊界處取值,因此,也有學者將其稱之為“邊界效應”問題[6–7],如圖1所示。信號曲線左端邊界處的下包絡曲線越過上包絡曲線,出現嚴重失真的“端點效應”問題,而且這種誤差會隨著分解層數增加不斷增加,將會使得分解結果誤差過大。由EMD算法的篩分過程可知,第1階IMF分量的準確與否能嚴重影響到后續IMF分量的準確性。如果第1階IMF分量有較大誤差,那么隨著篩分層數的增加,誤差會不斷累積傳播,將會導致用EMD算法分解信號后的結果失去意義。

圖1 端點效應
最近一次極值鏡像延拓法,可簡稱為“最近一次極值法”。與鏡像延拓法不同,它不是利用原始信號數據波形進行鏡像延拓,而是先求出原始信號的上下包絡曲線,然后通過一次極值法直接求解上下包絡曲線在邊界兩端的極值點,以此達到包絡曲線能準確到達邊界兩端的目的,從而避免端點效應問題的產生。由于該方法是利用原始信號極值序列中離左右兩端最近的一次極值作為拐點和鏡像對稱中心來延拓包絡曲線,故將其稱之為“最近一次極值鏡像延拓法”,“一次極值”的含義則是指在求包絡曲線時需對原始信號使用一次極值求解方法才能獲得該信號的極值序列。
設離散信號

其采樣步長為Δt,X(t)有M個極大值和N個極小值,對應的序列下標(Im,In)、時間(Tm,Tn)和函數值(U,V)記為


若要利用最近一次極值法來延拓原始信號的極值序列,一般有如下幾種情況需要討論。
(1)信號序列左端邊界




令左端需要延拓的極值點數為p,則
Tm(-p+1)、Tm(-p+2)、…、Tm(0)為左端延拓的極大值點所對應的時間,xm(-p+1)、xm(-p+2)、…、xm(0)為左端延拓的極大值點所對應的函數值;Tn(-p+1)、Tn(-p+2)、…、Tn(0)為左端延拓的極小值點所對應的時間,xn(-p+1)、xn(-p+2)、…、xn(0)為左端延拓的極小值點所對應的函數值。
(2)信號序列右端邊界




令右端需要延拓的極值點數為q,則為右端延拓的極大值點所對應的時間,為右端延拓的極大值點所對應的函數值;為右端延拓的極小值點所對應的時間為右端延拓的極小值點所對應的函數值。
綜上所述,最近一次極值法的主要思想和基本原理為:充分利用原始信號局部波形在其邊界附近拐點處的自然走向插補延展極值序列在邊界端點處的取值,旨在確保極值序列能準確到達邊界兩端。由于采用的是真實存在的波形拐點來順延極值序列,故延拓部分的波形趨勢能嚴格遵循信號序列的發展趨勢,從而避免端點效應的產生。
原始信號的極值序列到達兩端邊界處的延拓應當包括對其上下包絡曲線左右兩端共計四個極值點的求解,下面以一個非線性、周期性信號x(t)為例,對利用最近一次極值法求解其左端邊界值的具體過程進行介紹。
設原始信號為x(t):
(1)先求出原始信號x(t)的所有極小值點如圖2中的黑色“+”點所示,由于這些極小值點是第一次求出的,為了以示區別稱其為一次極小值點
(6)求出水平距離(橫坐標值)tm后,以中心對稱軸為鏡像軸找出水平距離(橫坐標值)tm在中心對稱軸右側的鏡像水平距離(橫坐標軸)tp,即tp=2tm;
(7)過鏡像水平距離(橫坐標值)tp做垂直于橫坐標軸并平行于縱坐標軸的垂線Lυ,且該垂線與已求出的下包絡曲線相交于一點P′(tp,xp),該點即為下包絡曲線左端點P(t1,xp)關于中心對稱軸的鏡像點P′(tp,xp);
(8)再過剛求出的鏡像點P′(tp,xp)做垂直于縱坐標軸并平行于橫坐標軸的水平線Lh,該水平線與縱坐標軸交于一點P(t1,xp),該點即為所求的信號下包絡曲線延拓至左端點處的端點值。
需要進一步指出的是:由于信號x(t)具有周期性,依據該信號的周期重復性,可知圖2中的Pr應該就是該信號x(t)下包絡曲線延拓至左邊界處的精確解,而從圖2上也能直觀地發現該精確解Pr與用最近一次極值法求出的端點值P恰好是重合的,這也就充分證明了利用最近一次極值法所求信號x(t)的包絡曲線在邊界處端點值的有效性和準確性。

圖2 根據最近一次極值法求端點值
對一仿真信號x(t)先用最近一次極值法進行延拓處理,然后再進行EMD分解,并與其它延拓方法的EMD分解結果相比較[8–13],以此來探討最近一次極值法能否有效抑制端點效應的產生。
設仿真信號x(t)由頻率為120 Hz的余弦信號和基頻為20 Hz、調頻為10 Hz的調幅調頻信號疊加而成,具體形式如式(39)所示。

其中,采樣頻率Fs=1 000Hz。原始信號及其組成分量如圖3所示。
下面依次使用偶延拓EMD算法(Even Extension of EMD,EE-EMD)、鏡像延拓EMD算法(Mirror Extension of EMD,ME-EMD)、周期延拓EMD算法(Periodic Extension of EMD,PE-EMD)、最近一次極值EMD算法(The Nearest Simple Extreme Mirror Extension of EMD,NSEME-EMD)分解仿真信號x(t),分解時間依次為0.172 s、0.219 s、0.344 s、0.547 s,結果見圖4和圖5。
圖4所示的IMF1對應于余弦信號cos(120 πt)分量,IMF2對應于調幅調頻信號(1+0.2sin(10 πt))·cos(40 πt+0.5sin(20 πt))分量,殘余分量Res可視為EMD算法篩分過程中所產生的分解誤差。
圖5所示的error1表示分解出的IMF1分量與余弦信號cos(120 πt)分量之間的誤差,error2表示分解出的IMF2分量與余弦信號分量之間的誤差。由圖4和圖5可知:EE-EMD算法對端點效應的抑制效果較差,ME-EMD算法對端點效應的抑制效果稍好,PE-EMD算法對端點效應的抑制效果較好,NSEME-EMD算法對端點效應的抑制效果最好。
為了準確評價上述四種EMD算法的分解精度,此處采用基于能量的評價指標θi進行計算分析[5,9,16,18]

圖3 原始信號及其組成分量

其中:θi為第i個IMF的評價指標;RMSi和RMSre為分解得到的第i個IMF和相應真實值的有效值。
RMS按下式計算

其中m為信號的采樣點數;s(j)為信號序列。
通過式(40)可知計算得出的θi≥0 ,如果θi=0,則說明端點效應對EMD分解無影響;如果θi>0,則說明端點效應對EMD分解有影響,并且θi值越大影響也越大。
四種EMD算法分解仿真信號的性能評價指標計算結果見表1。

表1 性能評價指標計算結果表/(%)
由表1可知:四種EMD算法的性能指標計算結果從大到小按序排為EE-EMD>ME-EMD>PE-EMD>NSEME-EMD,換句話說,EE-EMD算法在分解復雜信號x(t)的過程中造成的能量泄露最大,對信號序列x(t)的分解精度最低,ME-EMD算法和PEEMD算法在分解復雜信號x(t)的過程中造成的能量泄露較小,對信號序列x(t)的分解精度較低,而本文提出的NSEME-EMD算法在分解復雜信號x(t)的過程中造成的能量泄露最小,對信號序列x(t)的分解精度最高。

圖4 基于四種EMD算法的分解結果

圖5 根據四種EMD算法得到的IMF分量的誤差
通過仿真實例驗證了最近一次極值法能有效抑制端點效應,故將其應用到實際泵站廠房結構振動監測信號的分析中去。對寧夏新海水務有限公司同心分公司三泵站的廠房結構進行振動監測,從現場采集到的振動信號中選取一條采樣頻率為1 000 Hz、時長為2 s的豎向加速度響應信號的時程數據序列進行分析,如圖6所示。
分別采用Huang提出的傳統EMD方法和本文提出的最近一次極值法對該加速度響應信號進行分解,具體分解結果見圖7。
經計算,采用傳統EMD方法分解得到前5階各個IMF分量相對應的性能評價指標θ1~θ5分別為2.576%、3.864%、6.182%、10.510%、18.918%,而采用最近一次極值法分解得到的前5階各個IMF分量相對應的性能評價指標θ1~θ5分別為1.989%、2.486%、3.232%、4.363%、6.109%,由此可見相比于傳統EMD方法,最近一次極值法對振動信號的分解精度高,能量泄漏少,對端點效應的抑制效果好。

圖6 實測豎向加速度響應信號時程數據序列

圖7 實測豎向加速度響應信號分解結果
(1)通過最近一次極值法得到的邊界端點值能較好匹配原始信號在邊界附近的自然走向,實現原始信號在邊界端點處的光滑過渡,有效降低原始信號極值點在邊界處的波動,從原理上來說不失為一種解決端點效應的好方法。
(2)通過仿真算例和應用實例對最近一次極值法的分解性能進行對比分析,結果表明該方法能有效解決端點效應問題,同時具有簡單易行、適應性強、通用性好等特點。