DFT和FFT实验上传

1、装 订 线DFT和FFT实验一、实验目的和要求1、掌握DFT变换2、掌握DFT性质3、掌握快速傅立叶变换（FFT）二、实验内容和原理1、实验内容1)求有限长离散时间信号的离散时间傅立叶变换并绘图。· 已知 · 已知2)已知序列，绘制及其离散傅立叶变换的幅度、相位图。3)设，,其中，randn(n)为高斯白噪声。求出，m=2,3,4的matlab采用不同算法的执行时间。4)研究高密度频谱和高分辨率频谱。 设有连续信号· 以采样频率对信号x(t)采样，分析下列三种情况的幅频特性。· 采集数据长度N=16点，做N=16点的FFT，并画出幅频特性。·

2、采集数据长度N=16点，补零到256点，做N=256点的FFT，并画出幅频特性。· 采集数据长度N=256点，做N=256点的FFT，并画出幅频特性。观察三种不同频率特性图，分析和比较它们的特点以及形成的原因。2、实验原理1）DFT序列x(n)的离散时间傅里叶变换（DTFT）表示为，如果x(n)为因果有限长序列，n=0,1,.,N-1，则x(n)的DTFT表示为x(n)的离散傅里叶变换（DFT）表达式为序列的N点DFT是序列DTFT在频率区间0,2上的N点灯间隔采样，采样间隔为2/N。通过DFT，可以完成由一组有限个信号采样值x(n)直接计算得到一组有限个频谱采样值X(k)。X(k)

3、的幅度谱为，其中下标R和I分别表示取实部、虚部的运算。X(k)的相位谱为。离散傅里叶反变换（IDFT）定义为。2）FFT快速傅里叶变换（FFT）是DFT的快速算法，并不是一个新的映射。FFT利用了函数的周期性和对称性以及一些特殊值来减少DFT的运算量，可使DFT的运算量下降几个数量级，从而使数字信号处理的速度大大提高。若信号是连续信号，用FFT进行谱分析时，首先必须对信号进行采样，使之变成离散信号，然后就可以用FFT来对连续信号进行谱分析。为了满足采样定理，一般在采样之前要设置一个抗混叠低通滤波器，且抗混叠滤波器的截止频率不得高于与采样频率的一半。比较DFT和IDFT的定义，两者的区别仅在于指

4、数因子的指数部分的符号差异和幅度尺度变换，因此可用FFT算法来计算IDFT。三、主要仪器设备Matlab四、操作方法和实验步骤1、认真分析原函数，取点2、用matlab编写程序，运行程序得出结果五、实验数据记录、处理和分析1、求有限长离散时间信号的离散时间傅立叶变换并绘图。· 已知 · 已知【解答】思路：这是一道DFT的题，按照题目要求只需要取11个点即可。第（1）小题M文件源代码N=11; %取点个数为11个j=sqrt(-1); %定义j为复数单位f=inline('(0.9*exp(j*pi/3)n','n'); %定义一个函数f(n)

