CN103575984A - 基于凯塞窗双谱线插值fft的谐波分析方法 - Google Patents

Abstract

基于凯塞(Kaiser)窗函数双谱线插值FFT的谐波分析方法适用于对电力网电压、电流和谐波分析与监测。首先使用线性调频Z变换从含有谐波的电力信号中高精度提取基波信号参数(幅值、频率和相位)。然后对被分忻电力信号加Kaiser窗函数截断电力信号，并用FFT计算出剩余信号的频谱。再依据各谐波频率对Kaiser窗函数在频域内插值，精确计算出各电力谐波的参数。本发明能有效克服基波频率波动与白噪声对谐波分忻的影响，在非整数周期截断条件下，对21次谐波信号的频率计算相对误差为1.4×10％，幅值计算相对误差≤0.002％，初相位计算相对误差≤0.0001％。

Description

基于凯塞窗双谱线插值FFT的谐波分析方法

技术领域

本发明涉及一种电网电压和电流波形畸变的分析和自动监测算法，可用于各种电网电压和电流博畸变的分析仪器和自动监测装置。属于电力测量和自动化技术领域。

背景技术

随着各种非线性负荷特别是电力电子设备在电力系统中的广泛应用，电网谐波污染日益严重。FFT算法因其算法简单，有较好的实用性，因此成为当前谐波分析的主要算法。

FFT算法需要在完全同步采样的条件下，当信号的测量时间不等于信号周期的整数倍，或当信号含有非整数次谐波时，由传统FFT算法得到的各次谐波的频率和真实的谐波频率之间有较大的误差。实际电网频率总会在额定频率附近波动，消除频率同步误差一般采用两种方法：一是由硬件实现同步采样，在采集系统中加入锁相同步技术，优点是信号处理比较简单，但由于锁相环响应较慢，不能及时跟踪信号频率的快速变化，从而不能实现真正意义的同步采样；另一种方法是通过选择谱能量主要集中在主瓣，旁瓣谱能量小、且幅值衰减块的窗函数，以减小谱间干扰，即频谱的长范围泄露；通过双谱线间插值修正，仪减小栅栏效应，进而提高谐波估计精度。许多学者采用加窗插值法都有效地提高了谐波估计的精度，如专利《纳托尔窗函数连续频谱内插电力谐波参数获取方法》(专利号：201110154995.3)和专利《汉宁窗函数连续频谱内插电力谐波参数获取方法》(专利号：201110154993.4)。

由于Kaiser窗可定义一组可调的窗函数，其主瓣能量和旁瓣能量的比例近乎最大，且可自由选择主瓣宽度和旁瓣高度之间的比重，因此，对信号加权更加灵活，通过改变Kaiser窗的形状参数，可以满足不同的设计要求。通过实验证明基于Kaiser窗的双谱线插值FFT方法较上述两个专利更能克服基波频率波动与白噪声对谐波分析的影响。

发明内容

本发明的目的是提供一个基于凯塞(Kaiser)窗双谱线插值FFT的谐波分析方法，可用于各种电网电压和电流波形畸变的分析仪器和自动监测装置。

本发明的Kaiser窗双谱线插值FFT的谐波分析方法采用如下步骤：

步骤1：含谐波分量的信号x(t)以采样频率f_s均匀采样得到离散时间信号为：

公式1

其中：H为谐波项数；当h＝1时，r_h＝1，f₁，4，

分别为基波的频率、幅值和初相位。当h≠1时，r_h，A_h，

分别为第h项谐波的次数、幅值和初相位。

步骤2：用Kaiser窗函数W(·)对公式1中信号x(n)进行离散抽样，并忽略负频率点-r_hf₁处信号的旁瓣影响，得到加窗后的离散傅里叶变换的表达式为

公式2

式中，r_hf₁＝k_nΔf为第h项谐波的频率，Δf为离散抽样间隔；W(·)是Kaiser窗的连续频谱函数，表达式为

$W\left(w\right)=\frac{N}{I_0\left(\mathrm{\β}\right)}\frac{\sinh \left(\sqrt{{\mathrm{\β}}^2-{(\mathrm{Nw}/2)}^2}\right)}{\sqrt{{\mathrm{\β}}^2-{(\mathrm{Nw}/2)}^2}}$ 公式3

其中，I₀(β)是第1类变形零阶贝塞尔函数，β为窗函数的形状参数，N为窗函数的长度。β值越大，Kaiser窗函数频谱的旁瓣峰值越小，渐近衰减速率越大。

