李 明 王曉茹

（西南交通大學(xué)電氣工程學(xué)院 成都 610031）

1 引言

電力系統(tǒng)中不僅存在大量的整數(shù)次諧波，也存在著非整數(shù)次諧波（次諧波，間諧波）。間諧波對(duì)電力系統(tǒng)及設(shè)備的危害很大，會(huì)引起燈光閃爍、低頻繼電器異常運(yùn)行、無(wú)源濾波器過(guò)流跳閘、感應(yīng)電動(dòng)機(jī)噪聲和振動(dòng)等問(wèn)題。因此，需對(duì)電力系統(tǒng)間諧波進(jìn)行分析，從而為電力系統(tǒng)的監(jiān)控和保護(hù)提供依 據(jù)[1]。

快速傅里葉變換（Fast Fourier Transform，F(xiàn)FT）是最常用的電力系統(tǒng)諧波分析方法，但其頻率分辨率較低，并且，非同步采樣的情況下，還會(huì)導(dǎo)致長(zhǎng)范圍泄漏和柵欄效應(yīng)[2]。文獻(xiàn)[3]采用加窗插值的方式來(lái)對(duì)FFT 進(jìn)行改進(jìn)，可減少長(zhǎng)范圍泄漏和柵欄效應(yīng)的影響。但通過(guò)加窗插值的方式不能提高算法的頻率分辨率，并且，在時(shí)域內(nèi)加平滑窗還會(huì)進(jìn)一步降低頻率分辨率。為了檢測(cè)出信號(hào)中所有的諧波和間諧波分量，窗寬一般高達(dá)幾十個(gè)信號(hào)周期，不利于實(shí)時(shí)檢測(cè)。特征值算法[4-8]通過(guò)將自相關(guān)矩陣中的信息空間分解為信號(hào)子空間和噪聲子空間，從而達(dá)到諧波分析的目的，具有較高的頻率分辨率，故被廣泛應(yīng)用于電力系統(tǒng)間諧波檢測(cè)中。但其運(yùn)算量較大，硬件實(shí)現(xiàn)比較困難。小波變換[9-10]和時(shí)頻分析[11]也被用于分析電力系統(tǒng)間諧波，但頻率分辨率和估計(jì)精度較低，運(yùn)算量也較大。

自回歸（Autoregressive，AR）模型譜估計(jì)方法中的Burg 算法具有很高的頻率分辨率，并且由于采用Levinson 遞推，可進(jìn)行高效計(jì)算，故被應(yīng)用于電力系統(tǒng)間諧波譜估計(jì)中[12-13]。但Burg 算法對(duì)初始相位較為敏感，在對(duì)較短的間諧波信號(hào)進(jìn)行分析時(shí)尤為明顯，導(dǎo)致出現(xiàn)譜峰偏移，并在階數(shù)較大時(shí)易出現(xiàn)譜線(xiàn)分裂。針對(duì)Burg 算法的這一缺點(diǎn)，已有多種改進(jìn)算法，Marple 算法[14]解除了Levinson 遞推公式這一強(qiáng)約束條件，其頻率偏移程度較低，且不存在譜線(xiàn)分裂現(xiàn)象，譜估計(jì)性得到提高，但是，該算法不能保證得到一個(gè)穩(wěn)定的AR 模型，雖然可以設(shè)定遞推結(jié)束條件來(lái)保證算法的穩(wěn)定性，但會(huì)使譜估計(jì)結(jié)果趨于保守，可能會(huì)丟失一些譜峰。漢明窗Burg算法[15]和最優(yōu)窗Burg 算法[16]通過(guò)對(duì)預(yù)測(cè)誤差平均功率進(jìn)行加窗，可以降低頻率偏移和譜線(xiàn)分裂程度，并且可以得到一個(gè)穩(wěn)定的AR 模型，而最優(yōu)窗Burg算法通過(guò)使平均頻率誤差最小，其譜估計(jì)性能要進(jìn)一步優(yōu)于漢明窗Burg 算法。

本文采用最優(yōu)窗Burg 算法，將其應(yīng)用于電力系統(tǒng)間諧波譜估計(jì)。與Burg 算法相比，該算法對(duì)初始相位不敏感，頻率偏移和譜線(xiàn)分裂程度較低，譜估計(jì)性能較好，與各種特征值算法相比，其計(jì)算復(fù)雜度較低。仿真結(jié)果驗(yàn)證了該方法的有效性。

2 基于加窗Burg 算法的間諧波檢測(cè)原理

設(shè)電力系統(tǒng)間諧波信號(hào)為

