fft - 使用 FFTW 的零填充 FFT

Question

要在频域中插入信号，可以在时域中填充零并执行 FFT。

假设给定向量 X 中的元素数为 N 并且 Y 与 X 相同但在一侧用 N 零填充。然后下面给出相同的结果。

$$\hat{x}(k)=\sum_{n=0}^{2N-1} Y(n)e^{i2\pi k n/2N},\quad k=0,...,2N-1,$$
$$\hat{x}(k)=\sum_{n=0}^{ N-1} X(n)e^{i2\pi k n/2N},\quad k=0,...,2N-1.$$

现在如果我们使用FFTW 包，第一个方程需要2N 内存空间用于输入向量，而第二个只需要N 内存空间(我不知道是否可以在现有的 FFTW 包中完成)!计算复杂度也从2N^2log(2N)降低到2N^2log(N)。每当我们进行 2D FFT 或 3D FFT 时，问题就会变得更糟。是否可以使用 FFTW 包进行第二种方法？不过，这在 MATLAB 中相当容易做到。

Answer 1

如果 x 是一个 2N 信号，在 N 之上用零填充，它的 DFT 写成:

• 如果 k 是偶数:

因此，偶数频率的系数来自x(n)的N点离散傅立叶变换。

• 如果 k 是奇数:

因此，奇数频率的系数来自x(n)exp(i*M_PI*n/N)的N点离散傅立叶变换。

因此，零填充的 2N 信号的离散傅里叶变换恢复为长度为 N 的信号的两个 DFT，并且可以使用 fftw 来计算它们。

总计算时间为 2*c*N*ln(N)，其中 c 是常数。它预计比 DFT c*2*N*ln(2*N) 的直接计算更快。请记住 ln(2*N)=ln(2)+ln(N) :随着 N 变大，与 ln 相比，直接计算的额外工作可以忽略不计(N) :技巧变得无用， 即使维度大于一。它不影响复杂性。

此外，FFTW 非常高效，如果安装正确，它会使用您 PC 的许多功能，而且在任何情况下都很难做得比这更好，即使使用了所提供的技巧也是如此。 最后，如果输入信号是真实的，你可以使用 fftw_plan fftw_plan_dft_r2c_2d : 傅立叶空间中只有一半的系数被计算和存储。

关于内存要求，如果实在内存不足，可以使用FFTW_IN_PLACE标记并使用相同的数组进行输入和输出。然而，它稍微慢一些。

上面介绍的过程可以扩展为计算用 (L-1)N 个零填充的 N 点信号的 LN 信号的 DFT:它恢复到长度为 N 的 L 个 DFT 的计算。

与 FFTW 相比，您是否有任何引用资料显示 MATLAB 如何处理和优化填充信号的 DFT？

编辑:关于 3D 案例的进一步研究:

填充的 3D 信号 x(n,m,p) 的 3D DFT 是:

如果 k_n、k_m 和 k_p 是偶数:

如果 k_n 和 k_m 是偶数且 k_p 是奇数:

...有8个案例。

因此，将大小为 NxNxN 的 3D x 填充到 2Nx2Nx2N 的 3d dft 的计算恢复到大小为 NxNxN 的 8 个 3d dft 的计算。尺寸a 3d dft是3个1d dft的组合，尺寸N的dft总数为3x8xNxN，而直接计算需要尺寸为2N的3x(2N)*(2N)dft。计算时间是 24cN^3ln(N) 对比 24cN^3ln(2N) :小的增益是可能的...同样 fftw 很快...

然而，我们不使用黑盒 3d fft，而是通过在每个方向上执行 1d dft，一次计算大小为 N 的 8 个 dft。

• 沿 N 的 1d dft:2 种情况，NxN dfts => 2cN^3ln(N)
• 沿 M 的 1d dft:2 种情况，2NxN dfts => 4cN^3ln(N)
• 沿 P 的 1d dft:2 种情况，2Nx2N dfts => 8cN^3ln(N)

因此，总计算时间预计为 14cN^3ln(N) 而不是 24cN^3ln(2N):小增益是可能的...再次fftw 很快...

此外，计算

只需要调用一次 exp :首先计算 w=exp(I*M_PI/N) 然后更新 wn=wn*w; x(n)=x(n)*wn 或在精度成为问题时使用 pow。