实验7_dft_fft算法.docVIP

实验7 DFT/FFT算法 【实验目的】 1.加深对离散信号的DFT的理解 2.在MATLAB中实现FFT算法 3.熟练掌握DFS和DFT的表示和实现方法，以及两种变换的应用和物理意义 【实验原理】 周期序列傅立叶级数 DFS 周期为N的序列，傅立叶级数展开记为： DFS:， IDFS: 上式表明，周期为N 周期序列也是一个以N为周期的周期序列；频谱只有N个是独立的谐波成分。 MATLAB并没有关于离散傅里叶级数计算的公式，为了计算周期序列的离 散傅里叶级数方便，编写相应的函数d f s。 %实现周期序列的离散傅里叶级数 %dfs.m %xn是输入的时域序列，N为其长度 %Xk是输出的频域序列，因长度也为N，所以不用输出 function[Xk]=dfs(xn,N) n=[0:1:N-1]; k=[0:1:N-1]; WN=exp(-j*2*pi/N); nk=n*k; WNnk=WN.^nk; Xk=xn*WNnk; %逆离散傅里叶级数系数 %idfs.m function[xn]=idfs(Xk,N) n=[0:1:N-1]; k=[0:1:N-1]; WN=exp(-j*2*pi/N); nk=n*k; WNnk=WN.^(-nk); xn=(Xk*WNnk)/N; 离散傅立叶变换DFT 若设离散时间非周期信号为序列x(n)，则序列的傅里叶变换DTFT为： N点序列，傅立叶变换记为 IDFT: %离散傅里叶变换 %dft.m function[Xk]=dft(xn,N) n=[0:1:N-1]; k=[0:1:N-1]; WN=exp(-j*2*pi/N); nk=n’*k; WNnk=WN.^nk; Xk=xn*WNnk; %逆离散傅里叶变换 %idft.m function[xn]=idft(Xk,N) n=[0:1:N-1]; k=[0:1:N-1]; WN=exp(-j*2*pi/N); nk=n*k; WNnk=WN.^(-nk); xn=(Xk*WNnk)/N; 3．快速傅立叶算法 FFT 考虑DFT中WN的周期性和对称性，可以化简算法，得到FFT，分为两种： DIT（时间抽取的）－FFT FIT（频率抽取的）－FFT FFT函数 y=fft(x) 当x的个数是2的整数倍时，采用基2，其他采用混合基 y=fft(x,n) (n点计算，不够时补零) 举例：绘出有限长序列的频谱（有限长为5，样值为1，长度为20）。 L=5;N=20; n=1:N; xn=[ones(1,L),zeros(1,N-L)]; Xk=dft(xn,N); magXk=abs([Xk(N/2+1:N) Xk(1:N/2+1)]); k=[-N/2:N/2]; figure(1) subplot(2,1,1);stem(n,xn); xlabel(n);ylabel(x(n));title(DFT: L=5,N=20) subplot(2,1,2);stem(k,magXk); xlabel(k);ylabel(X(k));axis([-N/2,N/2,-0.5,5.5]) 举例2：x(n)=cos(0.48πn)（0=N=10），求X(K)； n1=[0:1:9]; x=cos(0.48*pi*n1); y1=x(1:1:10); subplot(2,1,1);stem(n1,y1); title(x(n) 0=n=9); xlabel(n),axis([0,10,-2.5,2.5]) Y1=fft(y1);magY1=abs(Y1(1:1:6)); k1=0:1:5;w1=2*pi/10*k1; subplot(2,1,2);stem(w1/pi,magY1); title(Samples of DTFT Magnitude); xlabel(frequency in pi unite);axis([0,1,0,10]) 【实验内容】 1.绘出有限长序列的频谱（有限长为L，样值为1，长度为N）。比较图形，说明了什么问题？ L＝5，N＝20； 程序：L=5;N=20; n=1:N; xn=[ones(1,L),zeros(1,N-L)]; Xk=dft(xn,N); magXk=abs([Xk(N/2+1:N) Xk(1:N/2+1)]); k=[-N/2:N/2]; figure(1) subplot(2,1,1);stem(n,xn); xlabel(n);ylabel(x(n));title(DFT: L=5,N=20) subplot(2,1,2);stem(k,magXk); xla