Diffraction pattern of plane wave through Circular Aperture

I have done diffraction pattern for infinite slit (in the y-dir.) with a =10 microns (total opening = 20 microns). but cannot do same thing for circular aperture with Diameter = 15 microns. Here is my code for slit aperture. Requesting for circular aperture-plane wave diffraction pattern solution.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% This program creates a plane wave function  U(r), the transmission            
% profile of a 20 um slit p(x), and a combination of the two, g(x),
% a distance away d computed from the user input Fresnel number.                                   
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

clear all

% cwd = pwd;
% cd(tempdir);
% pack
% cd(cwd)

posx=[];
posy=[];
% posy=[];
lambda=1.31; % wavelength
k=(2*pi)/lambda;
z=1; 
Nf=input('Enter Fresnel Number: ');
slit=10;
d=(slit^2)/(lambda*Nf);
%d=input('Enter distance of observation plane in microns: '); 

Ur=[];
aperture=[];

% 2-D Plane Wave Amplitude U(x)

for x=-30:1:30
    for y=-30:1:30
    posx=[posx x];
    posy=[posy y];
    u=exp(-j*k*z); % plane wave
    Ur=[Ur u];
    % begin slit definition 
    if (x>=-slit)&(x<=slit)
        amp=1;
    else
        amp=0;
    end %
    aperture=[aperture amp];
end
    
IU=abs(Ur).^2;

% Complex Wave after Aperture f(x)
fx=Ur.*aperture;
If=abs(fx).^2;

% g(x) and Intensity Determination |g(x)|^2

hx=[];
aa=[];
ho=(j/(lambda*d))*exp(-j*k*d);

for a=-70:0.025:70
    aa=[aa a];
    const=(-j*pi)/(lambda*d);
    h1=a^2;% Fresnel Appx. of transfer function of free space      
    h2=exp(const*h1);
    hx=[hx h2];
end
hhx=ho*hx;
gx=conv2(fx,hhx,'same'); % direct convolution
Ig=abs(gx).^2;

% intensity normalization
maxg1=max(Ig);

% Pupil Function p(x,y)
figure(1);
plot(posx, aperture);
title('Slit Transmission Profile');
xlabel('X-Distance (microns)');
ylabel('Transmission Value');

% Aperture Function f(x,y)
figure(2);
plot(posx,If);
title('Aperture Shadowed Beam Pattern (d=0)');
xlabel('X-Distance (microns)');
ylabel('Normalized Intensity');

% Observation Plane g(x)
figure(3);
plot(posx,Ig/maxg1);
t1=sprintf('Aperture Affected Beam Pattern (Nf=%f, d=%f microns)',Nf,d);
title(t1);
xlabel('X-Distance (microns)');
ylabel('Normalized Intensity');

Answer 1

[Bjorn Gustavsson]

Your first problem is surely that all your variables are 1-D arrays and it seems you still only model some kind of linear aperture. To model a circular aperture you will need to model a fully 2-D circular aperture. Start by making posx and posy 2-D arrays and then make aperture into a circular pupil.

[syed khaleduzzaman]

Thank you for the guidance. can you help to me to show how the code would look like for circular aperture?

[Bjorn Gustavsson]

Sure I can help you to get started.

%% Initialize the geometry of everythig, for convencience use wavelengths as length-unit
lambda=1;
%% First set up the screen where we should detect our diffraction-pattern
x_screen = linspace(-5e7,5e7,1001);
y_screen = linspace(-5e7,5e7,1001);
[x_screen,y_screen] = meshgrid(x_screen,y_screen);
z_screen = 5e7*ones(size(x_screen));
%% Insert plot of the screen
d_pupil = 50; 
x_pupil = linspace(-26,26,201);
y_pupil = linspace(-26,26,201);
[x_pupil,y_pupil] = meshgrid(x_pupil,y_pupil);
z_pupil = zeros(size(x_pupil));
iAperture = x_pupil.^2+y_pupil.^2 <= d_pupil^2/4;
%% Insert plot of the screen
E0 = z_pupil; % E-field amplitude at aperture
%% Diffraction
% the light from each point will expand as a spherical wave
% therefore the contribution dE_i2j at each point on the screen 
% from each point inside the aperture will be of the form:
%
% dE_i2j = E0_i/|r_j-r_i|*exp(-1i*(dot(k_i2j,r_j-r_i)-w*t))
%
% here the subscript i refers to points in the aperture, subscript j to
% points on the screen, r_i and r_j are the 3-D position-vectors of the
% corresponding point k_i2j are the wave-vector pointing from r_i to r_j
% and w is the wave-frequency (which we will obviously average away when
% integrating in time) This you can solve by looping over all points.

This will be a general solution that is very slow. However by being reasonably clever you can utilize the symmetries of this problem to significantly reduce the run-time. The reason I suggest you do it this way is that this brute-force approach allows you to vary the phase of the light at the aperture, or the orientation of the screen etc. This should give you a proper understanding of how the physics of diffraction work.

HTH

[syed khaleduzzaman]

Thank you !

[Bjorn Gustavsson]

You're welcome, happy that it helped.