Verifying the uncertainty principle with fft for Gaussian

Hi all,

I would like to verify that i indeed perform correctly the fourier transform using fft.

In particular, i would like to explicitly get the result using the definition in:

https://mathworld.wolfram.com/FourierTransformGaussian.html

which (i think) looks exactly like the continious version of the definition of fft in the documentation

https://www.mathworks.com/help/matlab/ref/fft.html

I.e, i would like to see that the fourier transform of: exp(- pi * x^2) will be exp(-pi*k^2), like in the link above.

Here is a very simple code:

>> v=[-100:0.01:100];

>> g=exp(-(pi)*v.^2); % Sigma=sqrt(0.5/pi)

>> f=fftshift(abs(fft(g)));

So i expect that when i plot f as a function of v i will see what i had for g (up to normalization which is less important at the moment)

i mean, the same standart deviation. Though, i am clearly getting: f= exp(-pi/4 .* v^2) and not : f= exp(-pi .* v^2) . So i am missing some factor.

Anyone can help?

Thanks

Answer 1

[David Goodmanson]

Hello Ori,

The problem stems from the following, where I am using t and f rather than x and k (see comment). For an fft, once you have a time array with a certain spacing delt, then you are not free to pick the frequency array arbitrarily. For an n-point fft, the freq array spacing, fdelt, must obey the golden rule for ffts (see code below). So you can't use t (or v) for both the time and frequency arrays. It turns out that in this case with the chosen parameters, delf differs from what it should be by a factor of approximately 2, leading to the factor of 4 that you appear to be off by. Once the frequency array is fixed up, the gaussians agree.

Note the use of ifftshift instead of fftshift in one spot in the code. Ifftshift and fftshift are not the same when n is odd

Comment -- on the Mathworld site I think they are using very poor notation. In configuration space vs momentum space, the the fourier transform exponential is usually denoted as exp(ikx), where k is the wave number, 2pi / lambda . They have exp(2pi ikx) which means that now k represents 1/lambda. Admittedly there is not a well accepted variable name for 1/lambda. I tend to call it ilam for lack of something better. (Spectroscopists use inverse centimeters cm^-1, but that is unit-dependent). To avoid confusion I will use time and frequency instead, where the exponential exp(2pi i f t) is unambiguous.

delt = .01;
t = -100:delt:100;
n = numel(t);
gt = exp(-pi*t.^2);
gf = fftshift(fft(ifftshift(gt)));
delf = 1/(n*delt);   % fft golden rule: delt*delf = 1/n
n2 = (n-1)/2;
f = delf*(-n2:n2);

figure(1)
subplot(2,1,1)
plot(t,gt)
xlim([-4 4])
subplot(2,1,2)
plot(f,real(gf))    % remove small nuisance imaginary part of gf
xlim([-4,4])

[ori grossman]

thanks!