張 璨，文福拴，王建軍，劉煥磊

(1.浙江大學 電氣工程學院，浙江 杭州310027;2.山西省電力公司 長治供電分公司，山西 長治046011)

0 引 言

中國臺灣中央大學的黃鍔教授(Norden E.Huang)在美國航空暨太空總署 (The National Aeronautics and Space Administration，NASA)工作期間提出了一種新的時間序列信號分析方法，即著名的Hilbert-黃變換(Hilbert-Huang Transform，HHT)［1，2］。HHT是對以傅里葉變換為基礎(chǔ)的線性和穩(wěn)態(tài)譜分析的一個重大突破，可以對非線性非平穩(wěn)信號進行分析［3］。HHT 包含兩個過程:a.通過經(jīng)驗?zāi)B(tài)分解 (empirical mode decomposition，EMD)將信號分解成有限個本征模態(tài)函數(shù)(intrinsic mode function，IMF)分量;b.對IMF進行Hilbert 變換得到瞬時頻率和瞬時幅值。EMD是HHT 的核心過程，因為只有對經(jīng)過EMD 之后的IMF 進行Hilbert 變換，得到的瞬時頻率和瞬時幅值才具有物理意義［4］。然而，在實際應(yīng)用時EMD 存在模態(tài)混疊問題，即當信號中的2 個組成分量的頻率在2 倍頻內(nèi)時，EMD 無法將其分解開［5，6］。為提高EMD 的頻率分辨率，解決2 倍頻內(nèi)不可分的問題，有些學者已經(jīng)做了一些研究工作。例如，文獻［7］ 提出了掩蔽信號EMD 法，通過采用比最高頻率分量頻率更高的掩蔽信號將該高頻分量分離出來，這在一定程度上解決了常規(guī)EMD 方法無法分離2 倍頻分量的問題，但并未給出如何系統(tǒng)地選取掩蔽信號方法;文獻［8］提出了利用快速傅里葉變換構(gòu)造掩蔽信號的方法;文獻［9］ 提出了基于快速傅里葉變換和信號能量的掩蔽信號構(gòu)造方法，較好的解決了文獻［7］ 中的方法不能有效分解低頻信號的問題;文獻［10］ 提出了頻率外差EMD 方法，通過頻率外差將高低頻信號翻轉(zhuǎn)，從而擴大了高頻和低頻信號間的差距，進而可采用常規(guī)EMD 方法分解，但由于EMD 分解的誤差和邊界效應(yīng)會導(dǎo)致在低頻段產(chǎn)生偽模態(tài)分量，經(jīng)頻率外差后這些偽模態(tài)分量將翻轉(zhuǎn)為高頻分量，與真實高頻分量混在一起而難以分辨。文獻［11，12］ 通過改進外差頻率選擇方法使得信號中各組成分量的頻率只發(fā)生偏移而不翻轉(zhuǎn)，這就避免了偽模態(tài)分量翻轉(zhuǎn)為高頻分量的問題。

針對EMD 的模態(tài)混疊問題，現(xiàn)有的研究工作大多集中于改進掩蔽信號EMD 方法和頻率外差EMD 方法，而很少有采用濾波方法來解決這一問題;對于是否發(fā)生模態(tài)混疊通常只是基于EMD 分解得到的各IMF 圖來直觀判斷，缺乏量化指標。

在上述背景下，本文將HHT 應(yīng)用于非線性信號分析之中，并研究了EMD 模態(tài)混疊問題的解決方法。首先，通過對非線性信號采用HHT 瞬時頻率圖分析，判斷是否發(fā)生模態(tài)混疊。之后，對于發(fā)生模態(tài)混疊的部分利用傅里葉變換確定發(fā)生模態(tài)混疊的信號頻率，并將高通低通濾波與HHT 結(jié)合，較好的解決了EMD 無法對2 倍頻內(nèi)信號分解的問題。最后，用算例說明了所提方法的有效性。

1 Hilbert-黃變換

1.1 經(jīng)驗?zāi)B(tài)分解 (EMD)

EMD 將信號分解為具有不同尺度特征的本征模態(tài)函數(shù)(IMF)分量，這些IMF 分量滿足以下兩個條件［7］:

(1)極值點個數(shù)和過零點個數(shù)相等或者最多相差1;

(2)在任意時刻，由極大值和極小值構(gòu)成的上下包絡(luò)線的平均值為0。

對于任意信號S(t)，將其分解為各本征模態(tài)函數(shù)分量的步驟如下［1，3，7］:

(1)對所有局部最大值點做3 次樣條插值，形成上包絡(luò)線;

(2)對所有局部最小值點做3 次樣條插值，形成下包絡(luò)線;

(3)計算上下包絡(luò)線的平均值曲線M1(t)，則信號S(t)與M1(t)之差即為P1(t):

(4)如果P1(t)滿足上述IMF 的兩個條件，則其為第一個IMF，否則將其作為原始信號重復(fù)步驟(1)到步驟(3)，得到P11(t):

式中:M11(t)為P1(t)的上下包絡(luò)線均值。

