一种加窗插值FFT谐波分析方法

一种加窗插值FFT谐波分析方法
摘要：由于很难实现同步采样和整周期截断，因此，利用fft算法分析电网谐波信号时存在频谱泄露和栅栏效应，影响算法的分析精度。

加窗插值fft是抑制频谱泄露和消除栅栏效应的有效方法，在此提出一种基于3项3阶nuttall窗插值fft的谐波分析方法，推导了插值系数公式以及各次谐波的频率、幅值和相位的修正公式。

对该算法与hanning窗、blackman窗插值fft算法进行matlab 仿真对比研究，验证了该算法具有更高的分析精度。

关键词：谐波； fft；窗函数；插值；电力系统
引言
随着大量电力电子装置和非线性负载在电力系统中的广泛应用，使电网中产生了大量的高次谐波［12］，严重威胁电网的电能质量和用户设备的安全运行，因此谐波的准确测量具有重要意义。

快速傅里叶变换(fft)是最主要的电力谐波分析方法，但由于实际工程应用中很难实现同步采样和整周期截断，因此，fft方法存在频谱泄露和栅栏效应，影响谐波分析精度。

加窗插值fft算法是抑制频谱泄漏和消除栅栏效应的有效方法，其原理是通过加窗运算抑制频谱泄露、通过插值运算消除栅栏效应。

常用的窗函数有hanning窗［35］、blackman窗［6］、blackman harris窗［7 8］等。

hanning窗的插值公式简单，计算量小，但分析精度较低；blackman和blackman harris窗插值fft算法的分析精度高，
但插值公式复杂，计算量大。

本文首先分析了fft算法频谱泄露的原因和3项3阶nuttall窗函数的特点，然后推导出其插值fft算法的计算公式，其插值系数具有简单的显式表达式，谐波的频率、幅值和相位的修正公式简单明了，易于实现。

同hanning窗和blackman窗插值fft算法的仿真对比研究结果表明，所提出算法更加适合于电力系统谐波的精确测量。

1fft频谱泄漏的原因
离散傅里叶变换(dft)变换是针对有限长序列信号进行傅里叶变换的一种数值分析方法。

为分析简便，设单一频率信号表达式
为:xm(t)=amej(2πfmt+φm)(1)式中am，fm，φm分别为幅值、频率和相位。

以采样频率fs对其进行均匀采样n点，得到离散序列x(n)：
x(n)=xm(n)·wr(n)，n=0,1,2,…,n-1(2)
式中wr(n)为矩形窗。

应用dft变换对采样信号x(n)进行分析时，隐含在时域上对其进行周期延拓。

在理想同步采样条件下，x(n)周期延拓后的序列与原连续信号x(t)的采样序列完全相同，如图1所示。

此时，dft变换能够精确分析原连续信号x(t)的谐波参数。

非同步采样时，x(n)周期延拓的序列不再是原连续信号x(t)的采样序列，如图2所示。

此时，dft变换将会产生频谱泄漏和栅栏效应，不能精确分析原连
续信号x(t)的谐波参数。

非同步采样序列x(n)周期延拓后在边界处产生的采样点跳变是产生频谱泄漏和栅栏效应的根本原因。

图1同步采样时的波形拓展图2非同步采样时的波形拓展根据调制定理，信号x(n)经离散时间傅里叶变换（dtft）后的频谱序列为［9］：x(λ)=amwr(λ－λm)ej［－n－1nπ(λ－λm)+φm］ (3) wr(λ)=sin(λπ)/sin(λπ/n) (4)式中：λ=nf/fs为整数；λm=n·fm/fs；wr(λ)是矩形窗的幅度谱。

设λm=km+δm，其中km为正整数，0≤δm＜1，则非同步采样时的幅度频谱如图3所示，可见，此时整数位置上的频谱x(km)与谐波的真实频谱x(λm)不一致，即发生了栅栏效应。

为消除dft算法的频谱泄漏和栅栏效应，需要选择合适的窗函数。

2nuttall窗插值fft算法
nuttall窗是一类余弦组合窗［10］，其时域表示为：w(n)=∑m
－1m=0(－1)mamcos(2πn·m/n)(5)式中：m为窗函数的项数，
n=0,1,2,…,n-1。

当窗函数为3项3阶nuttall窗时：a0=0.375，a1=0.5，a2=0.125，其旁瓣衰减为47 db，旁瓣衰减速度为30 db，适合电力谐波的准确分析。

