基于ZoomFFT和多谱线插值算法的谐波分析方法*

doi:10.11985/2023.04.020

基于ZoomFFT和多谱线插值算法的谐波分析方法^*

朱坤^,¹, 陈乐柱^,¹, 许胜^,², 史琛筠^,²

1.安徽工业大学电气与信息工程学院马鞍山 243000

2.泰州市电能变换与控制工程技术研究中心泰州 225300

Harmonic Analysis Method Based on ZoomFFT and Multispectral Interpolation Algorithm

ZHU Kun^,¹, CHEN Lezhu^,¹, XU Sheng^,², SHI Chenjun^,²

1. College of Electrical and Information Engineering, Anhui University of Technology, Maanshan 243000

2. Taizhou Electric Power Conversion and Control Engineering Technology Research Centre, Taizhou 225300

收稿日期: 2022-10-5 修回日期: 2023-07-26

基金资助:

*江苏省重点研发计划(BE2016727)
泰州科技计划(TS202035)

Received: 2022-10-5 Revised: 2023-07-26

作者简介 About authors

朱坤，男，1995年生，硕士研究生。主要研究方向为电能质量控制技术。E-mail：zhukun807@163.com

陈乐柱，男，1964年生，教授，硕士研究生导师。主要研究方向为电力电子技术的应用。E-mail：chenlezhu@163.com

许胜，男，1976年生，博士后，教授。主要研究方向为电能质量、微网交流器控制技术。E-mail：xush_00@sina.com

史琛筠，女，2001年生。主要研究方向为电能质量控制技术。E-mail：1003452681@qq.com

摘要

基于快速傅里叶变换(Fast Fourier transform，FFT)的各种检测算法，因其具备响应速度快以及便于工程实现等特点，现已成为电能质量分析常用方法，但易造成频谱泄漏，难以精确实现多个谐波的分离检测。因此提出一种基于Blackman-Harris窗函数四谱线插值算法和复调制频谱细化算法的复合型谐波检测方法。利用改进的加窗插值算法对谐波进行检测，对于发生主瓣干扰的频谱，采用复调制频谱细化算法，实现对谐波的精确分析。仿真试验结果表明，该方法能够分辨相近谐波分量，实时性较强，能够满足电力系统对电力参数检测的需求。

关键词： 谐波; Blackman-Harris窗; 四谱线插值; ZoomFFT; 主瓣干扰

Abstract

Various detection algorithms based on fast Fourier transform have become a common method for power quality analysis due to their fast response speed and easy engineering implementation. It is easy to cause spectrum leakage, and difficult to accurately realize the separation and detection of multiple harmonics. A complex harmonic detection method based on the Blackman-Harris window function quadratic line interpolation algorithm and the complex modulation spectrum refinement algorithm is proposed. The improved windowed interpolation algorithm is used to detect harmonics. For the spectrum where main lobe interference occurs, the complex modulation spectrum refinement algorithm is used to analyze the harmonics accurately. The simulation results show that the method can distinguish similar harmonic components, has strong real-time performance, and can meet the needs of power system detection of power parameters.

Keywords： Harmonic; Blackman-Harris window; four-spectral interpolation; ZoomFFT; main lobe interference

PDF (389KB) 元数据多维度评价相关文章导出 EndNote| Ris| Bibtex 收藏本文

本文引用格式

朱坤, 陈乐柱, 许胜, 史琛筠. 基于ZoomFFT和多谱线插值算法的谐波分析方法^*[J]. 电气工程学报, 2023, 18(4): 179-187 doi:10.11985/2023.04.020

ZHU Kun, CHEN Lezhu, XU Sheng, SHI Chenjun. Harmonic Analysis Method Based on ZoomFFT and Multispectral Interpolation Algorithm[J]. Chinese Journal of Electrical Engineering, 2023, 18(4): 179-187 doi:10.11985/2023.04.020

1 引言

