DFT和FFT实验(上传)

1、DFT和FFT实验一、实验目的和要求1、掌握DFT变换2、掌握DFT性质3、掌握快速傅立叶变换(FFT)二、实验容和原理1、实验容,1)求有限长离散时间信号的离散时间傅立叶变换X(ej )并绘图。装n订？已知 x(n) 0.9e 30 n 10线？已知 x(n) 2n10 n 10J 2)已知序列x(n) cos(0.82 n) 2sin( n), 0 n 50,绘制x(n)及其离散傅立叶变换 X(k)的幅度、：相位图。r:3)设 x(n) sin(0.2 n) randn(n) , 0 n N 1,其中，randn(n)为高斯白噪声。求出 N4m , m=2,3,4的matlab采用不同算法

2、的执行时间。94)研究高密度频谱和高分辨率频谱。：设有连续信号x(t) cos(2 6.5 103t) cos(2 7 103t) cos(2 9 103t)?以采样频率fs 32kHz对信号x(t)采样，分析下列三种情况的幅频特性。9?采集数据长度N=16点，做N=16点的FFT,并画出幅频特性。?采集数据长度N=16点，补零到256点，做N=256点的FFT,并画出幅频特性。:?采集数据长度N=256点，做N=256点的FFT,并画出幅频特性。:观察三种不同频率特性图，分析和比较它们的特点以及形成的原因。2、实验原理1) DFT序列x(n)的离散时间傅里叶变换(DTFT)表示为X(ej )

3、x( n)e jnN 1如果x(n)为因果有限长序列，n=0,1,.,N-1 ,则x(n)的DTFT表示为X(ej ) x(n)e n 0N 1j 2 nkx(n)的离散傅里叶变换(DFT)表达式为X(k) x(n)e N (k 0,1，N 1)n 0序列的N点DFT是序列DTFT在频率区间0,2可上的N点灯间隔采样，采样间隔为2”N。通过DFT ,可以完成由一组有限个信号采样值x(n)直接计算得到一组有限个频谱采样值X(k)o X(k)的幅度谱为X(k) VXR(k)X2(k),其中下标R和I分别表示取实部、虚部的运算。X(k)的相位谱为(k) arctanXR(k)一、4，1 N离散傅里叶

4、反变换(IDFT)定义为x(n) 一N n1j_2_nkX(k)e N (n 0,1,., N 1)。02) FFTj2-n快速傅里叶变换(FFT)是DFT的快速算法，并不是一个新的映射。FFT利用了 e N函数的周期性和对称性以及一些特殊值来减少DFT的运算量，可使 DFT的运算量下降几个数量级，从而使数字信号处理的速度大大提高。若信号是连续信号,用FFT进行谱分析时，首先必须对信号进行采样，使之变成离散信号， 然后就可以用FFT来对连续信号进行谱分析。为了满足采样定理，一般在采样之前要设置一个抗混叠低通滤波器, 且抗混叠滤波器的截止频率不得高于与采样频率的一半。比较DFT和IDFT的定义，

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

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

