我的问题是关于 Numpy 的 FFT 函数中使用的算法。
Numpy 的文档说它使用 Cooley-Tukey 算法。但是,您可能知道,该算法仅在点数 N 是 2 的幂时才有效。
numpy 是否填充我的输入向量 x[n] 以计算其 FFT X[k]?(我不这么认为,因为我在输出中的点数也是 N)。我如何才能真正“看到” numpy 用于其 FFT 功能的代码?
干杯!
文档说 numpy 的 FFT 基于FFTPACK。
在 FFTPACK 文档中,我发现以下内容:
子程序 rffti(n,wsave)
子程序 rffti 初始化用于 rfftf 和 rfftb 的数组 wsave。计算 n 的素数分解以及三角函数列表并将其存储在 wsave 中。
标准的 Cooley-Tukey 算法是“radix-2 with decimation in time”,它递归地将一个大小为 n 的 FFT 的计算减少为 2 个大小2*n
为 n 的 FFT,加上 n 个大小为 2 的 FFT。有一个相同的通用分解版本该算法将大小的 FFT 转换为大小m*n
为 m 的 n 个 FFT 加上大小为 n 的 m 个 FFT。FFTPACK 中的准备例程计算输入大小的素因数分解这一事实似乎表明这就是他们正在做的事情。因此,除非您选择素数的元素,或者您的元素数具有非常大的素数,否则您仍然应该获得相当不错的加速。
几年前,我在博客上写过 Cooley-Tukey 算法的radix-2和一般因式分解版本。阅读这些内容可能有助于理解 NumPy 内部发生的事情。下图是从那里拍摄的,描绘了 CT FFT:
根据我的经验,算法不会自动填充,或者至少其中一些不会。例如,对长度不长 == 2 的幂的信号运行 scipy.signal.hilbert 方法大约需要 45 秒。当我自己用零填充信号到这样的长度时,花了 100 毫秒。
YMMV,但基本上每次运行信号处理算法时都需要仔细检查。