Dolph-Chebyshev 窗插值FFT 的諧波參數(shù)估計(jì)

田文博，余健明，馬小津，朱 博

（西安理工大學(xué)自動(dòng)化學(xué)院，西安710048）

隨著高壓直流輸電系統(tǒng)快速發(fā)展以及非線性負(fù)荷在電力系統(tǒng)中的廣泛應(yīng)用，電力系統(tǒng)諧波和間諧波（包括次諧波）污染日益嚴(yán)重[1，2]，高精度諧波檢測(cè)手段對(duì)電力系統(tǒng)諧波治理有重要意義。

電力系統(tǒng)的諧波測(cè)量常常采用快速傅里葉變換FFT 實(shí)現(xiàn)。但是，由于電力系統(tǒng)的頻率在額定工頻左右波動(dòng)，這樣就無法確保這個(gè)時(shí)刻變化的頻率是采樣頻率分辨率的正整數(shù)倍，從而無法同步采樣，就產(chǎn)生了柵欄效應(yīng)和頻譜泄漏現(xiàn)象[3～5]，嚴(yán)重影響了檢測(cè)的精度。針對(duì)FFT 算法的缺點(diǎn)，國內(nèi)外學(xué)者提出了一系列加窗插值FFT 算法。常見的加窗插值FFT 算法有加矩形窗[6]、Bartlett、Hanning窗[7，8]、Blackman 窗[9，10]、Blackman-Harris 窗[11]和Rife-Vincent 窗[12，13]等插值FFT 算法。采用多項(xiàng)式逼近的有效形式計(jì)算頻率修正系數(shù)和振幅修正系數(shù)[14]，計(jì)算量小，實(shí)時(shí)性強(qiáng)，提高了加Blackman-Harris 窗插值FFT 算法的計(jì)算速度。

通常選擇窗函數(shù)應(yīng)使其頻譜：①主瓣寬度盡量小，以使過渡帶盡量陡；②旁瓣相對(duì)于主瓣越小越好，即能量盡可能集中于主瓣內(nèi)[15]。由于矩形窗，Bartlett 窗和廣義余弦窗的頻譜具有固定的主瓣寬度和固定的最大旁瓣電平，在對(duì)弱信號(hào)分量的提取時(shí)往往誤差很大。

Dolph-Chebyshev 是一種局部?jī)?yōu)化的時(shí)窗函數(shù)，它滿足窗函數(shù)的最大振幅比準(zhǔn)則，即在規(guī)定旁瓣電平的條件下，主瓣寬度最小。其主瓣峰值與旁瓣峰值之比可根據(jù)實(shí)際需要自由選擇，其可以通過調(diào)節(jié)Dolph-Chebyshev 形狀參數(shù)，減小旁瓣電平，從而減小不同諧波分量之間的泄漏的干擾，提高估算精度，尤其是對(duì)弱信號(hào)分量更為明顯。

本文討論Dolph-Chebyshev 窗的頻譜特性，提出了基于Dolph-Chebyshev 窗插值FFT 算法諧波參數(shù)估計(jì)方法，建立諧波參數(shù)估計(jì)的數(shù)學(xué)模型，推導(dǎo)出了信號(hào)各分量頻率、幅值和初相角的計(jì)算公式。利用該算法通過對(duì)弱信號(hào)分量的參數(shù)估計(jì)，證明該方法的優(yōu)越性。

1 Dolph-Chebyshev 窗

在給定旁瓣高度下，Dolph-Chebyshev 窗的主瓣寬度最小，具有等波動(dòng)性，也就是說，其所有的旁瓣都具有相等的高度，比通常的時(shí)窗函數(shù)具有明顯的優(yōu)點(diǎn)。其物理意義是切比雪夫多項(xiàng)式在單位圓上做N 點(diǎn)等間隔抽樣，然后再作離散傅里葉反變換得到的。Dolph-Chebyshev 窗的數(shù)學(xué)定義如下[16]。

時(shí)域表示為

離散頻譜為

式中：N 為窗函數(shù)的長(zhǎng)度；Rp為形狀參數(shù)；ΔF 由Rp決定的中間變量。參數(shù)N、Rp、ΔF、M 之間的關(guān)系為

