基于非迭代多點插值DFT的諧波電能計量算法

唐旭明，翁東波，宋雅楠，郝晶晶，張 明，李 冰

(1.國網淮南供電公司，安徽 淮南 232007；2.煙臺東方威思頓電氣有限公司，山東 煙臺 264000)

0 引言

伴隨國民經濟快速發展的需求，現代化電氣高速鐵路、特高壓輸電等加快投入使用。一方面，這促進了固態開關器件、非線性電力電子開關負載、工業設備整流器和逆變器新型電力電子器件在電網系統中大量應用，使得系統控制、保護等更為準確、穩定。另一方面，非線性電力電子器件應用產生大量諧波(穩態或動態)[1-2]所帶來的電能質量問題，不僅危及供電系統整體安全、穩定運行，而且直接影響電能計量的準確性，進而造成巨大的經濟損失。因此，在這種情況下，精確諧波電能計量顯得尤為必要和重要。這已成為電力公用事業及其客戶的重要研究課題。

在非正弦含諧波情況下，諧波標準IEEE STD 1459[3]有功功率定義為：

(1)

式(1)表明，準確估計測量直流、基波及各次諧波參數(幅值、相位及頻率)是確保諧波電能精確計量的前提。常用諧波估計算法分為時域和頻域兩大類。由于頻域方法具有高精度和較低計算負擔，在處理實時應用時它們通常比時域方法更受青睞。

插值離散傅里葉變換(interpolated discrete Fourier transform,IpDFT)是較廣泛采用的頻域方法。該類算法采用加窗和插值的方法，分別處理非相干采樣引起的頻譜泄漏，以及由頻譜中離散頻率引起的柵欄效應引起的估計不準確影響，并基于快速傅里葉變換(fast Fourier transform,FFT)使得算法實現簡單且可用于實時處理。典型IpDFT算法有非參數插值(補零、線性調頻和頻率翹曲[3-4])和特定時間窗函數方法。其中，通過剖析頻譜泄漏現象的本質，采用加余弦窗(包括多級卷積窗[5-6])的線性插值及多點譜線加權插值(multipoint weighted interpolated discrete fourier transform,MWIpDFT)方法是目前研究的重要方向[7-11]。該類方法通常分粗略估計和精細補償校正兩步實現。粗略估計通過基于頻率估計器對采樣信號使用加窗FFT(windowing FFT,WFFT)并搜索幅度譜全局最大值實現；精細補償校正通常通過內插與最大譜線相鄰譜線來實現，以補償校正柵欄效應導致的偏移誤差。

在插值補償校正時，由于偏移補償量與DFT譜線幅度譜的非線性關系通常采用多項式近似逼近，因而在一定程度上降低了補償精度。文獻[12]提出提升插值離散傅里葉變換(enhanced interpolated discrete Fourier transform,EIpDFT)算法。其核心在于采用迭代過程消除了與共軛分量及其他頻次諧波對參數估計結果的影響，大幅度提高了頻率校正及幅值相位估計精度。但迭代過程的高額計算成本使得其在實時計算應用存在困難。

因此，本文在EIpDFT算法的基礎上，采用非迭代方式擴展應用基于任意階最大衰減旁瓣窗(maximum decay sidelobe windows,MSDW)的多點插值方法，在短時間窗長內對多頻諧波信號的基本正弦分量進行參數估計。最后，根據諧波標準IEEE STD 1459功率定義，將得到的諧波參數用于諧波電能準確計量。一方面，非迭代方式有效降低了每次參數估計的計算負擔。另一方面，短時間窗截斷使得該類算法還可應用于快速波動情況下電網信號諧波參數估計，擴大了IpDFT的應用范圍。

1 通用加窗傅里葉變換

假設信號中僅含有單一頻率分量，其離散時間采樣表示為：

(2)

式中：An、fin和θin分別為信號幅值、頻率和起始相位；fs、N為采樣頻率和采樣點數。

實際情況下，考慮非同采樣情況，信號頻率可以表示為：

