毕业论文-快速傅里叶变换算法及其在信号处理中的应用

1、武汉工程大学毕业设计(论文)说明书 2015 届毕业设计(论文) 题 目快速傅里叶变换算法及其在信号处理中的应用专 业 班 级2011电子信息工程02 学 号1104030231 姓 名周汝耀 指 导 教 师华夏讲师 学 院 名 称电气信息学院 2011 年 6 月 9 日摘 要快速傅氏变换（FFT），是离散傅氏变换的快速算法，它是根据离散傅氏变换的奇、偶、虚、实等特性，对离散傅立叶变换的算法进行改进获得的。它对傅氏变换的理论并没有新的发现，但是对于在计算机系统或者说数字系统中应用离散傅立叶变换，可以说是进了一大步。傅里叶变换的理论与方法在“数理方程”、“线性系统分析”、“信号处理、仿真”等很

2、多学科领域都有着广泛应用,由于计算机只能处理有限长度的离散的序列,所以真正在计算机上运算的是一种离散傅里叶变换.虽然傅里叶运算在各方面计算中有着重要的作用，但是它的计算过于复杂，大量的计算对于系统的运算负担过于庞大，使得一些对于耗电量少，运算速度慢的系统对其敬而远之，然而，快速傅里叶变换的产生，使得傅里叶变换大为简化，在不牺牲耗电量的条件下提高了系统的运算速度，增强了系统的综合能力，提高了运算速度，因此快速傅里叶变换在生产和生活中都有着非常重要的作用，对于学习掌握都有着非常大的意义。 关键词：快速傅氏变换；快速算法；简化；广泛应用AbstractFast Fourier Transform (

3、FFT), is a discrete fast Fourier transform algorithm, which is based on the Discrete Fourier Transform of odd and even, false, false, and other characteristics of the Discrete Fourier Transform algorithms improvements obtained. Its Fourier transform theory has not found a new, but in the computer sy

4、stem or the application of digital systems Discrete Fourier Transform can be said to be a big step into. Fourier transform theory and methods in the mathematical equation and linear systems analysis and signal processing, simulation, and many other areas have a wide range of applications, as the com

5、puter can only handle a limited length of the sequence of discrete, so true On the computers operation is a discrete Fourier transform. Fourier Although all aspects of computing in the calculation has an important role, but its calculation was too complicated, a lot of computing system for calculati

6、ng the burden is too large for some Less power consumption, the slow speed of operation of its system at arms length, however, have the fast Fourier transform, Fourier transform greatly simplifying the making, not in power at the expense of the conditions to increase the speed of computing systems,

7、and enhance the system The comprehensive ability to improve the speed of operation, the Fast Fourier Transform in the production and life have a very important role in learning to master all have great significance.KeyWords：Fast Fourier Transform; fast algorithm; simplified; widely used 目 录摘要.IAbstr

8、act . II1.绪论1.1 选题背景11.2 课题研究的意义22.快速傅里叶变换原理及性质2.1快速傅里叶变换原理32.2快速傅里叶变换的优越性42.3快速傅里叶变换的意义43.快速傅里叶变换的算法3.1快速傅里叶变换算法63.2 Cooley=Tukey FFT算法83.3 Rader-Brennr FFT算法93.4 Goertsel 算法104.快速傅里叶变换在信号处理中的理论应用4.1利用FFT计算连续时间信号的傅里叶变换134.2利用FFT计算离散信号的线性卷积174.3利用FFT进行离散信号压缩194.4利用FFT对离散信号进行滤波224.5利用FFT提取离散信号中的最强正弦分

9、量245.快速傅里叶变换在数字信号分析与处理的实际应用5.1快速傅里叶变换在喇曼光谱信号噪声平滑中的应用295.2采用异步实现的快速傅里叶变换处理器315.3快速傅里叶算法在哈特曼夏克传感器波前重构算法中的应用33致谢36参考文献37 1 绪论傅立叶变换在生产生活中的重要性非常突出，它将原本难以处理的时域信号相比比较容易地转换成易于分析的频域信号，我们可以利用一些专业工具对这些频域信号进行处理、加工，使信号转化为可以对其进行各种数学变换的数学公式，对其进行处理。最后还可以根据傅立叶反变换将这些频域信号转换成原来的时域信号，这是一种特殊的积分变换。它能够将满足一定条件的某个函数表示成为正弦基函数

