基于加海明窗的FFT高精度谐波检测算法分析及其在测控设备低频干扰信号检测中的应用.doc

基于加海明窗的FFT高精度谐波检测算法分析及其在测控设备低频干扰信号检测中的应用技术室 马东阳摘要：大量电子元器件构成的测控系统，各分系统之间不可避免的会产生低频谐波干扰，采用加Hamming窗的FFT算法，在有效抑制频谱泄漏的情况下，对低频谐波信号进行检测。分析了算法的原理，设计了算法的应用流程，通过Matlab仿真测试，最终在测控设备低频谐波干扰信号检测中取得了较好的应用效果。关键词：低频干扰；快速傅立叶变换；海明窗；频谱泄漏1 引言现代测控系统功能日益复杂，所用电子元器件不断增多，随着设备列装时间增长，由大量电子元器件构成的各个分系统不可避免的会产生低频谐波干扰信号，对测控系统测量精度产生影响。本文分析了基于加海明窗的高精度快速傅立叶变换算法的物理意义及其在应用中需要注意的问题，并结合试验任务对测控设备低频干扰信号进行了检测，取得了较好的应用效果。2 FFT算法分析2.1 FFT算法结果的具体物理意义FFT是离散傅立叶变换的快速算法，可以将一个信号由时域变换到频域。时域上的任何周期信号或有限长信号y(t)，都可以用无数个余弦函数之和来表示： （1）式中：fk 为频率；A(fk)为幅值；(fk)为相位。基于FFT 获取y(t)频谱的基本过程如下：（1）对时域波形y(t)等间隔采样，使之变成适合计算机处理的N 点离散信号y(n)；（2）根据下式对y(n)进行离散快速傅立叶变换，即FFT，得：，其中 （2）（3）记Y(k)的实部为，虚部为，频率对应的幅度和相位为 （3）由上可知，假设采样得到N个采样点，经过FFT变换之后，就可以得到N个点的FFT结果。为了方便进行FFT运算，通常N取2的整数次方。假设采样频率为Fs，信号频率F，采样点数为N；那么FFT之后结果就是一个N点的复数。若原始信号的峰值为A，那么FFT的结果的每个点（除了第一个点直流分量之外）的模值就是A的N/2倍；而第一个点是直流分量，它的模值是直流分量的N倍。每个点的相位，就是在该频率下信号的相位；第一个点表示直流分量（即0Hz），而最后一个点N的再下一个点（实际上这个点是不存在的）则表示采样频率Fs，这中间被N-1个点平均分成N等份，每个点的频率依次增加。2.2 FFT算法在实际应用中需要注意的几个问题2.2.1栅栏效应由2.1可知，在采样点数为N，采样频率为Fs时，第n点所表示的频率为： （4）由公式（4）可以看出，Fn所能分辨到频率为Fs/N，如果采样频率Fs为1024Hz，采样点数为1024点，则可以分辨到1Hz。1024Hz的采样率采样1024点，刚好是1秒，也就是说采样1秒时间的信号并做FFT，则结果可以分析到1Hz，如果采样2秒时间的信号并做FFT，则结果可以分析到0.5Hz。这就相当于透过栅栏观赏风景,只能看到频谱的一部分,而其它频率点看不见,因此很可能使一部分有用的频率成分被漏掉,此种现象被称为栅栏效应。当时域采样满足采样定理时，栅栏效应不会有什么影响；而频域采样的栅栏效应则影响很大，“挡住”或丢失的频率成分有可能是重要的或具有特征的成分，使信号处理失去意义。减小栅栏效应可用提高采样间隔也就是频率分辨力的方法来解决。间隔小，频率分辨力高，被“挡住”或丢失的频率成分就会越少。2.2.2频谱泄漏对于频率为Fs的正弦序列，它的频谱应该只是在Fs处有离散谱。但是，在利用FFT求它的频谱时，对时域做了截短，结果使信号的频谱不只是在Fs处有离散谱，而是在以Fs为中心的频带范围内都有谱线出现，我们可以理解为是从Fs频率上“泄漏”出去的，这种现象称为频谱泄漏。同样在处理测控设备采集的离散数据x(n)时，通常需要将序列x(n)进行截短，取一部分数据进行分析，这就不可避免会出现频谱泄漏的情况。频谱泄漏引起的直接后果就是降低频谱分辨率，令主谱线旁边有很多旁瓣，这就会造成谱线间的干扰，严重时旁瓣的能量强到分不清是旁瓣还是信号本身。为了减少频谱能量泄漏，可采用不同的截取函数对信号进行截断。泄漏与窗函数频谱的两侧旁瓣有关，如果两侧瓣的高度趋于零，而使能量相对集中在主瓣，就可以较为接近于真实的频谱。为此，在时间域中可采用不同的窗函数来截断信号。典型的加窗序列有海明、汉宁、布莱克曼、高斯等窗序列。本文在处理实时测量数据时采用了加海明窗的方法。海明窗在时域的表达式如下： （5）图1 Matlab仿真的海明窗函数图形3 算法实现及在测控设备低频干扰信号检测中的应用3.1 FFT算法实现FFT算法的核心部分即快速傅立叶变换的蝶形运算，需要处理的是复数，故在参与运算前需要将测控设备采样的离散数据变换成复数。算法主要应用流程如下：（1）在时域上给采样得到的离散数据加海明窗，将加窗的结果序列做复数变换；（2）将复数序列做傅立叶变换；（3）变换后的结果根据第2.1节的分析，解算出相应频率对应的幅度值。算法总体处理流程如图2所示：图2 程序流程图算法核心部分是对复数进行的递归多级蝶形运算，在FFT函数中实现，主要代码如下： if (input.Length = 1) return new Complex input0 ; int length = input.Length; int half = length / 2; Complex output = new Complexlength; double fac = -2.0 * Math.PI / length; if (invert) fac = -fac; Complex evens = new Complexhalf; for (int i = 0; i half; i+) evensi = input2 * i; Complex evenResult = FFT(evens, invert); Complex odds = new Complexhalf; for (int i = 0; i half; i+) oddsi = input2 * i + 1; Complex oddResult = FFT(odds, invert); for (int k = 0; k half; k+) double fack = fac * k; Complex oddPart = oddResultk * new Complex(Math.Cos(fack), Math.Sin(fack); outputk = evenResultk + oddPart; outputk + half = evenResultk - oddPart; return output;3.2 算法处理结果与Matlab仿真结果对比将本文所设计算法的处理结果读入Matlab绘图，同时用Matlab自带的FFT算法进行仿真验证，结果如下。图3 本文所用算法和Matlab仿真结果4 结论本文在分析研究FFT算法的原理和应用重点基础上，对测控设备采集到的离散数据序列做了时域加Hamming窗的处理，进而采用FFT算法对其进行时-频变换；在频域上对信号进行频谱分析，检测测控系统内部低频谐波干扰，并用Matlab对原始数据进行了仿真分析。对比发现本文采用的数据处理流程满足测控系统低频谐波信号检测的要求，做成模块后，易于被各设备移植使用，方便编程扩展。参考文献1 胡广书数字信号处理M北京：清华大学出版社2 谢明，丁康频谱分析的校正方法J振动工程学报3 齐国清几种基于FFT的频率估计方法精度分析J振动工程学报4 齐国清，贾欣乐插值FF