下载文档原格式
快速傅里叶变换发展史快速傅立叶变换(FFT)个人日记2010-04-16 12:24:48 阅读163 评论0 字号：大中小订阅近十多年来数字信号处理技术同数字计算机、大规模集成电路等先进技术一样，有了突飞猛进的发展，日新月异，已经形成了一门具有强大生命力的技术科学。
由于它本身具有一系列的优点，所以能有效地促进各工程技术领域的技术改造和学科发展，应用领域也更加广泛、深入，越来越受到人们的重视。
在数字信号处理中，离散傅里叶变换(Discrete Fourier Transform,DFT)是常用的变换方法，它在各种数字信号处理系统中扮演着重要的角色。
傅里叶变换已有一百多年的历史了，我们知道频域分析常常比时域分析更优越，不仅简单，且易于分析复杂信号。
但用较精确的数字方法，即DFT进行谱分析，在FFT出现以前是不切实际的。
这是因为DFT计算量太大。
直到1965年出现了DFT]运算的一种快速方法以后，情况才发生了根本的变化。
快速傅里叶变换〔Fast Fourier Transfonn,FFT〕并不是与离散傅里叶变换不同的另一种变换，而是为了减少DFT计算次数的一种快速有效的算法。
当时Garwin在自己的研究中极需要一个计算傅立叶变换的快速方法，而正在写有关傅里叶变换的文章，Tukey概括地对Garwin介绍了一种方法，它实质上就是后来著名的Cooley-Tukey算法。
在Garwin的迫切要求下，1963年，IBM公司的Cooley根据Tukey的想法编写了第一个FFT算法程序。
在FFT算法中，Tukey主要利用了旋转因子的周期性和对称性。
这两个性质使DFT运算中的某些项可以合并，使DFT运算尽量分解为更少点数的DFT运算。
因为DFT的运算量与Pow(N,2)成比例，所以如果将一个大点数的DFT分解为若干个小点数的DFT 的组合，将有效地减少运算量。
Cooley在计算机上实现该算法时，为节省存储空间和减少寻址时间，采用了3维标号映射方法和在算法内部的循环结构，这些结构和技巧对后来的FFT算法研究及实现同样产生了很大影响。
快速傅⾥叶变换快速傅⾥叶变换快速傅⾥叶变换（FFT ）是根据计算量的最⼩化原理来设计和实施离散傅⾥叶变换(DFT)计算的⽅法。
1965年，库利（T.W.Cooley ）和图基（J.W.tukey ）发表了著名的《计算机计算傅⾥叶级数的⼀种算法》论⽂。
从此掀起了快速傅⾥叶变换计算⽅法研究的热潮。
快速傅⾥叶变换（FFT ）的出现，实现了快速、⾼效的信号分析和信号处理，为离散傅⾥叶变换(DFT)的⼴泛应⽤奠定了基础。
1.1离散傅⾥叶变换(DFT)的计算设x(n)是⼀个长度为M 的有限长序列，则定义x(n)的N 点离散傅⾥叶变换为∑-===10)()]([)(N n kn NW n x n x DFT k X 其中由于计算⼀个X(k)值需要N 次复乘法和(N-1)次复数加法，因⽽计算N 个X(k)值，共需N2次复乘法和N(N-1)次复加法。
每次复乘法包括4次实数乘法和2次实数加法，每次复加法包括2次实数加法，因此计算N 点的DFT 共需要4N2次实数乘法和(2N2+2N ·(N-1))次实数加法。
当N 很⼤时，这是⼀个⾮常⼤的计算量。
1.2减少DFT 计算量的⽅法减少DFT 的计算量的主要途径是利⽤k N W 的性质和计算表达式的组合使⽤，其本质是减少DFT 计算的点数N 以便减少DFT 的计算量。
k N W 的性质：（1）对称性： (2)周期性： (3) 可约性： (4) 特殊点：选择其中⼀个证明N N j k N j N k N j N k N e e e W 222)2(22πππ--+-+==ππj k N j e e --=2k N j e π2--=k N W -=FFT 算法是基于可以将⼀个长度为N 的序列的离散傅⾥叶变换逐次分解为较短的离散傅⾥叶变换来计算这⼀基本原理的。
这⼀原理产⽣了许多不同的算法，但它们在计算速度上均取得了⼤致相当的改善。
0,1,,1k N =-()*nk nk N N W W -=()()nk N n k n N k N N NW W W ++==nk mnk N mN W W =//nk nk m N N mW W =01N W =/21N N W =-(/2)k N k N NW W +=-在这⾥讨论两类基本的FFT 算法。
[实验2] 快速傅里叶变换 (FFT) 实现一、实验目的1、掌握FFT 算法和卷积运算的基本原理；2、掌握用C 语言编写DSP 程序的方法；3、了解利用FFT 算法在数字信号处理中的应用。
二、实验设备 1. 一台装有CCS 软件的计算机； 2. DSP 实验箱的TMS320C5410主控板； 3. DSP 硬件仿真器。
三、实验原理 （一）快速傅里叶变换傅里叶变换是一种将信号从时域变换到频域的变换形式，是信号处理的重要分析工具。
离散傅里叶变换（DFT ）是傅里叶变换在离散系统中的表示形式。
但是DFT 的计算量非常大， FFT 就是DFT 的一种快速算法, FFT 将DFT 的N 2 步运算减少至 ( N/2 )log 2N 步。
离散信号x(n)的傅里叶变换可以表示为∑=-=10][)(N N nk N W n x k X ， Nj N e W /2π-=式中的W N 称为蝶形因子，利用它的对称性和周期性可以减少运算量。
一般而言，FFT 算法分为时间抽取（DIT ）和频率抽取（DIF ）两大类。
两者的区别是蝶形因子出现的位置不同，前者中蝶形因子出现在输入端，后者中出现在输出端。
本实验以时间抽取方法为例。
时间抽取FFT 是将N 点输入序列x(n) 按照偶数项和奇数项分解为偶序列和奇序列。
偶序列为：x(0), x(2), x(4),…, x(N-2)；奇序列为：x(1), x(3), x(5),…, x(N-1)。
这样x(n) 的N 点DFT 可写成：()()∑++∑=-=+-=12/0)12(12/02122)(N n kn NN n nkNW n x Wn x k X考虑到W N 的性质，即2/)2//(22/)2(2][N N j N j N W e e W ===--ππ因此有：()()∑++∑=-=-=12/02/12/02/122)(N n nkN k NN n nkN W n x WWn x k X或者写成：()()12()kN X k X k W X k =+由于X 1(k) 与X 2(k) 的周期为N/2，并且利用W N 的对称性和周期性，即：k N N k N W W -=+2/可得：()()12(/2)kN X k N X k W X k +=-对X 1(k) 与X 2(k)继续以同样的方式分解下去，就可以使一个N 点的DFT 最终用一组2点的DFT 来计算。
4.2 课后习题详解4-1 如果一台通用计算机的速度为平均每次复乘40ns ，每次复加5ns ，用它来计算512点的DFT[x （n ）]，问直接计算需要多少时问?用FFT 运算需要多少时间?若做128点快速卷积运算，问最低抽样频率应是多少?解：①直接利用DFT 计算：复乘次数为N 2，复加次数为N （N-1）。
②利用FFT计算：复乘次数为，复加次数为N㏒2N 。
（1）直接计算复乘所需时间复加所需时间所以（2）用FFT 计算复乘所需时间复加所需时间所以4-2 N ＝16时，画出基-2按频率抽选法的FFT 流图采用输入自然顺序，输出倒位序），统计所需乘法次数（乘±1，乘±j 都不计在内）。
根据任一种流图确定序列x （n ）＝4cos （n π／2）（0≤n ≤15）的DFT 。
解：按频率抽取法的FFT 流图中的复数乘法出现在减法之后，其运算量为复数乘法：；复数加法：；由于N ＝16，有，，，不需要乘法。
按频率抽取，见图4-1（a ）。
图4-1（a ）运算量：复数乘法：由于，，，不需要乘法。
由图P4.2（a ）可知，共有的个数为1＋2＋4＋8＝15有的个数为1＋2＋4＝7所以总的乘法次数为32-15-7＝10（个）复数加法：举例：对序列x （n ）＝4cos （n π／2）（0≤n ≤15）可表示为由于N ＝16，可采用P4.2（b ）的流图。
设Xi （k ）＝（i ＝1，2，3，4）分别为第i 级蝶形结构的输出序列，则由P4.2（b ）的流图可知由于采用的是顺序输入、逆序输出的结构，因此输出X （k ）与X 4（k ）为逆序关系，即，为k 二进制逆序值由此可知，x （n ）的DFT 为X （4）＝X 4（2）＝32，X （12）＝X 4（3）＝12图4-1（b ）4-3 用MATLAB 或C 语言编制以下几个子程序。
（1）蝶形结运算子程序；（2）求二进制倒位序子程序；（3）基-2 DIT FFT 流程图，即迭代次数计算子程序。
快速傅立叶变换(FFT)算法实验摘要：FFT（Fast Fourier Transformation），即为快速傅里叶变换，是离散傅里叶变换的快速算法，它是根据离散傅里叶变换的奇、偶、虚、实等特性，对离散傅立叶变换的算法进行改进获得的。
这种算法大大减少了变换中的运算量，使得其在数字信号处理中有了广泛的运用。
本实验主要要求掌握在CCS环境下用窗函数法设计FFT快速傅里叶的原理和方法；并且熟悉FFT快速傅里叶特性；以及通过本次试验了解各种窗函数对快速傅里叶特性的影响等。
引言：快速傅里叶变换FFT是离散傅里叶变换DFT的一种快速算法。
起初DFT的计算在数字信号处理中就非常有用，但由于计算量太大，即使采用计算机也很难对问题进行实时处理，所以并没有得到真正的运用。
1965年J.W.库利和T.W.图基提出快速傅里叶变换，采用这种算法能使计算机计算离散傅里叶变换所需要的乘法次数大为减少,特别是被变换的抽样点数N越多,FFT算法计算量的节省就越显著。
从此，对快速傅里叶变换（FFT）算法的研究便不断深入，数字信号处理这门新兴学科也随FFT的出现和发展而迅速发展。
根据对序列分解与选取方法的不同而产生了FFT的多种算法，基本算法是基2DIT和基2DIF。
FFT 的出现，使信号分析从时域分析向频域分析成为可能，极大地推动了信号分析在各领域的实际应用。
FFT在离散傅里叶反变换、线性卷积和线性相关等方面也有重要应用。
一、 实验原理：FFT 并不是一种新的变换，它是离散傅立叶变换（DFT ）的一种快速算法。
由于我们在计算DFT 时一次复数乘法需用四次实数乘法和二次实数加法；一次复数加法则需二次实数加法。
每运算一个X （k ）需要4N 次复数乘法及2N+2（N-1）=2（2N-1）次实数加法。
所以整个DFT 运算总共需要4N^2次实数乘法和N*2(2N-1)=2N(2N-1)次实数加法。
如此一来，计算时乘法次数和加法次数都是和N^2成正比的，当N 很大时，运算量是可观的，因而需要改进对DFT 的算法减少运算速度。
function y=MyFFT_FB(x,n)%MYFFT_TB:My Fast Fourier Transform Frequency Based%按频率抽取基2-fft算法%input:% x -- 输入的一维样本% n -- 变换长度，缺省时 n=length(x) 当n小于x数据长度时，x数据被截断到第n个数据% 当n大于时，x数据在尾部补0直到 x 含n个数据%output:% y -- 1*n的向量，快速傅里叶变换结果%variable define:% N -- 一维数据x的长度% xtem -- 临时储存x数据用% m,M -- 对N进行分解 N=2^m*M,M为不能被2整除的整数% two_m -- 2^m% adr -- 变址，1*N的向量% l -- 当前蝶形运算的级数% W -- 长为 N/2的向量，记录 W(0,N),W(1,N),...W(N/2-1,N)% d -- 蝶形运算两点间距离% t -- 第l级蝶形运算含有的奇偶数组的个数% mul -- 标量，乘数% ind1,ind2 -- 标量，下标% tem -- 标量，用于临时储存%参考文献：% 81c 输入参数个数检查msg=nargchk(1,2,nargin);error(msg);%% 输入数据截断或加0N=length(x);if nargin==2if N<n % 加0xtem=x;x=zeros(1,n);x(1:N)=xtem;N=n;else % 截断xtem=x;x=xtem(1:n);N=n;endend%% 对N进行分解 N=2^m*M[m,M]=factorize(N);two_m=N/M;%% 变换if m~=0%% 如果N可以被2整除adr=address(m,M,two_m);y=x; % 蝶形运算级数 l=m 时%% 计算W向量W=exp(-2*pi*i* ( 0:N/2-1 ) /N);%% 蝶形运算d=N/2;t=1;for l=1:m% 加for ii=0:t-1ind1=ii*2*d+1;ind2=ind1+d;for r=0:d-1tem=y(ind1)+y(ind2);y(ind2)=y(ind1)-y(ind2);y(ind1)=tem;ind1=ind1+1;ind2=ind2+1;endend% 乘for r=0:d-1mul=W(r*t+1);for ii=0:t-1y(ii*2*d+d+1+r) = y(ii*2*d+d+1+r)*mul;endendd=d/2;t=t*2;end%% 直接傅里叶变换if M~=1 % N 分解含有非2因数M时，对y中每M个数据做直接傅里叶变换 for ii=1:two_my((ii-1)*M+1 : ii*M ) = DDFT( y((ii-1)*M+1 : ii*M ) );end%% 变址输出y=y(adr+1);else%% 如果N 不能被2整除y=DDFT(x);endend%% 内嵌函数 ====================================================== function y=DDFT(x)%% 直接离散傅里叶变换%input:% x -- 样本数据，N维向量%output:% y -- N维向量%参考文献：% 结构动力学，克拉夫，P82% variable define% s -- sum，用于求和N=length(x);y=zeros(size(x));for n=1:Ns=0;for m=1:Ns=s+x(m)*exp( -i*2*pi*(m-1)*(n-1)/N );endy(n)=s;endendfunction [m,M]=factorize(N)%% 对N分解m=0;while trueif mod(N,2)==0m=m+1;N=N/2;elsebreak;endendendfunction adr=address(m,M,two_m)%% 变址% b -- 2^m * m 的矩阵，用来存储二进制数据% ds -- 数，公差adr=zeros(two_m,M);b=de2bi(0:two_m-1,m);%转换为2进制注：matlab中二进制[0 1 1]=6 b=b(:,end:-1:1);% 逆序adr(:,1)=bi2de(b);%2进制转换为10进制if M~=1ds=two_m;adr=adr(:,1)*ones(1,M);adr=adr+ds*ones(size(adr,1),1)*(0:M-1);adr=reshape(adr',1,[]);endend。
快速傅里叶变换FFT的C语言算法彻底研究LED音乐频谱显示的核心算法就是快速傅里叶变换，FFT的理解和编程还是比较难的，特地撰写此文分享一下研究成果。
一、彻底理解傅里叶变换快速傅里叶变换（Fast Fourier Transform)是离散傅里叶变换的一种快速算法，简称FFT，通过FFT可以将一个信号从时域变换到频域。
模拟信号经过A/D转换变为数字信号的过程称为采样。
为保证采样后信号的频谱形状不失真，采样频率必须大于信号中最高频率成分的2倍，这称之为采样定理。
假设采样频率为fs，采样点数为N，那么FFT结果就是一个N点的复数，每一个点就对应着一个频率点，某一点n(n 从1开始)表示的频率为：fn=(n-1)*fs/N。
举例说明：用1kHz的采样频率采样128点，则FFT结果的128个数据即对应的频率点分别是0，1k/128，2k/128，3k/128，…，127k/128 Hz。
这个频率点的幅值为：该点复数的模值除以N/2（n=1时是直流分量，其幅值是该点的模值除以N）。
二、傅里叶变换的C语言编程1、对于快速傅里叶变换FFT，第一个要解决的问题就是码位倒序。
假设一个N 点的输入序列，那么它的序号二进制数位数就是t=log2N.码位倒序要解决两个问题：①将t位二进制数倒序；②将倒序后的两个存储单元进行交换。
如果输入序列的自然顺序号i用二进制数表示，例如若最大序号为15，即用4位就可表示n3n2n1n0，则其倒序后j对应的二进制数就是n0n1n2n3，那么怎样才能实现倒序呢？利用C语言的移位功能！程序如下，我不多说，看不懂者智商一定在180以下！复数类型定义及其运算#define N 64 //64点#define log2N 6 //log2N=6/*复数类型*/typedef struct{float real;float img;}complex;complex xdata x[N]; //输入序列/*复数加法*/complex add(complex a,complex b){complex c;c.real=a.real+b.real;c.img=a.img+b.img;return c;}/*复数减法*/complex sub(complex a,complex b){complex c;c.real=a.real-b.real;c.img=a.img-b.img;return c;}/*复数乘法*/complex mul(complex a,complex b){complex c;c.real=a.real*b.real - a.img*b.img;c.img=a.real*b.img + a.img*b.real;return c;}/***码位倒序函数***/void Reverse(void){unsigned int i,j,k;unsigned int t;complex temp;//临时交换变量for(i=0;i<N;i++)//从第0个序号到第N-1个序号{k=i;//当前第i个序号j=0;//存储倒序后的序号，先初始化为0for(t=0;t<log2N;t++)//共移位t次，其中log2N是事先宏定义算好的{j<<=1;j|=(k&1);//j左移一位然后加上k的最低位k>>=1;//k右移一位，次低位变为最低位}if(j>i)//如果倒序后大于原序数，就将两个存储单元进行交换（判断j>i是为了防止重复交换）{temp=x[i];x[i]=x[j];x[j]=temp;}}}2、第二个要解决的问题就是蝶形运算①第1级（第1列）每个蝶形的两节点“距离”为1，第2级每个蝶形的两节点“距离”为2，第3级每个蝶形的两节点“距离”为4，第4级每个蝶形的两节点“距离”为8。
《数字信号处理》课程设计报告按频率抽取的DFT快速算法分析及MATLAB实现专业：通信工程班级：组次：姓名：学号：目录摘要 (1)关键字 (1)0 引言 (1)1 按频率抽取的DFT快速算法原理 (1)2 DIF-FFT的运算规律及编程思想 (2)2.1 原位计算 (2)2.2 序列的倒序 (2)2.3 旋转因子的变换规律 (2)2.4 蝶形运算规律 (4)2.5 编程思想及程序框图 (4)3 DIF-FFT算法运算量分析 (5)4 MATLAB程序实现 (5)5 结束语 (7)参考文献 (7)按频率抽取的DFT 快速算法分析及MATLAB 实现摘要：DFT 是数字信号分析与处理中的一种重要变换。
但直接计算DFT 的计算量与变换区间长度N 的平方成正比，计算量非常大。
DFT 的快速算法使运算效率提高了1~2个数量级，为数字信号处理技术应用于各种信号的实时处理创造了条件。
为了对FFT 有更加深入的了解，本文对DIF-FFT 的原理进行了分析，并给出MATLAB 程序实现的方法与步骤。
关键词：DFT;DIF-FFT;MATLAB;0 引言DFT 是数字信号分析与处理中的一种重要变换。
但直接计算DFT 的计算量与变换区间长度N 的平方成正比，计算量非常大。
DFT 的快速算法使运算效率提高了1~2个数量级，为数字信号处理技术应用于各种信号的实时处理创造了条件。
本文通过对按频率抽取的DFT 快速算法原理介绍与MATLAB 实现以期使我们对傅里叶快速算法有更全面的理解，为我们以后更复杂的快速算法学习打下基础。
1 按频率抽取的DFT 快速算法原理设序列x(n)的长度为M N 2=,将序列前后对半分开，得到两个子序列，如下：式中：k kN NW )1(2/-=将x(k)分解成偶数组与奇数组，当k 取偶数（k=2m,m=0,1,…，N/2-1）时：∑∑-=-=++=++=12/02/12/02)]2()([)]2()([)2(N n mnN N n mn N W N n x n x W N n x n x m x （1）当k 取奇数（k=2m+1,m=0,1,…,N/2-1）时，∑∑-=-=+⋅+-=+-=+12/02/12/0)12()]2()([)]2()([)12(N n nmN n N N n m n N W W N n x n x W N n x n x m x （2）nk N N n Nk N kN n NN n nkN N n N N n nkN N n nk NN n nk NW W N n x n x W N n x W n x W n x Wn x Wn x k X ∑∑∑∑∑∑-=⎪⎭⎫ ⎝⎛+-=-=-=-=-=⎥⎦⎤⎢⎣⎡⎪⎭⎫ ⎝⎛++=⎪⎭⎫ ⎝⎛++=+==1202/212012012101202)(2)()()()()(令：⎭⎬⎫++=+-=)2()()()]2()([)(12Nn x n x n x W N n x n x n x n N 其中， n=0,1,2,…，N/2-1 将)()(21n x n x 和分别代入（1）、（2）式，可得：⎪⎭⎪⎬⎫∑∑-=-===+12/02/112/02/2)()2()()12(N n mnN N n nmN W n x m X W n x m X （3）（3）式表明，X(k)按奇偶k 值分为两组，其偶数组是)(1n x 的N/2点DFT ，奇数组则是)(2n x 的N/2点DFT 。
)(1n x 、)(2n x 和x （n ）之间的关系可以用图1所示的蝶形运算流图符号表示。
图2表示N=8的DIF-FFT 运算流图。
图1 DTF-FFT 蝶形运算流图符号图2 DIF-FFT 的运算流图（N=8）－1x (n )x (n ＋N / 2)nN W x (n )＋x (n ＋N / 2)[x (n )－x (n ＋N / 2)]n N W －1－10NW 2N W x (0)x (1)x (2)x (3)－1－1x (4)x (5)x (6)x (7)0NW 1N W 2N W 3N W X (0)X (4)X (2)X (6)X (1)X (5)X (3)X (7)N W 2N W －1－1－1－1－1－1－1－10NW 0N W 0NW 0NW2 DIF-FFT 的运算规律及编程思想 2.1 原位计算M N 2=点的FFT 共进行M 级运算，每级由N/2个蝶形运算组成。
同一级中，每个蝶形的两个输入数据只对计算本蝶形有用，而且每个蝶形的输入、输出数据结点又同在一条水平线上，这就意味着计算完一个蝶形后，所得输出数据可立即存入原输入数据所占用的存储单元。
这样，经过M 级运算后，原来存放输入序列数据的N 个存储单元中便依次存放X （k ）的N 个值。
原位计算可节省大量内存，从而使设备陈本降低。
2.2 序列的倒序由图2可知，DIF-FFT 算法输入序列为自然序列，而输出为倒序排列。
因此M 级运算完后，要对输出数据进行倒序才能得到自然顺序的X(k)。
图3为顺序与倒序二进制对照图。
图3 顺序与倒序二进制对照图2.3 旋转因子的变换规律N 点的DFT 快速傅里叶运算流图中，每级都有N/2个蝶形。
每个蝶形都要乘以因子PN W ，称其为旋转因子，P 为旋转因子的指数。
但各级的旋转因子和循环方式都有所不同。
为了编写计算程序，下面列出旋转因子P N W 与运算级数的关系。
用L 表示从左到右的运算级数（L=1,2,…，M ），第L 级共有12-L 个不同的旋转因子。
顺序倒序十进制数I二进制数 二进制数 十进制数J0 000 000 0 1 001 100 4 2 010 010 2 3 011 110 6 4 100 001 1 5 101 101 5 6 110 011 3 71111117对N=M 2的一般情况，第L 级的旋转因子为：J PN L W W 2= J=0,1,2，…，12-L -1因为 M L M L M LN --⋅=⨯=2222所以 LM M L J NJ N P N W W W --⋅⋅==22 J=0,1,2,12-L -1L M J P-⋅=22.4 蝶形运算规律对M N 2=点FFT ，输入倒位序，输出自然序，设第L 级运算每个蝶形的两节点距离为B行，则第L 级运算：{)()()(*)]()([)(1111B J A J A J A W B J A J A B J A L L L PNL L L ++⇐+-⇐+----2.5 编程思想及程序框图观察图2可以归纳出一些对编程有用的运算规律：第L 级中，每个蝶形的两个输入 数据相距B=12-L 个点；每级有B 个不同的因子；同一旋转因子对应着间隔为L2点的LM -2个蝶形。
频率抽取法输入为自然顺序，输出为倒序。
图4为大致流程图。
图5 为DIF-FFT 运算程序框图6倒序程序框图图4 大致流程图计算x 的长度n,不到2的数幂，补0开始读入x(n)蝶形运算序列倒序输出 结束DIF-FFT运算程序框图倒序程序框图图 5 图63 DIF-FFT算法运算量分析按频率抽取的 FFT算法是将频域信号序列X(K)分解为奇偶两部分，但算法仍是由时域信号序列开始逐级运算，把 N点分成N/2点计算FFT，可以把直接计算离散傅里叶变换所需的2N次乘法缩减到次。
4 MATLAB程序实现clc;clear all;close all;%关闭程序，清屏x=[1,2,3,4,5,6,7];x1=x;%本程序对输入序列实现DIF-FFT基2算法，点数取大于等于长度的2的幂次m=nextpow2(length(x)); %求的x长度对应的2的最低幂次mN=2^m;if length(x)<Nx=[x,zeros(1,N-length(x))]; %若的长度不是2的幂，补0到2的整数幂endfor L=m:-1:1 %将DFT做m次基2分解，从左到右，对每次分解作DFT运算D=2^L;u=1; %旋转因子u初始化WN=exp(-1i*2*pi/D); %本次分解的基本DFT因子WN=exp(-i*2*pi/D)for j=1:D/2 %本次跨越间隔内的各次蝶形运算for k=j:D:N %本次蝶形运算的跨越间隔为Dkp=k+D/2; %确定蝶形运算的对应单元下标temp=x(k); %保存x(k)的值x(k)=x(k)+x(kp); %加法运算x(kp)=(temp-x(kp))*u; %乘法运算endu=u*WN; %修改旋转因子，多乘一个基本DFT因子WNendendnxd=bin2dec(fliplr(dec2bin([1:N]-1,m)))+1; %求1：2^m数列的倒序disp('自编程序结果：');y=x(nxd)disp('系统自带函数结果：');y1=fft(x1,N)自编程序结果：y =Columns 1 through 528.0000 -9.6569 + 4.0000i -4.0000 - 4.0000i 1.6569 - 4.0000i 4.0000 Columns 6 through 81.6569 + 4.0000i -4.0000 + 4.0000i -9.6569 - 4.0000i系统自带函数结果：y1 =Columns 1 through 528.0000 -9.6569 + 4.0000i -4.0000 - 4.0000i 1.6569 - 4.0000i 4.0000Columns 6 through 81.6569 + 4.0000i -4.0000 + 4.0000i -9.6569 - 4.0000i5 结束语通过这次课程设计，能够提高我独立思考，解决学习问题的能力，并且重新温习了DIT-FFT运算，自学了DIF-FFT运算，对以前学过的知识掌握得更加牢固，同时也增加了MATLAB编程的信心。
参考文献[1] 唐向宏,岳恒立,郑雪峰.MATLAB及在电子信息类课程中的应用(第2版)[M].北京:电子工业出版社, 2009.6[2] 高西全,丁玉美.数字信号处理(第三版)[M].西安：西安电子科技大学出版社，2008。
