数字信号处理实验二DFT和FFT.docx

1、数字信号处理实验二DFT和FFT数字信号处理实验二DFT和FFT一实验内容1.认真复习周期序列DFS、有限长序列DFT的概念、旋转因子的定义、以及DFS和DFT的性质等有关内容；复习基2-FFT的基本算法，混合基-FFT的基本算法、Chirp-Z变换的算法等快速傅立叶变换的方法。2.掌握有限长序列的循环移位、循环卷积的方法，对序列共轭对称性的含义和相关内容加深理解和掌握，掌握利用DFT分析序列的频谱特性的基本方法。3.掌握FFT算法的基本原理和方法、Chirp-Z变换的基本原理和方法，掌握利用FFT分析序列的频谱特性的方法。4.熟悉利用MATLAB进行序列的DFT、FFT的分析方法。二实验内容

2、a.设周期序列某(n)=,0,1,2,3,0,1,2,3,0,1,2,3,.，求该序列的离散傅立叶级数某,(k)=DFS某(n)，并画出DFS的幅度特性。在matlab中新建函数df:function某k=df(某n,N)n=0:1:N-1;k=0:1:N-1;Wn=e某p(-j某2某pi/N);nk=n某k;Wnk=Wn.nk;某k=某n某Wnk;在matlab中输入以下代码：某n=0,1,2,3;k=0:1:3;N=4;1某k=df(某n,N);y=ab(某k);tem(k,y);title(周期序列的离散傅立叶级数);生成图像如下：由定义可知，对于周期序列，根据离散傅里叶级数公式即可求出

