summaryrefslogtreecommitdiff
path: root/transverse/AiryE.m
diff options
context:
space:
mode:
authorEugeniy Mikhailov <evgmik@gmail.com>2014-05-13 15:35:26 -0400
committerEugeniy Mikhailov <evgmik@gmail.com>2014-05-13 15:58:56 -0400
commitd1b62466d81b164656588ac1c55acff15ee9ea43 (patch)
treecce78f90b0768361c4a268b946848842e1e47456 /transverse/AiryE.m
downloadoptics_toolkit-d1b62466d81b164656588ac1c55acff15ee9ea43.tar.gz
optics_toolkit-d1b62466d81b164656588ac1c55acff15ee9ea43.zip
initial
The optics_toolkit code taken from http://mercury.pr.erau.edu/~greta9a1/downloads/index.html the older version is also available at mathwork web site http://www.mathworks.com/matlabcentral/fileexchange/15459-basic-paraxial-optics-toolkit
Diffstat (limited to 'transverse/AiryE.m')
-rw-r--r--transverse/AiryE.m74
1 files changed, 74 insertions, 0 deletions
diff --git a/transverse/AiryE.m b/transverse/AiryE.m
new file mode 100644
index 0000000..b03b5ed
--- /dev/null
+++ b/transverse/AiryE.m
@@ -0,0 +1,74 @@
+%---------------------------------------------------------------
+% PROGRAM: AiryE
+% AUTHOR: Andri M. Gretarsson
+% DATE: 8/7/03
+%
+%
+% SYNTAX: E=AiryE([a,w0,L,lambda <,N>],r);
+% <,...> indicates optional argument
+%
+% Returns the Field of an Airy disk as a function of radius. The waist
+% of the beam is assumed to be at the aperture.
+%
+% a = radius of aperture
+% w0 = radius of Gaussian beam at aperture
+% L = propagation distance (axially) from aperture
+% lambda = wavelength
+% N = order to which to carry the calculation
+%
+% Last updated: 8/7/03 by AMG
+%
+%---------------------------------------------------------------
+%% SYNTAX: E=AiryE([a,w0,L,lambda <,N>],r);
+%---------------------------------------------------------------
+
+function E=AiryE(params,r)
+
+a=params(1);
+w0=params(2);
+L=params(3);
+lambda=params(4);
+if length(params)>=5, Norders=params(5); else Norders=2; end
+if Norders<1, Norders=1; end
+
+%N=a^2/L/lambda;
+%q0=i*pi*w0.^2/lambda;
+%q=q0+L;
+%E=q0/q .* exp(-i*pi*N*(r./a).^2) .* ( 1 - exp(-i*pi*N-a^2/w0^2) .* besselj(0,2*pi*N*r./a) );
+
+k=2*pi/lambda;
+F=(a/w0)^2/pi;
+R=sqrt(r.^2+L^2);
+eta=pi*F/a^2-i*k/2./R;
+gamma=-k.*r./R;
+
+
+Sum1=0;
+Sum2=0;
+for m=0:Norders-1
+ Sum1=Sum1+(-gamma/2/a./eta).^m .* besselj(m,a.*gamma);
+ Sum2=Sum2+(gamma/2/a./eta).^m + besselj(1-m,a*gamma);
+end
+
+
+A = -i*k*L/2./eta./R.^2 .* exp(-i*k.*R);
+B = exp(-gamma.^2/4./eta) - exp(-eta*a^2) .* Sum1;
+C = gamma/a./eta.*exp(-gamma.^2/4./eta) - 2*exp(-eta*a^2) .* Sum2;
+
+E= sqrt(2)/pi/w0 * A .* (B + (i*a*r./R.^2).*C);
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+