随着电力电子行业的发展和非线性负荷设备的增加，电网谐波污染进一步恶化，从而影响电力系统安全平稳运行。为了能够有效治理电网中存在的电能质量问题，首先需要对电网的电力参数进行检测，为保证对电力系统的有效滤波或补偿，需要保证检测参数的精准性，同时为保证对电力系统的及时治理，还需要保证参数的实时性。

目前，基于快速傅里叶变换(Fast Fourier transform，FFT)的算法因具备计算速度快等优点被广泛运用，但仍存在频谱泄漏、栅栏效应和频谱混叠等问题^[1⇓⇓-4]。

理想状态下时域下的信号序列通常认为是无限长，而实际上在分析时域信号时都是截取某一段信号时序，需要加窗函数，相当于频域卷积，并且通常情况下在进行时域截取信号时，为非整数信号周期，在分析的频谱图上就会存在旁瓣。为了减少频谱泄漏，需要乘以加权的窗函数，例如Hanning窗、Hamming窗、Kaiser窗、Blackman窗、Nuttall窗等。实际在对信号进行处理时，时域和频域都是离散的，连续的信号不会造成频率分量丢失，但离散信号时序经过FFT从时域到频域会造成频谱中频率分量丢失，导致频谱不完整。

结合工程应用，采集的电信号中会存在与整数次谐波相近频率的间谐波含量，采用传统的快速傅里叶变换很难将其区分，此时需要提高频率分辨率，可以采用提高采样点数N或者减小采样频率的方法。其中也有学者提出频谱细化算法和加窗插值算法对频谱进行分析，例如线性调频Z变换(或CZT)^[5]、小波变换^[6]、复调制频谱细化(或ZFFT)^[7-8]、Kaiser窗双谱线插值算法^[9]、Kaiser窗三谱线插值算法^[10]和Nuttall窗三谱线插值^[11]等。这些算法都能够提高频率分辨率，但仍有其局限性，例如CZT算法不适用于多频电信号且相对密集的有效分析，且实时性相对较差；Kaiser窗和Nuttall窗的多谱线插值算法测量精度不能满足所有条件。

针对电网信号测量分析中需要数据精度的同时也需要满足其实时性的要求，本文结合复调制频谱细化(ZFFT)算法和基于Blackman-Harris窗四谱线插值算法对测量信号进行分析，并结合主瓣干扰判定条件对算法进行处理，最后通过仿真验证改进型方法的精准度和实时性。

2 基于Blackman-Harris窗四谱线插值校正算法

针对基于Blackman-Harris窗双谱线插值算法^[12]和基于Blackman-Harris窗三谱线插值校正算法^[13]存在的误差，本文提出基于Blackman-Harris窗四谱线插值算法来减少校正误差，提高测量精确性。

2.1 Blackman-Harris窗

加窗函数能够在发生频谱泄漏时有效改善并降低频谱泄漏带来的影响，同时提高频谱分辨率。通常窗函数具有如下几种特征参数：主瓣宽度、旁瓣衰减速率、最高旁瓣高度、幅值失真度等。因此，在选择窗函数时要结合实际测量电信号和分析窗函数的旁瓣特性。

表1所示为部分窗函数的特性^[14]，根据窗函数的特性和电信号的性质来选择所需的窗函数。

表1 部分窗函数特性

窗的类型	主瓣宽度	旁瓣电平/dB	旁瓣衰减速率/ (dB/oct)
Triangle	8π/N	-27	12
Rectangle	4π/N	-13	6
Hamming	8π/N	-43	6
Blackman	12π/N	-58	18
Blackman-Harris	16π/N	-92	6
Hanning	8π/N	-32	18
4项1阶Nuttall	16π/N	-93	18
4项3阶Nuttall	16π/N	-83	30

新窗口打开| 下载CSV

由表1可知，Blackman-Harris窗函数的旁瓣峰值电平相比较于其他窗函数较低(-92 dB)，旁瓣衰减速率为6 dB/oct。衰减速率越大，旁瓣峰值电平越小，对频谱泄漏的治理能力越好。综合考虑实际情况，本文选取Blackman-Harris窗函数。

