季 宇， 孫云蓮， 李 晶

(1.武漢大學(xué)電氣工程學(xué)院， 武漢 430072； 2.五環(huán)化工有限公司， 武漢 430072)

RELAX算法在間諧波檢測(cè)中的應(yīng)用①

季 宇1， 孫云蓮1， 李 晶2

(1.武漢大學(xué)電氣工程學(xué)院， 武漢 430072； 2.五環(huán)化工有限公司， 武漢 430072)

間諧波是頻率非基波整數(shù)倍的諧波成分，傳統(tǒng)諧波測(cè)量方法將會(huì)因?yàn)榉峭讲蓸佣a(chǎn)生泄漏誤差。為準(zhǔn)確測(cè)量間諧波，文中提出基于松弛譜估計(jì)(RELAX)的檢測(cè)方法。因?yàn)樵摲椒ㄡ槍?duì)復(fù)信號(hào)進(jìn)行分析，在檢測(cè)前需先對(duì)采樣信號(hào)進(jìn)行希爾伯特變換構(gòu)成復(fù)序列，計(jì)算結(jié)束后的相位測(cè)量結(jié)果也要進(jìn)行修正。由于RELAX算法的復(fù)雜度集中在快速傅里葉變換(FFT)處，故相比于其它方法其計(jì)算速度是有優(yōu)勢(shì)的。最后通過(guò)對(duì)具體間諧波信號(hào)的仿真測(cè)量，驗(yàn)證了該方法的有效性和實(shí)用性。

間諧波； 非同步采樣； RELAX算法； 希爾伯特變換； 快速傅里葉變換

隨著電力系統(tǒng)中非線性負(fù)荷的大量投用，電力系統(tǒng)中的諧波情況也愈發(fā)復(fù)雜，不僅存在頻率是基頻整數(shù)次的整數(shù)次諧波，而且存在非整數(shù)次諧波成份。這些對(duì)電力系統(tǒng)的安全、穩(wěn)定運(yùn)行造成了極大的危害。而對(duì)諧波間諧波進(jìn)行有效分析治理的基礎(chǔ)便是對(duì)其準(zhǔn)確的檢測(cè)。

