黃傳金，曹文思，陳鐵軍，邱道尹

（1.中州大學(xué) 工程技術(shù)學(xué)院，河南 鄭州 450044；2.華北水利水電大學(xué) 電力學(xué)院，河南 鄭州 450045；3.鄭州大學(xué) 電氣工程學(xué)院，河南 鄭州 450001）

0 引言

近年來，電力系統(tǒng)中的電能質(zhì)量擾動(dòng)信號(hào)檢測(cè)成為研究的熱點(diǎn)［1-7］，尤其是諧波擾動(dòng)信號(hào)和間諧波信號(hào)檢測(cè)備受關(guān)注［8-14］。快速傅里葉變換（FFT）不能處理非線性、非平穩(wěn)的信號(hào)，在處理間諧波時(shí)還存在頻譜泄漏和柵欄現(xiàn)象等缺點(diǎn)［15］。采用小波理論分析非線性、非平穩(wěn)的信號(hào)有很大的局限性，一旦分解層數(shù)確定，其頻率分辨率也就恒定，缺乏自適應(yīng)，同時(shí)還受Heisenberg測(cè)不準(zhǔn)原理的限制［16］?；贖ilbert-Huang變換（HHT）的電能質(zhì)量檢測(cè)方法雖取得了較好的效果，但在經(jīng)驗(yàn)?zāi)B(tài)分解中采用三次樣條插值擬合包絡(luò)信號(hào)時(shí)易出現(xiàn)過包絡(luò)、欠包絡(luò)現(xiàn)象；HHT中過多的“篩分”次數(shù)導(dǎo)致端點(diǎn)效應(yīng)污染整個(gè)數(shù)據(jù)段的程度較大，而且基于HHT的時(shí)頻分析方法得到的瞬時(shí)頻率常出現(xiàn)負(fù)值，這是難以解釋的物理現(xiàn)象［17］。

2005年，Smith在前人研究的基礎(chǔ)上提出了一種新的非線性、非平穩(wěn)信號(hào)的分解方法——局部均值分解 LMD（Local Mean Decomposition）［18］。LMD 將復(fù)雜信號(hào)分解為若干乘積函數(shù)PF（Product Function）之和，每個(gè)PF分量由一個(gè)純調(diào)頻函數(shù)和一個(gè)包絡(luò)函數(shù)之積組成，PF分量的頻率可由純調(diào)頻函數(shù)直接求取，包絡(luò)函數(shù)是PF分量的瞬時(shí)幅值。

LMD方法已成功應(yīng)用于腦電信號(hào)的檢測(cè)、瞬時(shí)頻率信號(hào)的提取和機(jī)械故障診斷［17-20］。由于LMD避免了HHT中過包絡(luò)、欠包絡(luò)的現(xiàn)象，具有端點(diǎn)效應(yīng)小，所求頻率是正的、連續(xù)的、有物理意義的時(shí)變頻率等優(yōu)點(diǎn)，基于LMD的時(shí)頻分析方法有望在電力系統(tǒng)故障檢測(cè)與信號(hào)分析領(lǐng)域獲得應(yīng)用。

1 LMD

LMD本質(zhì)上根據(jù)信號(hào)固有的包絡(luò)特征自適應(yīng)地將一個(gè)非線性、非平穩(wěn)信號(hào)按頻率遞減的順序逐級(jí)分離，得到若干個(gè)具有一定物理意義的PF分量，每個(gè)PF分量由純包絡(luò)信號(hào)和調(diào)頻信號(hào)之積組成。以信號(hào) x（t）為例，給出其分解步驟［19-20］。

a.找出信號(hào)x（t）的局部極值點(diǎn) ni，任意 2個(gè)相鄰的局部極值點(diǎn)平均值記為mi，則有：

將式（1）中所有相鄰的局部均值點(diǎn)mi和mi+1用折線連接，然后用滑動(dòng)平均法對(duì)其進(jìn)行平滑處理，得到局部均值函數(shù) m11（t）。

