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
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
|
% Coating reflectivity and transmittivity.
%--------------------------------------------------------------------------
%
% Front end code for "multidiel.m" by Orfanidis.
%
% Returns the reflectivity (fractional reflected power) and transmittivity
% (fractional transmitted power) of a periodic, multilayer dielectric
% coating as a function of the coating parameters and input beam properties
% (incident wavelength, polarization, angle of incidence, etc.)
%
% To use coatingRT.m, the function multidiel.m by Orfanidis must be in your
% matlab path. The function multidiel.m is available on the mathworks website
% in the package ewa.zip from Orfanidis:
% http://www.mathworks.com/matlabcentral/fileexchange/loadFile.do?objectId=4456
%
% AUTHOR: Andri M. Gretarsson, 6/2007.
% LAST UPDATED: 7/2009 by AMG
%
% SYNTAX: [R,T,anginc] = coatingRT(nH,nL,N,lambda,lambda0,pol,...
% <nreflecs,Lprotect,nprotect,nb,anginc,LH,LL>)
%
% OR (if the first two arg's are vectors):
%
% SYNTAX: [R,T,anginc] = coatingRT(n,L,[],lambda,lambda0,pol,...
% <nreflecs,Lprotect,nprotect,nb,anginc>)
%
% where n is a vector of [incident medium index, the coating layer indices of refraction, substrate index]
% and L is a vector of length length(n)-2 containing the coating layer OPTICAL thicknesses
% (in units of lambda0). In other words the arguments are in the same format as the first two argument
% of multidiel.m
%
% INPUT VARIABLES
% ---------------
% nH = index of refraction of the high-index coating layers (e.g. Ta2O5) OR vector
% of layer indices
% nL = index of refraction of the high-index coating layers (e.g. SiO2) OR vector
% of layer OPTICAL thicknesses
% N = number of coating layer pairs
% lambda = (1xN vector) incident wavelength in nanometers. This is the wavelength
% of the light that is actually hitting the coating. Note: Units are
% nanometers.
% lambda0 = center wavelength in nm. This is the wavelength of
% light in terms of which the layer thicknesses are specified. For
% a high-reflective coating with quarter-wave layers, the center
% wavelength is usually also the wavelength of maximum reflectivity. (Hence
% the term "center wavelength.") Note: Units are nanometers.
% pol = polarization of the incident light, should be 'te' for
% s-polarization or 'tm' for p-polarization. ('te' stands for
% "transverse-electric." In other words, the electric field is
% transverse to the plane of incidence.) It's also OK to specify 'p' or
% 's' as these are translated to 'tm' and 'te' respectively.
% nreflecs = This should usually be 0. Setting this to a higher number
% causes the function to add in that number of _secondary_
% reflections from the back surface of the sample (assumed to be
% uncoated).
% Lprotect = OPTICAL thickness (in terms of lambda0) of any protective coating layer laid
% down on top of the coating. Usually, the topmost periodic layer
% is a low index layer (usually SiO2). However, the topmost layer is
% usually a diffeent thickness layer (e.g. 1/2 wavelength thick) rather than a 1/4
% wavelength layer like the others.
% This is accounted for in this code by setting Lprotect to 1/4,
% and nprotect=nL. In other words the 1/2 wavelength topmost
% layer is constructed as two 1/4 wavelength layers (the periodic
% layer and the protective layer) each with the same index of
% refraction. And these are indeed the default values for Lprotect
% and nprotect. Default: Lprotect=1/4.
% nprotect = Index of refraction of any additional coating layer laid down
% to protect the periodic coating. Default: nprotect=nL.
% nb = Index of refraction of the substrate. Default: nb=1.46 (silica).
% anginc = 1xM vector of incidence angles in radians. Default is 500
% evenly spaced angles between zero and pi/2.
% LH = OPTICAL thickness in terms of lambda0 of the high index layers if not 1/4.
% Default: LH=1/4.
% LL = OPTICAL thickness in terms of lambda0 of the low index layers if not 1/4.
% Default: LL=1/4.
%
% OUTPUT VARIABLES
% ----------------
% R, T = MxN matrix of Reflectivity and Transmissivity (fractional reflected or
% transmitted _power_). R and T are vectors of size(anginc).
% lambda_mesh = This is the lambda mesh needed to plot R or T as a surface
% anginc_mesh = This is the anginc mesh needed to plot R or T as a surface
% r = The amplitude reflectance (R = r x r*).
%
%
% EXAMPLE 1 (R as a function of incident angle and wavelength)
% ---------
% lambda = [199:3:1201]*1e-9;
% lambda0 = 1064e-9;
% anginc = [0:0.5:90]*pi/180;
% [R,T,lambda_mesh,anginc_mesh,r] = coatingRT(2.02,1.46,15,lambda,lambda0,'tm',...
% 0,0,0,1.46,anginc);
% h=pcolor(lambda_mesh*1e9,anginc_mesh*180/pi,R);
% set(h,'edgecolor','none');
% colorbar;
% caxis([0 1]);
%
%
% EXAMPLE 2 (R as a function of incident angle at lambda=900 nm)
% ---------
% lambda = 900e-9;
% lambda0 = 1064e-9;
% anginc = [0:0.1:90]*pi/180;
% [R,T,lambda,anginc,r] = coatingRT(2.02,1.46,15,lambda,lambda0,'tm',...
% 0,0,0,1.46,anginc);
% plot(anginc*180/pi,R);
%
%
% EXAMPLE 3 (R as a function of wavelength at an incident angle of 70 degrees)
% ---------
% lambda_in = [500:0.1:1300]*1e-9;
% lambda0 = 1064e-9;
% anginc_in = 70*pi/180;
% [Rp,Tp,lambda,anginc] = coatingRT(2.02,1.46,15,lambda_in,lambda0,'tm',...
% 0,0,0,1.46,anginc_in);
% [Rs,Ts,lambda,anginc,r] = coatingRT(2.02,1.46,15,lambda_in,lambda0,'te',...
% 0,0,0,1.46,anginc_in);
% plot(lambda*1e9,Rp,lambda*1e9,Rs);
% axis([500 1300 0 1.1]);
%
%
% EXAMPLE 4 (R as a function of wavelength at normal incidence with an succession of irregular coatings)
% ---------
% for s=1:25
% n=[1,repmat([2.0+0.2*rand(1,1),1.5+0.2*rand(1,1)],1,15),1.46];
% L=0.25+[0.1*rand(1,30)];
% lambda_in = [500:0.1:1300]*1e-9;
% lambda0 = 900e-9;
% anginc_in = 0;
% [Rp,Tp,lambda,anginc] = coatingRT(n,L,[],lambda_in,lambda0,'tm',5,[],[],[],anginc_in);
% plot(lambda*1e9,Rp);
% axis([500 1300 0 1.1]);
% drawnow
% pause(0.1)
% end
%
% % EXAMPLE 5 (R from a simple uncoated glass surface as a function of incident angle at lambda=632 nm)
% ---------
% lambda = 632e-9;
% lambda0 = 1064e-9;
% anginc = [0:0.1:90]*pi/180;
% [R,T,lambda,anginc,r] = coatingRT([1 1.5],[],[],lambda,lambda0,'p',...
% [],[],[],[],anginc);
% plot(anginc*180/pi,R); set(gca,'xlim',[0 90]); shg;
%
%--------------------------------------------------------------------------
% SYNTAX: [R,T,lambda_mesh,anginc_mesh,r] = coatingRT(nH,nL,N,lambda,...
% lambda0,pol <,nreflecs,Lprotect,nprotect,nb,anginc,LH,LL>)
%--------------------------------------------------------------------------
function [R,T,lambda_mesh,anginc_mesh,r] = coatingRT(nH,nL,N,lambda,lambda0,pol,varargin)
m=6; n=1;
if nargin>=m+n, nreflecs=varargin{n}; else nreflecs=0; end; n=n+1;
if nargin>=m+n, Lprotect=varargin{n}; else Lprotect=0; end; n=n+1;
if nargin>=m+n, nprotect=varargin{n}; else nprotect=nL; end; n=n+1;
if nargin>=m+n, nb=varargin{n}; else nb=1.46; end; n=n+1;
if nargin>=m+n, anginc=varargin{n}; else anginc=linspace(0,90,500)*pi/180; end; n=n+1;
if nargin>=m+n, LH=varargin{n}; else LH=1/4; end; n=n+1; % OPTICAL thicknesses of the layers in units of lambda0,
if nargin>=m+n, LL=varargin{n}; else LL=1/4; end; n=n+1; % where lambda0 is the design wavelength in units of nm
na = 1; % refractive index of medium outside the sample (usually air or vacuum);
if strcmp(pol,'p')|strcmp(pol,'P')
pol='tm';
elseif strcmp(pol,'s')|strcmp(pol,'S')
pol='te';
elseif strcmp(pol,'tm')|strcmp(pol,'te')
else
error('Polarization string was not understood. Should be ''te'' or ''tm''.');
end
if length(nH)>1
n_in=nH;
n_out=reverse(n_in);
L_in=nL;
L_out=reverse(L_in);
na=n_in(1);
nb=n_in(end);
else
coating=repmat([nL,nH], 1, N);
L_in = (coating==nL)*LL + (coating==nH)*LH; % lengths of the layers in order as seen entering sample
if Lprotect~=0;
L_in=[Lprotect,L_in];
coating=[nprotect,coating];
end
L_out= reverse(L_in); % lengths of the layers in order as seen exiting sample
n_in = [na,coating,nb]; % indices for the layers in order as seen entering sample
n_out = reverse(n_in); % indices for the layers in order as seen exiting sample
end
angincg=asin(na*sin(anginc)/nb); %Snell's law
r_in=zeros(length(anginc),length(lambda));
r_out=zeros(length(anginc),length(lambda));
R_in=zeros(length(anginc),length(lambda));
R_out=zeros(length(anginc),length(lambda));
R_g=zeros(length(anginc),1); %Fresnel Eqn's don't depend on lambda explicitly, so neither does R_g
R=zeros(length(anginc),length(lambda));
[r_gs r_gp]=fresnel(nb,na,angincg*180/pi); % nb is the "left hand" medium in the call to fresnel because the beam is inside the sample
R_gs=r_gs.*conj(r_gs);
R_gp=r_gp.*conj(r_gp);
for s=1:length(anginc)
r_in(s,:)=multidiel(n_in,L_in,lambda/lambda0,anginc(s)*180/pi,pol);
r_out(s,:)=multidiel(n_out,L_out,lambda/lambda0,angincg(s)*180/pi,pol);
R_in(s,:)=(r_in(s,:).*conj(r_in(s,:)));
R_out(s,:)=(r_out(s,:).*conj(r_out(s,:)));
if strcmp(pol,'te'), R_g(s)=R_gs(s); else R_g(s)=R_gp(s); end
if strcmp(pol,'te'), r_g(s)=r_gs(s); else r_g(s)=r_gp(s); end
if nreflecs~=0
k=R_out(s,:)*R_g(s);
% next line is sum from 0 to nreflecs internally reflecting light beams
% NOTE: nreflecs = 2 is equivalent the having 3 beams in total reflected
% from the sample (primary reflection + 2 secondary reflections).
R(s,:)=R_in(s,:) + ( 1-R_in(s,:) ) * R_g(s) .* (1-k.^nreflecs)./( 1-k ) .* ( 1-R_out(s,:) );
if nargout==5,
warning('The amplitude reflection coefficient, r, does not include the effect of internal bounces. (Not yet implemented).');
end
r(s,:)=r_in(s,:);
else
R(s,:)=R_in(s,:);
r(s,:)=r_in(s,:);
end
end
T=1-R;
[lambda_mesh,anginc_mesh]=meshgrid(lambda,anginc);
|