目前電力系統(tǒng)中的間諧波源主要是以電弧爐為主的波動(dòng)性負(fù)荷及變頻調(diào)速裝置類(lèi)的非線性負(fù)荷。此外，高壓直流輸電技術(shù)從原理上分析等同于交直交的變流系統(tǒng)，因此也可劃為第二類(lèi)間諧波源。從頻域上分析，前者由于電弧阻抗的時(shí)變性導(dǎo)致頻譜呈連續(xù)譜，后者則是具有特征頻率的離散譜，且其頻率可通過(guò)互調(diào)理論進(jìn)行分析。無(wú)論是連續(xù)譜還是離散譜都具有間諧波的統(tǒng)一特征即頻率非工頻整數(shù)倍且幅值較小，這一特征決定了以往很多諧波檢測(cè)方法對(duì)于間諧波將不再適用。目前間諧波的檢測(cè)方法可以分為兩大類(lèi)，即非參數(shù)估計(jì)法和參數(shù)估計(jì)方法。前者包括快速傅里葉變換(FFT)、插值傅里葉變換[1～3]、小波分析[4]、支持向量機(jī)法[5]等；后者則主要是基于現(xiàn)代譜估計(jì)的方法[6]，包括奇異值最小二乘法、PRONY法、BURG法以及一大類(lèi)空間譜估計(jì)方法。這些方法中快速傅里葉變換作為最常用的方法在非同步采樣時(shí)會(huì)出現(xiàn)譜峰偏離，引起參數(shù)估計(jì)誤差。插值傅里葉方法針對(duì)傳統(tǒng)傅里葉變換進(jìn)行算法修正，通過(guò)加窗和插值的方法克服頻譜泄露及欄柵效應(yīng)。然而，這類(lèi)算法運(yùn)算量較大，一般需要求解高階方程。小波分析類(lèi)方法是將觀測(cè)信號(hào)利用子帶分離技術(shù)進(jìn)行頻帶劃分，得到各個(gè)頻率諧波。而小波濾波器的混疊性及噪聲敏感是這類(lèi)方法的缺點(diǎn)。支持向量機(jī)法是利用支持向量機(jī)的回歸特性，結(jié)合迭代變權(quán)方法估計(jì)諧波參數(shù)。這種方法能夠保證全局收斂性，且估計(jì)參數(shù)是全局最優(yōu)的，但參數(shù)選擇問(wèn)題限制了該方法的廣泛應(yīng)用。參數(shù)估計(jì)方法中奇異值最小二乘法和PRONY法可以直接測(cè)量間諧波的全部參數(shù)，但計(jì)算量較大，且由于在推導(dǎo)公式中存在假設(shè)條件(如白噪聲假設(shè))，因此實(shí)際檢測(cè)結(jié)果可能存在較大偏差。BURG法實(shí)質(zhì)是格型濾波算法，利用萊文森遞推公式計(jì)算模型參數(shù)。該方法避免了求解相關(guān)矩陣，但對(duì)于高階模型及大樣本情況分析效果較差。子空間法通過(guò)分離信號(hào)與噪聲空間給出觀測(cè)信號(hào)的偽譜，因此只能用于頻率檢測(cè)，且在信噪比較低時(shí)效果較差。本文中介紹的RELAX方法屬于非參數(shù)估計(jì)法，能夠準(zhǔn)確計(jì)算諧波的全部參數(shù)。同時(shí)，計(jì)算過(guò)程中利用了FFT，但克服了FFT非同步采樣時(shí)的估計(jì)誤差，且具有一定抗噪性能。

1 RELAX算法介紹[8，9]

RELAX算法是一種具有超分辨率的譜分析算法，該算法假設(shè)復(fù)諧波模型[3]為

n=0，1，…，N-1

(1)

其中αk是諧波幅值，此處可理解為帶相位信息，即αk=|αk|ejφ，en是觀測(cè)噪聲，K是諧波個(gè)數(shù)(實(shí)際運(yùn)算中可以先使用熵值估計(jì)方法或FFT方法對(duì)主要諧波個(gè)數(shù)進(jìn)行估計(jì))。在此基礎(chǔ)上將諧波參數(shù)估計(jì)問(wèn)題轉(zhuǎn)化為下面的優(yōu)化問(wèn)題：

J1(f1，α1，…，fk，αk)=

(2)

其中y=[y0y1…yN-1]T，ω(fk)=[1，ej2πfk，…，ej2πfk(N-1)]T，‖·‖2表示歐氏范數(shù)。設(shè)

P=[ω(f1)，ω(f2)，…，ω(fk)]

(3)

因?yàn)镻的列向量相互獨(dú)立，故可將范式中的表達(dá)式理解為一個(gè)向量投影問(wèn)題，則原優(yōu)化式可表示成

(4)

由泛函分析可知，當(dāng)投影向量作正交投影時(shí)可取得最小范數(shù)解，最小范數(shù)解(即最小二乘解)為

(5)

此時(shí)待估計(jì)的參數(shù)只有各個(gè)諧波頻率。設(shè)諧波頻率矢量f=[f1，f2，…，fk]T，及正交投影矩陣Ω=P(PHP)-1PH，優(yōu)化式最小范數(shù)值即為正交投影矩陣補(bǔ)空間上投影長(zhǎng)度?；诖?，優(yōu)化目標(biāo)簡(jiǎn)化為：

J2(f1，f2，…，fk)=‖Ω⊥y‖2=

‖(I-Ω)y‖2

(6)