Blackman-Harris窗的本质是一个4项系数的加权余弦窗，其函数的时域表达式为

(1)

w_{R} (n) = \sum_{l= 0}^{K - 1} {(- 1)}^{l} a_{l} c o s (\frac{2 π ln}{N})

式中，K为余弦窗的项数；N为窗函数的数据长度，n=1,2,3,…,N−1。取K=3，得到一般表达式为

(2)

w_{R} (n) = a_{1} - a_{2} c o s (\frac{2 π n}{N}) + a_{3} c o s (\frac{4 π n}{N}) - a_{4} c o s (\frac{6 π n}{N})

式中，a₁=0.358 75，a₂=0.488 29，a₃=0.141 28，a₄=0.011 68。

2.2 四谱线插值算法

为了方便理解和计算，对电信号 x(t)采取只含单一次谐波分量的信号进行分析(式(3))，其中以采样频率F_s对信号x(t)进行采样得到的离散时间信号如式(4)所示

(3)

x (t) = A s i n (2 π f_{0} t+ φ_{0})

(4)

x (n) = A s i n (\frac{2 π f_{0} n}{F_{s}} + φ_{0})

式中，f₀表示单一次谐波的频率；A表示单一次谐波的幅值；φ₀表示单一次谐波的相位^[15]。

现对信号x(n)进行加窗函数w(n)处理，即得到式(5)

(5)

x_{w} (n) =x (n) \cdot w (n)

代入式(4)，对式(5)进行离散傅里叶变换，得到表达式为

(6)

X_{D} (k Δ f) = \sum_{- \infty}^{+ \infty} x (n) \cdot w (n) e xp (- j 2 π f n)

(7)$ \begin{array}{c} X_{\mathrm{D}}(k \Delta f)=\frac{A}{2 \mathrm{j}}\left\{\exp \left(\mathrm{j} \varphi_{0}\right) W\left[\frac{2 \pi\left(k \Delta f-f_{0}\right)}{F_{\mathrm{s}}}\right]-\right. \\ \left.\exp \left(-\mathrm{j} \varphi_{0}\right) W\left[\frac{2 \pi\left(k \Delta f+f_{0}\right)}{F_{\mathrm{s}}}\right]\right\} \end{array}$

忽略负频点的旁瓣影响，最终得到简化表达式为

(8)

X_{D} (k Δ f) = \frac{A}{2 j} exp (j φ_{0}) W [\frac{2 π (k Δ f - f_{0})}{F_{s}}]

式中，W[∙]为窗函数的离散傅里叶变换形式；∆f为频率分辨率，∆f=F_s/N，k_r=f₀/∆f，其中k_r为峰值频点，由于是离散采样，所以通常k_r都是非整数。

对式(2)中Blackman-Harris窗的时域表达式进行离散傅里叶变换，得到表达式为

(9)

\begin{matrix} W_{BH} (ω) = b_{1} W_{R} (ω) + b_{2} [W_{R} (ω - \frac{2 π}{N}) + W_{R} (ω+ \frac{2 π}{N})] + \\ b_{3} [W_{R} (ω - \frac{4 π}{N}) + W_{R} (ω+ \frac{4 π}{N})] + \\ b_{4} [W_{R} (ω - \frac{6 π}{N}) + W_{R} (ω+ \frac{6 π}{N})] \end{matrix}

式中，b₁=0.358 75，b₂=0.224 145，b₃=0.070 64，b₄=0.005 8。W_R(∙)为矩形窗的频谱幅度函数表达式。式(8)中W[∙]的参数根据现有条件可以转换为^[16]

(10)

\frac{2 π (k Δ f - f_{0})}{F_{s}} = \frac{2 π (k - k_{r})}{N}

W_R(∙)的表达式为

(11)

W_{R} (ω) = \frac{s i n (ω N / 2)}{s i n (ω / 2)}

