FFTW介绍及FFTW库的使用.pdf

1、 1 / 5 第一 部分 、FFTW 介绍 一、FFTW 介绍 FFTW 由麻 省理 工学 院计 算 机科学 实验 室超 级计 算技 术组开 发的 一套 离散 傅立 叶变换 (DFT)的 计算 库， 开源 、高效 和 标准 C 语 言编 写的 代 码使其 得到 了非 常广 泛的 应用，Intel 的 数学库 和 Scilib(类 似于 Matlab 的 科学 计算 软件) 都 使用 FFTW 做 FFT 计算 。 FFTW 是 计算 离散 Fourier 变换(DFT) 的 快速 C 程序 的 一个完 整集 合。 1 、 它可 计算 一维 或多 维、 实 和 复数 据以 及任 意规 模的 DF

2、T ； 甚至 包括 正弦/ 余弦 变 换和 离散哈 特莱 变换(DHT) 。 2 、FFTW 输 入数 据长 度任意 。 3 、FFTW 支 持任 意多 维数据 。 4 、FFTW 支持 SSE 、SSE2 、 Altivec 和 MIPS 指令 集。 5 、FFTW 还包 含对 共享 和 分布式 存储 系统 的并 行变 换。 二、FFTW 的基 本结 构： FFTW 不是 采用 固定 算法 计 算变换 ， 而是 根据 具体 硬件 和变换 参数 来调 整使 用不 同算法 ， 以 期 达到 最佳 效果。 因此， 变换 被分 成两 个阶 段。 首 先，FFTW 规划 针对 目 标计算 机的 最快

3、变换的 计算 途径 ， 并 生成 一个包 含此 信息 的数 据结 构。 然 后， 对输 入数 据进 行变换 。 该 规划 可以被 多次 使用 。 在 一个 典型的 高性 能应 用中 ， 总 是在执 行相 同参 数条 件的 任务， 因而 ， 相 对复杂 但结 果可 被重 复使 用的初 始化 是值 得的 。另 一方面 ，当 你需 要某 一参 数的单 次变换 ， 初始化 就显 得过 于费 时。 基 于此 FFTW 提供 基于 启发式 和 先例 的快 速初 始化 。 总的 来 说， FFTW 的一个 显著 特点 就是 ， 对 某一参 数类 型的 单次 变换 优势不 大， 但对 于参 数相 同的多 次变

4、 换具 有更快 的平 均速 度。 FFTW 为了 加快 用户 的使 用 集成速 度， 提供 了三 种不 同层次 的接 口。 1. 连续 数据 的单 一变 换的基 本 接口 。 2. 计算 多重 和步 进阵 列数据 的 高级 接口 。 3. 支持 通用 数据 布局 、多重 和 步进 的顶 级接 口。 大部分 的用 户使 用基 本接 口就可 以满 足需 要， 顶级 接口需 要更 小心 使用 以避 免出错 ， 因 此需要 花更 多的 时间 去理 解掌握 。FFTW 不仅 提供 了 数据自 适应 ，还 提供 了高 级用户 定制 功 能。例 如， 由于 代码 空间 不足， 可以 去掉 用户 不需 要的功

5、 能代 码。 相反 ，FFTW 还 可以 拓展 数据结 构。 第 二部 分、FFTW 的数据 类型 2.1、 数据 任何调 用FFTW 的 程序 都要 包含它 的头 文件fftw3.h 及其库(在windows 下共 享库fftw- 3.2.1.dll) 。 2.1.1 、复数 类型 FFTW 使 用包 含两 个元 素的double 型数 据来 表示 一个复 数 ， 第 一个 元素 表示 实部 ， 第二 个元素 是虚 部。 typedef double fftw_complex2； 2.1.2 、精度 FFTW 默 认使 用双 精度 浮点进 行 运算 ，用 户可 以使用float 和 long

