close all; 
clc;
M=512;
L=10e-3;              % 5 mm; side length of the incident plane (m) 
dx=L/M;              % dx (sampling interval, m) 
rph = 0.2; % random phase magnitude=0.2*pi
sigma=0.4*L;
x=-L/2:dx:L/2-dx; 
y=x; 
[xx yy] = meshgrid(x,y);
Intensity=exp((-xx.^2-yy.^2)/(2*sigma^2));
imagesc(Intensity);
colorbar
pxs=dx ; % pixel size (m)
lambda=633*10^-9; % wavelength (m)
%% Prepare field at the input plane
k=2*pi/lambda;    % rad/m-1, wave_number
rng(1); % Create seed for random number generator

%% generate intensity and phase distribution of incident light
Int=zeros(512);
Phase=zeros(512);
for ii=1:512
    for jj=1:512
          Phase(ii,jj) = 2.*pi*rph*(randn()-0.5) ;
          Int(ii,jj)=exp((-x(ii)^2-y(jj)^2)/(2*sigma^2));
    end
end

%% numerical propagation to generate defocused images

field1=sqrt(Int).*exp(1i.*Phase); 
[V,sd] = Youngslits1D(field1,pxs) ;
figure;plot(sd,V,'bo'); hold on
xlabel('double slit separation [um]'); ylabel('Fringe Visibility')

return