试验四利用FFT对信号进行频谱分析1112

1、实验四 利用DFT对信号进行频谱分析、目的要求进一步加深对线性卷积的理解和分析能力；通过编程，上机调试程序，进一步增强使用计算机解决问题的能力;掌握线性卷积与循环卷积软件实现的方法，并验证二者之间的关系;掌握利用快速卷积的计算机实现方法；熟练掌握用MATLAB实现FFT及信号的频谱分析。二、实验环境PC机；MATLAB语言环境三、实验原理及实例分析FFT是DFT的快速算法，凡是可以利用离散傅里叶变换 DFT来进行计算的场合， 都可以用利用FFT算法及数字信号处理技术加以实验。FFT在数字通信、语言信号处理、数字图像处理、匹配滤波器以及功率谱计算、系统分析与仿真等诸多领域都得到了广泛的应用。(一

2、)利用MATLAB佥证两个序列的线性卷积和圆周卷积的关系线性卷积线性时不变系统(Linear Time-Invariant System, or L. T. I系统)输入、输出间的关系为:y(n),则系统输出为:当系统输入序列为x(n),系统的单位脉冲响应为h(n),输出序列为00y(n) = x(m)h(n -m) = x(n)* h(n)m 二.：二-bey(n) =、h(m)x(n-m) =h(n)*x(n)m=-：7上式称为离散卷积或线性卷积。图4.1示出线性时不变系统的输入、输出关系。6(n)- L. T. I 一 h(n)-bey(n);1 x(m)h(n - m)m 二：二图4.

3、1线性时不变系统的输入、输出关系循环卷积设两个有限长序列 X1(n)和x2 (n),均为N点长/、- D F T 八Xi (n) *Xi(k)D F TX2(n)X2(k)如果 X3(k) =Xi(k) X2(k)X3(n)= Xi(m)2(n-m) Rz(n) ,m=0N 1x、Xi(m) X2 (n - m)N 1m _0=Xi(n) N X2(n)0 _ n _ N -1上式称为循环卷积或圆周卷积 TOC o 1-5 h z 注：Xi(n)为Xi(n)序列的周期化序列；Xi(n)RN (n)为Xi(n)的主值序列。上机编程计算时，X3(n)可表示如下：nN i HYPERLINK l b

4、ookmark45 o Current Document X3(n)= x(m)X2(n-m)-二：x1 (m)X2(Nn-m)m =0m 印-d两个有限长序列的线性卷积序列Xi(n)为L点长，序列X2(n)为p点长，X3(n)为这两个序列的线性卷积，则X3(n) 为-bex3 (n) =xi (m)x2 (n - m)m 二二二且线性卷积X3(n)的最大长L + P i ,也就是说当n w -1和n之L + P i时 X3(n) =0。循环卷积与线性卷积的关系序列Xi(n)为L点长，序列X2(n)为P点长，若序列Xi (n)和X2(n)进彳f N点的循环卷 积，其结果是否等于该两序列的线性卷