6、 double 来改 变 运算 精度， 虽然 更高 精度 的数 据进行 运算 理论 上会 得到 更高精 度的 结果 ， 但 是由 于软件 系统 和硬 件的限 制往 往适 得其 反， 因此用 户要 结合 软硬 件 系统 来 选择 FFTW 的 数据 类型。 2 / 5 2.1.3 、分配 存储 空间 void *fftw_malloc(size_t n); void fftw_free(void *p); FFTW 一般 只是 简单 调用 系 统的函 数来 分配 存储 空间 ，通常 也不 必担 心有 重大 的内存 开 销，但 我们 强烈 建议 使用 该函数 来为 用户 的数 据分 配存储 空间

7、。 第 三部 分、FFTW 基本 接口 介绍 3.1 、基本 接口 3.1.1. 一维 复数 DFT 方案： 选择 最佳 的变 换算 法获得 变换 结果 。 使用 FFTW 计 算一 个大 小为 N 的 一维 DFT 很简 单， 步骤如 下： #include fftw_complex *in, *out; fftw_plan p; in = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * N); out = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * N); p = fftw_pl

8、an_dft_1d(N, in, out, FFTW_FORWARD, FFTW_ESTIMATE); fftw_execute(p); /* repeat as needed */ fftw_destroy_plan(p); fftw_free(in); fftw_free(out); 首先为 输入 输出 数组 分配 空间， 你可 以使 用任 何方 式进行 分配 ，但 我们 建议 使用 fftw_malloc,因之 不但 具有一 般 malloc 功能 还在 支持 SIMD 指令 的情 况下 能够 提 供数组 对齐 功能。 数据 是一 个 fftw_complex 型数 组。 分配 空间 后

9、，就 可以 创建 方案 了， 该方案 包括 了 执行变 换的 所有 需要 准备 的信息 。下 面的 函数 就是 用来创 建方 案的 ： fftw_plan fftw_plan_dft_1d(int n, fftw_complex *in, fftw_complex *out, int sign, unsigned flags); 第 一 个参数 n 是你 要进 行变 换 的数 据大 小， n 可 以是 任何大 小的 正整 数。 接下 来的两 个 参数是 输入 输出 数组 的指 针。这 两个 指针 可以 相同 ，表示 就地 转换 以节 省空 间。 第 四 个参数 sign，可 以 是 FFTW_

10、FORWARD( -1) 或 FFTW_BACKWARD (+1 ) ， 指出 变换的 方 向，从 技术 角度 讲， 它是 变换中 的指 数标 志。 3 / 5 参数 flags 通常 可以 是 FFTW_MEASURE 或 FFTW_ESTIMATE ， FFTW_MEASURE 指示 FFTW 计算几 种 FFT 算 法执 行这个 n 大 小的 变换 所花 费的 时间， 从中 找出 最优 算法 。处理 器所 花 费的时 间决 定于 硬件 性能 和 n 大小 。 FFTW_ESTIMATE 则 相反 ， 它只 是建 立一个 FFTW 认为 合 理的方 案（ 也许 不是 最优 的） 。 说白

11、了， 如果 你想 使 用同一 大小 相同 类型 的输 入数据 执行 多 次变换 ，那 么初 始化 建立 方案的 时间 被均 摊而 变得 不那么 重要 了， 选 择 FFTW_MEASURE ； 否则选 择 FFTW_ESTIMATE。在 选 择 FFTW_MEASURE 创 建 方案 的过 程中， 输入输 出 数据 的数 据将被 覆盖 ，因 此用 户应 该在初 始化 使用 变换 数据 前完成 变换 方案 的创 建。 一旦变 换方 案创 建完 成你 就可以 刷新 使用 你的 输入/ 输出数 组而 执行 多次 变换 ，实际 的 变 换 计算 工作由 fftw_execute(plan)来完 成：

12、void fftw_execute(const fftw_plan plan); 如果你 想对 相同 大小 不同 数据值 的数 组进 行变 换， 你可以 使 用 fftw_plan_dft_1d 创建 一 个新的 变换 方案 ，FFTW 将 会自动 使用 之前 的方 案数 据以节 约时 间。 完 成 变换 后调用 fftw_destroy_plan(plan)销毁 方案 ： void fftw_destroy_plan(fftw_plan plan); 使用 fftw_malloc 创建 的数组 需 要调用 fftw_free 释放 ， 而不 是使 用普 通的 内存 释 放函数 （或者 禁用

