1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
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);
|