步骤3：由于栅栏效应，第i个谐波的频率k_iΔf很难正巧位于抽样频点上，即k_i一般不是整数。设在峰值频点附近抽样得到的幅值最大和次最大谱线分别为k_i1和k_i2，k_i1≤k_i≤k_i2＝k_i1+1，这两条谱线的幅值分别为y₁＝|xw(k_i1Δf)|和y₂＝|xw(k_l2Δf)|。由于0≤k_i-k_i1＜1，并引入参数λ_i＝k_i-k_i1，β＝aπ，同时将公式3中频率信号平移(N-1)/2，得到关于λ_i的一元超越方程，如公式4所示：

$\frac{y_1}{y_2}=\frac{\left|\sinh \left(\mathrm{\π}\sqrt{a^2-{\mathrm{\λ}}_i^2}\right)\right|\left|\mathrm{\π}\sqrt{a^2-{(1-{\mathrm{\λ}}_i)}^2}\right|}{\left|\mathrm{\π}\sqrt{a^2-{\mathrm{\λ}}_i^2}\right|\left|\sinh \left(\mathrm{\π}\sqrt{a^2-{(1-{\mathrm{\λ}}_i)}^2}\right)\right|}$ 公式4

借助Matlab平台的solve(.)函数求解参数λ_j。

步骤4：根据参数λ_j，对j次奇次谐波的频率、幅值和相位分别由公式5，公式6，公式7得到。

f_j＝(k_j1+λ_j)f_s/N                                                                公式5

$A_j=\frac{2y_{j1}{I}_0\left(\mathrm{\β}\right)\sqrt{{\mathrm{\β}}^2-{\left(\mathrm{\π}{\mathrm{\λ}}_j\right)}^2}}{N\left|\sinh \left(\sqrt{{\mathrm{\β}}^2-{\left(\mathrm{\π}{\mathrm{\λ}}_j\right)}^2}\right)\right|}$ 公式6

公式7

第l次偶次谐波的频率、幅值、相位分别由公式8、公式9和公式10所得：

f_l＝(k_l2+(1-λ_l))f_s/N                                                            公式8

$A_l=\frac{2y_{l2}{I}_0\left(\mathrm{\β}\right)\sqrt{{\mathrm{\β}}^2-{\left(\mathrm{\π}{\mathrm{\λ}}_l\right)}^2}}{N\left|\sinh \left(\sqrt{{\mathrm{\β}}^2-{\left(\mathrm{\π}{\mathrm{\λ}}_l\right)}^2}\right)\right|}$ 公式9

公式10

有益效果：本发明提供了一种电力谐波分析的新方法。Kaiser窗可定义一组可调的窗函数，自由选择主瓣与旁瓣衰减之间的交换关系，结合FFT进行电力谐波分析，精确计算出各次谐波的参数。对包含21次谐波的动态信号仿真和三相多功能谐波电能表实际应用表明，采用本发明的电力谐波分析方法，具有较高的计算精度，设计灵活、实用价值高，且三相多功能谐波电能表准确度达到0.2S级，2～21次谐波分析满足GB/T14549-1993的A类谐波测量仪器要求。

具体实施方式

为了实现上述目的，本发明的实施可直接用分压器或从电压传感器P图而成侧取得电网的母线电压信号、从电流传感器CT取得电流信号，经过适当的信号调理后送达信号采用入口。

[0000]步骤1.采用被分析电力信号电压或电流，应用专利1中提到的CZT值计算含谐波分量的信号x(n)，n∈[0.N]。

其中N为线性调频Z变换时在频域内的抽样点数。

设被分析电力信号为：

其中：片为谐波项数。

