实验三-FFT算法的应用实验报告

成绩

数字信号处理

实验报告

实验名称： FFT算法的应用

实验班级：

姓名：

学号（后两位）：

指导教师：

实验日期： 2010.11.23

实验3 FFT算法的应用

一、实验目的
1、加深对离散信号的DFT的理解；

2、在MATLAB中实现FFT算法。

二、实验原理

N点序列的DFT和IDFT变换定义式如下：

利用旋转因子具有周期性，可以得到快速算法（FFT）。

在MATLAB中，可以用函数x=fft(x,N)和x=ifft(x,N)计算N点序列的DFT正、反变换。

三、预习要求

1、在MATLAB中，熟悉函数fft、ifft的使用；

2、阅读扩展练习中的实例，学习在MATLAB中的实现FFT算法的实现；

3、利用MATLAB编程完成计算，绘出相应图形。并与理论计算相比较，说明实验结果的原因。

例：对连续的单一频率周期信号按采样频率采样，截取长度N分别选N =20和N =16，观察其DFT结果的幅度谱。

解：此时离散序列，即k=8。用MATLAB计算并作图，函数fft用于计算离散傅里叶变换DFT，程序如下：

k=8;
n1=[0:1:19];
xa1=sin(2*pi*n1/k);
subplot(2,2,1)
plot(n1,xa1)
xlabel('t/T');ylabel('x(n)');
xk1=fft(xa1);xk1=abs(xk1);
subplot(2,2,2)
stem(n1,xk1)
xlabel('k');ylabel('X(k)');
n2=[0:1:15];
xa2=sin(2*pi*n2/k);
subplot(2,2,3)
plot(n2,xa2)
xlabel('t/T');ylabel('x(n)');
xk2=fft(xa2);xk2=abs(xk2);
subplot(2,2,4)
stem(n2,xk2)
xlabel('k');ylabel('X(k)');

图1 不同的截取长度的正弦信号及其DFT结果

计算结果示于图1，(a)和(b)分别是N=20时的截取信号和DFT结果，由于截取了两个半周期，频谱出现泄漏；(c) 和(d) 分别是N=16时的截取信号和DFT结果，由于截取了两个整周期，得到单一谱线的频谱。上述频谱的误差主要是由于时域中对信号的非整周期截断产生的频谱泄漏。

四、实验内容

1、2N点实数序列

N=64。用一个64点的复数FFT程序，一次算出，并绘出的图形。(按照基于2的蝶型结构的递推公式编程)

编程如下：

N=64;n=[0:1:N-1];

n1=2*n;n2=2*n+1;k=[0:1:N-1];

xn1=cos(2*pi/N*7*n1)+1/2*cos(2*pi/N*19*n1);

xn2=cos(2*pi/N*7*n2)+1/2*cos(2*pi/N*19*n2);

XK1=fft(xn1);XK2=fft(xn2);

X1=XK1+exp(-j*pi*k/N).*XK2;

X2=XK1-exp(-j*pi*k/N).*XK2;

X1=[X1 zeros(1,N)];X2=[zeros(1,N) X2];

XK=X1+X2;

k=[0:1:2*N-1];

XK=abs(XK); stem(k,XK);

xlabel('k');ylabel('|X(k)|');

title('X(k)=DFT[x(n)]2N')

2、已知某序列在单位圆上的N=64等分样点的Z变换为:

。

用N点IFFT程序计算出和。

MATLAB编程如下：

N=64;

k=[0:1:63];

xk=1./(1-0.8*exp(-j*2*pi*k/N));

xn=ifft(xk,64);

stem(k,xn)

xlabel('k');ylabel('x(n)');

disp('xn序列');disp(xn);

xn的图像为：

xn序列：

1.0000 - 0.0000i 0.8000 + 0.0000i 0.6400 + 0.0000i 0.5120 + 0.0000i

0.4096 + 0.0000i 0.3277 + 0.0000i 0.2621 + 0.0000i 0.2097 + 0.0000i

