数字信号处理：第7章快速傅里叶变换FFT

1、第七章 快速傅里叶变换（FFT）第七章学习目标理解按时间抽选的基-2FFT算法的算法原理、运算流图、所需计算量和算法特点理解按频率抽选的基-2FFT算法的算法原理、运算流图、所需计算量和算法特点理解IFFT算法了解CZT算法理解线性卷积的FFT算法引言FFT: Fast Fourier Transform有限长序列通过离散傅里叶变换（DFT）将其频域离散化成有限长序列，但其计算量太大，很难实时处理，因此引出了快速傅里叶变换。1965年，Cooley(库利)-Turky(图基) 发表文章机器计算傅里叶级数的一种算法，提出FFT算法，解决DFT运算量太大，在实际使用中受限制的问题。FFT并不是一种

2、新的变换形式，它只是DFT的一种快速算法，并且根据对序列分解与选取方法的不同产生了多种算法。FFT的应用：频谱分析、滤波器实现、实时信号处理、离散傅里叶反变换、线性卷积和线性相关等方面。DSP芯片实现。TI公司的TMS 320c30，10MHz时钟，基2-FFT1024点FFT时间15ms。7.1 直接计算DFT的问题及改进途径1. 运算量复数乘法复数加法一个X(k)NN 1N个X(k)(N点DFT)N 2N (N 1)实数乘法实数加法一次复乘42一次复加2一个X (k)4N2N+2 (N 1)=2 (2N 1)N个X (k)(N点DFT)4N 22N (2N 1)2.4. FFT算法分类:时

3、间抽选法 DIT: Decimation-In-Time频率抽选法 DIF: Decimation-In-Frequency N点DFTN/2点DFTN/4点 DFT2点 DFT1个 2个 4个 N/2个3. FFT算法的基本思想： 利用DFT系数 特性，合并DFT运算中的某些项，将大点数的DFT分解成若干个小点数的DFT。7.2 按时间抽取的基-2FFT算法（ DIT ）1、算法原理设序列点数 N = 2M，M 为整数。 若不满足，则补零将序列x(n)按n的奇偶分成两组：N为2的整数幂的FFT算法称基-2FFT算法。则x(n)的DFT:两个k一样吗？只是X(k)前半部分利用周期性求X(k)的

4、后半部分:图7-1 按时间抽取蝶形运算流图N点DFT的第一次按时间抽取过程：第一步：x(n)按n的奇偶分解成偶序号序列和奇序号序列第二步：对分解后的序列进行N/2点DFT第三步：对两序列的DFT进行蝶形运算得到X(k)按n的奇偶分解即图7-2 一个N点DFT分解为两个 点DFT（N=8）分解后的运算量：复数乘法复数加法N点DFT一个N / 2点DFTN 2(N / 2)2N (N 1)N / 2 *(N / 2 1)两个N / 2点DFTN 2 / 2N (N / 2 1)一个蝶形12N / 2个蝶形N / 2N总计运算量减少了近一半N / 2仍为偶数，进一步分解：N / 2 N / 4图7-

5、3 一个 点DFT分解为两个 点DFT（N=8）同理：其中：图7-4 一个N点DFT分解为四个 点DFT（N=8）图7-2 一个N点DFT分解为两个 点DFT（N=8）这样逐级分解，直到两点DFT，两点DFT不需要乘法运算。当N = 8时，即分解到X3(k)，X4(k)，X5(k)，X6(k)，k = 0, 1两点DFT可用一个蝶形运算实现例如：2点DFT的运算流图 当N = 2M时，N点DFT可完全分解为M 级蝶形运算，每级有N/2个蝶形结例如：N=8 = 23 M=3 图7-5 按时间抽取FFT信号流图（N=8）第一级第二级第三级2、运算量 当N = 2M时，共有M级蝶形， 每级有N/2个

6、如下的蝶形单元：总共复数乘法：复数加法：一个蝶形单元只需一次复数乘法，两次复数加法。可以共享直接计算DFT与FFT算法的计算量比较FFT复数乘法：DFT复数乘法：直接计算DFT与FFT算法所需乘法次数的比较图3、算法特点（ ）1）蝶形运算图7-6 按时间抽取蝶形运算结m表示第m级迭代，i和j表示数据所在的行数，其中：最后一级第M级共有N/2个系数，即 ，前推一级系数减少一半，只用后一级的偶序号上的系数，即每个蝶形的两节点距离2）同址运算(原位运算) 两点构成一个蝶形单元，并且这两点不再参与别的蝶形单元的运算。 当数据输入到存储器以后，每一级运算的结果仍然存储在这同一组存储器中，直到最后输出，中

7、间无需其他存储器。第一级第二级第三级蝶形运算单元 由于每一步分解都是按输入序列在时间上的奇偶次序，分解成两个半长的子序列，所以称“按时间抽取算法”。存储器3）码位倒序 输入序列x(n)的排列次序不符合自然顺序。此现象是由于按n的奇偶分组进行DFT运算而造成的，这种排列方式为“码位倒序”。 所谓“码位倒序”，是指按二进制表示的数字首尾位置颠倒，重新按十进制读数。表 顺序和倒序二进制数对照表例：N =8=234）存储单元输入序列x(n) : N个存储单元系数 ：N / 2个存储单元4.小结 全部计算分解为M级，或称为M次迭代。 输入倒序，输出正序。 每级都包含N/2个蝶形单元。 每级的若干蝶形单元