b.求出包絡(luò)估計(jì)值ai：

將式（2）中所有相鄰的局部均值點(diǎn)ai和ai+1用折線連接，然后用滑動(dòng)平均法對(duì)其進(jìn)行平滑處理，得到包絡(luò)估計(jì)函數(shù) a11（t）。

c.將局部均值函數(shù) m11（t）從原始信號(hào) x（t）中分離出來，得到信號(hào) h11（t）：

d.用 h11（t）除以包絡(luò)估計(jì)函數(shù) a11（t），得到調(diào)頻信號(hào) s11（t）：

對(duì) s11（t）重復(fù)上述步驟，得到 s11（t）的包絡(luò)估計(jì)函數(shù) a12（t）。若 a12（t）不等于 1，說明 s11（t）不是一個(gè)純調(diào)頻信號(hào)，需要重復(fù)上述迭代過程n次，直至s1n（t）為一個(gè)純調(diào)頻信號(hào)，即 s1n（t）的包絡(luò)估計(jì)函數(shù)a1（n+1）（t）=1，所以有：

迭代終止的條件為：

實(shí)際應(yīng)用中，為避免過多分解次數(shù)，設(shè)一個(gè)變動(dòng)量 Δ，可令 1-Δ≤a1n（t）≤1+Δ 時(shí)，迭代終止。

e.將迭代過程中產(chǎn)生的所有包絡(luò)估計(jì)函數(shù)做乘積，得到包絡(luò)信號(hào) a1（t）：

f.將式（8）得到的包絡(luò)信號(hào) a1（t）與純調(diào)頻信號(hào)s1n（t）做乘積，得到原始信號(hào) x（t）的第 1 個(gè) PF 分量：

g.將 PF1（t）分量從原始信號(hào) x（t）分離出來，得到信號(hào) u1（t），將 u1（t）作為新的數(shù)據(jù)重復(fù)以上步驟，循環(huán)k次，直到uk（t）為單調(diào)函數(shù)為止。

從以上步驟可以看出，原始信號(hào)可以由uk（t）和所有PF分量重構(gòu)，即：

2 基于局部均值的瞬時(shí)頻率求取

由式（11）可知信號(hào)可分解為若干PF分量之和，而每個(gè)PF分量由純包絡(luò)信號(hào)a（t）和純調(diào)頻函數(shù)s（t）=cos φ（t）組成，其頻率 f可由純調(diào)頻函數(shù) s（t）直接求解，即：

將式（12）展開并求導(dǎo)，則可以求出 s（t）的瞬時(shí)頻率，即相應(yīng)PF分量的瞬時(shí)頻率。因?yàn)槭峭ㄟ^對(duì)余弦函數(shù)求導(dǎo)得出PF的瞬時(shí)頻率，所以要求s（t）值在 ±1之間，如果 s（t）值近似等于 ±1，則用 ±1代替。這種求取頻率的方法直觀而且簡(jiǎn)單，簡(jiǎn)稱直接法。和采用HHT求取瞬時(shí)頻率的方法相比，直接法求取的瞬時(shí)頻率始終是正值，不會(huì)出現(xiàn)HHT中負(fù)頻率的現(xiàn)象。

3 基于LMD的仿真與實(shí)例分析

3.1 間諧波信號(hào)

設(shè)間諧波信號(hào)如下：

其中，A1=1 V，A2=0.3 V，A3=0.1 V；m1=1，m2=5.5，m3=11；f=50 Hz。

