应用于电力系统谐波分析的采样周期自适应方法

周 洲

(成都电业局，四川 成都 610021)

电网中大量非线性负荷的应用导致谐波污染日益严重，对谐波成分的准确分析有利于对电能质量进行合理评估。离散傅里叶变换(discrete Fourier transform，DFT)和快速傅里叶变换(fast Fourier transform，FFT)因适用于计算机计算而被广泛应用，但是电网频率波动会导致信号采样的不同步，此时会出现影响测量精度的频谱泄漏问题［1］。

很多种方法已经可以降低由于非同步采样中频谱泄漏引起的测量误差。在频域，可以选择适当的窗函数来降低旁瓣带来的长范围泄漏，并通过对频率、相位和幅值进行插值的手段减少短范围泄漏［2-10］。

在时域，可以在截断数据后补零来降低栅栏效应［1］，通过波峰检测分析谐波参数，但此方法克服不了谐波间的频谱干扰且计算量较大;时域同步插值通过对采样序列进行拉格朗日插值近似实现序列的同步采样［11-14］，此方法需要知道信号的实际频率。文献［15］提出了一种基于相关计算的采样窗口同步化方法，该方法无需知道信号的实际频率，通过序列迭代计算使得采样序列长度自动满足同步化收敛条件。但是在实际应用于稳态谐波计算时，由于该方法采样周期在迭代过程中保持恒定，且实际系统中电网频率在一定范围波动，这样会使得算法迭代时间长，同步化序列长，并且最终的同步化序列无法采用FFT算法。因此，提出了一种新的采样周期自适应方法用于谐波计算，通过计算固定长度相邻信号序列的相关系数并结合提出的二次修正方程，不断迭代调整采样周期，得到满足FFT算法要求的同步化采样序列。最后通过FFT计算得到各次谐波幅值相位。仿真算例证实了该方法的有效性。

1 基于序列相关计算的采样窗口同步化方法［15］

连续信号等间距离散化后得到有限长的离散时间序列，表示为信号向量 S={S0，S1，…SN-1}。两个长度相同的信号向量S1和S2的相关计算为这两个向量的内积。

其中，θ є［0，π］，表示两向量之间的夹角。

定义两离散信号序列的相关系数如下。

Kp的取值范围为［-1，1］，Kp值用来指示两个离散信号之间的相似程度。Kp=1表示两个离散信号相位一致;Kp=-1表示两个离散信号相位相反。

基于相关系数的采样窗口同步化计算步骤如下。

考虑一离散信号S

1)选择p的初始值，确定信号Xp为迭代的初始序列，其中 Xp={S0，S1，…Sp-1}。

2)定义一个与Xp长度相同的后继信号Yp={Sp，Sp+1，…S2p-1}。

3)计算Xp与Yp的相关系数Kp。

4)视察Kp是否满足收敛条件(Kp=1)，若不满足，令p=p+1，重复步骤1)、2)和3)直到序列满足收敛条件。

在此也可以称Kp为周期指示器。可以通过考察Kp与1的接近程度来估计采样窗口的同步化程度。当Kp=1时，表示序列Xp和Yp完全相同，采样已经实现了同步化;对序列Xp进行DFT便可以精确得到信号各次谐波参数，此时频谱泄漏的影响为零。Kp≠1时则继续进行迭代计算直到满足收敛条件为止。在实际应用中，若要求对各种信号都满足收敛条件Kp=1，则会使得采样序列非常长。一般情况下根据精度要求以ε=|Kp-1|≤εset作为收敛条件。

该方法的优点在于，不需要预先知道信号各成分的频率，采样窗口可以通过基于相关计算的迭代方法自适应求得，从而实现信号的整周期截断与同步采样。很显然，此方法在信号非平稳的情况下是不适用的。往往满足收敛条件的同步化序列非常长，在固定采样周期的非同步采样下尤为明显。

2 基于相关系数的采样周期自适应

上节所述的信号周期自适应方法可以在不知道信号频率的情况下对其实现同步采样，但是其存在下面两个问题:①要满足收敛条件ε=| Kp-1|≤εset的采样序列较长，迭代运算次数较多且DFT的运算量较大;②通过迭代计算得到采样序列的点数pmax无法保证满足pmax=2M(M为正整数)，以至在一般情况下对采样序列不能使用FFT算法。

现考虑固定每周期采样点数N=2M。初始采样周期为Ts0=Nf0(f0表示信号没有频偏时的频率，对于50 Hz系统，f0=50 Hz)。

通过计算XN与YN(YN为XN的后继N点序列)的相关系数，根据修正方程，不断地修正信号的采样周期，使得信号满足收敛条件实现同步采样。

计算中相关系数越是接近1，表示此时的XN与YN越接近，采样周期的修正量就应该越小;反之，相关系数越是远离1，表示XN与YN相似程度越差，采样周期的修正量就应该越大。这样可以加快信号的收敛速度。

设第i+1次采样周期修正为

ΔNi称为采样周期修正因子，i=0，1，2…。

