############################################################## # scriptfile: calculate_standard_image.lsf # # Description: Calculate image field from object field assuming # illumination with multiple plane waves focused at # a spot. # # Input parameters: # run_filename : the name of the file with the pre-calculated S parameter data # input_filter : illuminatino filter. The example shows setting a simple NA_in, but other filters could be created # Ein : the incident polarization vector # x0_illumination, y0_illumination, z0_illumination = position of illumination spot # # Ouput parameters: # output_filter: collection filter. The example shows simple NA_min to NA_max, but others can be created # # Imaging parameters: # M = demagnification factor # x_image, y_image, z_image : define imaging volume # # Notes on conventions: # The relationship between ES, EP and Ex,Ey is: # (EP) = ( cos(phi) sin(phi) ) * (Ex) # (ES) (-sin(phi) cos(phi) ) (Ey) # Circular polarization: # https://kb.lumerical.com/ref_sim_obj_planewave_circular_polarization.html # See also https://en.wikipedia.org/wiki/Jones_calculus # For light propagating upwards (along z axis) the relationship between ER, EL and Ex, Ey is: # (Ex) = (1/sqrt(2))(1 1) * (ER) # (Ey) (i -i) (EL) # and # (ER) = (1/sqrt(2))(1 -i) * (Ex) # (EL) (1 i) (Ey) # Note: if source is propagating DOWN, this is reversed. # Right polarized incident light will become Left on reflection from a flat surface # # Copyright 2018 Lumerical Solutions Inc. ######################################################### # clear existing workspace clear; ######################################################### # USER DEFINED PARAMETERS ######################################################### ######################################################### # run_filename : the name of the file with the pre-calculated S parameter data # This filename MUST be the same in the run and build S-matrix files run_filename = "metal_sim"; # name for ldf file with S-matrix, no extension # load the S matrix and unit vectors S_filename = pwd + "/" +run_filename + "_S_matrix"; loaddata(S_filename); ######################################################### ######################################################### # input_filter : illuminatino filter. The example shows setting a simple NA_in, but other filters could be created # Here we must construct the "input_filter" to define which input angles are used Ux_in = meshgridx(ux0,uy0); # do not modify this line Uy_in = meshgridy(ux0,uy0); # do not modify this line # construct the input illumination filter NA_in = 0.7; # define max NA for illumination Uxy_in = sqrt( Ux_in^2 + Uy_in^2 ); # do not modify this line input_filter = Uxy_in <= NA_in*1.001; #phi_in = atan2(Uy_in,Ux_in); #input_filter = input_filter & ((abs(phi_in)< pi/4) | (abs(phi_in-pi)< pi/4) | (abs(phi_in+pi)< pi/4)); polarimage(ux0,uy0,input_filter,"","","input filter"); ######################################################### ######################################################### # x0_illumination, y0_illumination, z0_illumination = position of illumination spot # these should be a single value x0_illumination = 0e-6; y0_illumination = 0e-6; z0_illumination = 0e-6; ######################################################### ######################################################### # input_polarization definition Ein = [1,0]; # linear polarized along x #Ein = [0,1]; # linear polarized along y #angle_in = -pi/4; Ein = [cos(angle_in), sin(angle_in) ]; # linearly polarized at angle angle_in (radians) #Ein = [1,1i]/sqrt(2); # circular RIGHT (because incident light is propagating down) #Ein = [1,-1i]/sqrt(2); # circular LEFT (because incident light is propagating down) ######################################################### ######################################################### # output polarizer projection (Jones matrix) # See also https://en.wikipedia.org/wiki/Jones_calculus output_polarizer = [ 1,0; 0,1]; # no polarizer #output_polarizer = [ 1,0; # 0,0]; # linear x #output_polarizer = [ 0,0; # 0,1]; # linear y #angle_out = pi/4; # linear at any angle, radians #output_polarizer = [ (cos(angle_out))^2,cos(angle_out)*sin(angle_out); # cos(angle_out)*sin(angle_out),(sin(angle_out))^2]; # linear at any angle #output_polarizer = [ 1,1i; # -1i,1]/2; # Right circular #output_polarizer = [ 1,-1i; # 1i,1]/2; # Left circular ######################################################### # output_filter: collection filter. The example shows simple NA_min to NA_max, but others can be created Ux_out = meshgridx(ux0,uy0); # do not modify this line Uy_out = meshgridy(ux0,uy0); # do not modify this line # construct the output collection filters NA_min = 0.7; NA_max = 0.9; output_filter = (Ux_out^2 + Uy_out^2 < NA_max^2) & (Ux_out^2 + Uy_out^2 >= NA_min^2); #phi_out = atan2(Uy_out,Ux_out); #output_filter = output_filter & ((abs(phi_out - pi/2)< pi/4) | (abs(phi_out+pi/2)< pi/4)); polarimage(ux0,uy0,output_filter,"","","output filter"); ######################################################### ######################################################### # M = demagnification factor M = 1; ######################################################### ######################################################### # x_image, y_image, z_image : define imaging volume x_image = linspace(-2.2e-6,2.2e-6,201); y_image = linspace(-2.1e-6,2.1e-6,201); z_image = 0e-6; #defocus ######################################################### ######################################################### # END USER DEFINED PARAMETERS ######################################################### # k vector lamda=c/f; # wavelength in microns k0=2*pi/lamda; ExS_FF = matrix(length(ux0),length(uy0)); ExP_FF = matrix(length(ux0),length(uy0)); EyS_FF = matrix(length(ux0),length(uy0)); EyP_FF = matrix(length(ux0),length(uy0)); # loop over all input angles for(i=1:length(ux0)) { ?"Calculating, " + num2str(round(i/length(ux0)*100)) + "% complete"; for(j=1:length(uy0)) { if(input_filter(i,j)) { uz0 = sqrt( 1-ux0(i)^2-uy0(j)^2 ); uxy0 = sqrt( ux0(i)^2 + uy0(j)^2 ); theta = acos(uz0); phi = atan2(uy0(j),ux0(i)); kx0 = -ux0(i)*k0; ky0 = -uy0(j)*k0; kz0 = -uz0*k0; input_factor0_xS = -sin(phi); input_factor0_xP = cos(phi); input_factor0_yS = cos(phi); input_factor0_yP = sin(phi); ESS = pinch(S(i,j,1,:,:,1)); ESP = pinch(S(i,j,1,:,:,2)); EPS = pinch(S(i,j,2,:,:,1)); EPP = pinch(S(i,j,2,:,:,2)); input_factor = exp( (-1i*kx0)*x0_illumination + (-1i*ky0)*y0_illumination + (-1i*kz0*z0_illumination) ); ExS_FF = ExS_FF + input_factor*(input_factor0_xS*ESS + input_factor0_xP*EPS); ExP_FF = ExP_FF + input_factor*(input_factor0_xS*ESP + input_factor0_xP*EPP); EyS_FF = EyS_FF + input_factor*(input_factor0_yS*ESS + input_factor0_yP*EPS); EyP_FF = EyP_FF + input_factor*(input_factor0_yS*ESP + input_factor0_yP*EPP); } } } # use sphereical coordinate system of theta and phi # combine contributions from both input polarizations EP_FF = Ein(1)*ExP_FF + Ein(2)*EyP_FF; ES_FF = Ein(1)*ExS_FF + Ein(2)*EyS_FF; # apply output polarizer # calculate 2D mesh of unit vectors Ux = meshgridx(ux0,uy0); Uy = meshgridy(ux0,uy0); sin_phi = sin(atan2(Uy,Ux)); cos_phi = cos(atan2(Uy,Ux)); Ex_FF = cos_phi*EP_FF - sin_phi*ES_FF; Ey_FF = sin_phi*EP_FF + cos_phi*ES_FF; tempx = output_polarizer(1,1)*Ex_FF + output_polarizer(1,2)*Ey_FF; tempy = output_polarizer(2,1)*Ex_FF + output_polarizer(2,2)*Ey_FF; Ex_FF = tempx; Ey_FF = tempy; EP_FF = cos_phi*Ex_FF + sin_phi*Ey_FF; ES_FF =-sin_phi*Ex_FF + cos_phi*Ey_FF; # change to sphereical coordinates to apply magnification Et = EP_FF; Ep = ES_FF; # calculate 2D mesh of unit vectors Ux = meshgridx(ux0,uy0); Uy = meshgridy(ux0,uy0); ## apply demagnification factor Ux=(Ux)*M; Uy=(Uy)*M; ux0 = ux0 * M; uy0 = uy0 * M; # assume magnification cannot lead to imaging # with NA > 1 output_filter = output_filter & (Ux^2 + Uy^2 < 1); Uxy = real(sqrt(Ux^2+Uy^2)); Uz = real(sqrt(1-Uxy^2)); sin_phi = sin(atan2(Uy,Ux)); cos_phi = cos(atan2(Uy,Ux)); cos_theta = Uz; sin_theta = Uxy; # calculate demagnified field components in image plane Exf = output_filter * (-Ep*sin_phi + Et*cos_phi*cos_theta); Eyf = output_filter * ( Ep*cos_phi + Et*sin_phi*cos_theta); Ezf = output_filter * ( -Et * sin_theta); # calculate the image, uzing chirped z-transform kx = ux0*k0; ky = uy0*k0; dkx = kx(2)-kx(1); dky = ky(2)-ky(1); Kz = Uz*k0; px = 1:length(x_image); py = 1:length(y_image); Ex_image = matrix(length(x_image),length(y_image),length(z_image)); Ey_image = matrix(length(x_image),length(y_image),length(z_image)); Ez_image = matrix(length(x_image),length(y_image),length(z_image)); Uz = Uz + 1e-5; for(pz=1:length(z_image)) { z0 = z_image(pz); Ex_image(px,py,pz) = czt(Exf*exp(1i*Kz*z0)*sqrt(1/Uz/k0^2),kx,ky,x_image,y_image)*dkx*dky/(2*pi); Ey_image(px,py,pz) = czt(Eyf*exp(1i*Kz*z0)*sqrt(1/Uz/k0^2),kx,ky,x_image,y_image)*dkx*dky/(2*pi); Ez_image(px,py,pz) = czt(Ezf*exp(1i*Kz*z0)*sqrt(1/Uz/k0^2),kx,ky,x_image,y_image)*dkx*dky/(2*pi); } E2_image = abs(Ex_image)^2 + abs(Ey_image)^2 + abs(Ez_image)^2; #Plot image: image(x_image/1e-6, y_image/1e-6 ,E2_image,"x (um)","y (um)","|E|^2 image");