(5)重復(fù)篩選直到第k 次篩選時由式(3)得到的P1k(t)滿足IMF 的兩個條件:

在實際計算時，可利用由式(4)求得的VD值判斷每次篩選結(jié)果是否為IMF:

式中:r 為信號采樣點數(shù);VD的值通常在0.2 ～0.3之間。

(6)令C1(t)= P1k(t)，則C1(t)即為第一個IMF，其包含了原信號中周期最短的分量。將C1(t)從S(t)中分離出來:

(7)將R1(t)作為原始信號重復(fù)以上步驟n次，可獲得信號S(t)的n 個IMF:

(8)當Rn為單調(diào)函數(shù)而不能從中再分解出IMF 時，分解過程結(jié)束。最后可以得到:

1.2 Hilbert 變換

對經(jīng)過EMD 得到的各IMF 進行Hilbert 變換［1，3，7］。定義C(t)的Hilbert 形式為

構(gòu)造解析信號

幅值函數(shù)為

瞬時相位為

瞬時頻率為

2 采用HHT 的非線性信號分析

2.1 EMD 模態(tài)混疊評判指標

采用HHT 分析非線性信號，若不存在模態(tài)混疊則可直接得到各次組成分量的信息;否則，就需要采取進一步的措施。對于是否發(fā)生模態(tài)混疊，目前常用的方法是通過EMD 分解得到的各IMF 圖直觀判斷，缺乏量化指標。有鑒于此，這里提出一種通過HHT 瞬時頻率判斷是否發(fā)生模態(tài)混疊的指標。定義頻率落差為

式中:fl和fr-m 分別為HHT 瞬時頻率圖中第l 個和第r-m 個頻率點。由于EMD 存在邊界效應(yīng)，所以在式 (13)中求取最大值和最小值時去除了邊界的部分點。

用dset表示閾值，則可通過下式判斷是否發(fā)生模態(tài)混疊:

dset由專家根據(jù)所分析的具體信號以及對誤差的容忍情況來確定。當某個本征模態(tài)函數(shù)分量的瞬時頻率不滿足式 (14)時，則判定存在模態(tài)混疊，需采取措施消除。

2.2 基于高通低通濾波的EMD

可以首先采用傅里葉變換確定信號中的組成分量，之后應(yīng)用高通低通濾波將位于2 倍頻內(nèi)的高頻信號和低頻信號分別分離出來，最后采用EMD 方法對分離開的信號進行分解。

假設(shè)對發(fā)生模態(tài)混疊的信號X(t)采用傅里葉變換后確定的各組成分量的頻率為fa，fb和fc(fa＜fb＜fc)，且fc/fb＜2，fb/fa＞2。由于fb和fc的頻率位于2 倍頻內(nèi)，所以采用EMD 時會發(fā)生模態(tài)混疊。采用高通低通濾波可得到如下2 個信號:

式中:FH(·)和FL(·)分別為高通濾波變換和低通濾波變換。

頻率為fc的分量位于XH(t)中，頻率為fa和fb的分量位于XL(t)中。對XH(t)和XL(t)采用EMD 方法就不會發(fā)生模態(tài)混疊。

2.3 非線性信號分析步驟

采用HHT 分析非線性信號的步驟如下:

(1)對所要分析的信號采用HHT 得到各模態(tài)分量的瞬時頻率圖;

(2)計算各模態(tài)分量的頻率落差，找到開始發(fā)生模態(tài)混疊的模態(tài)函數(shù)分量;

(3)從原始信號中去除發(fā)生模態(tài)混疊之前的各個分量，得到分析信號;

(4)對分析信號進行傅里葉變換，并對位于2 倍頻內(nèi)的分量分別進行高通和低通濾波;

(5)對濾波后的信號采用HHT 得到各模態(tài)分量圖及其瞬時頻率圖和瞬時幅值圖;

(6)將步驟(5)中的分量與步驟(2)中未發(fā)生模態(tài)混疊的分量整合，得到完整的信號各模態(tài)函數(shù)分量。

采用HHT 分析非線性信號的流程如圖1所示。

3 算例分析

圖1 采用HHT 分析非線性信號的流程Fig.1 The flowchart of the HHT based nonlinear signal analysis

為驗證所提方法的有效性，本節(jié)給出算例結(jié)果。假設(shè)信號 S (t)= 10 sin (2π50t)+ 5 sin(2π100t)+4 sin (2π150t)+6 sin(2π350t)，dset=50 Hz。首先利用常規(guī)EMD 方法對原始信號s(t)進行分解，得到各IMF 的結(jié)果如圖2 所示，IMF 的瞬時頻率如圖3 所示。

從圖2 的EMD 分解結(jié)果中可直觀看出IMF2發(fā)生了模態(tài)混疊。由圖3 可得到前2 個IMF 的頻率落差，即 df1= 28 Hz，df2= 84 Hz;由式(14)可判斷出從第2 個IMF 開始出現(xiàn)模態(tài)混疊。從原始信號S(t)中除去IMF1 后可得到分析信號X(t)。對X(t)進行傅里葉變換，可得到圖4所示的結(jié)果。

