Dolph-Chebyshev 窗插值FFT 的谐波参数估计

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

（西安理工大学自动化学院，西安710048）

随着高压直流输电系统快速发展以及非线性负荷在电力系统中的广泛应用，电力系统谐波和间谐波（包括次谐波）污染日益严重[1，2]，高精度谐波检测手段对电力系统谐波治理有重要意义。

电力系统的谐波测量常常采用快速傅里叶变换FFT 实现。但是，由于电力系统的频率在额定工频左右波动，这样就无法确保这个时刻变化的频率是采样频率分辨率的正整数倍，从而无法同步采样，就产生了栅栏效应和频谱泄漏现象[3～5]，严重影响了检测的精度。针对FFT 算法的缺点，国内外学者提出了一系列加窗插值FFT 算法。常见的加窗插值FFT 算法有加矩形窗[6]、Bartlett、Hanning窗[7，8]、Blackman 窗[9，10]、Blackman-Harris 窗[11]和Rife-Vincent 窗[12，13]等插值FFT 算法。采用多项式逼近的有效形式计算频率修正系数和振幅修正系数[14]，计算量小，实时性强，提高了加Blackman-Harris 窗插值FFT 算法的计算速度。

通常选择窗函数应使其频谱：①主瓣宽度尽量小，以使过渡带尽量陡；②旁瓣相对于主瓣越小越好，即能量尽可能集中于主瓣内[15]。由于矩形窗，Bartlett 窗和广义余弦窗的频谱具有固定的主瓣宽度和固定的最大旁瓣电平，在对弱信号分量的提取时往往误差很大。

Dolph-Chebyshev 是一种局部优化的时窗函数，它满足窗函数的最大振幅比准则，即在规定旁瓣电平的条件下，主瓣宽度最小。其主瓣峰值与旁瓣峰值之比可根据实际需要自由选择，其可以通过调节Dolph-Chebyshev 形状参数，减小旁瓣电平，从而减小不同谐波分量之间的泄漏的干扰，提高估算精度，尤其是对弱信号分量更为明显。

本文讨论Dolph-Chebyshev 窗的频谱特性，提出了基于Dolph-Chebyshev 窗插值FFT 算法谐波参数估计方法，建立谐波参数估计的数学模型，推导出了信号各分量频率、幅值和初相角的计算公式。利用该算法通过对弱信号分量的参数估计，证明该方法的优越性。

1 Dolph-Chebyshev 窗

在给定旁瓣高度下，Dolph-Chebyshev 窗的主瓣宽度最小，具有等波动性，也就是说，其所有的旁瓣都具有相等的高度，比通常的时窗函数具有明显的优点。其物理意义是切比雪夫多项式在单位圆上做N 点等间隔抽样，然后再作离散傅里叶反变换得到的。Dolph-Chebyshev 窗的数学定义如下[16]。

时域表示为

离散频谱为

式中：N 为窗函数的长度；Rp为形状参数；ΔF 由Rp决定的中间变量。参数N、Rp、ΔF、M 之间的关系为

图1（a）给出了当Rp=0.312 5 固定不变（此时旁瓣峰值约为-200 dB），N = [129，257，513]时，Dolph-Chebyshev 窗的归一化对数频谱。可看到当变化时该窗函数旁瓣电平为恒定不变。

图1 Dolph-Chebyshev 窗频谱Fig.1 Spectra of Dolph-Chebyshev windows

图1（b） 给出了当N = 129 固定不变，Rp=[0.457 2，0.351 2，0.307 6]时Dolph-Chebshev 窗的归一化对数频谱。形状参数Rp越小，该窗函数旁瓣电平越小即旁瓣相对于主瓣越小。因此在进行插值FFT 时，通过选择合适的形状参数Rp，增大窗谱旁瓣衰减，减小泄漏影响，提高精度。

利用拉格朗日插值法得该窗最大旁瓣电平γ与形状参数Rp的数学关系为

当用式（3）估计最大衰减电平时，绝对误差不超过1.536 0 dB。利用拉格朗日插值法得该窗主瓣宽度w 与最大旁瓣电平γ 的关系为

当用式（4）估计最大衰减电平时，绝对误差不超过0.418 2。

实际应用中，确定所要达到的第一旁瓣相对主瓣衰减的分贝数γ，然后根据式（3）确定窗函数的形状参数，从而在已知N 情况下确定窗函数。并通过式（4）估计此时主瓣的宽度w。

2 加Dolph-Chebyshev 窗插值FFT 算法

含谐波分量的信号x（t）以采样频率fs均匀采样获得的离散时间信号为

式中：H 为谐波的次数；Ah为第h 次谐波的幅值；fh为第h 次谐波的频率值；φh为第h 次谐波初相角。