根据上述要求，采用二次抛物线的修正方程为

其中，AN为迭代加速因子;N为采样点数;KNi表示第i-1次序列XNi与YNi的相关系数。

当N=64时修正方程的函数图形如图1所示。

图1 采样周期修正因子ΔNi的修正方程函数图

从图1可以看出，相关系数越是接近1，修正量就越小，这样使得收敛更快且保证信号序列由快到慢接近收敛点。在某种程度上来讲，加速因子AN的取值越大，收敛速度越快。

基于相关系数的采样周期自适应具体步骤如下。

1)根据Ts0采样得到2N点的信号序列，前N个点序列为XN0，后N个点序列为YN0，根据式(2)计算得到相关系数KN0。

2)根据修正方程ΔN0=AN和式(3)分别计算

3)根据步骤2)分别得到两个采样周期值，采样并计算KN1和KN2。若KN1＞KN2，说明信号频率偏大，式(3)分母中的符号在迭代过程中定为“+”号;若KN1＜KN2，说明信号频率偏小，定为“-”号。

4)根据步骤2)得到的采样周期对信号继续采样，计算 KNi。

5)根据修正方程ΔNi=AN，和式(3)计算得到第i+1次的采样周期。

6)根据新的采样周期继续采样2N点的信号序列，计算KNi，判断是否满足收敛条件。满足则结束计算，不满足则返回步骤5)继续迭代计算，直到满足为止。

最后，对满足收敛条件的最新采样序列做FFT计算得到信号各次谐波的幅值和相对相位(由于采样周期的不断调整使得各次谐波相位发生变化，相位计算采用以基波相位作为零度基准值)。

3 数字仿真计算

含有各次谐波的电力系统信号表示为

其中，m为正整数;Am和φm分别表示m次谐波的幅值和相角，具体参数如表1所示。

表1 谐波参数设置

计算1:对两种方法的计算量进行比较。

首先，对信号进行基于采样窗口同步化的迭代计算［15］，采样频率为固定值3200 Hz，初始迭代序列Xp的p=60，收敛标准Kp＞0.99999。在不同频率偏差的情况下得到收敛所需的迭代次数、加法和乘法运算次数如表2所示。再对信号进行基于采样周期自适应的迭代计算，此时，初始采样频率fs0=3200 Hz，收敛标准 Kp＞0.9999999999。每周期采样点数恒是N=64，通过迭代计算不断地修正采样周期最终使得采样序列达到收敛条件，具体的运算量亦如表2所示。

从表中看出，后者的迭代次数及运算量要明显小于前者。采用基于采样周期自适应的迭代计算时，随着加速因子的增加可以有效地提高收敛速度。

计算2:对谐波计算精度进行数字比较。

假设信号实际频率f=49.5 Hz，对式(4)信号分别采用以下方法进行谐波计算。方法1，直接采用FFT算法计算;方法2，采用基于信号周期自适应的方法计算［15］;方法3，采用基于采样周期自适应的方法计算。方法1的采样频率fs=3200 Hz，窗口长度为64点。方法2和方法3的采样频率、收敛条件等与计算1时相同。

计算结果如表3所示，从各次谐波的分布可以看出，信号中的偶次谐波含量要远比奇次谐波小。

现令

作为测量精度指示因子。ΔAm和Δφm为m次谐波幅值和相位的绝对误差，由于基波和8次谐波存在零值，在此不参与α因子计算。最后单独做绝对误差比较。

表2 两种方法运算量的比较

根据表中数据计算，采用方法1时，奇次谐波的α因子为

偶次谐波α因子为

由于非同步采样，尤其是对于含有量较小的偶次谐波而言测量误差较大。

采用方法2时，奇次、偶次谐波精度指示因子分别为 α奇次=0.000862221和 α偶次=0.0014。其精度和FFT相比有着明显地提高。

采用方法3时，根据最后一列数据，得到奇次、偶次谐波精度指示因子分别为 α奇次=0和α偶次=0.000016766。并且从表中可以看出收敛标准的提高导致测量精度的提高，但在相同加速因子的条件下会增加迭代次数。与方法2相比，增加后的运算量仍甚小，并且精度要远高于方法2。

对于基波，方法3的测量误差为零值。方法1、2都有一定程度的误差存在。对于8次谐波，方法3计算得到的幅值误差要远远小于其他方法，但是其相位误差较大。

在采用方法2时，由于信号的频率未知，导致满足收敛标准的信号序列很长，且不能保证满足FFT的采样点数要求(如表中的点数为1099点)，这样通过DFT计算谐波参数时计算量就会更加大。

4 应用讨论

4.1 加速因子的取值范围

加速因子AN的大小决定了计算迭代的初始步长。从这里的仿真结果可以看出，加速因子AN的取值越大，迭代计算的初始步长越长，算法收敛速度也就越快，谐波计算的实时性也就越高。但是，对于一个工频周波信号，若盲目提高加速因子取值，可能会引发计算收敛序列变成若干个周期，从而影响计算实时性和DFT计算量。因此，加速因子AN应根据实际信号的频偏特性在一起范围内取值。这个取值范围取决于信号频率偏差范围。