從圖4 可看出X(t)中包含了150 Hz，100 Hz和50 Hz 這3 個分量，150 Hz 與100 Hz 以及100 Hz與50 Hz 位于2 倍頻內(nèi)，這就解釋了圖2中IMF 發(fā)生模態(tài)混疊的原因。對X(t)經(jīng)過高通低通濾波后再采用EMD 方法分解，并與圖2 中未發(fā)生模態(tài)混疊的分量整合，結(jié)果如圖5 所示。各IMF 的瞬時頻率和瞬時幅值分別如圖6 和圖7所示，端點處的異常是由濾波器的延遲造成的。從圖5 ～圖7 可看出采用所提方法可有效解決模態(tài)混疊問題，最后所得到的結(jié)果與原始信號S(t)一致。

圖2 采用常規(guī)EMD 方法對信號的分解結(jié)果Fig.2 Signal decomposition results by the conventional EMD method

圖3 各IMF 的瞬時頻率Fig.3 The instantaneous frequency of each IMF

圖4 X(t)的傅里葉譜分析Fig.4 Fourier spectrum analysis results of X(t)

4 結(jié) 論

在對Hilbert －黃變換的基本原理和存在問題進行分析的基礎(chǔ)上，提出了基于HHT 瞬時頻率來判斷是否發(fā)生模態(tài)混疊的指標。之后，通過順序采用傅里葉變換和高通低通濾波這兩種方法，來解決常規(guī)經(jīng)驗?zāi)B(tài)分析方法無法分解2 倍頻內(nèi)信號的問題。將濾波后再采用EMD 獲得的各模態(tài)函數(shù)分量與之前未發(fā)生模態(tài)混疊的分量整合，獲得了原始信號完整的各模態(tài)函數(shù)分量及相應(yīng)的頻率和幅值。最后，用算例說明了所提方法的正確性。

圖5 濾波后EMD 分解結(jié)果Fig.5 The signal decomposition results by employing the EMD method after filtering

圖6 濾波后應(yīng)用HHT 得到的各IMF 的瞬時頻率Fig.6 The instantaneous frequency of the IMF by applying the HHT after filtering

圖7 濾波后應(yīng)用HHT 得到的各IMF 的瞬時幅值Fig.7 The instantaneous amplitude of the IMF by applying the HHT after filtering

［1］Huang N E，Shen Z，Long S R，et al.The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis［J］.Proceedings of the Royal Society of London Series A，1998，45(4):903 －995.

［2］Huang N E，Shen Z，Long S R.A new view of non-linear water waves:the Hilbert spectrum［J］.Annual Review of Fluid Mechanics，1999，31(1):417 －457.

［3］蘇玉香，劉志剛，李科亮，等.Hilbert-Huang 變換在電氣化鐵路諧波檢測中的應(yīng)用［J］.電網(wǎng)技術(shù)，2008，32(18):30 －35.

［4］胡維平，莫家玲，杜明輝.經(jīng)驗?zāi)B(tài)分解中的頻域分辨率及其改進方法［J］.華南理工大學學報(自然科學版)，2007，35(5):15 －19.

［5］李成鑫，劉俊勇，楊嘉湜，等.移頻經(jīng)驗?zāi)B(tài)分解在低頻振蕩參數(shù)提取中的應(yīng)用［J］.電力系統(tǒng)自動化，2011，35(20):1 －6.

［6］楊德昌，REHTANZ C，李 勇，等.基于改進希爾伯特－黃變換算法的電力系統(tǒng)低頻振蕩分析［J］.中國電機工程學報，2011，31(10):102 －108.

［7］Deering R，Kaiser J F.The use of a masking signal to improve empirical mode decomposition［C］// Proceedings of IEEE International Conference on Acoustic，Speech，and Signal Processing.NC，USA:Durham，2005:485 －488.

［8］Senroy N，Suryanarayanan S，Ribeiro P F.An improved Hilbert-Huang method for analysis of time-varying waveforms in power quality ［J］.IEEE Trans on Power Systems，2007，22(4):1843 －1850.

［9］Laila D S，Messina A R，Pal B C.A refined Hilbert-Huang transform with applications to interarea oscillation monitoring［J］.IEEE Trans on Power Systems，2009，24(2):610 －620.

［10］Senroy N，Suryanarayanan S.Two techniques to enhance empirical mode decomposition for power quality applications［C］// Proceedings of IEEE Power Engineering Society General Meeting，June 24 －28，2007，Tampa，F(xiàn)L，USA:Tampa，2007:1 －6.

［11］胡紅英，殷福亮.頻率外差經(jīng)驗?zāi)Ｊ椒纸飧倪M算法［J］.農(nóng)業(yè)機械學報，2010，41(10):209 －213.

［12］李成鑫，劉俊勇，姚良忠，等.基于改進頻移經(jīng)驗?zāi)B(tài)分解的低頻振蕩參數(shù)提取［J］.電力系統(tǒng)自動化，2012，36(15):8 －13.