图3非同步采样时的幅度频谱余弦组合窗的dtft表达式为：w(ej ω)=∑m－1m=0(－1)mam2［wr(ej(ω－2πmn))+wr(ej(ω+2πmn))］
≈∑m－1m=0am2(wr(ω－2πmn)+wr(ω+2πmn))·e－jn－12ω=w(ω)e－jn－12ω (6)式中:wr(ejω)=wr(ω)e－jn－12ω是矩形窗的dtft频谱;wr(ω)=sin(nω/2)/sin(ω/2)是矩形窗的幅度谱;w(ω)为余弦组合窗的幅度谱。

信号x(n)加余弦组合窗的频谱序列为：x(λ)=amw(λ－λm)ej ［－n－1nπ(λ－λm)+φm］ (7)
w(λ)=∑m－1m=0am2［wr(λ－m)+wr(λ+m)］ (8)定义x(n)幅度谱线上的两个相邻峰值之比为：βm=|x(km+1)|/|x(km)|(9)由于n 比较大，而且0≤δm＜1，对于3项3阶nuttall窗，可得：βm=(2+δm)/(3－δm) (10)即：δm=(3βm－2)/(1+βm)(11)根据δm可以估计出谐波的频率、幅值和相位：fm=(km+δm)fs/n(12)
am=1n|x(km)|·πδm(1－δ2m)(4－δ2m)1.5sin(δmπ)(13) φm=angle［x(km)］－δmπn－1n(14)3仿真结果
为验证提出的3项3阶nuttall窗插值fft算法的有效性，将该算法与hanning窗及blackman窗插值fft算法进行matlab仿真对比研究。

设信号表达式为：x(t)=∑9m=1amcos(2πmf1t+φm)(15)式中各谐波的具体参数如表1～表3所示。

取采样频率为10 khz，采样点数2 048点，约10个基波周期。

采用hanning窗、blackman窗及提出的3项3阶nuttall窗插值fft算法得到的频率、幅值和相位分析结果如表1～表3所示。

由
仿真结果可见：hanning窗的分析精度最低，blackman窗次之，而3项3阶nuttall窗的分析精度最高。

4结语
本文分析了fft算法产生频谱泄露的根本原因和3项3阶nuttall 窗的旁瓣特性，并且推导了其插值fft算法的谐波频率、幅值和相位的计算公式。

提出的谐波分析方法具有插值修正公式简单，计算量小，分析精度高的优点。

仿真结果表明，该算法适合电力系统谐波的准确测量。

表3相位参数的仿真结果
谐波
次数设定值
/(°)估计值 /(°)hanningblackmannuttall12020.020 720.019 120.017 524042.030 341.186 039.900 635049.976 050.020 050.053 447069.977 370.030 970.064 05100100.052 8100.069 8100.086 76130130.060 4130.043 1130.014 07200200.005 1200.020 4200.035 58250249.903 4249.991 4250.058 29290290.061 4290.066 9290.069 6误差标准差0.637 60.363 20.052 0
参考文献
［1］李红，杨善水.电力系统谐波检测的现状与发展［j］.现代电子技术，2004，27(9):8184.
［2］李圣清，朱英浩，周有庆，等.电网谐波检测方法的综述［j］.高电压技术，2004，30(3):3942.
［3］祁才君，陈隆道，王小海.应用插值fft算法精确估计电网谐波参数［j］.浙江大学学报：工学版，2003,37(1):112116. ［4］熊杰锋，王柏林，孙艳.电力系统间谐波和谐波分析的海宁窗插值算法［j］.自动化仪表，2010，31(4):2526.
［5］高索丹，巴鹏.hanning窗在电力系统谐波分析中的应用［j］.自动化技术与仪表，2008，27(11):124127.
［6］周俊，王小海，祁才君.基于blackman窗函数的插值fft
在电网谐波信号分析中的应用［j］.浙江大学学报：理学版，2006,33(6):650653.
［7］zhang f, geng z, yuan w. the algorithm of interpolating windowed fft for harmonic analysis of electric power system \[j\]. ieee transactions on power delivery, 2001, 16(2): 160 164.
［8］赵文春，马伟明，胡安.电机测试中谐波分析的高精度fft 算法［j］.中国电机工程学报，2001,21(12):8387.
［9］grandke t. interpolation algorithms for discrete fourier transforms of weighted signals \[j\]. ieee transactions on instrument and measurement, 1983, 32(2): 350 355.
［10］nuttall a h. some windows with very good sidelobe behabior \[j\]. ieee transaction on acoustics, speech, and signal processing, 1981, assp29(1): 8491.。