matlabfftifft

Matlab: for even real functions, FFT complex result, IFFT real result

I am testing the validity of the FFT and IFFT functions in Matlab.

I can compare the outputs of these functions to a well-known mathematical fact: the Fourier transform of an even, real function (like a Gaussian centered at 0), is another even, real function (FFT[real, 0-centered Gaussian] = real, 0-centered Gaussian). This fact should hold both for FFT and IFFT.

First I make my grid:

nx = 256; % grid total pixel count
X = 500; % grid size (um)
dx = X/nx; %  grid spacing (um)
x = linspace(-nx/2,nx/2-1,nx)*dx; % x grid (um)

df = 1/(nx*dx); % spectral grid spacing (1/um)
f = linspace(-nx/2,nx/2-1,nx)*df; % f grid (1/um)

And I make my Gaussian: 

A = 1; % magnitude (arbitrary units) 
x_fwhm = 7; % Full width at half maximum diameter (um)

x0 = x_fwhm/sqrt(2*log(2)); % 1/e^2 radius (um)
y = A*exp(-2*x.^2./(x0)^2); % Gaussian (arbitrary units)

And apply the Fourier transform, with FFT: 

yFFT = fftshift(fft(fftshift(y)));

Alternatively, with IFFT: 

yIFFT = fftshift(ifft(fftshift(y)));

Plotting the results: Shown hare are (a) Gaussian defined in space: y, (b) real part of FFT(y), (c) real part if IFFT(y), (d) imaginary part of FFT(y), and (e) imaginary part of IFFT(y).

The IFFT does a perfect job: yIFFT is a purely real Gaussian. However, FFT yields a complex number: a very small imaginary part exists. This is fine, since an error should be expected in the fourier transform algorithms, and it's negligible anyway. What confuses me is why there is no error at all in IFFT? Are the FFT and IFFT algorithms that much different?

*** Note: fftshift and ifftshift are equivalent here, since my array has an even number of elements.

Solution

  • Handling of real-valued time-domain signal is a fairly common occurrence. So much so that the ifft function has built-in handling for the corresponding symmetry that arises in the frequency domain, as described under the "Algorithm" section of the documentation:

    The ifft function tests whether the vectors in Y are conjugate symmetric. A vector v is conjugate symmetric when the ith element satisfies v(i) = conj(v([1,end:-1:2])). If the vectors in Y are conjugate symmetric, then the inverse transform computation is faster and the output is real.

    In other words, ifft constructs the imaginary part of yIFFT to be exactly 0 since it detects that your input has conjugate symmetry.

    On the other hand, even time-domain signals are relatively less common and Mathworks didn't deem it necessary to perform a similar test in the fft function. That said, you could still take advantage of the conjugate symmetry test by computing the FFT using the ifft function through

    % compute fft(x,[],dim) using ifft:
size(x,dim) * conj(ifft(conj(x),[],dim))