莫旭陽,錢承山,王志偉,陶海燕

(南京信息工程大學信息與控制學院,南京 210044)

?

基于諧波測量的改進余弦窗函數(shù)*

莫旭陽,錢承山*,王志偉,陶海燕

(南京信息工程大學信息與控制學院,南京 210044)

摘要:電網(wǎng)電壓監(jiān)測設備在監(jiān)測運行設備的電流電壓波形時通常采用諧波分析法,針對傳統(tǒng)FFT諧波分析法在電網(wǎng)頻率波動時諧波測量會產(chǎn)生柵欄效應和頻譜泄漏的問題,在研究余弦窗特性的基礎上提出了一種改進3項余弦窗。通過頻譜衰減曲線與傳統(tǒng)余弦窗對比,采用基于數(shù)據(jù)樣本的非線性擬合方法,擬合出了改進窗插值算法的諧波幅值相位修正系數(shù),使改進窗插值算法能實際應用于諧波測量。通過諧波測量仿真對比,驗證了加改進余弦窗的諧波測量準確度更高。

關(guān)鍵詞:諧波測量;FFT分析法;頻譜泄漏;余弦窗;加窗插值算法

電網(wǎng)監(jiān)測設備是保證電網(wǎng)安全穩(wěn)定運行的重要電子設備,電網(wǎng)監(jiān)測設備一般通過監(jiān)測運行設備的諧波電流來判斷電網(wǎng)設備的故障和老化。傳統(tǒng)諧波測量采用FFT快速傅里葉變換法,在電網(wǎng)頻率穩(wěn)定時諧波幅值和頻率測量能滿足測量要求。然而電網(wǎng)頻率具有時變性,在電網(wǎng)頻率波動時,非同步測量會產(chǎn)生柵欄效應和頻譜泄漏[1],采用傳統(tǒng)傅里葉變換測量諧波會產(chǎn)生較大偏差,降低諧波測量的準確性。

柵欄效應和頻譜泄漏主要是由于測量設備的采樣頻率和電網(wǎng)頻率非同步造成的,非同步采樣下信號的截取長度與信號周期之比是非整數(shù)關(guān)系,造成諧波分析傅里葉變換后的頻譜線與信號頻率不重合從而產(chǎn)生偏差。對采樣信號進行加窗插值是目前降低柵欄效應和頻譜泄漏的主要方法[2]。傳統(tǒng)加窗算法主要采用余弦窗,經(jīng)典余弦窗有漢寧窗、海明窗、布萊克曼窗以及B-H窗[3]。

本文在研究余弦窗特性的基礎上提出了一種基于布萊克曼窗的3項余弦改進窗,該窗較傳統(tǒng)余弦窗在旁瓣峰值和旁瓣衰減速度均有提升,通過仿真驗證了改進窗的有效性,并通過非線性擬合方法擬合出了改進窗與插值算法結(jié)合后的幅值相位修正系數(shù),完善了窗函數(shù),增加了改進窗的實用性。

1　余弦窗函數(shù)的特性

余弦窗是定義在[-L/2,L/2]關(guān)于X=0軸對稱的余弦函數(shù),其在時域中的表達式為:

(1)

對式(1)進行歸一化后余弦窗可化為:

(2)

其中

(3)

取窗函數(shù)w(t)在t=L/2處的左右極限得:

由函數(shù)的連續(xù)性和對稱性可知,若

(4)

則窗函數(shù)在t=±L/2處連續(xù),同時在整個定義域中連續(xù)。

對式(2)窗函數(shù)w(t)分別求1階導數(shù)、2階導數(shù)后可得:

同理可證明1階導數(shù)w′(t)在t=±L/2處連續(xù),同時在整個定義域中連續(xù)。當式(5)成立時,w″(t)在t=±L/2處連續(xù),同時在整個定義域中連續(xù)。

(5)

式(3)～式(5)為余弦窗及其1階2階導數(shù)連續(xù)的約束條件[5],其中式(3)為余弦窗的歸一化特性。

