EMD 算法在檢測電力故障信號中的應(yīng)用*

楊 彥 ，王 浩 ，曹德健，趙 力*

(1.鹽城紡織職業(yè)技術(shù)學(xué)院，江蘇 鹽城 650500;2.東南大學(xué)信息科學(xué)與工程學(xué)院，南京 210096)

由于電力系統(tǒng)的復(fù)雜性，在遠(yuǎn)程傳輸以及運(yùn)載中會產(chǎn)生的不平穩(wěn)的故障信號，這也就使得信號突變的檢測在電力系統(tǒng)中占有重要的地位。其常用的傳統(tǒng)數(shù)學(xué)工具為短時傅里葉變換，短時傅里葉變換將信號劃分為許多小的時間間隔，用傅里葉變化分析每個時間間隔，一邊確定該時間間隔的頻率。然而對于短時傅里葉變換，時間—頻率窗口的寬度對于所有頻率的譜具有不變特性，這一點(diǎn)不適合于非平穩(wěn)信號的高頻與低頻部分的特性分析。后來很多工程應(yīng)用者嘗試將小波變換理論引入電力系統(tǒng)故障分析中，因為小波變換能夠在時域和頻域同時得到比較高的分辨率，然而對信號進(jìn)行小波變換必須預(yù)先確定小波的基函數(shù)，且信號分解的效果取決于小波基函數(shù)的選擇，并不能夠保證最優(yōu)的小波分解效果［1］。1998年，美國國家宇航局(NASA)的Norden和E.Huang 等人在對瞬時頻率的概念進(jìn)行深入研究的基礎(chǔ)上，創(chuàng)立了希爾伯特－黃變換，該方法被認(rèn)為是近年來對已傅里葉變換為基礎(chǔ)的線性和穩(wěn)態(tài)譜分析的一個重大突破［1］。它是一種基于經(jīng)驗?zāi)Ｊ椒纸釫MD(Empirical Mode Decomposition)的時頻分析處理方法。它將非平穩(wěn)信號按不同尺度的波動或趨勢逐級分解成若干個本征模式分量IMF(Intrinsic Mode Function)，然后對IMF 進(jìn)行Hilbert 變換，得到信號的時頻譜圖具有良好的時頻分辨率，從而具有良好的局部視頻特性表現(xiàn)能力，因此能夠更加準(zhǔn)確地反映原信號的特性，可以用來檢測突變信號。該方法吸取了小波變換的多分辨率的優(yōu)勢，同時克服了小波變換中選取小波基的困難。能夠更加準(zhǔn)確地檢測出突變信號，從而正確定位出信號的故障點(diǎn)。

1 電力系統(tǒng)中的故障信號

電力系統(tǒng)的不斷發(fā)展，使得輸送電壓和容量的逐步提高，有輸電線故障造成的損失也越來越大。電力系統(tǒng)主要有短路故障、斷路故障、復(fù)雜故障。用Simulink 中的Power system 仿真工具包來進(jìn)行故障信號仿真，通過仿真可以得到如圖1 幾種故障的A相電壓信號波形圖。

圖1 幾種故障的電壓波形圖

2 EMD 用于故障信號處理的基本原理

經(jīng)驗?zāi)Ｊ椒纸馐?998年Huang 提出來的一種全新的信號處理方法［2－3］。它通過對信號的篩選將信號分解成不同頻率的IMF，IMF 具有以下特點(diǎn):第1 IMF 的極值(包含極大值和極小值)數(shù)與過零點(diǎn)的數(shù)目相等或最多相差一個;第2 在任意時刻，IMF 的上包絡(luò)線與下包絡(luò)線的均值必須是零。具體的分解過程如下［4］:

第1 步:首先找到待分析信號x(t)所有的極大值點(diǎn)并將其使用3 次樣條函數(shù)插值(在下面做了描述)成為原數(shù)據(jù)序列的上包絡(luò)線;找到x(t)所有的極小值點(diǎn)將其使用3 次樣條函數(shù)插值成為原數(shù)據(jù)序列的下包絡(luò)線，然后求出上包絡(luò)線和下包絡(luò)線的平均包絡(luò)線m1(t);

第2 步:將原數(shù)據(jù)序列減去平均包絡(luò)后即可得到一個去掉低頻的新數(shù)據(jù)序列h1(t)=x(t)－m1(t);然后判斷h1(t)是否滿足IMF 的條件，若不滿足，將h1(t)看作新的x(t)，重復(fù)上述處理過程，直到h1(t)滿足IMF 條件時，記c1(t)=h1(t)，視為IMF1;

第3 步:將r(t)=x(t)－h(huán)1(t)看作新的x(t)，重復(fù)以上(1)和(2)的步驟，即可以依次得到IMF2，IMF3，…，直到cn(t)或r(t)滿足給定的終止條件時篩選結(jié)束。

按照上面的分解步驟，原始數(shù)據(jù)序列x(t)可以由這些IMF 分量及一個均值或者趨勢項r(t)來表示:

其中每一個IMF 分量都代表一個特征尺度的數(shù)據(jù)序列，所以可以說EMD 過程實際上是將原數(shù)據(jù)序列分解成為不同特征波動的序列。

其中Re 表示取實部，ai(t)表示每一個本征模函數(shù)經(jīng)希爾伯特變換后的瞬時振幅，fi(t)表示每一個本征模函數(shù)經(jīng)希爾伯特變換后的瞬時頻率。若以時間、瞬時頻率為自變量，幅度以等高線的形式表示，則可以將三者之間的關(guān)系表示為一個三維時頻譜，又稱為Hilbert 幅度譜，簡稱為Hilbert 譜，記為H(f，t)，則有

下面進(jìn)一步來分析Hilbert 譜，定義邊緣譜為h(f)，有

其中T 表示信號序列的時間長度。根據(jù)邊緣譜的表達(dá)式，可以看出，邊緣譜是對三維時頻譜在時間上進(jìn)行積分，便形成了只有頻率和幅值的二維的關(guān)系圖，邊緣譜表征了整組數(shù)據(jù)在每個頻率點(diǎn)的積累的幅度值分布，而傅里葉頻譜在每一點(diǎn)的幅度值表征了信號里有一個含有此頻率的三角函數(shù)的組分，傅里葉頻譜和邊緣譜是不同意義的兩種頻譜與幅值關(guān)系表達(dá)式。

從上述理論，我們可以看出，EMD 算法是具有自適應(yīng)性的，主要表現(xiàn)在:生成的基函數(shù)具有自適應(yīng)性，IMF 的頻率分辨率也是具有自適應(yīng)性的，且同時具有濾波特性［5］。所以可以把EMD 分解看做是一組具有自適應(yīng)特性的帶通濾波器，其截至頻率和帶寬隨被分解信號的不同而不同［6－7］。同小波變換相比，EMD 是有自適應(yīng)的，而小波分析則沒有，因為小波的分解尺度被選定，分解得到的時域波形就是固定頻率段的，且頻率段只與分析頻率有關(guān)，與信號本身無關(guān)［8］。這就是EMD 相比于小波變換的優(yōu)點(diǎn)。

3 EMD 方法存在的缺點(diǎn)與改進(jìn)

EMD 分解的一個重要步驟是構(gòu)造信號的上下包絡(luò)線，從而得到信號的瞬時平均包絡(luò)。Huang 等提出的以信號的極大值點(diǎn)來擬合上包絡(luò)線，極小值點(diǎn)來擬合下包絡(luò)線，其中擬合算法采用的是3 次樣條插值方法。然而，我們知道，在信號的開始端和結(jié)束端往往不是極大值或者極小值點(diǎn)，因此3 次樣條曲線容易在數(shù)據(jù)的兩個端點(diǎn)出現(xiàn)發(fā)散現(xiàn)象，且EMD的過程也使得這種發(fā)散現(xiàn)象不只局限于端點(diǎn)處的小范圍，而是伴隨著每次的篩選過程逐漸往內(nèi)傳播，會引起最終結(jié)果較大的失真，這種現(xiàn)象稱為“端點(diǎn)效應(yīng)”。此外，我們知道，EMD 的過程就是對信號進(jìn)行篩選的過程，篩選對信號有兩個結(jié)果:第1，駝峰減少，使得極大值為正數(shù)或零，極小值為負(fù)數(shù)或零;第2，光滑了非對稱幅度。然而過度的篩選會使得信號在幅度上失去物理意義，且當(dāng)前的電力設(shè)備很多都不能滿足這種過度的篩選迭代，在這里將這種現(xiàn)象稱為過度篩選現(xiàn)象。

在電力系統(tǒng)的復(fù)雜環(huán)境下，信號波形可能含有多個諧波和噪聲，端點(diǎn)效應(yīng)影響程度可能會加大，所以需要采取一定的方法來改善這種端點(diǎn)效應(yīng)造成的誤差。目前抑制這種效應(yīng)主要采用兩種方法:第1，尋找更加精確的求解上下包絡(luò)的方法來替代已有的3 次樣條插值方法，比如有約束3 次樣條插值法，分段冪函數(shù)法，高次樣條插值法;第2，在信號的兩端選取合適的極大值和極小值點(diǎn)，使上下包絡(luò)線能夠完整的包絡(luò)整個信號。在這里，我們采用第2種方法來減少端點(diǎn)效應(yīng)所帶來的誤差，由于電力信號在一般情況下非故障點(diǎn)表現(xiàn)為一系列正弦信號之和，故電力信號在非故障處可以看作周期信號，本文將端點(diǎn)處的值作為極大值或極小值進(jìn)行處理，當(dāng)端點(diǎn)值為正數(shù)時，將其作為極大值，對其取反作為極小值;當(dāng)端點(diǎn)值為負(fù)數(shù)時，將其作為極小值，對其取反作為極大值。采用這種方法信號的包絡(luò)線得到了很大的改善，能夠減少端點(diǎn)效應(yīng)所帶來的誤差，能夠滿足大多數(shù)電力設(shè)備的要求。

