<p align="left">快速傅里叶变化(FFT)是离散傅里叶变化(DFT)的快速算法。一个长度为N的有限长序列的DFT为</p>
<p align="center"></p>
<p align="left">式中的 为 。离散傅里叶反变换为</p>
<p align="center"></p>
<p align="left">FFT算法的基本思想是利用 的对称性和周期性,即 =- 和 (r为任意整数),通过将长序列分解为短序列,再由短序列的DFT组合来实现整个长序列的DFT,从而减少DFT的复数乘法和复数加法次数,提高DFT的运算速度和效率。将N点的长序列分解为短序列可以有多种方式,常用的有基2时间抽取FFT算法,基2频率抽取FFT算法,基4抽取FFT算法和分裂基FFT算法。基2抽取算法简单有效,基4抽取算法和分裂基抽取算法与基2抽取算法原理类似,但基4算法和分裂基算法效率更高,运算量更少。本文主要介绍基2抽取FFT算法。</p>
<p align="left"><strong> 1、基2时间抽取FFT算法:<sup></sup></strong></p>
<p align="left">基2时间抽取FFT算法是将N点序列x[n]按序号奇偶分解得到两个 点子序列 和 ,通过计算这两个子序列的DFT 和 实现 的DFT 。算法要点为N点序列x[n]对时域下标按奇数和偶数分解,即</p>
<p align="center">(1)<br> </p>
<p align="left">对频域X[k]按前一半和后一半分解,即 。频域的前一半可由下式求得</p>
<p align="center"></p>
<p align="left">由DFT的周期性和 的对称性有</p>
<p align="center"></p>
<p align="center"></p>
<p align="center">-=</p>
<p align="left">故频域的后一半可由下式求得</p>
<p align="left"></p>
<p align="left">将前一半和后一半组合起来,即得原来的长序列x[n]的DFT 。</p>
<p align="left">可对 和 继续按式(1)分解下去,直到每个子序列的长度为1为止。任何一个N= 点的DFT,可由通过M次分解,得到长度为1的子序列。对于N不是2的整数次幂的情况,可在序列后面添0补齐到2的整数次幂。每级分解运算可用一个蝶形运算流程图来表示,如图</p>
<p align="center"></p>
<p align="left">下图是对一个长度为8的序列进行基2时间抽取FFT运算的流程图</p>
<p align="center"></p>
<p align="left">直接计算N点序列的DFT需要N<sup>2</sup>次复数乘法,需要N(N-1)次复数加法。而用基2时间抽取FFT算法,则乘法次数为 ,加法次数 ,运算次数大大减少。</p>
<p align="left"><strong>2、基2频率抽取FFT算法</strong></p>
<p align="left">基2频率抽取FFT算法则是对频域下标按奇数和偶数分解,对时域下标按前一半和后一半分解。 的DFT 可表示为</p>
<p align="center"></p>
<p align="left">将 分解成 和 两个长度为 的子序列,如下:</p>
<p align="center">(2)</p>
<p align="left">则 的DFT 可按频率下标分奇数和偶数表示为</p>
<p align="center"></p>
<p align="center"></p>
<p align="left">可继续将 和 按式(2)分解下去,直到每个子序列长度为1。任何一个N= 点的DFT,可由通过M次分解,得到长度为1的子序列。对于N不是2的整数次幂的情况,可在序列后面添0补齐到2的整数次幂。每级分解运算也可用一个蝶形运算流程图来表示,如图</p>
<p align="left"></p>
<p align="left">下图是对一个长度为8的序列进行基2频率抽取FFT运算的流程图:</p>
<p align="left"></p>
<p align="left"><strong>3、编程实现</strong></p>
<p align="left">Matlab有一个内置的函数fft()可以实现快速傅里叶变换,由于它是built-in函数,所以无法查看源码来判断它是应用的那种分解方式来实现的快速傅里叶变换。下面将用Matlab语言来编写实现类似fft()的快速傅里叶变换函数myfft()。本函数采用的FFT算法为基2时间抽取FFT算法,函数定义源码如下:</p>
<p align="left">%----------------------------------</p>
<p align="left">%myfft(A,M)有两个形参,A为要进行快速傅里叶变换的序列,M为序列长度对2的对数,</p>
<p align="left">%即如果序列A[n]的长度为8,则M=3</p>
<p align="left">function [A] = myfft(A,M)</p>
<p align="left">N=2^M; %序列的长度为N</p>
<p align="left">LH=N/2; %序列长度的一半</p>
<p align="left">N1=N-2;</p>
<p align="left">J=LH</p>
<p align="left">for I=1:1:N1 %将输入序列A[n]按运算流程图第一列的顺序排好</p>
<p align="left"> if I<J</p>
<p align="left"> temp=A(I+1);</p>
<p align="left"> A(I+1)=A(J+1);</p>
<p align="left"> A(J+1)=temp;</p>
<p align="left"> end</p>
<p align="left"> K=LH;</p>
<p align="left"> while J>=K</p>
<p align="left"> J=J-K;</p>
<p align="left"> K=K/2;</p>
<p align="left"> end</p>
<p align="left"> J=J+K;</p>
<p align="left">end</p>
<p align="left">for L=1:1:M |
|