Error in FFT?

Brief problem here, I will try to add graphs later. I can create a function that returns a gaussian, discretize it, and take the FFT with fft() or rfft() from the FFTW module. This does not return a gaussian as it should, altertating coefficients have opposite signs. If I take absolute values these appear to follow a gaussian, as I would expect.

Anyone else have this problem?

The Fourier transform of a shited gaussian is a modulated gaussian: center your input at 0 to get rid of the modulation.

1 Like

As in center it in the discretization mesh? IF so, I did before asking. Otherwise, I’m not sure what you mean because arrays don’t have any inherent notion of 0.

In particular, note that the “origin” for a discrete Fourier transform (DFT) is at the first index (with periodic boundaries so that it “wraps around”), which is not the middle of your input array. If you center the Gaussian at the middle of your input array, for an even length, then by the shift theorem it multiplies the k-th output by (-1)ᵏ, which explains the sign alternation that you observed. (The fftshift function can shift the origin from the center to the beginning of the array and back.)

(If you google “fft gaussian” you will see that this is a common confusion for people using FFT functions in many languages.)

Also, you should technically understand that an FFT computes the DFT, not a continuous Fourier transform. The DFT of a sampled Gaussian is not quite a sampled Gaussian, because of the truncated tails. This is a small effect if the width of the Gaussian is much smaller than your DFT windows, however. (To get an exact relation, you would compute the DFT of an infinite series of periodic Gaussians, one per sampling window, adding up all of the tails.)

In general, if you want to use an FFT you need to understand what it computes and how it relates to the continuous Fourier transform, Fourier series, discrete-time Fourier transform, or whatever else you want to compute. There are lots of books and other materials on this subject, e.g. Brigham’s book or Oppenheim and Schafer.

4 Likes

I am well aware of the effects of shifting, and also well aware that the DFT and CFT are not the same, I studied the subject from Brad Osgood’s excellent introductory manusrcipt on the topic. However, the CFT is a limiting case of the DFT so this cannot account for the difference. I’ve worked quite a bit with Fourier transforms in the past, but it’s all been with acoustic analysis where the only important thing was the power spectrum and I never bothered about phase… until now. I’ll try the fftshift function and see.

Thanks @stevengj

It also might be helpful to note here that DSP.jl recently acquired the zerophase::Bool keyword argument for all the window functions, so gaussian(1024, 0.1; zerophase=true) will give you a gaussian window that is centered at index 1.

Docs here (though I just realized the unicode rendering is pretty ugly in the browser. I should fix that at some point.)

2 Likes