Home > SPlaC v1_0 > Plane > PlnMultiRef.m

PlnMultiRef

PURPOSE ^

Computes the reflection/refraction problem at a multilayer interface for TM or TE polarizations

SYNOPSIS ^

function stResMulti = PlnMultiRef (sPol, nNbSurf, lambda, Cepsilon, CL, aideg, d)

DESCRIPTION ^

 Computes the reflection/refraction problem at a multilayer interface for TM or TE polarizations
 The problem is defined by j=1..nNbSurf+1 regions separated by
 i=1..nNbSurf interfaces at z=Z_i (see Fig. F.3). The first interface is at z=0.
 Region 1 is z<0 and contains the incident wave (incident from z=-infty)
 The plane of incidence is by convention xOz.
 
 Parameters:
 - sPol:   string
           sPol='TE' for TE or sPol='TM' for TM (default value)
           polarizations
 - nNbSurf: scalar integer [1x1]
            number of interfaces (i.e. there are nNbSurf+1 regions)
            see Sec. F.4.1 and Fig. F.3 for more details.
 - lambda: column vector [Lx1]
           wavelengths in NANOMETERS (nm)
 - Cepsilon: cell of nNbSurf+1 column vectors (same size as lambda):
             i.e.: {(nNbSurf+1) * [Lx1]}
             lambda-dependent epsilons in each region.
 - CL:     cell of nNbSurf scalar real numbers: {nNbSurf * [1x1]}
           region thicknesses (L_j=Z_{j+1}-Z_j) in NANOMETERS (i.e. distance between two interfaces).
           Note that CL{1}=L_1 must be zero by convention (see Fig. F.3)
 - aideg:  row vector [1xA]
           angles of incidence in DEGREES.
 - d:      optional scalar real number [1x1] (d=0 if not supplied)
           distance in NANOMETERS from interface where local field EFs are estimated
      
 Returns: a stResMulti structure containing fields as listed below:
 Note that EM fields are defined as in Eq. F.11 (TM) and Eq. F.18 (TE)

 For TM waves, we have:
 ** The Fresnel coefficients, which can be used for general properties:
 - rP: matrix [LxA]
       Fresnel reflection coefficients (Eqs. F.45, F.79)
       for each wavelength and angle of incidence
 - tP: matrix [LxA]
       Fresnel transmission coefficients (Eqs. F.46, F.80)
       for each wavelength and angle of incidence

 ** The field amplitudes, which can be used to calculate the EM field at
 any points:
 - HiyOvH1y: Cell of nNbSurf+1 matrices, {(nNbsurf+1)*[LxA]}
             field amplitude of first wave in each region i
             characterized by the ratio \bar{H}_{iy}/H_{1y}
             (see Eq. F.74 for a definition of \bar{H}_{iy})
             given for each wavelength and angle of incidence
 - HipyOvH1y: Cell of nNbSurf+1 matrices, {(nNbsurf+1)*[LxA]}
              same as HiyOvH1y but for second wave: \bar{H}'_{iy}/H_{1y}
 - EixOvE1: Cell of nNbSurf+1 matrices, {(nNbsurf+1)*[LxA]}
            \bar{E}_{ix}/E_1
            where E_1=H_{1y}/(\epsilon_0 c \sqrt{\epsilon_1})
 - EizOvE1: {(nNbsurf+1)*[LxA]}
            \bar{E}_{iz}/E_1
 - EipxOvE1: {(nNbsurf+1)*[LxA]}
            \bar{E}'_{ix}/E_1
 - EipzOvE1: {(nNbsurf+1)*[LxA]}
            \bar{E}'_{iz}/E_1

 ** The local field intensity enhancement factors (LFIEF) at a distance d on either side of
 each of the nNbSurf interface.
 For interface i, "outside" means towards region i and inside means
 towards region i+1.
 - MoutPerp: {nNbSurf * [LxA]}
             perpendicular LFIEF, i.e. |E_z|^2/|E_1|^2 at a distance d
             on the "outside" of each interface
 - MoutPara: {nNbSurf * [LxA]}
             parallel LFIEF, i.e. |E_x|^2/|E_1|^2 at a distance d
             on the "outside" of each interface
 - MinPerp: {nNbSurf * [LxA]}
             perpendicular LFIEF, i.e. |E_z|^2/|E_1|^2 at a distance d
             on the "inside" of each interface
 - MinPara: {nNbSurf * [LxA]}
             parallel LFIEF, i.e. |E_x|^2/|E_1|^2 at a distance d
             on the "inside" of each interface


 For TE waves, the output is similar but the name are changed to reflect
 the different nature of the fiels:
 - rS: matrix [LxA]
       Fresnel reflection coefficients (Eqs. F.62, F.79)
       for each wavelength and angle of incidence
 - tS: matrix [LxA]
       Fresnel transmission coefficients (Eqs. F.63, F.80)
       for each wavelength and angle of incidence
 - EiyOvE1y: {(nNbsurf+1)*[LxA]}
             field amplitude of first wave in each region i
             characterized by the ratio \bar{E}_{iy}/E_{1y}
             (see Eq. F.74 for a definition of \bar{E}_{iy})
             given for each wavelength and angle of incidence
 - EipyOvE1y: {(nNbsurf+1)*[LxA]}
              same as EiyOvE1y but for second wave: \bar{E}'_{iy}/E_{1y}
 - HixOvH1: {(nNbsurf+1)*[LxA]}
            \bar{H}_{ix}/H_1
            where H_1=E_{1y}*(\epsilon_0 c \sqrt{\epsilon_1})
 - HizOvH1: {(nNbsurf+1)*[LxA]}
            \bar{H}_{iz}/H_1
 - HipxOvH1: {(nNbsurf+1)*[LxA]}
            \bar{H}'_{ix}/H_1
 - HipzOvH1: {(nNbsurf+1)*[LxA]}
            \bar{H}'_{iz}/H_1
 - MoutPerp: {nNbSurf * [LxA]}
             perpendicular LFIEF, is zero for TE waves
 - MoutPara: {nNbSurf * [LxA]}
             parallel LFIEF, i.e. |E_x|^2/|E_1|^2 at a distance d
             on the "outside" of each interface
 - MinPerp: {nNbSurf * [LxA]}
             perpendicular LFIEF is zero for TE waves
 - MinPara: {nNbSurf * [LxA]}
             parallel LFIEF, i.e. |E_x|^2/|E_1|^2 at a distance d
             on the "inside" of each interface
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 This file is part of the SPlaC v1.0 package (copyright 2008)
 Check the README file for further information