对x（n）加Dolph-Chebyshev 窗，得到

对xw（n）做DTFT 变换得

对上式进行频率离散化，并考虑到当估计第h次谐波分量的参数时，f 处于fh附近的频谱，故忽略负频率点-fh对第h 次谐波影响和其他分量对第h 次谐波的能量泄漏的影响，并可近似得

式中：离散频率间隔Δf=fs/N，N 是数据截断的长度；W[·]为Dolph-Chebyshev 窗的连续频谱。

对频率为fh信号，若fs不是fh正整数倍，即不能对该频率信号的同步采样，则第h 次谐波的实际频率fh=kΔf 将不能正好处于Xw（kΔf）的离散频点上，此时k 不为整数。对此，设峰值点k 左右两根谱线分别是第k1和k2条谱线，这两条谱线应该是k 附近的幅值最大和次最大谱线，显然有k1≤k≤k2（其中k2=k1+1）。

令y1和y2为k1和k2处对应谱线幅值，即y1=|Xw（k1Δf）|，y2=|Xw（k2Δf）|，引入参数δ=k1-k（其中0≤δ≤1），则有

将式（9）用Dolph-Chebshev 窗离散频谱Wc[·]表示得

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

式中，M、ΔF 分别为确定的常数。

式（9）为δ 的一元超越方程，可利用Matlab 平台中fzero（）函数来求出此次方程在闭区间[0，1]内的正解。然后根据

求得第h 次谐波频率。利用下式求得第h 次谐波幅值和初相角分别为

3 仿真分析

3.1 与其他窗插值FFT 谐波参数估计算法比较

采用文献[6，8，17]给出的测试信号模型进行仿真实验分析。

式中：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 窗的形状参数取Rp=0.351 2（相当窗谱旁瓣衰减200 dB）。由于Dolph-Chebshev窗的长度必须是奇数，这里将Dolph-Chebshev 窗的窗长度取2 049，fs=1 000 Hz。表1 列出了采用加矩形窗[6]、Hanning 窗[8]、Kaiser 窗[17]的插值FFT 算法与本文提出的算法对比，各项数据均以绝对误差列出。表中将Dolph-Cheybshev 简称Chebyshev。

从表1 中可以看出，Dolph-Chebshev 窗具有良好的频谱泄漏抑制能力，与经典窗算法比较，证明算法的正确性，并且全部计算结果精度获得了大幅度的提高。

表1 不同加窗插值FFT 算法仿真结果对比Tab.1 Comparison of simulation results from different windowed interpolation FFT algorithms

3.2 对弱信号分量的参数估计

在某个幅值较小的谐波信号（弱信号分量，如2 次谐波）附近存在一个幅值较大的谐波（如基波）信号时，较大幅值信号的频谱泄漏将导致弱信号分量被干扰甚至淹没，弱信号分量的检测误差往往很大。利用文献[18，19]给出的信号进行测试，该信号的基波及谐波分量参数如表2 所示。这里将文献[18，19]的计算结果一并与本文提出的加Dolph-Chebshev 窗插值FFT 算法做比较。采样频率fs=1 500 Hz，数据截断长度为N=512。这里Dolph-Chebshev 窗的形状参数取Rp=0.351 2（相当窗谱旁瓣衰减200 dB）。（文献[18]仅给出了基波频率和各次谐波幅值仿真结果，而文献[19]未给出各次谐波频率，故表2 只给出部分比较结果）。

表2 仿真信号的谐波成分Tab.2 Harmonic components of the simulated signal

仿真结果如表3 所示，其中Ef、EA、Eφ是本文提出的算法分析得到的频率、幅值、初相角的相对误差；EA1、Eφ1是文献[19]提出的算法分析得到的幅值、初相角的相对误差；EA2是文献[18]提出的算法幅值的相对误差。

表3 不同算法下的计算结果比较Tab.3 Comparisons of the proposed algorithm with others%

采用本文算法的基波频率误差为3.1×10-9%，优于文献[18]最优误差2×10-5%，略逊于文献[19]的1.8×10-9%；从表3 中可以看出，本文方法对初相角检测的多数结果优于文献[19]。对于幅值的结果则是优于文献[19]4 个数量级，优于文献[18]5个数量级。总体上，对弱分量参数估计精度比文献[18，19]有一定的提高。

尤其是在基波附近，受到基波影响最为强烈的2 次谐波，检测精度都高于文献[18，19]。