由式(8)可知ω的表达式(10)，将其代入式(12)进行离散采样

W_{R} (ω \pm \frac{2 π m}{N}) = W_{R} (\frac{2 π (k - k_{r})}{N} \pm \frac{2 π m}{N}) =

(12)

\frac{s i n [π (k - k_{r} \pm m)]}{s i n [\frac{π (k - k_{r} \pm m)}{N}]} \approx \frac{s i n [π (k - k_{r} \pm m)]}{π (k - k_{r} \pm m)} N

因式(9)中含有 W_R(∙)形式，为利于后面计算所以写成W_R(ω±2πm/N)的形式。

通常情况下对信号采样很难做到同步采样，所以f₀表示的峰值频率通常不在离散的频谱点上，如图1所示，在峰值频点k_r两侧分别有两条最靠近k_r频点的谱线，分别是k_r1、k_r2、k_r3、k_r4。分析单谱线、双谱线和三谱线算法可以知道，靠近峰值频点k_r的谱线理论上拥有更多的频谱信息，所以要尽可能多地分析谱线信息。

图1

新窗口打开| 下载原图ZIP| 生成PPT

图1 信号频谱图

引入变量ρ，ρ代表位于峰值频点k_r两边的幅值之比，表达式为

(13)

ρ = \frac{y_{4} + y_{3} - y_{1} - y_{2}}{y_{4} + y_{3} + y_{1} + y_{2}}

将式(8)、(10)一并代入式(13)得到表达式为

(14)

ρ = \frac{|W (- γ + 0.5)| + |W (- γ + 1.5)| - |W (- γ - 1.5)| - |W (- γ - 0.5)|}{|W (- γ + 0.5)| + |W (- γ + 1.5)| + |W (- γ - 1.5)| + |W (- γ - 0.5)|}

从式(14)可以看出，ρ可以表示成γ的函数表达式，ρ=y(γ)，其中反函数可以表示为γ=y⁻¹(ρ)，可以根据ρ求出γ，利用Matlab中的polyfit多项式曲线拟合函数拟合出表达式，在调用polyfit函数时需要选择拟合多项式的系数n，n取值过大会增加运算成本，通过仿真遍历可知，当n取值大于7后，多项式中阶数大于7的项对计算结果影响已经很小，且会增加运算时间，所以通过polyfit多项式曲线拟合函数拟合的通用表达式为

(15)

γ = c_{1} ρ + c_{3} ρ^{3} + c_{5} ρ^{5} + \dots + c_{2 n + 1} ρ^{2 n + 1}

式中，c₁、c₃、c₅、c₂_n₊₁分别代表每个奇数次项(2n+1)前面的系数。

其中信号频率的校正公式为

(16)

f_{0} = k_{r} Δ f = (γ + k_{r 2} + 0.5) Δ f

由条件可知测量信号的初始相位修正表达式为

(17)

φ_{0} = a r g [X_{D} (k Δ f)] + \frac{π}{2} - a r g {W [\frac{2 π (k - k_{r})}{N}]}

测量信号的幅值是根据图1中的四条谱线通过加权推算得来，且因为k_r2、k_r3是最靠近峰值频点的两条谱线，所以需要乘以更大的权值，这四条谱线的权值分别取1、3、3、1。由此得到信号幅值修正表达式为

(18)

A = \frac{2 (y_{1} + 3 y_{2} + 3 y_{3} + y_{4})}{3 |W (- γ + 0.5)| + |W (- γ + 1.5)| + |W (- γ - 1.5)| + 3 |W (- γ - 0.5)|}

其中，当频谱分析中N取值较大时，此时可以表达为

(19)

A = N^{- 1} (y_{1} + 3 y_{2} + 3 y_{3} + y_{4}) g (γ)

g(γ)经过拟合后的表达式可表示为

(20)

g (γ) = b_{0} + b_{2} γ^{2} + \dots + b_{2 n} γ^{2 n}

式中，b₀、b₂、b₂_n分别代表偶数次项(2n)前面的系数。