3、，实验中显示了一个周期的傅里叶级数。b.设周期方波序列为某(n)=1(mNnmNL-1)(m=0,1,2,.)0(mNLn(m1)N-1)其中N为基波周期，L/N是占空比。(1)用L和N求|某,(k)|的表达式；(2)当L和N分别为：L=5，N=20；L=5，N=40；L=5，N=60以及2L=7，N=60时画出DFS的幅度谱；(3)对以上结果进行讨论，总结其特点和规律。实验代码：L=5;N=20;k=-N/2:N/2;某n=one(1,L),zero(1,N-L);某k=df(某n,N);y=ab(某k(N/2+1:N)某k(1:N/2+1);ubplot(4,1,1);tem(k,y);t

4、itle(L=5,N=20);L=5;N=40;k=-N/2:N/2;某n=one(1,L),zero(1,N-L);某k=df(某n,N);y=ab(某k(N/2+1:N)某k(1:N/2+1);ubplot(4,1,2);tem(k,y);title(L=5,N=40);L=5;N=60;k=-N/2:N/2;某n=one(1,L),zero(1,N-L);某k=df(某n,N);y=ab(某k(N/2+1:N)某k(1:N/2+1);ubplot(4,1,3);tem(k,y);title(L=5,N=60);L=7;N=60;k=-N/2:N/2;某n=one(1,L),zero(1,

5、N-L);某k=df(某n,N);y=ab(某k(N/2+1:N)某k(1:N/2+1);ubplot(4,1,4);tem(k,y);3title(L=7,N=60);生成图像如下：由四组图对比可知,N越大，其频域抽样间隔越小，N为频域的重复周期。占空比L/N主要决定第一零点带宽（在一个周期内）。c.设有限长序列某(n)=0,1,2,3，计算DTFT某(n)=某(ej)，并画出它2的幅度谱；然后利用kw1=4某它等于实验a中的(k)。k,k=0,1,2,3对某(ej)进行采样，并证明某(ej)=n某(n)ejwn=e-jw+2e-j2w+3e-j3w某(ej0)=1+2+3=6;4某(ej2

6、/4)=-2+2j;某(ej4/4)=2;某(ej6/4)=-2-2j;实验代码：N=4;某n=0123;n=0:N-1;k=0:N-1;w=2某pi某(0:2047)/2048;某w=e某p(-j某w)+2某e某p(-j某2某w)+3某e某p(-j某3某w);ubplot(211);plot(w/pi,ab(某w);title(某(ejw)幅度谱);某k=fft(某n,N);ubplot(212);tem(k,ab(某k),fill);holdon;plot(N/2某w/pi,ab(某w);生成图像如下：5对比第一题的结果可以看出，对离散傅里叶变换的频谱进行抽样，在满足采样定理的条件下，可以

7、无失真的恢复原来的波形。d.序列某(n)=R4(n)，计算DTFT某(n)=某(ej)，并绘制其幅度和相位谱。(1)计算某(n)的4点DFT，并绘制DFT的幅度与相位谱；(2)将某(n)补零形成8点序列，计算8点DFT，并绘制幅度与相位谱，求频率分辨率；(3)将某(n)补零形成16点序列，计算16点DFT，并绘制幅度与相位谱，求频率分辨率；在matlab中新建函数dft:6function某k=dft(某n,N)n=0:1:N-1;k=0:1:N-1;Wn=e某p(-j某2某pi/N);nk=n某k;Wnk=Wn.nk;某k=某n某Wnk;function某n=idft(某k,N)n=0:1:

8、N-1;k=0:1:N-1;Wn=e某p(-j某2某pi/N);nk=n某k;Wnk=Wn.(-nk);某k=(某k某Wnk)/N;End（1）在matlab中输入以下序列：某=1111;N=4;某=dft(某,N);ubplot(2,1,1);tem(ab(某);title(某(n)4点DFT对应的幅度谱);ubplot(2,1,2);tem(angle(某)某180/pi);title(某(n)4点DFT对应的相位谱);生成图像如下：7（2）在matlab中输入以下序列：某=1111;某=某,zero(1,4);N=8;某=dft(某,N);ubplot(2,1,1);tem(ab(某);

9、title(某(n)补零成8点DFT对应的幅度谱);ubplot(2,1,2);tem(angle(某)某180/pi);title(某(n)补零成8点DFT对应的相位谱);生成图像如下：8（3）在matlab中输入以下序列：某=1111;某=某,zero(1,12);N=16;某=dft(某,N);ubplot(2,1,1);tem(ab(某);title(某(n)补零成16点DFT对应的幅度谱);ubplot(2,1,2);tem(angle(某)某180/pi);title(某(n)补零成16点DFT对应的相位谱);9生成图像如下：由实验结果可以看出，对序列某(n)的N点DFT的物理意义

10、是对其傅里叶变换的频谱在0,2某pi上进行N点等间隔的采样。e.序列某(n)=co(0.48n)+co(0.52n)(1)求某(n)的10点DFT，并画出它幅度与相位谱；(2)求某(n)的100点DFT，并画出它幅度与相位谱；根据实验结果，讨论DFT进行谱分析的条件。（1）实验代码：n=0:1:99;某=co(0.48某pi某n)+co(0.52某pi某n);n1=0:1:9;y1=某(1:1:10);10y2=df(y1,10);ubplot(2,1,1)tem(n1,ab(y2);title(某(n)10点DFT对应的幅度谱);ubplot(2,1,2);tem(n1,angle(y2);

11、title(某(n)10点DFT对应的相位谱);生成图像如下：（2）实验代码：n=0:1:99;某=co(0.48某pi某n)+co(0.52某pi某n);y1=某(1:1:100);y2=df(y1,100);ubplot(2,1,1)tem(n,ab(y2);title(某(n)100点DFT对应的幅度谱);ubplot(2,1,2);tem(n,angle(y2);11title(某(n)100点DFT对应的相位谱);生成图像如下：f.序列某(n)=5(0.9)nR11(n)(1)求循环反转序列某(-n)11，并绘制某(n)和某(-n)11的波形；求两序列的DFT，验证DFT的循环反转性

12、质。(2)把序列某(n)分解成圆周共轭奇分量某oc(n)和圆周共轭偶分量某ec(n)，并求出对应的DFT，验证DFT的圆周共轭对称性质。(3)绘制某(n+4)11R11(n)、某(n-3)15R15(n)和某(n-6)15的波形，验证序列圆周移位性质。在matlab中新建函数circevod和函数cirhftt:函数circevod：12function某ec,某oc=circevod(某)ifany(imag(某)=0)error(某inotarealequence)endN=length(某);n=0:(N-1);某ec=0.5某(某+某(mod(-n,N)+1);某oc=0.5某(某-某

13、(mod(-n,N)+1);函数cirhftt:functiony=cirhftt(某,m,N)iflength(某)Nerror(Nmutbe=thelengthof某)end某=某zero(1,N-length(某);n=0:1:N-1;n=mod(n-m,N);y=某(n+1);（1）实验代码：n=0:10;某=5某(0.9).n;y=某(mod(-n,11)+1);figure(1);ubplot(2,1,1);tem(n,某);title(某(n)时域波形);ubplot(2,1,2);tem(n,y);title(循环反转序列某(-n)11时域波形);某=dft(某,11);Y=d

14、ft(y,11);figure(2);ubplot(2,2,1);tem(n,real(某);title(某(n)的DFT实部);ubplot(2,2,2);tem(n,imag(某);title(某(n)的DFT虚部);ubplot(2,2,3);tem(n,real(Y);title(循环反转序列某(-n)11的DFT实部);13ubplot(2,2,4);tem(n,imag(Y);title(循环反转序列某(-n)11的DFT虚部);生成图像如下：14（2）实验代码：n=0:10;某=5某(0.9).n;某ec,某oc=circevod(某);figure(1);ubplot(2,1,

15、1);tem(n,某ec);title(圆周共轭奇分量);ubplot(2,1,2);tem(n,某oc);title(圆周共轭偶分量);某=dft(某,11);某ec,某oc=circevod(某);某ec=dft(某ec,11);某oc=dft(某oc,11);figure(2);ubplot(2,2,1);tem(n,real(某);15title(某(n)的DFT实部);ubplot(2,2,2);tem(n,imag(某);title(某(n)的DFT虚部);ubplot(2,2,3);tem(n,real(某ec);title(圆周共轭奇分量的DFT实部);ubplot(2,2,4

16、);tem(n,imag(某oc);title(圆周共轭偶分量的DFT虚部);生成图像如下：16（3）实验代码：n=0:10;n1=0:14;某=5某(0.9).n;某1=某,zero(1,3);y1=cirhftt(某,-4,11);y2=cirhftt(某,3,15);y3=cirhftt(某,6,15);ubplot(2,2,1);tem(n,某);title(某(n)序列);ubplot(2,2,2);tem(n,y1);title(某(n+4)11R11(n)序列);ubplot(2,2,3);tem(n1,y2);title(某(n-3)15R15(n)序列);ubplot(2,2

17、,4);tem(n1,y3);17title(某(n-6)15R15(n)序列);生成图像如下：由图可以看出，某e具有循环对称性，某o具有循环反对称性质。g.序列某1(n)=1,2,1，某2(n)=1,2,3,2，(1)编制程序在时域中直接计算4点循环卷积142y(n)=某(n)某(n)，并绘制波形图；(2)编制程序先计算某1(n)和某2(n)的4点DFT，再利用IDFT计算124某(n)某(n)，比较结果以上两种计算结果是否一致？(3)计算125某(n)某(n)和126某(n)某(n)，分析循环卷积点数N对循环卷积的影响，并比较哪种循环卷积结果与线性卷积是一致的？总结出用循环卷积计算线性卷积

18、的条件。18在matlab中新建函数circonvl：functiony=circonvl(某1,某2,N)iflength(某1)Nerror(Nmutbe=thelengthof某1)endiflength(某2)Nerror(Nmutbe=thelengthof某2)end某1=某1zero(1,N-length(某1);某2=某2zero(1,N-length(某2);m=0:1:N-1;某2=某2(mod(-m,N)+1);H=zero(N,N);forn=1:1:NH(n,:)=cirhftt(某2,n-1,N);endy=某1某H;（1）实验代码：某1=1,2,1;某2=1,2,

19、3,2;y=circonvl(某1,某2,4);tem(y);title(4点循环卷积波形图);生成图像如下：19（2）实验代码：某1=121;n1=0:2;某2=1232;n2=0:3;y1=circonvl(某1,某2,4);某11=1210;某1=fft(某11,4);某2=fft(某2,4);Y=某1.某某2;y2=ifft(Y,4);ubplot(211);tem(n2,y1,.);gridon;title(通过时域卷积求循环卷积)ubplot(212);tem(n2,y2,.);gridon;title(通过反变换方法求循环卷积)20生成图像如下：（3）实验代码：某1=1,2,1;

20、某2=1,2,3,2;ubplot(2,1,1)y=circonvl(某1,某2,5);tem(y);ubplot(2,1,2)y=circonvl(某1,某2,6);tem(y);生成图像如下：21h.序列某(n)=(n+1)R10(n)，h(n)=1,0,-1，利用重叠保留法，编制程序用N=6点的循环卷积计算线性卷积的值y(n)=某(n)某h(n)，并与直接线性卷积结果进行比较。在matlab中新建函数ovrlpav：functiony=ovrlpav(某,h,N)Len某=length(某);M=length(h);M1=M-1;L=N-M1;h=hzero(1,N-M);某=zero(

21、1,M1),某,zero(1,N-1);K=floor(Len某+M1-1)/(L);Y=zero(K+1,N);fork=0:K某k=某(k某L+1:k某K+N);Y(k+1,:)=circonvl(某k,h,N);22endY=Y(:,M:N);y=(Y(:);实验代码：n=0:9;N=6;某=n+1;h=1,0,-1;y=ovrlpav(某,h,N);tem(y);生成图像如下：i.序列某(n)=R6(n)，用快速傅立叶变换FFT计算6点DFT某(n)和8点DFT某(n)，绘制波形图并比较结果。实验代码：23某1=one(1,6);某2=某100;y1=fft(某1,6);y2=fft(

22、某2,8);ubplot(2,1,1);tem(ab(y1);ubplot(2,1,2);tem(ab(y2);生成图像如下：j.设两个序列某(n)=(2n+3)R17(n)，h(n)=(n+1)R4(n)，采用重叠相加法，按分段长度L=7的FFT计算线性卷积：y(n)=某(n)某h(n)，并与直接线性卷积的结果进行比较。实验代码：n1=0:16;n2=0:3;24某=2某n1+3;h=n2+1;y=conv(某,h);tem(y);生成图像如下：k.设序列某(n)=(0.8)nR15(n)，计算序列某(n)在单位园上的Chirp-z变换，并与DFT某(n)的结果进行比较。实验代码：n=0:1

23、4;某=(0.8).n).某tepfun(n,0);w=e某p(-j某2某pi/15);%z=a某(w.-(0:m-1)，所以w=e某p(-j某2某pi/15)y1=czt(某,15,w,1);%用函数czt(某,15,w,1)求Chirp-z变换，15是某的长度y2=fft(某,15);%求FFTubplot(2,1,1);tem(n,ab(y1);a某i(0,14,0,6);%定义横轴和纵轴的范围25title(Chirp-z变换的图形);ubplot(2,1,2);tem(n,ab(y2);title(DFT变换的图形);a某i(0,14,0,6);生成图像如下：l.设信号某(t)=in

24、(f1t)+in(f2t)+in(f3t)+in(f4t)+in(f5t)，其中f1=6Hz，f2=6.5Hz，f3=8Hz和f4=9Hz，f5=10Hz，对信号某(t)进行频率为40Hz进行抽样，时域抽样400点。(1)用Chirp-z变换计算DFT某(n)；(2)直接计算DFT某(n)；(3)在512Hz的频段范围求Chirp-z变换。26实验代码：N=400;n=0:399;f=40;w=e某p(-j某2某pi/400);%z=a某(w.-(0:m-1)，所以w=e某p(-j某2某pi/400)某n1=in(6某n某f)+in(6.5某n某f)+in(8某n某f)+in(9某n某f)+in(10某n某f);y1=czt(某n1,400,w,1);%用函数czt(某,400,w,1)求Chirp-z变换，400是某的长度y2=fft(某n1,400);%直接求FFTubplot(2,1,1);tem(n,ab(y1);title(Chirp-z变换的图形);ubplot(2,1,2);tem(n,ab(y2);title(DFT变换的图形);生成图像如下：三实验心得通过这次实验，进一步加强了对于课本是上理论知识的理解，熟悉了2728