############################################
# file: PSL_Cu_scattering.lsf 
# Reproduces the results in Kim et al. 
# figures number 6 to 9.
# 
# Copyright 2008 Lumerical Solutions Inc
############################################

#closeall;
#clear;


###########################################
# setup parameter sweep
radius=matrix(3);
wavelength=matrix(3);

fig8=0;   # if =1, make fig 8. if =0, make fig7.

if (fig8) {  # PSL Fig 8
  radius(1)=155e-9/2;
  radius(2)=155e-9/2;
  radius(3)=155e-9/2;
  wavelength(1)=442e-9;
  wavelength(2)=532e-9;
  wavelength(3)=633e-9;
  mat="PSL";
} else {  # Cu Fig 7
  radius(1)=92e-9/2;
  radius(2)=123e-9/2;
  radius(3)=155e-9/2;
  wavelength(1)=442e-9;
  wavelength(2)=442e-9;
  wavelength(3)=442e-9;
  mat="Cu";
}


###########################################
# analysis variable initialization
mname="R"; # set monitor name
resX=401;  # resolution of far field projections
resY=001;    
farfieldfilter(0.5);
theta=matrix(resX);
dS_dO_all= matrix(resX,3);
Pc_all =   matrix(resX,3);
eta_all =  matrix(resX,3);


###########################################
# do paramater sweep
for (loop=1:3) {
  rad=radius(loop);
  lambda=wavelength(loop);
  file="Mat"+mat+"_Dia"+num2str(2*rad*1e9)+"_Lambda"+num2str(lambda*1e9);
 
  PSL_Cu_scattering_helper;

  theta= theta;
  dS_dO_all(1:resX,loop)=dS_dO;
  Pc_all(1:resX,loop)=Pc;
  eta_all(1:resX,loop)=eta;
}


###########################################
# plot results
plot(theta,log10(dS_dO_all),"theta (deg)","dS/dO (log scale)","Diff Scat cross section");
if (fig8) { legend(num2str(wavelength(1)*1e9),num2str(wavelength(2)*1e9),num2str(wavelength(3)*1e9));
}  else   { legend(num2str(radius(1)*2*1e9),num2str(radius(2)*2*1e9),num2str(radius(3)*2*1e9));
}

plot(theta,Pc_all,"theta (deg)","Pc","Degree of Circular Polarization");
if (fig8) { legend(num2str(wavelength(1)*1e9),num2str(wavelength(2)*1e9),num2str(wavelength(3)*1e9));
}  else   { legend(num2str(radius(1)*2*1e9),num2str(radius(2)*2*1e9),num2str(radius(3)*2*1e9));
}
plot(theta,eta_all,"theta (deg)","Polarization angle","Polarization angle");
if (fig8) { legend(num2str(wavelength(1)*1e9),num2str(wavelength(2)*1e9),num2str(wavelength(3)*1e9));
}  else   { legend(num2str(radius(1)*2*1e9),num2str(radius(2)*2*1e9),num2str(radius(3)*2*1e9));
}
load("PSL_Cu_scattering.lsf");