1、第三章离散傅里叶变换 DFT 及其快速算法 FFT 3 1离散傅里叶变换的定义及物理意义 3 2DFT的主要性质 3 3频域采样 3 5DFT FFT 应用举例 3 4DFT的快速算法 快速傅里叶变换 FFT 本章主要讲述 傅里叶变换 建立以时间t为自变量的 信号 与以频率f为自变量的 频率函数 频谱 之间的某种变换关系 所以 时间 或 频率 取连续还是离散值 就形成各种不同形式的傅里叶变换对 傅里叶变换的几种形式 3 1离散傅里叶变换的定义及物理意义 模拟域FT LT数字域FT ZT数字域DFT 时间域t 连续 频率域 s 连续 时间域n 离散 频率域k 离散 频率域 z 连续 返回 DFT
2、的图形解释 Z变换 DTFT DFT的取值范围 四种傅立叶变换 离散傅立叶变换 DFT 实现了信号首次在频域表示的离散化 使得频域也能够用计算机进行处理 并且这种DFT变换可以有多种实用的快速算法 使信号处理在时 频域的处理和转换均可离散化和快速化 因而具有重要的理论意义和应用价值 是本课程学习的一大重点 本节主要介绍 返回 3 1 1DFT定义 3 1 2DFT与ZT FT DFS的关系 3 1 3DFT的矩阵表示 3 1 1DFT定义 设序列x n 长度为M 定义x n 的N点DFT为式中 N称为离散傅里叶变换区间长度 要求N M 为书写简单 令 因此通常将N点DFT表示为定义X k 的N
3、点离散傅里叶逆变换 IDFT 为 长度为N的离散序列 返回 回到本节 例3 1 分别计算x n 的8点 16点DFT 解 x n 的8点DFT为x n 的16点DFT为 返回 回到本节 是在频率区间上的等间隔采样 返回 回到本节 3 1 2DFT与ZT FT DFS的关系 DFT有明确的物理意义 我们可以通过比较序列的DFT FT ZT 并将DFT与周期序列的DFS联系起来 得到DFT的物理意义 DFT和FT ZT之间的关系假设序列的长度为M N M将N点DFT和FT ZT的定义重写如下 返回 回到本节 比较前面三式 得到 k 0 1 2 N 1 k 0 1 2 N 1结论 1 序列的N点DF
4、T是序列傅里叶变换在频率区间 0 2 上的N点等间隔采样 采样间隔为2 N 2 序列的N点DFT是序列的Z变换在单位圆上的N点等间隔采样 频率采样间隔为2 N 返回 回到本节 DFT与z变换 DFT与DTFT变换 序列x n 的N点DFT是x n 的Z变换在单位圆上的N点等间隔采样 X k 为x n 的傅立叶变换在区间上的N点等间隔采样 这就是DFT的物理意义 变量 周期 分辨率 返回 回到本节 DFT和DFS之间的关系 周期延拓取主值有限长序列周期序列主值区序列有限长序列周期序列主值区间序列 返回 回到本节 返回 回到本节 周期序列DFS 有限长序列的DFT 对比二者发现 是的主值区序列 条
5、件N M 返回 回到本节 DFS DFT 返回 回到本节 DFT与DFS之间的关系 有限长序列x n 的DFT变换X k 就是x n 的周期延拓序列的DFS系数的主值序列 返回 回到本节 DFS与FT之间的关系 周期延拓序列的频谱特性由其傅里叶级数的系数确定 幅度相差一个常数因子 DFT的是的主值区序列 所以x n 的DFT表示的是周期序列的频谱特性 返回 回到本节 3 1 3DFT的矩阵表示 周期序列的DFS以及有限长序列x n 的DFT如下可以发现它们右边的函数形式一样 但k的定义域不同 X k 只是的主值区序列 或者说X k 以N为周期进行周期延拓即是 用后面两式表示二者的关系 返回 回
6、到本节 式 3 1 5 3 1 8 说明了DFT和DFS之间的关系 这些关系式成立的条件是N M 即DFT的变换区间N不能小于x n 的长度M 如果该条件不满足 按照式 3 1 5 将x n 进行延拓时 中将发生时域混叠 由式 3 1 8 得到的X k 不再是x n 的DFT 这时以上讲的DFS和DFT之间的关系不再成立 3 1 7 3 1 8 返回 回到本节 也可以表示成矩阵形式式中 X是N点DFT频域序列向量 x是时域序列向量 DN称为N点DFT矩阵 定义为 3 1 12 返回 回到本节 也可以表示为矩阵形式 称为N点IDFT矩阵 定义为从式 3 1 12 和式 3 1 14 我们可以发现
7、 3 1 14 返回 回到本节 3 2DFT的主要性质 与序列的FT类似 DFT也有许多重要的性质 其中一些性质本质上与FT的相应性质相同 但是某些其他性质稍微有些差别 返回 线性性质 DFT的隐含周期性 循环移位性质 复共轭序列的DFT DFT的共轭对称性 循环卷积定理 离散巴塞伐尔定理 线性性质设有限长序列的长度分别为 a和b为常数 则式中 返回 回到本节 DFT的隐含周期性在第一节中 DFT和IDFT只定义了X k 和x n 在变换区间上的N个值 如果使DFT中k的取值域为 就会发现X k 是以N为周期的 即X k mN X k 称X k 的这一特性为DFT的隐含周期性 物理意义 X k
8、 为在区间上的N点等间隔采样 以2 为周期 X k 以N为周期 返回 回到本节 循环移位性质 有限长序列的循环移位设序列x n 的长度为M 对x n 以N N M 为周期进行周期延拓 得到定义x n 的循环移位序列为上式表示将序列x n 以N为周期进行周期延拓 再左移m个单位并取主值序列 就得到x n 的循环移位序列y n 下图中 a b c 和 d 分别描述了x n 和y n 图中M 6 N 8 m 2 返回 回到本节 序列的循环移位过程示意图 返回 回到本节 循环移位性质设序列x n 长度为M x n 的循环移位序列为 N M则复共轭序列的DFT假设用表示x n 的复共轭序列 长度为N 且
9、 则 k 0 1 2 N 1式中 同样 返回 回到本节 DFT的共轭对称性上一章介绍了序列FT的共轭对称性 DFT也有类似的共轭对称性质 但FT中的共轭对称是指对坐标原点的共轭对称 在DFT中指的是对变换区间的中心 即N 2点的共轭对称 有限长共轭对称序列和共轭反对称序列假设有限长序列满足下式 n 0 1 2 N 1则称为共轭对称序列 假设有限长序列满足下式 n 0 1 2 N 1则称其为共轭反对称序列 返回 回到本节 任一有限长序列x n 都可以用它的共轭对称分量和共轭反对称分量之和表示 即将上式中的n用N n代替 并两边取共轭 得到由上面两式得到和与原序列x n 的关系为 返回 回到本节
10、DFT的共轭对称性质假设序列x n 长度为N 其N点DFT用X k 表示 将序列x n 分成实部和虚部 相应x n 的DFT分成共轭对称和共轭反对称两部分 即式中 则 返回 回到本节 将序列x n 分成共轭对称和共轭反对称两部分 相应x n 的DFT分成实部和虚部两部分 即式中 则 返回 回到本节 实信号DFT的特点设x n 是长度为N的实序列 其N点DFT用X k 表示 我们从 的结论可知道X k 具有共轭对称性质 即如果将X k 写成极坐标形式 由共轭对称性质 说明X k 的模关于k N 2点偶对称 利用DFT的共轭对称性质可以减小实序列的DFT计算量 a 利用计算一个复序列的N点DFT
11、很容易求得两个不同的实序列的N点DFT b 实序列的2N点DFT 可以用复序列的N点DFT得到 返回 回到本节 a 设是实序列 长度均为N 用它们构成一个复序列对上式进行N点DFT 得到利用 的结论可以得到这样只计算一个N点DFT 得到X k 用上面两式容易得到两个实序列的N点DFT 3 2 18 3 2 19 返回 回到本节 b 通过复序列的N点DFT得到实序列的2N点DFT 设是一个长度为2N的实序列 首先分别用中的偶数点和奇数点形成两个长度为N的新序列 即再由构造长度为N的复序列x n 即计算x n 的N点DFT 因为均是实序列 利用式 3 2 18 和式 3 2 19 得到 最后由以得
12、到实序列v n 的2N点DFT 即V k 返回 回到本节 实序列2N点的DFT 拆分重组为N点复序列DFT 例如是实序列 长度为2N拆分重组N点DFT 返回 回到本节 循环卷积定理时域循环卷积定理是DFT中很重要的定理 具有很强的实用性 已知系统输入和系统的单位脉冲响应 计算系统的输出 以及FIR滤波器用FFT实现等 都是该定理的重要应用 1 两个有限长序列的循环卷积设序列h n 和x n 的长度分别为N和M h n 与x n 的L点循环卷积定义为式中 L称为循环卷积的长度 L max N M 为了区别线性卷积 用表示循环卷积 用表示L点循环卷积 即x n 返回 回到本节 有限长序列循环卷积的
13、矩阵形式上式中右边第一个矩阵称为x n 的L点循环矩阵 它的特点是 a 第一行是x n 的L点循环倒相 x 0 不动 后面其它反转180 放在他的后面 b 第二行是第一行向右循环移一位 c 第三行是第二行向右循环移一位 依次类推 返回 回到本节 例3 2 计算下面给出的两个长度为4的序列h n 与x n 的4点和8点循环卷积 解 按照式 3 2 21 写出h n 与x n 的4点循环卷积矩阵形式为 返回 回到本节 h n 与x n 的8点循环卷积矩阵形式为 返回 回到本节 2 DFT的时域循环卷积定理设h n 和x n 长度分别为N和M yc n 为序列h n 和x n 的L点循环卷积 即x
14、n 则式中时域循环卷积定理表明 DFT将时域循环卷积关系 变换成频域的相乘关系 用时域循环卷积定理计算两个序列循环卷积运算的方框图如下图所示 图3 2 3用DFT计算两个有限长序列L点循环卷积运算的方框图 返回 回到本节 3 DFT的频域循环卷积定理设h n 和x n 长度分别为N和M 并且 则式中 L max N M DFT 时域循环卷积频域乘积离散巴塞伐尔定理设长度为N的序列x n 的N点DFT为X k 则 返回 回到本节 3 3频域采样 时域和频域的对偶原理对时间序列x n 的连续频谱函数在频域等间隔采样 则采样得到的离散频谱对应的时域序列必然是原时间序列x n 的周期延拓序列 而且仅对
15、时域有限长序列 当满足频域采用定理时 才能由频域离散采样恢复原来的连续频谱函数 或原时间序列 时域采样频域周期延拓时域周期延拓频域采样本节讨论 频域采样定理 频率采样条件 频域内插公式 返回 频域采样与频域采样定理设任意序列x n 的Z变换为而且X z 的收敛域包含单位圆 以2 N为采样间隔 在单位圆上对X z 进行等间隔采样得到实质上 是对x n 的频谱函数的等间隔采样 因为以2 为周期 所以是以N为周期的频域序列 返回 回到本节 根据离散傅里叶级数理论 必然是一个周期序列的DFS系数 经推导 我们能够得到上式说明频域采样所对应的时域周期序列是原序列x n 的周期延拓序列 延拓周期为N 根据
16、DFT与DFS之间的关系知道 分别截取和的主值序列则和构成一对DFT 3 3 3 3 3 4 3 3 5 3 3 6 返回 回到本节 3 3 3 式表明是对X z 在单位圆上的N点等间隔采样 即对在频率区间 0 2 上的N点等间隔采样 3 3 3 式 3 3 6 式说明 对应的时域有限长序列就是原序列x n 以N为周期的周期延拓序列的主值序列 综上所述 可以总结出频域采样定理 如果原序列x n 长度为M 对在频率区间 0 2 上等间隔采样N点 得到 则仅当采样点数N M时 才能由频域采样恢复 否则将产生时域混叠失真 不能由恢复原序列x n 定理告诫我们 只有当时域序列x n 为有限长时 以适当
17、的采样间隔对其频谱函数采样 才不会丢失信息 返回 回到本节 返回 回到本节 2 频域内插公式所谓频域内插公式 就是用频域采样表示X z 和 频域内插公式是FIR数字滤波器的频率采样结构和频率采样设计法的理论依据 设序列x n 的长度为M 在Z平面单位圆上对X z 的采样点数为N 且满足频域采样定理 N M 则有 3 3 7 3 3 8 返回 回到本节 将式 3 3 8 代入式 3 3 7 得到式中 令则 3 3 9a 3 3 9b 3 3 11 所以 返回 回到本节 式 3 3 9a 和式 3 3 11 称为用表示X z 的z域内插公式 称为z域内插函数 将带入 c 式并化简 得到用表示的内插
18、公式和内插函数 式 3 3 12 和式 3 3 13 将用于FIR数字滤波器的频率采样设计法的误差分析 3 3 12 3 3 13 返回 回到本节 保证了各采样点上的值与原序列的频谱相同 采样点之间为采样值与对应点的内插公式相乘 并叠加而成 返回 回到本节 3 4DFT的快速算法 快速傅里叶变换 FFT DFT使计算机在频域处理信号成为可能 但是当N很大时 直接计算N点DFT的计算量非常大 快速傅里叶变换 FFT FastFourierTransform 可使实现DFT的运算量下降几个数量级 从而使数字信号处理的速度大大提高 工程应用成为可能 人们已经研究出多种FFT算法 它们的复杂度和运算效
19、率各不相同 本章主要介绍最基本的基2FFT算法及其编程方法 返回 3 4 1直接计算DFT的特点及减少运算量的基本途径 DFT计算量 长度为N的序列x n 的N点DFT为计算X k 的每一个值需要计算N次复数乘法和N 1次复数加法 那么计算X k 的N个值需要N2次复数乘法和 N 1 N次复数加法运算 当N 1时 N点DFT的复数乘法和复数加法运算次数均与N2成正比 返回 回到本节 减少运算量方法 长序列分为短序列 由于N点DFT的运算量随N2增长 因此 当N较大时 减少运算量的途径之一就是将N点DFT分解为几个较短的DFT进行计算 则可大大减少其运算量 旋转因子的周期性和对称性 的周期性 的
20、对称性 返回 回到本节 快速傅里叶变换就是不断地将长序列的DFT分解为短序列的DFT 并利用的周期性和对称性及其一些特殊值来减少DFT运算量的快速算法 时间域抽取 频率域抽取 基2时间抽取 Decimationintime DIT FFT算法 基2频率抽取 Decimationinfrequency DIF FFT算法 返回 回到本节 3 4 2基2DIT FFT算法 基2FFT要求DFT变换区间长度N 2M M为自然数 DIT FFT算法序列x n 的N点DFT为将上面的和式按n的奇偶性分解为 返回 回到本节 令x1 l x 2l x2 l x 2l 1 因为 所以上式可写成上式说明 按n的
21、奇偶性将x n 分解为两个N 2长的序列x1 l 和x2 l 则N点DFT可分解为两个N 2点DFT来计算 用X1 k 和X2 k 分别表示x1 l 和x2 l 的N 2点DFT 即 3 4 4 3 4 5 返回 回到本节 将式 3 4 5 和式 3 4 6 代入式 3 4 4 并利用X1 k X2 k 的隐含周期性可得到这样 就将N点DFT的计算分解为计算两个N 2点离散傅里叶变换X1 k 和X2 k 再计算式 3 4 7 3 4 6 3 4 7 返回 回到本节 蝶形图 蝶形图及运算功能 8点DFT一次时域抽取分解运算流图 返回 回到本节 8点DIT FFT运算流图 返回 回到本节 2 DI
22、T FFT的运算效率 DIT FFT与DFT所需乘法次数比较曲线 返回 回到本节 由8点DIT FFT运算流图可见 N 2M时 其DIT FFT运算流图由M级蝶形构成 每级有N 2个蝶形 因此 每级需要N 2次复数乘法运算和N次复数加法运算 M级蝶形所需复数乘法次数CM 2 和复数加法次数CA 2 分别为直接计算N点DFT的复数乘法次数为N2 复数加法次数为 N 1 N 当N 1时 N2 CM 2 1 所以N越大 DIT FFT运算效率越高 DIT FFT算法与DFT所需乘法次数与N的关系曲线见前页图示 返回 回到本节 例3 3 N 210 1024时 DIT FFT的运算效率为而当N 211
23、 2048时有 DFT的乘法次数DIT FFT的乘法次数 返回 回到本节 3 DIT FFT的运算规律及编程思想运算规律 1 原位计算 2 旋转因子的变化规律 3 蝶形运算规律 4 序列的倒序 5 编程思想及程序框图 返回 回到本节 1 原位计算观察每个蝶形的两个输入和两个输出蝶形的输出可存入原输入数据所占存储单元利用同一组存储单元存储输入 输出数据的方法 称为原位 址 计算 节约内存节省寻址的时间 返回 回到本节 2 旋转因子的变化规律N点DIT FFT运算流图中 每级都有N 2个蝶形 每个蝶形都要乘以因子 称其为旋转因子 p称为旋转因子的指数 由于各级的旋转因子和循环方式都有所不同 为了编
24、写计算程序 应先找出旋转因子与运算级数的关系 用L表示从左到右的运算级数 L 1 2 M 返回 回到本节 由8点DIT FFT运算流图可以发现 第L级共有2L 1个不同的旋转因子 N 23 8时的各级旋转因子表示如下 L 1时L 2时L 3时 对N 2M的一般情况 第L级的旋转因子为 返回 回到本节 由于所以这样 就可按上面两个式子确定第L级运算的旋转因子 实际编程序时 L为最外层循环变量 返回 回到本节 3 蝶形运算规律设序列x n 经时域抽选 倒序 后 存入数组A中 如果蝶形运算的两个输入数据相距B个点 应用原位计算 则蝶形运算可表示成如下形式 式中下标L表示第L级运算 AL J 则表示第
25、L级运算后数组元素A J 的值 即第L级蝶形的输出数据 而AL 1 J 表示第L级运算前A J 的值 即第L级蝶形的输入数据 返回 回到本节 用实数运算完成上述蝶形运算 可按下面的算法进行 设式中 下标R表示取实部 I表示取虚部 则设则 返回 回到本节 4 序列的倒序DIT FFT算法的输出X k 为自然顺序 但为了适应原位计算 其输入序列不是按x n 的自然顺序排序 这种经过M 1次偶奇抽选后的排序称为序列x n 的倒序 倒位 因此 在运算之前应先对序列x n 进行倒序 由于N 2M 因此顺序数可用M位二进制数 nM 1nM 2 n1n0 表示 M次偶奇时域抽选过程如右图所示 第一次按最低位
26、n0的0和1将x n 分解为偶奇两组 第二次又按次低位n1的0 1值分别对偶奇组分解 依次类推 第M次按nM 1位分解 最后得到二进制倒序数 形成例序的树状图 N 23 返回 回到本节 不按顺序排列 按顺序输出 返回 回到本节 5 编程思想及程序框图先从输入端 第1级 开始 逐级进行 共进行M级运算 在进行第L级运算时 依次求出B 2L 1个不同的旋转因子 每求出一个旋转因子 就计算完它对应的所有2M L个蝶形 这样 我们可用三重循环程序实现DIT FFT运算 程序框图如右程序运行后 数组A中存放的是x n 的N点DFT 即X k A k 图所示 返回 回到本节 4 IDFT的高效算法上述FF
27、T算法流图也可以用于IDFT 比较DFT和IDFT的运算公式 只要将DFT运算式中的因子改为 最后乘以1 N 就是IDFT的运算公式 所以 只要将上述的DIT FFT算法中的旋转因子改为 最后的输出再乘以1 N 就可以用来计算IDFT 如果流图的输入是X k 则输出就是x n 返回 回到本节 我们还可以直接调用FFT子程序计算IFFT 由于对上式两边同时取共轭 得到这样 可以先将X k 取共轭 然后直接调用FFT子程序 或者送入FFT专用硬件设备进行FFT运算 最后对FFT结果取共轭并乘以1 N得到序列x n 返回 回到本节 3 5DFT FFT 应用举例 DFT因为具有快速算法FFT 应用非
28、常广泛 限于篇幅 这里仅介绍DFT在线性卷积和频谱分析两方面的应用 本节主要论述 返回 3 5 1用DFT FFT 计算两个有限长序列的线性卷积 3 5 2用DFT计算有限长序列与无限长序列的线性卷积 3 5 3用DFT对序列进行谱分析 3 5 1用DFT FFT 计算两个有限长序列的线性卷积 当h n 或x n 序列较长时 直接计算线性卷积的时间会很长 满足不了实时处理的要求 可用FFT将时域卷积变换为频域相乘来计算线性卷积 达到快速实时要求 下面求出线性卷积结果和循环卷积结果相等的条件 设h n 长度为N x n 长度为M yc n 和y n 分别表示h n 与x n 的L点循环卷积和线性
29、卷积 有 3 5 1 3 5 2 返回 回到本节 在式 3 5 1 中因此将上式与式 3 5 2 对比 方括号部分就是移位iL的线性卷积 因此得到 3 5 5 3 5 3 3 5 4 返回 回到本节 式 3 5 5 说明 L点循环卷积yc n 等于线性卷积y n 以L为周期的周期延拓序列的主值序列 由于y n 长度为N M 1 所以 只有当L N M 1时 式 3 5 5 给出的周期延拓无混叠 才能使yc n y n L N M 1是循环卷积结果和线性卷积结果相等的重要条件 只要满足该条件 就可以用下图计算线性卷积 返回 回到本节 循环卷积计算线性卷积的运算量 小于直接计算线性卷积的运算量可预
30、先计算并存储 乘法的运算次数可降低 又称快速卷积法 返回 回到本节 例3 4 求 返回 回到本节 3 5 2用DFT计算有限长序列与无限长序列的线性卷积 用FFT计算有限长序列与无限长序列的线性卷积问题 h n 为某滤波器的单位脉冲响应 长度有限 输入信号x n 很长 h n 要补许多零再进行计算 计算量有很大的浪费解决方法 重叠相加法重叠保留法重叠相加法 将x n 分段 每段长度为M 然后依次计算各段与h n 的卷积 再由各段的卷积结果得到y n 返回 回到本节 假设h n 和x n 是因果序列 对x n 进行分段 每段长为M 则x n 可以表示为如下形式 的长度为L N M 1 按照图3
31、2 3用L点FFT和IFFT计算 再按式 3 5 8 求和 得到 返回 回到本节 计算步骤 1 计算并保存H k DFT h n L L N M 1 i 0 2 读入并计算 3 4 n 0 1 2 L 1 5 计算 6 i i 1 返回 2 说明 一般x n 是因果序列 假设初始条件 返回 回到本节 1 重叠相加法 由分段卷积的各段相加构成总的卷积输出 返回 回到本节 回到本节 返回 3 5 3用DFT对序列进行谱分析 本节只讨论序列的频谱分析 连续信号的谱分析在第5章介绍 理论依据 序列的N点DFT就是序列的频谱函数在频率区间上的N点等间隔采样 具体方法 1 根据频率分辨率要求确定DFT变换
32、区间长度N 在数字频率域 如果要求频率分辨率为D弧度 而N点DFT意味着频谱采样间隔为2 N 即DFT能够实现的频率分辨率为2 N 因此要求2 N D 即 另外为满足基2FFT对点数N的要求 一般取 M为正整数 2 计算x n 的N点DFT 并绘制频谱图 返回 回到本节 例3 5 假设 用DFT分析x n 的频谱 频率分辨率为0 02 要求画出幅频特性曲线和相频特性曲线 解 1 根据频率分辨率求N 因为2 N 0 02 即N 100 取N 100 2 计算x n 的N点DFT 返回 回到本节 的幅频特性和相频特性曲线如下图 a 和 b 所示 由上图可见 x n 的频谱变化很缓慢 所以用DFT对该信号进行谱分析时 频率分辨率可以再低一些 即N小一些 也可以得到正确的频谱 当N 16时 频率分辨率为 8 幅频特性和相频特性曲线如图3 5 4 c 和 d 所示 与N 100的情况相差很小 返回 回到本节