式中，M 為信號(hào)中所含諧波和間諧波的個(gè)數(shù)；Ai，fi，?i分別為第i 次諧波的幅值、頻率和初始相位；fs為采樣頻率；η(n)為白噪聲序列。

式（1）可轉(zhuǎn)化為[17]

由式（2）可見(jiàn)，電力系統(tǒng)間諧波信號(hào)y(n)可看作為AR 模型。另外，由于事先并不知道信號(hào)所含諧波和間諧波的個(gè)數(shù)，所以AR 模型階數(shù)m 需要預(yù)估，現(xiàn)有的AR 模型階次準(zhǔn)則（如最終預(yù)測(cè)誤差準(zhǔn)則，Akaike 信息準(zhǔn)則等）所預(yù)估的結(jié)果都不理想，一個(gè)經(jīng)驗(yàn)法則是，對(duì)于一個(gè)較短的間諧波信號(hào)，可令m 的取值范圍為N/3～N/2（N 為信號(hào)采樣點(diǎn)數(shù)），這樣得到的譜估計(jì)效果較好[18]。

AR 模型中較常用的方法為Burg 算法，該方法利用Levinson 遞推公式由低階到高階來(lái)求取AR 模型參數(shù)，Levinson 遞推公式為[16,19]

式中，i=1,2,…,m?1。

由式（3）可以看出，在m?1 階解已知的情況下，欲求出m 階解，僅有am(m)未知，所以，關(guān)鍵的問(wèn)題是如何求取反射系數(shù)am(m)。Burg 算法通過(guò)在前、后向預(yù)測(cè)誤差平均功率最小的意義下直接求解am(m)。首先定義前、后向預(yù)測(cè)誤差分別為

將式（3）代入到式（4）和式（5）中，可得

式（6）和式（7）分別為前、后項(xiàng)預(yù)測(cè)誤差的階數(shù)遞推公式，則m 階預(yù)測(cè)誤差平均功率為

對(duì)預(yù)測(cè)誤差平均功率em進(jìn)行加窗后可得[18]

式中，wm(n)為m 階窗函數(shù)。如果wm(n)選取合適，可明顯降低譜峰偏移程度，具體原因詳見(jiàn)第3 節(jié)。

由式（10）可見(jiàn)，|am(m)|＜1，因此，加窗Burg 算法可以得到一個(gè)穩(wěn)定的AR 模型[18]。將am(m)代入式（3）中，可求出AR 模型的m 階解am(i)。根據(jù)所求得的am(i)和預(yù)測(cè)誤差em可進(jìn)一步求得信號(hào)的功率譜密度

由于數(shù)字頻率ω為連續(xù)函數(shù)，需要對(duì)其進(jìn)行數(shù)字離散化。將式（11）轉(zhuǎn)化為

數(shù)字頻率間隔Δω 可以根據(jù)需要設(shè)定，Δω 越小，頻率估計(jì)精度越高，同時(shí)計(jì)算量也越大。為了兼顧頻率估計(jì)精度和計(jì)算復(fù)雜度，本文假定Δω為π/fs，則模擬頻率間隔為0.5Hz。在式（12）所對(duì)應(yīng)的功率譜圖上，各峰值處所對(duì)應(yīng)的頻率即為各次諧波和間諧波的頻率。在此基礎(chǔ)上，可用非線(xiàn)性最小二乘法進(jìn)一步求取各次諧波和間諧波的幅值和相位信息，推導(dǎo)過(guò)程詳見(jiàn)文獻(xiàn)[8]。

3 誤差分析和最優(yōu)窗的推導(dǎo)

由式（1）可知，間諧波信號(hào)可以看作多個(gè)正弦波和噪聲疊加之和。為簡(jiǎn)單分析起見(jiàn)，只選取其中一 個(gè) 正 弦 波 x(n)=Aisin(2πfin/fs+?i)， 令 Ai=1，ω=2πfi/fs，? =?i?π/2，則x(n)=cos(ωn+?)。由于前、后向預(yù)測(cè)誤差的初值均為x(n)，將其代入到式（10），可得一階系數(shù)[16]

式（13）可表示為

式中

由式（14）可知，由于存在θ 項(xiàng)，從而使數(shù)字頻率ω 偏離其實(shí)際值，導(dǎo)致譜峰偏移現(xiàn)象。由式（15）可知，誤差項(xiàng)θ 的大小與信號(hào)初始相位? 密切相關(guān)?？赏ㄟ^(guò)加合適的窗函數(shù)w1(n)來(lái)減小誤差項(xiàng)θ 的影響。首先令頻率誤差Δf=θ/2π，可得平均頻率誤差的方差為

