3D Fourier Transform

Question

For my work, I need to perform discrete Fourier transforms (DFTs) on large images. In the current example, I need a 3D FT for a 1921 x 512 x 512 image (along with a 2D FFT of 512 x 512 images). Right now I am using the numpy package and its related function np.fft.fftn () . The example code snippet below shows 2D and 3D FFT times on a grid with a uniform / slightly smaller 2D / 3D size with an arbitrary number as follows:

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

Conclusion:

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

, , , , . ( , 2 , (- > )), , 3D- ~ 5 ( ). , , , /grid- . 2D .

1921x512x512, np.fft.fftn() ~ 5 . , scipy , , MATLAB FFT ~ 5 , , MATLAB . , , -, MATLAB FFTW, python. , pyFFTW ? , 1921 , 2 (17, 113), , . 512 - . MATLAB 2048?

, ( , !), python, , .

Answer 1

, , FFTW pyfftw numpy.fft scipy.fftpack. DFT , : FFT Python

:

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 documentation. , DFT .

Answer 2

FFT Intel MKL (Math Kernel Library), , FFTW. Intel mkl-fft Python, numpy.fft. , , :

pip install mkl-fft

, - .

, numpy 1.17 ( ) FFT:

FFT- fftpack pocketfft

(Fortran77 FFTPACK Paul N. Swarztrauber), pocketfft , , . , , pocketfft , O (N log N) O (N * N) . , , , .