(3)

式中：Δf為N點采樣下的頻率分辨率；η為信號截斷周期數；l(整數)為η的四舍五入值；τ∈[-0.5,0.5]為信號頻率在離散譜線間的所處位置，即歸一化情況的頻率偏差。

對N點離散信號進行加窗w(n)傅里葉變換，得：

λ-η)ejθin+W(l+λ+η)e-jθin

(4)

式中：λ用于確定在插值中使用的多條譜線位置[13]；W()為窗函數的離散頻譜表達式。

根據H階最大衰減旁瓣余弦窗的定義：

(5)

對應權重系數分別為：

(6)

(7)

則當H>1、ζ=N及N>1時，H階最大衰減旁瓣余弦窗的頻譜表達式近似為：

(8)

根據式(8)，H≤6的窗函數系數及重要性能指標，包括主瓣寬度(MLBW)、旁瓣衰減速度(SLdec,單位dB/oct)及等效噪聲帶寬(ENBW)。最大衰減旁瓣余弦窗參數如表1所示。

表1 最大衰減旁瓣余弦窗參數

2 插值校正

多點插值校正方法在前述加窗DFT基礎上，通常分兩步實現：粗略估計和精細補償校正。以常用三譜線插值修正算法為例[15]。首先，利用局部峰值搜索(周期圖)，在等間隔離散頻譜上尋找實際峰值頻率點附近最大譜線，將其作為粗略頻率估計值。其次，聯立峰值譜線左右兩個次最大譜線幅值，采用非線性擬合方式得到頻率修正系數，完成精細補償校正過程(若需進一步計算幅值和相位，則可再次通過非線性擬合方式得到頻率修正系數與幅值及相位的關系，將其代入得到幅值和相位的補償校正)。

值得注意的有兩點：①非線性擬合方式精度與其采用的計算方法有著直接的關系，同時也額外增加了計算成本；②在利用三譜線非線性擬合求頻率修正的過程中，常用方法往往忽略了共軛負頻率部分旁瓣泄漏的影響(如文獻[5])。插值校正流程如圖1所示。

圖1 插值校正流程圖

為解決上述兩個問題，文獻[16-18]提出了采用迭代的提升插值離散傅里葉變換算法，對應流程如圖1所示。該方法優點在于充分考慮了共軛負頻率旁瓣衰減對校正方法的影響，即圖1中迭代過程窗譜函數完整包含了負頻率部分，對應公式為：

(9)

基于上述特點，迭代算法可以在很短信號周期內提供準確的諧波信號參數估計。但迭代過程的高額計算成本，使得其應用于實時計算時存在困難。

3 非迭代插值算法

①頻率估計。

假設信號頻率處于其對應離散傅里葉變換的第kp根譜線附近，聯合峰值譜線左右兩根次最大譜線(kp±1)，進而根據加窗傅里葉變換式(4)，可得：

(10)

式(10)可進一步化簡為：

AB=0

(11)

式中：0為零矩陣。

(12)

由于幅值和相位組成的變量矩陣不可能為0，因此方程系數矩陣的行列式ΔA等于0，即：

(13)

進一步根據式(8)，可得下述比例關系：

(14)

將式(14)代入式(13)，利用行列式性質，化簡為：

(15)

通過對上述行列式分解，可得關于η的代數方程為：

Γ1-η2Γ2=0

(16)

(17)

則可得頻率估計為：

(18)

②幅值和相位估計。

假設：

(19)

則根據式(4)可得：

(20)

同理，聯合左右兩根次最大譜線(kp±1)，以矩陣方式重寫式(20)，可得：

CD=E

(21)

矩陣C、D和E分別為：

(22)

D=[ab]T

(23)

E=[XR(kp)XR(kp)]T

(24)

聯立式(21)～式(24)，采用最小二乘法，可得：

(25)

則幅值和相位分別為：

(26)