式中，δ 為脈沖函數(shù)。為了盡量減小由θ 項(xiàng)引起的誤差，需使平均頻率誤差方差〈var(Δf)〉最小，可求得一階最優(yōu)窗

同理，可求得m 階最優(yōu)窗

最優(yōu)窗的具體推導(dǎo)見(jiàn)文獻(xiàn)[16]。

為了減少運(yùn)算量，可將式（18）轉(zhuǎn)化為如下遞推形式：

式中

由式（18）可以看出，窗函數(shù) wm(n)只與信號(hào)長(zhǎng)度N 和階數(shù)m 相關(guān)，而與信號(hào)的頻率、幅值和相位無(wú)關(guān)。故用該窗分析式（1）中的任意次諧波和間諧波時(shí)，均可有效減少由誤差項(xiàng)θ 引起的頻率偏差，從而減少間諧波信號(hào)的譜峰偏移，并可有效抑制譜線(xiàn)分裂。

4 算法復(fù)雜度比較

對(duì)最優(yōu)窗 Burg 算法分析可知，由于采用Levinson 遞推求取AR 模型系數(shù)，從而避免了自相關(guān)矩陣估計(jì)，相對(duì)于各種特征值法（如 Pisarenko算法、Music 算法、Esprit 算法，Min-Norm 算法等），其算法復(fù)雜度較低。為了進(jìn)一步減少最優(yōu)窗 Burg算法運(yùn)算量，可令式（12）中的比例系數(shù)em為一非零常數(shù)。另外，雖然各種特征值法的運(yùn)算量和存儲(chǔ)量不同，但均需進(jìn)行自相關(guān)矩陣估計(jì)，故可用自相關(guān)矩陣估計(jì)來(lái)代替各種特征值法。具體算法復(fù)雜度比較見(jiàn)表1 和表2。

表1 各種算法復(fù)雜度比較 Tab.1 The comparison of computational complexity for various algorithms

表2 各種算法復(fù)雜度定量比較 Tab.2 The quantitative comparison of computational complexity for various algorithms

由表1 和表2 可見(jiàn)，特征值算法僅自相關(guān)矩陣估計(jì)這一步，所需的運(yùn)算量和存儲(chǔ)量均已明顯高于最優(yōu)窗Burg 算法，如果再將矩陣的特征值分解等運(yùn)算包含在內(nèi)，其計(jì)算復(fù)雜度將會(huì)更高。并且，兩者計(jì)算復(fù)雜度的差距還會(huì)隨著信號(hào)采樣點(diǎn)數(shù)N 和AR模型階數(shù)m 的增加而進(jìn)一步增加。

5 仿真算例

5.1 最優(yōu)窗Burg 算法、漢明窗Burg 算法和原Burg算法譜估計(jì)性能比較

設(shè)分析信號(hào)為：x(t)=sin(2π×50t+5π/4)，采樣頻率為1000Hz，采樣點(diǎn)數(shù)為45。當(dāng)AR 模型階數(shù)為2 時(shí)，結(jié)果如圖1a 所示。當(dāng)AR 模型階數(shù)為4 時(shí)，結(jié)果如圖1b 所示。

從圖1a 可以看出，由于Burg 算法檢測(cè)到的譜峰位置與相位密切相關(guān)。故受到信號(hào)初始相位的影響，出現(xiàn)明顯的譜峰偏移。圖上檢測(cè)到的譜峰為47Hz，與實(shí)際譜峰的偏差為3Hz。偏移程度較大。加窗Burg 算法通過(guò)對(duì)預(yù)測(cè)誤差平均功率進(jìn)行加窗處理，從而降低了譜峰偏移程度。漢明窗Burg 算法和最優(yōu)窗Burg 算法檢測(cè)的譜峰分別為49.5Hz 和50Hz，譜估計(jì)性能明顯優(yōu)于Burg 算法。而最優(yōu)窗Burg 算法通過(guò)使平均頻率誤差方差最小化，從而使式（14）中的誤差項(xiàng)達(dá)到最小，譜估計(jì)性能比漢明 窗Burg 算法進(jìn)一步提高。

圖1 功率譜估計(jì)結(jié)果比較 Fig.1 Comparison of the results of power spectral estimation

從圖1b 可以看出，當(dāng)AR 模型階數(shù)增大為4 時(shí)，原Burg 算法出現(xiàn)明顯的譜線(xiàn)分裂，檢測(cè)到40Hz 和48.5Hz 兩個(gè)譜峰，而漢明窗Burg 算法和最優(yōu)窗Burg算法均檢測(cè)到一個(gè)譜峰，說(shuō)明加窗可以有效抑制譜線(xiàn)分裂。但漢明窗Burg 算法檢測(cè)到的譜峰為49Hz。頻率偏差為1Hz，而最優(yōu)窗Burg 算法檢測(cè)到的譜峰為50Hz。說(shuō)明加最優(yōu)窗不僅可以抑制譜線(xiàn)分裂，而且可以最大程度減少譜峰偏移。