10、的线性组合或者积分。然而，它在运算上过于复杂，过于宏大的运算过程，对于一些相对简单的低功耗处理器来说，难以自如应对，因此，快速傅里叶变换则显出了它的优越性。快速傅氏变换（FFT），即离散傅氏变换的快速算法，它是根据离散傅氏变换的奇、偶、虚、实等特性，对离散傅立叶变换的算法进行改进获得的1。对于计算机处理信号方面上是一大进步。系统的速度不但取决于其本身的速度，而且还在相当大的程度上取决于运用的算法，算法运算量的大小直接影响到对设备的控制质量。通过傅立叶变换（DFT），运用测试软件进行检测，我们可以看出，快速傅里叶变换大大的提高了运算速度，它为各系统的设计方案提供了简单算法，有着非常重要的意义。1

11、. 1 选题背景近十多年来，数字信号处理技术同大规模集成电路、数字计算机等，都有了突飞猛进的发展，日新月异，早已成为了一门具有强大生命力的技术科学。因为它本身就具有一系列的优点，所以能够有效地促进工程技术领域的技术改造和科学发展，应用领域也更加地广泛、深入，越来越受到人们的重视。在信号处理中，离散傅里叶变换（Discrete Fourier Transform，DFT）是比较常用的变换方法之一，它在各种数字信号处理系统中扮演着及其重要的角色。由于离散傅里叶变换（DFT）而发现了频率离散化，可以直接用它来分析信号的频谱、计算滤波器的频率响应、以及实现信号通过线系统的卷积运算等，因而在信号的谱分析

12、等方面有着非常大的作用。傅里叶变换已经有一百多年的历史了，我们熟知频域分析往往比时域分析更优越，不仅简单明了，而且易于分析较为复杂的信号。但需要用较为精准的数字方法，即DFT，进行谱分析，在快速傅氏变换（FFT）出现以前是不切实际的。由于DFT的计算量太大，即使运用计算机也很难对问题进行实时的有效处理，所以DFT并没有得到真正的应用。直到1965年库利（J.W.Cooly）和图基（J.W.Tukey）首次发现DFT的一种快速算法，局面才发生根本性的变化。继库利和图基算法出来之后，桑德（G.Sander）等快速算法也相继出现，又经过其他学者一步步改进，很快就出现了通用型的快速傅里叶变换，简称FF

13、T。快速傅里叶变换（Fast Fourier Transform，FFT）并非与离散傅里叶变换完全不同的另一种变换，而是为了减少DFT计算次数而诞生的一种快速、有效的算法。应当指出的是，也是因为当时电子数字计算机的“落后”条件也促成了这个算法的提出。它使得DFT的运算量大大的缩小简化，它推动了近30年来信号处理技术止步不前的前进发展，成为了数字信号处理应用领域里强有力的工具，为DFT乃至数字信号处理技术的实际应用创造了良好的条件，从而使DFT在实际使用中得以广泛的应用2。数字信号处理器（DSP），是一种可编程的高性能处理器。近年来发展尤为迅速，它不仅应用于数字信号处理方面，而且在图像处理、语音

14、处理、通信等领域得到广泛的应用。之前通用的微处理器在运算速度上已经很难适应信号实时处理的高要求。DSP处理器中集成了高速的乘法硬件，能快速、准确地进行大量数据的乘法以及加法的运算。数字信号处理区别于普通的科学计算与分析，它强调运算的实时性。除了需要普通微处理器所强调的高速运算和控制能力之外，鉴于实时数字信号处理的特点，在处理器结构、指令系统、指令流程上做了很大程度上的改进。1. 2 课题研究的意义如上所述，基于对DSP的快速傅里叶变换算法的研究，从而使FFT算法能够有效地在DSP芯片上实现。DSP芯片的出现，使FFT的实现更加方便。多数的DSP芯片都能够在一个指令周期内完成一次乘法和加法，并且

