CN103995181A - 一种数字化变电站电能质量谐波分析方法 - Google Patents

Abstract

本发明涉及电力系统数字化变电站电能质量技术领域，提出一种数字化变电站电能质量谐波分析方法，包括：步骤1：接收数字化变电站合并单元MU发送的实时采样数据x[n]；步骤2：对接收的实时采样数据x[n]进行软件数字滤波，滤除高次谐波及扰动信号，从而测得系统频率；步骤3：获取所述数字化变电站采样率，并根据所获取的采样率，相应进行谐波、间谐波分析，得到分析结果；步骤4：根据步骤3的分析结果，进行电能质量的后续指标运算、界面显示和上送至电能质量分析主站。本发明方法满足各种采样率下谐波、间谐波的计算需求，且能大幅度减少计算量，满足实际工程应用需求。

Description

一种数字化变电站电能质量谐波分析方法

技术领域

本发明涉及电力系统数字化变电站电能质量技术领域，具体涉及一种基于数字化变电站的电能质量谐波分析方法。 

背景技术

随着数字化变电站的发展，新型电子式互感器、智能断路器的实用化，以及DL/T860标准的实施，通信网络正逐渐取代变电站间隔层和过程层的并行电缆，并成为保护和控制系统的重要组成部分，合并单元MU以等间隔采样的方式把采集数据发送至电能质量监测装置，采样率一般为4K、8K、12.8K、51.2K等。但是由于电网频率是不断波动的，加上电能质量标准规定频率跟踪范围为±7.5Hz，在这么大的频率偏差范围内，每周波采样点并不是理想的80个点(4K采样率下)；而常规电能质量分析方法都是实时跟踪采样频率进行同步采样，即每周波采样点数是固定的，因此采用常规的电能质量谐波分析算法完全不适用于数字化变电站电能质量分析，简单的DFT运算会产生严重的频谱混叠、栅栏效应及频谱泄露等问题，这将使信号参数的分析(频率、幅值和相位)不准确，无法满足数字化变电站电能质量谐波分析要求。 

针对这一问题，国内外专家做了各种研究，其中，申请号201310145750.3的发明申请提出了一种测量电力系统信号频率及谐波参数的方法，通过对采样数据加7项Rife-Vincent窗的FFT数据中 真实频率点附近的2根最大谱线进行修正，从而得到谐波的频率，幅值和相位。但该发明存在以下不足之处： 

(1)计算量过大，不能满足工程实际需求。MU可能输出10K的采样率供电能质量监测装置分析使用，即每周波200个点，而由于该发明采用双峰插值算法，每个谐波点需要计算3次DFT运算，若系统同时接入100个通道，每秒计算量为2000*1000*5*100*3≈3G次浮点数乘法(2000表示10周波数据窗，每个数据窗2000个点，1000表示计算0-50次间谐波，即需要计算500次谐波实部、虚部，5表示1秒有5个数据窗，100表示100个通道)，这对现有CPU是无法承受的，而且随着数字化变电站电能质量需求的发展，MU输出高采样率给电能质量质量监测装置是必然趋势，按照目前传统电能质量的运行情况，MU今后可能输出51.2K采样率供电能质量分析使用，因此，该方法不能满足工程应用需求。 

(2)该发明并没有明确间谐波(小数次谐波)计算方法。随着电能质量监测技术的发展，间谐波越来越受到电力行业的重视，而该方法中只提及整数次谐波，并未提及间谐波的计算方法。 

发明内容

有鉴于此，本发明提出一种数字化变电站电能质量谐波分析方法，以满足各种采样率下谐波、间谐波的计算需求，且能大幅度减少计算量，满足实际工程应用需求。 

本发明提供的一种数字化变电站电能质量谐波分析方法，包括： 

步骤1：接收数字化变电站合并单元MU发送的实时采样数据x[n]； 