5.2 多個(gè)諧波和間諧波分量檢測(cè)

設(shè)信號(hào)除基波（50Hz）外還含有1.3、3、5、5.2、7、9 次諧波分量。幅值分別為基波的7%、4%、1.3%、2%、1%、0.3%，相位分別為π/8、π/4、π/8、π/7、π/4、0，基波相位為π/3。信號(hào)的采樣頻率為1000Hz，采樣點(diǎn)數(shù)N 為51。AR 模型預(yù)估階數(shù)為20，信號(hào)中加入了60dB 的高斯白噪聲，采用最優(yōu)窗Burg算法與FFT 算法得到的結(jié)果如圖2 所示。

圖2 功率譜估計(jì)結(jié)果 Fig.2 The results of power spectral estimation

從圖2b 可以看出，圖上共有4 個(gè)譜峰，分別在基波、3 次諧波、5 次諧波和7 次諧波附近。9 次諧波幅值較低，在非同步采樣的情況下，受到頻譜泄漏的影響從而檢測(cè)不到。由于觀(guān)測(cè)窗較短，頻率分辨率較低，導(dǎo)致基波和1.3 次諧波，5 次諧波和5.2次諧波無(wú)法區(qū)分開(kāi)來(lái)。3 次諧波和7 次諧波雖然可以檢測(cè)出來(lái)，但檢測(cè)到的頻率為 156.8627Hz、352.9412Hz，與實(shí)際頻率偏差較大。從圖2a 可以看出，最優(yōu)窗Burg 算法可以檢測(cè)到頻率為50、64.5、150、249.5、260.5、350、450.5Hz 的7 個(gè)譜峰，所有諧波和間諧波分量都被檢測(cè)出來(lái)。與實(shí)際譜峰相比，平均頻率偏差僅為0.2857Hz，檢測(cè)結(jié)果較為精確，譜估計(jì)性能明顯優(yōu)于FFT。

6 結(jié)論

本文提出了最優(yōu)窗Burg 算法的間諧波譜估計(jì)方法。仿真結(jié)果表明，對(duì)于較短的間諧波信號(hào)，采用該算法可以精確檢測(cè)出各次諧波和間諧波的頻率信息，譜估計(jì)性能較好。與Burg 算法相比，加窗Burg 算法可以有效降低譜峰偏移程度、抑制譜線(xiàn)分裂。而選擇最優(yōu)窗獲得的譜估計(jì)效果要優(yōu)于漢明窗。因最優(yōu)窗Burg 算法利用了Levinson 遞推，計(jì)算效率高，相對(duì)于各種特征值算法，更有利于硬件實(shí)現(xiàn)。

[1] 林海雪.電力系統(tǒng)中的間諧波問(wèn)題[J].供用電,2001,18(3):6-9.

Lin Haixue.Interharmonics in electrical power system[J].Distribution and Utilization,2001,18(3):6-9.

[2] 龐浩,李東霞,俎云霄,等.應(yīng)用 FFT 進(jìn)行電力系統(tǒng)諧波分析的改進(jìn)算法[J].中國(guó)電機(jī)工程學(xué)報(bào),2003,23(6):50-54.

Pang Hao,Li Dongxia,Zu Yunxiao,et al.An improved algorithm for harmonic analysis of power system using FFT technique[J].Proceedings of the CSEE,2003,23(6):50-54.

[3] 祁才君,王小海.基于插值FFT 算法的間諧波參數(shù)估計(jì)[J].電工技術(shù)學(xué)報(bào),2003,18(1):92-95.

Qi Caijun,Wang Xiaohai.Interharmonics estimation based on interpolation FFT algorithm[J].Transactions of China Electrotechnical Society,2003,18(1):92-95.

[4] 王志群,朱守真,周雙喜.基于Pisarenko 諧波分解的間諧波估算方法[J].電網(wǎng)技術(shù),2004,28(15):72-77.

Wang Zhiqun,Zhu Shouzhen,Zhou Shuangxi.Inter-harmonics estimation by Pisarenko harmonic decomposition method[J].Power System Technology,2004,28(15):72-77.

[5] 高培生,谷湘文,吳為麟.基于空間譜和支持向量回歸機(jī)的間諧波分析[J].電力系統(tǒng)自動(dòng)化,2007,31(24):67-70.

