% File: http://web.mit.edu/8.13/matlab/Examples/antenna.m
% Author: Scott Sewell
% Date: 2008-June-11

% This example Matlab script calculates the Airy pattern as described in
% Section 10.2.5 "The Circular Aperture" within "OPTICS: 4th Edition" by
% Eugene Hecht.
% It further explores the effect of convolving this irradiance distribution
% with the finite size of the source (sun).
% Finally, it demonstrates some of the nifty built in grahpics routines
% like the 'surf' function.
% This script uses the built-in function BESSELJ(NU,Z). This is the Matlab function for Bessel
% functions of the first kind.  NU = 1 for our case since we want J_1.  
% Theta = Z. 
% Type 'help besselj' at the Matlab command prompt for more information.

sundiam = 0.5 *1.20% Angular diameter of sun at 21cm
lambda = .21; % 21cm radiation expressed in [m]
k = 2*pi/lambda; % wavenumber of the 21cm radiation expressed in [m-1]
a = 7/2 * .3048 ; % radius of the SRT expressed in [m] (30.48cm per foot)
thetax = -20:.1:20; % extent of region to model in [Degrees]
arg = k * a* sin(thetax*pi/180);
i = real((besselj(1, arg) ./ arg ).^2);
I = i./max(i);
I(201)=1.0;
finitesun=rectpuls(thetax,sundiam);
finitesun=double(finitesun); % If the code doesn't run then try commenting out this line, you 
                             % are probably running an old version of MatLab. The conv
                             % function was modified in MatLab version 7.8
o=conv(I,finitesun);

cut=length(o)/4;
o=o(cut:length(o)-cut);
observed=o./max(o);

figure(1);
clf;
plot(thetax,I,thetax,finitesun,thetax,observed);
grid on; 
set(gca,'FontSize',14);
xlabel('\theta (degrees)', 'fontsize', 16);
ylabel('I_0/I', 'fontsize', 16);
title('1-D SRT Irradiance Distributions (Airy pattern)', 'fontsize', 16);
legend('Diffraction Theory',['Finite Sun =',num2str(sundiam)],'Observed Irradiance');
legend;

figure(2);
clf;
[X,Y]=meshgrid(-20:.1:20);
R=sqrt(X.^2+Y.^2);
z=(besselj(1,R)./R).^2;
Z=z./max(max(z));
surf(X,Y,Z,'FaceColor','red','EdgeColor','none');
camlight left;lighting phong;
set(gca,'FontSize',14);
title({'Irradiance at SRT resulting from','Fraunhofer Diffraction at a circular aperture'},'fontsize',16);