2

I have a problem with the inverse fourier transform of a function that is supposed to be smooth, but practically it becomes saw-shaped. The function in Fourier domain is: f(m,n)=im(heavisede(m)+heaviside(-m))/(m^2+n^2). Applying the ifft2 to this function must result in a smooth function as F(x,y)=x/(x^2+y^2). However, doing this leads to a saw-shaped function in x-direction.

I would appreciate any suggestion in this regard

Here is the code:

clc
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= i *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')
Bentoy13
  • 4,886
  • 1
  • 20
  • 33
Mohammad Bazrafshan
  • 51
  • 4
  • Did you do `fftshift`/`ifftshift`? – user3528438 Aug 09 '15 at 20:14
  • Please post runnable code. What are `m`, `im` etc? – Luis Mendo Aug 10 '15 at 00:10
  • yes I did fftshift (i think in a correct way). here is the code: – Mohammad Bazrafshan Aug 10 '15 at 07:31
  • here is the code: 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 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 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') – Mohammad Bazrafshan Aug 10 '15 at 07:41
  • m and n are the variables in fourier domain corresponding to x and y in space domain, respectively. and i is sqrt(-1). – Mohammad Bazrafshan Aug 10 '15 at 07:44
  • @MohammadBazrafshan, please include that code in the question itself, by editing it. – A. Donda Aug 10 '15 at 20:22
  • @ A. Donda I just added the code to the question. Many thanks in advance for your help :-) – Mohammad Bazrafshan Aug 11 '15 at 07:37
  • Is that correct that your Fourier signal is constant phase? Moreover, for me, `m.*(heaviside(m)+heaviside(-m))` == `m`, so is it right or is there a typo? – Bentoy13 Aug 11 '15 at 07:52
  • @ Bentoy13 Actually, the only thing that matters to me is the amplitude and not the phase, that's why I am comparing the abs values. The fourier transform of this Fxy=x/(x^2+y^2) , performed by Mathematica, results in this fmn. – Mohammad Bazrafshan Aug 11 '15 at 08:21

0 Answers0