相对较慢的Python NumPy三维傅里叶变换

Question

相对较慢的Python NumPy三维傅里叶变换

5

为了我的工作，我需要对大型图像执行离散傅里叶变换（DFT）。在当前示例中，我需要对一个1921 x 512 x 512的图像执行3D FT（以及512 x 512图像的2D FFT）。目前，我正在使用numpy包和相关函数np.fft.fftn()。下面的代码片段示例性地展示了在等大小/略小的2D/3D随机数生成网格上进行2D和3D FFT的时间：

import sys
import numpy as np
import time

tas = time.time()
a = np.random.rand(512, 512)
tab = time.time()
b = np.random.rand(100, 512, 512)

tbfa = time.time()

fa = np.fft.fft2(a)
tfafb = time.time()
fb = np.fft.fftn(b)
tfbe = time.time()

print "initializing 512 x 512 grid:", tab - tas
print "initializing 100 x 512 x 512 grid:", tbfa - tab
print "2D FFT on 512 x 512 grid:", tfafb - tbfa
print "3D FFT on 100 x 512 x 512 grid:", tfbe - tfafb

输出：

initializing 512 x 512 grid: 0.00305700302124
initializing 100 x 512 x 512 grid: 0.301637887955
2D FFT on 512 x 512 grid: 0.0122730731964
3D FFT on 100 x 512 x 512 grid: 3.88418793678

我遇到的问题是需要经常进行这个过程，因此每张图片的处理时间应该很短。在我的电脑上测试（中端笔记本电脑，为虚拟机分配了2GB RAM (--> 因此测试网格较小)），如您所见，3D FFT需要约5秒钟（数量级）。现在，在工作中，机器要好得多，是集群/网格架构系统，FFT速度更快。在两种情况下，2D FFT都可以几乎瞬间完成。

然而，使用1921x512x512，np.fft.fftn()需要大约5分钟。考虑到我猜测scipy的实现速度不会快多少，并且在MATLAB中同样大小的网格的FFT完成时间约为5秒，我的问题是是否有一种方法可以将此过程加速到或几乎达到MATLAB的速度。我的FFT知识有限，但显然MATLAB使用FFTW算法，而Python没有。通过一些pyFFTW包，是否有合理的机会获得类似的速度？另外，1921似乎是个不幸的选择，只有2个质因数（17、113），所以我认为这也起了一定作用。另一方面，512是一个适合的二次幂。如果不用用0填充到2048，是否可以实现类似于MATLAB的时间？我之所以问，是因为我将不得不经常使用FFT（在这种情况下，这些差异将对结果产生巨大影响！），如果在Python中无法减少计算时间，我将不得不转向其他更快的实现。

- bproxauf

如果pyfftw失败，请尝试与R或Octave的fft实现进行比较。如果其中任何一个运行更快，您可以从Python中调用这些实现（不知道惩罚有多大）。 - xvan

2个回答

0

你可以尝试使用英特尔MKL（数学核心库）的FFT，它比FFTW更快。英特尔为Python提供了mkl-fft，可以替换numpy.fft。你需要做的就是输入以下命令：

pip install mkl-fft

然后，再次运行您的程序，无需进行任何更改。

此外，numpy 1.17（即将发布）将具有新的FFT实现：

用pocketfft库替换基于fftpack的FFT模块

两种实现都有相同的祖先（Paul N. Swarztrauber的Fortran77 FFTPACK），但pocketfft包含了额外的修改，可以在某些情况下提高精度和性能。对于包含大质因数的FFT长度，pocketfft使用Bluestein算法，保持O（N log N）的运行时间复杂度，而不是对于质数长度恶化为O（N * N）。此外，接近质数长度的实值FFT的准确性已经得到改善，并且与复值FFT相当。

- marcin

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- francis · Accepted Answer

是的，通过接口 pyfftw 使用FFTW可能会比 numpy.fft 或

 scipy.fftpack 减少计算时间。这些DFT算法的实现性能可以在基准测试中进行比较，例如此处：一些有趣的结果在Python中提高FFT性能中报告。
我建议使用以下代码进行测试：
import pyfftw
import numpy
import time
import scipy

f = pyfftw.n_byte_align_empty((127,512,512),16, dtype='complex128')
#f = pyfftw.empty_aligned((33,128,128), dtype='complex128', n=16)
f[:] = numpy.random.randn(*f.shape)

# first call requires more time for plan creation
# by default, pyfftw use FFTW_MEASURE for the plan creation, which means that many 3D dft are computed so as to choose the fastest algorithm.
fftf=pyfftw.interfaces.numpy_fft.fftn(f)

#help(pyfftw.interfaces)
tas = time.time()
fftf=pyfftw.interfaces.numpy_fft.fftn(f) # here the plan is applied, nothing else.
tas = time.time()-tas
print "3D FFT, pyfftw:", tas

f = pyfftw.n_byte_align_empty((127,512,512),16, dtype='complex128')
#f = pyfftw.empty_aligned((33,128,128), dtype='complex128', n=16)
f[:] = numpy.random.randn(*f.shape)


tas = time.time()
fftf=numpy.fft.fftn(f)
tas = time.time()-tas
print "3D FFT, numpy:", tas

tas = time.time()
fftf=scipy.fftpack.fftn(f)
tas = time.time()-tas
print "3D FFT, scipy/fftpack:", tas

# first call requires more time for plan creation
# by default, pyfftw use FFTW_MEASURE for the plan creation, which means that many 3D dft are computed so as to choose the fastest algorithm.
f = pyfftw.n_byte_align_empty((128,512,512),16, dtype='complex128')
fftf=pyfftw.interfaces.numpy_fft.fftn(f)

tas = time.time()
fftf=pyfftw.interfaces.numpy_fft.fftn(f) # here the plan is applied, nothing else.
tas = time.time()-tas
print "3D padded FFT, pyfftw:", tas


对于一个尺寸为127*512*512的图像，在我的普通电脑上，我得到了以下结果：

3D FFT, pyfftw: 3.94130897522
3D FFT, numpy: 16.0487070084
3D FFT, scipy/fftpack: 19.001199007
3D padded FFT, pyfftw: 2.55221295357

所以，pyfftw比numpy.fft和scipy.fftpack快得多。使用填充甚至更快，但计算的内容是不同的。
最后，根据文档，pyfftw在第一次运行时似乎较慢，因为它使用标志FFTW_MEASURE。只有当需要连续计算相同大小的DFT时才是好事。