對于過度篩選現(xiàn)象，本文降低了篩選的要求，即降低了IMF 的要求。我們知道，當(dāng)一個信號完全滿足IMF 的要求時，其上下包絡(luò)線是關(guān)于時間軸完全對稱的，包絡(luò)線均值是處處為零的。一般情況，即使在誤差范圍內(nèi)，也需要進(jìn)行很多次的篩選，在這里，Huang等人使用esd標(biāo)準(zhǔn)來降低篩選條件。esd的表達(dá)式為

其中hi－1(t)和hi(t)分別代表經(jīng)過第i－1 次和第i次篩選后的信號。一般情況下，當(dāng)esd下降到0.2或者0.3 時，結(jié)束篩選過程。然而在esd中分母h2i－1(t)可能會出現(xiàn)零的情況，一般情況下這種為零現(xiàn)象是不可避免的，因此，esd的計算可能是不準(zhǔn)確的。

從上面所述，當(dāng)信號滿足IMF 的要求時，其上下包絡(luò)線是關(guān)于時間軸完全對稱的，所以在篩選的過程中，信號的上下包絡(luò)線是關(guān)于時間軸趨向?qū)ΨQ的，那么我們可以根據(jù)包絡(luò)線的對稱程度來決定是否結(jié)束篩選。在這里，我們設(shè)第i 次的上下包絡(luò)線數(shù)據(jù)序列分別為mup(n)和mdown(n)，為了表示兩個數(shù)據(jù)變量的近似程度，我們采用了相關(guān)系數(shù)標(biāo)準(zhǔn)。根據(jù)相關(guān)系數(shù)的定義，設(shè)ri為第i 次時的上下包絡(luò)線相關(guān)系數(shù)，表達(dá)式為

在一般情況下，我們將相關(guān)系數(shù)ri設(shè)置為大于等于0.985 時結(jié)束篩選。這樣可以有效降低EMD算法的循環(huán)次數(shù)。

4 故障信號的EMD 仿真分析

圖2為三相短路時的電壓信號波形，我們在0.03 s 時設(shè)置故障，在0.065 s 時切除。

圖2 三相短路時的電壓信號波形

應(yīng)用前兩節(jié)方法對該信號進(jìn)行EMD 分解，使用Matlab 仿真，結(jié)果為圖3。

圖3 信號的EMD 分解

對IMF 做希爾伯特變換，可得到圖4。

圖4 希爾伯特變換后的時頻關(guān)系圖

通過在Matlab 中的放大觀察，我們看到信號在0.030161和0.065125 出的頻率很高，也就是故障發(fā)生和切除的時間，其誤差0.125 ms～0.161 ms，故障時間定為準(zhǔn)確，且運(yùn)行時間在許可的時間內(nèi)。

5 結(jié)論

電力系統(tǒng)對于繼電保護(hù)裝置的要求是快速，準(zhǔn)確。本文采用了EMD 方法將信號分解，根據(jù)其IMF分量的希爾伯特變換得到信號的時頻特性，從而能夠準(zhǔn)確定位故障發(fā)生的時間，在下一步的工作中，我們應(yīng)該尋求對故障信號進(jìn)行處理分類，使得繼電保護(hù)器能夠迅速的分辨出故障的類型，以便做出相應(yīng)的動作。

［1］李天云，趙妍，李楠.基于EMD 的Hilbert 變換應(yīng)用于暫態(tài)信號分析［J］.電力系統(tǒng)自動化，2005，29(4):49－52.

［2］Huang N E，Shen Z，Long S R.The Empirical Mode Decomposition and the Hilbert Spectrum for Nonlinear and Nonstationary Time Series Analysis［C］//Proc.R.Soc.London A，1998，4:903－995.

［3］Huang N E，Stever R L.A New View of Nonlinear Water Waves:the Hilbert Spectrum［C］//Annual Review of Fluid Mechanics.1999，1:417－457.

［4］Balocchi R，Menicucci D，Varanini M.Empirical Mode Decomposition to a Approach the Problem of Detecting Sources from a Reduced Number of Mixtures［C］//Cancun，Mexico:Proceedings of the 25th Annual International Conference of the IEEE EMBS，2003.

［5］Raia V K，Mohanty A R.Bearing Fault Diagnosis Using FFT of Intrinsic Mode Functions in Hilbert-Huang Transform［J］.Mechanical Systems and Signal Processing，2007，21(6):2607－2615.

［6］Wu Z H，Huang N E.Ensemble Empirical Mode Decomposition:A Nosie Assisted Data Analysis Method［R］.Calverton Center for Ocean-Land-Atmosphere Studies，2005:855－895.

［7］Gabriel R，Patrick F.One or Two Frequencies? The Empirical Mode Decomposition Answers［J］.IEEE Transactions on Signal Processing，2008，56(1):85－96.

［8］Zheng G T，McFadden P D.A Time-Frequency Distribution for Analysis of Signals with Transient Components and Its Application to Vibration Analysis Transaction of the ASME［C］//Proc.R.Soc.London A，July 1999，2:328－333.