第四章快速傅里叶变换FFT.ppt

第四章第四章快速傅里叶变换快速傅里叶变换（（FFT））2021/8/61主要内容主要内容qDIT-FFT算法算法 qDIF-FFT算法算法qIFFT算法算法qChirp-FFT算法算法q线性卷积的线性卷积的FFT算法算法§4.1 引言引言qFFT: Fast Fourier Transformq1965年，年，Cooley-Turky 发表文章发表文章《机器计算傅《机器计算傅里叶级数的一种算法》，提里叶级数的一种算法》，提出出FFT算法，解决算法，解决DFT运算量太大，在实际使用中受限制的问题运算量太大，在实际使用中受限制的问题qFFT的应用频谱分析、滤波器实现、实时信号的应用频谱分析、滤波器实现、实时信号处理等qDSP芯片实现芯片实现TI公司的公司的TMS 320c30，，10MHz时钟，基时钟，基2-FFT1024点点FFT时间时间15msq典型应用：信号频谱计算、系统分析等典型应用：信号频谱计算、系统分析等 系统分析系统分析 频谱分析与功率谱计算频谱分析与功率谱计算§4.2 直接计算直接计算DFT的问题及改进途径的问题及改进途径1、、 DFT与与IDFT2、、DFT与与IDFT运算特点运算特点复数乘法复数乘法复数加法复数加法一个一个X(k)NN – 1N个个X(k)(N点点DFT)N 2N (N – 1)同理：同理：IDFT运算量与运算量与DFT相同。

相同实数乘法实数乘法实数加法实数加法一次复乘一次复乘42一次复加一次复加2一个一个X (k)4N2N+2 (N – 1)=2 (2N – 1)N个个X (k)(N点点DFT)4N 22N (2N – 1)3、降低、降低DFT运算量的考虑运算量的考虑FFT算法分类算法分类:q时间抽选法时间抽选法DIT: Decimation-In-Timeq频率抽选法频率抽选法DIF: Decimation-In-Frequency§4.3 按时间抽取（按时间抽取（DIT）的）的FFT算法算法（（Decimation  In  Time））1、算法原理、算法原理设序列点数设序列点数 N = 2L，，L 为整数 若不满足，则补零若不满足，则补零将序列将序列x(n)按按n的奇偶分成两组：的奇偶分成两组：N为为2的整数幂的的整数幂的FFT算法称算法称基基-2FFT算法算法将将N点点DFT定义式分解为两个长度为定义式分解为两个长度为N/2的的DFT记：记： ……… ………（（1 1））（这一步利用：（这一步利用： ））再利用周期性求再利用周期性求X(k)的后半部分的后半部分将上式表达的运算用一个专用将上式表达的运算用一个专用“蝶形蝶形”信流图表示。

信流图表示注：注：a. 上支路为加法，下支路为减法；上支路为加法，下支路为减法；        b. 乘法运算的支路标箭头和系数乘法运算的支路标箭头和系数用用““蝶形结蝶形结””表示上面运算的分解：表示上面运算的分解： 分解后的运算量：分解后的运算量：复数乘法复数乘法复数加法复数加法一个一个N/2点点DFT(N/2)2N/2 (N/2 –1)两个两个N/2点点DFTN2/2N (N/2 –1)一个蝶形一个蝶形12N/2个蝶形个蝶形N/2N总计总计运算量减少了近一半运算量减少了近一半进一步分解进一步分解由于由于 ，， 仍为偶数，因此，两个仍为偶数，因此，两个 点点DFTDFT又可同样进一步分解为又可同样进一步分解为4 4个个 点的点的DFTDFT““蝶形蝶形””信流图表示信流图表示 N点点DFT分解为四个分解为四个N/4点的点的DFTq类似的分解一直继续下去，直到分解为最后类似的分解一直继续下去，直到分解为最后的两类蝶形运算为止的两类蝶形运算为止(2点点DFT).q如上述如上述N=8=23，，N/4=2点中：点中： 类似进一步分解类似进一步分解1点点DFTx(0)1点点DFTx(4)X3(0)X3(1)进一步简化为蝶形流图：进一步简化为蝶形流图：X3(0)X3(1)x(0)x(4)因此因此8 8点点FFTFFT时间抽取方法的信流图如下时间抽取方法的信流图如下————FFT运算量与运算特点运算量与运算特点 1．． N=2L时，共有时，共有L=log2N级运算；每一级有级运算；每一级有N/2个蝶形结。