上面分析的過(guò)程將多參數(shù)估計(jì)轉(zhuǎn)化為單一的頻率參數(shù)的估計(jì)，然而計(jì)算最小范數(shù)解相對(duì)復(fù)雜，但可以基于此得到一個(gè)很重要的結(jié)論，就是當(dāng)諧波個(gè)數(shù)K小于實(shí)際個(gè)數(shù)時(shí)， 已測(cè)出的部分諧波參數(shù)仍是可靠的，這是由正交理論決定的。

J3(fk)=Min‖yk-αkω(fk)‖2

(7)

求解該優(yōu)化式，可得

Maxfk|ωH(fk)yk|

(8)

(9)

由此可以看出，fk的估計(jì)可以利用FFT方法進(jìn)行計(jì)算，且幅值最大處的頻率，即為所估計(jì)頻率值，而該頻率處FFT變換幅值與估計(jì)幅值成比例。實(shí)際應(yīng)用RELAX方法計(jì)算時(shí)，在計(jì)算出某一對(duì)參數(shù)后，應(yīng)利用優(yōu)化式J1判斷參數(shù)估計(jì)是否合適，如果未達(dá)到誤差要求則循環(huán)迭代，直到參數(shù)估計(jì)的足夠好后再計(jì)算后面的參數(shù)。這樣每增加一對(duì)參數(shù)的估計(jì)，就需對(duì)前面已估計(jì)出的參數(shù)進(jìn)行更新。

需要說(shuō)明的是上面介紹的RELAX算法不能直接用來(lái)分析實(shí)信號(hào)，但通過(guò)以下步驟修正后可以準(zhǔn)確檢測(cè)諧波(間諧波)信號(hào)參數(shù)。

步驟1X(t)=x(t)+jH[x(t)]

在這樣反差顯著且形式緊迫的大環(huán)境下， 結(jié)合計(jì)算機(jī)信息和3D圖形處理分析技術(shù)以及互聯(lián)網(wǎng)技術(shù)的發(fā)展契機(jī)，BIM應(yīng)運(yùn)而生, 并迅速引起行業(yè)內(nèi)領(lǐng)袖企業(yè)的關(guān)注和重視。 BIM作為一個(gè)完整的概念被明確提出可以回溯到20世紀(jì)90年代 [4]。而其真正興起是在新千年后，在軟件行業(yè)巨頭Autodesk， Bentley Systems和Graphisoft等公司的共同推動(dòng)下，一系列以數(shù)字化建造為核心概念的新一代軟件產(chǎn)品(例如AECOsim以及稍后占北美市場(chǎng)主導(dǎo)地位的Autodesk Revit)被推出，以取代傳統(tǒng)CAD產(chǎn)品作為行業(yè)變革的技術(shù)基礎(chǔ)并支持建筑行業(yè)內(nèi)數(shù)字信息交互和協(xié)同合作。

(10)

步驟2φ=φ+π/2

(11)

(注：H[·]表示Hilbert變換，且Hilbert變換在負(fù)頻處有π/2的相位延遲)

至此，整個(gè)算法流程如下：

(1)估計(jì)主諧波個(gè)數(shù)K。

(2)初始化αk及fk，對(duì)觀測(cè)數(shù)據(jù)進(jìn)行HILBERT變換，取變換后復(fù)值觀測(cè)數(shù)據(jù)作待測(cè)數(shù)據(jù)。

(3)計(jì)算待估計(jì)分量yk。

(5)根據(jù)優(yōu)化目標(biāo)J1判斷參數(shù)估計(jì)是否準(zhǔn)確，滿足條件則依次估計(jì)其它分量參數(shù)否則重復(fù)估計(jì)已測(cè)參數(shù)。

(6)對(duì)估計(jì)相位進(jìn)行修正。

2 仿真實(shí)驗(yàn)

構(gòu)造如下含間諧波的虛擬電壓信號(hào)，

V(t)=10sin(2π50t+π/4)+

2sin(2π73t+π/3)+