15、提供了专门的FFT指令，完成一次指令的周期只需10ns，使得FFT算法在DSP芯片上实现的速度更加快速。快速傅里叶变换为频谱分析、卷积、相关数字滤波器设计与实现与功率谱计算、传递函数建模、图像处理等，提供了快速有效的运算方法。FFT技术应用DSP芯片，从而可以提供使调制、解调、压缩、解压缩的数据传输更为高效的信号处理解决方案，因而广泛应用于雷达、图像处理、通信、生物医学和声纳领域。2.快速傅里叶变换原理及性质 数字信号中的傅里叶变换,通常是采用离散型傅里叶变换(DFT)。DFT 存在的缺点就是计算量太大,不易进行实时处理。比如，计算一个N 点的DFT ,一般需要次复数乘法和N(N-1)次复数加

16、法运算.因此,当N较大或要求对信号进行实时处理时,往往很难实现达到所需的运算速度。1965年,J.W.Cooly和J.W.Tukey发现了DFT的一种快速算法,经过后来学者的进一步改进, 很快便形成了一套高效的运算方法,即现在通用的快速傅里叶变换, 简称FFT( The Fast Fourier Transform)。快速傅里叶变换的实质是利用式(3-1)中的权函数的对称性和周期性,把N点DFT进行一系列分解和组合,使整个DFT的计算过程变成一系列叠代运算过程,使DFT的运算量大大简化,为DFT及数字信号的实时处理和应用创造了非常良好的条件3。2. 1 快速傅里叶变换原理快速傅里叶变换原理：1

17、. 将长序列DFT分解为短序列的DFT2. 利用旋转因子的对称性、周期性、可约性。将时域序列逐次分解为一组子序列，依据旋转因子的特性，由子序列的DFT来实现整个序列的DFT4。其中：快速傅里叶变换分为两种，分为基2时间抽取算法和基2频率抽取算法基2-时间抽取(Decimation in time)FFT算法其中:r=0,1,2 （2-1）基2-频率抽取(Decimation in frequency)FFT算法 （2-2）2. 2 快速傅里叶变换的优越性设为项的复数序列，依据DFT变换，任一的计算都需要有次复数乘法和（）次复数加法，而且一次复数乘法等同于四次实数乘法和两次实数加法，同样的，一次

18、复数加法等同于两次实数加法，即使我们把一次复数乘法和一次复数加法定义成一次“运算”（四次实数乘法和四次实数加法），那么求出项复数序列的, 即N点DFT变换大约就需要次运算。当点甚至更多的时候，需要N2=1048576次运算。在FFT中，利用WN的对称性和周期性，把一个N项序列（设,为正整数），分为两个项的子序列，而且每个点的DFT变换需要次运算，再运用N次运算把两个点的DFT 变换重新组合成一个N点的DFT变换。如此变换以后，总的运算次数就变成了。承接上面的例子，当时，总的运算次数就变成了525312次，这样看来，节省了大约50%的运算量。而如果我们将这种“一分为二”的思想不断进行下去，直到分

19、成两两一组的DFT运算单元，那么N点的 DFT变换就只需要次的运算，在1024点时，运算量仅有10240次，是先前的直接算法的1%，点数越多，运算量的节约就越大，这就是 FFT的优越性.当然，FFT提高了运算速度,但是,也对参与运算的样本序列作出了限制,即要求样本数为2N点.离散傅里叶变换DFT则无上述限制5。2. 3 快速傅里叶变换的意义 傅立叶变换的物理意义：傅立叶变换是数字信号处理技术领域一项很重要的算法。要知道傅立叶变换算法的意义，首先要了解傅立叶变换原理的意义。傅立叶变换原理表明：任何连续测量的时序或信号，都能够表示成为不同频率的正弦波信号的无限叠加。而利用该原理而创立的傅立叶变换算