个蝶形结2．每一级有．每一级有N个数据（中间数据），且每级只用到本个数据（中间数据），且每级只用到本级的转入中间数据，适合于迭代运算级的转入中间数据，适合于迭代运算3．计算量：．计算量： 每级每级N/2次复乘法，次复乘法，N次复加每蝶形只乘一次，次复加每蝶形只乘一次，加减各一次）共有加减各一次）共有L*N/2=N/2log2N 次复乘法；次复乘法；复加法复加法L*N=Nlog2N 次与直接次与直接DFT定义式运算定义式运算量相比量相比(倍数倍数) N2/(Nlog2N) 当 N大时，此倍数大时，此倍数很大比较比较DFT 参考参考P150 表表4-1 图图4-6可以直观看出，当点数可以直观看出，当点数N越大时，越大时，FFT的优点更突出的优点更突出按时间抽取按时间抽取FFT蝶形运算特点蝶形运算特点 1、关于、关于FFT运算的混序与顺序处理（位倒序处理）运算的混序与顺序处理（位倒序处理） 由于输入序列按时间序位的奇偶抽取，故输入序由于输入序列按时间序位的奇偶抽取，故输入序列是混序的，为此需要先进行混序处理列是混序的，为此需要先进行混序处理混序规律：混序规律： x(n)按按n位置进行码位（二进制）倒置位置进行码位（二进制）倒置规律输入，而非自然排序，即得到混序排列。

所规律输入，而非自然排序，即得到混序排列所以称为位倒序处理以称为位倒序处理位倒序实现：位倒序实现：（（1））DSP实现采用位倒序寻址实现采用位倒序寻址（（2）通用计算机实现可以有两个方法：一是严格按）通用计算机实现可以有两个方法：一是严格按照位倒序含义进行；二是倒进位的加照位倒序含义进行；二是倒进位的加N/2倒位序倒位序自然序自然序0000000010041001010220101106301100114100101551010113611011177111倒位序倒位序例例 计算计算 ，， 计算 点点FFTFFT用时间抽取输入倒序算法，用时间抽取输入倒序算法，问倒序前寄存器的数问倒序前寄存器的数 和倒序后和倒序后 的的数据值？数据值？解：倒序前解：倒序前 倒序倒序 倒序为倒序为 倒序后倒序后 DIT FFT中最主要的蝶形运算实现中最主要的蝶形运算实现（（1）参与蝶形运算的两类结点（信号）间）参与蝶形运算的两类结点（信号）间“距离距离”（码地址）与其所处的第几级蝶形有关；第（码地址）与其所处的第几级蝶形有关；第m级的级的“结距离结距离”为为 (即原位计算迭代）即原位计算迭代）（（2）每级迭形结构为）每级迭形结构为q 蝶形运算两节点的第一个节点为蝶形运算两节点的第一个节点为k值，表值，表示成示成L位二进制数，左移位二进制数，左移L – m位，把右位，把右边空出的位置补零，结果为边空出的位置补零，结果为r的二进制数。

的二进制数3）） 的确定：的确定： 第第m级的级的r取值：取值：DIT算法的其他形式流图算法的其他形式流图q输入倒位序输出自然序输入倒位序输出自然序q输入自然序输出倒位序输入自然序输出倒位序q输入输出均自然序输入输出均自然序q相同几何形状相同几何形状q输入倒位序输出自然序输入倒位序输出自然序q输入自然序输出倒位序输入自然序输出倒位序参考参考P154-155时间抽取、时间抽取、 输入自然顺序、输入自然顺序、 输出倒位序的输出倒位序的FFTFFT流图流图         例例  用用FFT算算法法处处理理一一幅幅N×N点点的的二二维维图图像像，，如如用用每每秒秒可可做做10万万次次复复数数乘乘法法的的计计算算机机，，当当N=1024时时，，问问需需要要多多少少时时间间（（不不考虑加法运算时间）？考虑加法运算时间）？        解解 当当N=1024点点时时，，FFT算算法法处处理理一一幅幅二二维维图图像像所所需需复复数数乘乘法约为法约为      次次，，仅仅为为直直接接计计算算DFT所所需需时时间间的的10万万分之一 即原需要即原需要3000小时，现在只需要小时，现在只需要2 分钟。

分钟 §4.4 按频率抽取（按频率抽取（DIF）的）的FFT算法算法q与与DIT-FFT算法类似分解，但是抽取的是算法类似分解，但是抽取的是X(k)即分解即分解X(k)成奇数与偶数序号的两个序列成奇数与偶数序号的两个序列q设：设： N = 2L，，L 为整数将为整数将X(k)按按k的奇偶分的奇偶分组前，先将输入组前，先将输入x(n)按按n的顺序分成前后两半：的顺序分成前后两半：（（Decimation In Frequency)一、算法原理一、算法原理下面讨论下面讨论按按k k的奇偶将的奇偶将X(k)X(k)分成两部分：分成两部分：显然：显然：令：令：用蝶型结构图表示为：用蝶型结构图表示为：x1(0)x1(1)-1x1(2)x1(3)-1x2(0)x2(1)-1x2(2)x2(3)-1N/2点DFTN/2点DFTx(0)x(7)x(1)x(2)x(3)x(4)x(5)x(6)X1(0)=X(0)X2(0)=X(1)X1(1)=X(2)X1(2)=X(4)X1(3)=X(6)X2(1)=X(3)X2(2)=X(5)X2(3)=X(7)N/2仍为偶数，进一步分解：仍为偶数，进一步分解：N/2 →→ N/4x3(0)x3(1)-1-1x4(0)x4(1)N/4点DFTN/4点DFTx1(0)x1(1)x1(2)x1(3)X3(0)=X1(0)=X(0)X4(0)=X1(1)=X(2)X3(1)=X1(2)=X(4)X4(1)=X1(3)=X(6)q按照以上思路继续分解，即一个按照以上思路继续分解，即一个N/2的的DFT分解成两个分解成两个N/4点点DFT，，直到只计算直到只计算2点的点的DFT，，这就是这就是DIF-FFT算法。

