串行FFT递归算法(蝶式递归计算原理)求傅里叶变换

串行 FFT 递归算法 蝶式递归计算原理 求傅里叶变换 摘要 FFT 即为快速傅氏变换 是离散傅氏变换的快速算法 它是根据离散傅氏变换的 奇 偶 虚 实等特性 对离散傅立叶变换的算法进行改进获得的 它对傅氏变换的 理论并没有新的发现 但是对于在计算机系统或者说数字系统中应用离散傅立叶变换 可以说是进了一大步 设 x n 为 N 项的复数序列 由 DFT 变换 任一 X m 的计算都需要 N 次复数乘法 和 N 1 次复数加法 而一次复数乘法等于四次实数乘法和两次实数加法 一次复数加 法等于两次实数加法 即使把一次复数乘法和一次复数加法定义成一次 运算 四次 实数乘法和四次实数加法 那么求出 N 项复数序列的 X m 即 N 点 DFT 变换大约就 需要 N 2 次运算 当 N 1024 点甚至更多的时候 需要 N2 1048576 次运算 在 FFT 中 利用 WN 的周期性和对称性 把一个 N 项序列 设 N 2k k 为正整数 分为两个 N 2 项 的子序列 每个 N 2 点 DFT 变换需要 N 2 2 次运算 再用 N 次运算把两个 N 2 点的 DFT 变换组合成一个 N 点的 DFT 变换 这样变换以后 总的运算次数就变成 N 2 N 2 2 N N 2 2 继续上面的例子 N 1024 时 总的运算次数就变成了 525312 次 节省了大约 50 的运算量 而如果我们将这种 一分为二 的思想不断进行下去 直到分成两两一组的 DFT 运算单元 那么 N 点的 DFT 变换就只需要 Nlog 2 N 次的运 算 N 在 1024 点时 运算量仅有 10240 次 是先前的直接算法的 1 点数越多 运算 量的节约就越大 这就是 FFT 的优越性 关键字 FFT 蝶式计算 傅里叶变换 目录 一 题目及要求 1 1 1 题目 1 二 设计算法 算法原理 1 2 1 算法原理与设计 1 2 2 设计步骤 2 三 算法描述 设计流程 4 3 1 算法描述 4 3 2 流程图 5 四 源程序代码及运行结果 7 4 1 源程序代码 7 4 2 运行结果 12 五 算法分析 优缺点 13 5 1 算法分析 13 5 2 优缺点 14 六 总结 15 七 参考文献 16 1 一 题目及要求 1 1 题目 对给定的 利用串行 FFT 递归算法 蝶式递归计算 23 27 16 8 30 74 22 21 原理 计算其傅里叶变换的结果 1 2 要求 利用串行递归与蝶式递归原理 对给定的向量求解傅里叶变换的结果 二 设计算法 算法原理 2 1 算法原理与设计 蝶式递归计算原理 令 为 n 2 次单位元根 则有 将 b 向量的偶数项 和奇数项 分别记 为 和 注意推导中反复使用 图 2 1 公式图形 2 2 ni e 2 T n bbb 220 T n bbb 131 T n bbb 1 10 2 T n bbb 1 10 2 ppsnnn 1 1 1 ln2 2 2 2 设计步骤 DFTaaaaaabbb T n T n nnn 的是因此 向量 1 11 10220 222 DFTaaaaaabbb T n T n n n nn 的是因此 向量 1 1 1 1 10131 2 2 22 对于以上的分析可画出如图 2 2 所示的离散傅里叶变换递归计算流图 图 2 3 就 是一个按此递归方法计算的 n 8 的 FFT 蝶式计算图 1 1 0 2 1 0 1 1 1 2 2 2 1 10 1 1 12 2 2 4 1 1 2 0 1 12 2 4 1 2 1 12 2 4 1 2 0 1 0 2 2 2 2 2 2 222 2 2 222 2 222 2 2 n k k k lk n l ll n l ll n l ll l ll n k k lk ll laa aa aaaaaa aa aaaaaa aaaa aaaa abb n n n n nnn n n nnn n nnn n n 偶数时 1 1 0 2 1 0 1 1 11 2 2 22 1 10 1 1 112 2 2 24 1 1 2 0 1 112 1 2 1 112 2 24 1 2 0 1 12 1 1 12 1 12 1 12 1 2 12 2 1 12 0 1 0 12 12 2 2 2 22 222 2 22 222 22 22 2 22 2 2 2 2 2 2 n k k k klk n l ll n l ll n l l l ll n ln ll l ll n k k kl ll laa aaaaaaaa aaaaaaaa aaa aaaa aaa aaaa abb n n n nn nnn n nn nnn nn nn n nn n n n n n n 奇数时 3 FFT 的蝶式递归计算图 有计算原理推出 图 2 2 递归计算流图 特别的 的 FFT 蝶式计算图 展开的 图 2 3 蝶式计算图 按输入元素展开 前面将输出序列之元素 按其偶下标 和 k a j blj2 展开 导出 和递归计算式 按此构12 lj 1 0 2 2 n n k k k lk aa 1 0 2 2 n n k k k klk aa 造出了如图 1 所示的 FFT 递归计算流程图 事实上 我们也可以将输入序列之元素 a 0 a 1 a n 1 DFT n 2 DFT n 2 1 a a n 2 1 n 2 1n 1 a a n 2 1 1 aa n 2 0 aan 2 1n 1 aa n 2 0 aa n 2 1 1 an 2 1 an 2 an 2 1 n 2 1 4 按其偶下标 和几下标 展开 则导出另一种形式的 FFT 递归计算式 k ak212 k 1 0 n k k kj j awb 三 算法描述 设计流程 3 1 算法描述 SISD 上的 FFT 分治递归算法 输入 a a0 a1 an 1 输出 B b0 b1 bn 1 Procedure RFFT a b begin if n 1 then b0 a0 else 1 RFFT a0 a2 an 2 u0 u1 un 2 1 2 RFFT a1 a3 an 1 v0 v1 vn 2 1 3 z 1 4 for j 0 to n 1 do 4 1 bj uj mod n 2 zvj mod n 2 4 2 z z endfor endif end 注 1 算法时间复杂度 t n 2t n 2 O n t n O nlogn 5 n 8 的 FFT 蝶式计算图 图 3 1 FFT 蝶式计算图 n 6 的 FFT 递归计算流程图 图 3 2 FFT 递归计算流程图 3 2 流程图 6 飞飞 是 否 是 否 是 否 输入序列长度 size x 输入序列对应值 例 如 5 j3 输入 5 3 计算出前 size x 2 个 exp j 2 k size x 个值 即 W 的值 级数 i 输出 fft 结果序列 开始 结束 计算出该级需要的 W 的个数 l 该级该组起始下标 j 级数 i 加 1 组起始下标加 2 l 该级该组元素序数 k X j k X j k l X j k l W size x 2 l k X j k l 1 K 加 1 7 四 源程序代码及运行结果 4 1 源程序代码 FFT include 整个程序输入和输出利用同一个空间 x N 节约空间 include include define N 1000 定义输入或者输出空间的最大长度 typedefstruct double real doubleimg complex 定义复数型变量的结构体 void fft 快速傅里叶变换函数声明 void initW 计算 W 0 W size x 1 的值函数声明 void change 码元位置倒置函数函数声明 void add complex complex complex 复数加法 void mul complex complex complex 复数乘法 void sub complex complex complex 复数减法 void divi complex complex complex 复数除法 void output 输出结果 complex x N W 输出序列的值 intsize x 0 输入序列的长度 只限 2 的 N 次方 double PI pi 的值 int main inti system cls PI atan 1 4 8 printf Please input the size of x n 输入序列的长度 scanf d printf Please input the data in x N such as 5 6 n for i 0 i size x i 输入序列对应的值 scanf lf lf initW 计算 W 0 W size x 1 的值 fft 利用 fft 快速算法进行 DFT 变化 output 顺序输出 size x 个 fft 的结果 return 0 进行基 2 FFT 运算 蝶形算法 这个算法的思路就是 先把计算过程分为 log size x log 2 1 级 用 i 控制级数 然后把每一级蝶形单元分组 用 j 控制组的第一个元素起始下标 最后算出某一级某一组每一个蝶形单元 用 k 控制个数 共 l 个 voidfft inti 0 j 0 k 0 l 0 complexup down product change 实现对码位的倒置 for i 0 i log size x log 2 i 循环算出 fft 的结果 l 1 i for j 0 j size x j 2 l for k 0 k l k 算出第 i 级内 j 组蝶形单元的结果 算出 j 组中第 k 个蝶形单元 mul x j k l W size x 2 l k size 2 l 是该级 W 的相邻上标差 l 是该级该组取的 W 总个数 add x j k product sub x j k product 9 x j k up up 为蝶形单元右上方的值 x j k l down down 为蝶形单元右下方的值 void initW 计算 W 的实现函数 inti W complex malloc sizeof complex size x 申请 size x 个复数 W 的空间 这部申请的空间有点多 实际上只要申请 size x 2 个 即可 for i 0 i size x 2 i 预先计算出 size x 2 个 W 的值 存放 由于蝶形算法只需要前 size x 2 个值即可 W i real cos 2 PI size x i 计算 W 的实部 W i img 1 sin 2 PI size x i 计算 W 的虚部 void change 输入的码组码位倒置实现函数 complex temp unsigned short i 0 j 0 k 0 double t for i 0 i0 10 j j 1 if j i temp x i x i x j x j temp void output 输出结果实现函数 inti printf The result are as follows n for i 0 i 0 0001 如果虚部的值大于 0 0001 输出 jx img 的形式 printf j 4f n x i img else if fabs x i img real a real b real 复数实部相加 c img a img b img 复数虚部相加 void mul complex a complexb complex c 复数乘法实现函数 c real a real b real a img b img 获取相乘结果的实部 c img a real b img a img b real 获取相乘结果的虚部 void sub complex a complexb complex c 复数减法实现函数 c real a real b real 复数实部相减 c img a img b img 复数虚部相减 void divi complex a complexb complex c 复数除法实现函数 c real a real b real a img b img b real b real b img b img 获取相除结果的实部 c img a img b real a real b img b real b real b img b img 获取相除结果的虚部 12 4 2 运行结果 1 处理器 p 8 图 4 1 当时串行 FFT 输出结果 8 7 6 5 4 3 2 1 2 处理器 p 8 当时输出结果与计算结果相符如图 4 2 所示 23 27 16 8 30 74 22 21 图 4 2 运行图 13 五 算法分析 优缺点 5 1 算法分析 1 FFT 算法的基本原理是把长序列的 DFT 逐次分解为较短序列的 DFT 按照抽 取方式的不同可分为 DIT FFT 按时间抽取 和 DIF FFT 按频率抽取 算法 按照蝶 形运算的构成不同可分为基 2 基 4 基 8 以及任意因子 2n n 为大于 1 的整数 基 2 基 4 算法较为常用 2 总体结构说明 输入数据为串行的数据流 故在第一级蝶形运算模块前加入串并转换模块 将串 行数据流转换为并行的两列数据流以适应基 2 蝶形运算模块的输入信号要求 由于每级蝶形运算一次处理的两个输入数据不能直接由前一级蝶形运算一次性输 出 故在两个蝶形运算单元之间插入延时对齐模块 将前一级蝶形运算的结果 两列 并行的数据流 作适当的延时并通过转接器对齐 形成后一级蝶形运算模块所需要的 2 列输入序列 在最后一级蝶形运算后加入串并转换模块 将 2 列并行的数据流合成为 1 列 最 后加入倒序模块将 DIF FFT 得到的倒序输出序列整理为顺序输出 旋转因子产生模块产生各级基 2 蝶形运算所需的旋转因子 由运算流图可以看出 最后一级的旋转因子其实是 1 故可省略最后一级蝶形运算单元中的旋转因子乘法器 因此用一个双口 ROM 将两组数据分别输出到第一级和第二级的蝶形运算单元即可 基 2 蝶形运算模块由两个复数加法器和一个复数乘法器构成 旋转因子由 ROM 产 生后 作为复数乘法器的输入之一 与前面复数加法器得到的结果相乘完成一次蝶形 运算 为提高系统的运行速度可在蝶形运算单元中插入流水线寄存中间结果 3 蝶形运算单元如下所示 14 图 5 1 运算单元 5 2 优缺点 优点 结构清晰 可读性强 而且容易用数学归纳法来证明算法的正确性 因此 它为设计算法 调试程序带来很大方便 缺点 递归算法的运行效率较低 无论是耗费的计算时间还是占用的存储空间都 比非递归算法要多 15 六 总结 经过一周的课程设计 使我更加深