20、法则利用直接能测量到的原始信号，并以累加方式来计算该信号中不同正弦波信号的频率、相位和振幅。与傅立叶变换算法对应的是反傅立叶变换算法。该反变换从本质上说也就是一种累加处理，这样便可以将单独改变的正弦波信号转换成一个信号。因此，也可以说，傅立叶变换将原来难以处理的时域信号转换成了易于分析的频域信号（信号的频谱），我们可以利用一些专业工具对这些频域信号进行加工、处理。最后还可以利用傅立叶反变换将这些频域信号转换成原来的时域信号。从现代数学的眼光来看，傅里叶变换其实就是一种特殊的积分变换。它能够将满足一定条件下的某个函数表示成为正弦基函数的线性组合或者积分形式。在不同的研究领域里，傅里叶变换具有多种

21、形式各异的变体形式，如连续傅里叶变换和离散傅里叶变换。在数学领域，尽管最初傅立叶分析是作为热过程的解析分析的工具，但是其思想方法仍然具有典型的还原论和分析主义的特征。任意的函数通过一定的分解，都能够表示为正弦函数的线性组合的形式，而正弦函数在物理上是被充分研究而相对简单的函数类：1. 傅立叶变换是线性算子,若赋予适当的范数,它还是酉算子;2. 傅立叶变换的逆变换容易求出,而且形式与正变换非常类似;3. 正弦基函数是微分运算的本征函数,从而使得线性微分方程的求解可以转化为常系数的代数方程的求解.在线性时不变杂的卷积运算为简单的乘积运算,从而提供了计算卷积的一种简单手段;4. 离散形式的傅立叶的物