步骤2：用Kaiser窗函数w((n)对信号x(n)进行处理，得到加窗后信号的离散傅里叶变换表达式：

式中，W(·)是Kaiser窗的连续频谱函数，其频域表示为：

$W\left(w\right)=\frac{N}{I_0\left(\mathrm{\β}\right)}\frac{\sinh \left(\sqrt{{\mathrm{\β}}^2-{(\mathrm{Nw}/2)}^2}\right)}{\sqrt{{\mathrm{\β}}^2-{(\mathrm{Nw}/2)}^2}}$

步骤3：在峰值频点附近抽样得到的幅值最大和次最大谱线，分别记为k_i1和k_i2，则这两条谱线的幅值分别为y₁＝|xw(k_i1Δf)|和y₂＝|xw(k_i2Δf)|。令λ_i＝k_i-k_i1，则有

$\frac{y_1}{y_2}=\frac{\left|W(-2\mathrm{\π}{\mathrm{\λ}}_i/N)\right|}{\left|W\left(2\mathrm{\π}(1-{\mathrm{\λ}}_i)N\right)\right|}$

[0000]步骤4：设Kaiser窗函数中β＝πa，将Kaiser函数中频域信号平移(N-1)/2，有：

$W\left(w\right)=\frac{N-1}{I_0\left(\mathrm{\πa}\right)}\frac{\sinh \left(\sqrt{{\left(\mathrm{\πa}\right)}^2-{((N-1)w/2)}^2}\right)}{\sqrt{{\left(\mathrm{\πa}\right)}^2-{((N-1)w/2)}^2}}{e}^{-j\frac{N-1}{2}w}$

设w＝2kπ/N，则有

$\frac{y_1}{y_2}=\frac{\left|\sinh \left(\mathrm{\π}\sqrt{a^2-{\mathrm{\λ}}_i^2}\right)\right|\left|\mathrm{\π}\sqrt{a^2-{(1-{\mathrm{\λ}}_i)}^2}\right|}{\left|\mathrm{\π}\sqrt{a^2-{\mathrm{\λ}}_i^2}\right|\left|\sinh \left(\mathrm{\π}\sqrt{a^2-{(1-{\mathrm{\λ}}_i)}^2}\right)\right|}$

步骤5：根据参数λ_j，对j次奇次谐波的频率、幅值和相位分别为：

f_j＝(k_j1+λ_j)f_s/N

$A_j=\frac{2y_{j1}{I}_0\left(\mathrm{\β}\right)\sqrt{{\mathrm{\β}}^2-{\left(\mathrm{\π}{\mathrm{\λ}}_j\right)}^2}}{N\left|\sinh \left(\sqrt{{\mathrm{\β}}^2-{\left(\mathrm{\π}{\mathrm{\λ}}_j\right)}^2}\right)\right|}$

第l次偶次谐波的频率、幅值、相位分别为：

f_l＝(k_l2+(1-λ_l))f_s/N

$A_l=\frac{2y_{l2}{I}_0\left(\mathrm{\β}\right)\sqrt{{\mathrm{\β}}^2-{\left(\mathrm{\π}{\mathrm{\λ}}_l\right)}^2}}{N\left|\sinh \left(\sqrt{{\mathrm{\β}}^2-{\left(\mathrm{\π}{\mathrm{\λ}}_l\right)}^2}\right)\right|}$

Claims (1)

1.一种基于凯塞窗双谱线插值FFT的谐波分析方法，其特征是该方法采用如下步骤：

1)步骤1：含谐波分量的信号x(t)以采样频率f_s均匀采样得到的离散时间信号为：

公式1

其中，H为谐波项数；当h＝1时，r_h＝1，f_l，A_l， 分别为基波的频率、幅值和初相位。当h≠1时，r_h，A_h， 

分别为第h项谐波的次数、幅值和初相位。

2)步骤2：用Kaiser窗函数W(·)对公式1中信号x(n)进行离散抽样，并忽略负频率点-r_hf_l处信号的旁瓣影响，得到加窗后的离散傅里叶变换的表达式为：

公式2

式中，r_hf_l＝k_nΔf为第h项谐波的频率，Δf为离散抽样间隔；W(·)是Kaiser窗的连续频谱函数，表达式为：

公式3

其中，I₀(β)是第1类变形零阶贝塞尔函数，β为窗函数的形状参数。

3)步骤3：由于栅栏效应，第i个谐波的频率k_lVf很难正巧位于抽样频点上，即k_i一般不是整数。设在峰值频点附近抽样得到的幅值最大和次最大谱线分别为k_i1和k_i2，k_i1≤k_i≤k_i2＝k_i1+1，这两条普线的幅值分别为y₁＝|xw(k_i1 Vf)|和y₂＝|xw(k_i2Δf)|。由于0≤k_i-k_i1＜1，并引入参数λ_i＝k_i-k_i1，β＝aπ，同时将公式3中频率信号平移(N-1)/2，得到关于λ_i的一元超越方程，如公式4所示：

公式4

借助Matlab平台的solve(·)函数求解参数λ_j。

4)步骤4：根据参数λ_j，对j次奇次谐波的频率、幅值和相位分别由公式5，公式6，公式7得到。

f_j＝(k_j1+λ_j)f_s/N                                                         公式5 

公式6

公式7

第l次偶次谐波的频率、幅值、相位分别由公式8、公式9和公式10所得：

f_l＝(k_l2+(1-λ_l))f_s/N                                                 公式8

公式9

公式10 。