采樣頻率為6400 Hz，采樣時(shí)長(zhǎng)為0.4 s，變動(dòng)量Δ=0.001，為驗(yàn)證LMD的優(yōu)點(diǎn)，端點(diǎn)未作處理。運(yùn)用LMD分離出的各PF分量和相應(yīng)的純調(diào)頻函數(shù)如圖1（a）所示，可知LMD按頻率遞減的順序逐級(jí)分離出PF1、PF2和PF3分量；由純調(diào)頻函數(shù)求得的各PF分量的瞬時(shí)頻率函數(shù)和采用HHT方法求取的瞬時(shí)頻率函數(shù)如圖1（b）所示；分別運(yùn)用LMD和HHT求取的瞬時(shí)幅值函數(shù)如圖 1（c）所示。由圖 1（b）、（c）可知，用LMD分析間諧波時(shí)即使端點(diǎn)不作處理，分解得到的波形在端點(diǎn)處的效果可以與采用鏡像拓延后的HHT方法相比擬；另外，這2種方法求取的瞬時(shí)頻率和幅值都存在微小的波動(dòng)，但根據(jù)LMD方法求取的瞬時(shí)特征參數(shù)的波動(dòng)幅度明顯較小。究其原因，一是LMD分解信號(hào)時(shí)采用平均滑動(dòng)法擬合包絡(luò)線，避免了HHT中采用三次插值時(shí)易出現(xiàn)過包絡(luò)和欠包絡(luò)現(xiàn)象，其端點(diǎn)效應(yīng)較小；二是LMD通過除以包絡(luò)均值函數(shù)得到PF分量，相對(duì)于HHT中采用減法的篩分過程，LMD得到一個(gè)PF分量時(shí)篩分次數(shù)較少，而篩分次數(shù)越少，端點(diǎn)效應(yīng)向內(nèi)污染的程度就越輕。

對(duì)HHT和LMD中所求瞬時(shí)參數(shù)波動(dòng)問題，通常有2種方法：一種是取數(shù)據(jù)擬合后的穩(wěn)態(tài)值作為檢測(cè)值；另一種是將除端點(diǎn)部分的波峰值和波谷值相加，然后取平均值作為檢測(cè)值。本算例采用最小二乘擬合，限于篇幅，本文僅顯示了LMD所求間諧波頻率f2的經(jīng)數(shù)據(jù)擬合后的結(jié)果，如圖1（d）所示。LMD、HHT經(jīng)擬合后求取瞬時(shí)征參數(shù)以及加Blackman-Harries窗插值離散傅里葉算法求取的瞬時(shí)特征參數(shù)如表1所示。由表1可知，與HHT方法相比，LMD檢測(cè)的頻率和幅值準(zhǔn)確性較高；與加窗FFT相比，用LMD方法獲取的間諧波頻率準(zhǔn)確性略遜一籌，但就檢測(cè)的幅值而言，LMD方法比加窗FFT方法的檢測(cè)效果更好。

在裝有Win7的神舟K480P-i5G筆記本（可用內(nèi)存2.92 G）和MATLAB2011a仿真平臺(tái)上，從分離原始信號(hào)直至獲得信號(hào)的頻率和幅值，LMD方法耗時(shí)1.148 769 s，HHT方法（HHT篩分停止條件采用Rilling 提出的準(zhǔn)則［21］）運(yùn)行了 2.203376 s，加窗 FFT花費(fèi)了0.370305 s。與HHT相比，LMD方法的運(yùn)行時(shí)間減少了1.054607 s。一方面，LMD篩分次數(shù)較少，另一方面，LMD通過對(duì)調(diào)頻函數(shù)求導(dǎo)獲得瞬時(shí)頻率，而HHT方法求取信號(hào)的瞬時(shí)特征參數(shù)時(shí)需用Hilbert變換及其逆變換，所以HHT耗時(shí)較多。但與加窗FFT相比，LMD和HHT耗費(fèi)的時(shí)間都較多。3種方法采樣頻率為6400 Hz。

圖1 間諧波信號(hào)分析結(jié)果Fig.1 Results of interharmonic signal analysis

表1 檢測(cè)結(jié)果Tab.1 Results of detection

3.2 基于LMD的諧波失真信號(hào)分析［11］

