第03章离散傅里叶变换及其快 速算法 邹江 zhujiang@public.wh.hb.cn
第03章 离散傅里叶变换及其快 速算法 邹江 zoujiang@public.wh.hb.cn
3.5快速傅里叶变换(FFT) 351DFT的计算量 离散傅里叶变换在实际应用中是非常重要的,利用它可以计算信 号的频谱、功率谱和线性卷积等。但是,如果使用定义式(3.20)来 直接计算DFT,当N很大时,即使使用高速计算机,所花的时间也 太多。因此,如何提高计算DFT的速度,便成了重要的研究课题 1965年库利( Cooley和图基( Tukey)在前人的研究成果的基础上提出了 快速计算DFT的算法,之后,又出现了各种各样快速计算DFT的方 法,这些方法统称为快速傅里叶变换( Fast Fourier Transform),简称为 FFT. FFT的出现,使计算DFT的计算量减少了两个数量级,从而成 为数字信号处理强有力的工具。 快速傅里叶变换(F「T)是离散傅里叶变换(OFT)的快速算法。它是 DSP领域中的一项重大突破,它考虑了计算机和数字硬件实现的约 束条件、研究了有利于机器操作的运算结构,使DFT的计算时间缩 短了1~2个数量级,还有效地减少了计算所需的存储容量,FT技术 的应用极大地推动了DSP的理论和技术的发展
3. 5 快速傅里叶变换(FFT) 3.5.1 DFT的计算量 离散傅里叶变换在实际应用中是非常重要的，利用它可以计算信 号的频谱、功率谱和线性卷积等。但是，如果使用定义式(3.20)来 直接计算DFT，当N很大时，即使使用高速计算机，所花的时间也 太多。因此，如何提高计算DFT的速度，便成了重要的研究课题。 1965年库利 (Cooley)和图基(Tukey)在前人的研究成果的基础上提出了 快速计算DFT的算法，之后，又出现了各种各样快速计算DFT的方 法，这些方法统称为快速傅里叶变换(Fast Fourier Transform)，简称为 FFT。FFT的出现，使计算DFT的计算量减少了两个数量级，从而成 为数字信号处理强有力的工具。 快速傅里叶变换(FFT)是离散傅里叶变换(DFT)的快速算法。它是 DSP领域中的一项重大突破，它考虑了计算机和数字硬件实现的约 束条件、研究了有利于机器操作的运算结构，使DFT的计算时间缩 短了1~2个数量级，还有效地减少了计算所需的存储容量，FFT技术 的应用极大地推动了DSP的理论和技术的发展
在导出FFT算法之前,首先来估计一下直接计算DFT所需的计算量 DFT的定义 ∑x(n)W,0≤k≤N-1 X(R)=DFT[x(n)]=3 m= 0, 其它 「∑x()W*,0≤n≤N-1 x(n)= IDFT[X(k)] 其 其中 e =cos(2π/N)-jsin(2π/N
DFT的定义 在导出FFT算法之前，首先来估计一下直接计算DFT所需的计算量。 其中
将DFT定义式展开成方程组 X(0)=x(0)w°+x(1)Wg1+…+x(N-1)W-1 X(1)=x(0)Wk°+x(1)Wk+…+x(N-1)Wk-1 X(2)=x(0)W3+x(1)W3+…+x(N-1)W3w-1) ⅹ(N-1)=x(0)WN-10+x(1)W-11+…+x(N-1)W-1) 将方程组写成矩阵形式 X(0) WWl W(N-1) 「x(0) X(1) W 1·0 l·1 W 1·(N-1) x(1) X(N-1 W-1)0WA-1)1…W-1)( x(N-1) 用向量表示 X-W
将DFT定义式展开成方程组 将方程组写成矩阵形式 用向量表示
用复数表示 Wx={RewW]·Rex-Im[W]·Im[x]} j{Re[W]·Im[x]+Re[x]·Im[w] 从矩阵形式表示可以看出,由于计算一个X(k)值需要N次复乘法和 (N-1)次复数加法,因而计算N个X(k)值,共需№次复乘法和N(N-1)次 复加法。每次复乘法包括4次实数乘法和2次实数加法,每次复加 法包括2次实数加法,因此计算N点的DFT共需要4N2次实数乘法和 (2N2+2N(N-1)次实数加法。当N很大时,这是一个非常大的计算量。 N=8 FF算法主要利用了W的两个性质: 2x (1)对称性,即 N Re 2)周期性,即(m)*=∥n (N-nk k=0 N k=1 k=2 W=Wr为任意整数 图3.14W的特性
从矩阵形式表示可以看出，由于计算一个X(k)值需要N次复乘法和 (N-1)次复数加法，因而计算N个X(k)值，共需N 2次复乘法和N(N-1)次 复加法。每次复乘法包括4次实数乘法和2次实数加法，每次复加 法包括2次实数加法，因此计算N点的DFT共需要4N2次实数乘法和 (2N2+2N·(N-1))次实数加法。当N很大时，这是一个非常大的计算量。 FFT算法主要利用了WN k的两个性质： (1)对称性，即 (2)周期性，即 用复数表示： r为任意整数
FFT算法是基于可以将一个长度为N的序列的离散傅里叶变换 逐次分解为较短的离散傅里叶变换来计算这一基本原理的。这 原理产生了许多不同的算法,但它们在计算速度上均取得了 大致相当的改善。 在本章中我们集中讨论两类基本的FFT算法。 第一类称为按时间抽取( Decimation -in-Tme)的基2FFT算法, 的命名来自如下事实:在把原计算安排成较短变换的过程中 序列x(n)通常看作是一个时间序列)可逐次分解为较短的子序列 第二类称为按频率抽取( Decimation-in-Frequency)的基2FFT算法 在这类算法中是将离散傅里叶变换系数序列X(k)分解为较短的 子序列。 前面两种算法特别适用于N等于2的幂的情况。 对于N为合数的情况,本章也将介绍两种处理方法
FFT算法是基于可以将一个长度为N的序列的离散傅里叶变换 逐次分解为较短的离散傅里叶变换来计算这一基本原理的。这 一原理产生了许多不同的算法，但它们在计算速度上均取得了 大致相当的改善。 在本章中我们集中讨论两类基本的FFT算法。 第一类 称为按时间抽取(Decimation-in-Time)的基2FFT算法，它 的命名来自如下事实：在把原计算安排成较短变换的过程中， 序列x(n)(通常看作是一个时间序列)可逐次分解为较短的子序列。 第二类称为按频率抽取(Decimation-in-Frequency)的基2FFT算法， 在这类算法中是将离散傅里叶变换系数序列X(k)分解为较短的 子序列。 前面两种算法特别适用于N等于2的幂的情况。 对于N为合数的情况，本章也将介绍两种处理方法
3.5.2时间抽选基2F算法(库里一图基算法) 这种算法简称为时间抽选FFT算法,其基本出发点是,利用旋 转因子WN的对称性和周期性,将一个大的DFT分解成一些逐次 变小的DFT来计算。 分解过程遵循两条规则: ①对时间进行偶奇分解; ②对频率进行前后分解。 设N=2M,M为正整数。为了推导方便,取N=23=8,即离散时间 信号为 x(n)={x(0),x(1),x(2),x(3),x(4),x(5),x(6),x(7)} 按照规则(1),将序列x(n)分为奇偶两组,一组序号为偶数,另一 组序号为奇数,即 x(0),x(2),x(4),x(6)|x(1),x(3),x(5),x(7) 分别表示为: g(r)=x(2r)偶数项 N 0,1 h2(r)=x(2r+1)—奇数项 2
3. 5. 2 时间抽选基2FFT算法(库里—图基算法) 这种算法简称为时间抽选FFT算法，其基本出发点是，利用旋 转因子WN k的对称性和周期性，将一个大的DFT分解成一些逐次 变小的DFT来计算。 分解过程遵循两条规则： ①对时间进行偶奇分解； ②对频率进行前后分解。 设N＝2 M ，M为正整数。为了推导方便，取N＝2 3＝8，即离散时间 信号为 按照规则(1)，将序列x(n)分为奇偶两组，一组序号为偶数，另一 组序号为奇数，即 分别表示为：
根据DFT的定义 N x(k)=∑x(n)W=∑x(n)+>(m)W对 偶数n 奇数n N/2-1 N/2-1 ∑x(2W+∑x(2r+1)W (2r+1)k N/2-1 N/2-1 ∑g()W)*+W∑h()(W系) 〃=0 0 因为WM2=W2,所以上式变为 N/2-1 N/2-1 x(k)=∑g()W对n+W∑h()W对a G(k)+WH(),k=0,1,…N-1364) 上式中的G(k)和H(k)都是N2点的DFT
根据DFT的定义 因为 WN 2=WN/2 1，所以上式变为 上式中的G(k)和H(k)都是N/2点的DFT。 (3.64)
按照规则(2),将Xk)分成前后两组,即 X(k)={X(0),X(1),X(2),X(3)|x(4),X(5),X(6),X(7) 由(364)表示的是N2点DFT,前4个k值的DFT可表示为 X(k)=G(k)+WH(k),k=0,1,2,2~1365) 后4个k值的X(k)表示为: N Xk+ >g(r)w? (+2) W ∑h(r)W (艹+2) 2 N/2 N/2 r=0 因为/(+学) N w,Wz=1所以 N/2 N N/2-1 Xk+ ∑g()W2-W∑h()W N =G(k)-WMH(k),k=0,1,2, 3.66
按照规则(2)，将X(k)分成前后两组，即 由(3.64)表示的是N/2点DFT，前4个k值的DFT可表示为 后4个k值的X(k)表示为： 因为 所以 (3.66) (3.65)
按照式(365)和式(366)可画出图3.15所示的信号流程图。 c(0) X(O)=G(0)+WH(0) x(2) N c(1) 点 X(1)=G(1)+WH(1) G(2) DFT X(2)=c(2)+W(2 G(3) x(6) X(3)=C(3)+W(3) H(0) x X(4)=G(0)-WH(0) x(3) H(1) 点 X(5)=G(1)-WH(1) x(5) H(2) FT X(6)=G(2)-WH(2) H(3) oX(7)=G(3)-WAH(3) 图315N点的DFT计算分解成两个N/2点的DFT计算的时间抽选流程图(N=8)
按照式(3.65)和式(3.66)可画出图3．15所示的信号流程图