圖1（a）給出了當(dāng)Rp=0.312 5 固定不變（此時(shí)旁瓣峰值約為-200 dB），N = [129，257，513]時(shí)，Dolph-Chebyshev 窗的歸一化對(duì)數(shù)頻譜。可看到當(dāng)變化時(shí)該窗函數(shù)旁瓣電平為恒定不變。

圖1 Dolph-Chebyshev 窗頻譜Fig.1 Spectra of Dolph-Chebyshev windows

圖1（b） 給出了當(dāng)N = 129 固定不變，Rp=[0.457 2，0.351 2，0.307 6]時(shí)Dolph-Chebshev 窗的歸一化對(duì)數(shù)頻譜。形狀參數(shù)Rp越小，該窗函數(shù)旁瓣電平越小即旁瓣相對(duì)于主瓣越小。因此在進(jìn)行插值FFT 時(shí)，通過選擇合適的形狀參數(shù)Rp，增大窗譜旁瓣衰減，減小泄漏影響，提高精度。

利用拉格朗日插值法得該窗最大旁瓣電平γ與形狀參數(shù)Rp的數(shù)學(xué)關(guān)系為

當(dāng)用式（3）估計(jì)最大衰減電平時(shí)，絕對(duì)誤差不超過1.536 0 dB。利用拉格朗日插值法得該窗主瓣寬度w 與最大旁瓣電平γ 的關(guān)系為

當(dāng)用式（4）估計(jì)最大衰減電平時(shí)，絕對(duì)誤差不超過0.418 2。

實(shí)際應(yīng)用中，確定所要達(dá)到的第一旁瓣相對(duì)主瓣衰減的分貝數(shù)γ，然后根據(jù)式（3）確定窗函數(shù)的形狀參數(shù)，從而在已知N 情況下確定窗函數(shù)。并通過式（4）估計(jì)此時(shí)主瓣的寬度w。

2 加Dolph-Chebyshev 窗插值FFT 算法

含諧波分量的信號(hào)x（t）以采樣頻率fs均勻采樣獲得的離散時(shí)間信號(hào)為

式中：H 為諧波的次數(shù)；Ah為第h 次諧波的幅值；fh為第h 次諧波的頻率值；φh為第h 次諧波初相角。

對(duì)x（n）加Dolph-Chebyshev 窗，得到

對(duì)xw（n）做DTFT 變換得

對(duì)上式進(jìn)行頻率離散化，并考慮到當(dāng)估計(jì)第h次諧波分量的參數(shù)時(shí)，f 處于fh附近的頻譜，故忽略負(fù)頻率點(diǎn)-fh對(duì)第h 次諧波影響和其他分量對(duì)第h 次諧波的能量泄漏的影響，并可近似得

式中：離散頻率間隔Δf=fs/N，N 是數(shù)據(jù)截?cái)嗟拈L(zhǎng)度；W[·]為Dolph-Chebyshev 窗的連續(xù)頻譜。

對(duì)頻率為fh信號(hào)，若fs不是fh正整數(shù)倍，即不能對(duì)該頻率信號(hào)的同步采樣，則第h 次諧波的實(shí)際頻率fh=kΔf 將不能正好處于Xw（kΔf）的離散頻點(diǎn)上，此時(shí)k 不為整數(shù)。對(duì)此，設(shè)峰值點(diǎn)k 左右兩根譜線分別是第k1和k2條譜線，這兩條譜線應(yīng)該是k 附近的幅值最大和次最大譜線，顯然有k1≤k≤k2（其中k2=k1+1）。

令y1和y2為k1和k2處對(duì)應(yīng)譜線幅值，即y1=|Xw（k1Δf）|，y2=|Xw（k2Δf）|，引入?yún)?shù)δ=k1-k（其中0≤δ≤1），則有

將式（9）用Dolph-Chebshev 窗離散頻譜Wc[·]表示得

將式（2）代入式（10）得

式中，M、ΔF 分別為確定的常數(shù)。

式（9）為δ 的一元超越方程，可利用Matlab 平臺(tái)中fzero（）函數(shù)來求出此次方程在閉區(qū)間[0，1]內(nèi)的正解。然后根據(jù)

