Why do scipy and numpy fft plots look different?

Question

I am currently doing some spectrum analysis for a piece of coursework, though we haven't been explicitly taught Fourier transforms yet. I have been playing around with the various fft algorithms in scipy and numpy on some data that I know what the answer should look like

In this case its an AM signal at 8kHz carrier frequency and 1kHz modulated sine wave on top so should have 3 clear peaks on an fft

When applying scipy.fftpack.rfft and numpy.fft.rfft I get the following plots respectively:

Scipy:

Numpy:

While the shape of the 2 FFTs are roughly the same with the correct ratios between the peaks, the numpy one looks much smoother, whereas the scipy one has slightly smaller max peaks, and has much more noise.

I'm assuming this is largely down to different applications of a Discrete Fourier Transform algorithm, and have seen other articles about how the scipy implementation is faster in run time. But I was wandering what it is that specifically causes the difference, and which one is actually more accurate?

EDIT: Code used to generate plots:

data = pd.read_csv("./Waveforms/AM waveform Sine.csv", sep = ',', dtype = float)

data = data.as_matrix()
time = data[:,0]
voltage = data[:,1]/data[:,1].max() # normalise the values

#scipy plot:
plt.figure()
magnitude =  scipy.fftpack.rfft(voltage)
freq = scipy.fftpack.rfftfreq(len(time),np.diff(time)[0])
plt.figure()
plt.plot(freq, np.absolute(magnitude), lw = 1)
plt.ylim(0,2500)
plt.xlim(0,15)

#numpy plot
magnitude = np.fft.rfft(voltage)
freq = np.fft.rfftfreq(len(time),np.diff(time)[0])

plt.figure()
plt.plot(freq, np.absolute(magnitude), lw = 1)
plt.ylim(0,2500)
plt.xlim(0,15)

Pierre de Buyl Pierre de Buyl · Accepted Answer · 2017-11-15T20:57:48

From NumPy's doc for rfft:

Returns:

out : complex ndarray

The truncated or zero-padded input, transformed along the axis indicated by axis, or the last one if axis is not specified. If n is even, the length of the transformed axis is (n/2)+1. If n is odd, the length is (n+1)/2.

It is not written explicitly but the "transformed data" is here complex.

From SciPy's doc for rfft

z : real ndarray

The returned real array contains:

[y(0),Re(y(1)),Im(y(1)),...,Re(y(n/2))]              if n is even
[y(0),Re(y(1)),Im(y(1)),...,Re(y(n/2)),Im(y(n/2))]   if n is odd

Conclusion: the storage is different.

For a starter, look at the length of magnitude, it will be different in both cases. I give an example below for clarity:

In [33]: data = np.random.random(size=8)

In [34]: np.fft.rfft(data)
Out[34]: 
array([ 3.33822983+0.j        ,  0.15879369+0.48542266j,
        0.00614876+0.03590621j, -0.67376592-0.69793372j,  1.51730861+0.j        ])

In [35]: scipy.fftpack.rfft(data)
Out[35]: 
array([ 3.33822983,  0.15879369,  0.48542266,  0.00614876,  0.03590621,
       -0.67376592, -0.69793372,  1.51730861])

The first element in both cases is the so-called "DC component" (the mean of the signal).

Then, you can recognize in the SciPy version the succession of real and imaginary parts of the NumPy version.

Why do scipy and numpy fft plots look different?

1 Answers