0.1678 + 0.0000i 0.1342 + 0.0000i 0.1074 + 0.0000i 0.0859 + 0.0000i

0.0687 - 0.0000i 0.0550 - 0.0000i 0.0440 - 0.0000i 0.0352 - 0.0000i

0.0281 - 0.0000i 0.0225 - 0.0000i 0.0180 - 0.0000i 0.0144 - 0.0000i

0.0115 - 0.0000i 0.0092 - 0.0000i 0.0074 - 0.0000i 0.0059 - 0.0000i

0.0047 - 0.0000i 0.0038 - 0.0000i 0.0030 - 0.0000i 0.0024 - 0.0000i

0.0019 - 0.0000i 0.0015 - 0.0000i 0.0012 - 0.0000i 0.0010 + 0.0000i

0.0008 + 0.0000i 0.0006 + 0.0000i 0.0005 + 0.0000i 0.0004 + 0.0000i

0.0003 + 0.0000i 0.0003 + 0.0000i 0.0002 + 0.0000i 0.0002 + 0.0000i

0.0001 + 0.0000i 0.0001 + 0.0000i 0.0001 + 0.0000i 0.0001 + 0.0000i

0.0001 + 0.0000i 0.0000 + 0.0000i 0.0000 + 0.0000i 0.0000 - 0.0000i

0.0000 + 0.0000i 0.0000 + 0.0000i 0.0000 + 0.0000i 0.0000 - 0.0000i

0.0000 + 0.0000i 0.0000 - 0.0000i 0.0000 - 0.0000i 0.0000 - 0.0000i

0.0000 - 0.0000i 0.0000 - 0.0000i 0.0000 - 0.0000i 0.0000 - 0.0000i

五、结果分析

理论值计算：

因，k=0,1，……，63。所以，其中|z|>0.8，所以，n=0,1,2，……63。当n=0时，x（n）=1.0，当n=1时，x（n）=0.8，等等，

六、心得体会

这次试验加深了离散信号的DFT的理解，了解和掌握了利用MATLAB进行FFT算法，经过前几次的试验，感觉比前几次试验做的稍微好了点，但是还是依然出现些问题，还是需要在课下多加练习。

第二篇：FFT算法实现实验报告

FFT算法实现实验报告

辛旸 PB10210006

实验目的

1、加深对快速傅里叶变换的理解。

2、掌握FFT 算法及其程序的编写。

3、掌握算法性能评测的方法。

实验内容

1.编写自己的FFT算法：

代码如下：

function [ X ] = Sampling( x,N )

%myFFT 实现FFT时域取样算法

% 输出：生成FFT序列X(k),输入：欲变换序列x（n），FFT变换长度N（可缺省）

（1） if ~exist('N','var'); %检查是否有变换长度N输入

（2） N=length(x); %若无，则令N等于序列长度

（3） end

（4） if N<length(x); %如果N小于序列长度，则对序列进行截短

（5） x=x(:,1:N);

（6） else

（7） x=[x,zeros(1,N-length(x))]; %如果N大于序列长度，对序列补零进行延长

（8） end;

（9） for i=1:1:length(x)/2+1; %判断N是2的多少次方

（10） if 2^i>=length(x); %若N不是2的整数幂

（11） N=2^i; %增大N为2的整数幂

（12） break;

（13） end

（14） end

（15） x=[x,zeros(1,N-length(x))]; %确保要变换的序列长度为2^i

（16） k1=zeros(1,N);

（17） X=zeros(1,N);

（18） w=zeros(1,N);

（19） for m=0:1:N-1; %确定反序序列k1和正序序列k的关系

（20） k=m;

（21） for n=i-1:-1:0; %从高位开始依次将各位移至反序位

（22） k1(m+1)=k1(m+1)+fix(k/(2^n))*(2^(i-1-n));

（23） k=rem(k,2^n);

（24） end;

（25） end

（26） for l=1:1:N;

（27） X(k1(l)+1)=x(l); %生成反序序列

（28） w(l)=exp(-1i*2*pi/N*(l-1)); %生成旋转因子