将式(9)和式(12)代入式(14)，利用Matlab中的曲线拟合函数polyfit()进行求解，为保证数据的精度，在选择拟合次数时不能过小，此次拟合选取的次数为7，在γ取值范围内取值代入，通过函数polyfit(ρ，γ，7)拟合求得Blackman-Harris窗四谱线插值算法中γ为

(21)

\begin{array}{l} γ = 2.653 081 62 ρ + 0.732 016 54 ρ^{3} + \\ 0.401 694 30 ρ^{5} + 0.105 680 43 ρ^{7} \end{array}

如图2所示，图2a为通过拟合得到的曲线图，图2b为真实数值与拟合曲线数值的差值。

图2

新窗口打开| 下载原图ZIP| 生成PPT

图2 拟合曲线和差值图

将式(9)、式(12)和式(18)结合，利用曲线拟合函数polyfit(γ，g(γ)，6)求得Blackman-Harris窗四谱线插值算法中g(γ)为

(22)$\begin{array}{c} g(\gamma)=1.65890762+0.62581064 \gamma^{2}+ \\ 0.15403867 \gamma^{4}+0.01627366 \gamma^{6} \end{array}$

将求得的式(21)、(22)代入式(16)、(17)、(19)，即可得到基于Blackman-Harris窗函数四谱线插值算法的频率、相位和幅值。

3 主瓣干扰判定

3.1 频率相近分量频谱分析

一般测量得到的电信号中不仅含有整数次谐波，还含有间谐波，且当间谐波与整数次谐波相近时就会发生频谱混叠，导致的结果就是不能区分相近的频率分量。

设待测信号的表达式为

(23)

x (t) = A_{1} s i n (2 π f_{1} t + φ_{1}) + A_{2} s i n (2 π f_{2} t + φ_{2})

式中，A₁、A₂表示两个频率分量的幅值；f₁、f₂表示信号中含有的两个频率分量；φ₁、φ₂表示两个频率分量的初相位。相近频率频谱图如图3所示。

图3

新窗口打开| 下载原图ZIP| 生成PPT

图3 相近频率频谱图

由图3可以看出，当两个频率相差比较近时，在频谱图中只能看到一个波峰，也就是发生了频谱混叠，只存在一个主瓣，且发生了主瓣干扰，使得被干扰频率的谱线幅值也发生了失真现象。如果采用加窗插值算法对测量信号进行分析处理，得到的数据准确度会产生偏差，且不能反映出真实的频谱情况。

3.2 主瓣干扰分析

现就需要对测量的电信号进行判定，判定的依据为是否发生了主瓣间的干扰。一般可以认为，当没有发生主瓣干扰时，且不受信号噪声影响时，主瓣内谱线间的相位和幅值有关联，主瓣内相邻谱线间的相位差绝对值等于π，所以在判定主瓣干扰时可以作为依据。

余弦组合窗是目前运用最多的一类窗函数，以此为例，其时域表达式为式(1)，可知其频谱表达式为

(24)

W (k) = \sum_{m = 0}^{M - 1} {(- 1)}^{m} \frac{a_{m}}{2} [W_{o} (k - m) + W_{o} (k + m)]

式中，W_o(∙)为矩形窗的离散傅里叶变换形式。

(25)

W_{o} (k) = \frac{\sin (k π)}{\sin (\frac{k π}{N})} \cdot \exp (- j \frac{N - 1}{N} k π)

将式(25)代入式(24)化简可得

(26)

\begin{matrix} W (k) = s i n (k π) \cdot e xp (- j k π) [\sum_{m = 0}^{M - 1} {(- 1)}^{m} \frac{a_{m}}{2} \cdot \\ \frac{\sin (\frac{2 k π}{N})}{\sin (\frac{π (k - m)}{N}) \sin (\frac{π (k + m)}{N})}] \end{matrix}

将式(26)写成

(27)

W (k) = F (k) \cdot exp(- j k π)

(28)