理想余弦窗應具有較小的旁瓣和較快的旁瓣衰減速度,從而有效抑制信號發(fā)生頻譜泄漏[6-7]。文獻[5]指出余弦窗的旁瓣衰減速度及其最大旁瓣大小與其窗函數(shù)及其導數(shù)在定義域內(nèi)的連續(xù)性有關(guān)。對于2項余弦窗而言,當余弦窗在定義域內(nèi)連續(xù),即余弦窗系數(shù)ak同時滿足式(3)和式(4)時,傅里葉變換后的余弦窗旁瓣衰減速度為-18dB/Octave,最大旁瓣為-32dB,特征余弦窗為漢寧窗,2項系數(shù)為a0=0.5,a1=0.5,頻譜衰減曲線如圖1所示。當余弦窗在定義域內(nèi)不連續(xù),即滿足式(3),不滿足式(4)時,旁瓣衰減速度降為-6dB/Octave,最大旁瓣減小為-43dB,特征余弦窗為海明窗,2項系數(shù)為a0=0.54,a1=0.46,頻譜衰減曲線如圖2所示。對于3項余弦窗而言,當余弦窗在定義域內(nèi)連續(xù),二階導數(shù)不連續(xù),即系數(shù)ak滿足式(3)、式(4)不滿足式(5)時,旁瓣衰減速度為-18dB/Octave,最大旁瓣為-58dB。特征窗余弦窗為布萊克曼窗,3項系數(shù)為a0=0.42,a1=0.5,a2=0.08,頻譜衰減曲線如圖3所示。當余弦窗2階導數(shù)在定義域內(nèi)連續(xù),即系數(shù)ak同時滿足式(3)～式(5)時,旁瓣衰減速度達到3項余弦窗極大值,平均衰減達到-30dB/Octave,同時最大旁瓣大小達到3項余弦窗極大值,大小為-46dB,特征余弦窗為exact-布萊克曼窗,3項系數(shù)為a0=0.375,a1=0.5,a2=0.125,頻譜衰減曲線如圖4所示。由此推論當3項余弦窗系數(shù)滿足式(3)、式(4),不滿足式(5)時,可以獲得較小的旁瓣和較快的旁瓣衰減速度。

圖1　漢寧窗幅頻衰減曲線

圖2　海明窗幅頻衰減曲線

圖3　布萊克曼窗幅頻衰減曲線

圖4　exact布萊德曼窗幅頻衰減曲線

2　改進的3項余弦窗函數(shù)

2項漢寧窗和3項布萊克曼窗由于窗項數(shù)少,計算量小,旁瓣衰減速度快被廣泛應用于實際測量中,本文提出了基于3項布萊克曼窗的改進余弦窗,窗函數(shù)系數(shù)滿足前文所述的式(3)、式(4),同時不滿足式(5),窗函數(shù)時域表達式為:

t≤|L/2|

(6)

其離散表達式和頻譜函數(shù)分別為式(7)、式(8),其中N為采樣點數(shù)

n=0,1,…,N-1

(7)

(8)

其中

(9)

為了便于觀察改進余弦窗頻譜的衰減特性,定義了式(10)的頻譜衰減函數(shù),并將頻率進行歸一化處理

(10)

通過MATLAB軟件編程得到改進余弦窗的頻譜衰減函數(shù)曲線,如圖5所示,可以看出改進窗的頻譜衰減明顯優(yōu)于其他傳統(tǒng)余弦窗,改進窗的最大旁瓣達到-64dB,優(yōu)于最大旁瓣衰減為-32dB的漢寧窗和最大旁瓣衰減為-58dB的布萊克曼窗;改進窗旁瓣衰減速度為-18dB/Octave,與漢寧窗和布萊克曼窗基本持平。相比窗系數(shù)滿足式(5)的exact布萊克曼窗,雖然其旁瓣衰減速度達到了-30dB/Octave,但是其最大旁瓣只有-46dB的衰減,遠低于-64dB的改進窗,改進窗具有更小的最大旁瓣和較快的旁瓣衰減速度,旁瓣對主瓣影響更小,能更好抑制頻譜泄漏。

圖5　改進窗幅頻衰減曲線

3　基于改進余弦窗插值公式修正

由于電網(wǎng)頻率具有時變性,實際測量中很難做到同步采樣,非同步采樣會使得實際諧波頻率難以落在頻譜線上,使頻譜產(chǎn)生柵欄效應,插值算法能通過插值修正減小柵欄效應偏差,還原真實頻率,插值算法包括單峰插值和雙峰插值[8]。單峰插值算法是搜索最大峰值,雙峰插值算法則搜索最大和次大2個峰值。單峰插值易受到頻譜泄漏和噪音干擾,本文采用雙峰插值對改進的余弦窗進行插值修正。

假設存在一個幅值為A,頻率為f0,初始相位為φ待測諧波信號,以fs采樣率對其進行采樣N點后得到如下離散信號[9]:

(11)

對其所加時域表達式為w(n)的窗函數(shù)后進行傅里葉變換得到如下形式的離散信號:

(12)

由于負頻點-f0離正頻點f0較遠,其頻峰處的旁瓣影響可以忽略,則正頻點f0處的頻譜函數(shù)表達式為:

(13)

對上式頻率進行歸一化處理后頻譜函數(shù)變換為:

(14)