13、、删 除 ）。 DFT 结果存 放在 输出 数组 里， 零频率 直流 响应 放 在 out0 中。 如果 输 入和输 出不 是同 一数 组 ， 则 输出不 会引 起输 入数 据的 修改 （反 之变 换结 果将 覆盖 输入数 据 ） 。 用户需 要注 意的 是 FFTW 不 是 使用 通常的 DFT 算 法。 而 是计 算一 个正 向接 着一个 逆 向变 换(反之 亦然) 而得 出结 果。 可以参 考 What FFTW Really Computes 。 如果 你有 一个 支持 C99 的 编译器 如 GCC ，你 可以 包含#include ， 使用双 精度 复数 类型 来进 行运算。 另外

14、 ， 除了 FFTW_ESTIMATE 和 FFTW_MEASURE 标 志， 还有许 多 其它 的标 志存 在， 例如 ， 如果 可以 忍受更 长的 时间 ，你 可以 使用 FFTW_PATIENT 标志 来生成 一个 更有 效率 的方 案（ 见 FFTW Reference ） 。 你甚 至可 以将 方案保 存起 来留 待以 后的 应用， 见 Words of Wisdom-Saving Plans。 3.1.2. 一维 实数 DFT 在许多 应用 场合， 输入 数 据是实 数序 列， 而 输出 满 足 “Hermitian ” 冗余： outi 是 outn-i 的共轭 。可 以利 用这

15、 个特 点来改 善速 度和 内存 使用 量。 为了换 取速 度和 空间 的优 势， 用户 需要 牺牲 一些 简 单的 FFTW 复数 变换 内容。首 先 ，输 入和输 出的 数组 大小 会不 同：输 入是 大小 为 n 的实 数组， 输出 是 n/2+1 的复 数组( 去掉 了冗 余的输 出) ； 此外 还要 求输 入数据 的就 地转 化需 要填 充。 第二 点 ， 逆 变换(复数到 实 数)默认 情 况下会 破坏 输入 数组 。 这 些不便 对于 用户 来说 可能 不是个 严重 的问 题， 但是 ， 了 解这 些问 题 很有必 要。 实数的 变换 程序 和复 数几 乎是相 同的 ：分 配空

16、间( 最好 使用 fftw_malloc) ，创 建一 个 fftw_plan 并可 使用 fftw_execute(plan)多次 执行 该方案 ， 使用 fftw_destroy_plan(plan) (and fftw_free) 清理 和销 毁方 案 。唯一 的不 同是 输入 的数 据类型 为双 精度 实数 型， 并且使 用一 个 不同的 程序 来创 建方 案。 一维： fftw_plan fftw_plan_dft_r2c_1d(int n, double *in, fftw_complex *out, unsigned flags); fftw_plan fftw_plan_dft

17、_c2r_1d(int n, fftw_complex *in, double *out, unsigned flags); 实数输 入复 数输 出(r2c)和复 数 输入 实数 输出(c2r) 不同 于 复数 DFT ， 没有 sign 参数。r2c DFTs 总是 FFTW_FORWARD，而 c2r DFTs 总是 FFTW_BACKWARD.。这 里 n 是逻 辑 大小， 而不 需要数 组的 物理 大小。 需 要特别 说明 的是 实数 组 有 n 个 元素 ， 复数 组有 n/2+1 个 元 素。 对于 4 / 5 就 地 转换 来说， 输 入和 输出 数 组放 在一 起， 需要 保证

18、 数 组足 够大， 因 此实 数组 的 长度 因该是 2(n /2+1) 。 综上，c2r 变 换即 使在 非就 地转换 的情 况下 也会 破坏 输入数 组的 数据 内容 ，如 果有必 要 可以使 用 FFTW_PRESERVE_INPUT 这个 flag 来避 免， 缺 陷 是会 牺牲 一定 的性 能。 该 标志 目前 还不支 持多 维实 数 DFT 。 熟悉实 数 DFT 的读 者知 道， 输 出 复数 组的第 0 个 元素( 直流) 和第 n/2 个元 素( 当 n 为偶数 时 叫“ Nyquist”频 率 ) 完 全是实 数。 因此 一些 实现 为了 使输 入 和输 出数 据长 度一

