a saw-shaped IFFT in matlab
3 visualizaciones (últimos 30 días)
Mostrar comentarios más antiguos
bazrafshan88@gmail.com
el 10 de Ag. de 2015
Respondida: Bjorn Gustavsson
el 11 de Ag. de 2015
I have a problem with the inverse Fourier transform of a function that is supposed to be smooth, but in practice, it becomes saw-shaped. The function in the Fourier domain is:
f(m,n)=i*m*(heaviside(m)+heaviside(-m))/(m^2+n^2)
Applying ifft2 to this function must result in a smooth function as:
F(x,y)=x/(x^2+y^2)
where m and n are the variables in fourier domain corresponding to x and y in space domain, respectively and i is sqrt(-1).
However, doing this leads to a saw-shaped function in the x-direction.
Here is the code I have written:
clear all
N=pow2(8);
Lx=0.01;
XX=linspace(-Lx,Lx,N);
[X,Y]=meshgrid(XX,XX);
dx=XX(2)-XX(1);
% the grid in fourier domain is chosen in such a way to avoid the
% singularity at the origin
mm=2*pi/N/dx*linspace(-N/2,N/2,N);
[m,n]=meshgrid(mm);
fmn=1i*m.*(heaviside(m)+heaviside(-m))./(m.^2+n.^2);
% the numerical ifft2 of the spectrum
fxy=1/dx^2*fftshift(ifft2(ifftshift(fmn)));
% the analytical ifft2 of the spectrum
Fxy=X./(X.^2+Y.^2);
subplot(2,1,1)
mesh(X,Y,abs(fxy))
title('the output by ifft2')
subplot(2,1,2)
mesh(X,Y,abs(Fxy))
title('the expected (analytical) inverse transform')
I would appreciate any suggestion in this regard.
Mohammad
2 comentarios
Walter Roberson
el 10 de Ag. de 2015
Formally speaking your function is not defined at m=0, but if you assign any non-infinite value to heaviside(0) it becomes -i*m/(m^2+n^2) . It is not obvious to me that the 2D inverse fourier of that would be of the form you indicate.
Respuestas (1)
Bjorn Gustavsson
el 11 de Ag. de 2015
I think that the problem you have is that your frequency coordinates are defined without 0 when you calculate fmn - this means that your fmn doesn't have any DC-component, then you send that into the ifft2 function after properly doing ifftshift. However the ifft2 function will treat first row and column as the DC/0-frequency component of an FFT. So when you send ifft2 something you intend to be symmetric, ifft2 doesn't agree with you so you'll get something other than you expect. To give ifft2 something it will transform right you'd have to give it the 0-frequency components, and to do that you'll have to decide how to handle the divergence to infinity...
HTH
0 comentarios
Ver también
Categorías
Más información sobre Transforms en Help Center y File Exchange.
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!