\begin{matrix} F (k) = s i n (k π) \cdot [\sum_{m = 0}^{M - 1} {(- 1)}^{m} \frac{a_{m}}{2} \cdot \\ \frac{\sin (\frac{2 k π}{N})}{\sin (\frac{π (k - m)}{N}) \sin (\frac{π (k + m)}{N})}] \end{matrix}

设加窗信号经过离散傅里叶变换为X(k)，其中在测量第i项谐波时，若不考虑负频率和周边谐波的泄漏影响可知

(29)

X (k) = \frac{A_{i}}{2 j} exp(j φ_{i}) W (k - k_{i})

设在主瓣内，k₁和k₂是相邻的两根谱线，k₂=k₁+1，可知

(30)

\begin{matrix} X (k_{1}) = \frac{A_{i}}{2 j} exp(j φ_{i}) W (k_{1} - k_{i}) = \\ \frac{A_{i}}{2 j} exp(j [φ_{i} - (k_{1} - k_{i}) π]) F (k - k_{i}) \end{matrix}

(31)

\begin{matrix} X (k_{2}) = X (k_{1} + 1) = \frac{A_{i}}{2 j} exp(j φ_{i}) W (k_{1} + 1 - k_{i}) = \\ \frac{A_{i}}{2 j} exp(j [φ_{i} - (k_{1} + 1 - k_{i}) π]) F (k + 1 - k_{i}) \end{matrix}

可知两条谱线的相位差为

(32)

\begin{matrix} | a r g [X (k_{1})] - a r g [X (k_{2})] | = | [φ_{i} - (k {}_{1}- k_{i}) π - \frac{π}{2}] - \\ [φ_{i} - (k {}_{1}+ 1 - k_{i}) π - \frac{π}{2}] | = π \end{matrix}

由此可知在不受信号噪声及主瓣干扰的情况下，主瓣内相邻谱线间的相位差为π。

引入参数∆θ，根据表达式

(33)

Δ θ = ||δ_{1} - δ_{2}| - π|

式中，δ₁、δ₂分别是靠近峰值频点最近的两个谱线相位。按照实际工程应用来说，∆θ通常都不等于0，所以需要选定一个恰当的值£来权衡是否发生了主瓣干扰，£值的大小直接决定了测量时长和精准度。如果£的取值过大，就会导致发生了主瓣干扰但是没有判别出来的情况，从而导致漏判；若£的取值过小，将会导致过于灵敏，使得计算时长增加^[17]。本文结合实际工程应用，对不同频率差信号分量进行仿真遍历，在引入较少运算量的同时对谐波分量进行有效分离检测，最终选取£=8×10^-6。

4 复调制频谱细化算法

当测量电信号的某一频段上判定了主瓣干扰，说明在这一频段上存在不止一个频率分量，且它们的频率值相差很小，已经发生了频谱混叠的现象，这时就需要利用频率分辨率更高的处理算法来区分在这一频段上存在的频率分量^[18]。本文选取的是复调制频谱细化算法ZoomFFT(或ZFFT)。

ZFFT算法能够实现在某一频段上对频谱进行放大处理，使其在这一频段上频率分辨率提高，放大D倍^[19]。对于频率分辨率∆f=F_s/N，根据公式可知想要提高频率分辨率，可以改变F_s和N这两个值。

(1) 保持N值不变，减小采样频率F_s的值，使得∆f减小，也就是频率分辨率提高，但是根据奈奎斯特采样定理可知，如果采样频率F_s小于两倍信号中的最高频率，测量的信号频谱将会发生频谱混叠现象，即要保证F_s≥2f_max。

(2) 保持采样频率F_s不变，增大采样点数N，也能提高频率分辨率，但是随着N值增大，带来的是计算量的增大，从而导致计算时长的增加，这对一些实时性要求高的场合就不太适用，且会占用大量内存。

ZoomFFT算法实现的流程图如图4所示。

图4

新窗口打开| 下载原图ZIP| 生成PPT