EXAMPLE OF OUTPUT ^

Example figure output

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 
0002 %%%%%%%%%%%%%%%%%%%%%%%
0003 % Check parameters
0004 %%%%%%%%%%%%%%%%%%%%%%%
0005 % optional parameters
0006 if nargin < 7, d = 0; end
0007 
0008 if ~((strcmp(sPol,'TM')) || (strcmp(sPol,'TE')))
0009     disp 'Warning (PlnMultiRef): Polarization not properly specified. Assuming TM polarization.';
0010     sPol='TM';
0011 end
0012 
0013 if (CL{1}~=0)
0014     disp 'Warning (PlnMultiRef): L_1 must be zero. Check CL parameter. Changing to CL{1}=0';
0015     % enforce L_1=0 (by convention)
0016     CL{1}=0;
0017 end
0018 
0019 %%%%%%%%%%%%%%%%%%%%%%%
0020 % Computations
0021 %%%%%%%%%%%%%%%%%%%%%%%
0022 
0023 nNbAngle=length(aideg); % nb of columns (angles), i.e. A
0024 nNbLambda=length(lambda); % nb of rows (wavelengths), i.e. L
0025 
0026 k0=2*pi./lambda; % column [Lx1]
0027 
0028 % angle of incidence in radian
0029 ai=aideg/180*pi; % row [1xA]
0030 
0031 % calculate k_x (same for all regions), Eq. F.42
0032 k1=k0.*sqrt(Cepsilon{1}); % column [Lx1]
0033 kx=k1 * sin(ai); % matrix [LxA] (lambda in rows and angles in columns)
0034 
0035 % calculates k_{iz} in each region (matrices) and place in a cell
0036 Ckiz=cell(1,nNbSurf+1);
0037 for ii=1:(nNbSurf+1)
0038     % each CKiz{ii} is a matrix [LxA]
0039     % Computed from Eq. F.71
0040     Ckiz{ii}=sqrt(repmat(Cepsilon{ii}.*k0.^2,1,nNbAngle)-kx.^2); 
0041     % Note that k_{iz} is assumed non-zero
0042     % except possibly in last region i=nNbSurf+1
0043     % Check this and issue a warning if not the case.
0044     % Assume k_{iz} is zero if k_{iz}/k0 smaller than 1e-6
0045     if (ii<nNbSurf+1)
0046         indPb=find(abs(Ckiz{ii}./repmat(k0,1,nNbAngle))<1e-6);
0047         for jj=1:length(indPb)
0048             % issue warning for each problem
0049             disp (['Warning: kiz is close to zero in region ' ...
0050                 num2str(ii) ' at wavelength ' num2str(lambda(mod(indPb(jj),nNbLambda))) ...
0051                 ' nm and angle ' num2str(aideg(floor(indPb(jj)-1)/nNbLambda)) ' degrees.']);
0052             disp (['kiz = ' num2str(Ckiz{indPb(jj)}/k0(mod(indPb(jj),nNbLambda))) ' * k0']);
0053         end
0054     end
0055 end
0056 
0057 
0058 % initialize M matrix to identity
0059 M11=1;
0060 M21=0;
0061 M12=0;
0062 M22=1;
0063 
0064 % loop through the interfaces to compute Mi (Eq. F.76) and M (Eq. F.78)
0065 Mi11=cell(1,nNbSurf);
0066 Mi21=cell(1,nNbSurf);
0067 Mi12=cell(1,nNbSurf);
0068 Mi22=cell(1,nNbSurf);
0069 for ii=1:nNbSurf
0070     % calculate Ki [LxA]
0071     if (strcmp(sPol,'TM')) % TM wave, see Eq. F.16 and Eq. F.72
0072         Ki=repmat(Cepsilon{ii}./Cepsilon{ii+1},1,nNbAngle).*Ckiz{ii+1}./Ckiz{ii};
0073     else % TE wave, see Eq. F.22
0074         Ki=Ckiz{ii+1}./Ckiz{ii};
0075     end
0076     
0077     % calculate matrix Mi components (all are matrices [LxA])
0078     % using Eq. F.76
0079     AuxMi=exp(i*CL{ii}*Ckiz{ii});
0080     Mi11{ii}=0.5*(1+Ki)./ AuxMi;
0081     Mi21{ii}=0.5*(1-Ki).* AuxMi;
0082     Mi12{ii}=0.5*(1-Ki)./ AuxMi;
0083     Mi22{ii}=0.5*(1+Ki).* AuxMi;
0084 
0085     % update M matrix (Eq. F.78)
0086     M11new=M11.*Mi11{ii} + M12.*Mi21{ii};
0087     M21new=M21.*Mi11{ii} + M22.*Mi21{ii};
0088     M12new=M11.*Mi12{ii} + M12.*Mi22{ii};
0089     M22new=M21.*Mi12{ii} + M22.*Mi22{ii};
0090     M11=M11new;
0091     M21=M21new;
0092     M12=M12new;
0093     M22=M22new;
0094 end
0095 
0096 %%%%%%%%%%%%%%%%%%%%%%%
0097 % Outputs
0098 %%%%%%%%%%%%%%%%%%%%%%%
0099 
0100 % calculate Fresnel coefficients
0101 % see Eq. F.79 and F.80 (for TM)
0102 % same expressions for TE but relate to electric field
0103 tF=1./M11;  
0104 rF=M21.*tF; % Eq. F.80 (for TM)
0105 
0106 % initialize local field intensity EFs at each interface
0107 MoutPerp=cell(1,nNbSurf);
0108 MoutPara=cell(1,nNbSurf);
0109 MinPerp=cell(1,nNbSurf);
0110 MinPara=cell(1,nNbSurf);
0111 
0112 
0113 % treat separately TM and TE from now on
0114 % All field outputs are cells of matrices [LxA]
0115 if (strcmp(sPol,'TM'))
0116 %%%%%%%%%%%%%%%%%
0117 % TM waves
0118 %%%%%%%%%%%%%%%%%
0119     % initialize field amplitudes in each region
0120     HiyOvH1y=cell(1,nNbSurf+1); 
0121     HpiyOvH1y=cell(1,nNbSurf+1); 
0122     EixOvE1=cell(1,nNbSurf+1);
0123     EpixOvE1=cell(1,nNbSurf+1);
0124     EizOvE1=cell(1,nNbSurf+1);
0125     EpizOvE1=cell(1,nNbSurf+1);
0126 
0127     % Fields in last region (nNbSurf+1)
0128     HiyOvH1y{nNbSurf+1}=tF; % H_{iy} / H_{1y}
0129     HpiyOvH1y{nNbSurf+1}=0; % H'_{iy} / H_{1y}
0130     % use for reference (possibly complex) E field: E1= H_{1y}/|H_{1y}| |E_1|
0131     % E1=H_{1y}/ [(epsilon_0 c) sqrt{epsilon_1}]
0132     AuxE1=repmat(sqrt(Cepsilon{1})./k0./Cepsilon{nNbSurf+1},1,nNbAngle);
0133     EixOvE1{nNbSurf+1}=HiyOvH1y{nNbSurf+1}.* Ckiz{nNbSurf+1} .* AuxE1;
0134     EpixOvE1{nNbSurf+1}=0;
0135     EizOvE1{nNbSurf+1}=-HiyOvH1y{nNbSurf+1}.* ( (Cepsilon{1}./Cepsilon{nNbSurf+1}) * sin(ai) );
0136     EpizOvE1{nNbSurf+1}=0;
0137 
0138     % loop (downwards) to compute field amplitudes in other regions
0139     for ii=nNbSurf:-1:1
0140         HiyOvH1y{ii}=Mi11{ii}.*HiyOvH1y{ii+1}+Mi12{ii}.*HpiyOvH1y{ii+1};
0141         HpiyOvH1y{ii}=Mi21{ii}.*HiyOvH1y{ii+1}+Mi22{ii}.*HpiyOvH1y{ii+1};
0142     
0143         AuxE1=repmat(sqrt(Cepsilon{1})./k0./Cepsilon{ii},1,nNbAngle);
0144         EixOvE1{ii}=HiyOvH1y{ii}.* Ckiz{ii} .* AuxE1;
0145         EpixOvE1{ii}=-HpiyOvH1y{ii}.* Ckiz{ii} .* AuxE1;
0146         AuxE2=(Cepsilon{1}./Cepsilon{ii}) * sin(ai) ; % matrix product
0147         EizOvE1{ii}=-HiyOvH1y{ii}.* AuxE2;
0148         EpizOvE1{ii}=-HpiyOvH1y{ii}.* AuxE2;
0149     end
0150 
0151     % loop again to compute local field EFs on each of the nNbSurf
0152     % interfaces
0153     for ii=1:nNbSurf
0154         MoutPerp{ii}=abs( EizOvE1{ii}.*exp(i*(CL{ii}-d)*Ckiz{ii})+ ...
0155             EpizOvE1{ii}.*exp(-i*(CL{ii}-d)*Ckiz{ii}) ).^2;
0156         MoutPara{ii}=abs( EixOvE1{ii}.*exp(i*(CL{ii}-d)*Ckiz{ii})+ ...
0157             EpixOvE1{ii}.*exp(-i*(CL{ii}-d)*Ckiz{ii}) ).^2;
0158 
0159         MinPerp{ii}=abs( EizOvE1{ii+1}.*exp(i*d*Ckiz{ii+1})+ ...
0160             EpizOvE1{ii+1}.*exp(-i*d*Ckiz{ii+1}) ).^2;
0161         MinPara{ii}=abs( EixOvE1{ii+1}.*exp(i*d*Ckiz{ii+1})+ ...
0162             EpixOvE1{ii+1}.*exp(-i*d*Ckiz{ii+1}) ).^2;
0163     end
0164 
0165 
0166 else 
0167 %%%%%%%%%%%%%%%%%
0168 % TE waves
0169 %%%%%%%%%%%%%%%%%
0170     
0171     % initialize field amplitudes in each region
0172     EiyOvE1y=cell(1,nNbSurf+1); 
0173     EpiyOvE1y=cell(1,nNbSurf+1); 
0174     HixOvH1=cell(1,nNbSurf+1);
0175     HpixOvH1=cell(1,nNbSurf+1);
0176     HizOvH1=cell(1,nNbSurf+1);
0177     HpizOvH1=cell(1,nNbSurf+1);
0178 
0179     % Fields in last region (nNbSurf+1)
0180     EiyOvE1y{nNbSurf+1}=tF; % E_{iy} / E_{1y}
0181     EpiyOvE1y{nNbSurf+1}=0; % E'_{iy} / E_{1y}
0182     % use for reference (possibly complex) H field: H1= E_{1y}/|E_{1y}| |H_1|
0183     % H1=E_{1y}* [(epsilon_0 c) sqrt{epsilon_1}]
0184     AuxH1=repmat(1./(sqrt(Cepsilon{1}).*k0),1,nNbAngle);
0185     HixOvH1{nNbSurf+1}=-EiyOvE1y{nNbSurf+1}.* Ckiz{nNbSurf+1} .* AuxH1;
0186     HpixOvH1{nNbSurf+1}=0;
0187     AuxH2=repmat(sin(ai),nNbLambda,1);
0188     HizOvH1{nNbSurf+1}=EiyOvE1y{nNbSurf+1}.* AuxH2;
0189     HpizOvH1{nNbSurf+1}=0;
0190 
0191     % loop (downwards) to compute field amplitudes in other regions
0192     for ii=nNbSurf:-1:1
0193         EiyOvE1y{ii}=Mi11{ii}.*EiyOvE1y{ii+1}+Mi12{ii}.*EpiyOvE1y{ii+1};
0194         EpiyOvE1y{ii}=Mi21{ii}.*EiyOvE1y{ii+1}+Mi22{ii}.*EpiyOvE1y{ii+1};
0195     
0196         HixOvH1{ii}=-EiyOvE1y{ii}.* Ckiz{ii} .* AuxH1;
0197         HpixOvH1{ii}=EpiyOvE1y{ii}.* Ckiz{ii} .* AuxH1;
0198         HizOvH1{ii}=EiyOvE1y{ii}.* AuxH2;
0199         HpizOvH1{ii}=EpiyOvE1y{ii}.* AuxH2;
0200     end
0201 
0202     % loop again to compute local field EFs on each of the nNbSurf
0203     % interfaces
0204     for ii=1:nNbSurf
0205         MoutPara{ii}=abs( EiyOvE1y{ii}.*exp(i*(CL{ii}-d)*Ckiz{ii})+ ...
0206             EpiyOvE1y{ii}.*exp(-i*(CL{ii}-d)*Ckiz{ii}) ).^2;
0207 
0208         MinPara{ii}=abs( EiyOvE1y{ii+1}.*exp(i*d*Ckiz{ii+1})+ ...
0209             EpiyOvE1y{ii+1}.*exp(-i*d*Ckiz{ii+1}) ).^2;
0210         MoutPerp{ii}=zeros(nNbLambda,nNbAngle);
0211         MinPerp{ii}=zeros(nNbLambda,nNbAngle);
0212     end
0213 end
0214 
0215 % Finally,
0216 % create structure for results
0217 
0218 % results common to TM and TE waves
0219 stResMulti.sPol=sPol;
0220 stResMulti.nNbSurf=nNbSurf;
0221 stResMulti.lambda=lambda;
0222 stResMulti.Cepsilon=Cepsilon;
0223 stResMulti.CL=CL;
0224 stResMulti.aideg=aideg;
0225 stResMulti.d=d;
0226 
0227 stResMulti.MoutPerp=MoutPerp;
0228 stResMulti.MoutPara=MoutPara;
0229 stResMulti.MinPerp=MinPerp;
0230 stResMulti.MinPara=MinPara;
0231 
0232 
0233 if (strcmp(sPol,'TM')) % TM waves
0234     stResMulti.rP=rF;
0235     stResMulti.tP=tF;
0236     stResMulti.HiyOvH1y=HiyOvH1y;
0237     stResMulti.HpiyOvH1y=HpiyOvH1y;
0238     stResMulti.EixOvE1=EixOvE1;
0239     stResMulti.EizOvE1=EizOvE1;
0240     stResMulti.EpixOvE1=EpixOvE1;
0241     stResMulti.EpizOvE1=EpizOvE1;
0242 else % TE waves
0243     stResMulti.rS=rF;
0244     stResMulti.tS=tF;
0245     stResMulti.EiyOvE1y=EiyOvE1y;
0246     stResMulti.EpiyOvE1y=EpiyOvE1y;
0247     stResMulti.HixOvH1=HixOvH1;
0248     stResMulti.HizOvH1=HizOvH1;
0249     stResMulti.HpixOvH1=HpixOvH1;
0250     stResMulti.HpizOvH1=HpizOvH1;
0251 end
0252 
0253 
0254 
0255 
0256 
0257 
0258

This web page is part of the SPlaC package © 2008. Contact: Eric Le Ru
Generated on Wed 03-Dec-2008 11:10:14 by m2html © 2003 (adapted)