為了和基于HHT的時(shí)頻分析方法相對(duì)比，本文所取諧波失真信號(hào)和文獻(xiàn)［11］中的前3個(gè)算例信號(hào)相類似。

3.2.1 短時(shí)諧波失真信號(hào)

設(shè)原始信號(hào)x（t）由基波信號(hào)和3次短時(shí)（0.08 s≤t≤0.2 s）諧波失真信號(hào)組成，基波信號(hào)的振幅為1 V，3次諧波信號(hào)的振幅為0.5 V；其波形如圖2中的x（t）所示。采樣頻率為3200 Hz，運(yùn)用LMD處理x（t）得到PF分量如圖2中PF1、PF2所示；由LMD得到的PF分量按局部頻率從高到低的順序排列，故高頻諧波信號(hào)在PF1分量中，由于PF分量由調(diào)幅函數(shù)和調(diào)頻函數(shù)組成，故分析PF1分量可知基波和擾動(dòng)諧波的幅值、頻率以及擾動(dòng)發(fā)生和恢復(fù)的時(shí)刻，PF1分量的幅值函數(shù)和瞬時(shí)頻率分別如圖2中的a1和f1所示。

圖2 基于LMD的短時(shí)諧波信號(hào)分析Fig.2 Short-time harmonic signal analysis based on LMD

從圖2中的PF1分量可以看出在0.08～0.2 s時(shí)段內(nèi)的信號(hào)成分是擾動(dòng)的諧波分量，其余時(shí)間的信號(hào)為基波分量；由PF2分量可知，在0.08～0.2 s時(shí)段內(nèi)的信號(hào)成分是基波信號(hào)；圖2中的a1和f1存在不同程度的波動(dòng)，經(jīng)數(shù)據(jù)擬合后得知，本文檢測(cè)的基波信號(hào)幅值為0.9990 V，3次諧波幅值為0.4996 V，基波頻率為50.07 Hz，3次短時(shí)諧波頻率為150.01 Hz，另外，由瞬時(shí)頻率函數(shù)f1還可知擾動(dòng)信號(hào)發(fā)生時(shí)刻 0.08002 s，恢復(fù)時(shí)刻在0.2060 s。

和文獻(xiàn)［11］采用HHT方法相比，LMD方法得到的瞬時(shí)頻率更準(zhǔn)確，而且求取的頻率是正的、連續(xù)的、具有物理意義的時(shí)變頻率，沒有出現(xiàn)無法解釋的負(fù)頻率。

3.2.2 暫態(tài)諧波失真信號(hào)

圖3中x（t）所示的暫態(tài)諧波失真信號(hào)表達(dá)式為：

其中，當(dāng) 0.12 s≤t≤0.24 s時(shí) p（t）=1，其他時(shí)刻 p（t）=0；A1=1 V，A2=0.5 V，A3=0.3 V；m2=3，m3=7；f=50 Hz。采樣頻率為3200 Hz。用LMD求取的暫態(tài)諧波失真信號(hào)時(shí)頻分析結(jié)果如圖3所示。

圖3 基于LMD的暫態(tài)諧波信號(hào)時(shí)頻分析Fig.3 Transient harmonic signal time-frequency analysis based on LMD

從圖3中的PF1和PF2分量可知，高次的諧波擾動(dòng)分量包含在PF1分量中，低頻次的諧波擾動(dòng)信號(hào)在PF2分量中；由圖3中PF1的瞬時(shí)頻率函數(shù)f1可知高次諧波擾動(dòng)發(fā)生時(shí)刻為0.1205 s，在0.2418 s恢復(fù)；由圖3中PF2的瞬時(shí)頻率函數(shù)f2知低頻諧波擾動(dòng)信號(hào)產(chǎn)生于0.1206 s，終止于0.2417 s。擬合PF1和PF2的瞬時(shí)特征函數(shù)可知高頻擾動(dòng)信號(hào)頻率為349.7 Hz，振幅為0.2993 V；低頻擾動(dòng)信號(hào)的頻率為150.1 Hz，振幅為0.4985 V。