设信号的频率偏差为Δf'，理想的采样周期应为

表3 3种方法谐波参数计算结果

根据式(8)可以计算得到

在采样周期自适应的迭代过程中，步长ΔNi应小于等于ΔN'。另由于迭代步长是单调递减的，因此第一次迭代计算的值应为最大ΔNmax，必须满足

得到AN的取值范围为

根据电力系统的运行规范和标准，电力系统的频率往往被控制在50 Hz±0.5 Hz之间，因此信号的频率偏差是有上限的，根据这个上限，可以确定加速因子AN的取值上限。取值为0.5 Hz。

4.2 收敛条件

计算结果的精度取决于用户所定的计算收敛条件，即为:ε=| Kp-1|≤εset。若Kp取值越高，则计算误差越小;反之，则计算误差越大。但提高计算精度的代价是影响计算的实时性，因此，在实际应用过程中，Kp的取值可以根据测量的精度要求、计算实时性和环境的噪声水平等综合考虑后取值，在先保证信号实时性的条件下，即采样的时间必须小于某个数tmax的条件下，取相关系数最高的那组周波数据进行DFT计算，得到相应的谐波值。

5 结论

采用采样周期自适应的方法可以有效解决由于非同步采样引起的频谱泄露问题。具体有以下几个特点。

1)该方法无需知道电网信号的实际频率，随着电网频率的偏移，迭代计算可以不断适应频率的动态偏移实现同步采样。

2)迭代过程中通过对采样周期进行修正，并引入修正加速因子可以提高信号序列收敛速度，即使在频率偏差较大时也能即时使序列很快收敛。

3)在迭代过程中采样点数恒定，每次迭代的计算量相同;而采样周期固定时，每次迭代的计算量随着迭代次数的增加而增加。当迭代次数较多时计算量很大。

4)采样点数满足N=2M，因此可以采用FFT计算，提高运算速度。

5)实际分析时也可用时域多项式插值手段对采样信号进行采样周期调整的重采样。

［1］Ignacio，Santamaria-Caballero，Carlos J，Pantaleon-Prieto，Jesus Ibanez，et al.Improved Procedures for Estimating Amplitudes and Phases of Harmonics with Application to Vibration Analysis［J］.IEEE Transactions on Instrumentation and Mesurement，1998，47(1):209-2144.

［2］Dusan Agrez.Weighted Multipoint Interpolated DFT to Improve Amplitude Estimation of Multifrequency Signal［J］.IEEE Transactions on Instrumentation and Mesurement，2002，51(2):287-292.

［3］D.Belega，D.Dallet.Frequency Estimation via Weighted Multipoint Interpolated DFT［J］.IET Science，Measurement＆ Technology，2008，2(1):1-8.

［4］潘文，钱俞寿，周鹗.基于加窗差值FFT的电力谐波测量理论(I)窗函数研究［J］.电工技术学报，1994，9(1):50-54.

［5］潘文，钱俞寿，周鹗.基于加窗差值FFT的电力谐波测量理论(II)双插值 FFT理论［J］.电工技术学报，1994，9(2):53-56.

［6］张伏生，耿中行，葛耀中.电力系统谐波分析的高精度FFT 算法［J］.中国电机工程学报，1999，19(3):63-66.

［7］孙红伟，李梅，袁健华，等.用于电力系统谐波分析的加窗FFT 算法研究［J］. 高电压技术，2004，30(8):52-55.

［8］钱昊，赵荣祥.基于插值FFT算法的间谐波分析［J］.中国电机工程学报，2005，25(21):87-91.

［9］赵文春，马伟明，胡安.电机测试中谐波分析的高精度FFT算法［J］.中国电机工程学报2001，21(12):83-87.

［10］柴旭峥，文习山，关根志，等.一种高精度的电力系统谐波分析算法［J］.中国电机工程学报，2003，23(9):67-70.

［11］陈隆道，钱照明，张圣训.周期域分析中非同步取样数据的同步化［J］.电子学报，2001，29(7):1-4.

［12］曾泽昊，余有灵，许维胜.一种减小泄漏的同步化算法［J］.电测与仪表，2005(479):12-14.

［13］曾泽昊，余有灵，许维胜.基于插值同步算法的频谱泄漏分析与仿真［J］.系统仿真技术，2005，1(3):148-152.

［14］惠锦，杨洪耕.改进的时域同步插值算法在基于TDA谐波和间谐波测量中的应用［J］.电工技术学报，2008，23(7):124-129.

［15］T.X.Zhu.Exact Harmonics/Imterharmonics Calculation Using Adaptive Window Width［J］.IEEE Transactions on Power Delivery，2007，47(1):2279-2288.

［16］郑恩让，杨润贤，高森.关于电力系统FFT谐波检测存在问题的研究［J］. 继电器，2006，34(18):52-57.

［17］布赖母，柳群译.快速傅立叶变换［M］.上海:上海科学技术出版社，1979.