19、致， 将 “Nyquist” 元素放 到了 直流 元素 的虚 部里。 然而 这种 实现 不能 够很好 的概 括多 维变 换， 而且节 省下 来的 空间有 限， 因此 FFTW 不 支持。 一维 r2c 和 c2r DFTs 的 另一 种 接口在 r2r 中 体现( 见:The Halfcomplex-format DFT) ， 那种 半复数 格式 的输 出与 输入 数据的 类型 和长 度相 同。 那种接 口对 于多 维实 数变 换也不 常用 ， 某 些情况 下可 能会 产生 较好 的效果 。 第 四部 分、FFTW 创建 方案( 计划) 函 数的 flags 参数 说明 4.1 、方案 的 f

20、lags 参数 FFTW 中 所有 的方 案程 序都支 持 flags 参数， 该参 数采 用位或 的方 法提 供 0 到多参 数 的配 置。 这些 参数 严格( 且 有时 效限制)的控 制方 案的 生成过 程， 并可 以加 入和 取消一 些 支持 的算 法限制 。 注意： 在创 建方 案之 后一 定要初 始话 输入 数组 ， 除 非使用 一个 之前 保存 的方 案， 因 为在 方案的 生成 过程 中， 输入 数组有 被改 写。 唯一 的例 外是 FFTW_ESTIMATE 和 FFTW_WISDOM_ONLY flag ，下面 会提 到。 在所有 情况 下， 如果 之前 保存过 相同 或者

21、更大 规模 的方案 ，可 以导 入之 。例 如在 FFTW_ESTIMATE 模式 下可以 导 入所 有的 已存 方案 ，而 在 FFTW_PATIENT 模式下 只有 patient 或 exhaustive 模式 下创 建的 方 案可 以使 用。 详见 Words of Wisdom-Saving Plans。 4.2 、方案 设计 方式 flags FTW_ESTIMATE ： 规定 根据 输入数 据长 度等 信息 直接 挑选一 个算 法( 可 能不 是最优 的) 来， 而不是 去计 算比 较各 种算 法得到 ，因 而方 案创 建的 比较快 ，且 不会 覆盖 输入 输出阵 列。 FFTW

22、_MEASURE：告 诉 FFTW 通过 几种 算法 计算 输入 阵列实 际所 用的 时间 来确 定最优 算 法。根 据你 的机 器硬 件配 置情况 ，可 能需 要点 时间( 通常为 几秒 钟) ，该 参数 是 默认选 项。 FFTW_PATIENT ：类似 FFTW_MEASURE ，但 进行 更广泛 的 算法 筛选 来确 定更 理想 的 方案 （特别 是对 于大 规模 的变 换任务 ） ， 同时 也会 带来 方 案的创 建时 间成 倍增 加。 FFTW_EXHAUSTIVE ：类似 FFTW_PATIENT ， 但进 行更为 广 泛的 算法 筛选 ，即 使是 通 常认 为并不 理想 的算

23、法也 不放 过，其 结果 相较 于后 者可 能更优 秀但 耗时 无疑 会更 长。 FFTW_WISDOM_ONLY ： 一 个特殊 的方 案设 计模 式， 该模式 的方 案仅 当任 务支 持 wisdom 时 才 会被 创建， 否则 返回 NULL 。 该 模式 通常 与其 他模式 一 起工 作， 如 “FFTW_WISDOM_ONLY |FFTW_PATIENT” 将会 在支持 wisdom 的情 况下 创建 一个 基于 FFTW_PATIENT 模式 的方案 。 在 用户需 要判 断任 务是 否支 持 wisdom 时 使用 FFTW_WISDOM_ONLY， 比 如， 用户 的 任务 在