22、理系统内,频率是个不变的性质,从而系统对于复杂激励的响应可以通过组合其对不同频率正弦信号的响应来获取;5. 著名的卷积定理指出:傅立叶变换可以化复变换可以利用数字计算机快速的算出(其算法称为快速傅立叶变换算法(FFT)。正是由于上述的良好性质,傅里叶变换在物理学、数论、组合数学、信号处理、概率、统计、密码学、声学、光学等领域都有着广泛的应用6。3. 快速傅里叶变换的算法3.1 快速傅里叶变换算法 快速傅里叶变换FFT 是离散傅里叶变换DFT 的一种快速算法，只有FFT 才能在现实中有实际应用的意义。 (3-1)由(3-1)式可知，对每一个n，计算X()须作N次复数乘法及N-1次复数加法，要完成

23、这组变换共需次乘法及N(N-1)次复数加法。但以下介绍的快速傅里叶变换的算法，可大大减少运算次数，提高工作效率。当时，n和k可用二进制数表示： (3-2) （3-3）又记,则（3-1）式可改写为 （3-4）式中： (3-5)因为所以（3-2）可改成 （3-6） (3-7)则式（3-7）即为式（3-6）的分解形式。将初始数据代入式（3-7）的第一个等式，可得每一组计算数据，一般将痗L-1组计算数据代入式（3-7）的第L个等式，计算后可得第L组计算数据（L，），计算公式也可表示为 = （3-8）式中 （3-9） 根据式（3-8），第L个数组中每个 的计算只依赖于上一个数组的两个数据这两个数据的标号

24、相差，即，而且这两个数据只用于计算第L个数组中标号的数据（等号右端为二进制数）。当分别取和时，分别有。因此，用上一组的两个数据计算所得的两个新数据仍可储存在原来位置，计算过程中只需要N个存储器。将与称为第L个数组中的对偶结点对。计算每个对偶结点对只需一次乘法，事实上由式（3-8）可得 （3-10）式中： ；别为式（3-9）中取，时对应的P值。因，于是对偶结点的有如下关系： ,因此式（3-8）可表示为 （3-11） P的求法：在中，i写成二进制数右移位，就成为颠倒位序得式（3-7）中，前面的r个等式，每个等式均对应一组数据进行计算，每组数据都有N/对结点，根据式（3-11），每对结点只需作次乘法

25、和次加法，因此，每组数据只需N/2次乘法和N次加法，因而完成r组数据的计算共需Nr/2次乘法和Nr次加法。3. 2 Cooley-Tukey FFT算法的核心是将一层运算映射为二层运算。作点时，若不是素数，则可分解为，那么由的 （3-12）通过映射： （3-13）可得到 而，可化简为 （3-14）从而式（3-12）转化为 （3-15）其中。以点为例，映射方式为：，则计算流图如图3-1所示。 n1 k2f0 0 W0 0 F0f4 1 W0 1 F5f8 2 W0 2 F10f12 3 W0 3 F15f16 4 W0 0 F1f1 0 W1 1 F6f5 1 W2 2 F11f9 2 W3 3

26、 F16f13 3 f17 4 W0 0 F2 W2 1 F7f2 0 W4 2 F12f6 1 W6 3 F17f10 2f14 3 W0 0 F3f18 4 W3 1 F8 W6 2 F13f3 0 W9 3 F18f7 1f11 2 W0 0 F4f15 3 W4 1 F9f19 4 W8 2 F14 W12 3 F19k1=0n2=0n2=1n2=2n2=3k1=1k1=2k1=3k1=4图3-1 Cooley-Tukey 20 点FFT算法3. 3 Rader-Brenner FFT算法Rader-Brenner 算法是类似于 Cooley-Tukey 算法，但是采用的旋转因子都是纯

27、虚数，以增加加法运算和降低了数值稳定性为代价减少了乘法运算。这方法之后被split-radix variant of Cooley-Tukey所取代，与Rader-Brenner算法相比，有一样多的乘法量，却有较少的加法量，且不牺牲数值的准确性7。 Bruun以及QFT算法是不断的把DFT分成许多较小的DFT运算。(Rader-Brenner以及QFT算法是为了2的指数所设计的算法，但依然可以适用在可分解的整数上。Bruun算法则可以运用在可被分成偶数个运算的数字)。尤其是Bruun算法，把FFT看成是 ，并把它分解成 与 的形式。 另一个从多项式观点的快速傅立叶变换法是Winograd 算法

28、。此算法把 分 解成cyclotomic多项式，而这些多项式的系数通常为1，0，-1。 这样只需要很少的乘法量(如果有需要的话)，所以winograd是可以得到最少乘法量的快速傅立叶算法，对于较小的数字，可以找出有效率的算方式。更精确地说，winograd算法让DFT可以用2K点的DFT来简化，但减少乘法量的同时，也增加了非常多的加法量。Winograd也可以利用剩余值定理来简化DFT。 Rader算法提出了利用点数为N(N为质数)的DFT进行长度为N-1的回旋折积来表示原本的DFT，如此就可利用折积用一对基本的FFT来计算 DFT。另一个prime-size的FFT算法为chirp-Z算法。

29、此法也是将DFT用折积来表示，此法与Rader算法相比，能运用在更一般的转换 上，其转换的基础为Z转换8。3.4 Goertsel 算法如前所述, 点时域序列的离散付里叶变换式为， （3-16）这点频域序列是同时被算出的，不可能只计算其中某一个或几个指定点。Goetzel 算法是为了解决这个问题而提出的。这个算法把离散付里叶变换看作一组滤波器，将输入端的时域序列与其中一个滤波器的冲激响应序列进行卷积运算，求滤波器的输出序列，即得序列的一点。这种算法利用旋转因子的周期性，使DFT运算化为线性滤波运算9。 由于 故式（3-16）可化为 （3-17）定义序列为 （3-18）可见是由两个序列卷积而得到

30、的序列。 （3-19）其中，是输入的点序列，另一个序列被看作滤波器的冲激响应序列。 （3-20）对比式(3-17）和式(3-18），可知：按式（3-19）进行卷积运算，当时，滤波器的输出就是： （3-21）对式(3-20)进行Z变换，可得滤波器的系统函数 （3-22）这是一个一阶系统。图3-2示出这个系统的信号流图，相应的差分方程为 ， (3-23) 按照此式进行递推运算，到了 时刻，即可依据式（3-21）得到。 按照式(3-22)进行运算时，可先算好旋转因子，储存起来。每次递推包含一次复数乘法。按式(3-16)直接计算点离散付里叶变换,需要次实数乘法和次实数加法。按照上述Goertzel算法

31、，所需的实数乘法和实数加法都是次。所以当不大时，上述算法的效率稍差10。下面介绍改进的Goertzel 算法，这种算法所需的实数乘法次数约为直接方法的一半。图3-2 用一阶系统实现Goertzel 算法 图 3-3 用二阶系统实现Goert算法把式（3-22）的分子和分母都乘上因子，就得到第个滤波器的系统函数为 （3-24）与此相应的信号流图示于图3-3。由式（3-24)可见，滤波器是一个二阶系统，有一对复数共轭极点和一个复数零点。为了便于运算，在图3-3所示的流图中，设立状态变量 和。按照图3-3计算时，步骤有二，即： 1.实现一对复数极点输入点依次取 ，进行递推运算。每次运算中，更新状态变

32、量和。作次迭代所需的计算量是次实数乘法和次实数加法。 2.实现复数零点。是一个点序列， 。在点上。计算状态变量和。这时，按照图3-2算出滤波器的输出，此即。所需的计算量是4次实数乘法和次实数加法11。 综上所述，计算一点需要进行次实数乘法和次实数加法。这种算法要求的乘法次数约为直接算法的一半。在这种较为有效的方案中，仍具有这样的优点，即必须计算和存储的系数只有和。还要说明图3-3所示的算法的另一个优点。当输入序列为实序列时，离散付里叶变换序列是对称的，即。容易证明，图3-3的网络形式在计算时和计算时具有完全相同的极点，但前者的零点系数与后者的零点系数成复共轭关系。由于零点仅在最后的迭代中实现，

33、所以诸极点要求的次乘法和 次加法可以用来计算离散付里叶变换的两个值。因此，若用Goertzel 算法计算离散付里叶变换的所有个点的值，需要的乘法次数近似为，加法次数近似为。然而，它同直接计算离散付里叶变换一样，计算量仍然正比于。4.快速傅里叶变换在信号处理中的理论应用4.1 连续时间信号的快速傅里叶变换设是连续时间信号，并假设时，则其傅里叶变换由下式给出： （4-1） 令是一固定的正实数，是一个固定的正整数。当时，利用FFT算法可计算。已知一个固定的时间间隔，选择让T足够小，使得每一个T秒的间隔内，的变化很小，则式中积分可近似为 （4-2）假设足够大，对于所有的整数，幅值很小，则式（4-2）变

34、为 （4-3）当时，式（4-2）两边的值为 （4-4）其中代表抽样信号的点。最后令，则上式变为 （4-5）首先用FFT算法求出，然后可用上式求出时的。应该强调的是，式（4-3）只是一个近似表示，计算得到的只是一个近似值。通过取更小的抽样间隔，或者增加点数，可以得到更精确的值。如果时，幅度谱很小，对应于奈奎斯特抽样频率，抽样间隔选择比较合适。如果已知信号只在时间区间内存在，可以通过对时的抽样信号补零，使足够大12。将连续时间傅立叶变换进行数字近似，用函数fft（快速傅立叶算法）高效地计算这个近似值。很多信号都能用（4-1）式连续时间傅立叶变换（CTFT）来表示。利用MATLAB可以计算（CTFT

35、）积分的数值近似。利用在密集的等间隔t的样本上的求和来近似这个积分，就可以用函数fft高效地计算这个近似值。所用的近似式是根据积分的定义得到的，即 （4-6）对于一般信号，在足够小的下，上式右边的和式是对于CTFT积分的一个好的近似。若信号对于和为零，那么这个近似式就能写成 （4-7）式中，N为一整数。可以利用函数fft对一组离散的频率计算上式中的和式。如果N个样本是存在向量x内的话，那么调用函数X=tau*fft(x)就可以计算出 （4-8）式中 以及N假设为偶数。为了计算高效，fft在负的频率样本之前先产生正频率样本。为了将频率样本置于上升的顺序，能用函数fftshift。为了将存入X中的

36、的样本排列成使就是对于，在上求得的CTFT,可用X=fftshift(tau*fft(x)。例4-1 利用FFT计算傅里叶变换如图4-1所示的信号其傅里叶变换为：利用下面的命令，可得到的近似值和准确值。 图4-1 连续时间信号x(t)N=input(Input N:);T=input(Input T:);%计算X(w)近似值t=0:T:2;x=t-1 zeros(1,N-length(t);X=fft(x);gamma=2*pi/(N*T);k=0:10/gamma;Xapp=(1-exp(-i*k*gamma*T)/(i*k*gamma)*X;%计算真实值X(w)w=0.05:0.05:10

37、;Xact=exp(-i*w)*2*i.*(w.*cos(w)-sin(w)./(w.*w);plot(k*gamma,abs(Xapp(1:length(k),o,w,abs(Xact);legend(近似值,真实值);xlabel(频率(rad/s);ylabel(|X|)运行程序后输入N=128，T=0.1，此时，得到实际的和近似的傅里叶变换的幅度谱如图4-2所示，此时近似值已经相当准确。通过增加NT可以增加更多的细节，减少T使得到的值更精确。再次运行程序后输入N=512，T=0.05，此时，得到实际的和近似的傅里叶变换的幅度谱如图4-3所示。图4-2 N=128，T=0.1时的幅度谱

38、图4-3 N=512，T=0.05时的幅度谱4.2 FFT计算离散信号的线性卷积已知两个离散时间信号与，取，对和右端补零，使得 （4-9）利用FFT算法可以求得和的L点DFT，分别是和，利用DTFT卷积性质，卷积等于乘积的L点DFT反变换，这也可以通过FFT算法得到。例4-2 利用FFT计算线性卷积已知，其中为单位阶跃序列，信号如图4-4所示。由于当时，很小，故可以取为17；N取10，。利用下面的Matlab命令，可得到、的卷积图形如图4-4所示。subplot(3,1,1);n=0:16;x=0.8.n;stem(n,x);xlabel(n);ylabel(xn);subplot(3,1,2

39、);n=0:15;y=ones(1,10) zeros(1,6); stem(n,y);xlabel(n);ylabel(yn)subplot(3,1,3);L=26;n=0:L-1;X=fft(x,L);Y=fft(y,L);Z=X.*Y;z=ifft(Z,L);stem(n,z);xlabel(n);ylabel(zn)图4-4 信号xn、yn及其卷积zn=xn*yn利用下面的Matlab命令，可得到信号xn、yn的幅度谱与相位谱如图4-5所示。subplot(2,2,1);L=26;k=0:L-1;n=0:16;x=0.8.n;X=fft(x,L);stem(k,abs(X);axis(

40、0 25 0 5);xlabel(k);ylabel(|Xk|)subplot(2,2,2);stem(k,angle(X);axis(0 25 -1 1);xlabel(k);ylabel(Angle(Xk)(弧度)subplot(2,2,3);y=ones(1,10);Y=fft(y,L);stem(k,abs(Y);axis(0 25 0 10);xlabel(k);ylabel(|Yk|)subplot(2,2,4);stem(k,angle(Y);axis(0 25 -3 3);xlabel(k);ylabel(Angle(Yk)(弧度)图4-5 信号xn、yn的幅度谱与相位谱4.3

41、 FFT进行离散信号压缩利用FFT算法对离散信号进行压缩的步骤如下：1）通过采样将信号离散化；2）对离散化信号进行傅里叶变换；3）对变换后的系数进行处理，将绝对值小于某一阈值的系数置为0，保留剩余的系数；4）利用IFFT算法对处理后的信号进行逆傅里叶变换13。 例4-3 对单位区间上的下列连续信号以采样频率进行采样，将其离散化为个采样值用FFT分解信号，对信号进行小波压缩，然后重构信号。令绝对值最小的80%系数为0，得到重构信号图形如图4-6a)所示，均方差为0.0429，相对误差为0.0449；令绝对值最小的90%系数为0，得到重构信号图形如图4-6b)所示，均方差为0.0610，相对误差为

42、0.0638。 a) 绝对值最小的80%系数为0的重构信号（FFT） b) 绝对值最小的90%系数为0的重构信号（FFT）图4-6 用FFT压缩后的重构信号相关Matlab程序如下:function wc=compress(w,r)%压缩函数compress.m%输入信号数据w,压缩率r%输出压缩后的信号数据if(r1) error(r 应该介于0和1之间!);end;N=length(w);Nr=floor(N*r);ww=sort(abs(w);tol=abs(ww(Nr+1);wc=(abs(w)=tol).*w; function unbiased_variance,error=fft

43、comp(t,y,r)%利用FFT做离散信号压缩%输入时间t,原信号y,以及压缩率rif(r1) error(r 应该介于0和1之间!);end;fy=fft(y);fyc=compress(fy,r); %调用压缩函数compress.myc=ifft(fyc);plot(t,y,r,t,yc,b);legend(原信号,重构信号);unbiased_variance=norm(y-yc)/sqrt(length(t);error=norm(y-yc)/norm(y);输入以下Matlab命令：t=(0:255)/256;f=t+cos(4*pi*t)+1/2*sin(8*pi*t);unb

44、iased_variance,error=fftcomp(t,f,0.8) unbiased_variance = 0.0429error =0.0449如果用Harr尺度函数和Harr小波分解信号，对信号进行小波压缩，然后重构信号。令绝对值最小的80%系数为0，得到重构信号图形如图4-7a)所示，均方差为0.0584，相对误差为0.0611；令绝对值最小的90%系数为0，得到重构信号图形如图4-7b)所示，均方差为0.1136，相对误差为0.1190。 a) 绝对值最小的80%系数为0的重构信号（Harr） b) 绝对值最小的90%系数为0的重构信号（Harr）图4-7 用Harr小波压缩后

45、的重构信号相关Matlab程序如下function unbiased_variance,error=daubcomp(t,y,n,r)%利用Daubechies系列小波做离散信号压缩%输入时间t,原信号y,分解层数n,以及压缩率r%输出原信号和压缩后重构信号的图像,以及重构均方差和相对l2误差if(r1) error(r应该介于0和1之间!);end;c,l=wavedec(y,n,db1);cc=compress(c,r); %调用压缩函数compress.myc=waverec(cc,l,db1);plot(t,y,r,t,yc,b);legend(原信号,重构信号);unbiased_v

46、ariance=norm(y-yc)/sqrt(length(t);error=norm(y-yc)/norm(y);输入以下Matlab命令：=(0:255)/256;f=t+cos(4*pi*t)+1/2*sin(8*pi*t);unbiased_variance,error=daubcomp(t,f,8,0.8) unbiased_variance = 0.0584error = 0.0611结论：在信号没有突变、快变化或者大致上具有周期性的信号，用FFT可以处理得很好（甚至比小波还要好）。4. 4 利用FFT对离散信号进行滤波利用FFT算法对信号进行滤波的步骤如下：1）通过采样将信号离

47、散化；2）对离散化信号进行傅里叶变换；3）对变换后的系数进行处理，将绝对值大于某一阈值的系数置为0，保留剩余的系数；4）利用IFFT算法对处理后的信号进行逆傅里叶变换14。例4-4 对被白噪声污染的信进行频谱分析，从中鉴别出有用的信号。要求：将信号的幅度换算成实际的幅度，信号的频率换算成实际的频率。相关Matlab程序如下fs=1000;%采样频率N=512;%数据点数n=0:N-1;randn(1,N);t=0:1/fs:(N-1)/fs;%采样时间序列f0=100;%信号频率x=2+3*cos(100*pi*t-pi/6)+(3/2)*cos(150*pi*t+pi/2);subplot(3,1,1);plot(t,x);xlabel(t);ylabel(sin(2+3*cos(100*pi*t-pi/6)+(3/2)*cos(150*pi*t+pi/2);title(时域信号);Y=fft(x,N);%对信号进行FFT变换magY=abs(Y);%求得FFT