Differences between MATLAB and Numpy/Scipy FFT

Question

EDIT: As it turns out this is still a question of floating point rounding error like others. The asymmetry in fft vs ifft absolute error comes from the difference in the magnitudes of the numbers (1e10 vs 1e8).

---

So there are many questions about the differences between Numpy/Scipy and MATLAB FFT's; however, most of these come down to floating point rounding errors and the fact that MATLAB will make elements on the order of 1e-15 into true 0's which is not what I'm after.

I am seeing a totally different issue where for identical inputs the Numpy/Scipy FFT's produce differences on the order of 1e-6 from MATLAB. At the same time for identical inputs the Numpy/Scipy IFFT's produce differences on the order or 1e-9. My data is a complex 1D vector of length 2^14 with the zero point in the middle of the array (If you know how to share this let me know). As such for both languages I am calling fftshift before and after the fft (ifft) operation.

My question is where is this difference coming from and, more importantly, why is it asymmetric with the fft and ifft? I can live with a small difference but 1e-6 is large when it accumulates over a large number of fft's.

The functional form of the fft (I'm not doing anything else to it) for either language is:

def myfft
    return fftshift(fft(fftshift(myData)))
def myifft
    return fftshift(ifft(fftshift(myData)))

I have the data saved in a .mat file and load it with scipy.io.loadmat into python. The data is a (2**14,) numpy array

The fft differences are calculated and plotted with

myData = loadmat('mydata.mat',squeeze_me=True)
plt.figure(1)
py = myfft(myData['fft_IN'])
mat = myData['fft_OUT']
plt.plot(py.real-mat.real)
plt.plot(py.imag-mat.imag)
plt.title('FFT Difference')
plt.legend(['real','imaginary'],loc=3)
plt.savefig('fft_diff')

and the ifft differences are calculated with

myData = loadmat('mydata.mat',squeeze_me=True)
plt.figure(1)
py = myifft(myData['ifft_IN'])
mat = myData['ifft_OUT']
plt.plot(py.real-mat.real)
plt.plot(py.imag-mat.imag)
plt.title('FFT Difference')
plt.legend(['real','imaginary'],loc=3)
plt.savefig('fft_diff')

Versions: Python:3.7 MATLAB:R2019a Scipy:1.4.1 Numpy:1.18.5

Answer 1

As it turns out this is still a question of floating point rounding error like all the other MATLAB vs numpy fft questions.

For my data the output of the fft function has numbers on the order of 1e10. This means that a precision of around 1e-16 on a float of this size is an absolute error less than or equal to 1e-6. The asymmetry in fft vs ifft absolute error comes from the output of the ifft being around 1e8. As such, this absolute error would then be less than or equal to 1e-8 which is exactly what we see.

Credit for this goes to @CrisLuengo who also helpfully pointed out that the ordering of fftshift and ifftshift for proper handing of odd length arrays.

Answer 2

You'll have to come up with a better workable example to show what you're after (also I don't have MATLAB, just Octave, and likely many others). I ran a quick code of fft and back with no issues. Be aware, generally DFTs (FFTs) are extremely nuanced to work with. You need to consider sampling, windowing, etc. very carefully.

Also, why the comparison to MATLAB to begin with, are you trusting it more, or just want to learn more about why one package produces an answer vs another? MATLAB uses fftw under the hood, which is very well tested and documented, but it doesn't mean that all the above nuances aren't coming into play in a different way.

import numpy as np
import matplotlib.pyplot as plt

fft = np.fft.fft
ifft = np.fft.ifft
def myfft(myData):
    return fft(myData)
def myifft(myData):
    return ifft(myData)

myData = np.exp(-np.linspace(-1, 1, 256)**2 / (2 * .25**2))
plt.figure(1)
fft_python = myifft(myfft(myData))
plt.plot(myData - fft_python.real)
plt.plot(fft_python.imag)
plt.title('FFT Difference')
plt.legend(['real','imaginary'],loc=3)
plt.savefig('fft_diff')