文献[18]最优结果是采用双峰谱线修正加Blackman 窗的实验结果，Blackman 最大旁瓣峰值电平为-57 dB。而文献[19]采用的Rife-Vincent（I）窗，旁瓣峰值电平为-74.5 dB。而本采用的Dolph-Chebshev 窗通过调节形状参数，使最大旁瓣衰减为-200 dB，极大地抑制了基波频谱泄漏对邻近弱分量的影响。从上述仿真实验结果来看，证明了这一点。

4 结语

Dolph-Chebshev 窗是一种具有最大振幅比的时窗函数，其频谱旁瓣衰减的分贝数可以通过调整形状参数进行自由选择。选择合适的形状参数，使窗谱具有较低的旁瓣电平，这时Dolph-Chebshev窗具有良好的频率泄漏抑制特性。文章将该窗应用于插值FFT 谐波参数估计，提出了基于Dolph-Chebyshev 窗插值FFT 算法的谐波参数估计方法，与已有的加窗FFT 进行比较，能有效地消除各次谐波相互干扰，具有很高的计算精度，对于信号中弱分量的参数估计效果尤为明显，适合于高准确度谐波分析。

[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）.电力系统中的间谐波问题（Interharmonics in electrical power system）[J].供用电（Distribution&Utilization），2001，18（3）：6-9.

[3]李庚银，陈志业，宁宇（Li Gengyin，Chen Zhiye，Ning Yu）.快速傅里叶变换的两种改进算法（Two fast Fourier transform algorithms）[J].电力系统自动化（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）.电力系统谐波和间谐波检测方法综述（Detection methods of harmonics and inter-harmonics in power system）[J]. 电力系统及其自动化学报（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]柴旭峥，文习山，关根志，等（Chai Xuzheng，Wen Xishan，Guan Genzhi，et al）. 一种高精度的电力系统谐波分析算法（An algorithm with high accuracy for analysis of power system harmonics）[J]. 中国电机工程学报（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 窗函数的插值FFT 在电网谐波信号分析中的应用（Estimation of electrical harmonic parameters by using the interpolated FFT algorithm based on Blackman window）[J].浙江大学学报：理学版（Journal of Zhejiang University Science Edition），2006，33（6）：650-653.

[10]许珉，杨阳，章梦哲，等（Xu Min，Yang Yang，Zhang Mengzhe，et al）.一种加三项余弦窗的加窗插值FFT 算法（An interpolated FFT algorithm based on three-item cosine window）[J].电力系统保护与控制（Power System Protection and Control），2010，38（17）：11-15.

[11]许珉，刘凌波（Xu Min，Liu Lingbo）.基于三次样条函数的加Blackman - harris 窗插值FFT 算法（Blackmanharris window interpolated FFT algorithm based on cubic spline function）[J].电力自动化设备（Electric Power Automation Equipment），2009，29（2）：59-63.

[12]曾博，滕召胜，高云鹏，等（Zeng Bo，Teng Zhaosheng，Gao Yunpeng，et al）. 基于Rife-Vincent 窗的高准确度电力谐波相量计算方法（An accurate approach for power harmonic phasor calculation based on Rife-Vincent window）[J]. 电工技术学报（Transactions of China Electrotechnical Society），2009，24（8）：154-159.

[13]张文强，杨耀民，许珉，等（Zhang Wenqiang，Yang Yaomin，Xu Min，et al）.基于三次样条函数的加Rife—Vincent（Ⅲ）窗FFT 插值算法（The Rife—Vincent（Ⅲ）window interpolation FFT algorithm by using cubic spline function）[J].电力系统保护与控制（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].电力系统保护与控制（Power System Protection and Control），2010，38（3）：11-14，19.

[15]潘文，钱俞寿，周鹗（Pan Wen，Qian Yushou，Zhou E）.基于加窗插值FFT 的电力谐波测量理论——（Ⅰ）窗函数研究（Power harmonics measurement based on windows and interpolated FFT（I）Study of windows）[J].电工技术学报（Transactions of China Electrotechnical Society），1994，24（1）：50-54.

[16]《数学辞海》编辑委员会.数学辞海[M].5 卷.北京：中国科学技术出版社，2002.

[17]高云鹏，滕召胜，温和，等（Gao Yunpeng，Teng Zhaosheng，Wen He，et al）. 凯塞窗插值FFT 的电力谐波分析与应用（Harmonic analysis based on Kaiser window interpolation FFT and its application）[J].中国电机工程学报（Proceedings of the CSEE），2010，30（4）：43-48.

[18]庞浩，李东霞，俎云霄，等（Pang Hao，Li Dongxia，Zu Yunxiao，et al）.应用FFT 进行电力系统谐波分析的改进算法（An improved algorithm for harmonic analysis of power system using FFT technique）[J].中国电机工程学报（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]. 中国电机工程学报（Proceedings of the CSEE），2009，29（10）：115-120.