從以上檢測(cè)結(jié)果可知本文所提方法可準(zhǔn)確確定多諧波擾動(dòng)信號(hào)中每個(gè)諧波的頻率、振幅和發(fā)生及恢復(fù)的時(shí)刻，并且頻率沒有發(fā)生突變，是正的、連續(xù)的、具有物理意義的時(shí)變頻率。

3.3 時(shí)變諧波信號(hào)

時(shí)變諧波信號(hào)如圖4中的x（t）所示。基波信號(hào)在0.12 s幅值從1 V驟升為2 V，在0.24 s幅值又驟降為1 V；在0.12 s加入幅值為1.6 V的3次諧波和幅值為1.2 V的7次諧波信號(hào)；0.24 s時(shí)3次諧波和7次諧波幅值分別驟降為0.8 V和0.6 V；時(shí)間長(zhǎng)度為0.4 s，采樣頻率為3200 Hz，用LMD求取的時(shí)變諧波信號(hào)時(shí)頻分析結(jié)果如圖4所示。

圖4 基于LMD的時(shí)變諧波信號(hào)分析結(jié)果Fig.4 Results of time-varying harmonic signal time-frequency analysis based on LMD

由圖4中PF1和PF2可知在0.12～0.4 s時(shí)段，PF1分量對(duì)應(yīng)7次諧波信號(hào)，PF2分量對(duì)應(yīng)3次諧波信號(hào)。對(duì)圖4中的中PF1的瞬時(shí)特征函數(shù)a1和f1進(jìn)行數(shù)據(jù)擬合，可知7次諧波信號(hào)的幅值在0.12～0.24 s時(shí)段為 1.198 6 V，在0.24～0.4 s時(shí)段為0.6002 V；0.12～0.4 s時(shí)段檢測(cè)的7次諧波頻率為349.7 Hz。運(yùn)用相同的方法，可知在0.12～0.4 s時(shí)段3次諧波信號(hào)的幅值為1.5994 V，0.24～0.4s時(shí)段為0.7992V；檢測(cè)的3次諧波頻率為149.8 Hz。從而可知，本文所提方法還可以確定多頻諧波時(shí)變信號(hào)中每個(gè)時(shí)變信號(hào)的幅值、頻率以及擾動(dòng)發(fā)生和終止的時(shí)刻。

3.4 實(shí)例信號(hào)分析

某500 kV變電站采用電容式電壓互感器對(duì)電壓進(jìn)行測(cè)量，電容式電壓互感器二次側(cè)裝有諧振性阻尼器，35 kV母線分南母、北母2段。在進(jìn)行35 kV南母電容器組投運(yùn)操作時(shí)發(fā)生35 kV南母中間電壓互感器高壓保險(xiǎn)B相熔斷事故。圖5給出了35 kV南母進(jìn)行電容器組投切試驗(yàn)時(shí)采用尼高力波形記錄儀（Nicolet vision）錄取的35 kV側(cè)B相電壓暫態(tài)波形（采樣頻率為2 kHz，本文取最初的10個(gè)周期）及其對(duì)應(yīng)的頻譜圖?；贚MD的時(shí)頻分析結(jié)果如圖6所示。

