Low-level function to calculate VSH expansions for PWE for a fixed lambda, but many r and theta The VSHs used can either be of type 1 (sBessel='j') or type 3 (sBessel='h1') For PWE, we have: |m|=1, c_{n,1}=c_{n,-1}, d_{n,1}=-d_{n,-1}). The fields Ecr, Ect, Esf given in the results are discussed in the supplementary information. lambda0 and epsilon0 MUST BE A SCALAR here (one wavelength only) Parameters: - lambda0: scalar [1 x 1] wavelengths in nm - epsilon0: scalar [1 x 1] epsilon(lambda0) of dielectric where field is evaluated - cn1, dn1: 2 row vectors [1 x nNmax] containing the Mie coefficients c_{n,1} and d_{n,1} for the field expansion of M^(i)_{n,1} and N^(i)_{n,1}, respectively - r: column [R x 1], possibly with zero-component (if sBessel='j') spherical coordinate r (in nm) of points - theta: possibly row vector [1 x T] with spherical coordinate theta of points - sBessel: string defining the Bessel function to be used sBessel='j' for or sBessel='h1' - stPinTaun: structure (optional) with functions of theta pi_n and tau_n if omitted, then the functions are computed from scrath it is faster to pass this structure as argument if these functions have already been calculated Returns: stEAllPhi, structure with 3 fields containing matrices [R x T] representing the three components E_{cr}, E_{ct}, E_{sf} such as: of E = E_{cr} cos(phi) e_r + E_{ct} cos(phi) e_theta + E_{sf} sin(phi) e_phi - stEAllPhi.Ecr is E_{cr} - stEAllPhi.Ect is E_{ct} - stEAllPhi.Esf is E_{sf} This file is part of the SPlaC v1.0 package (copyright 2008) Check the README file for further information
0001 nNmax=size(cn1,2); 0002 % check that only one lambda 0003 if ~((size(lambda0,1)==1) && (size(epsilon0,1)==1) && (size(cn1,1)==1) && (size(dn1,1)==1)) 0004 disp 'Error in PweEgenRThetaAllPhi: only one wavelength allowed, lambda0 and epsilon0 should be scalar.') 0005 end; 0006 0007 nNbR=length(r); 0008 nNbTheta=length(theta); 0009 0010 % get Zn(rho) for radial dependence 0011 rho=2*pi*sqrt(epsilon0)/lambda0 * r; % column [R x 1] 0012 % exclude rho=0 (treat separately) and replace by rho=1 0013 indRho0=find(rho==0); 0014 rho(indRho0)=1; 0015 stZnAll=GenZnAll(nNmax,rho,sBessel); % fields are [R x nNmax] 0016 0017 % get theta dependence if not provided 0018 if nargin < 8 0019 stPinTaun=PwePinTaun(nNmax,transpose(theta)); % fields are [T x nNmax] 0020 end 0021 0022 n=1:nNmax; %[1 x nNmax] 0023 % need n-dep coeff for series (mu_n *n *(n+1) and mu_n) 0024 % can be multiplied already by c_{n,1} and d_{n,1} to improve speed 0025 cffnr=sqrt((2*n+1)/(4*pi)); %[1 x nNmax] 0026 mun=cffnr./(n.*(n+1)); 0027 cn1mun=cn1.*mun; % cn1*mu_n [1 x nNmax] for Et and Ef 0028 dn1cffnr=dn1.*cffnr; %[1 x nNmax] for Er 0029 dn1mun=dn1.*mun; % cn1*mu_n [1 x nNmax] for Et and Ef 0030 clear cffnr mun; 0031 0032 % loop over rho's to avoid using 3-dimensional matrices: 0033 % At a fixed rho, the sum over n for all theta can be carried out 0034 % using a matrix product of the theta-and-n-dependent [T x nNmax] matrix 0035 % by a n-dependent column vector [nNmax x 1] 0036 % the result, a [T x 1] matrix is then transposed to a [1 x T] line 0037 Ersum=zeros(nNbR,nNbTheta); 0038 Etsum=zeros(nNbR,nNbTheta); 0039 Efsum=zeros(nNbR,nNbTheta); 0040 for ll=1:nNbR 0041 if isempty(find(indRho0==ll,1)) % will treat case rho(ll)==0 separately 0042 % for Er, vecN=d_{n,1} * Z_n^1(rho) * mu_n * n *(n+1) 0043 vecNdep=transpose(stZnAll.Z1(ll,:).*dn1cffnr); % [nNmax x 1] 0044 % Ersum=sum_n(pi_n(t) * vecN_n) 0045 % Do the sum as a matrix product (and check for loss of precision) 0046 Ersum(ll,:)=transpose(GenCheckSumMatVec(stPinTaun.pin, vecNdep, ... 0047 'Ersum','PweEgenRThetaAllPhi')); % [1 x T] 0048 0049 % for Et and Ef 0050 vecNdep=transpose(i*stZnAll.Z0(ll,:).*cn1mun); % [nNmax x 1] 0051 vecNdep2=transpose(stZnAll.Z2(ll,:).*dn1mun); % [nNmax x 1] 0052 % Do the sums and check for loss of precision 0053 tmp1=GenCheckSumMatVec(stPinTaun.pin, vecNdep, ... 0054 'tmp1_1','PweEgenRThetaAllPhi'); 0055 tmp2=GenCheckSumMatVec(stPinTaun.taun, vecNdep2, ... 0056 'tmp2_1','PweEgenRThetaAllPhi'); 0057 Etsum(ll,:)=transpose(GenCheckSum2Mat(tmp1,tmp2, ... 0058 'Etsum','PweEgenRThetaAllPhi')); % [1 x T] 0059 0060 tmp1=GenCheckSumMatVec(stPinTaun.taun, vecNdep, ... 0061 'tmp1_2','PweEgenRThetaAllPhi'); 0062 tmp2=GenCheckSumMatVec(stPinTaun.pin, vecNdep2, ... 0063 'tmp2_2','PweEgenRThetaAllPhi'); 0064 Efsum(ll,:)=transpose(GenCheckSum2Mat(tmp1,tmp2, ... 0065 'Efsum','PweEgenRThetaAllPhi')); % [1 x T] 0066 else 0067 % special case where rho=0 0068 % the following is equivalent to \mathbf{E}=-4/3 * mu_1 *dn1 \mathbf{e}_x 0069 Ersum(ll,:)=repmat(2/3*dn1mun(1),1,nNbTheta); 0070 Etsum(ll,:)=2/3*dn1mun(1)* cos(theta); 0071 Efsum(ll,:)=Ersum(ll,:); 0072 disp 'r=0' 0073 end 0074 end 0075 0076 % Results are all [R x T] matrices 0077 stEAllPhi.Ecr=-2*repmat(sin(theta),nNbR,1).*Ersum; 0078 stEAllPhi.Ect=-2*Etsum; 0079 stEAllPhi.Esf=2*Efsum; 0080 0081