5、W=0:2*pi/1000:2*pi; %定义离散域的基本频率W为数组，间距为2*pi/1000Xw=zeros(size(W); %定义一个与W位数相等的数组 for n=0:N-1 Xw=Xw+f(n)*exp(-j*W*n); end %对f(n)函数做DFT变换 xn=;for n=0:N-1 xn(n+1)=f(n); end %将f(n)的值放进数组xn里面，便于最后画出xn的图像magXw=abs(Xw); %定义一个数组magXw,将abw(Xw)的值赋给它angleXw= angle(Xw); %定义数组angleXw,将angle(Xw)的值赋给它figure(1); pl

6、ot(xn,'.-');xlabel('n'); ylabel('x(n)'); %画出xn的图figure(2);k=0:1:N-1;plot(W,magXw,'-') xlabel('W'); ylabel('|X(W)|'); %画出magXw的图像figure(3);plot(W,angleXw,'-');xlabel('W'); ylabel('angle(X(W)'); %画出angleXw的图像运行结果xn图像：X(W)的幅度图X(W)的

7、相位图 【分析】X=n=010(0.9ej3)n·e-jn=1-(0.9ej3-)111-0.9ej3-，X=|1-0.9ej3-11|1-0.9ej3-|可见X的幅度频谱有11-1=10个极大，11-1=10个极小。而X的相位则有11-1=10个极大，11-1=10个极小，并且相位在-2和2之间摆动。第（2）小题M文件源代码N=10; %取点个数为11个j=sqrt(-1); %定义j为复数单位f=inline('2n','n'); %定义一个函数f(n)W=0:2*pi/1000:2*pi; %定义离散域的基本频率,将其设置为间距为2*pi/100

8、0的数组Xw=zeros(size(W); %定义一个数组Xw，位数与W相等 for n=-N:N Xw=Xw+f(n)*exp(-j*W*n); end %对f(n)函数做DTFT变换 xn=; %定义一个数组xnfor n=-N:N xn(n+1+N)=f(n);end %将f(n)的值放进数组xn里面，便于最后画出xn的图像 magXw=abs(Xw); %将X(W)的模值放进数组magXw angleXw=angle(Xw); %将X(W)的相位放进数组magXwfigure(1); plot(xn,'.-');xlabel('n'); ylabel(&

9、#39;x(n)'); %画出xn的图figure(2);plot(W,magXw,'-') xlabel('W'); ylabel('|X(w)|'); %画出magXw的图像figure(3);plot(W,angleXw,'-');xlabel('W'); ylabel('angle(X(w)'); %画出angleXw的图像xn图像X(W)的幅度图X(W)的相位图【分析】X=n=-10102n·e-jn=2e-j-101-2e-j211-2e-j，|X|=2-10|1-2e

10、-j21|1-2e-j|可见X的幅度有一个极大值，一个极小值。X的相位在-到之间来回振动，并且中间出现突变的情况。2)已知序列，绘制及其离散傅立叶变换的幅度、相位图。【解答】思路：这是一道DFT的题，按照题目要求只需要取51个点即可。M文件源代码N=51; %取点个数为50个j=sqrt(-1); %定义j为复数单位f=inline('cos(0.82*pi*n)+2*sin(pi*n)','n'); %定义一个函数f(n)Xk=; %定义一个数组XkW=2*pi/N; %定义离散域的基本频率for k=1:N Xk(k)=0; for n=0:N-1 Xk(k

11、)=Xk(k)+f(n)*exp(-j*(k-1)*W*n); endend %对f(n)函数做DFT变换xn=; %定义一个数组xnfor n=0:N-1 xn(n+1)=f(n);end %将f(n)的值放进数组xn里面，便于最后画出xn的图像magXk=; %定义一个数组magXkfor k=1:N magXk(k)=abs(Xk(k);end %将X(kW)的模值放进数组magXkangleXk=; %定义数组angleXkfor k=1:N angleXk(k)=angle(Xk(k);end %将X(kW)的相位放进数组magXkfigure(1);plot(xn,'.-&

12、#39;);xlabel('n'); ylabel('x(n)'); %画出xn的图figure(2);k=0:1:N-1;plot(k,magXk,'+-') xlabel('k'); ylabel('|X(k)|'); %画出magXk的图像figure(3);plot(k,angleXk,'x-');xlabel('k'); ylabel('angle(X(k)'); %画出angleXk的图像命令窗口中的运行及其结果：xn图像Xk的幅度图Xk的相位图【分析】X

13、k=n=050cos0.82n+2sinne-jk251n=n=050ej0.82n+e-j0.82n2e-jk251n=12n=050(ej0.82-k251n+e-j0.82+k251n)=12(1-(ej0.82-k251)511-ej0.82-k251+1-e-j0.82+k251511-e-j0.82+k251)可见Xk的幅度频谱拥有两个峰值，Xk的相位频谱在-到之间来回振动，且中间存在3个台阶式的向下跳变，一个台阶式的向上跳变。3设，,其中，randn(n)为高斯白噪声。求出，m=2,3,4的matlab采用不同算法的执行时间。【解答】思路：计算DFT算法和FFT算法的运行时间可以

14、使用的etime函数M文件源代码m=input('m=:'); %输入m值N=4m; %求出所取得x(n)的点数j=sqrt(-1); %定义j为复数单位arr=; %定义一个数组arrW=2*pi/N; %定义离散域的基本频率dft_time=0; %定义dft_time为0t1=clock; %此处为dft计算的时间起点for k=1:N arr(k)=0; for n=0:N-1 arr(k)=arr(k)+(sin(0.2*pi*n)+rand(1)*exp(-j*(k-1)*W*n); endend %对x(n)做dft变换dft_time=etime(clock,t

15、1) %得出dft变换所花的时间Wn=exp(-j*2*pi/N); %求出旋转因子 xn=zeros(1,N); %定义一个N位数组fft_time=0; %定义fft_time为0t2=clock; %此处为fft计算的时间起点for n=0:N-1 xn(n+1)=sin(0.2*pi*n)+randn(1); %将x(n)的值放入数组xn中 end n1=fliplr(dec2bin(0:N-1); %码位倒置步骤1：将码位转换为二进制，再进行倒序n2=bin2dec(n1); %码位倒置步骤2：将码位转换为十进制后翻转x=zeros(2*m+1,N); %定义一个（2m+1)

16、5;N的矩阵for i=1:N x(1,i)=xn(n2(i)+1);end %将码位倒序后的值重新赋值进入矩阵第一列for i=1:2*m %进行第v级计算 Number=2(2*m-i); Interval_of_Unit=2(i-1); %每组中每个计算单元的间距 Interval_of_Group=2i; %每组之间的间距 Wnr=; %定义一个新数组Wnr for r=1:2(i-1) Wnr(r)=Wn(r-1)*N/2i); end %将每一级运算的指数因子赋值给Wnr for k=0:Number-1 for l=1:2(i-1) x(i+1,l+k*2i)=x(i,l+k*2

17、i)+Wnr(l)*x(i,l+k*2i+2(i-1); x(i+1,l+k*2i+2(i-1)=x(i,l+k*2i)-Wnr(l)*x(i,l+k*2i+2(i-1); end end end %对x(n)做fft变换fft_time=etime(clock,t2) %得出dft变换所需要的时间命令窗口中的运行及其结果m=:2dft_time = 0fft_time = 0m=:3dft_time = 0.02200000000000fft_time = 0.00400000000000m=:4dft_time =fft_time =【分析】从实验结果可以看出，对于同样的x(n)，FFT变

18、换的计算时间小于DFT变换的计算时间，并且随着取样点数的增多，这种差距越来越明显。这是由于FFT算法利用了旋转因子的对称性和周期性，将长序列分解为短序列，分级进行蝶形运算使得复数乘法的次数减少到分解前的一半。通过对短序列的计算进行适当的组合，从而达到了删除重复运算，提供运算速度的目的。4)研究高密度频谱和高分辨率频谱。 设有连续信号· 以采样频率对信号x(t)采样，分析下列三种情况的幅频特性。· 采集数据长度N=16点，做N=16点的FFT，并画出幅频特性。· 采集数据长度N=16点，补零到256点，做N=256点的FFT，并画出幅频特性。· 采集数据长

19、度N=256点，做N=256点的FFT，并画出幅频特性。观察三种不同频率特性图，分析和比较它们的特点以及形成的原因。【解答】思路：此处为FFT变换，要注意采样数据长度N和FFT变换点数M不同，因此当M>N时，应该将x(n)中非0-N-1范围内的点补零。M文件源代码N=input('N=:'); %输入采样数据长度NM=input('M=:'); %输入FFT变换点数Mfs=32000;v=log2(M); %给采样频率fs赋值为32000，v为FFT变换的级数Wn=exp(-j*2*pi/M); %定义旋转因子Wnj=sqrt(-1); %定义复数单位jx

20、n=zeros(1,M); %定义一个M位数组xnfor n=0:N-1xn(n+1)=cos(2*pi*6.5*1000*n/fs)+cos(2*pi*7*1000*n/fs)+cos(2*pi*9*1000*n/fs);end %将x(n)的值赋值进入xn中n1=fliplr(dec2bin(0:M-1);%码位倒置步骤1：将码位转换为二进制，再进行倒序n2=bin2dec(n1);%码位倒置步骤2：将码位转换为十进制后翻转x=zeros(v+1,M); %定义一个（v+1)×M的矩阵xfor i=1:M x(1,i)=xn(n2(i)+1);end %将码位倒序后的x(n)值赋

21、值进入矩阵x的第一行for i=1:v %进行第v级计算 Number=2(v-i); %每一级计算中的“群”数 Interval_of_Unit=2(i-1); %每组中每个计算单元的间距 Interval_of_Group=2i; %每组之间的间距 Wnr=; %定义一个数组Wnr for r=1:2(i-1) Wnr(r)=Wn(r-1)*M/2i); end %将每一级运算的指数因子赋值给Wnr for k=0:Number-1 for l=1:2(i-1) x(i+1,l+k*2i)=x(i,l+k*2i)+Wnr(l)*x(i,l+k*2i+2(i-1); x(i+1,l+k*2i+2(i-1)=x(i,l+k*2i)-Wnr(l)*x(i,l+k*2i+2(i-1); end end end %对x(n)做fft变换Xk=; %定义一个数组Xkfor k=