FFT nomalization problem ?

Values of my intensity profile are much bigger after I used fftshift(fft2()) and values should be equal. My question is : is there a normalization problem when I do the fft() of my field U2 and what to do to have the right result ?

Here is my code so far :

 clear all;
 close all;
 clc;

 %%constants declaration and parameter calculation.

 E0 = 10;                                  %maximum field amplitude.
 z = 6.8e-3;                               %distance between fiber and lens in m.
 w0 = 15e-6;                               %waist of beamlets in m.
 lambda = 0.98e-6;                         %wavelength in m.
 Nx = 3;                                   %number of fiber over x axis.
 Ny = 3;                                   %number of fiber over y axis.
 dx = .4e-3;                               %pitche over x axis in m.
 dy = .4e-3;                               %pitche over y axis in m. 

 zr = pi*w0^2/lambda;                      %Rayleigh range in m.
 k = 2*pi/lambda;                          %wavenumber in m-1.
 wz = w0*sqrt(1+(z^2)/(zr^2));             %w(z) in m.
 Rz = z*(1+(zr^2)/(z^2));                  %R(z) in m.
 phiz = atan(z/zr);                        %Gouy phase.
 R = 1.4*wz;                               %lens radius in m. 

 %%[X,Y] domain creation.

 Npts = 2000;
 Nspots = 9;

 res1 = 1e-6;
 [X1,Y1] = meshgrid((-Npts/2:Npts/2-1).*res1,(-Npts/2:Npts/2-1).*res1);
 res2 = 50e-6;
 [X2,Y2] = meshgrid((-Npts/2:Npts/2-1).*res2,(-Npts/2:Npts/2-1).*res2);

 %%aperture center vectors.

 X0 = [-dx 0 dx -dx 0 dx -dx 0 dx];
 Y0 = [-dy -dy -dy 0 0 0 dy dy dy];

 %%field/intensity calculation and plot.

 U1 = zeros(Npts,Npts,Nspots);
 U2 = zeros(Npts,Npts,Nspots);
 tmp = E0*(w0/wz);
 for j = 1:length(X0)
     U1(:,:,j) = tmp*exp(-((X1-X0(j)).^2+(Y1-Y0(j)).^2)/(wz.^2));
     U2(:,:,j) = tmp*exp(-((X2-X0(j)).^2+(Y2-Y0(j)).^2)/(wz.^2));
 end
 U1 = sum(U1,3);
 U2 = sum(U2,3);

 UCC = conj(U1);
 I = U1.*UCC;
 figure(1);
 imagesc(X1(1,:),Y1(:,1),I);
 title('Real intensity profile');
 xlabel('x in m');
 ylabel('y in m');

 %%FFT calculation and plot.

 Uf = fftshift(fft2(U2));
 UfCC = conj(Uf);
 If = Uf.*UfCC;
 figure(2);
 imagesc(X2(1,:),Y2(:,1),If);
 title('FFT intensity profile');
 xlabel('x in m');
 ylabel('y in m');

Answer 1

[Star Strider]

I rarely use fft2, but I know that with fft it is necessary to normalise it by dividing it by the length of the time-domain signal.

See if:

Uf = fftshift(fft2(U2))/numel(U2);

does what you want.

[Regis Robles]

for :

 Uf = fftshift(fft2(U2))/numel(U2);

It changes the intensity by a e-13 factor. The result should be a sort of sum of intensities from the near field pic.

Here the pics of my tests :

[Star Strider]

I have no experience with the physics you are simulating, so I have no idea what the results should look like. I thought you were having problems with the magnitude of the fft2 output, so I suggested a way to scale it.

[Regis Robles]

You thought well. I am having problem with the magnitude because a should have the same scaling before and after I used fftshift(fft2(U2));. I think I am doing wrong with fft thus I will check how FT really works. Thanks for your help.

[Star Strider]

My pleasure.