电力信号同步采样算法

张 红

（桂林师范高等专科学校数学与计算机科学系，桂林 541001）

利用离散傅里叶变换DFT（discrete Fourier transform）对信号进行分析，具有正交、完备、快速等许多优点，它可以精确地确定出平稳波形中各次谐波的幅值。它的理论基础是假定输入信号是周期信号，可以分解为恒定直流分量与整数倍基频周期分量之和，因此利用正交三角函数基或正交指数基把时域信号变换成频域信号，分离出电参量信号的基波及各次谐波分量，从而可得到信号各分量的幅值、频率和相位，然后再由相应公式可以很方便的求出其余的参数值。在短时傅里叶变换、小波变换、现代谱估计和二次变换等理论和方法还没有得到充分完善之前，基于DFT的谐波分析方法在谐波分析和检测领域中仍将发挥其重要的作用［1］。

但DFT存在一些固有的限制，如：取样信号必须是周期性；取样信号必须为基波信号的整数倍；取样速率必须高于信号最高频率的2倍以上。当所处理的信号不符合这些条件时，则将产生栅栏效应、泄漏效应、混叠效应［2，3］。若采样频率与信号基频不同步时，截取并周期延拓而成的信号将不等于原信号，会造成波形间断，这样必然会使变换结果偏离真实值，出现频谱泄漏现象。减少频谱泄漏的方法有加窗法［5～10］、非同步采样法［11，12］，频谱校正［13，14］等。

在观察电力系统的电压电流信号过程中发现其频谱具有以下特点：属于离散频谱，频率在50 Hz左右［4］的基波及频率为基波整数倍的谐波和频率为基波非整数倍的间谐波构成，基波占主要成分，一般奇次谐波分量较小，偶次谐波分量和间谐波分量的含量更少。电力信号采集过程中往往是设定一个固定的采样速率在对其进行A／D采样，这样得到的采样值就会有可能与电力信号的基频不成整倍数关系。

因此，为了保证离散电力信号与基波信号严格的成整倍数关系，需要对A／D采样得到的离散信号用一个与基频成整倍数的采样速率进行转换。而电力信号的基频应尽可能的精确估计，才能对其进行采样速率的转换。

本文对A／D后的离散电力信号运用软件无线电技术中精细频率估计的方法来估计和跟踪其频率。实验发现利用这个算法估计离散电力信号的频率可以达到0.001数量级的精度。估计得出电力信号的真实信号频率后需要对原来的离散信号进行采样率转换。本文对离散电力信号采样速率转换运用了改进Farrow算法查表滤波器系数的方法来实现。其中Farrow滤波器的单位脉冲响应系数运用了MATLAB中的Fdatool工具设计。本文仿真验证了该算法的正确性和有效性，该方法不仅能准确地估算出电力信号的频率，且精度和稳定度好，计算量较小、速度快。

1 电力信号频率的估计

1.1 电力信号及其傅里叶变换

设稳态电力信号如式子（1）所示

式中，fi、Ai、φi分别为信号各分量的频率、幅值和相位参数。

国家标准的电力信号频率是50 Hz，而谐波干扰和其他因素干扰所造成的频率偏差使电力信号频率出现波动，其中正常范围是（50±0.2 Hz）［4］。所以假设电力信号的工频频率f＝50 Hz，那么一个信号周期T＝1／f＝20 ms。以采样率fs＝6400 Hz对信号x（t）采样，每个周期采样长度是N＝128点。

采样信号的时域表示为

式中：n＝0，1，2，…，N；Ts为采样时间间隔，Ts＝1／f s。

分别连续取两段长度都为N＝128点的数据得到xm（n）和xn（n），n＝0，1，2，…，N－1。然后对xm（n）和xn（n）做离散傅里叶变换，计算变换公式为

对xm（n）和xn（n）分别做DFT变换后得到Xm（k）和Xn（k）。

1.2 精细频率分辨率

通过分析信号相位关系可以得到精细频率分辨率。如果一段数据中的最高频率分量在m时刻为Xm（k），k为输入信号的频率分量。从DFT输出中得到输入信号的初始相位为

式中，Im和Re分别表示复数的实部和虚部。

