第3章(2)快速傅里叶变换

3 3快速傅里叶变换 FFT并不是一种新的变换形式 它只是DFT的一种快速算法 并且根据对序列分解与选取方法的不同而产生了FFT的多种算法 计算DFT复数运算量 计算一个X k 值 运算量为例k 1则要进行N次复数乘法 N 1 次复数加法所以 要完成整个DFT运算 其计算量为 N N次复数相乘 N N 1 次复数加法 利用的固有对称性和周期性改善DFT的运算效率 的对称性 的周期性 因为 由此可得出 时间抽取基 2FFT算法Decimation in Time DIT 一 算法原理 设输入序列长度为N 2M M为正整数 将该序列按时间顺序的奇偶分解为越来越短的子序列 称为基2按时间抽取的FFT算法 也称为Coolkey Tukey算法 其中基数2 N 2M M为整数 若不满足这个条件 可以人为地加上若干零值 加零补长 使其达到N 2M 二 按时间抽选的基 2FFT算法 1 算法推导设序列点数N 2M M为整数 若不满足 则补零 将序列x n 按n的奇偶分成两组 N为2的整数幂的FFT算法称基 2FFT算法 则x n 的DFT 一个N点的DFT被分解为两个N 2点DFT X1 k X2 k 这两个N 2点的DFT按照 再应用W系数的周期性 求出用X1 k X2 k 表达的后半部的X k N 2 的值 再利用周期性求X k 的后半部分 求出后半部的表示式 后半部的k值所对应的X1 k X2 k 则完全重复了前半部分的k值所对应的X1 k X2 k 的值 频域中的N个点频率成分为 结论 只要求出 0 N 2 1 区间内的各个整数k值所对应的X1 k X2 k 值 即可以求出 0 N 1 整个区间内全部X k 值 这就是FFT能大量节省计算的关键 由于N 2L 因此N 2仍为偶数 可以依照上面方法进一步把每个N 2点子序列 再按输入n的奇偶分解为两个N 4点的子序列 按这种方法不断划分下去 直到最后剩下的是2点DFT 两点DFT实际上只是加减运算 2 蝶形结 即蝶式计算结构也即为蝶式信号流图上面频域中前 后半部分表示式可以用蝶形信号流图表示 X1 k X2 k 作图要素 1 左边两路为输入 2 右边两路为输出 3 中间以一个小圆表示加 减运算 右上路为相加输出 右下路为相减输出 4 如果在某一支路上信号需要进行相乘运算 则在该支路上标以箭头 将相乘的系数标在箭头旁 5 当支路上没有箭头及系数时 则该支路的传输比为1 例子 求N 23 8点FFT变换 1 先按N 8 N 2 4 做4点的DFT 先将N 8DFT分解成2个4点DFT 可知 时域上 x 0 x 2 x 4 x 6 为偶子序列x 1 x 3 x 5 x 7 为奇子序列频域上 X 0 X 3 由X k 给出X 4 X 7 由X k N 2 给出 将N 8点分解成2个4点的DFT的信号流图 X 4 X 7 同学们自已写 2 N 2 4点 N 4 2点 FFTa先将4点分解成2点的DFT 因为4点DFT还是比较麻烦 所以再继续分解 若将N 2 4点 子序列按奇 偶分解成两个N 4点 2点 子序列 即对将x1 r 和x2 r 分解成奇 偶两个N 4点 2点 点的子序列 b求2点的DFT c一个2点的DFT蝶形流图 d另一个2点的DFT蝶形流图 同理 3 将N 4 2点 DFT再分解成2个1点的DFTa求2个一点的DFT 最后剩下两点DFT 它可分解成两个一点DFT 但一点DFT就等于输入信号本身 所以两点DFT可用一个蝶形结表示 取x 0 x 4 为例 b2个1点的DFT蝶形流图 进一步简化为蝶形流图 4 一个完整N 8的按时间抽取FFT的运算流图 3按时间抽取FFT算法的特点 原位运算 in place 码位倒序规则 算法特点 原位计算 m表示第m级迭代 k j表示数据所在的行数 例 N 8FFT运算 输入 看出 用原位运算结构运算后 A 0 A 7 正好顺序存放X 0 X 7 可以直接顺序输出 码位倒序规则 我们从输入序列的序号及整序规律得到码位倒读规则 由N 8蝶形图看出 原位计算时 FFT输出的X k 的次序正好是顺序排列的 即X 0 X 7 但输入x n 都不能按自然顺序存入到存储单元中 而是按x 0 x 4 x 2 x 6 的顺序存入存储单元即为乱序输入 顺序输出 这种顺序看起来相当杂乱 然而它是有规律的 即码位倒读规则 例子 以N 8为例 01234567 000001010011100101110111 自然顺序 二进制码表示 码位倒读 码位倒置顺序 000100010110001101011111 04261537 看出 码位倒读后的顺序刚好是数据送入计算机内的顺序 倒位序 频率抽取基 2FFT算法Decimation in Frequency DIF 算法原理 设序列点数N 2M M为整数 将X k 按k的奇偶分组前 先将输入x n 按n的顺序分成前后两半 按k的奇偶将X k 分成两部分 令 则X 2r 和X 2r 1 分别是x1 n 和x2 n 的N 2点DFT 记为X1 k 和X2 k N 2仍为偶数 进一步分解 N 2N 4 同理 其中 逐级分解 直到2点DFT 当N 8时 即分解到x3 n x4 n x5 n x6 n n 0 1 DIF与DIT比较1 相同之处DIF与DIT两种算法均为原位运算 DIF与DIT运算量相同 运算量 不同之处DIF与DIT两种算法结构倒过来DIF为输入顺序 输出乱序 运算完毕再运行 二进制倒读 程序 DIT为输入乱序 输出顺序 先运行 二进制倒读 程序 再进行求DFT 蝶形结不同DIT的复数相乘出现在减法之前 DIF的复数相乘出现在减法之后 DIF与DIT比较2 快速傅里叶逆变换 IFFT 以上所讨论的FFT的运算方法同样可用于IDFT的运算 简称为IFFT 即快速付里叶反变换 从IDFT的定义出发 可以导出下列二种利用FFT来计算IFFT的方法 利用FFT计算IFFT的思路1 把FFT的时间抽取法 用于IDFT运算时 由于输入变量由时间序列x n 改成频率序列X k 原来按x n 的奇 偶次序分组的时间抽取法FFT 现在就变成了按X k 的奇偶次序抽取了 同样 频率抽取的FFT运算用于IDFT运算时 也应改变为时间抽取的IFFT 利用FFT计算IFFT的思路2 只须将频域成份一个求共轭变换 即首先将X k 的虚部乘以 1 先取X k 的共轭 得X k 然后将X k 直接送入FFT程序即可得出Nx n 最后再对运算结果取一次共轭变换 并乘以常数1 N 即可以求出IFFT变换的x n 的值 此为DFT可用FFT程序 利用FFT计算IFFT的思路3 直接调用FFT子程序计算IFFT的方法 快速傅里叶的应用 快速卷积快速相关 快速卷积 重叠相加法 快速相关 fftfiltfftifft 3 4与本章内容有关的MATLAB函数 例已知信号 求N点DFT的幅值谱和相位谱 M文件如下 N 64 fs 100 dt 1 fs n 0 N 1 x 2 sin 5 pi n dt 5 cos 18 pi n dt y fft x N mag 2 abs y N pha angle