数字信号处理_第四章

第4章快速傅里叶变换(FFT)第4章快速傅里叶变换(FFT)4.1引言4.2基2FFT算法4.3进一步减少运算量的措施第4章快速傅里叶变换(FFT)§4.1引言FFT:FastFourierTransform1965年，Cooley-Turky发表文章《机器计算傅里叶级数的一种算法》，提出FFT算法，解决DFT运算量太大，在实际使用中受限制的问题。FFT的应用频谱分析、滤波器实现、实时信号处理等。DSP芯片实现TI公司的TMS320c30，10MHz时钟，基2-FFT1024点FFT时间15ms。第4章快速傅里叶变换(FFT)4.2基2FFT算法一般情况下，x(n)为复数序列，对某一个k值，直接按(4.2.1)式计算一个系数X(k)值需要：复数乘法：N次复数加法：(N-1)次计算全部X(k)值需要：N2次N(N-1)次当N1时，N(N-1)N2；随着N增大复乘、复加次数非线性增大。10()(),0,1,,1NknNnXkxnWkN(4.2.1)4.2.1直接计算DFT的特点及减少运算量的基本途径一、直接计算DFT的运算量长度为N的有限长序列x(n)的DFT为：第4章快速傅里叶变换(FFT)二、减少运算量的基本途径在DFT定义式中只有两种运算：x(n)与WmN的乘和加，WmN的特性对运算量必有影响。（1）根据旋转因子WmN的周期性、对称性和特殊值减少乘法运算次数。WmN的可约性为：WmN的对称性表现为：WmN的周期性表现为：22()jmlNjmmlNmNNNNWeeW(4.2.2)2[]mNmNmmNNNNNmmNN或者(4.2.3)WmN的特殊值为：j420;1;1nNrnrNWW第4章快速傅里叶变换(FFT)（2）把N点DFT分解为几个较短DFT的组合，可使乘法次数大大减少。FFT算法就是不断地把长序列的DFT分解成几个短序列的DFT，并利用旋转因子的周期性和对称性来减少DFT的运算次数。因为，DFT的运算次数N2若序列长度N，则运算次数。最常用的算法是基2-FFT(N＝2M)。第4章快速傅里叶变换(FFT)4.2.2时域抽取法（DIT）基2-FFT基本原理根据减少运算量的途径，巧妙地在时域或频域进行不同的抽取分解与组合，可得到不同的快速算法。FFT算法基本上分为两大类：时域抽取法FFT(DecimationInTime-FFT，DIT-FFT)频域抽取法FFT(DecimationInFrequency-FFT，DIF-FFT)我们介绍基2DIT-FFT算法。所谓基2是序列x(n)的长度N满足：2,MNM为自然数（1）算法思想：进行时域M级奇偶抽取，并利用：将N点DFT变成M级蝶形运算。（2）特点：运算流程图结构规则，可原位计算，程序简单。第4章快速傅里叶变换(FFT)设序列x(n)的长度为N，且满足：2,MNM为自然数按n的奇偶把x(n)分解为两个N/2点的子序列：12()(2),0,1,12()(21),0,1,12NxrxrrNxrxrr则x(n)的DFT为：第4章快速傅里叶变换(FFT)所以/21/211/22/21200()()()()()NNkrkkrkNNNNrrXkxrWWxrWXkWXk由于222222/2jkrNjkrkrkrNNNWeeW其中X1(k)和X2(k)分别为x1(r)和x2(r)的N/2点DFTk=0,1,…,N-1由于X1(k)和X2(k)均以N/2为周期，则由复指数序列的周期性可得:)2()()2()(2211NkXkXNkXkX虽然k=0,1,…,N-1，但从此式只能方便地得到N/2个DFT系数，那么，还有N/2个DFT系数如何方便地得到呢？第4章快速傅里叶变换(FFT)所以:2NkkNNWW另外由旋转因子WmN的对称性：1212()()()0,1,12()()()0,1,122kNkNNXkXkWXkkNNXkXkWXkk(4.2.7)(4.2.8))()()2()2()2(21221kXWkXNkXWNkXNkXkNNkN这样将N点DFT分解为两个N/2点的DFT运算，用如图所示的流图符号（蝶形）表示。要完成一个蝶形运算，需一次复数乘和两次复数加法运算。仅经过一次分解，就使运算量减少近一半：N2(N2/4)+(N2/4)+N/2第4章快速傅里叶变换(FFT)图4.2.2N点DFT的一次时域抽取分解图(N=8)N/2点DFTWN0N/2点DFTWN1WN2WN3x(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章快速傅里叶变换(FFT)由于N=2M，与第一次分解相同，将x1(r)按奇偶分解成两个N/4长的子序列x3(l)和x4(l)，即：那么，X1(k)又可表示为：式中：N/4点DFT13/2413/24()()(),0,1,,/41(/4)()()kNkNXkXkWXkkNXkNXkWXk(4.2.10)3241()(2),0,1,,1()(21)4xlxlNlxlxl第4章快速傅里叶变换(FFT)用同样的方法，将x2(r)按奇偶分解，可得：其中/4155/450/4166/4605262()()[()]()()[()]()(2),0,1,/41()(21)NklNiNklNiXkxlWDFTxlXkxlWDFTxlxlxllNxlxll=0l=0N/4点DFT第4章快速傅里叶变换(FFT)图4.2.3N点DFT的第二次时域抽取分解图(N=8)N/4点DFTWN12WN12WN0WN1WN2WN3X1(0)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)x(0)X3(0)X3(1)X4(0)X4(1)x(4)x(2)x(6)x(1)x(5)x(3)x(7)N/4点DFTN/4点DFTN/4点DFTWN02WN0213/2413/24()()(),0,1,,/41(/4)()()kNkNXkXkWXkkNXkNXkWXk(4.2.10)第4章快速傅里叶变换(FFT)图4.2.4N点DIT-FFT运算流图(N=8)WN0WN1WN2WN3WN0WN2WN0WN2WN0WN0WN0WN0x(0)x(4)x(2)x(6)x(1)x(5)x(3)x(7)A(0)A(1)A(2)A(3)A(4)A(5)A(6)A(7)A(0)A(1)A(2)A(3)A(4)A(5)A(6)A(7)A(0)A(7)X(0)X(1)X(2)X(3)X(4)X(5)X(6)X(7)A(0)A(1)A(2)A(3)A(4)A(5)A(6)A(7)N=8时，x3(l)和x4(l)已是两个长度为2的子序列：经过M-1次分解，将N点DFT分解成了N/2个2点的DFT，因此N=2M时，FFT总共需要M级运算，每级N/2个蝶形。WN/21=WN2第4章快速傅里叶变换(FFT)4.2.3DIT―FFT算法与直接计算DFT运算量的比较每一级运算都需要N/2次复数乘和N次复数加(每个蝶形需要两次复数加法)。所以，M级运算共需运算量为：复数乘法：m(M)=（N/2）M=（N/2）log2N复数加法：a(M)=NM=Nlog2N例如，N=210=1024时221048576204.8(/2)log5120NNN直接计算FFT算法在实际运算中，FFT的运算量还要减少，因为WNN=1，WN0=1，WNN/2=-1，WNN/4=-j，此时不需进行复数乘法。第4章快速傅里叶变换(FFT)4.2.4DIT-FFT的实现1.原位计算DIT-FFT的运算过程很有规律。N=2M点的FFT共进行M级运算，每级由N/2个蝶形运算组成。同一级中，每个蝶形的两个输入数据只对计算本蝶形有用，而且每个蝶形的输入、输出数据结点又同在一条水平线上，即计算完一个蝶形后，所得输出数据可立即存入原输入数据所占用的存贮单元。这种利用同一存贮单元存贮蝶形计算输入、输出数据的方法称为原位(址)计算。原位计算可节省大量内存，使设备成本降低。第4章快速傅里叶变换(FFT)2．旋转因子的变化规律如上所述，N点DIT-FFT运算流图中，每级都有N/2个蝶形。每个蝶形都要乘以因子，称其为旋转因子，p为旋转因子的指数。但各级的旋转因子和循环方式都有所不同。pNW用L表示从左到右的运算级数(L=1，2，…，M)。观察图4.2.4不难发现，第L级共有2L－1个不同的旋转因子。第4章快速傅里叶变换(FFT)对N=2M的一般情况，第L级的旋转因子为第4章快速傅里叶变换(FFT)第4章快速傅里叶变换(FFT)第4章快速傅里叶变换(FFT)5、序列的倒序DIT-FFT算法的输入序列的排序看起来似乎很乱，但实际排序是很有规律的。实现按时间抽取FFT算法的关键是将输入数据排列成满足连续运用奇偶分解所需的次序。N=2M，原输入序列的顺序表示为：n=(nM-1nM-2…n1n0)2对输入数据次序的变化可根据一个简单的位对换规则进行（称为倒位序）：n’=(n0n1…nM-2nM-1)2当把输入数据进行了重新排序，则输出结果是正确的次序。例：n=23，输入序“6”重新排序后为“3”，因为6=1102倒位序后n’=0112=3。第4章快速傅里叶变换(FFT)01234567奇偶抽取02461357奇偶抽取04261537倒位序表4.2.1顺序和倒序二进制数对照表第4章快速傅里叶变换(FFT)4.2.5频域抽取基－2FFT算法算法的推导频域抽取算法是把时间序列前后对半分解为两个长为N/2点的序列，则：)(nx第4章快速傅里叶变换(FFT)∴的N点DFT按k的奇偶分组可分为两个N/2的DFT)(kX)(nx为奇数为偶数kkWkkNN11)1(2/rnNNnWNnxnxrX2120)]2()([)2(rnNNnWNnxnx2120)]2()([rnNNnWNnxnxrX2120)]2()([)12(nNrnNNnWWNnxnx2120)]2()([nNNNWnxnxnxnxnxnx)]()([)()()()(2221令：第4章快速傅里叶变换(FFT)12/02/1)()2(NnrnNWnxrX则12/02/2)()12(NnrnNWnxrX这一结论表明：求的N点DFT再次分解成求两个N/2点DFT)(nx)(kXDIF-FFT的蝶式运算流图)(nx)2/(Nnx)2/()(NnxnxnNWNnxnx)]2/()([1nNW第4章快速傅里叶变换(FFT)DIF-FFT的一次分解运算流图先蝶式运算，后DFT。例如：N=8时第4章快速傅里叶变换(FFT)DIF-FFT的二次分解运算流图通常N/2仍然为2的整数幂，继续将N/2点DFT分成偶数组和奇数组，这样每个N/2点DFT又可分解成两个N/4点DFT，其输入序列分别是和按上下对半分开后通过蝶式运算构成的4个子序列，如下图所示：)(1nx)(2nx第4章快速傅里叶变换(FFT)按照以上方法继续分解下去，经过M-1次分解，最后分解为N/2个两点DFT，这N/2个2点DFT的输出就是N点DFT的结果X(k)，如下图所示：第4章快速傅里叶变换(FFT)有关说明以上给出了N=8时完整的DIF-FFT的运算流图。由于这种方法是按在频域进行奇偶分解，因此称之为频域抽取基－2ＦＦＴ运算。比较ＤＩＦ－ＦＦＴ与ＤＩＴ－ＦＦＴ相同点：运算次数与存储量相同