因为在短时间内输入频率变化不大，假设在m时刻之后间隔不长的n时刻，1 ms长数据的DFT分量Xn（k）也是最大分量。在n时刻、k频率分量时，输入信号的初始相角为

利用这两个相位角可得到精细频率为

由式（6）得到的频率分辨率比DFT得到的频率分辨率更精确。为了保持这个频率不模糊，相位差θn（k）－θm（k）必须小于2π。当相位差为2π时，不模糊带宽是1／（n－m），其中，n－m是两组连续数据之间的延迟时间［15］。

1.3 电力信号频率校正

用工频频率f＝50 Hz，与精细频率fΔ相加得到电力信号的新的频率为

这样就可以得到一个更为逼近电力信号频率的频率，可以对接着进来的电力信号x（t）再取两个周期使用采样率

式中，N＝128。

进行新一次采样得到两段新的信号，重新做DFT变换，利用式（6）～式（8）得到新的采样率，这样不断重复这个过程就可以不断逼近实际的电力信号频率，就可达到很高的采样同步。

2 电力信号的同步采样

估计得到新的采样率后，就对之前由A／D得到的原采样率的离散信号进行采样率转换。这里运用对Farrow滤波器改进的方法实现离散信号的采样率转换。

2.1 改进的Farrow滤波器

假设一组数字信号序列x（nTs），Ts为采样间隔，发送方采样率为fs＝1／Ts。若接收方采样率为fy＝1／Ty与发送方不一致，则发送方与接收方必须进行采样率转换，转换公式［16，17］为

对于任意因子的采样率转换，设mTy＝（km＋Δ）Ts，0≤Δ＜1，km为整数，则式（9）可以写成

简化为

令k＝km－n，可以得到

2.2 滤波器的插值运算

用P阶多项式，如泰勒展开式、拉格朗日多项式来逼近h（k＋Δ）［17］，其展开式为

本文滤波器设计为P＝7，采样得到N＝8个系数，用h（n）表示，对这8个系数做DFT变换得到

为了得到任意的采样率下的离散信号，对H（k）做插值处理，让它逼近连续信号，以方便实现得到任意数值下的采样率。根据离散逆傅里叶变换公式

在这里为了得到插值的效果，n的取值可以是任意步长的小数，也可以是整数，把式（15）展开得

插值后需要考虑h（n）的误差［18］，即

在这里从500个点到10000个点每隔100点做插值实验发现采用2000个点的就可以达到所要求的信号恢复效果，其中插值2000点和10000点相减的差值，见图1。

从图1中可以发现这个差值非常小，说明差值2000点和10000点相差微小。把2000点的插值保存为文件形式，以作为信号重构时候查表用。

图1 h（n）插值2000点和10000点相减的结果示意Fig.1 Subtraction result of h（n）interpolation of 2000 points and 10000 points

2.3 信号的采样速率转换

在得到新的比较合适的采样率h（n）的系数表之后，从硬件A／D得到的离散的电力信号x（t）的采样率按照为工频频率f＝50 Hz的整数倍来采样得到离散信号x（n），x（n）对h（n）进行查表卷积，就可以得到新的采样率的电力信号，即

3 仿真分析及实践

在PC机上进行实验，仿真软件为MATLAB7.0。在编程过程中把采样率归一化为fs＝1，固定每个电力信号周期采样点数N＝128，采样点间隔Ts＝1／fs，归一化因子M＝128×50，基波频率f0＝50／M Hz，其波动范围最大是0.2 Hz，所以构造仿真的电力信号频率f1＝50.2／M，编程的数据精度用long。

实验考虑3种信号环境：没有谐波和噪声、有3次、5次和7次谐波以及添加了谐波和20 dB的高斯白噪声的频率估计，运用精细频率估计方法，得出结果如表1所示。

表1 精细频率的估计仿真结果Tab.1 Simulation result of fine frequency Hz