步骤2：对接收的实时采样数据x[n]进行软件数字滤波，滤除高 次谐波及扰动信号，从而测得系统频率； 

步骤3：获取所述数字化变电站采样率，并根据所获取的采样率，相应进行谐波、间谐波分析，得到分析结果； 

步骤4：根据步骤3的分析结果，进行电能质量的后续指标运算、界面显示和上送至电能质量分析主站。 

所述根据所获取的采样率，相应进行谐波、间谐波分析，包括： 

当采样率大于等于10K时，进行如下步骤： 

步骤3.1.1：采样数据插值，将非同步采样转换成同步采样：把采样数据插值成临近的2的N次方次采样点； 

步骤3.1.2：根据插值后的采样率，判断采样率是否过高，过高时进行抽样采样，如当采样率大于等于25.6K时，抽样成12.8K采样率数据，降低运算复杂度 

步骤3.1.3：谐波计算：利用FFT算法计算整数次谐波的谐波值及谐波相角； 

步骤3.1.4：间谐波计算：将数据窗重采样转变成多个FFT计算，将x[n]进行重采样，得到5个采样序列x[5n]、x[5n+1]、x[5n+2]、x[5n+3]、x[5n+4]，并对所述5个采样序列进行FFT运算，得到间谐波的谐波值及谐波相角； 

当采样率小于10K时，进行如下步骤： 

步骤3.2.1：根据步骤2测得的系统频率f₀，计算出每周波实际采样点数； 

步骤3.2.2：假设所求通道幅值为A，初相位为θ，采样后设定为 如下波形： 

$X_n=A\sin (2\mathrm{\π}\frac{f_0}{f_s}n+\mathrm{\θ})$ 式一 

步骤3.2.3：取10周波采样数据，加Hanning窗，假定时域形式为ω(n)，其连续频谱为W(2πf),则加窗后采样信号的连续傅里叶变换忽略负频点-f₀处频峰的旁瓣影响，在正频点f₀附近的离散频谱函数可以表达为： 