求得第h 次諧波頻率。利用下式求得第h 次諧波幅值和初相角分別為

3 仿真分析

3.1 與其他窗插值FFT 諧波參數(shù)估計(jì)算法比較

采用文獻(xiàn)[6，8，17]給出的測(cè)試信號(hào)模型進(jìn)行仿真實(shí)驗(yàn)分析。

式中：A0=0.2 V；A1=6 V；f1=20.2 Hz；φ1=0.1°；A3=1 V；f3=60.6 Hz；φ3=0°。

Dolph-Chebshev 窗的形狀參數(shù)取Rp=0.351 2（相當(dāng)窗譜旁瓣衰減200 dB）。由于Dolph-Chebshev窗的長(zhǎng)度必須是奇數(shù)，這里將Dolph-Chebshev 窗的窗長(zhǎng)度取2 049，fs=1 000 Hz。表1 列出了采用加矩形窗[6]、Hanning 窗[8]、Kaiser 窗[17]的插值FFT 算法與本文提出的算法對(duì)比，各項(xiàng)數(shù)據(jù)均以絕對(duì)誤差列出。表中將Dolph-Cheybshev 簡(jiǎn)稱Chebyshev。

從表1 中可以看出，Dolph-Chebshev 窗具有良好的頻譜泄漏抑制能力，與經(jīng)典窗算法比較，證明算法的正確性，并且全部計(jì)算結(jié)果精度獲得了大幅度的提高。

表1 不同加窗插值FFT 算法仿真結(jié)果對(duì)比Tab.1 Comparison of simulation results from different windowed interpolation FFT algorithms

3.2 對(duì)弱信號(hào)分量的參數(shù)估計(jì)

在某個(gè)幅值較小的諧波信號(hào)（弱信號(hào)分量，如2 次諧波）附近存在一個(gè)幅值較大的諧波（如基波）信號(hào)時(shí)，較大幅值信號(hào)的頻譜泄漏將導(dǎo)致弱信號(hào)分量被干擾甚至淹沒，弱信號(hào)分量的檢測(cè)誤差往往很大。利用文獻(xiàn)[18，19]給出的信號(hào)進(jìn)行測(cè)試，該信號(hào)的基波及諧波分量參數(shù)如表2 所示。這里將文獻(xiàn)[18，19]的計(jì)算結(jié)果一并與本文提出的加Dolph-Chebshev 窗插值FFT 算法做比較。采樣頻率fs=1 500 Hz，數(shù)據(jù)截?cái)嚅L(zhǎng)度為N=512。這里Dolph-Chebshev 窗的形狀參數(shù)取Rp=0.351 2（相當(dāng)窗譜旁瓣衰減200 dB）。（文獻(xiàn)[18]僅給出了基波頻率和各次諧波幅值仿真結(jié)果，而文獻(xiàn)[19]未給出各次諧波頻率，故表2 只給出部分比較結(jié)果）。

表2 仿真信號(hào)的諧波成分Tab.2 Harmonic components of the simulated signal

仿真結(jié)果如表3 所示，其中Ef、EA、Eφ是本文提出的算法分析得到的頻率、幅值、初相角的相對(duì)誤差；EA1、Eφ1是文獻(xiàn)[19]提出的算法分析得到的幅值、初相角的相對(duì)誤差；EA2是文獻(xiàn)[18]提出的算法幅值的相對(duì)誤差。

表3 不同算法下的計(jì)算結(jié)果比較Tab.3 Comparisons of the proposed algorithm with others%

采用本文算法的基波頻率誤差為3.1×10-9%，優(yōu)于文獻(xiàn)[18]最優(yōu)誤差2×10-5%，略遜于文獻(xiàn)[19]的1.8×10-9%；從表3 中可以看出，本文方法對(duì)初相角檢測(cè)的多數(shù)結(jié)果優(yōu)于文獻(xiàn)[19]。對(duì)于幅值的結(jié)果則是優(yōu)于文獻(xiàn)[19]4 個(gè)數(shù)量級(jí)，優(yōu)于文獻(xiàn)[18]5個(gè)數(shù)量級(jí)。總體上，對(duì)弱分量參數(shù)估計(jì)精度比文獻(xiàn)[18，19]有一定的提高。