Gao Peisheng,Gu Xiangwen,Wu Weilin.Interharmonic analysis based on spatial spectrum and support vector regression machine[J].Automation of Electric Power Systems,2007,31(24):67-70.

[6] 沈睿佼,楊洪耕.基于特征空間求根法的非整數(shù)次諧波估計(jì)方法[J].中國(guó)電機(jī)工程學(xué)報(bào),2006,26(24):72-76.

Shen Ruijiao,Yang Honggeng.A new method for non-integer harmonics measurement based on root-eigenspace method[J].Proceedings of the CSEE,2006,26(24):72-76.

[7] 金國(guó)彬,李玲,李天云,等.間諧波高精度檢測(cè)新方法[J].電力系統(tǒng)及其自動(dòng)化學(xué)報(bào),2009,21(2):25-30.

Jin Guobin,Li Ling,Li Tianyun,et al.Novel method of high-accuracy detection for interharmonics[J].Proceedings of the CSU-EPSA,2009,21(2):25-30.

[8] 李晶,裴亮,郁道銀,等.一種用于電力系統(tǒng)諧波與間諧波分析的超分辨率算法[J].中國(guó)電機(jī)工程學(xué)報(bào),2006,26(15):35-39.

Li Jing,Pei Liang,Yu Daoyin,et al.A super- resolution algorithm for harmonic and inter-harmonic analysis in power system[J].Proceedings of the CSEE,2006,26(15):35-39.

[9] 任震,黃群古,黃雯瑩,等.基于多頻帶小波變換的電力系統(tǒng)諧波分析新方法[J].中國(guó)電機(jī)工程學(xué)報(bào),2000,20(12):38-41.

Ren Zhen,Huang Qungu,Huang Wenying,et al.New methods of power system harmonics analysis based on wavelet transform with multi frequency band[J].Proceedings of the CSEE,2000,20(12):38-41.

[10] 薛蕙,楊仁剛.基于連續(xù)小波變換的非整數(shù)次諧波測(cè)量方法[J].電力系統(tǒng)自動(dòng)化,2003,27(5):49-53.

Xue Hui,Yang Rengang.A novel method for non-integer hamonics measurement using continuous wavelet transform[J].Automation of Electric Power Systems,2003,27(5):49-53.

[11] 張宇輝,金國(guó)彬,李天云.基于自適應(yīng)最優(yōu)核時(shí)頻分布理論的間諧波分析新方法[J].中國(guó)電機(jī)工程學(xué)報(bào),2006,26(18):84-89.

Zhang Yuhui,Jin Guobin,Li Tianyun.A novel approach to interharmonics analysis based on adaptive optimal kernel time-frequency distribution[J].Proceedings of the CSEE,2006,26(18):84-89.

[12] 蔡忠法,陳隆道.基于AR 譜估計(jì)和Adaline 神經(jīng)元的間諧波分析[J].電力系統(tǒng)自動(dòng)化,2007,31(17):78-82.

Cai Zhongfa,Chen Longdao.Inter-harmonic analysis based on AR spectral estimation and Adaline neural network[J].Automation of Electric Power Systems,2007,31(17):78-82.

[13] 馬秉偉,劉會(huì)金,周莉,等.一種基于自回歸模型的間諧波譜估計(jì)的改進(jìn)算法[J].中國(guó)電機(jī)工程學(xué)報(bào),2005,25(15):79-83.

Ma Bingwei,Liu Huijin,Zhou Li,et al.An improved algorithm of interharmonics spectral estimation based on AR mode[J].Proceedings of the CSEE,2005,25(15):79-83.

[14] Marple L.A new autoregressive spectrum analysis algorithm[J].IEEE Transactions on Acoustics,Speech,and Signal Processing,1980,ASSP-28(4):441-454.

[15] Swingler D N.A modified Burg algorithm for maximum entropy spectral analysis[J].Proc.IEEE,1979,67(9):1368-1369.

[16] Kaveh M,Lippert G.An optimum tapered Burg algorithm for linear prediction and spectral analysis[J].IEEE Transactions on Acoustics,Speech,and Signal Processing,1983,ASSP-31(2):438-444.

[17] 張賢達(dá).現(xiàn)代信號(hào)處理[M].2 版.北京:清華大學(xué)出版社,1999.

[18] Proakis J G.統(tǒng)計(jì)信號(hào)處理算法[M].湯俊,譯.北京:清華大學(xué)出版社,2006.

[19] 皇甫堪,陳建文,樓生強(qiáng).現(xiàn)代數(shù)字信號(hào)處理[M].北京:電子工業(yè)出版社,2003.