离散傅里叶变换及其快速算法

1、离散傅里叶变换及其快速算法概述 傅里叶变换实现了时域到频域的转换，在连续信号和离散信号处理技术领域有广泛的应用。 为利用计算机计算傅里叶变换,对信号与频谱有如下要求：1. 信号与频谱应是离散的2. 数据长度都是有限的 本节介绍如何将傅里叶级数和傅里叶变换的分析方法应用于离散时间信号，它们是由傅里叶变换发展而来的一种变换方法。 离散傅里叶变换(DFT)和快速傅里叶变换(FFT)在理论上解决了利用计算机进行傅里叶分析的问题。离散时间序列 数字计算机只能处理有限长的数字信号。因此，必须把一个连续的变化的模拟信号转换成有限长的离散时间序列，才能由计算机来处理。这一转换称模拟信号数字化。 对x(t)进行

2、抽样，抽样间隔为， x(t)的离散值在时间t=r，写为xr。xr,r= ,-1,0,1,2,3, 叫做。周期序列的离散分析 连续周期函数x(t)，抽样间隔t=T/N0/2j/21( )eTktkTcx tdtT0j( )e 0, 1, 2,ktkkx tck 离散傅里叶级数(DFS)的性质 离散傅里叶系数 也是离散周期离散周期的，周期为N。21j01erkNNkrrcxN2 ()211jj(2)0011eerN krkNNkNNkNrrkrrcxxcNNkc21j01erkNNkrrcxN0, 1, 2,k 傅里叶变换对小结 傅里叶级数(FS)（时域：连续连续周期周期；频域：非周期非周期离散离

3、散） 傅里叶变换(FT)（时域：连续连续非周期非周期；频域：非周期非周期连续连续） 离散傅里叶级数(DFS)（时域：离散离散周期周期；频域：周期周期离散离散）j2( ) ( )( )edftX fF x tx tt1j2( )( )( )edftx tFX fX ff0/2j/21( )eTktkTcx tdtT0j( )e 0, 1, 2,ktkkx tck 21j01erkNNkrrcxN21j0erkNNrkkxc0, 1, 2,k 傅里叶变换对小结离散傅里叶变换(DFT) 工作中经常要对有限长序列有限长序列进行频谱分析,这就是我们这里要谈到的。对有限长数字序列进行分析处理，即对有限长数