尤其是在基波附近，受到基波影響最為強(qiáng)烈的2 次諧波，檢測(cè)精度都高于文獻(xiàn)[18，19]。

文獻(xiàn)[18]最優(yōu)結(jié)果是采用雙峰譜線修正加Blackman 窗的實(shí)驗(yàn)結(jié)果，Blackman 最大旁瓣峰值電平為-57 dB。而文獻(xiàn)[19]采用的Rife-Vincent（I）窗，旁瓣峰值電平為-74.5 dB。而本采用的Dolph-Chebshev 窗通過調(diào)節(jié)形狀參數(shù)，使最大旁瓣衰減為-200 dB，極大地抑制了基波頻譜泄漏對(duì)鄰近弱分量的影響。從上述仿真實(shí)驗(yàn)結(jié)果來看，證明了這一點(diǎn)。

4 結(jié)語

Dolph-Chebshev 窗是一種具有最大振幅比的時(shí)窗函數(shù)，其頻譜旁瓣衰減的分貝數(shù)可以通過調(diào)整形狀參數(shù)進(jìn)行自由選擇。選擇合適的形狀參數(shù)，使窗譜具有較低的旁瓣電平，這時(shí)Dolph-Chebshev窗具有良好的頻率泄漏抑制特性。文章將該窗應(yīng)用于插值FFT 諧波參數(shù)估計(jì)，提出了基于Dolph-Chebyshev 窗插值FFT 算法的諧波參數(shù)估計(jì)方法，與已有的加窗FFT 進(jìn)行比較，能有效地消除各次諧波相互干擾，具有很高的計(jì)算精度，對(duì)于信號(hào)中弱分量的參數(shù)估計(jì)效果尤為明顯，適合于高準(zhǔn)確度諧波分析。

[1]Ortmeyer T H，Chakravarthi K R，Mahmoud A A. The effects of power system harmonics on power system equipment and loads [J]. IEEE Trans on Power Apparatus and Systems，1985，104（9）：2555-2563.

[2]林海雪（Lin Haixue）.電力系統(tǒng)中的間諧波問題（Interharmonics in electrical power system）[J].供用電（Distribution&Utilization），2001，18（3）：6-9.

[3]李庚銀，陳志業(yè)，寧宇（Li Gengyin，Chen Zhiye，Ning Yu）.快速傅里葉變換的兩種改進(jìn)算法（Two fast Fourier transform algorithms）[J].電力系統(tǒng)自動(dòng)化（Automation of Electric Power Systems），1997，21（12）：37-40.

[4]Harris F J. On the use of windows for harmonic analysis with the discrete Fourier transform[J]. Proceedings of the IEEE，1978，66（1）：51-83.

[5]楊洪耕，惠錦，侯鵬（Yang Honggeng，Hui Jin，Hou Peng）.電力系統(tǒng)諧波和間諧波檢測(cè)方法綜述（Detection methods of harmonics and inter-harmonics in power system）[J]. 電力系統(tǒng)及其自動(dòng)化學(xué)報(bào)（Proceedings of the CSU-EPSA），2010，22（2）：65-69.

[6]Jain V K，Collins W L Jr，Davis D C.High-accuracy analog measurements via interpolated FFT[J]. IEEE Trans on Instrumentation and Measurement，1979，28 （2）：113-122.

[7]柴旭崢，文習(xí)山，關(guān)根志，等（Chai Xuzheng，Wen Xishan，Guan Genzhi，et al）. 一種高精度的電力系統(tǒng)諧波分析算法（An algorithm with high accuracy for analysis of power system harmonics）[J]. 中國電機(jī)工程學(xué)報(bào)（Proceedings of the CSEE），2003，23（9）：67-70.

[8]Grandke T. Interpolation algorithms for discrete Fourier transform of weighted signals[J]. IEEE Trans on Instrumentation and Measurement，1983，32（2）：350-355.