式中λ0是信號頻率f0以Δf=fs/N歸一化后對應的值,實際測量中信號頻率f0很難正好位于頻率分辨點上,即f0/Δf不是整數(shù),f0=k0Δf(k0不為整數(shù))。根據(jù)雙峰插值算法取信號f0對應的頻率點k0左右兩側(cè)譜線,分別記為第k1和第k2條譜線,顯然這兩條譜線為k0附近的最大和次大譜線,設兩條譜線幅值分別為y1=|Xλ(k1)|和y2=|Xλ(k2)|,k0k1k2滿足k1≤k0≤k2=k1+1。

引入輔助參數(shù)α[6],令α=k0-k1-0.5,

可推出

f0=k0Δf=(k1+α+0.5)·fs/N

(15)

(16)

(17)

將本文提出的改進窗對應的頻譜函數(shù)W(n)代入上式可得β關(guān)于α的函數(shù),記為β=g(α),其對應反函數(shù)記為α=g-1(β)。因為窗函數(shù)|W(n)|是偶對稱的,所以函數(shù)g(α)及其反函數(shù)g-1(β)都是奇函數(shù)。為了計算方便可將反函數(shù)逼近成以下奇次多項式形式:

α=a0β+a3β3+a5β5+a7β7

然而對于3項以上的余弦窗,其對應的函數(shù)g(α)很復雜,難以直接獲得其對應的反函數(shù)g-1(β)。本文采用了基于數(shù)據(jù)樣本的非線性擬合方法,通過MATLAB編程對設定好的諧波信號進行加窗和FFT變換得到變換后的數(shù)據(jù)樣本,所加窗為本文提出的3項改進余弦窗,通過數(shù)據(jù)樣本計算出各次諧波對應的α,β值,再利用MATLAB非線性擬合函數(shù)nlinfit()擬合出了奇函數(shù)α(β)對應多項式系數(shù),擬合后的α(β)多項式為:

α=1.9β+14.1β3-429.7β5+4191.6β7

雙峰插值后信號幅值修正采用了對k1和k2兩根譜線的幅值進行加權(quán)平均的方法,其幅值修正計算公式為:

(18)

當采樣點數(shù)N較大時,式(18)可以進一步化簡為:

A=N-1·2(y1+y2)·v(α)

(19)

其中v(α)為偶函數(shù),可以采用多項式逼近偶函數(shù)v(α),逼近多項式可表示為:

v(α)=b0+b2α2+b4α4+b6α6

(20)

由于α也為關(guān)于β的變量函數(shù),難以直接利用多項式進行函數(shù)逼近,同樣采用基于數(shù)據(jù)樣本的非線性擬合,通過樣本數(shù)據(jù)進行非線性擬合,擬合后的v(α)函數(shù)多項式為:

v(α)=2.7583+1.2513α2-2.3784α4

+8.0707α6

(21)

將式(21)代入式(19)即可得到插值后幅值修正式(21):

A=N-1·2(y1+y2)·(2.7583+1.2513α2

-2.3784α4+8.0707α6)

(22)

插值相位修正公式為式(23):

φ=arg[Xλ(k1)]+π/2-π[α+(-1)i·0.5]

(23)

經(jīng)過驗證擬合后的α值誤差為2%左右,屬于可以接受的誤差范圍。本文提出的基于數(shù)據(jù)樣本的方法可適應于其他諧波信號,數(shù)據(jù)樣本采用了所含諧波最高次數(shù)為7次的諧波信號,對于其他次數(shù)的諧波,修改對應的樣本信號函數(shù),通過本文方法可重新擬合出對應系數(shù),便于軟件實現(xiàn)。

4　仿真實驗

為了驗證所提出的改進3項余弦窗,本文通過MATLAB進行仿真研究,信號模擬了7次以內(nèi)的電網(wǎng)諧波信號。由于我國電網(wǎng)頻率波動范圍一般為(50±0.5)Hz[10],為了模擬電網(wǎng)頻率波動下加窗插值算法的作用,基波頻率f1取49.5Hz,采樣頻率fs為1 600Hz,采樣長度為1 024點,基波和各次諧波原始幅值和相位分別如表1和表2所示。

表1　諧波幅值仿真結(jié)果對比　單位:V

表2　諧波相位仿真結(jié)果對比　單位:(°)

仿真過程為先對信號進行采樣,對采樣后的數(shù)據(jù)進行加窗處理,然后進行快速傅里葉變換得到頻譜數(shù)組,通過雙峰插值算法,根據(jù)幅值和相位的修正式(22)、式(23)進行修正,最后得出修正后諧波的真實幅值和相位,并與傳統(tǒng)余弦窗進行了仿真對比,表1和表2分別為不加窗FFT測量、加漢寧窗、布萊克曼窗和改進余弦窗下各次諧波的幅值和相位仿真結(jié)果。DU和DΦ分別為測量結(jié)果與原始信號幅值和相位的絕對誤差。

