The Fourier transform of a Gaussian is not Gaussian, but it is wrong! - Python

Question

I am trying to use the Nump fft function, however, when I provide the function with a simple Gaussian function, the fft of this Gaussian function is not Gaussian, its proximity, but its half, so that each half is at both ends of the x axis.

Gaussian function I compute y = exp (-x ^ 2)

Here is my code:

from cmath import * from numpy import multiply from numpy.fft import fft from pylab import plot, show """ Basically the standard range() function but with float support """ def frange (min_value, max_value, step): value = float(min_value) array = [] while value < float(max_value): array.append(value) value += float(step) return array N = 256.0 # number of steps y = [] x = frange(-5, 5, 10/N) # fill array y with values of the Gaussian function cache = -multiply(x, x) for i in cache: y.append(exp(i)) Y = fft(y) # plot the fft of the gausian function plot(x, abs(Y)) show() 

The result is not quite right, because the FFT of a Gaussian function must be the Gaussian function itself ...

Answer 1

np.fft.fft returns the result in the so-called "standard order": ( from documents )

If A = fft(a, n) , then A[0] contains a term of zero frequency (the average value of the signal), which is always purely real for real inputs. Then A[1:n/2] contains terms with a positive frequency, and A[n/2+1:] contains terms with a negative frequency in decreasing order of negative frequency.

The np.fft.fftshift function arranges the result in the order expected by most people (and which is good for graphing):

The np.fft.fftshift(A) routine shifts the transforms and their frequencies to place the zero-frequency components in the middle ...

So using np.fft.fftshift :

 import matplotlib.pyplot as plt import numpy as np N = 128 x = np.arange(-5, 5, 10./(2 * N)) y = np.exp(-x * x) y_fft = np.fft.fftshift(np.abs(np.fft.fft(y))) / np.sqrt(len(y)) plt.plot(x,y) plt.plot(x,y_fft) plt.show() 

Answer 2

It is displayed with the center (i.e. average) at zero ratio. This is why the right half seems to be on the left, and vice versa.

EDIT: Examine the following code:

 import scipy import scipy.signal as sig import pylab x = sig.gaussian(2048, 10) X = scipy.absolute(scipy.fft(x)) pylab.plot(x) pylab.plot(X) pylab.plot(X[range(1024, 2048)+range(0, 1024)]) 

The last line will display X , starting at the center of the vector, then wrap it at the beginning.

Answer 3

Your result is not even close to Gaussian, not even divided into two halves.

To get the expected result, you will have to position your own Gaussian center with index 0, and the result will also be set this way. Try the following code:

 from pylab import * N = 128 x = r_[arange(0, 5, 5./N), arange(-5, 0, 5./N)] y = exp(-x*x) y_fft = fft(y) / sqrt(2 * N) plot(r_[y[N:], y[:N]]) plot(r_[y_fft[N:], y_fft[:N]]) show() 

The plot commands split the arrays into two halves and change them to get a nicer image.

Answer 4

The Fourier transform is implicitly repeated endlessly, since it is a transformation of a signal that is implicitly repeated endlessly. Note that when you pass the y you want to convert, the x values ​​are not provided, so the actually converted Gaussian is the center with a median value from 0 to 256, therefore 128.

Remember also that the translation f (x) is a phase change of F (x).

Answer 5

Following from Sven Marnach's answer, a simpler version would be as follows:

 from pylab import * N = 128 x = ifftshift(arange(-5,5,5./N)) y = exp(-x*x) y_fft = fft(y) / sqrt(2 * N) plot(fftshift(y)) plot(fftshift(y_fft)) show() 

This gives a graph similar to the one above.

The key (and this seems odd to me) is that NumPy assumed that the order of the data — both in the frequency and time domains — should be set to zero. This is not what I would expect from other FFT implementations, such as the FFTW3 libraries in C.

This was slightly distorted in the answers from unutbu and Steve Tjoa above, because they take the absolute value of the FFT in front of its graphics, thereby eliminating the phase problems that arise because they do not use the "standard order" in time.