about FFT and wave packet
2 次查看(过去 30 天)
显示 更早的评论
Hi there, I am trying to verify the fourier transformation on the Gaussian in momentum and position space. The analytical Gaussian in momentum space is
g = sqrt(1/(sqrt(pi)*s))*exp(-p^2/(4*s^2)), where s is the standard deviation. If we integrate g^2, we see that g^2 is normalized to 1.
The inverse fourier transformation of g is the Gaussian in position, which is given analytically as f f = sqrt(2*s/sqrt(pi))*exp(-2s^2*x^2)
Again, if we integrate f^2, we got the normalization to 1. I try the Fourier transformation of f in mathematica and I reproduce g. Now I need to do that numerically in matlab, I got the following code
N = 2048;
dx = 0.1963;
x = dx*[-N/2:(N/2-1)];
p = 0.0156*[-N/2+1:N/2];
dp = abs(p(2)-p(1));
s = 0.2;
f = sqrt(2*s/sqrt(pi))*exp(-2*s^2.*x.^2);
g = fftshift(fft(f))/sqrt(N); % numerical solution
gt = 1/(sqrt(2)*pi^(1/4)*sqrt(s))*exp(-p.^2/(8*s^2)); % analytical expression
sum(abs(f).^2*dx) % we see |f|^2 is normalized to one
sum(abs(g).^2*dp) % |g|^2 is not normalized to one
subplot(2,1,1); plot(x, abs(f).^2);
subplot(2,1,2); plot(p, abs(g).^2); hold on; plot(p, abs(gt).^2, 'r:');
I don't understand why g is not same as gt, which is the analytical solution.
0 个评论
回答(1 个)
Chris Turnes
2015-9-2
You are using different scaling factors for the two expressions; the FFT version is scaled with dx = 0.1963 and the "analytical expression" with dp = 0.0156.
First, you'll notice that while sum(abs(g).^2*dp) is not equal to 1, if you correct the scaling factor to use dx (which appears to define the grid over which g is calculated), you actually will get the proper normalization:
sum(abs(g).^2*dx)
ans =
1.0000e+00
Second, the analytical grid is defined over -N/2+1:N/2 while the numerical grid is -N/2:(N/2-1). If you're expecting the same results, the grids should be defined over the same point range.
To get the correct normalization then, you just need to re-scale the output of the FFT to be relative to dx, not dp:
g = g*sqrt(dx/dp);
Once you do that, you'll see
sum(abs(g).^2*dp)
ans =
1
and the plots will look correct.
2 个评论
Chris Turnes
2015-9-4
I suppose the normalization that you need will depend on exactly how the problem is set up and what you're trying to compare. I think the main issue here is just making sure that the time grid of the analytical solution and the time grid of the numerical solution that you are comparing against are the same.
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Fourier Analysis and Filtering 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!