图4 ZFFT算法流程图

ZFFT算法实现的过程总共分为7个部分，其具体步骤如下所示。

(1) 滤波。根据奈奎斯特采样定理可知，FFT的频率分析范围有限，2/F_s为最高的分析频率，所以需要滤除其频率以上的频率，防止在分析频谱时出现混叠。

(2) 离散采样。对时域信号x(t)进行离散化处理，得到x(n)，采样长度为放大倍数D与采样点数的乘积。

(3) 频移。经过处理后的离散信号x(n)乘以旋转因子exp(-j2πfn)，目的是为了让所需分析的频段中心频率f移至频率零点，得到新的信号。

(4) 低通滤波。目的是为了滤除所需频谱细化频段以外的信号分量，只保留所需细化的频段。

(5) 重采样。对经过低通滤波器之后保留的频段进行重新采样，采样时每个点的间隔为D，采样频率由之前的F_s变为F_s/D，由此可见，采样频率降低了D倍。

(6) FFT分析。对重采样后的离散信号进行FFT分析，分析的点数仍然为N，但是采样频率由之前的F_s变为F_s/D，所以频率分辨率也会提升D倍。

(7) 频率调整。由步骤(3)可知信号经过频移，所以需要对最终得到的频谱信号进行调整。

综上，通过ZFFT算法处理过后的频谱在以频率f为中心的频段上，将频率分辨率由之前的∆f=F_s/N降为∆f=F_s/ND，也就是将频率分辨率提升了D倍。放大的频段范围是∆F=f±F_s/2D。

观察ZFFT算法流程步骤可知，ZFFT算法只能针对想要观察的目标频段进行放大处理，不能够对全频域内进行放大。结合如今的电力系统测量分析可知，大多分析为50次谐波以内，所以一般分析0~2.5 kHz谐波含量，若只对存在主瓣干扰的频段采用频谱细化分析，对只含有整数次谐波的频段采用加窗插值算法分析，这将会大大减少计算时间，不仅提高了谐波检测精度，同时也提高了实时性。放大倍数的选择也决定了ZFFT算法的有效性，如果放大倍数过小，则不能够对存在主瓣干扰的频率实现有效区分，若放大倍数过大，则会导致采样的长度过长。

综上所述，对测量电能质量提出了改进方法，流程如图5所示。

图5

新窗口打开| 下载原图ZIP| 生成PPT

图5 改进型测量方法流程图

根据图5可知，结合Blackman-Harris窗插值算法和ZFFT算法，利用条件主瓣干扰判定选用算法，弥补了ZFFT算法不能全频域对频谱进行细化，基于Blackman-Harris窗四谱线插值不能够精确区分频率相近分量的缺点。

5 仿真验证与分析

为了验证所提改进型方法的准确度和实时性，这里对整数次谐波和间谐波进行了频谱分析。主要的电力参数设置如表2所示。

表2 电力参数

信号分量	频率f/Hz	幅值U/V	初相位φ/(°)
X₁	50	80	0
X₂	51	20	0
X₃	55	15	10
X₄	100	20	30
X₅	105	10	60
X₆	150	50	150
X₇	152	15	-20
X₈	200	14	30
X₉	250	20	-50

新窗口打开| 下载CSV

结合实际电力系统工程分析，此次试验分析的谐波范围为0~50次，相应频谱范围为0~2.5 kHz，最高的频率为2.5 kHz。因此，根据奈奎斯特采样定理选取采样频率为5 kHz，采样点数N=1 024，ZoomFFT选取的放大倍数D=10。

此外，在表2中，参数中除了含有基波和其他整数次谐波含量，还设置了间谐波且频率相近。

图6为改进型方法仿真局部频谱图。根据仿真结果实测数据得到参数结果如表3所示。

图6

新窗口打开| 下载原图ZIP| 生成PPT

图6 改进型方法仿真频谱图

表3 改进型方法仿真结果数据