8、组成“群”。第1级群数为N/2，第2级群数为N/4，第i级群数为N/2i，最后一级的群数为1。 每个蝶形单元都包含乘WNk与-WNk的运算（可简化为乘WNk与加、减法各一次）。 同一级中，各个群的W分布规律完全相同。7.3 按频率抽选的基-2FFT算法（ DIF ）1、算法原理 将X(k)按k的奇偶分组前，先将输入x(n) 按n的顺序分成前后两半：设序列点数 N = 2M，M 为整数。 若不满足，则补零。前一半x(n)：n=0,1,2,N/2-1后一半x(n)：n=N/2,N/2+1,N-1 按k的奇偶将X(k)分成两部分：当k=2r时：当k=2r+1时：其中： 则图7-8 按频率抽取蝶形运算

9、流图N点DFT的第一次按频率抽取过程：第一步：x(n)按n的顺序分成前一半序列和后一半序列第三步：对两序列进行N/2点DFT得到X(k)第二步：前一半和后一半序列进行蝶形运算得即例如：N=8 = 23 M=3r=0,1,2,3图7-9 一个N点DFT分解为两个 点DFT（N=8）分解后的运算量：N点DFT复数乘法N 2复数加法N (N 1)一个N / 2点DFT(N / 2)2N / 2 *(N / 2 1)两个N / 2点DFTN 2 / 2N (N / 2 1)一个蝶形12N / 2个蝶形N / 2N总计运算量减少了近一半 按r的奇偶将X(2r)分成两部分：当r=2l时：当r=2l+1时：

10、 按r的奇偶将X(2r)分成两部分：N /2仍为偶数，进一步分解：N /2 N /4图7-10 一个 点DFT分解为两个 点DFT（N=8）图7-10 一个 点DFT分解为两个 点DFT（N=8）图7-10 一个 点DFT分解为两个 点DFT（N=8）同理：图7-11 一个8点DFT分解成四个两点DFTN点DFT经过两次按频率抽取后：逐级分解，直到两点DFT，两点DFT不需要乘法运算。当N=8时，分解到x3(n)，x4(n)，x5(n)，x6(n)，n=0,1结论：两点DFT可用一个蝶形运算来实现蝶形运算蝶形运算两点DFT可用一个蝶形运算来实现例如：2点DFT的运算流图 当N = 2M时，N点

11、DFT可完全分解为M 级蝶形运算，每级有N/2个蝶形结例如：N=8 = 23 M=3 图7-12 按频率抽取FFT信号流图（N=8）2、运算量当N = 2M时，共有M级蝶形，每级N/2个蝶形，每个蝶形有1次复数乘法2次复数加法。比较DFT： 复数乘法：复数加法：图 直接计算DFT与FFT算法所需乘法次数的比较3、算法特点（ ）1）蝶形运算图7-13 按频率抽取蝶形运算结m表示第m级迭代，i和j表示数据所在的行数，其中：最后一级第M级共有N/2个系数，即 ，前推一级系数减少一半，只用后一级的偶序号上的系数，即每个蝶形的两节点距离2）同址运算(原位运算) 当数据输入到存储器以后，每一级运算的结果仍

12、然存储在这同一组存储器中，直到最后输出，中间无需其他存储器。3）码位倒序X(k)按k的码位倒序输出，x(n)按n的自然顺序输入4、DIT与DIF的异同基本蝶形不同DIT: 先复乘后加减DIF: 先减后复乘运算量相同都可原位运算DIT和DIF的基本蝶形互为转置复乘数复加数输入输出顺序不同DIT: 输入码位倒序，输出自然顺序DIF: 输出码位倒序，输如自然顺序DITDIF 将流图的所有支路方向都反向，并且交换输入与输出，但节点变量值不交换，这样即可从一种流图得到另一种。 这样，对每一种按时间抽取的FFT流图都存在一个按频率抽取的FFT流图，反之亦然。 频率抽取法与时间抽取法实质上是两种等价的FFT

13、运算。转置定理直接DFT方法 / CZT方法：当要求准确的N点DFT，且N是素数时7.4 N为复合数的FFT算法混合基算法基-2FFT算法：补零使满足混合基FFT算法：N是复合数7.5 快速傅里叶反变换（IFFT）算法比较：IDFT:DFT:改变FFT的IFFT程序 图7-15 按时间抽取IFFT信号流图（N=8）共轭FFT共轭乘1/ N直接调用FFT子程序计算IFFT的方法：直接用FFT程序的IFFT程序7.6 线性调频 z变换（CZT)算法FFT不适用于：只研究信号的某一频段，要求对该频段抽样密集，提高分辨率；研究非单位圆上的抽样值；需要准确计算N点DFT，且N为大的素数；等等。CZT算法

14、：对z变换采用螺线抽样，chirp-z变换线性调频 z变换1、算法原理 N点有限长序列，其z变换：沿z平面上的一段螺线作等分角抽样，抽样点zk：其中：M为要分析的复频谱点数 则图7-16 CZT计算路径 螺旋采样 抽样点：:起始抽样点的相角:相邻抽样点的角度差:分析路径的盘旋趋势(螺线的伸展率)W01:向内盘旋(内缩) W01:向外盘旋(外伸) :逆时针 :顺时针图7-16 CZT计算路径:起始抽样点的矢量半径长度谱分析的起始点位置求抽样点处的z变换：令令图7-17 CZT运算流图其中： 卷积可以通过圆周卷积来实现，这样可借用FFT快速算法。2、CZT的实现步骤(自学)图7-17 CZT运算流图取前M个点图7-18 CZT变换的圆周卷积1) 选择 ，且2) 形成L点序列f(n)：递归法求解求f(n) 的L点DFT：3）形成L点序列h(n):或求h(n) 的L点DFT：4）求乘积5）取g(n)的前M个点得g(k)6）求抽样点的z变换:求G(k) 的L点IDFT：说明令3、CZT算法的特点 可计算单位圆上任一段曲线上的Z变换，可任意给定起止频率； 周线不必是z平面上的圆，在语音分析中螺旋周线具有