用DFT对信号进行谱分析实验报告

用DFT(FFT)对信号进行谱分析2015年4月1日课程名称：数字信号处理实验名称：DFT对信号进行分析学号：姓名：______指导老师评定：签名：__________________一、实验目的1、在理论学习的基础上，通过本次实验加深对DFT的理解。2、熟悉应用FFT对典型信号进行频谱分析的方法。3、了解应用FFT进行信号频谱分析过程中可能出现的各种误差，以便在实际中正确应用FFT。二、实验原理在运用DFT进行频谱分析的时候可能会产生三种误差，现分析如下：（一）截断效应实际中的信号序列往往很长，甚至是无限长序列。为了方便，我们往往只取实际序列的一部分来近似它们。这种截短等价于给原信号序列乘以一个矩形窗函数。根据卷积定理，最终信号的频谱等于原信号的谱和矩形窗的谱的卷积，从而造成谱线加宽或称为频谱泄漏。矩形窗时间取得越长，矩形窗的频谱变窄，由截断引起的效应会减小。例如50Hz正弦波xa(t)=sin(2π·50t)，它的幅度曲线是线状谱，如图3.1(a)所示。如果将它截取0.09s的一段，相当于将它乘一长度为0.09s矩形窗函数，即xa(t)RTp(t)，Tp=0.09s,该信号的谱等于原信号的谱和矩形窗的谱的卷积，如图1（b）所示。矩形窗长度扩大Tp=0.18s,后，频谱泄漏会变小，如图1（c）。10.50－250－200－150－100－50050100150200250幅度f/Hz(a)10.50幅度－250－200－150－100－50050100150200250f/Hz(b)10.50幅度－250－200－150－100－50050100150200250f/Hz(c)图3.1用DFT对正弦波进行谱分析(a)50Hz正弦波的幅频曲线；(b)50Hz正弦波加窗后的幅频曲线(Tp=0.09s);(c)50Hz正弦波加窗后的幅频曲线(Tp=0.18s)同时，由于频谱泄漏，还会造成靠得很近的两个谱峰混淆为一个谱峰，或是强的谱线的旁瓣掩盖弱的谱线，称为谱间干扰，导致频谱分辨率降低。矩形窗时间取得越长，矩形窗的频谱变窄，由截断引起的效应会减小。泄漏和谱间干扰都会小。（二）频谱分析时种因素的综合考虑在实际用FFT对模拟信号进行谱分析时，首先要把模拟信号转换成数字信号，转换时要求知道模拟信号的最高截止频率，以便选择满足采样定理的采样频率，避免混叠效应。一般选择采样频率是模拟信号中最高频率的3~4倍。如果模拟信号不是严格的带限信号，会因为频谱混叠现象引起谱分析的误差，这种情况下可以预先将模拟信号进行预滤，或者尽量将采样频率取高一些。除选择采样频率外，还需要确定分辨率F，并由此进一步确定观测时间，避免栅栏效应。同时，取样时间加长，可以减小泄漏，从而使得相近频率的谱间干扰减小。最小的观测时间Tpmin和分辨率成倒数关系，一般用教材(6.7.8)式确定。最小的采样点数用教材(6.7.12)式确定。要求选择的采样点数和观测时间大于它的最小值。fsmin=2fc用FFT作谱分析时，一般取FFT的点数服从2的整数幂，这一点在上面选择采样点数时可以考虑满足，有时可以通过在序列尾部加0完成。如果要进行谱分析的模拟信号是周期信号，最好选择观测时间是信号周期的整数倍。如果不知道信号的周期，要尽量选择观测时间长一些，以减少截断效应的影响。举一个极端的例子，一个周期性正弦波，如果所取观察时间太短，例如取小于一个周期，它的波形和正弦波相差太大，肯定误差很大，但如果取得长一些，即使不是周期的倍数，这种截断效应也会小一些。三、实验内容及步骤1.初步练习：用信号x1(t)=cos20πt观察频谱泄漏现象，栅栏效应及正弦波抽样时的规律。信号频率为f=10Hz，周期为0.1秒。令采样频率为fs=100Hz，采样周期Ts=0.01s，满足采样定理要求及一般正弦信号抽样习惯。令截取长度为N=100点，共取10个整周期，根据fk=kfs/N计算可得，模拟信号x1(t)对应的k值为10,在k=10处观察到谱线。周围未观察到频谱泄漏。这时候因为截断的关系，泄漏依旧存在，其实质在于泄漏处的不为零的频谱都处在k值的中间，被栅栏挡住未显示（FFT只计算并显示了整数k对应的频谱点）。整数k值刚好处在主谱线以及泄漏频谱的零minmin2sCffNFFmin1pTF点上。主谱线被显示，泄漏的非零点未显示。假如抽样点较少，对于正弦函数，抽样有没有位于整周期处，则由于栅栏效应，k值实际取到的位置距离信号的实际频谱的位置就会太远，不能代表信号的实际频率。因此信号中若含有单个的正弦成分，在不知道正弦成分的准确周期的情况下，抽样点数要增加，直到前后两次做出的频谱接近到一定精度。由本实验应该看出，单个频率的信号，由于截断效应的影响，向两边泄漏，但存在一个主峰。如果k值抽样对应的频率刚好取在主峰未知，就能避开栅栏效应的影响，较好地反映信号包含的实际频率。2．进一步练习观察：由包含两个频率的信号x2(t)=cos20πt＋2cos55πt。分别取N=100,512,1024。可以看出，只有取到较大的值，即时间取长之后，才能使两个信号均有效避免栅栏效应的过度影响，反映两个信号的真实相对幅度大小（2倍关系）。如果取样时间过短，那么k值抽样时，就不一定都会抽取在主峰或其附近，造成频谱分析的误差。3．较复杂的练习：用信号x3(t)=cos0.48πt＋cos0.52πt，其频率分别为0.24Hz,0.26Hz。采样周期选1秒，满足采样定理要求。然后分别取（1）N=10；（2）N=10然后补90个零。四．实验结果：参考程序functionexperiment31%观察单个正弦波泄漏、栅栏、正弦波的整周期取样closeall%关闭所有图形窗口Ts=0.01;%取样周期fs=1/Ts;N=68;%取样点数68,取样长度为0.68秒，共6.8个信号周期n=0:N-1;t=n*Ts;%计算取样时间xn=cos(20*pi*t);%信号频率10Hz，周期0.1秒subplot(221);stem(t,xn);%绘制原信号波形axis([0,0.5,-10,10]);title(['T=0.01,N=',num2str(N)]);ylabel('x1(n)');Xk=fft(xn);MAGXk=abs(Xk);%进行FFT变换并求幅频特性fn=fs/N*n%对应各个n值的实际模拟频率值subplot(222);stem(n,MAGXk);%绘制幅频特性axis([0,20,0,50]);ylabel('X(k)');title('x(n)的频谱');N=100;%取样点数100,取样长度为1秒，共10个信号周期n=0:N-1;t=n*Ts;xn=cos(20*pi*t);subplot(223);stem(t,xn);axis([0,0.5,-10,10]);title(['T=0.01,N=',num2str(N)]);ylabel('x2(n)');Xk=fft(xn);MAGXk=abs(Xk);fn=fs/N*n%对应各个n值的实际模拟频率值subplot(224);stem(n,MAGXk);ylabel('X(k)');title('x(n)的频谱');functionexperiment32%进一步观察泄漏及栅栏效应closeall%关闭所有图形窗口Ts=0.01;%取样周期fs=1/Ts;N=100;%取样点数100,n=0:N-1;t=n*Ts;%计算取样时间xn=cos(20*pi*t)+2*cos(55*pi*t);%信号频率10Hz，27.5Hz,2倍幅值关系subplot(421);stem(t,xn);%绘制原信号波形ylim([-3,3]);title(['T=0.005,N=',num2str(N)]);ylabel('x1(n)');Xk=fft(xn);MAGXk=abs(Xk);%进行FFT变换并求幅频特性fn=fs/N*n%对应各个n值的实际模拟频率值subplot(422);stem(n,MAGXk);%绘制幅频特性axis([0,50,0,100]);ylabel('X1(k)');title('x1(n)的频谱');N=512;%取样点数扩大,取样长度增加n=0:N-1;t=n*Ts;%计算取样时间xn=cos(20*pi*t)+2*cos(55*pi*t);%信号频率10Hz，27.5Hz,2倍幅值关系subplot(423);stem(t,xn);%绘制原信号波形ylim([-3,3]);title(['T=0.005,N=',num2str(N)]);ylabel('x1(n)');Xk=fft(xn);MAGXk=abs(Xk);%进行FFT变换并求幅频特性fn=fs/N*n%对应各个n值的实际模拟频率值subplot(424);stem(n,MAGXk);%绘制幅频特性axis([0,256,0,512]);ylabel('X(k)');title('x(n)的频谱');N=1024;%取样点数扩大,取样长度增加n=0:N-1;t=n*Ts;%计算取样时间xn=cos(20*pi*t)+2*cos(55*pi*t);subplot(425);stem(t,xn);%绘制原信号波形ylim([-3,3]);title(['T=0.005,N=',num2str(N)]);ylabel('x1(n)');Xk=fft(xn);MAGXk=abs(Xk);%进行FFT变换并求幅频特性fn=fs/N*n%对应各个n值的实际模拟频率值subplot(426);stem(n,MAGXk);%绘制幅频特性axis([0,512,0,1024]);ylabel('X(k)');title('x(n)的频谱');N=2048;%取样点数扩大,取样长度增加n=0:N-1;t=n*Ts;%计算取样时间xn=cos(20*pi*t)+2*cos(55*pi*t);subplot(427);stem(t,xn);%绘制原信号波形ylim([-3,3]);title(['T=0.005,N=',num2str(N)]);ylabel('x1(n)');Xk=fft(xn);MAGXk=abs(Xk);%进行FFT变换并求幅频特性fn=fs/N*n%对应各个n值的实际模拟频率值subplot(428);stem(n,MAGXk);%绘制幅频特性axis([0,1024,0,2048]);ylabel('X(k)');title('x(n)的频谱');functionexperiment33%观察谱间干扰；观察截断、栅栏效应的综合影响。%本实验事先已知频率。%但是：若事先不知频率，则N值需要取大一些，当改变N值时，若频谱变化不大，则可认为误差小到一定范围，截断、栅栏等影响减小到一定范围。closeall%关闭所有图形窗口Ts=1;%取样周期1秒，频率1Hzfs=1/Ts;N=10;%取样点数,取样长度为NTs秒，频率分辨率为1/NTsHzn=0:N-1;t=n*Ts;%计算取样时间xn=cos(0.48*pi*t)+cos(0.52*pi*t);%信号频率0.24Hz,0.26Hzsubplot(221);stem(t,xn);%绘制原信号波形title(['T=1,N=',num2str(N)]);ylabel('x1(n)');Xk=fft(xn);MAGXk=abs(Xk);%对x(n)进行FFT变换并求幅频特性fn=fs/N*n%对应各个n值的实际模拟频率值subplot(222);stem(n,MAGXk);%绘制x1(n)幅频特性axis([0,N/2,0,max(MAGXk)]);fk=n/N/Ts%观察相应的k值对应的实际模拟频率值ylabel('X(k)');title('x(n)的频谱');M=100;%该部分观察补零对于频谱的影响，频谱会更密，但无助消除频谱间的混淆。m=0:M-1;xn2=[xn,zeros(1,90)];title(['T=1,N=',nu