第四章快速傅里叶变换 41引言 DFT的计算量与N的平方成正比 DFT是信号分析与处理中的一种重要变换。但因直接计算DFT的计算量与 变换区间长度N的平方成正比,当N较大时,计算量很大,所以在快速傅里叶 变换出现以前,直接用DFT算法进行谱分析和信号实时处理是不切实际的。 直到1965年发现了DFT的一种快速算法,情况才发生了根本变化。 快速傅立叶变换发展 令1965年,图基(JW.Tuky)和库力(TW. Coody)在 Math. Computation 杂志发表著名的“机器计算傅立叶级数的一种算法”一文。桑德图 基快速算法出现,经改进,形成一套高效运算方法,这就是快速傅里 叶变换。 ◇1984年,杜哈梅尔和霍尔曼提出分裂基快速算法,使运算效率进 步提高。 42基2FFT算法 本章主要讨论基2FFT算法及其编程思想 42.1直接计算DFT的特点及减少运算量的基本途径 1、按DFT定义计算DFT的计算量 根据前面知识,长度为N的有限长序列x(m)的DFT为 (k)=∑x(n)W,k=01…,N 考虑x(n)为复数序列的一般情况,对某一个k值,直接按(1)式计算X(k)值 需要N次复数乘法、(N-1)次复数加法。因此剧所有N个k值共需N2次复数
第四章 快速傅里叶变换 4.1 引言 DFT 的计算量与 N 的平方成正比 DFT 是信号分析与处理中的一种重要变换。但因直接计算 DFT 的计算量与 变换区间长度 N 的平方成正比,当 N 较大时,计算量很大,所以在快速傅里叶 变换出现以前,直接用 DFT 算法进行谱分析和信号实时处理是不切实际的。 直到 1965 年发现了 DFT 的一种快速算法,情况才发生了根本变化。 快速傅立叶变换发展 1965 年,图基(J. W. Tuky)和库力(T.W. Coody)在 Math.Computation 杂志发表著名的“机器计算傅立叶级数的一种算法”一文。桑德-图 基快速算法出现,经改进,形成一套高效运算方法,这就是快速傅里 叶变换。 1984 年,杜哈梅尔和霍尔曼提出分裂基快速算法,使运算效率进一 步提高。 4.2 基 2FFT 算法 ——本章主要讨论基 2FFT 算法及其编程思想 4.2.1 直接计算 DFT 的特点及减少运算量的基本途径 1、按 DFT 定义计算 DFT 的计算量 根据前面知识,长度为 N 的有限长序列 x n 的 DFT 为 1 0 , 0,1, , 1 N kn N n X k x n W k N (1) 考虑 x n 为复数序列的一般情况,对某一个 k 值,直接按(1)式计算 X k 值 需要 N 次复数乘法、(N-1)次复数加法。因此对所有 N 个 k 值,共需 2 N 次复数
乘法及NN1次复数加法运算。当N很大时,运算量相当可。 2、减少运算量的基本途径 因为(a+b)≥a2+b2,所以,把N点DFT分解为几个较短的DFT,可 是乘法次数大大减少。 旋转因子Wm具有明显的周期性和对称性,即 m+IN Wm+=eM =W Wm=WXN-m或者「WA-m=W 因此,可纵通过不断地把长序列的DFT分解成几个短序列的DFT,并利 W如的周期性和对称性来减少运算次数。 42.2时域抽取法基2FFT基本原理 FFT算法大致可以分为的域抽取法FFT(DIFT, Decimation-In- Time fft) 和颜域抽取法( DIF-FFT, Decimation- n-Frequency FFT)两大类。首先介绍 DIT-FFT算法。 1、序列分解 设序列x(n)的长度为N,且满足 N=2M,M为自然数 按n的奇偶把x(n)分解为两个点的子序列 x()=x(2).r=01-…N-1 x()2=x(2x+1,r=0.1…,1 则x(m)的DFT为
乘法及 N(N-1)次复数加法运算。当 N 很大时,运算量相当可观。 2、减少运算量的基本途径 因为 2 2 2 a b a b ,所以,把 N 点 DFT 分解为几个较短的 DFT,可 是乘法次数大大减少。 旋转因子 m WN 具有明显的周期性和对称性,即 2 2 j m lN j m m lN m N N W e e W N N (2) * m N m N m m W W W W N N N N 或者 (3) 2 N m m W W N N 因此,可以通过不断地把长序列的 DFT 分解成几个短序列的 DFT,并利用 kn WN 的周期性和对称性来减少运算次数。 4.2.2 时域抽取法基 2FFT 基本原理 FFT 算法大致可以分为时域抽取法 FFT(DIT-FFT, Decimation-In-Time FFT) 和频域抽取法(DIF-FFT, Decimation-In-Frequency FFT)两大类。首先介绍 DIT-FFT 算法。 1、序列分解 设序列 x n 的长度为 N,且满足 2 , M N M 为自然数 按 n 的奇偶把 x n 分解为两个 2 N 点的子序列 1 2 , 0,1, , 1 2 N x r x r r 2 2 1 , 0,1, , 1 2 N x r x r r 则 x n 的 DFT 为
(k)=∑x(m)+∑x(m) ∑x(2r)W+∑x(2r+1)x =∑x(r)形2+W∑x1()W 由于 所以 X(k)=X1(k)+X2(k),k=01,2…,N-1 其中x(k)和x(k)分别是x()和x:()的点DF,即 X()=∑x()w=DT[x(r) x2(k)=23()W*=DFT[x2(r)I ·由于x()和x()均以为周期,且=-形,所以X(k)又可表示为 X(k)=X1(k)+WX2(k),k=0, Xk+2=X{(4)-X2(k),k=0.1 这样,就将N点DFT分解为两个点的DFT和(7式以及(8)式的运算 2、蝶形运算 为了形象直观地描述(7)和(8)式的运算,引入以下符号
1 1 2 2 2 2 1 0 0 1 1 2 2 2 2 1 2 0 0 2 2 1 kn kn N N n n N N kr kr N N r r N N kr k kr N N N r r X k x n W x n W x r W x r W x r W W x r W 偶数 奇数 由于 2 2 2 2 2 2 j kr N j kr kr kr N W e e W N N 所以 1 2 , 0,1,2, , 1 k X k X k W X k k N N (4) 其中 X k 1 和 X k 2 分别是 x r 1 和 x r 2 的 2 N 点 DFT,即 1 2 1 1 1 0 2 N kr N r X k x r W DFT x r (5) 1 2 2 2 2 0 2 N kr N r X k x r W DFT x r (6) 由于 X k 1 和 X k 2 均以 2 N 为周期,且 2 N k k W W N N ,所以 X k 又可表示为 1 2 , 0,1, , 1 2 k N N X k X k W X k k (7) 1 2 , 0,1, , 1 2 2 k N N N X k X k W X k k (8) 这样,就将 N 点 DFT 分解为两个 2 N 点的 DFT 和(7)式以及(8)式的运算。 2、蝶形运算 为了形象直观地描述(7)和(8)式的运算,引入以下符号 C A B A+ BC A-BC
图421蝶形运算符号 采用该符号上述分解算法可表示为 x(1) M2点x(1) (3) 图422N点DFT的一次时域抽取分解图(N=8) 3、分解运算对提高运算效率的影响 由图42.1可见,要完成一次蝶形运算,需要一次复数乘和两次复数加运算。 由图422可见,经过一次分解后,计算一个N点DFT共需要计算两个一点DFT 和个蝶形运算。而计算一个点DFT需要 2丿次复数乘和(N 次复 数加运算。所以,按上述方法,N点DFT共需2 ≈(N>1时)次 复数乘法和N(2-1)+2=2次复数加法运算。由此可见,仅经过次分 就使运算量减少近一半 4、进一步处理 进一步,如果N=2m,还可以再分解。将x()按奇偶分解成两个长的子 序列x()和x(),即 -X N l=0,1, x x 那么,X(k)又可以表示为
图 4.2.1 蝶形运算符号 采用该符号上述分解算法可表示为 N/2点 DFT W N 0 N/2点 DFT W N 1 W N 2 W N 3 x(0) X1 (0) x(2) x(4) x(6) x(1) x(3) x(5) x(7) X1 (1) X1 (2) X1 (3) X2 (0) X2 (1) X2 (2) X2 (3) X(0) X(1) X(2) X(3) X(4) X(5) X(6) X(7) 图 4.2.2 N 点 DFT 的一次时域抽取分解图(N=8) 3、分解运算对提高运算效率的影响 由图 4.2.1 可见,要完成一次蝶形运算,需要一次复数乘和两次复数加运算。 由图 4.2.2 可见,经过一次分解后,计算一个 N 点 DFT 共需要计算两个 2 N 点 DFT 和 2 N 个蝶形运算。而计算一个 2 N 点 DFT 需要 2 2 N 次复数乘和 1 2 2 N N 次复 数加运算。所以,按上述方法,N 点 DFT 共需 2 2 2 2 2 2 N N N (N>>1 时)次 复数乘法和 2 1 2 2 2 2 N N N N 次复数加法运算。由此可见,仅经过一次分解, 就使运算量减少近一半。 4、进一步处理 进一步,如果 2 m N ,还可以再分解。将 x r 1 按奇偶分解成两个 4 N 长的子 序列 x l 3 和 x l 4 ,即 3 1 4 1 2 , 0,1, , 1 2 1 4 x l x l N l x l x l 那么, X k 1 又可以表示为
()=x(2)形2+x(2+)形 ∑x:()形4+W∑x() 2l= (k)+WX(k),0,1 N 式中 X()=∑x()=DF[x()] X()=∑x(H=DFT[x:(O] 同理,由X(k)和x4(k)的周期性和W的对称性WN4=-W,最后得到 X1(k)=X2(k)+WX4(k) N 0,1, 4 x(4)-2x() 用同样的方法可得 0.1 N Xs (k)-WNX(k) 其中 X(k)=∑x(0)=DF X()=∑x()w=DFT[x()] x5()=x2(2) x()=x1(21+1) 这样,经过第二次分解,又将N2点DFT分解为两个N4点DFT式。依此类推
1 1 4 4 2 2 1 1 1 / 2 1 / 2 0 0 1 1 4 4 3 / 4 4 / 4 0 0 2 3 4 2 2 2 1 , 0,1, , 1 2 N N M k l N N l l N N kl k kl N N N l l k N X k x l W x l W x l W W x l W N X k W X k (9) 式中 1 4 3 3 3 0 4 N kl N l X k x l W DFT x l 1 4 4 4 4 0 4 N kl N l X k x l W DFT x l 同理,由 X k 3 和 X k 4 的周期性和 2 m WN 的对称性 4 2 2 N k k W W N N ,最后得到 1 3 4 2 1 3 4 2 , 0,1, , 1 4 4 k N k N X k X k W X k N k N X k X k W X k (10) 用同样的方法可得 2 5 6 2 2 5 6 2 , 0,1, , 1 4 4 k N k N X k X k W X k N k N X k X k W X k (11) 其中 1 4 5 5 5 0 4 N kl N l X k x l W DFT x l 1 4 6 6 6 0 4 N kl N l X k x l W DFT x l 5 2 6 2 2 , 0,1, , 1 2 1 4 x l x l N l x l x l 这样,经过第二次分解,又将 N/2 点 DFT 分解为两个 N/4 点 DFT 式。依此类推