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

（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）时，系统自带函数无需采取补零操作，而自己编写的函数会先补零再变换，改变了频域取样密度，会得到与系统自