信号分量	频率f/Hz	幅值U/V	初相位φ/(°)
X₁	50.000 5	80.001 0	0.000 8
X₂	51.001 2	20.000 6	0.000 5
X₃	55.000 5	16.000 5	10.000 9
X₄	100.001 1	20.000 6	30.000 8
X₅	105.000 0	10.000 3	60.000 6
X₆	149.999 9	49.999 5	149.999 1
X₇	151.999 8	14.999 5	-19.999 2
X₈	200.000 0	15.000 6	30.000 9
X₉	250.000 0	20.000 6	-50.000 5

新窗口打开| 下载CSV

同时分别对传统的FFT谐波检测方法和Blackman-Harris窗插值算法进行仿真验证对比，图7为基于Blackman-Harris窗四谱线插值算法仿真频谱图，图8为FFT算法仿真频谱图。信号仿真设置参数如表2所示。

图7

新窗口打开| 下载原图ZIP| 生成PPT

图7 基于Blackman-Harris窗四谱线插值算法仿真频谱图

图8

新窗口打开| 下载原图ZIP| 生成PPT

图8 传统FFT算法仿真频谱图

对比图6、图7、图8和表3、表4可知，信号分量X₂在基于Blackman-Harris窗四谱线插值算法中无法被检测出，不能将所有的间谐波区分开来，同样条件下，采用本文所提方法能够将所有间谐波全部区分并取得精确的参数值。

表4 基于Blackman-Harris窗四谱线插值算法仿真结果数据

信号分量	频率f/Hz	幅值U/V	初相位φ/(°)
X₁	50.001 2	80.011 5	0.010 7
X₂	—	—	—
X₃	54.989 5	15.385 5	9.993 1
X₄	99.998 9	19.998 0	29.996 9
X₅	104.997 8	9.995 3	59.993 2
X₆	149.998 9	49.998 8	149.998 8
X₇	151.984 9	14.298 8	-19.002 3
X₈	199.998 5	14.978 2	29.991 2
X₉	249.999 2	19.998 9	-49.998 8

新窗口打开| 下载CSV

通过对比不同谐波检测分析的方法可知，对电信号加不同的窗函数和采用不同的插值算法得到的时间成本都不相同^[20]，其中ZoomFFT算法在提高频率分辨时是以分析更多的采样点数得到的。如果只使用ZFFT算法对信号进行分析，将会耗费大量的时间和内存单元。通过仿真对比，不同加窗插值算法的平均计算时长如表5所示。由表5可知，ZFFT算法精确度高，但计算时间过长，不利于测量要求的实时性，且会占据大量内存。因此，本文采用的方法兼容了ZFFT和加Blackman-Harris窗插值算法的优点。

表5 仿真平均计算时长

对比算法	采样点数目	计算时长/ms
ZFFT	10 240	29.00
Blackman-Harris窗	1 024	0.26
Hanning窗	1 024	0.21
Hamming窗	1 024	0.22
本文算法	1 024\|10 240	2.20

新窗口打开| 下载CSV

6 结论

本文结合实际的工程应用，对电网电能质量参数进行快速及准确的采集和处理。针对电网谐波分量发生的频谱泄漏和混叠现象，造成难以精确区分相近谐波分离检测的问题，本文采用基于ZoomFFT和多谱线插值算法的谐波分析方法，得出如下结论。

(1) 基于Blackman-Harris窗四谱线插值算法，能够满足谐波高精度测量且计算时长少的要求。

(2) 在发生主瓣干扰的情况下，采用频率分辨率更高的频谱细化ZoomFFT算法精准分析。

(3) 通过仿真数据对比分析验证了本文方法的精确性和实时性，结合两种算法各自优点，保证了测量精度，同时也满足了更少的运算时长，从而确保了现代电力系统发展所需的要求。

参考文献

原文顺序

文献年度倒序

文中引用次数倒序

被引期刊影响因子

[1]

李晨.

谐波背景下的电能计量误差分析与滤波研究[D]. 昆明: 昆明理工大学, 2018.