7、(Xw)的值赋给它figure。)；plot(xn,'.-');xlabel('n'); ylabel('x(n)'); %figure(2);k=0:1:N-1;plot(W,magXw,'-') xlabel('W); ylabel('|X(W)|'); figure(3);plot(W,angleXw,'-');xlabel(W); ylabel('angle(X(W)'); %运行结果xn图像：0.8-0.6画出xn的图画出magXw的图像画出angleXw 的图像0

8、.60.40.20-0.2-0.4-0.8-0.6-0.4-0.200.20.40.60.81nX(W)的幅度图X(W)的相位图)a【分析】11-1=10可见整的幅度频谱有ii-i=i0个极大，ii-i=i0个极小。而的相位则有ii-i=i0个极大, 个极小，并且相位在-区和区之间摆动。第(2)小题M文件源代码N=10;%j=sqrt(-1);%f=inline('2An','n'); %W=0:2*pi/1000:2*pi; %Xw=zeros(size(W); % for n=-N:N取点个数为11个定义j为复数单位定义一个函数f(n)定义离散域的基本频率，

9、将其设置为间距为2*pi/1000 的数组定义一个数组Xw,位数与W相等Xw=Xw+f(n)*exp(-j*W*n);end %对f(n) 函数做 DTFT变换xn=;%for n=-N:Nxn(n+1+N)=f(n);end%magXw=abs(Xw);angleXw=angle(Xw);figure。);plot(xn,'.-');xlabel('n'); ylabel('x(n)');figure(2);plot(W,magXw,'-')定义一个数组xn将f(n)的值放进数组xn里面，便于最后画出 xn的图像%将X(W)的模

10、值放进数组magXw%将X(W)的相位放进数组 magXw%画出xn的图xlabel('W); ylabel('|X(w)|'); % figure(3);画出magXw的图像plot(W,angleXw,'-');xlabel('W); ylabel('angle(X(w)'); %画出angleXw 的图像xn图像1200100080060040020010152025nX(W)的幅度图2200200018001600140012001000800600X(W)的相位图4-401234567【分析】曲卬） 1110x（n）=

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

12、i/N;%定义离散域的基本频率for k=1:NXk(k)=0;for n=0:N-1Xk(k)=Xk(k)+f(n)*exp(-j*(k-1)*W*n);endend%xn=;%for n=0:N-1xn(n+1)=f(n);对 f(n) 函数做DFT 变换定义一个数组xnend%magXk=;%for k=1:NmagXk(k)=abs(Xk(k);end%angleXk=;%for k=1:NangleXk(k)=angle(Xk(k);end%figure(1);plot(xn,'.-');xlabel('n'); ylabel('x(n)

13、9;); % figure(2);k=0:1:N-1;plot(k,magXk,'+-') xlabel('k'); ylabel('|X(k)|'); % figure(3);plot(k,angleXk,'x-');将 f(n) 的值放进数组xn 里面，便于最后画出定义一个数组magXk将 X(kW) 的模值放进数组magXk定义数组angleXk将 X(kW) 的相位放进数组magXk画出xn 的图画出magXk 的图像xnxlabel('k'); ylabel('angle(X(k)');

14、命令窗口中的运行及其结果：%画出angleXk 的图像xn 图像的图像10.80.60.40.20-0.2-0.4-0.6-0.8-1102030405060Xk的幅度图30252015100T11 l lin l【l I1_l I U L U L 05101520253035404550Xk的相位图【分析】ejt).B2Jm + e-p.82mi 一四色/c 51部)=y(情等cosO.827m + 2simm)e = J小勃+ 卦1 1 (J,必/3)>1:小如十碟5171。必+或)可见的幅度频谱拥有两个峰值，尺(的的相位频谱在一汽到冗之间来回振动，且中间存在3个台阶式的向 下跳变，

15、一个台阶式的向上跳变。3 .设 x(n) sin(0.2 n) randn(n) , 0 n N 1,其中，randn(n)为高斯白噪声。求出 N4m , m=2,3,4的matlab采用不同算法的执行时间。【解答】思路：计算DFT算法和FFT算法的运行时间可以使用的etime函数M文件源代码m=input('m=:');%输入m值N=4Am;%求出所取得x(n)的点数j=sqrt(-1);%定义j为复数单位arr=;%定义一个数组arrW=2*pi/N;%定义离散域的基本频率dft_time=0;%定义dft_time为0t1=clock;%此处为dft计算的时间起点for

16、k=1:Narr(k)=0;for n=0:N-1arr(k)=arr(k)+(sin(0.2*pi*n)+rand(1)*exp(-j*(k-1)*W*n); endend%dft_time=etime(clock,t1) %Wn=exp(-j*2*pi/N); %xn=zeros(1,N); %fft_time=0;%t2=clock;%for n=0:N-1对 x(n) 做 dft 变换得出 dft 变换所花的时间求出旋转因子定义一个N 位数组定义 fft_time为 0此处为 fft 计算的时间起点xn(n+1)=sin(0.2*pi*n)+randn(1);%end将 x(n) 的值

17、放入数组xn 中n1=fliplr(dec2bin(0:N-1);%n2=bin2dec(n1);%x=zeros(2*m+1,N);%for i=1:Nx(1,i)=xn(n2(i)+1);end%for i=1:2*m%Number=2A(2*m-i);Interval_of_Unit=2A(i-1); % Interval_of_Group=2Ai; % Wnr=;%for r=1:2A(i-1)Wnr(r)=WnA(r-1)*N/2Ai);end%码位倒置步骤1 ：将码位转换为二进制，再进行倒序码位倒置步骤2：将码位转换为十进制后翻转定义一个(2m+1)XN的矩阵将码位倒序后的值重新赋

18、值进入矩阵第一列进行第 v 级计算每组中每个计算单元的间距每组之间的间距定义一个新数组Wnr将每一级运算的指数因子赋值给Wnrfor k=0:Number-1for l=1:2A(i-1)x(i+1,l+k*2Ai)=x(i,l+k*2Ai)+Wnr(l)*x(i,l+k*2Ai+2A(i-1);x(i+1,l+k*2Ai+2A(i-1)=x(i,l+k*2Ai)-Wnr(l)*x(i,l+k*2Ai+2A(i-1);end对 x(n) 做 fft 变换得出 dft 变换所需要的时间endend%fft_time=etime(clock,t2) %命令窗口中的运行及其结果m=:2dft_tim

19、e =0fft_time =0 m=:3dft_time = 0.000ffttime =0.000 m=:4dft_time =0.000fft_time =0.000【分析】从实验结果可以看出，对于同样的x(n) ， FFT 变换的计算时间小于DFT 变换的计算时间，并且随着取样点数的增多，这种差距越来越明显。这是由于FFT 算法利用了旋转因子的对称性和周期性，将长序列分解为短序列，分级进行蝶形运算使得复数乘法的次数减少到分解前的一半。通过对短序列的计算进行适当的组合，从而达到了删除重复运算，提供运算速度的目的。4 ) 研究高密度频谱和高分辨率频谱。设有连续信号x(t) cos(2 6.5

20、 103t) cos(2 7 103t) cos(2 9 103t)?以采样频率fs 32kHz对信号x(t)采样，分析下列三种情况的幅频特性。?采集数据长度N=16点，做N=16点的FFT,并画出幅频特性。?采集数据长度N=16点，补零到256点，做N=256点的FFT,并画出幅频特性。?采集数据长度N=256点，做N=256点的FFT,并画出幅频特性。观察三种不同频率特性图，分析和比较它们的特点以及形成的原因。【解答】 思路：此处为FFT变换，要注意采样数据长度 N和FFT变换点数M不同，因此当M>N时，应该将x(n)中非 0-N-1 围的点补零。M文件源代码N=input('

21、;N=:');%M=input('M=:');%fs=32000;v=log2(M); %输入采样数据长度N输入 FFT 变换点数M给采样频率fs 赋值为 32000 ， v 为 FFT 变换的级数Wn=exp(-j*2*pi/M); %j=sqrt(-1);%xn=zeros(1,M); %定义旋转因子Wn定义复数单位j定义一个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=flip

22、lr(dec2bin(0:M-1);%码位倒置步骤1 ：将码位转换为二进制，再进行倒序n2=bin2dec(n1);%x=zeros(v+1,M); %for i=1:M码位倒置步骤2 ：将码位转换为十进制后翻转定义一个(v+1) XM的矩阵xx(1,i)=xn(n2(i)+1);end %将码位倒序后的for i=1:v % 进行第 v 级计算x(n) 值赋值进入矩阵x 的第一行Number=2A(v-i);%Interval_of_Unit=2A(i-1); %Interval_of_Group=2Ai; %Wnr=;%每一级计算中的“ 群 ” 数每组中每个计算单元的间距每组之间的间距定义一个数组Wnrfor r=1:2A(i-1)Wnr(r)=WnA(r-1)*M/2Ai);将每一级运算的指数因子赋值给Wnrend%for k=0:Number-1for l=1:2A(i-1)x(i+1,l+k*2Ai)=x(i,l+k*2Ai)+Wnr(l)*x(i,l+k*2Ai+2A(i-1)；x(i+1,l+k*2Ai+2A(i-1)=x(i,l+k*2Ai)-Wnr(l)*x(i,l+k*2Ai+2A(i-1)； endendend%Xk=;%for k=1:MXk(k)=x(v+1,k);end%figure。);k=