$\stackrel{\‾}{X}\left(\mathrm{k\Δf}\right)=\left(\frac{A}{2j}\right)e^{\mathrm{j\θ}}W(2\mathrm{\π}(\mathrm{k\Δf}-f_0)/f_s$ 式二 

式二中离散频率间隔为N是截断数据长度,

步骤3.2.4：非同步采样下，峰值频率f₀＝k0*Δf，设峰值左右两侧的谱线分别为第K1和K2条谱线，令由于0≤K0-K1≤1，所以引入一个辅助参数α＝K0-K1-0.5。 

则有： 

$\frac{y2-y1}{y2+y1}=\frac{\left|W(2\mathrm{\π}(-\mathrm{\α}+0.5)/N|-|W(2\mathrm{\π}(-\mathrm{\α}-0.5)/N\right|}{|W(2\mathrm{\π}(-\mathrm{\α}+0.5)/N|+|W(2\mathrm{\π}(-\mathrm{\α}-0.5)/N|}$ 式三 

令通过多项式逼近的方法，得出Hanning窗下有如下修正公式： 

α＝1.5β 

步骤3.2.5：根据公式四、五，计算实际的基波幅值、相位： 

A＝N^-1(y1+y2)(2.356194+1.155436α²+0.326078α³+0.078914α⁶式四 

$\mathrm{\θ}=\arg \left[X(k1*\mathrm{\Δf})\right]+\frac{\mathrm{\π}}{2}-\arg \left[W(2*\mathrm{\π}*\frac{k1-k0}{N})\right]$ 式五 

步骤3.2.6：根据步骤3.2.1所求得的每周波采样点数，取10周波采样数据窗； 

步骤3.2.7：对获取的10周波采样数据取四舍五入值，然后进行DFT运算； 

步骤3.2.8：根据式六求得间谐波子组值： 

$C_{\mathrm{isg},n}=\sqrt{{\mathrm{\Σ}}_{i=2}^8{C}_{k+i}^2}$ 式六 

式六中： 

表示第n次谐波频谱C_K与第n+1次谐波频谱之间不与其直接相连的连续7个间谐波频谱分量。 

附图说明

为了更清楚地说明本发明实施例或现有技术中的方案，下面将对实施例中所需要使用的附图作一简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施例，对于本领域普通技术人员而言，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。 

图1为本发明一实施例给出的数字化变电站电能质量谐波分析方法的流程示意图； 

图2为本发明一实施例给出的数字化变电站电能质量谐波分析方法中高采样率下的谐波计算流程示意图； 

图3为本发明一实施例给出的数字化变电站电能质量谐波分析方法中低采样率下的谐波计算流程示意图。 

具体实施方式

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整的描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例，都属于本发明保护的范围。 

如图1所示，本发明一实施例提供的一种数字化变电站电能质量谐波分析方法，包括： 

步骤1：接收数字化变电站合并单元MU发送的实时采样数据x[n]； 

步骤2：对接收的实时采样数据x[n]进行软件数字滤波，滤除高次谐波及扰动信号，从而测得系统频率；由于电力系统扰动以及高次谐波的影响，传统的过零点判断法无法准确测出系统频率。本实施例采用软件数字滤波技术滤除高次谐波及扰动信号，其中频率数据窗采用大于10周波数据窗，根据输出系统频率＝54.3Hz，判断出10个上升沿。 

步骤3：获取所述数字化变电站采样率，并根据所获取的采样率，相应进行谐波、间谐波分析，得到分析结果； 

步骤4：根据步骤3的分析结果，进行电能质量的后续指标运算、界面显示和上送至电能质量分析主站。 

在采样率大于10K的情况下(每周波200个点)，利用插值法完全满足谐波、间谐波的国标A级精度要求。而在采样率小于10K情况下，尤其是在目前MU输出4K采样率数据的环境下，线性插值法在高次谐波计算时不能满足精度要求，此时必须保证采样数据的原始性。因此，本方法根据所获取的采样率，分别对高采样率数据及低采样率数据进行谐波、间谐波分析： 

如图2所示，当采样率大于等于10K时，进行如下步骤： 

步骤3.1.1：采样数据插值，将非同步采样转换成同步采样：把采样数据插值成临近的2的N次方次采样点；MU目前常用采样率为10K、12.8K、20K、25.6K等，可分别插值成12.8K、12.8K、12.8K、25.6K，以此规律自适应采样率插值。可以采用Lagrange二次插值法，即选用周边3个点插值成目标采样点。若原始每周波采样点为2的N次方，则不需要重新插值。 

步骤3.1.2：根据插值后的采样率，判断是否需要进行抽样采样，当采样率大于12.8K时，可抽样成12.8K采样率数据，进而降低运算复杂度；本方法中，当采样率为51.2K、25.6K采样率时，可抽样成采样率为12.8K的序列。经验证明，重抽样后完全不影响计算精度。 

步骤3.1.3：谐波计算：利用FFT算法计算整数次谐波的谐波值及谐波相角；此时由于采样点为2的N次方，符合快速傅里叶运算要求，即可求出整数次谐波的谐波值，谐波相角。 

步骤3.1.4：间谐波计算：将数据窗重采样转变成多个FFT计算，将x[n]进行重采样，得到5个采样序列x[5n]、x[5n+1]、x[5n+2]、 x[5n+3]、x[5n+4]，并对所述5个采样序列进行FFT运算，得到间谐波的谐波值及谐波相角；由于间谐波分辨率为0.1Hz，需采用10周波数据窗，而10周波数据窗为非2的N次方，无法直接采用快速傅里叶计算。因此，本发明把数据窗重采样转变成多个FFT计算，将x[n]进行重新采样成x[5n]、x[5n+1]、x[5n+2]、x[5n+3]、x[5n+4]，这样每个序列可重新组合成2的N次方，可适用于快速傅里叶运算。 

$\begin{array}{c}X\left[K\right]=\underset{n=0}{\overset{\frac{N}{5}-1}{\mathrm{\Σ}}}X\left[5n\right]W_N^{5\mathrm{nK}}+\underset{n=0}{\overset{\frac{N}{5}-1}{\mathrm{\Σ}}}X[5n+1]W_N^{(5n+1)K}\\ +{\mathrm{\Σ}}_{n=0}^{\frac{N}{5}-1}X[5n+2]W_N^{(5n+2)K}+{\mathrm{\Σ}}_{n=0}^{\frac{N}{5}-1}X[5n+3]W_N^{(5n+3)K}\\ +{\mathrm{\Σ}}_{n=0}^{\frac{N}{5}-1}X[5n+4]W_N^{(5n+4)K}\end{array}$

根据上式可推出： 

$\begin{array}{c}X\left[K\right]=\underset{n=0}{\overset{\frac{N}{5}-1}{\mathrm{\Σ}}}X\left[5n\right]W_{N/5}^{\mathrm{nK}}+W_N^K\underset{n=0}{\overset{\frac{N}{5}-1}{\mathrm{\Σ}}}X[5n+1]W_{N/5}^{\mathrm{nK}}\\ +W_N^{2K}{\mathrm{\Σ}}_{n=0}^{\frac{N}{5}-1}X[5n+2]W_{N/5}^{\mathrm{nK}}+W_N^{3K}{\mathrm{\Σ}}_{n=0}^{\frac{N}{5}-1}X[5n+3]W_{N/5}^{\mathrm{nK}}\\ +W_N^{K}4{\mathrm{\Σ}}_{n=0}^{\frac{N}{5}-1}X[5n+4]W_{N/5}^{\mathrm{nK}}.\end{array}$

进而还可以推出： 

$\begin{array}{c}X\left[K\right]=\underset{n=0}{\overset{\frac{N}{5}-1}{\mathrm{\Σ}}}X1\left[n\right]W_{N/5}^{\mathrm{nK}}+W_N^K\underset{n=0}{\overset{\frac{N}{5}-1}{\mathrm{\Σ}}}X2\left[n\right]W_{N/5}^{\mathrm{nK}}\\ +W_N^{2K}{\mathrm{\Σ}}_{n=0}^{\frac{N}{5}-1}X3\left[n\right]W_{N/5}^{\mathrm{nK}}+W_N^{3K}{\mathrm{\Σ}}_{n=0}^{\frac{N}{5}-1}X4\left[n\right]W_{N/5}^{\mathrm{nK}}\\ +W_N^{K}4{\mathrm{\Σ}}_{n=0}^{\frac{N}{5}-1}X5\left[n\right]W_{N/5}^{\mathrm{nK}}.\end{array}$

如图3所示，当采样率小于10K时，进行如下步骤： 

步骤3.2.1：根据步骤2测得的系统频率f₀，计算出每周波实际采样点数； 

步骤3.2.2：假设所求通道幅值为A，初相位为θ，采样后设定为如下波形： 

$X_n=A\sin (2\mathrm{\π}\frac{f_0}{f_s}n+\mathrm{\θ})$ 式一 

步骤3.2.3：取10周波采样数据，加Hanning窗，假定时域形式为ω(n)，其连续频谱为W(2πf),则加窗后采样信号的连续傅里叶变换忽略负频点-f₀处频峰的旁瓣影响，在正频点f₀附近的离散频谱函数可以表达为： 

$\stackrel{\‾}{X}\left(\mathrm{k\Δf}\right)=\left(\frac{A}{2j}\right)e^{\mathrm{j\θ}}W(2\mathrm{\π}(\mathrm{k\Δf}-f_0)/f_s$ 式二 

式二中离散频率间隔为N是截断数据长度,

步骤3.2.4：非同步采样下，峰值频率f₀＝k0*Δf，设峰值左右两侧的谱线分别为第K1和K2条谱线，令由于0≤K0-K1≤1，所以引入一个辅助参数α＝K0-K1-0.5。 

则有： 

$\frac{y2-y1}{y2+y1}=\frac{\left|W(2\mathrm{\π}(-\mathrm{\α}+0.5)/N|-|W(2\mathrm{\π}(-\mathrm{\α}-0.5)/N\right|}{|W(2\mathrm{\π}(-\mathrm{\α}+0.5)/N|+|W(2\mathrm{\π}(-\mathrm{\α}-0.5)/N|}$ 式三 

令通过多项式逼近的方法，得出Hanning窗下有如下修正公式： 

α＝1.5β 

步骤3.2.5：根据公式四、五，计算实际的基波幅值、相位： 

A＝N^-1(y1+y2)(2.356194+1.155436α²+0.326078α³+0.078914α⁶式四 

$\mathrm{\θ}=\arg \left[X(k1*\mathrm{\Δf})\right]+\frac{\mathrm{\π}}{2}-\arg \left[W(2*\mathrm{\π}*\frac{k1-k0}{N})\right]$ 式五 

步骤3.2.6：根据步骤3.2.1所求得的每周波采样点数，取10周波采样数据窗； 

步骤3.2.7：对获取的10周波采样数据取四舍五入值，然后进行DFT运算； 

步骤3.2.8：根据式六求得间谐波子组值： 

$C_{\mathrm{isg},n}=\sqrt{{\mathrm{\Σ}}_{i=2}^8{C}_{k+i}^2}$ 式六 

式六中： 

表示第n次谐波频谱C_K与第n+1次谐波频谱之间不与其直接相连的连续7个间谐波频谱分量。 

综上所述，本发明能够解决数字变电站中由于非同步采样导致电能质量监测装置高次谐波分析不准确的技术难题，为数字化变电站电能质量装置接入提供技术支撑；此外，本发明可以自适应MU输出采样率，在不同采样率下均能计算谐波、间谐波参数，可扩展性强；另 外，本发明通过线性插值、分裂基FFT及重抽样算法，把大规模的DFT运算转成多个小规模FFT计算，使计算量大规模降低，极大的降低了硬件平台要求。如51.2K采样率时(每周波1024点)，利用普通的DFT运算，计算0-50次间谐波需要100*(10*1024)*(50*10*2)*5＝2457600000(4.76G次)，其中，各个数据的含义分别为：100为通道个数；10*1024为数据窗200ms，即10个周波；50*10*2表示间谐波的个数0.1-50.0次，分别计算实部虚部；5表示1秒钟计算5次。而采取本发明方法后，计算量直接下降至：100*256*8*5*5＝5120000(0.005G次)其中，各个数据的含义分别为：100为通道个数；256*8*5表示数据窗200ms，重抽样总计1280点，分为5个序列，每个序列256点,256为2的8次方；5表示1秒钟计算5次。从上述数据计算量可以看出，采用本发明方法后，计算高采样率下谐波时计算量为经典算法的950倍，硬件成本大幅度降低。 

Claims (2)

1.一种数字化变电站电能质量谐波分析方法，其特征在于，所述方法包括：

步骤1：接收数字化变电站合并单元MU发送的实时采样数据x[n]；

步骤2：对接收的实时采样数据x[n]进行软件数字滤波，滤除高次谐波及扰动信号，从而测得系统频率；

步骤3：获取所述数字化变电站采样率，并根据所获取的采样率，相应进行谐波、间谐波分析，得到分析结果；

步骤4：根据步骤3的分析结果，进行电能质量的后续指标运算、界面显示和上送至电能质量分析主站。

2.根据权利要求1所述的一种数字化变电站电能质量谐波分析方法，其特征在于，所述根据所获取的采样率，相应进行谐波、间谐波分析，包括：

当采样率大于等于10K时，进行如下步骤：

步骤3.1.1：采样数据插值，将非同步采样转换成同步采样：把采样数据插值成临近的2的N次方次采样点；

步骤3.1.2：根据插值后的采样率，判断采样率是否过高，过高时进行抽样采样，如当采样率大于等于25.6K时，抽样成12.8K采样率数据，降低运算复杂度；

步骤3.1.3：谐波计算：利用FFT算法计算整数次谐波的谐波值及谐波相角；

步骤3.1.4：间谐波计算：将数据窗重采样转变成多个FFT计算，将x[n]进行重采样，得到5个采样序列x[5n]、x[5n+1]、x[5n+2]、x[5n+3]、x[5n+4]，并对所述5个采样序列进行FFT运算，得到间谐波的谐波值及谐波相角；

当采样率小于10K时，进行如下步骤：

步骤3.2.1：根据步骤2测得的系统频率f₀，计算出每周波实际采样点数；

步骤3.2.2：假设所求通道幅值为A，初相位为θ的单一频率信号X(t)，在经过了采样率为f_s的模数变换后得到如下形式的离散信号：

$X_n=A\sin (2\mathrm{\π}\frac{f_0}{f_s}n+\mathrm{\θ})$ 式一

步骤3.2.3：取10周波采样数据，加Hanning窗，假定时域形式为ω(n)，其连续频谱为W(2πf),则加窗后采样信号的连续傅里叶变换忽略负频点-f₀处频峰的旁瓣影响，在正频点f₀附近的离散频谱函数可以表达为：

$\stackrel{\‾}{X}\left(\mathrm{k\Δf}\right)=\left(\frac{A}{2j}\right)e^{\mathrm{j\θ}}W(2\mathrm{\π}(\mathrm{k\Δf}-f_0)/f_s$ 式二

式二中离散频率间隔为N是截断数据长度,

步骤3.2.4：非同步采样下，理论峰值谱线位置K0＝f₀/Δf，其中Δf为DFT变换时的频谱分辨率，由于f₀并不是Δf的整数倍，因此实际K0谱线是不存在的，设峰值左右两侧的谱线分别为第K1和K2条谱线，令由于0≤K0-K1≤1，所以引入一个辅助参数α＝K0-K1-0.5。

则有：

$\frac{y2-y1}{y2+y1}=\frac{\left|W(2\mathrm{\π}(-\mathrm{\α}+0.5)/N|-|W(2\mathrm{\π}(-\mathrm{\α}-0.5)/N\right|}{|W(2\mathrm{\π}(-\mathrm{\α}+0.5)/N|+|W(2\mathrm{\π}(-\mathrm{\α}-0.5)/N|}$ 式三

令通过多项式逼近的方法，得出Hanning窗下有如下修正公式：

α＝1.5β

步骤3.2.5：根据公式四、五，计算实际的基波幅值、相位：

A＝N^-1(y1+y2)(2.356194+1.155436α²+0.326078α³+0.078914α⁶式四

$\mathrm{\θ}=\arg \left[X(k1*\mathrm{\Δf})\right]+\frac{\mathrm{\π}}{2}-\arg \left[W(2*\mathrm{\π}*\frac{k1-k0}{N})\right]$ 式五

步骤3.2.6：根据步骤3.2.1所求得的每周波采样点数，取10周波采样数据窗；

步骤3.2.7：对获取的10周波采样数据取四舍五入值，然后进行DFT运算；

步骤3.2.8：根据式六求得间谐波子组值：

$C_{\mathrm{isg},n}=\sqrt{{\mathrm{\Σ}}_{i=2}^8{C}_{k+i}^2}$ 式六

式六中：

表示第n次谐波频谱C_K与第n+1次谐波频谱之间不与其直接相连的连续7个间谐波频谱分量。