24、不 支持 wisdom 的情 况下 可以 重新创 建一 个方 案以 避免 用户的 数据 被覆 盖。 4.3 、算法 限制 flags FFTW_DESTROY_INPUT：允 许 就地 转换(out-of-place) 的变换 可以 覆盖 输入 阵列 为任意 内 容，该 模式 有时 候会 得到 更优秀 的方 案。 FFTW_PRESERVE_INPUT ： 规定就 地转 换(out-of-place) 的变换 不能 更改 输入 阵列 。 通常 情 5 / 5 况下， 特别 是 r2c 和 hc2r( 复数 到 实数) ，FFTW_DESTROY_INPUT 是默 认模 式。 其余的 情况 ，

25、将 通 过使用 FFTW_PRESERVE_INPUT 使 得变 换方 案放弃 那 些会 破坏 输入 阵列 的 算法。 对于 二 维 c2r 变换 ，方 案将 返回 NULL ，即 该任 务不 支持 FFTW_PRESERVE_INPUT 模 式。 FFTW_UNALIGNED ：禁止 算法发 布任 何对 输入 输出 阵列进 行对 齐的 请求(SIMD 禁用) 。 通常情 况下 该模 式很 少使 用， 因 为方 案会 自动 检测 失调阵 列。 唯一 使用 此模 式的是 当你 用一 个老的 方案 去对 新的 不同 长度的 阵列 进行 变换 时的 用来对 齐( 用 fftw_malloc 分配空

26、间的 数 据不必 使用 该模 式) 。 4.4 、方案 限时 extern void fftw_set_timelimit(double seconds); 该函数 在方 案中 规定 大约 最大用 时多 少秒 。 如 果参 数 seconds =FFTW_NO_TIMELIMIT( 默 认值， 负数) ， 则不 设限 时 。反之 ，创 建方 案过 程中 会对不 同的 方案 模式 挨个 测试， 直到 测 试完成 或超 时， 然后 从中 挑出 最 佳方 案。 例如 ， 使用 FFTW_PATIENT 模 式， 如果 加 上限 时， FFTW 会先 测试 FTW_ESTIMATE 模式 ，然后 FF

27、TW_MEASURE 模 式， 如果 时间 允许 最后是 FFTW_PATIENT 模式 。如 果设 定 的是 FFTW_EXHAUSTIVE 模 式，FFTW 最 终还 要测试 FFTW_EXHAUSTIVE 模式 。 注意 ， seconds 参 数只 是一 个粗略 的限 制； 实际 上， 可能由 于方 案创 建过 程被 放在了 一个 不能被 打断 的操 作中 间， 导致花 费更 多的 时间 而产 生超时 。不 过起 码方 案会 工作在 FFTW_ESTIMATE 模 式下( 相当 于 超时 时间 设定为 0) 。 第 五部 分、FFTW 使用 方案( 创建 计划) 5.1 、使用 方案

28、所有的 变换 类型 方案 信息 存放 在 fftw_plan(一个 不透明 的 指针 类型) ， 该类 型存放 着 包括 输入输 出数 据指 针的 所有 变换所 需信 息。 创建方 案不 仅仅 收集 数据 信息， 还包 括对 输入 数据 的分析 ， 从 多种 不同 特点 的算法 中选 择最适 合的 一种 或几 种并 进行数 据预 处理 ， 稍后 的执 行变换 仅按 照方 案提 供的 信息执 行变 换 即可，FFTW 将 变换 分成 方 案和执 行两 部分 对于 相同 类型方 案， 不同 输入 数据 的多次 变换 ， 可以省 略前 面的 步骤 ， 重 复调用 执行 变换 即可 ， 从 而在不 损失 变换 精度 和速 度的前 提下 令平 均变换 效率 大大 提高 。 void fftw_execute(const fftw_plan plan); 在不改 变 fftw_plan 的前 提 下，上 面的 函数 可以 多次 执行， 因而 可以 处理 连续 数据。 fftw_execute 是 FFTW 中唯一 线 程安 全的 函数 。 void fftw_destroy_plan(fftw_plan plan); 上面的 函数 用来 释放 方案 包括其 相关 的数 据。