5、积，完全取决于循环卷积的长度：当N之L + P -1时循环(圆周)卷积等于线性卷积，即Xi(n) N X2 (n) = Xi(n)*X2(n)当N L + P -1时，循环卷积等于两个序列的线性卷积加上相当于下式的时间混叠，即X3N(n)匕 X3(n rN)r=o00 n Nerror(N must not be less than length of x1);endif length(x2)Nerror(N must not be less than length of x2);end%以上语句判断两个序列的长度是否小于Nx1=x1,zeros(1,N-length(x1);%填充序列 x1

6、(n)使其长度为 N1+N2-1(序列%h(n)的长度为N1 ,序列x(n)的长度为N2)x2=x2,zeros(1,N-length(x2);%填充序列x2(n)使其长度为 N1+N2-1n=0:1:N-1;x2=x2(mod(-n,N)+1);%生成序列x2(-n)NH=zeros(N,N);for n=1:1:NH(n,:)=cirshiftd(x2,n-1,N);%该矩B的 k 行为 x2(k-1-n)Nendyc=x1*H;%计算循环卷function y=cirshiftd(x,m,N)%directly realize circular shift for sequence x%

7、y=cirshiftd(x,m,N);%x:input sequence whose length is less than N%m:how much to shift%N:circular length%y:output shifted sequenceif length(x)Nerror(length of x must be less than N);endx=x,zeros(1,N-length(x);n=0:1:N-1;y=x(mod(n-m,N)+1);%研究两者之间的关系clear all;n=0:1:11;m=0:1:5;N1=length(n);%生成x(n)%生成h(n)%

8、直接用函数conv计算线性卷积%用函数circonv计算N1点循环卷积%画图N2=length(m);xn=0.8.An;hn=ones(1,N2);yln=conv(xn,hn);ycn=circonv(xn,hn,N1);ny1=0:1:length(yln)-1;ny2=0:1:length(ycn)-1;subplot(2,1,1);stem(ny1,yln);subplot(2,1,2);stem(ny2,ycn);axis(0,16,0,4);运行结果:积卷性线积卷环循图4.1 线性卷积和循环卷积的比较【任务一】请同学修改程序，使得循环卷积可以代替线性卷积(二)利用快速卷积法计算两

9、个序列的卷积快速卷积法的计算框图如下所示：%产生序列x(n)%产生序列h(n)%序列x(n)和h(n)卷积后得到序列 y(n)的长度%对序列x(n)作L点DFT%对序列h(n)作L点DFT【实例4-2J已知序列x(n) =sin(0.4n)R15(n),h(n) =0.9nR2O(n),试利用快速卷积法计算这两个序列的卷积 y(n) =x(n)* h(n)。MAT LAB程序清单如下：Nx=15;Nh=20;n1=1:Nx-1;n2=0:Nh-1;xn=sin(0.4*n1).*(n1=0 &n1=0 &n2Nh);L=pow2(nextpow2(Nx+Nh-1);Xk=fft(xn,L);H

10、k=fft(hn,L);yn=ifft(Xk.*Hk) subplot(3,1,1);stem(xn,.);title(x(n);subplot(3,1,2);stem(hn,.);title(h(n);subplot(3,1,3);stem(yn,.);title(y(n);【任务二】修改程序，比较用直接调用conv函数计算结果与FFT函数计算的结果是否一致。(三)利用FFT对信号进行频谱分析【实例4-3】已知连续信号，x(t) =cos(8nt) + cos(16nt) + cos(20nt)选择采样频率fs =64Hz,变换区间在N=16, 32, 64三种情况下进行谱分析，分别打印其幅

11、频特性，并进行分析和讨论。%程序清单Fs=64;T=1/Fs;N=16;n=0:N-1;%FFT 的变换区间 N=16x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);%对 x6(t)16 点采样X6k16=fft(x6nT);%计算 x6nT 的 16 点 DFTX6k16=fftshift(X6k16);%将零频率移到频谱中心Tp=N*T;F=1/Tp;%频率分辨率Fk=-N/2:N/2-1;fk=k*F;%产生16点DFT对应的采样点频率（以零频率为中心）subplot(3,1,1);stem(fk,abs(X6k16),.); % 绘制

12、16 点 DFT 的幅频特性图box on; %坐标系右边和上边有边框,off 了就没有了title(16 点 |DFTx(nT)|);xlabel(f(Hz);ylabel(幅度);axis(-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k16)N=32;n=0:N-1;%FFT的变换区间 N=16x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);%对x6（t）32点采样X6k32=fft(x6nT);%计算x6nT的32点DFTX6k32=fftshift(X6k32);%将零频率移到频谱中心Tp=N*T;F=1/Tp;%

13、频率分辨率Fk=-N/2:N/2-1;fk=k*F;%产生32点DFT对应的采样点频率（以零频率为中心）subplot(3,1,2);stem(fk,abs(X6k32),.);box on % 绘制 32 点 DFT 的幅频特性图title(32 点 |DFTx(nT)|);xlabel(f(Hz);ylabel(幅度);axis(-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k32);N=64;n=0:N-1;%FFT的变换区间 N=16%对x6(t)64点采x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);X6k64=f

14、ft(x6nT);%计算x6nT的64点DFTX6k64=fftshift(X6k64);%将零频率移到频谱中心Tp=N*T;F=1/Tp;%频率分辨率Fk=-N/2:N/2-1;fk=k*F;%产生64点DFT对应的采样点频率（以零频率为中心）subplot(3,1,3);stem(fk,abs(X6k64),.); box on%绘制 64 点 DFT 的幅频特性图title( 64 点 |DFTx(nT)|);xlabel(f(Hz);ylabel( 幅度);axis(-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k64)【任务二】若分析修改fs =40,70,128,分析运行结果。【任务三】若分析修改N= 40,70,128,分析运行结果。四、实验内容1.有限长序列x（n）=2,1,0,1,3;h（n）=1,3,2,1），试利用FFT实现由DFT计算线性卷积，并与线性卷积直接计算（c