4sin(2π77t+π/5)+

6sin(2π121t+π/6)+e(t)

(12)

其中加入了信噪比為30dB的高斯白噪聲。下圖分別是觀測(cè)信號(hào)波形及使用FFT對(duì)觀測(cè)信號(hào)進(jìn)行頻域分析結(jié)果(頻域圖中幅值已歸一化)。

為使用RELAX方法，需先確定K值?，F(xiàn)定義譜峰分離判別式

(13)

當(dāng)兩個(gè)頻率處的譜峰滿足上述關(guān)系，即可認(rèn)為被正確識(shí)別。本例中經(jīng)計(jì)算K為3，這是由于73 Hz與77Hz頻率接近而無(wú)法分辨。為測(cè)試算法容錯(cuò)能力，實(shí)驗(yàn)中分別取K為3、4、5、6，測(cè)試算法的。圖3表示RELAX算法在不同K值下的分析結(jié)果，表1是具體的誤差分析。

圖1 混噪待測(cè)信號(hào)Fig.1 Signal mixed with noise

圖2 FFT頻域分析Fig.2 Frequency analysis via FFT表1 不同K值下RELAX算法誤差分析Tab.1 Error analysis of RELAX algorithm in different case of K

諧波數(shù)K實(shí)測(cè)頻率/Hz誤差/%實(shí)測(cè)幅值/V誤差/%實(shí)測(cè)相位/rad誤差/%349.980.049.891.10.754.5176.770.033.8830.596.1120.860.126.1220.512.6449.890.229.910.90.790.5972.560.62.073.51.050.2777.120.164.1230.630.27120.520.45.921.330.512.6549.790.429.782.20.772.9672.650.481.8861.031.6476.890.143.892.80.612.92120.390.51.8768.830.2356.07120.270.62.0366.170.6320.32649.850.39.594.11.2356.6172.730.371.12441.022.676.490.661.3765.750.2166.58120.390.55.783.670.488.3349.890.220.9890.20.2864.3576.520.621.2369.255.21729.2

結(jié)合圖3和表1可以看出，當(dāng)諧波估計(jì)數(shù)K小于實(shí)際諧波個(gè)數(shù)時(shí)，RELAX算法仍可準(zhǔn)確檢測(cè)諧波參數(shù)。當(dāng)K大于實(shí)際個(gè)數(shù)時(shí)，會(huì)產(chǎn)生某些頻率諧波的重復(fù)測(cè)量，并導(dǎo)致這些頻率處諧波參數(shù)估計(jì)不準(zhǔn)確。因此只需要將這些頻率諧波合并，并相應(yīng)減少K值重復(fù)計(jì)算便可將參數(shù)準(zhǔn)確測(cè)量。

圖3 RELAX在不同K值下頻率估計(jì)Fig.3 Frequency estimation via RELAX by different K

3 結(jié)語(yǔ)

RELAX算法的計(jì)算主要使用FFT，因此在實(shí)時(shí)性上得到保障，經(jīng)實(shí)際測(cè)試雖然計(jì)算精度不如PRONY及SVD-TLS方法，但測(cè)量時(shí)間遠(yuǎn)小于二者。與傳統(tǒng)FFT方法相比較，該方法不存在頻譜泄露問(wèn)題，可準(zhǔn)確計(jì)算頻率接近的諧波參數(shù)，因此可以應(yīng)用于間諧波的實(shí)時(shí)測(cè)量分析中。

[1] 王波, 楊洪耕(Wang Bo, Yang Honggeng).基于AR譜估計(jì)和插值FFT的間諧波檢測(cè)方法(A method for Interharmonic measurement based on AR model spectral estimation and interpolation FFT algorithm)[J].繼電器(Relay), 2006, 34(4):49-52,56.