[9]周俊，王小海，祁才君，等（Zhou Jun，Wang Xiaohai，Qi Caijun，et al）. 基于Blackman 窗函數(shù)的插值FFT 在電網(wǎng)諧波信號(hào)分析中的應(yīng)用（Estimation of electrical harmonic parameters by using the interpolated FFT algorithm based on Blackman window）[J].浙江大學(xué)學(xué)報(bào)：理學(xué)版（Journal of Zhejiang University Science Edition），2006，33（6）：650-653.

[10]許珉，楊陽，章夢(mèng)哲，等（Xu Min，Yang Yang，Zhang Mengzhe，et al）.一種加三項(xiàng)余弦窗的加窗插值FFT 算法（An interpolated FFT algorithm based on three-item cosine window）[J].電力系統(tǒng)保護(hù)與控制（Power System Protection and Control），2010，38（17）：11-15.

[11]許珉，劉凌波（Xu Min，Liu Lingbo）.基于三次樣條函數(shù)的加Blackman - harris 窗插值FFT 算法（Blackmanharris window interpolated FFT algorithm based on cubic spline function）[J].電力自動(dòng)化設(shè)備（Electric Power Automation Equipment），2009，29（2）：59-63.

[12]曾博，滕召勝，高云鵬，等（Zeng Bo，Teng Zhaosheng，Gao Yunpeng，et al）. 基于Rife-Vincent 窗的高準(zhǔn)確度電力諧波相量計(jì)算方法（An accurate approach for power harmonic phasor calculation based on Rife-Vincent window）[J]. 電工技術(shù)學(xué)報(bào)（Transactions of China Electrotechnical Society），2009，24（8）：154-159.

[13]張文強(qiáng)，楊耀民，許珉，等（Zhang Wenqiang，Yang Yaomin，Xu Min，et al）.基于三次樣條函數(shù)的加Rife—Vincent（Ⅲ）窗FFT 插值算法（The Rife—Vincent（Ⅲ）window interpolation FFT algorithm by using cubic spline function）[J].電力系統(tǒng)保護(hù)與控制（Power System Protection and Control），2009，37（12）：36-39.

[14]蔣春芳，劉敏（Jiang Chunfang，Liu Min）.基于雙插值FFT 算法的間諧波分析（Inter-harmonics analysis based on double interpolation FFT algorithm）[J].電力系統(tǒng)保護(hù)與控制（Power System Protection and Control），2010，38（3）：11-14，19.

[15]潘文，錢俞壽，周鶚（Pan Wen，Qian Yushou，Zhou E）.基于加窗插值FFT 的電力諧波測(cè)量理論——（Ⅰ）窗函數(shù)研究（Power harmonics measurement based on windows and interpolated FFT（I）Study of windows）[J].電工技術(shù)學(xué)報(bào)（Transactions of China Electrotechnical Society），1994，24（1）：50-54.

[16]《數(shù)學(xué)辭海》編輯委員會(huì).數(shù)學(xué)辭海[M].5 卷.北京：中國科學(xué)技術(shù)出版社，2002.

[17]高云鵬，滕召勝，溫和，等（Gao Yunpeng，Teng Zhaosheng，Wen He，et al）. 凱塞窗插值FFT 的電力諧波分析與應(yīng)用（Harmonic analysis based on Kaiser window interpolation FFT and its application）[J].中國電機(jī)工程學(xué)報(bào)（Proceedings of the CSEE），2010，30（4）：43-48.

[18]龐浩，李東霞，俎云霄，等（Pang Hao，Li Dongxia，Zu Yunxiao，et al）.應(yīng)用FFT 進(jìn)行電力系統(tǒng)諧波分析的改進(jìn)算法（An improved algorithm for harmonic analysis of power system using FFT technique）[J].中國電機(jī)工程學(xué)報(bào)（Proceedings of the CSEE），2003，23（6）：50-54.

[19]曾博，滕召勝，溫和，等（Zeng Bo，Teng Zhaosheng，Wen He，et al）.萊夫-文森特窗插值FFT 諧波分析方法（An approach for harmonic analysis based on Rife-Vincent window interpolation FFT）[J]. 中國電機(jī)工程學(xué)報(bào)（Proceedings of the CSEE），2009，29（10）：115-120.