13第十三讲-FFT的应用

第3章快速傅里叶变换第十三讲FFT的应用3.7利用FFT分析时域连续信号频谱3.8FFT的其他应用第3章快速傅里叶变换3.7利用FFT分析时域连续信号频谱3.7.1基本步骤图3-21时域连续信号离散傅里叶分析的处理步骤LPFsc(t)Ha(j)xc(t)A/Dx(n)w(n)v(n)DFTV(k)在实际应用中，前置低通滤波器的阻带不可能是无限衰减的，故由Xc(jΩ)周期延拓得到的X(ejω)有非零重叠，即出现频谱混叠现象。由于进行FFT的需要，必须对序列x(n)进行加窗处理，即v(n)=x(n)w(n)，加窗对频域的影响，是导致频谱泄漏。第3章快速傅里叶变换因为DFT对应的数字域频率间隔为Δω=2π/N，且模拟频率Ω和数字频率ω间的关系为ω=ΩΤ，其中Ω=2πf。所以，离散的频率函数第k点对应的模拟频率为:NTkTk2NTkfk数字域频率间隔Δω=2π/N对应的模拟域谱线间距为21sfFTNTNTN谱线间距，又称频谱分辨率（单位：Hz）。所谓频谱分辨率是指可分辨两频率的最小间距。它的意思是，如设某频谱分析的F=5Hz，那么信号中频率相差小于5Hz的两个频率分量在此频谱图中就分辨不出来。第3章快速傅里叶变换图3-22022F44F066F88F1010F1212F1414F1616F18fs＝1/TkfF＝1/tp010203040|V(k)|022T44T066T88T1010T1212T1414T1616T18tp＝1/FntT＝1/fs0246v(n)(b)(a)T为采样时间间隔（单位：s）；fs为采样频率（单位：Hz）；tp为截取连续时间信号的样本长度（又称记录长度，单位：s）；F为谱线间距，又称频谱分辨率（单位：Hz）。psptNTNfFNTt11第3章快速傅里叶变换在实际应用中,要根据信号最高频率fh和频谱分辨率F的要求，来确定T、tp和N的大小。（1）首先，由采样定理，为保证采样信号不失真，fs≥2fh(fh为信号频率的最高频率分量，也就是前置低通滤波器阻带的截止频率)，即应使采样周期T满足hfT21（2）由频谱分辨率F和T确定N,FTFfNs1为了使用FFT运算，这里选择N为2的幂次即N=2m，N大，分辨率好，但会增加样本记录时间tp。（3）最后由N,T确定最小记录长度，tp=NT。第3章快速傅里叶变换例3-3有一频谱分析用的FFT处理器，其采样点数必须是2的整数幂，假定没有采用任何特殊的数据处理措施，已给条件为:①频率分辨率≤10Hz;②信号最高频率≤4kHz。试确定以下参量：(1)最小记录长度tp(2)最大采样间隔T（即最小采样频率）(3)在一个记录中的最少点数N。第3章快速傅里叶变换解(1)由分辨率的要求确定最小长度tpsF1.01011所以记录长度为sFtp1.01(2)从信号的最高频率确定最大可能的采样间隔T（即最小采样频率fs=1/T）。按采样定理hsff2即sfTh3310125.01042121(3)最小记录点数N应满足80010104223FfNh取80010242210mN第3章快速傅里叶变换如果我们事先不知道信号的最高频率，可以根据信号的时域波形图来估计它。例如，某信号的波形如图3-23所示。先找出相邻的波峰与波谷之间的距离，如图中t1，t2，t3，t4。然后，选出其中最小的一个如t4。这里,t4可能就是由信号的最高频率分量形成的。峰与谷之间的距离就是周期的一半。因此，最高频率为)(214Hztfh知道fh后就能确定采样频率hsff2tt3t4ot1t2x(t)图3-23估算信号最高频率fh第3章快速傅里叶变换3.7.2可能出现的误差利用FFT对连续信号进行傅里叶分析时可能造成的误差如下。1.频谱混叠失真为了不产生混叠，必须满足hsff2也就是采样间隔T满足hsffT211对于FFT来说，频率函数也要采样，变成离散的序列，其采样间隔为F(即频率分辨率)。Ftp1psptNTNfFNTt11信号的最高频率分量fh与频率分辨率F存在矛盾关系，要想fh增加，则时域采样间隔T就一定减小，而fs就增加，此时若是固定N，必然要增加F，即分辨率下降。第3章快速傅里叶变换反之，要提高分辨率（减小F），就要增加tp，当N给定时，必然导致T的增加（fs减小）。要不产生混叠失真，则必然会减小高频容量（信号的最高频率分量）fh。要想兼顾高频容量fh与频率分辨率F，即一个性能提高而另一个性能不变（或也得以提高）的惟一办法就是增加记录长度的点数N，即要满足FfFfNhs2这个公式是未采用任何特殊数据处理（例如加窗处理）的情况下，为实现基本FFT算法所必须满足的最低条件。如果加窗处理，相当于时域相乘，则频域周期卷积，必然加宽频谱分量，频率分辨率就可能变差，为了保证频率分辨率不变，则须增加数据长度tp。第3章快速傅里叶变换2.栅栏效应利用FFT计算频谱，只给出离散点ωk=2πk/N或Ωk=2πk/(NT)上的频谱采样值，而不可能得到连续频谱函数，这就像通过一个“栅栏”观看信号频谱，只能在离散点上看到信号频谱，称之为“栅栏效应”,如果在两个离散的谱线之间有一个特别大的频谱分量，就无法检测出来。第3章快速傅里叶变换减小栅栏效应的一个方法就是要使频域采样更密，即增加频域采样点数N，在不改变时域数据的情况下，必然是在数据末端添加一些零值点，使一个周期内的点数增加，但并不改变原有的记录数据。频谱采样为2πk/N，N增加，必然使样点间距更近（单位圆上样点更多），谱线更密，谱线变密后原来看不到的谱分量就有可能看到了。必须指出，补零以改变计算FFT的周期时，所用窗函数的宽度不能改变。换句话说，必须按照数据记录的原来的实际长度选择窗函数，而不能按照补了零值点后的长度来选择窗函数。补零不能提高频率分辨率，这是因为数据的实际长度仍为补零前的数据长度。第3章快速傅里叶变换3.频谱泄漏与谱间干扰对信号进行FFT计算，首先必须使其变成有限时宽的信号，这就相当于信号在时域乘一个窗函数如矩形窗，窗内数据并不改变。时域相乘即v(n)=x(n)·w(n)，加窗对频域的影响，可用式卷积公式表示deWeXeVjjj)()(21)()(卷积的结果，造成所得到的频谱V(ejω)与原来的频谱X(ejω)不相同，有失真。这种失真最主要的是造成频谱的“扩散”（拖尾、变宽），这就是所谓的“频谱泄漏”。第3章快速傅里叶变换由上可知，泄漏是由于我们截取有限长信号所造成的。对具有单一谱线的正弦波来说，它必须是无限长的。也就是说，如果我们输入信号是无限长的，那么FFT就能计算出完全正确的单一线频谱。可是我们不可能这么做，而只能取有限长记录样本。如果在该有限长记录样本中，正弦信号又不是整数个周期时，就会产生泄漏。例如，一个周期为N=16的余弦信号x(n)=cos(2πn／１６）截取一个周期长度的信号即x1(n)=cos(2πn/16)R16(n),其16点FFT的谱图见3-24上图，若截取的长度为13，则其16点FFT的谱图见3-24下图。由此可见，频谱不再是单一的谱线，它的能量散布到整个频谱的各处。这种能量散布到其他谱线位置的现象即为“频谱泄漏”。第3章快速傅里叶变换图3-24余弦信号频谱泄漏示例图00246851015|X2(k)|00246851015kk|X1(k)|第3章快速傅里叶变换应该说明，泄漏也会造成混叠，因为泄漏将会导致频谱的扩展，从而使最高频率有可能超过折叠频率（fs/2），造成频率响应的混叠失真。泄漏造成的后果是降低频谱的分辨率。此外，由于在主谱线两边形成很多旁瓣，引起不同频率分量间的干扰(简称谱间干扰)，特别是强信号谱的旁瓣可能湮没弱信号的主谱线，或者把强信号谱的旁瓣误认为是另一信号的谱线，从而造成假信号，这样就会使谱分析产生较大偏差。在进行FFT运算时，时域截断是必然的，从而频谱泄漏和谱间干扰也是不可避免的。为尽量减小泄漏和谱间干扰的影响，需增加窗的时域宽度(频域主瓣变窄)，但这又导致运算量及存储量的增加；其次，数据不要突然截断，也就是不要加矩形窗，而是加各种缓变的窗（例如：三角形窗、升余弦窗、改进的升余弦窗等），使得窗谱的旁瓣能量更小，卷积后造成的泄漏减小。第3章快速傅里叶变换3.7.3应用实例一个长度为N的时域离散序列x(n)，其离散傅里叶变换X(k)（离散频谱）是由实部和虚部组成的复数，即)()()(1kjXkXkXR对实信号x(n)，其频谱是共轭偶对称的，故只要求出k在0，1，2，…,N/2上的X(k)即可。将X(k)写成极坐标形式)](arg[|)(|)(kXjekXkX式中，|X(k)|称为幅频谱，arg［X(k)］称为相频谱。第3章快速傅里叶变换实际中也常常用信号的功率谱表示，功率谱是幅频谱的平方，功率谱PSD定义为NkXkPSD2|)(|)(功率谱具有突出主频率的特性，在分析带有噪声干扰的信号时特别有用。由频谱图可以知道信号存在哪些频率分量，它们就是谱图中峰值对应的点。谱图中最低频率为k=0，对应实际频率为0（即直流）；最高频率为k=N/2，对应实际频率为f=fs/2；第3章快速傅里叶变换对处于0，1，2，…N/2上的任意点k，对应的实际频率为f=kF=kfs/N。由于所取单位不同，频率轴有几种定标方式。图3-25列出频率轴几种定标方式的对应关系。上图中f′为归一化频率，定义为sfff'f′无量纲，在归一化频率谱图中，最高频率为0.5。专用频谱分析仪器常用归一化频率表示。第3章快速傅里叶变换oooooN/4N/2k/(rad/s)f/Hzf′/rad0.5fs/2s/2图3-25第3章快速傅里叶变换例3-4已知信号x(t)=0.15sin(2πf1t)+sin(2πf2t)-0.1sin(2πf3t)，其中,f1=1Hz，f2=2Hz，f3=3Hz。x(t)包含三个正弦波，但从时域波形图3-26（a）来看，似乎是一个正弦信号，很难看到小信号的存在，因为它被大信号所掩盖。取fs=32Hz作频谱分析。解因fs=32Hz，故nnnnTxnx326sin1.0324sin322sin15.0)()(第3章快速傅里叶变换该信号为周期信号，其周期为N=32。现对x(n)作32点的离散傅里叶变换(DFT)，其幅度特性|X(k)|如图3-26(b)所示。图中仅给出了k=0,1,…,15的结果。k=16,17,…,32的结果可由|X(N-k)|=|X(k)|得出。因N=32,故频率分辨率F=fs/N=1Hz。由图3-26(b)可知,k=1,2,3所对应的频谱即为频率f1=1Hz，f2=2Hz,f3=3Hz的正弦波所对应的频谱。幅频图如图3-26（b）,该图中小信号成分清楚显示出来。可见小信号成分在时域中很难辨识而在频域中容易识别。第3章快速傅里叶变换图3-26已知信号的时域图和幅频图－2－10120102030(a)(b)0051015816nkx(n)|X(k)|第3章快速傅里叶变换3.8FFT的其他应用3.8.1线性卷积的FFT算法——快速卷积1.FFT的快速卷积法以FIR滤波器为例，因为它的输出等于有限长单位脉冲响应h(n)与有限长输入信号x(n)的离散线性卷积。设x(n)为L点，h(n)为M点，输出y(n)为10)()()(Mmmnxmhny第3章快速傅里叶变换y(n)也是有限长序列，其点数为L+M-1点。下面首先讨论直接计算线性卷积的运算量。由于每一个x(n)的输入值都必须和全部的h(n)值相乘一次，因而总共需要LM次乘法，这就是直接计算的乘法次数，以md表示为md=LM用FFT算法也就是用圆周卷积来代替这一线性卷积时，为了不产生混叠，其必要条件是使x(n)，h(n)都补零值点，补到至少N=M+L-1，即:0)()(nxnx0≤n≤L-1L≤n≤N-10)()(nhnh0≤n≤M-1M≤n≤N-1第