算法2个个1点的点的DFT蝶形流图蝶形流图 进一步简化为蝶形流图：进一步简化为蝶形流图：1点点DFTx3(0)1点点DFTx3(1)X(0)X(4)X(0)X(4)x3(0)x3(1)二、按频率抽取二、按频率抽取FFT蝶形运算特点蝶形运算特点1）原位计算）原位计算-1L级蝶形运算，每级级蝶形运算，每级N/2个蝶形，每个蝶形结构：个蝶形，每个蝶形结构：m表示第表示第m级迭代，级迭代，k，，j表示数据所在的行数表示数据所在的行数2）蝶形运算）蝶形运算对对N=2L点点FFT，输入自然序，输出倒位序，，输入自然序，输出倒位序，两节点距离：两节点距离：2L-m=N / 2m第第m级运算：级运算：q 蝶形运算两节点的第一个节点为蝶形运算两节点的第一个节点为k值，表示成值，表示成L位二进制数，左移位二进制数，左移m-1位，把右边空出的位位，把右边空出的位置补零，结果为置补零，结果为r的二进制数的二进制数存储单元存储单元输入序列输入序列x(n) : N个存储单元个存储单元系数系数 ：：N / 2个存储单元个存储单元三、三、DIT与与DIF的异同的异同q基本蝶形不同基本蝶形不同qDIT: 先复乘后加减先复乘后加减qDIF: 先减后复乘先减后复乘q运算量相同运算量相同q都可原位运算都可原位运算qDIT和和DIF的基本蝶形互为转置的基本蝶形互为转置§4.5 IDFT的的FFT算法算法（（FFT应用一）应用一） 一、从定义比较分析一、从定义比较分析与与DFT的比较：的比较： 1）、旋转因子）、旋转因子WN-kn 的不同；的不同； 2）、结果还要乘）、结果还要乘 1/N。

二、实现算法二、实现算法——直接使用直接使用FFT程序的算法程序的算法共轭共轭FFT共轭共轭乘乘1/ N直接调用直接调用FFT子程序计算子程序计算IFFT的方法：的方法：§4.6 线性调频线性调频Z变换（变换（Chirp-Z变换）算法变换）算法 （（FFT应用二）应用二）单位圆与非单位圆采样单位圆与非单位圆采样（（a）） 沿单位圆采样沿单位圆采样;    (b) 沿沿AB弧采样弧采样  螺线采样螺线采样 zk=AW-k k=0, 1, …, M-1 Chirp-Z变换的线性系统表示变换的线性系统表示 由于系统的单位脉冲响应由于系统的单位脉冲响应 可以想象为频率可以想象为频率随时间随时间(n)(n)呈线性增长的复指数序列在雷达系统中呈线性增长的复指数序列在雷达系统中, ,这这种信号称为种信号称为线性调频信号（线性调频信号（Chirp SignalChirp Signal）），因此，这，因此，这里的变换称为里的变换称为线性调频线性调频Z Z变换变换一、基本算法思路一、基本算法思路§4.7 线性卷积的线性卷积的FFT算法算法（（FFT应用三）应用三）若若L点点x(n)，，M点点h(n)，，则直接计算其线性卷积则直接计算其线性卷积y(n)需运算量：需运算量：若系统满足线性相位，即：若系统满足线性相位，即：则需运算量：则需运算量：FFT法：以圆周卷积代替线性卷积法：以圆周卷积代替线性卷积N总运算量：总运算量： 次乘法次乘法比较直接计算和比较直接计算和FFT法计算的运算量法计算的运算量讨论：讨论：1）当）当2）当）当见见P183 表表 x(n)长长度度很很长长时时，，将将x(n)分分为为L长长的的若若干干小的片段，小的片段，L与与M可比拟。