从表1结果可以看出，在没有谐波和噪声的时候信号经过第一次的估计，频率估计值就达到了0.001数量级精度，经过3次估计就完全收敛到50.200000 Hz；同时添加了谐波和噪声对频率估计的精度影响不大，这两种情况都是经过3次循环估计就完全收敛。验证了该算法应用在电力信号的估计以及跟踪中的有效性和可靠性，为下一步的信号采样率转换提供了准确的与基波成整数倍的采样频率，从而为运用FFT做谐波分析提供可靠的采样信号。信号的采样速率转换没有直接利用h（n）的表达式，而是利用MATLAB7.0中的工具Fdatool设计h（n）的系数，h（n）的参数为截止频率fc＝6400 Hz，采样率fs＝6400×4 Hz，阶数为7阶，选择Bartlett－hanning窗函数设计。

为了在查表实现对信号采样率转换过程当中达到更好的效果，对h（n）的参数利用逆DFT公式插值达到7×300个点，根据DFT的对称性在利用IDFT做插值的过程中取变换的长度为N＝3.5，插值效果如图2所示。

图2 插值成7×300个点的效果Fig.2 Diagram of 7×300 points

利用这个插值成7×300个点的h（n）对一个信号做D／A转换，其中0至100点的效果如图3所示，具有非常理想的采样速率转换效果。

图3 对信号重新采样后的效果Fig.3 Diagram of the re－sampling signal

设原始信号的数学式子为

信号中含有的最高谐波成分的频率为350 Hz，对以上信号设置采样频率fs＝128×50 Hz进行采样，并进行傅里叶变换得到如图4所示，其中只显示了原始信号的0～500点以及一个周波的频谱。

图4 对信号重新采样后的FFT频谱Fig.4 FFT of the re－sampling signal

对x（t）设置一个采样率f′s为f1的整数倍，取一个周期信号做FFT变换得到结果，以及让x（t）经过算法做频率估计和采样速率转换后也取一个周期信号做FFT变换得到的结果于表2中所示。因为重采样过程中存在泄漏，因此频率估计到第5次才收敛于50.20066477092068。

表2 理论整数倍采样率与同步采样后的信号FFT变换结果Tab.2 FFT comparison of the signals with the theoretical integer sampling rate and the synchronous sampling rate

从表2中可以看到，如果采样频率严格为电力信号的整数倍的时候，在频域中除了基波和谐波项，其他的幅度都为零。

采用本文的频率估计和重采样算法后信号的频谱仍然存在少量的泄漏，这是因为在采样速率转换过程中不可避免地会有泄漏，使得采样频率没有严格地等于信号频率的整数倍引起。但是基波幅度和谐波幅度最大比值仅仅为7.761520465715398×10－5，说明频谱泄漏很小，此精度完全能满足后面对电力信号的谐波分析，完全能准确地找出谐波项。

4 结论

（1）本文提出了运用精细频率估计来估算和跟踪电力信号的频率，这种方法可以无穷逼近电力信号的真实频率，一次估计据可以达到0.001的精度，仿真发现当信号中无论是否存在谐波和噪声，经过3次估算就可以完全得到电力信号的频率，当重采样中存在泄漏时经过5次循环估计的频率就收敛。

（2）利用MATLAB7.0自带的工具Fdatool设计好低通滤波器的脉冲响应的系数保存成表格形式，再运用Farrow滤波器的算法查表实现信号的采样率转换。通过查表而不用脉冲响应的表达式，在保证采样频率转换效果良好的同时可以很大程度上减少了信号采样速率转换的运算量。

（3）本文算法经仿真和实践检验证明可行有效，该算法可以应用于各种电力配电系统装置的前端处理中的电力信号的同步采样。

［1］ 曹国剑（Cao Guojian）.频谱校正理论在电参量检测中的应用研究（Research on the Application of the Spectrum Correction Theory in the Detection of the Electric Parameters）［D］.长沙：湖南大学电气与信息工程学院（Changsha：Institute of Electrical and Information Engineering，Hunan University），2004.

［2］ 陈春玲（Chen Chunling）.电能质量扰动分析与监测研究（Study on Disturbance Analysis Method and Monitoring of Power Quality）［D］.沈阳：沈阳农业大学信息与电气工程学院（Shenyang：Institute of Information and Electrical Engineering，Shenyang Agricultural University），2009.