（29） end

（30） for l=0:1:i-1; %控制FFT运算级数

（31） for m=1:1:N; %每一级中有N/2个蝶形运算

（32） if rem((m-1),2^(l+1))<2^l; %找到蝶形运算的上半部分

（33） b=X(m)+X(m+2^l)*w(2^(i-1-l)*rem((m-1),2^l)+1); %将结果暂存至b

（34） X(m+2^l)=a(m)-X(m+2^l)*w(2^(i-1-l)*rem((m-1),2^l)+1);

（35） X(m)=b; %实现原位运算

（36） end

（37） end

（38） end

2.选择实验1中的典型信号序列验证算法的有效性：

为方便比较两个算法，编写了myCompare函数计算两种算法的运行时间，并绘制频谱曲线

代码如下：

function [ t1,t2,e ] = myCompare( x,N )

%myCompare函数：比较自己编写的算法与系统自带算法的差异

% 输入：与变换信号序列x（n）和欲变换长度N

% 输出：自己编写的函数的执行时间t1，系统自带函数的执行时间t2，两者计算序列的差异平方和e tic;

X1=myFFT(x,N);

t1=toc;

tic;

X2=fft(x,N);

t2=toc;

subplot(1,2,1);plot(abs(X1));xlabel('k');ylabel('X(k');title('ÓÃ×Ô¼º±àÐ´µÄº¯ÊýµÃµ½µÄ±ä»»ÐòÁÐÆµÆ×');

subplot(1,2,2);plot(abs(X2));xlabel('k');ylabel('X(k');title('ÏµÍ³×Ô´øFFTº¯ÊýµÃµ½µÄ±ä»»ÐòÁÐÆµÆ×');

e=sum((X1-X2).^2);

end

对理想采样信号A=444.128，α=50*2^(1/2)*π，Ω=50*2^(1/2)*π，T=1/1000，序列长度50，用自己编写的FFT算法和系统自带算法做64点FFT变换后绘制频域序列，如下：

对高斯序列，p=8，q=8，序列长度16，用自己编写的FFT算法和系统自带算法做16点FFT变换后绘制频域序列，如下：

对衰减正弦序列α=0.01，f=0.05，序列长度100，用自己编写的FFT算法和系统自带算法做128点FFT变换后绘制频域序列，如下：

由以上结果可知，自己编写的算法运行结果与系统自带算法一致，且可以对信号进行截断或补零后再做变换。

3.对所编制FFT算法进行性能评估：

与自己编写的DFT算法进行性能比较：

对N点序列进行DFT变换需要N²/2次复乘，而对N点序列做基2-FFT只需N/2*log2（N）次复乘，因此运算量减少了很多，且随着序列长度增加，运算量差异变大。

与系统自带FFT算法进行性能比较：

由于系统自带FFT函数用C语言实现，无法查看源代码，只知道效率更高，而且在计算任意点的DFT（不指定变换长度N）时，系统自带函数无需采取补零操作，而自己编写的函数会先补零再变换，改变了频域取样密度，会得到与系统自带函数不同的结果。

比较自己编写的DFT算法、自己编写的FFT算法和系统自带算法三者运算时间，得到下表：

由此绘制曲线如图：

由比较结果可知，虽然运算时间曲线和理想曲线不完全相同，但三种算法的相对运算时间与理论预期一致。

实验总结及个人结论：

从对实验的比较结果中可知，自己编写的FFT算法有效且比自己编写的DFT算法快很多，但却始终比系统自带的算法慢，一方面是因为系统自带算法实现效率高，另一方面，在观察了自己编写的算法中各步执行时间后，发现函数用在生成反序序列的时间和实际进行FFT运算的时间相当，相当于多用了一倍的时间。另外，对任意长序列，自己编写的FFT算法只能补零后计算，而不能像系统自带函数一样算出实际的N点DFT。

相关推荐

实验三-FFT算法的应用实验报告

第二篇：FFT算法实现实验报告

FFT算法实现实验报告

专栏推荐