4、字序列进行分析处理，即必须对信号进行截断必须对信号进行截断。截取一段长为截取一段长为T的信号的信号p 截断截断的过程是，对实测连续信号的过程是，对实测连续信号 )(tx( )x ntNT 称为样本，称为样本，称为样本长度，称为样本长度，N为采样点。为采样点。离散傅里叶变换(DFT)的定义和基本概念 设x(n)为有限长序列： 正变换： 逆变换：21j0( ) ( )( )e 01nkNNnX kDFT x nx nkN21j01( )( )( )e 01nkNNkx nIDFT X kX knNN( ) 01( )0 x nnNx nn其他值DFT与DFS的区别 如何区分DFT与DFS这两个变换

5、对？ 在意义上有差别,在形式上相同。 DFS描述的是周期离散序列xr与其频谱ck的关系，它表明时域中的周期序列得到的频谱也是周期离散的。它是严格按照傅里叶分析的概念得来的。 由前面的分析已知，有限长序列是非周期性的，故其傅里叶变换应当是连续、周期性的频谱。DFTDFT是是DFSDFS的主值序列的主值序列，只是一种借用形式，一种算法。 注意：-离散傅里叶变换关系中，有限长序列都作为周期序列的一个周期来表示，都隐含有周期性意义。rxrN0N2Nkc02NNkN01N 0( )x r)(kXDFSDFT2Nrk非周期函数的离散傅里叶变换的物理逻辑过程 (a) 原函数，(b)截断后保留部分，(c)周期

6、拓广，(d)离散采样，(e)离散傅里叶变化后的离散谱 DFT表达式： 通常记 ，则DFT简化为 上两式可写为矩阵形式2 /ejNNW10( )( ) NnkNnX kx n W101( )( ) NnkNkx nX k WN21j0( )( )enkNNnX kx n21j01( )( )enkNNkx nX kN 例：利用DFT的矩阵表达式求4点序列 的DFT 解：由N=4得( )cos(2/)x nn N2 /ejNNWj 10( )( ) NrkNrX kx r W101( )( ) NrkNkx rX k WN x(n)和X(k)的波形图如下所示采样、采样定理和混频现象采样、采样定理和

7、混频现象1、采样、采样 数字信号处理的第一步是对连续模拟信号的数字信号处理的第一步是对连续模拟信号的离散化离散化( (即采样即采样) )。离散化离散化本身一般会带来本身一般会带来频率混叠误差频率混叠误差。其原理图如下：其原理图如下： kktktt)()( )(tx)(ktx数学描述为：理想抽样就是以数学描述为：理想抽样就是以周期性冲激串周期性冲激串来对连续时间信来对连续时间信号进行抽样。号进行抽样。周期单位脉冲周期单位脉冲 kkktkttxttxtx)()()( 周期单位冲激序列： 复指数形式FS：周期单位冲激序列的傅里叶变换( )()sns ttnTj2( )esnf tnns tS1/ss

8、fTj2j2222211( )ed() edssssssTTnf tnf tTTnsnssSs tttnTtTTt tT 1 1 1 1 1 1T1T 12T12T o 积分限在-Ts/2和Ts/2之间： 周期单位冲激序列的傅里叶变换：周期单位冲激序列的傅里叶变换j2221() edsssTnf tTnsnsStnTtT()(0)( )sstnTtTtj22211( )edsssTnf tTnssSttTTj211( ) ( )e()snf tsnnssS fF s tFfnfTT周期单位冲激序列的傅里叶变换各点处的单位脉冲按各点处的单位脉冲按等间隔排列所组成的序列等间隔排列所组成的序列 ,

9、2, 1, 0 kt kt 单位脉冲序列单位脉冲序列的傅式谱的傅式谱仍为脉冲序列，但其谱值为仍为脉冲序列，但其谱值为如上图（如上图（b）。）。 n s s 如上图（如上图（a）。）。谱线间隔为谱线间隔为其原理图如下：其原理图如下： kktktt)()( )(tx)(ktx数学描述为：理想抽样就是以数学描述为：理想抽样就是以周期性冲激串周期性冲激串来对连续时间信来对连续时间信号进行抽样。号进行抽样。周期单位脉冲周期单位脉冲 kkktkttxttxtx)()()( 单位脉冲序列单位脉冲序列的傅式谱可写成的傅式谱可写成 kt ()()nssnfffnf)(ktx根据根据卷积定理卷积定理，的傅式谱的傅

10、式谱为 Xf nssnssnssnXfXfffXfnfdfXfnfdfXfnf kkktkttxttxtx)()()( （1）幅值发生了变化（）幅值发生了变化（fs倍）倍）（2）周期延拓，周期为）周期延拓，周期为T=1/ t=fs2、采样定理和频混现象采样定理和频混现象则在离散信号谱则在离散信号谱这种现象称为这种现象称为频率混叠或频混频率混叠或频混。max2sff Xf如果原信号如果原信号x(t)中包含的最高频率成分中包含的最高频率成分中相应周期的谱会出现重叠，中相应周期的谱会出现重叠，为为混叠频率或混叠频率或Nyquist频率频率。即不产生频率混淆。即不产生频率混淆现象的临界条件。现象的临界

11、条件。采样定理又可称为采样定理又可称为：如果分析信号中最高频率成分不：如果分析信号中最高频率成分不超过混叠频率，则不出现频率混叠。超过混叠频率，则不出现频率混叠。2sNff反之，如果反之，如果频率成分的两倍，则采样后离散信号频谱中不会出现频率频率成分的两倍，则采样后离散信号频谱中不会出现频率混叠。这就是混叠。这就是采样定理采样定理。max2sff即采样频率大于等于分析信号中最高即采样频率大于等于分析信号中最高如何避免频混 消除混叠的方法有两种: 1.提高采样频率提高采样频率fs，即缩小采样时间间隔。缺点缺点: （1）采样频率高，采样频率高，内存占用量和计算量内存占用量和计算量 。 （2）许多）

12、许多信号本身可能含有全频带的频率成分信号本身可能含有全频带的频率成分； （3）信号采集系统包括采样频率上限。）信号采集系统包括采样频率上限。2.采用抗混叠滤波器。采用抗混叠滤波器。可先使信号通过一个低通滤波器，滤除高于fmax /2的信号频率成分，然后再进行采样。 实际上，由于信号频率都实际上，由于信号频率都不是严格有限的，而且，实际不是严格有限的，而且，实际使用的滤波器也都不具有理想使用的滤波器也都不具有理想滤波器在截止频率处的垂直截滤波器在截止频率处的垂直截止特性，故不足以把稍高于截止特性，故不足以把稍高于截止频率的频率分量衰减掉。止频率的频率分量衰减掉。 所以实际处理时一般应使采样频率满

13、足max2.5 4.0sff抗混滤波抗混滤波(anti-aliasing)离散傅里叶变换的泄漏问题(Leakage)在实际应用中，通常将所观测的信号限制在一定的时间间隔内，也就是说，在时域对信号进行截断操作，或称作加时间窗，亦即用时间窗函数乘以信号，即由卷积定理可知，时域相乘，频域为卷积，则有有时会造成能量分散现象，称之为。( )( ) ( )x tb t x t1( )()( )X fB fX fp 从数学意义上讲，无限长连续信号的截断相当于用一高度为从数学意义上讲，无限长连续信号的截断相当于用一高度为1 1，宽度为宽度为T T的矩形窗函数的矩形窗函数w(t)乘原信号乘原信号x(t) ，则，

14、则)()()(twtxtxT)()()(fWfXfXT截断信号时域截断信号时域截断信号频域截断信号频域离散傅里叶变换的泄漏问题(Leakage)2cos()(0tfAtx)()(2)(00ffffAfX2 02 1)(TtTttTfTfTfW)sin()(2 02 )2cos()(0TtTttfAtxT)()(sin)()(sin21)(0000ffTffTffTffTATfXT余弦信号余弦信号矩形窗函数矩形窗函数频域频域乘积乘积时域时域 余弦信号被矩形窗截断形成的泄漏余弦信号被矩形窗截断形成的泄漏 对于连续周期函数，在符合采样定理的条件下，保证窗函数b(t)的时段等于被截函数的周期T的数，可

15、以保证逆变换后准确地恢复原波形，不产生泄漏。 对于随机振动信号（非周期函数），控制泄漏的方法是采用特定的窗函数特定的窗函数，以达到控制旁瓣的效果。为了减少泄漏应该尽量寻找频域中接近为了减少泄漏应该尽量寻找频域中接近 (f)的窗函数的窗函数W(f)，即，即主瓣窄旁瓣小主瓣窄旁瓣小的窗函数。的窗函数。 频域中频域中|f|1/T的部分称为的部分称为W(f)的主瓣主瓣窄，分辨率高的主瓣主瓣窄，分辨率高 其余两旁的部分即附加频率分量称为旁瓣。其余两旁的部分即附加频率分量称为旁瓣。 旁瓣低，减少泄露旁瓣低，减少泄露1/T1/T 2/T 3/T2/T常用窗函数(1) 矩形窗(Rectangular)w(t)

16、=1(2) 汉宁窗(Hanning)w(t)=1-cos(2t/T),(3) 凯塞窗(Kaiser-Bessel)w(t)=1-2.4 cos(2t/T)+0.244 cos(4t/T)-0.00305 cos(6t/T)(4) 平顶窗（FlatTop）w(t)=1-1.93 cos(2t/T)+1.29 cos(4t/T)-0.388 cos(6t/T) +0.0322 cos(8t/T)窗函数用法 矩形窗：瞬态信号、伪随机或周期随机、窗长等于周期信号整周期时 汉宁窗：纯随机 平顶窗：周期或准周期信号 力窗或指数衰减窗：锤击法测频响函数时的力信号和脉冲响应信号快速傅里叶变换(FFT) 直接利

17、用DFT进行谱分析时，存在一个突出矛盾，即当序列长度N较大时计算量大、计算时间长、数据占用内存多，难以利用DFT进行实时处理，其应用受到很大的限制。 1965年库利(Cooley)和图基(Tukey)提出了一种DFT的快速算法，这就是FFT。FFT算法使计算量大大降低，计算时间减少，特别是当序列长度N较大时，效果更为显著。 FFT并不是一种新的变换形式，它只是DFT的一种快速算法。并且根据对序列分解与选取方法的不同而产生了FFT的多种算法。 FFT在离散傅里叶反变换、线性卷积和线性相关等方面也有重要应用。 DFT表达式： 通常记 ，则DFT简化为 上两式可写为矩阵形式2 /ejNNW10( )

18、( ) NnkNnX kx n W101( )( ) NnkNkx nX k WN21j0( )( )enkNNnX kx n21j01( )( )enkNNkx nX kNDFT的计算量DFT的计算量 有限长序列x(n)的DFT为： 将DFT定义式展开成方程组写为矩阵形式2 /ejNNW10( )( ) NnkNnX kx n W101( )( ) NnkNkx nX k WN 用向量表示 用复数表示 计算一个X(k)值需要N次复数乘法和(N-1)次复数加法，那么N个X(k)共需N2次复数乘法和N(N-1)次复数加法。每次复数乘法包括4次实数乘法和2次实数加法，每次复加法包括2次实数加法，因

19、此计算N点的DFT共需要4N2次实数乘法和(2N2+2N(N-1)次实数加法，如N=2048时，计算量为419万次。Re Re Im Im Re Im Re ImWxWxjWxxWWxXWx减少运算量方法：减少运算量方法：1. 长长序列分为短序列：序列分为短序列： 由于N点DFT的运算量随N2增长，因此，当N较大时，减少运算量的途径之一就是将N点DFT分解为几个较短的DFT进行计算，则可大大减少其运算量。2. 的的周期性和对称性：周期性和对称性： 周期性：周期性： 对称性：对称性：nkNW()2NnknkWW mod()nknkNWWmod()nkN表示用N除nk之后的余数FFT算法就是不断地

20、将长序列的DFT分解为短序列的DFT，利用 的对称性和周期性，将一个大的DFT分解成一些逐次变小的DFT，来减少DFT运算量的快速算法。 分解过程遵循两条规则： 对 时 间 进 行 偶 奇 分 解 按 时 间 抽 取 的 F F T 算 法(Decimation in time)DIT-FFT算法算法 对 频 率 进 行 前 后 分 解 按 频 率 抽 取 的 F F T 算 法(Decimation in frequency)DIF-FFT算法算法nkNW(2 )( )0,1,1(21)2xrNx nrxr(2 )( )0,1,1(21)2XmNX kmXm按时间抽取的FFT算法 要求DFT变换区间长度N=2M，M为自然数。 序列x(n)的N点DFT为 将上面的和式按n的奇偶性分解为 10( )DET ( )( ) 0,1,1Nk nNNnX kx nx n WkN /2 1/2 12(21)00( )( )( )(2 )(21)k nk nNNNNk rkrNNrrX kx n Wx n Wxr WxrW按时间抽取的FFT算法 如取N=23=8，即离散时间信号为 按照规则将序列x(n)