[2] 張伏生, 耿中行, 葛耀中(Zhang Fusheng,Geng Zhongxing,Ge Yaozhong).電力系統(tǒng)諧波分析的高精度FFT算法(FFT algorithm with high accuracy for harmonic analysis in power system)[J].中國(guó)電機(jī)工程學(xué)報(bào)(Proceedings of the CSEE), 1999 ,19(3):63-66.

[3] 易立強(qiáng), 鄺繼順(Yi Liqiang, Kuang Jishun). 一種基于FFT的實(shí)時(shí)諧波分析算法(FFT-based algorithm for real-time harmonic analysis)[J].電力系統(tǒng)及其自動(dòng)化學(xué)報(bào)(Proceedings of the CSU-EPSA), 2007, 19(2):98-102.

[4] 薛惠, 楊仁剛(Xue Hui, Yang Rengang).利用Morlet連續(xù)小波變換實(shí)現(xiàn)非整次諧波的檢測(cè)(Morlet wavelet based detection of noninteger harmonics)[J].電網(wǎng)技術(shù)(Power System Technology), 2002, 26(12):41-44.

[5] 占勇, 丁屹峰(Zhan Yong, Ding Yifeng). 電力系統(tǒng)諧波分析的穩(wěn)健支持向量機(jī)方法研究(A robust support vector algorithm for harmonics analysis of electric power system)[J].中國(guó)電機(jī)工程學(xué)報(bào)(Proceedings of the CSEE), 2004, 24(12):43-47.

[6] 石敏, 吳正國(guó), 伊為民(Shi Min, Wu Zhengguo, Yi Weimin). 基于多信號(hào)分類(lèi)法和普羅尼法的間諧波參數(shù)估計(jì)(Inter-harmonic parameter estimation based on multi-signal classification and Prony method)[J].電網(wǎng)技術(shù)(Power System Technology), 2005, 29(15):81-84.

[7] Ling J, Stoica P, Li J，etal.On using cyclic algorithms for sinusoidal parameter estimation[J].Electronics Letters, 2008, 44(19):1160-1161.

[8] Li Jian , Stoica Petre. Efficient mixed-spectrum estimation with applications to target feature extraction[J].IEEE Trans on Signal Processing, 1996, 44(2):281-295.

[9] IEC 61000-4-30,Testing and measurement techniques-Power quality measurements methods[S].

ApplicationofRELAXAlgorithmonInterharmonicDetection

JI Yu1， SUN Yun-lian1， LI Jing2

(1.School of Electrical Engineering, Wuhan University, Wuhan 430072, China；2.Wuhuan Chemical Co.,Ltd., Wuhan 430072, China)

Interharmonic is the special harmonic element which is not an integer multiple of the fundamental. Due to the asynchronous sampling problem, serious leakage error may arise with the measurement method for normal harmonics. In this paper, a new approach for sinusoidal parameter estimation based on RELAX(Relaxation algorithm)is presented. To get the correct result, the proposed method depends on two additional steps: the Hilbert Transform of the samples and the phase revise. The method can make the faster measurement in comparison with others since the Fast Fourier Transform(FFT) takes up the main calculation time of the RELAX. Finally, the simulation results have verified the effectiveness and practicability of the algorithm.

interharmonic； asynchronous sampling； relaxation algorithm； Hilbert transform； fast Fourier transform

2009-11-09

2010-01-18

TM93

A

1003-8930(2011)03-0135-04

季 宇(1982-)，男，博士，研究方向?yàn)殡娔苜|(zhì)量分析與控制。Email:jiyu820507@gmail.com

孫云蓮(1962-)，女，博士生導(dǎo)師，研究方向?yàn)殡娔苜|(zhì)量分析與控制，分布式發(fā)電并網(wǎng)控制技術(shù)。Email:ylsun@whu.edu.cn

李 晶(1980-)，女，碩士，研究方向?yàn)榛ぴO(shè)計(jì)中的配電技術(shù)。Email:lijingdk@cwcec.com