上述仿真結(jié)果可以看出當電網(wǎng)頻率發(fā)生波動時,傳統(tǒng)FFT法存在柵欄效應和頻譜泄漏,直接測量諧波會產(chǎn)生很大偏差,所測量的結(jié)果無法反應真實諧波信號,采用本文提出的改進余弦窗結(jié)合雙峰插值算法能有效地減少柵欄效應和頻譜泄漏,準確提取諧波信號,對頻率波動具有較好的魯棒性。

5　結(jié)論

加窗插值法能有效提高頻率波動下諧波測量的準確性,本文提出的改進余弦窗具有更低的旁瓣和旁瓣衰減速度,能有效減小非同步測量帶來的柵欄效應和頻譜泄漏,結(jié)合雙峰插值算法可以有效應用于諧波測量中。本文提出的基于數(shù)據(jù)樣本的插值系數(shù)擬合方法相比傳統(tǒng)求解復雜反函數(shù)方法更易實現(xiàn),稍有不足的是該方法依賴于擬合函數(shù)的準確度,若能進一步提高擬合準確度,測量結(jié)果將更準確。

參考文獻:

[1]林青松,陳清華.諧波檢測中頻譜泄漏問題的研究[J].工礦自動化,2010(2):53-56.

[2]蔣慶斌,錢金法.繼電保護器中的諧波分析[J].電子器件,2012,35(5):562-566.

[3]郭桂香.改進的FFT在電能質(zhì)量監(jiān)測系統(tǒng)中的應用[J].自動化與儀器儀表,2011(2):94-95.

[4]許珉,楊陽,章夢哲,等.一種加三項余弦窗的加窗插值FFT算法[J].電力系統(tǒng)保護與控制,2010,38(17):11-14.

[5]Albert H.Nuttall Some Windows with Very Good Sidelobe Behavior[J].IEEE Transactions on Acoustics,Speech,and Signal Processing,1981,29(1):84-91.

[6]Qian Hao,Zhao Rongxiang,Chen Tong.Inter-Harmonics Analysis Based on Interpolating Windowed FFT Algorithm[J].IEEE Transactions on Power Delivery,2007,22(2):1064-1069.

[7]熊杰鋒,王柏林,孫艷.電力系統(tǒng)間諧波和諧波分析的海寧窗插值算法[J].自動化儀表,2010,31(4):25-27.

[8]王巍,陳慧慧,劉偉偉,等.基于二階海寧卷積窗和雙峰插值傅里葉變換的諧波檢測方法[J].工礦自動化,2010(10):49-52.

[9]牛勝鎖,梁志瑞,張建華.基于三譜線插值FFT的電力諧波分析算法[J].中國電機工程學報,2012,32(16):130-135.

[10]楊拴科,何衛(wèi)鋒.非同步采樣對電力系統(tǒng)諧波分析精度影響的仿真研究[J].電工技術(shù)雜志,2003(5):36-38.

莫旭陽(1988-),男,漢族,江蘇蘇州人,南京信息工程大學信息與控制學院碩士研究生,主要研究方向為避雷器在線監(jiān)測、嵌入式系統(tǒng)、物聯(lián)網(wǎng)應用,moxuyang163@163.com;

錢承山(1971-),男,漢族,山東泰安人,南京信息工程大學信息與控制學院教授,博士,碩士生導師,主要研究方向為非線性系統(tǒng)控制、自動檢測技術(shù)、智能終端與物聯(lián)網(wǎng)應用等,qianchengshan@163.com。

ImprovedCosineWindowFunctionBasedonHamonicMeasure*

MOXuyang,QIANChenshan*,WANGZhiwei,TAOHaiyan

(School of Information and Control,Nanjing University of Information Science and Technology,Nanjing 210044,China)

Abstract:Grid voltage monitoring equipments in the monitoring current and voltage waveforms of operating equipment commonly use the method of harmonic analysis.To solve the problem that traditional FFT harmonic analytical algorithm produces fence effect and spectrum leakage when the grid frequency shifts,an improved method of three cosine windows is proposed on the characteristic research of cosine window.Comparison to traditional windows through the spectrum attenuation curve shows superiority,a nonlinear fitting method based on the data sample is used to calculate correction factors for harmonic amplitude and frequency of the improved cosine window interpolated algorithm.With the corrected factors,the improved algorithm can be applied to real harmonic measure.Contrast to harmonic simulation results,the improved cosine window shows higher accuracy.

Key words:harmonic measure;FFT analytical algorithm;spectrum leakage;cosine window;windowed interpolation algorithm

doi:EEACC:815010.3969/j.issn.1005-9490.2014.04.042

中圖分類號:TM714

文獻標識碼:A

文章編號:1005-9490(2014)04-0777-05

收稿日期:2013-04-09修改日期:2013-12-28

項目來源:企事業(yè)委托項目(2013h066);南京信息工程大學科研啟動基金項目(20100307)