由圖6（a）可知運(yùn)用LMD可將35 kV側(cè)B相電壓分解為 PF1和 PF2分量；由圖 6（b）中 PF1瞬時(shí)頻率函數(shù)f1可以看出投切電容器前的頻率為基波分量，投切后PF1的瞬時(shí)頻率增大，為3次諧波成分，除端點(diǎn)外將波動(dòng)的3次諧波頻率函數(shù)的波峰值與波谷值相加求取的平均值為150.15 Hz，由此可知，運(yùn)用LMD方法檢測(cè)信號(hào)頻率時(shí)有較高的精度；由PF2的瞬時(shí)頻率函數(shù)f2可知電容器投切后B相電壓有基波分量；從圖6（c）PF1的瞬時(shí)幅值函數(shù)可知沒有連接電容器時(shí)的基波幅值為30 kV，這與圖5（a）中基波幅值很吻合，接入電容器后的3次諧波幅值較大，可知電容器投切后B相電壓含有較多的3次諧波成分，而3次諧波成分較多是鐵磁元件嚴(yán)重飽和的一個(gè)重要特征；另外，圖6（c）從幅值曲線可以看出在投切電容器時(shí)引起了基波電壓下降，這和理論情況相符，采用平均值法求取的電壓下降后PF2的瞬時(shí)幅值為21.8 kV，可知基波電壓從30 kV下降到21.8 kV；另外，從圖 5（a）可知電容器在采樣點(diǎn) 129（對(duì)應(yīng)0.0645 s）時(shí)投運(yùn)，從圖 6（b）PF1的瞬時(shí)頻率函數(shù) f1可知在133采樣點(diǎn)（對(duì)應(yīng)0.0665 s）時(shí)產(chǎn)生了3次諧波，從圖6（c）PF1的瞬時(shí)幅值函數(shù)也可知過132采樣點(diǎn)（0.066 s）后的幅值在21.8 kV小范圍內(nèi)波動(dòng)，由此可知，根據(jù)PF1瞬時(shí)頻率f1或瞬時(shí)幅值函數(shù)可以確定基波下降的時(shí)刻和產(chǎn)生3次諧波的時(shí)刻。

另外，從圖 6（b）可直觀看出，與 HHT方法相比，本文所提方法求取的瞬時(shí)頻率具有端點(diǎn)效應(yīng)小、波動(dòng)程度輕和更接近真實(shí)值等優(yōu)點(diǎn)。

圖5 35 kV側(cè)B相電壓波形及其頻譜圖Fig.5 Waveform of 35 kV-side phase-B voltage and its spectrum

圖6 35 kV側(cè)B相電壓的LMD時(shí)頻分析結(jié)果Fig.6 Results of time-frequency analysis based on LMD for 35 kV-side phase-B voltage

4 結(jié)論

本文首次利用LMD方法對(duì)間諧波信號(hào)、單頻諧波信號(hào)、多頻諧波信號(hào)、多頻時(shí)變諧波信號(hào)和電容器組投運(yùn)時(shí)的35 kV側(cè)B相電壓信號(hào)進(jìn)行分析，從仿真及實(shí)驗(yàn)結(jié)果看出該方法可準(zhǔn)確判斷諧波擾動(dòng)發(fā)生和恢復(fù)的時(shí)刻以及諧波信號(hào)的幅值和頻率。間諧波信號(hào)仿真分析結(jié)果表明，與HHT相比，LMD方法具有運(yùn)算量小、速度快、求取的瞬時(shí)特征參數(shù)波動(dòng)幅度小、檢測(cè)精度高等優(yōu)點(diǎn)；從實(shí)際信號(hào)的頻率分析對(duì)比結(jié)果可知，和HHT方法相比，本文所提方法具有很大的優(yōu)越性。更重要的是用LMD方法所求取的瞬時(shí)頻率都是正值，而且所得波形受端點(diǎn)效應(yīng)的影響較小，為電力系統(tǒng)電能質(zhì)量檢測(cè)提供了一種新的途徑。值得注意的是，作為一種新提出的非線性信號(hào)分析方法，LMD還有迭代終止條件、滑動(dòng)平均跨度選取（本文取相鄰極值點(diǎn)最長(zhǎng)距離的1/3）和端點(diǎn)效應(yīng)等理論問題需進(jìn)一步研究。雖然LMD的端點(diǎn)效應(yīng)較小，但諧波相位分析是建立在準(zhǔn)確獲取信號(hào)初值之上，導(dǎo)致該方法暫不能用于分析諧波相位。隨著這些理論問題的解決，LMD將有良好的應(yīng)用前景。