基于非迭代散點插值校正可精確得到單一頻率信號參數估計。通過前述推導，上述估計方法可進一步擴展應用于諧波情況下各頻次參數估計，具體差異在于式(11)和式(20)所對應的系數矩陣構造。

4 諧波電能計量

根據前述得到的被測信號基波和各次諧波分量的幅值、相位和頻率，可完成基波和諧波電能精確計算。假設電網電壓電流信號分別為：

(27)

根據IEC STD 61000-4-7，在k個基波周期T內，單一頻率分量的標準諧波電能公式為：

(28)

則在計及基波及各次諧波正反向電能情況下，可得總電能為：

(29)

基于非迭代插值校正諧波電能計量流程如圖2所示。

圖2 基于非迭代插值校正諧波電能計量流程圖

5 算法仿真分析

5.1 窗譜近似誤差

窗頻譜及對應近似誤差如圖3所示。

圖3 窗頻譜及對應近似誤差圖

在前述頻率估計過程中，為滿足方程解析過程需要，根據文獻[19]近似得到最大衰減旁瓣余弦窗的頻譜，因此需要分析近似可能帶來的估計誤差。其中，根據H階最大衰減旁瓣余弦窗定義，其準確頻譜根據離散時間傅里葉變換得到，為：

(30)

窗頻譜近似首先會造成式(4)頻域離散采樣誤差。而從圖3可明顯看出，在主瓣內W(ζ)與W″(ζ)近似誤差非常小，基本在10-9范圍以下，表明以式(17)和式(18)連續三頻域離散采樣點線性組合的頻率估計是非常精確的。同時，用于方程求解的式(14)線性變換準確度也可以得到保證。

5.2 單一頻率信號參數估計

假設單一頻率離散采樣信號為：

s(n)=Aincos(2πfinnT+θin)

(31)

此處參數設置分別為：幅值Ain=5;頻率fin=49.52 Hz;初始相位θin=π/3;采樣周期T=1/5 000 s。

上述信號充分考慮了實際電網信號頻率波動范圍，因而具有典型性。分別采用本文設計的非迭代三譜線插值和典型加漢寧窗三點權重插值方法。從結果可以明顯看出，本文設計的方法較加漢寧窗三點權重插值方法精度高出約4個數量級，整體相對誤差零級在10-6左右；值得特別注意的是在短時截斷情況下，非迭代方法也可以準確得到信號參數，因而可應用于短時波動信號分析。

5.3 諧波電能計量

假設含諧波離散采樣信號為：

(32)

諧波分量參數如表2所示。基波頻率為49.8 Hz，采樣頻率為5 000 Hz。

表2 諧波分量參數

首先，以電壓作為分析對象，分別采用H=2、3、4、5所對應的4個不同離散時間窗，得到的不同階窗函數下周期數與估計誤差如圖4所示。

圖4 不同階窗函數下周期數與估計誤差

最后以5個基波周期截斷，分別計算電壓、電流各次諧波幅值、相位，得到功率及誤差，如表3所示。

表3 諧波分量誤差

6 結束語

常用多點線性插值方法未考慮負頻率共軛部分帶來的長泄漏干擾，對于高精度諧波電能計量，會帶來累積誤差。本文提出的基于非迭代多點插值DFT精確諧波電能計量算法，充分考慮并消除了與共軛分量相關的參數估計結果的影響，利用多階最大衰減旁瓣窗函數高精度近似，進一步消除近似誤差對窗譜函數的影響。而采用穩健的三點線性組合方式，有效避免了幅值和相位波動對頻率估計的影響，理論分析及仿真結果表明，該方法在短時窗情況也具有很高的計算精度，因而基于該非迭代插值方法的電能計量算法理論的可用于動態負荷諧波電能計量，較以往其他計量算法具有更廣的適用性。

鑒于篇幅有限，本文未進一步就算法抗噪特性進行分析和仿真，因而下一步將重點研究算法在含噪聲情況下性能評估和提升；同時，還將該算法進一步推廣應用于非平穩信號參數估計，以解決復雜工況下諧波電能計量問題。