［3］ 武艳华（Wu Yanhua）.频谱分析理论在谐波及间谐波检测中的应用研究（Research on the Application of the Spectrum Analysis Theory in the Detection of Harmonics and Inter－harmonics）［D］.长沙：湖南大学电气与信息工程学院（Changsha：Institute of Electrical and Information Engineering，Hunan University），2006.

［4］ 黄纯（Huang Chun）.谐波分析的加窗插值改进算法（Improved window and interpolation algorithm for analysis of power system harmonics）［J］.中国电机工程学报（Proceedings of the CSEE），2005，25（15）：26－32.

［5］ Andria G，Savino M，Trotta A.Windows and interpolation algorithms to improve electrical measurement accuracy［J］.IEEE Trans on Instrumentation and Measurement，1989，38（4）：856－863.

［6］ Heydt G T，Fjeld P S，Liu C C，et al.Applications of windowed FFT to electric power quality assessment［J］.IEEE Trans on Power Delivery，1999，l4（4）：1411－1416.

［7］ 潘文，钱俞寿，周鹗（Pan Wen，Qian Yushou，Zhou E）.基于加窗插值FFT的电力谐波测量理论（Ⅰ）──窗函数研究（Power harmonic measurement based on windows and ieterpolated FFT（I）：Study of win－dows）［J］.电工技术学报（Transactions of China Electrotechnical Society），1994，2（1）：50－54.

［8］ 张伏生，耿中行，葛耀中（Zhang Fusheng，Geng Zhongxing，Ge Yaozhong）.电力系统谐波分析的高精度FFT算法（FFT algorithm with high accuracy for harmonic analysis in power system）［J］.中国电机工程学报（Proceedings of the CSEE），1999，l9（3）：63－66.

［9］ 杨洪耕，惠锦，侯鹏（Yang Honggeng，Hui Jin，Hou Peng）.电力系统谐波和间谐波检测方法综述（Detection methods of harmonics and inter－harmonics in power system）［J］.电力系统及其自动化学报（Proceedings of the CSU－EPSA），2010，22（2）：65－69.

［10］赵文春，马伟明，胡安（Zhao Wenchun，Ma Weiming，Hu An）.电机测试中谐波分析的高精度FFT算法（FFT algorithm with high accuracy for harmonic analysis in the electric machine）［J］.中国电机工程学报（Proceedings of the CSEE），2001，21（12）：83－87，92.

［11］柴旭峥，文习山，关根志，等（Chai Xuzheng，Wen Xishan，Guan Genzhi，et al）.一种高精度的电力系统谐波分析算法（An algorithm with high accuracy for analysis of power system harmonics）［J］.中国电机工程学报（Proceedings of the CSEE），2003，23（9）：67－70.

［12］戴先中，唐统一（Dai Xianzhong，Tang Tongyi）.周期信号谐波分析的一种新方法（A novel method for spectrum analysis of periodic signals）［J］.仪器仪表学报（Chinese Journal of Scientific Instrument），1989，10（3）：248－255.

［13］丁康，张晓飞（Ding Kang，Zhang Xiaofei）.频谱校正理论的发展（Advances in spectrum correction theory）［J］.振动工程学报（Journal of Vibration Engineering），2002，13（1）：14－22.

［14］刘渝（Liu Yu）.快速高精度正弦波频率估计综合算法（A fast and accurate single frequency estimator synthetic approach）［J］.电子学报（Acta Electronica Sinica），1999，27（6）：126－128.

［15］崔保延.GPS软件接收机基础［M］.陈军，译.北京：电子工业出版社，2001.

［16］张盛耀（Zhang Shengyao）.软件无线电中用于采样速率转换的Farrow结构滤波器设计（Design of Farrow structure filter based sample rate conversion in software radio）［J］.广西通信技术（Guangxi Communication Technology），2007，3（1）：20－22，42.

［17］陈彩莲，于宏毅，沈彩耀，等（Chen Cailian，Yu Hongyi，Shen Caiyao，et al）.采样率转换中Farrow滤波器实现结构研究（Efficient implementation for sample rate conversion using Farrow structure）［J］.信息工程大学学报（Journal of Information Engineering University），2009，10（3）：329－332.

［18］Harris F.Performance and design of Farrow filter used for arbitrary resampling［C］∥International Conference on Digital Signal Processing.Santorini，Greece：1997.