庞吉耀
（南京磐能电力科技股份有限公司，江苏南京210032）
摘要：为寻求一种频率无关的实时相量计算方法，通过分析传统算法，提出一种可变数据窗正弦相量插值计算方法。该算法通过等间隔采样获得信号一个周期左右的样本数据，采样数据计算出信号的基波频率并确定信号的一个基频周期所需的最大采样点数N，计算信号N 点DFT和N+1点DFT，再利用两次DFT计算结果进行线性插值得到信号真实频点上各次谐波的实部和虚部，计算出信号各次谐波的幅值和瞬时相位。仿真分析和实际测量表明，该算法的计算精度和实时性较高，能够满足频率随机变化的应用要求。
教育期刊网 http://www.jyqkw.com
关键词 ：离散傅里叶变换；电力系统；相量计算；插值计算
中图分类号： TN911.23-34；TM930.1 文献标识码：A 文章编号：1004-373X（2015）12-0005-03
收稿日期：2014-12-14
0 引言
当电网处于工频时，基于定间隔采样技术的DFT算法具有良好的谐波滤波特性，测量结果十分精确，但当电网频率偏离50 Hz时，由于非同步采样带来频域泄漏导致传统测量算法难以同时满足计算量小、跟踪速度快和计算精度高等要求。就目前智能变电站的数据体系结构而言，保护装置和测控装置的数据都源于合并单元，且合并单元向数据订阅方输出的数据都是定间隔采样的数据，因此以往通过调整本地ADC 采样间隔的频率跟踪算法不适用智能变电站，针对全数字化的保护测控设备需要采用新的计算方法。
为提高传统DFT算法精度，目前已提出了多种改进方法来修正传统的DFT算法，减少频谱泄漏的影响。如使用加窗函数和内插方法来减小泄漏误差[1-4]，非同步采样点的同步化[5]和动态调整采样率实现同步采样的方法等，其中加窗和内插法通过选择合适的窗函数对采样数据进行加权，再对加窗后的序列进行FFT运算并进行双谱线插值，非同步采样点同步化法则是通过对原始采样序列进行同步化插值得到新的同步采样序列，而动态调整采样率方法需要采用锁相环实现硬件同步采样来减小由非同步采样所导致的截断误差，提高计算精度。
以上几种方法都可以在一定程度上减小频谱泄漏，但除去在中间的步骤推导时可能会出现一些舍入误差，计算量比较大外，还因为使用多个信号周期的数据导致算法实时性下降，限制了其应用。而无论是同步插值还是基于硬件同步采样的方法都导致采样间隔的不稳，无法利用现有的数字滤波技术对信号进一步处理。
本文根据传统傅里叶变换的原理，提出了一种快速频率跟踪算法，该算法首先测得信号的基波频率，据此确定信号的一个基频周期所需的最大采样点数N，然后计算信号N 点DFT和N+1点DFT，再利用两次DFT计算结果进行线性插值得到信号真实频点上各次谐波的实部和虚部，最后计算出信号各次谐波的幅值和瞬时相位。仿真分析和实际测量表明，该算法的计算精度和实时性较高，能够兼顾动态响应和精度，是同类算法中数据窗长度最短的，能够满足频率随机变化的应用要求。
1 算法的理论基础
依傅里叶变换原理，若信号的实际周期为T ，则其实部和虚部可表示如下：
X Tr = 1/T ∫0T x(t)cos ωtdt （1）
X Ti = 1/T ∫0T x(t)sin ωtdt （2）
式中：X Tr 是信号周期为T 时傅里叶算法所得实部；X Ti是信号周期为T 时傅里叶算法所得虚部。若以间隔Ts对信号进行等间隔采样，则必然存在整数N 使式（3）成立：
NTs≤T≤(N + 1)Ts （3）
按式（1）和式（2），当信号周期在[NTs,(N + 1)Ts] 之间连续变化时，算法的实部X Tr 和虚部X Ti 将分别构成一段连续曲线。根据拉格朗日插值定理，在误差允许的情况下，基频周期为Tx 信号的实部和虚部实际值可以分别由其端点处的值通过线性内插获得，以实部为例推导如下：
依离散傅里叶变换（DFT）的定义，将式（5）~式（8）离散化，使用离散求和替代积分，可以分别得到N 点DFT 的实部和虚部。
式中：X Nr 为N 点数据窗的实部；X Ni 为N 点数据窗的虚部。
类似地，使用信号的N + 1 个点的数据，估计实部和虚部可得如下结果：
根据式（13）和式（14）计算出的实部和虚部可以进一步利用复数取模等方法获得幅值和相角，而闭区间[NTs,(N + 1)Ts] 的大小将决定插值的误差。
当满足同步采样条件时，Tx Ts 必为整数，此时VTx= 0 ，则X Txr = X NTsr ，X Txi = X NTsi ，即满足同步采样条件时新算法等同于传统DFT 。不难验证，对于直流分量式（14）的计算结果为0，而式（13）的计算结果将等于直流分量且不损失精度。
2 算法实现
2.1 最大周期采样点数N 计算
本文算法需要先测得出信号周期，可以先计算信号频率后再导出周期，也可以通过采样序列中过零点的位置来获得信号周期[4，6-7]。由于在小信号时会导致周期或频率测不准，为保证算法稳定应根据实际情况限定频率或周期的范围，针对电力应用可限定在45~55 Hz。设信号基频周期为Tx（对应频率fx = 1/Tx），采样间隔为Ts（对应采样率fs = 1/Tx），则每周期得到的采样数据点数为：
Nx = Tx/Ts = fs /fx （15）
算法按下面的方法确定每个基频周期最大整数采样点数N ：
N = INT(Tx /Ts) = INT( fs /fx), N ≤Nx ≤N + 1 （16）
式中函数INT(x) 表示向下取整。显然满足同步采样条件时Nx = N 为整数。
2.2 应用举例
以市电交流信号为例，以采样率fs =3 200 次/s（Ts = 312.5 μs）对信号进行等间隔采样时，算法操作步骤如下：
（1）计算信号的当前频率，此处假设所得频率为fx =49.5 Hz；
（2）获得信号每周期最大整数采样点个数及所需序列长度fs /fx = 3 200/ 49.5 = 64.646 4 ，向下取整得N = floor( fs fx) = 64；
（3）从当前采样点向前，分别取64 点和65 点数据构建数据窗，并计算64点和65点的DFT；
（4）按式（13）和式（14）计算信号基波分量的实部和虚部，并由此计算幅值和相角；
（5）用类似的方法计算所需要的谐波实部及虚部，并导出谐波的幅值和相角。
3 应用效果分析
3.1 频域插值间距
根据数字信号处理原理，取采样序列N 点数据计算DFT时获得的基波谱线位于频点1 NTs 处，而当取采样序列N + 1 点数据计算DFT时获得的基波谱线位于频点1 (N + 1)Ts 处，以采样率fs = 3 200 为例，本文用来插值的两个谱线之间间距为：
可见本算法的谱线间距要比现有算法[8]小一个数量级，而式（17）表明随着采样率的增加，谱线间距会更小，因此仅采用的最简单的线性内插就可以获得满意的算法精度。而目前的算法由于频域插值间距大，所以不能采用线性插值，需要更复杂的插值方法。
3.2 开销分析
目前的同类算法为了达到所需要的精度至少需要两个周期以上的数据[1，8]，也即动态响应能力更差。与目前的频率修正DFT 算法相比，本文需要的存储开销最小，在同步采样条件下只需要保留信号的一个周期采样数据，即使在非同步采样条件下也只需要额外增加一个采样数据。因此本算法在精度和动态响应方面取得比较好的平衡。此外本文采用的是矩形窗直接截断，同目前的采用其他窗函数的算法相比没有加窗计算过程，而且采用简单线性插值较采用窗函数的频谱函数和频率偏移校正量来校正信号相量的开销小。
3.3 误差分析
为分析该算法的误差，构建幅度为1且含有基波和2，3，5谐波正弦序列模型如式（18），并以频率和初相角为变量，采样间隔Ts 为参变量，在Matlab中利用上述算法计算信号的幅度和相位与理论值进行误差：
x(n) = sin(2πfx nTs + -1) + sin(4πfx nTs + -2 ) +sin(6πfx nTs + -3) + sin(10πfx nTs + -5) （18）
式中：Ts 为采样间隔；fx 为信号频率；-1,-2 ,-3,-5 为谐波初相角。
首先固定采样率为3 200，然后以频率和初相角为参变量进行误差扫描（其中频率步距0.1 Hz，相角步距1°），结果以频率为x 轴，初相角为y 轴，而幅值和相位误差为z 轴绘制三维误差图形，如图1所示，算法本身精度能够满足测控装置要求。
随后改变采样率，再次以同样的频率和相位步距在Matlab中进行误差扫描。在电力系统常用的几种采样率下的最大误差数据如表1所示，结果表明即使在采样率3 200次/s的条件下依然可以满足电力系统对测控装置的精度要求。
4 结语
本文从满足保护的快速性和测控精度的行业要求出发，采用最短矩形窗截取采样序列的一个周期数据，通过变长度DFT和频域线性插值实现频域泄漏修正，与同类方法相比具有明显的优势：计算仅需要一个周期的采样数据，内存消耗少；单周期即可获得信号的相量值，比同类算法动态响应特性好；通过增加采样率（减小采样间隔），可以增加计算精度但不影响动态特性；计算结果不受信号频率的影响；谱线fxl 和fxh 的间距（对应频域分辨率）在同等采样间隔下是同类算法中最短的。
通过仿真和实际测量表明，本文给出的算法精度和快速性完全可以满足继电保护要求，而算法精度也满足测量要求。目前该算法已成功应用在我公司的DMP5500智能变电站保护测控设备和DMP5000微电网保护控制系统中，且获国家知识产权局的发明专利授权[9]。
作者简介：庞吉耀（1972—），男，江苏淮安人，高级工程师。研究方向为电力用嵌入式平台和嵌入式系统。
教育期刊网 http://www.jyqkw.com
参考文献
[1] 黄纯，江亚群.谐波分析的加窗插值改进算法[J].中国电机工程学报，2005，25（15）：26-32.
[2] 李新福，刘教民，崔玉龙.电器检测的谐波分析算法[J].河北工业大学学报，2003（5）：23-26.
[3] 庞浩，李东霞，俎云霄，等.应用FFT 进行电力系统谐波分析的改进算法[J].中国电机工程学报，2003（6）：50-54.
[4] 梁志国，孙璟宇.信号周期的一种数字化测量方法[J].仪器仪表学报，2003，24（z1）：195-198.
[5] 曾泽昊，余有灵，许维胜.基于插值同步算法的频谱泄漏分析与仿真[J].系统仿真技术，2005（3）：147-152.
[6] 刘建林，陈兵.一种高精度单频信号频率估计算法[J].无线电工程，2011（4）：26-28.
[7] 吴峻，吴崇昊.基于Kaiser 窗的改进傅氏测频算法[J].江苏电机工程，2008，27（4）：17-19.
[8] 高云鹏，滕召胜，卿柏元.基于Kaiser 窗双谱线插值FFT 的谐波分析方法[J].电子学报，2000（2）：287-292.
[9] 庞吉耀.一种频率无关的正弦信号相量计算方法：中国，201310055886[P].2013-05-09.
[10] 闵勇，丁仁杰，韩英铎，等.自适应调整采样率的相量在线测量算法研究[J].电力系统自动化，1998（10）：10-13.
[11] 杨洪耕，惠锦，侯鹏.电力系统谐波和间谐波检测方法综述[J].电力系统自动化，2010（2）：65-69.