可比拟1 1、重叠相加法、重叠相加法 则：则： 输出：输出：其中：其中：可以用圆周卷积计算：可以用圆周卷积计算： 选选 ，上面圆周卷积可用，上面圆周卷积可用FFTFFT计算 N    由于由于yi(n)长度为长度为N，而，而xi(n)长度长度L ，必有，必有M-1点重叠，点重叠， yi(n)应相加才能构成最后应相加才能构成最后y(n)的重叠相加法图形重叠相加法图形和上面的讨论一样，和上面的讨论一样， 用用FFT法实现重叠相加法的步骤如下法实现重叠相加法的步骤如下: ①① 计算计算N点点FFT，， H(k)=DFT［［h(n)］；］；②② 计算计算N点点FFT，，Xi(k)=DFT［［xi(n)］；］；③③ 相乘，相乘，Yi(k)=Xi(k)H(k)；；④④ 计算计算N点点IFFT，，yi(n)=IDFT［［Yi(k)］；］；⑤⑤  将各段将各段 yi(n)（包括重叠部分）相加，（包括重叠部分）相加，  重叠相加的名称是由于各输出段的重叠部分相加而得名的重叠相加的名称是由于各输出段的重叠部分相加而得名的 例例 已知序列已知序列x[n]=n+2,0 n 12, h[n]={1,2,1}试试利用重叠相加法计算线性卷积利用重叠相加法计算线性卷积, 取取L=5 。

y[n]={2, 7, 12, 16, 20, 24, 28, 32, 36, 40, 44, 48, 52, 41, 14}解解: 重叠相加法重叠相加法x1[n]={2, 3, 4, 5, 6}x2[n]={7, 8, 9, 10, 11}x3[n]={12,13, 14}y1[n]={2, 7, 12, 16, 20, 17, 6}y2[n]={ 7, 22, 32, 36, 40, 32, 11}y3[n]={12, 37, 52, 41, 14}2 2、重叠保留法、重叠保留法q此方法与上述方法稍有不同此方法与上述方法稍有不同q先将先将x(n)分段，每段分段，每段L=N-M+1个点，这是相同的个点，这是相同的q不同之处是，序列中补零处不补零，而在每一段不同之处是，序列中补零处不补零，而在每一段的前边补上前一段保留下来的（的前边补上前一段保留下来的（M-1）个输入序列）个输入序列值，值， 组成组成L+M-1点序列点序列xi(n) q如果如果L+M-1<2m，， 则可在每段序列末端补零值点，则可在每段序列末端补零值点，补到长度为补到长度为2m，这时如果用，这时如果用DFT实现实现h(n)和和xi(n)圆周卷积，则其每段圆周卷积结果的前（圆周卷积，则其每段圆周卷积结果的前（M-1）个点的值不等于线性卷积值，必须舍去。

个点的值不等于线性卷积值，必须舍去重叠保留法示意图重叠保留法示意图重叠保留法示意图重叠保留法示意图y[k]={2, 7, 12, 16, 20, 24, 28, 32, 36, 40, 44, 48, 52, 41, 14}x1[k]={0, 0, 2, 3, 4}x2[k]={3, 4, 5, 6 ,7}x3[k]={6 ,7 , 8, 9, 10}y1[k]= x1[k] h[k]= {11,  4,  2,  7,  12}x4[k]={9, 10 , 11, 12,13}y2[k]= x2[k] h[k]= {23,  17,  16,  20,  24}y3[k]= x3[k] h[k]= {35,  29,  28,  32,  36}y4[k]= x4[k] h[k]= {47,  41,  40,  44,  48}x5[k]={12,13, 14, 0, 0}y5[k]= x5[k] h[k]= {12,   37,   52,   41,   14}解解: : 重叠保留法重叠保留法例例 已知序列已知序列x[n]=n+2,0 n 12, h[n]={1,2,1}试试利用重叠利用重叠保留法保留法计算线性卷积计算线性卷积, 取取L=5 。

语音信号消噪过程语音信号消噪过程(a)信号淹没在啸叫噪声中；信号淹没在啸叫噪声中；      (b) 信号与噪声的功率谱；信号与噪声的功率谱； (b)(c) 去噪后的功率谱；去噪后的功率谱； (d) 重构原语音信号重构原语音信号FFT应用举例应用举例。