Welcome to PyFocus’s documentation!¶
Citing¶
You can find detailed information about the theorical calculations used in this package and a full description of many common use cases in the following paper: https://doi.org/10.1016/j.cpc.2022.108315. If use this package for your research, please cite this DOI.
Installation¶
At https://github.com/StefaniLab/PyFocus we provide a Windows executable version of PyFocus’s GUI, called “PyFocus.exe”, that can be downloaded and executed without any installations or Python enviroment.
For instalation of the highlevel functions and the GUI class for use in custom scripts, first we install a new enviroment in wich we use Python 3.6. This version is needed for using the GUI, but is not required for the highlevel functions we provide. From an anaconda prompt with admin privileges:
conda create name pyfocus python=3.6
conda activate pyfocus
To install from pip, for package name issues the package name was chosen as PyCustomFocus:
pip install PyCustomFocus
Dependencies¶
PyFocus uses various packages. Here we show an example of the needed installation, some packages are preferably installed using conda:
conda install numpy
conda install scipy
conda install qdarkstyle
pip install PyCustomFocus
pip install config==0.5.1
pip install tqdm==4.62.3
pip install matplotlib==3.3.4
pip install PyQt5==5.15.4
pip install qtpy==1.11.2
pip install pyqtgraph==0.11.1
pip install configparser==5.0.2
Package versions are specified for consistency. If a package does not have the specified version available (represented by the ==version part of pip install), you can try removing this restrain (by deleting the ==version part of the command).
After this, we can use a python script.
Importing¶
The modules used in the examples are imported with:
from PyFocus import sim, plot
import numpy as np
Examples using PyFocus User Interface¶
First we import the user interface class and define an instance:
from PyFocus import user_interface
ui=user_interface.UI()
To show the interface we use the function “show”:
ui.show()
If this example does not work, please contact us. You can continue with the rest of the examples using PyFocus’s highlevel functions, or use the previously mentioned “PyFocus.exe” windows executable.
Examples using PyFocus in custom codes¶
First, we define the parameters of the simulation:
 NA (float):
Numerical aperture
 n (float or array):
Type depends on interface: float for no interface, array for interface. Refraction index for the medium of the optical system.
 h (float):
Radius of aperture of the objective lens(mm)
 w0 (float):
Radius of the incident gaussian beam (mm)
 wavelength (float):
Wavelength in vacuum
 gamma (float):
Parameter that determines the polarization, arctan(ey/ex) (gamma=45, beta=90 gives right circular polarization and a donut shape)
 beta (float):
parameter that determines the polarization, phase difference between ex and ey (gamma=45, beta=90 gives right circular polarization and a donut shape)
 z (float):
Axial position for the XY plane (nm)
 x_steps (float):
Resolution in the X or Y coordinate for the focused field (nm)
 z_steps (float):
Resolution in the axial coordinate (Z) for the focused field (nm)
 x_range (float):
Field of view in the X or Y coordinate in which the focused field is calculated (nm)
 z_range (float):
field of view in the axial coordinate (z) in which the focused field is calculated (nm)
 I_0 (float):
Incident field intensity (mW/cm^2)
 L (float, optional):
Distance between phase mask and objective lens (mm), only used if propagation=True
 R (float, optional):
Phase mask radius (mm), only used if propagation=True
 ds (array, optional):
Thickness of each interface in the multilayer system (nm), must be given as a numpy array with the first and last values a np.inf. Only used if interface=True
 z_int (float, optional):
Axial position of the interphase. Only used if interface=True
 figure_name (string, optional):
Name for the images of the field. Also used as saving name if using the UI
NA=1.4
n=1.5
h=3
w0=5
wavelength=640
I0=1
We define an incident beam of left circular polarization:
gamma=45
beta=90
Next, we define the range in which we realize the simulation, the resolution and the name of the resulting figures:
z=0
x_steps=5
z_steps=10
x_range=1000
z_range=2000
figure_name='Example'
For the first examples we will neglect the propagation of the incident field and won’t simulate an interface and so the parameters L, R, ds and z_int won’t be used by the code. For completion, we will define them as ‘’:
L=''
R=''
ds=''
z_int=''
To simplify the code, we define the “parameters” array:
parameters=np.array((NA, n, h, w0, wavelength, gamma, beta, z, x_steps, z_steps, x_range, z_range, I0, L, R, ds, z_int, figure_name), dtype=object)
Gaussian beam without modulation¶
Simulation of a focused gaussian beam without phase modulation by using the “no_mask” function. We obtain the “fields” tuple, which contains 6 arrays with the resulting field, and then we plot them using the “plot_XZ_XY” function. This function returns a figure of the total intensity and polarization, here called fig1 and a figure of the intensity of each cartesian component:
fields=sim.no_mask(False,False,*parameters)
fig1,fig2=plot.plot_XZ_XY(*fields,x_range,z_range,figure_name)
VP mask modulation¶
Simulation of a focused gaussian beam with VP modulation by using the “VP” function. Remember that the 2 first ariables are boolean parameters that define if we simulate the propagation of the incident field and if there is an interface.
fields=sim.VP(False,False,*parameters)
fig1,fig2=plot.plot_XZ_XY(*fields,x_range,z_range,figure_name)
Propagation and modulation by a VP mask¶
To simulate the propagation of the incident field from the phase mask to the objective lens, we redefine the parameters L and R and set the propagation variable of VP to True:
L=1000
R=5
parameters=np.array((NA, n, h, w0, wavelength, gamma, beta, z, x_steps, z_steps, x_range, z_range, I0, L, R, ds, z_int, figure_name), dtype=object)
fields=sim.VP(True,False,*parameters)
fig1,fig2=plot.plot_XZ_XY(*fields,x_range,z_range,figure_name)
Interface and modulation by a VP mask¶
To simulate a glasswater interface located at the focal plane, we redefine the parameters n, ds and z_int and set the interface variable of VP to True:
n=np.array((1.5,1.33))
ds= np.array((np.inf,np.inf))
z_int=0
parameters=np.array((NA, n, h, w0, wavelength, gamma, beta, z, x_steps, z_steps, x_range, z_range, I0, L, R, ds, z_int, figure_name), dtype=object)
fields=sim.VP(False,True,*parameters)
fig1,fig2=plot.plot_XZ_XY(*fields,x_range,z_range,figure_name)
Multilayer system¶
For a multilayer system, we add 2 more layers of refraction index 0.14+3.55j and 1.54, and thicknesses 44 and 24:
n=np.array((1.5,0.14+3.55j,1.54,1.33))
ds= np.array((np.inf,44,24,np.inf))
z_int=0
parameters=np.array((NA, n, h, w0, wavelength, gamma, beta, z, x_steps, z_steps, x_range, z_range, I0, L, R, ds, z_int, figure_name), dtype=object)
fields=sim.VP(False,True,*parameters)
fig1,fig2=plot.plot_XZ_XY(*fields,x_range,z_range,figure_name)
Custom mask examples¶
Missalignment¶
First, we simulate the displacement of the gaussian beam and the VP mask. We define the auxiliar variables rho2 and phi2, which correspond to the displaced coordinates. We set the displacement to 0.4 times the aperture radius of the objective lens.
d=0.4*h
rho2=lambda rho,phi:(rho**2+d**22*rho*d*np.cos(phi))**0.5
phi2=lambda rho,phi:np.arctan2(rho*np.sin(phi),rho*np.cos(phi)d)
Then we define the mask function and begin the simulation. We set the 2D integration divisions to be 200 for both variables:
entrance_field=lambda rho,phi,w0,f,k: np.exp((rho2(rho,phi)/w0)**2)
custom_mask=lambda rho,phi,w0,f,k: np.exp(1j*phi2(rho,phi))
divisions_theta=200
divisions_phi=200
fields=sim.custom(entrance_field,custom_mask,False,False,*parameters,divisions_theta,divisions_phi)
fig1,fig2=plot.plot_XZ_XY(*fields,x_range,z_range,figure_name)
Inclination (tilt)¶
For inclination in an angle of 3.11*10**5 radians:
angle=3.11*10**5
entrance_field=lambda rho,phi,w0,f,k: np.exp((rho/w0)**2)
custom_mask=lambda rho,phi,w0,f,k: np.exp(1j*(phi+k*np.sin(ang)*rho*np.cos(phi)))
divisions_theta=200
divisions_phi=200
fields=sim.custom(entrance_field,custom_mask,False,False,*parameters,divisions_theta,divisions_phi)
fig1,fig2=plot.plot_XZ_XY(*fields,x_range,z_range,figure_name)
TIRF¶
After setting a waterglass interface, to simulate TIRF with modulation of a VP mask we define the mask function the following way:
n=np.array((1.5,1.33))
ds= np.array((np.inf,np.inf))
z_int=0
parameters=np.array((NA, n, h, w0, wavelength, gamma, beta, z, x_steps, z_steps, x_range, z_range, I0, L, R, ds, z_int, figure_name), dtype=object)
theta_crit=np.arcsin(n[1]/n[0])
entrance_field=lambda rho,phi,w0,f,k: np.exp((rho/w0)**2)
def custom_mask(rho,phi,w0,f,k):
if rho>f*np.sin(theta_crit):
return np.exp(1j*phi)
else:
return 0
fields=sim.custom(entrance_field,custom_mask,False,True,*parameters,200,200)
fig1,fig2=plot.plot_XZ_XY(*fields,x_range,z_range,figure_name)
Classes for the user interface¶
Highlevel simulation functions¶
Highlevel functions to simulate the foci obtained when: Focusing a gaussian beam, focusing a gaussian beam modulated with a VP mask and focusing a gaussian beam modulated with a custom phase mask
The functions use the ‘parameters’ array, defined as:
parameters=np.array((NA, n, h, w0, wavelength, gamma, beta, z, x_steps, z_steps, x_range, z_range, I0, L, R, ds, z_int, figure_name), dtype=object)
With:
 NA (float):
Numerical aperture
 n (float or array):
Type depends on presence of multilayer: float for no multilayer, array for multilayer. Refraction index for the medium of the optical system.
 h (float):
Radius of aperture of the objective lens(mm)
 w0 (float):
Radius of the incident gaussian beam (mm)
 wavelength (float):
Wavelength in vacuum
 gamma (float):
Parameter that determines the polarization, arctan(ey/ex) (gamma=45, beta=90 gives right circular polarization and a donut shape)
 beta (float):
parameter that determines the polarization, phase difference between ex and ey (gamma=45, beta=90 gives right circular polarization and a donut shape)
 z (float):
Axial position for the XY plane (nm)
 x_steps (float):
Resolution in the X or Y coordinate for the focused field (nm)
 z_steps (float):
Resolution in the axial coordinate (Z) for the focused field (nm)
 x_range (float):
Field of view in the X or Y coordinate in which the focused field is calculated (nm)
 z_range (float):
field of view in the axial coordinate (z) in which the focused field is calculated (nm)
 I_0 (float):
Incident field intensity (mW/cm^2)
 L (float, optional):
Distance between phase mask and objective lens (mm), only used if propagation=True
 R (float, optional):
Phase mask radius (mm), only used if propagation=True
 ds (array, optional):
Thickness of each interface in the multilayer system (nm), must be given as a numpy array with the first and last values a np.inf. Only used if multilayer=True
 z_int (float, optional):
Axial position of the interphase. Only used if multilayer=True
 figure_name (string, optional):
Name for the images of the field. Also used as saving name if using the UI

sim.
VP
(propagation=False, multilayer=False, NA=1.4, n=1.5, h=3, w0=5, wavelength=640, gamma=45, beta=90, z=0, x_steps=5, z_steps=8, x_range=1000, z_range=2000, I0=1, L='', R='', ds='', z_int='', figure_name='')¶ Simulate the field obtained by focusing a gaussian beam modulated by a VP mask
 Args:
 propagation (bool):
True for calculating and ploting the field incident on the lens by fraunhofer’s difractin formula, in which case R and L are needed; False for calculating the field incident on the lens while neglecting the propagation
 multilayer (bool):
True for calculating the focused field with an multilayer, in which case ds and z_int are needed; Flase for calculating the field obtained without an multilayer
 Parameters:
NA,n,h,w0,wavelength,gamma,beta,z,x_steps,z_steps,x_range,z_range,I0,L,R,ds,z_int,figure_name: Simulation parameters, defined as part of the ‘parameters’ array
 Returns:
 arrays:
ex_XZ,ey_XZ,ez_XZ,ex_XY,ey_XY,ez_XY, each one is a matrix with the amplitude of each cartesian component on the XZ plane (ex_XZ,ey_XZ,ez_XZ) or on the XY plane (ex_XY,ey_XY,ez_XY)
Each index of the matrixes corresponds to a different pair of coordinates, for example:
ex_XZ[z,x] with z each index of the coordinates np.linspace(z_range/2,z_range/2,2*int(z_range/z_steps/2)) and x each index for np.linspace(x_range/2**0.5,x_range/2**0.5,2*int(x_range/x_steps/2**0.5)) in which the field is calculated.
IMPORTANT: It is worth noting the x positions are sqrt(2) times higher to allow a square field of view for the XY plane (because the maximum radial position to be calculated for the XY plane is higher than the maximum x or y position)
ex_XZ[y,x2] with y each index of the coordinates np.linspace(x_range/2,x_range/2,2*int(x_range/x_steps/2)) and x each index for np.linspace(x_range/2,x_range/2,2*int(x_range/x_steps/2)) in which the field is calculated
The XZ plane is given by y=0 and the XY plane by the z parameter
The unit of the intensity is the same as the initial intensity (mW/cm^2)

sim.
no_mask
(propagation=False, multilayer=False, NA=1.4, n=1.5, h=3, w0=5, wavelength=640, gamma=45, beta=90, z=0, x_steps=5, z_steps=8, x_range=1000, z_range=2000, I0=1, L='', R='', ds='', z_int='', figure_name='')¶ Simulate the field obtained by focusing a gaussian beam without being modulated in phase
 Args:
 propagation (bool):
Kept for homogeneity of the functions, True yields no diference from False, where the field incident on the lens is calculated by neglecting the propagation towards the objective lens
 multilayer (bool):
True for calculating the focused field with an multilayer, in which case ds and z_int are needed; Flase for calculating the field obtained without an multilayer
 Parameters:
NA,n,h,w0,wavelength,gamma,beta,z,x_steps,z_steps,x_range,z_range,I0,L,R,ds,z_int,figure_name: Simulation parameters, defined as part of the ‘parameters’ array
 Returns:
 arrays:
ex_XZ,ey_XZ,ez_XZ,ex_XY,ey_XY,ez_XY, each one is a matrix with the amplitude of each cartesian component on the XZ plane (ex_XZ,ey_XZ,ez_XZ) or on the XY plane (ex_XY,ey_XY,ez_XY)
Each index of the matrixes corresponds to a different pair of coordinates, for example:
ex_XZ[z,x] with z each index of the coordinates np.linspace(z_range/2,z_range/2,2*int(z_range/z_steps/2)) and x each index for np.linspace(x_range/2**0.5,x_range/2**0.5,2*int(x_range/x_steps/2**0.5)) in which the field is calculated.
IMPORTANT: It is worth noting the x positions are sqrt(2) times higher to allow a square field of view for the XY plane (because the maximum radial position to be calculated for the XY plane is higher than the maximum x or y position)
ex_XZ[y,x2] with y each index of the coordinates np.linspace(x_range/2,x_range/2,2*int(x_range/x_steps/2)) and x each index for np.linspace(x_range/2,x_range/2,2*int(x_range/x_steps/2)) in which the field is calculated
The XZ plane is given by y=0 and the XY plane by the z parameter
The unit of the intensity is the same as the initial intensity (mW/cm^2)

sim.
custom
(entrance_field, custom_mask, propagation=False, multilayer=False, NA=1.4, n=1.5, h=3, w0=5, wavelength=640, gamma=45, beta=90, z=0, x_steps=5, z_steps=8, x_range=1000, z_range=2000, I0=1, L='', R='', ds='', z_int='', figure_name='', divisions_theta=200, divisions_phi=200, plot_Ei=True)¶ Simulate the field obtained by focusing a gaussian beam modulated by a custom phase mask
 Args:
 entrance_field (function):
Analytical function that defines the entrance beam, must be a function of the 5 internal variables: rho, phi, w0, f and k, with:
rho: Radial coordinate from 0 to the aperture radius of the objective.
phi: Azimutal coordinate from 0 to 2pi.
w0: Radius of the incident gaussian beam.
f: Focal distane of the objective lens
k: Wavenumber in the objective lens medium IMPORTANT:(mm)
 custom_mask (2D array or function):
Custom phase mask to be simulated
If 2D array: Complex amplitude of the custom phase mask for each value of the spherical coordinates theta and phi: custom_mask[phi_position,theta_position] for phi_position a value in np.linspace(0,2*np.pi,divisions_phi) and theta_position a value in np.linspace(0,alpha,divisions_theta)
If function: Analytical function that defines the phase mask, must be a function of the 5 internal variables: rho, phi, w0, f and k defined previously.
The real part defines the amplitude modulation of the entrance beam
 propagation (bool):
True for calculating and ploting the field incident on the lens by fraunhofer’s difractin formula, in which case R and L are needed; False for calculating the field incident on the lens while neglecting the propagation
 multilayer (bool):
True for calculating the focused field with an multilayer, in which case ds and z_int are needed; Flase for calculating the field obtained without an multilayer
 Parameters:
NA,n,h,w0,wavelength,gamma,beta,z,x_steps,z_steps,x_range,z_range,I0,L,R,ds,z_int,figure_name: Simulation parameters, defined as part of the ‘parameters’ array
 Returns:
 arrays:
ex_XZ,ey_XZ,ez_XZ,ex_XY,ey_XY,ez_XY, each one is a matrix with the amplitude of each cartesian component on the XZ plane (ex_XZ,ey_XZ,ez_XZ) or on the XY plane (ex_XY,ey_XY,ez_XY)
Each index of the matrixes corresponds to a different pair of coordinates, for example:
ex_XZ[z,x] with z each index of the coordinates np.linspace(z_range/2,z_range/2,2*int(z_range/z_steps/2)) and x each index for np.linspace(x_range/2**0.5,x_range/2**0.5,2*int(x_range/x_steps/2**0.5)) in which the field is calculated.
IMPORTANT: It is worth noting the x positions are sqrt(2) times higher to allow a square field of view for the XY plane (because the maximum radial position to be calculated for the XY plane is higher than the maximum x or y position)
ex_XZ[y,x2] with y each index of the coordinates np.linspace(x_range/2,x_range/2,2*int(x_range/x_steps/2)) and x each index for np.linspace(x_range/2,x_range/2,2*int(x_range/x_steps/2)) in which the field is calculated
The XZ plane is given by y=0 and the XY plane by the z parameter
The unit of the intensity is the same as the initial intensity (mW/cm^2)
Internal functions¶
IMPORTANT: In these functions the wavelength corresponds to the wavelegth on the lens’s medium unless stated otherwise
Functions for a VP mask¶
Functions for the simulation of the foci obtained by a VP mask

VP_functions.
VP_integration
(alpha, n, f, w0, wavelength, x_steps, z_steps, x_range, z_range)¶ Generate the II arrays, which are the result of the integration for different positions along the radius and z
This matrixes are later used to calculate the focused field
 Args:
 alpha:
semiangle of aperture
 wavelength:
wavelength in the medium (equals wavelength in vacuum/n)
 x_steps:
resolution in the x or y coordinate (nm)
 z_steps:
resolution in the axial coordinate,z (nm)
 x_range:
field of view in the x or y coordinate in which the field is calculated (nm)
 z_range:
field of view in the axial coordinate, z, in which the field is calculated (nm)
The other parameters are specified in sim
 Returns:
 (arrays):
II1,II2,II3

VP_functions.
VP_fields
(II1, II2, II3, II4, II5, wavelength, I0, gamma, beta, x_steps, z_steps, x_range, z_range, phip0, n, f, zp0)¶ Given the II matrixes calculate the field on the focus
 Args:
 phip0:
Gives an azimutal offset for the XZ plane calculus
 wavelength:
wavelength in the medium (equals wavelength in vacuum/n)
 zp0:
axial position of the XY plane
The other parameters are specified in sim.py
 Returns:
 arrays:
Ex,Ey,Ez,Ex2,Ey2,Ez2, each one is a matrix with the amplitude of each cartesian component on the XZ plane (Ex,Ey,Ez) or on the XY plane (Ex2,Ey2,Ez2)
Each index of the matrixes corresponds to a different pair of coordinates, for example:
ex[z,x] with z each index of the coordinates np.linspace(z_range/2,z_range/2,2*int(z_range/z_steps/2)) and x each index for np.linspace(x_range/2**0.5,x_range/2**0.5,2*int(x_range/x_steps/2**0.5)) in which the field is calculated
ex2[y,x2] with y each index of the coordinates np.linspace(x_range/2,x_range/2,2*int(x_range/x_steps/2)) and x each index for np.linspace(x_range/2,x_range/2,2*int(x_range/x_steps/2)) in which the field is calculated
The XZ plane is given by y=0 and the XZ plane by z=zp0
The radial field of view in the XZ plane is sqrt(2) times bigger to allow a square field of view for the XY plane (the maximum radial position is higher than the maximum x or y position)

VP_functions.
VP_fraunhofer
(gamma=45, beta= 90, steps=500, R=5, L=100, I_0=1, wavelength=640, FOV=11, w0=5, limit=2000, div=1, plot=True, folder='', figure_name='')¶ Calculate and plot the field inciding on the lens by Fraunhofer’s difraction formula
 Args:
 limit:
Ammount of iterations the scipy.quad command can do
 div:
Ammount of divisions in which the integration is divided in order to avoid the scipy.quad function from failing to converge
 plot (bool):
True plots the inciding field’s intensity and amplitude
 wavelength:
wavelength given in the medium (equals wavelength in vacuum/n)
The other parameters are specified in sim.py
 Returns:
 array:
E_rho: the inciding amplitude along the radial coordinate for an x polarized beam
 arrays:
Ex and Ey, the x and y components of this amplitude, in a matrix over the x and y coordinates so it can be ploted easily

VP_functions.
VP_integration_with_propagation
(alpha, n, f, radius_VP, wavelength, zp0, z_steps, x_steps, x_range, laser_width, E_rho, div)¶ Given the inciding field E_rho, which only depends on the radial coordinate, generate the I matrixes, which are the same as in VP_integration
Since the calculus takes a long time, only the field along the XY plane is calculated
wavelength is given in the medium (equals wavelength in vacuum/n)
The other parameters are specified in sim.py

VP_functions.
VP_fields_with_propagation
(I1, I2, I3, I4, I5, wavelength, I0, gamma, beta, x_steps, z_steps, x_range, phip0, n, f, zp0)¶ Given the I matrixes calculate the field on the focus Since the calculus takes a long time, only the field along the XY plane is calculated parameter phip0 has no purpose, is only left to have the same variables for the functions wavelength is given in the medium (equals wavelength in vacuum/n)
The other parameters are specified in sim.py
 Returns:
 arrays:
Ex,Ey,Ez,Ex2,Ey2,Ez2, each one is a matrix with the amplitude of each cartesian component on the XZ plane (Ex,Ey,Ez) or on the XY plane (Ex2,Ey2,Ez2)
Each index of the matrixes corresponds to a different pair of coordinates, for example:
ex[z,x] with z each index of the coordinates np.linspace(z_range/2,z_range/2,2*int(z_range/z_steps/2)) and x each index for np.linspace(x_range/2**0.5,x_range/2**0.5,2*int(x_range/x_steps/2**0.5)) in which the field is calculated
ex2[y,x2] with y each index of the coordinates np.linspace(x_range/2,x_range/2,2*int(x_range/x_steps/2)) and x each index for np.linspace(x_range/2,x_range/2,2*int(x_range/x_steps/2)) in which the field is calculated
The XZ plane is given by y=0 and the XZ plane by z=zp0
The radial field of view in the XZ plane is sqrt(2) times bigger to allow a square field of view for the XY plane (the maximum radial position is higher than the maximum x or y position)
Functions for no mask¶
Functions for the simulation of the field obtained by focuisng a gaussian beam

no_mask_functions.
no_mask_integration
(alpha, n, f, w0, wavelength, x_range, z_range, z_steps, r_steps)¶ Generate the II arrays, which are the result of the integration for different positions along the radius and z
This matrixes are later used to calculate the focused field
 Args:
 alpha:
semiangle of aperture
 wavelength:
wavelength in the medium (equals wavelength in vacuum/n)
 r_steps:
resolution in the x or y coordinate (nm)
 z_steps:
resolution in the axial coordinate,z (nm)
 x_range:
field of view in the x or y coordinate in which the field is calculated (nm)
 z_range:
field of view in the axial coordinate, z, in which the field is calculated (nm)
The other parameters are specified in sim
 Returns:
 (arrays):
II1,II2,II3

no_mask_functions.
no_mask_fields
(II1, II2, II3, wavelength, I0, beta, gamma, z_steps, r_steps, x_range, z_range, phip0, n, f, zp0)¶ Given the II matrixes calculate the field on the focus
 Args:
 phip0:
Azimutal offset for the XZ plane calculus
 wavelength:
wavelength given in the medium (equals wavelength in vacuum/n)
 zp0:
axial position of the XY plane
The other parameters are specified in sim
 Returns:
 arrays:
Ex,Ey,Ez,Ex2,Ey2,Ez2, each one is a matrix with the amplitude of each cartesian component on the XZ plane (Ex,Ey,Ez) or on the XY plane (Ex2,Ey2,Ez2)
Each index of the matrixes corresponds to a different pair of coordinates, for example:
ex[z,x] with z each index of the coordinates np.linspace(z_range/2,z_range/2,2*int(z_range/z_steps/2)) and x each index for np.linspace(x_range/2**0.5,x_range/2**0.5,2*int(x_range/r_steps/2**0.5)) in which the field is calculated
ex2[y,x2] with y each index of the coordinates np.linspace(x_range/2,x_range/2,2*int(x_range/r_steps/2)) and x each index for np.linspace(x_range/2,x_range/2,2*int(x_range/r_steps/2)) in which the field is calculated
The XZ plane is given by y=0 and the XZ plane by z=zp0
The radial field of view in the XZ plane is sqrt(2) times bigger to allow a square field of view for the XY plane (the maximum radial position is higher than the maximum x or y position)
Functions for a custom phase mask¶
Functions for the simulation of the foci obtained by an arbitrary phase mask

custom_mask_functions.
generate_incident_field
(maskfunction, alpha, f, divisions_phi, divisions_theta, gamma, beta, w0, I0, wavelength)¶ Generate a matrix for the field X and Y direction of the incident field on the lens, given the respective maskfunction
 Args:
 maskfunction (function):
Analytical function that defines the phase mask, must be a function of the 5 internal variables: rho, phi, w0, f and k, with:
rho: Radial coordinate from 0 to the aperture radius of the objective.
phi: Azimutal coordinate from 0 to 2pi.
w0: Radius of the incident gaussian beam.
f: Focal distane of the objective lens (mm)
k: Wavenumber in the objective lens medium (mm)
The real part defines the amplitude of the incident field
 divisions_phi,divisions_theta:
Number of divisions in the phi and theta coordinates to use the 2D integration for the calculation of the focused field
The rest of the parameters are specified in sim
 Returns:
 arrays:
ex_lens,ey_lens
This arrays have the value of the amplitude of the incident field for each value of theta and phi: ex_lens[phi_position,theta_position] for phi_position a value in np.linspace(0,2*np.pi,divisions_phi) and theta_position a value in np.linspace(0,alpha,divisions_theta)

custom_mask_functions.
generate_rotated_incident_field
(maskfunction, alpha, f, divisions_phi, divisions_theta, gamma, beta, w0, I0, wavelength)¶ Generate a matrix for the field X and Y direction of the incident field on the lens evaluated at phi180º for ex_lens and at phi270º for ey_lens, given the respective maskfunction
 Args:
 maskfunction (function):
Analytical function that defines the phase mask, must be a function of the 5 internal variables: rho, phi, w0, f and k, with:
rho: Radial coordinate from 0 to the aperture radius of the objective.
phi: Azimutal coordinate from 0 to 2pi.
w0: Radius of the incident gaussian beam.
f: Focal distane of the objective lens (mm)
k: Wavenumber in the objective lens medium (mm)
The real part defines the amplitude of the incident field
 divisions_phi,divisions_theta:
Number of divisions in the phi and theta coordinates to use the 2D integration for the calculation of the focused field
The rest of the parameters are specified in sim
 Returns:
 arrays:
ex_lens,ey_lens
This arrays have the value of the amplitude of the incident field for each value of theta and phi: ex_lens[phi_position,theta_position] for phi_position a value in np.linspace(0,2*np.pi,divisions_phi) and theta_position a value in np.linspace(0,alpha,divisions_theta)

custom_mask_functions.
plot_in_cartesian
(Ex, Ey, r_range, alpha, f, figure_name)¶ Plot the fields Ex and Ey, who are described in the same coordinates as ex_lens and ey_lens. To do so the field in the closest cartesian coordinates for each position is calculated
 Args:
 r_range:
Radial distance in which to plot the field (total distance in x and y is 2*r_range)
 alpha:
Semiangle of aperture of the objective lens, given by alpha=np.arcsin(NA / n)
 f:
Focal distance of the objective lens, given by the sine’s law: f=h*n/NA
 figure_name:
Name for the ploted figures
 Returns:
I_cartesian, Ex_cartesian, Ey_cartesian: Intensity and amplitude of the incident field calculated in cartesian coordinates

custom_mask_functions.
custom_mask_objective_field
(h, gamma, beta, divisions_theta, divisions_phi, N_rho, N_phi, alpha, focus, custom_field_function, R, L, I0, wavelength, w0, fig_name, plot=True)¶ Calculate the incident field on the objective by fraunhofer’s difraction formula for a custom phase mask
The resultant matrix Ex and Ey are returned in polar coordinates (each row is a different value of phi and each column a different rho)
 Args:
 N_rho and N_phi:
Number of divisions for the calclation of the 2D integral in rho and phi respectively (this are not the coordinates in which the field is calculated)
 divisions_theta,divisions_phi:
Number of divisions for the field indicent on the lens on the theta and phi coordinates respectively. This field is later used to calculate the focused field.
The rest of the parameters are specified in sim
 Returns:
 arrays:
Ex, Ey: Amplitude of the incident field calculated on the pair of coordinates: theta_values=np.linspace(0,alpha,divisions_theta), phip_values=np.linspace(0,2*np.pi,divisions_phi)
by using sine’s law this corresponds to the radial positions rhop_values=np.sin(theta_values)*focus
 arrays:
I_cartesian, Ex_cartesian, Ey_cartesian: Intensity and amplitude of the incident field calculated in cartesian coordinates

custom_mask_functions.
custom_mask_focus_field_XY
(ex_lens, ey_lens, alpha, h, wavelength, zp0, resolution_x, divisions_theta, divisions_phi, x_range, countdown=True, x0=0, y0=0, field_is_already_rotated=False)¶ 2D integration to calculate the field focused by a high aperture lens on the XY plane
 Args:
 ex_lens,ey_lens:
X and Y component of the incident field. This arrays have the value of the amplitude of the incident field for each value of theta and phi: ex_lens[phi_position,theta_position] for phi_position a value in np.linspace(0,2*np.pi,divisions_phi) and theta_position a value in np.linspace(0,alpha,divisions_theta)
 zp0:
Axial position for the XY plane (given by z=zp0)
 resolution_x:
Resolution for the field at the focus, the same for x and y
 divisions_theta,divisions_phi:
Resolution for the 2D calculus (must be the same as the size of ex_lens and ey_lens)
 wavelength:
Wavelength (nm) in the medium (equals wavelength in vacuum/n)
 x0, y0:
Used for centering the XY field at an x0, y0 position (nm)
 countdown (bool):
True for shoing a progress bar with the time elapsed and expected to finish the calculation, True not recomended if this function is to be used many times
 field_is_already_rotated (bool):
Internal variable used to tell if the field is already evaluated at phi180º and phi270º
The rest of the parameters are specified in sim
 Returns:
 arrays:
ex,ey,ez: Cartesian components of the focused field on the XY plane, given by z=zp0

custom_mask_functions.
custom_mask_focus_field_XZ_XY
(ex_lens, ey_lens, alpha, h, wavelength, z_range, resolution_z, zp0, resolution_x, divisions_theta, divisions_phi, x_range, x0=0, y0=0, plot_plane='X', field_is_already_rotated=False)¶ 2D integration to calculate the field focused by a high aperture lens on the XY and XZ plane
 Args:
 ex_lens,ey_lens:
X and Y component of the incident field. This arrays have the value of the amplitude of the incident field for each value of theta and phi: ex_lens[phi_position,theta_position] for phi_position a value in np.linspace(0,2*np.pi,divisions_phi) and theta_position a value in np.linspace(0,alpha,divisions_theta)
 zp0:
Axial position for the XY plane (given by z=zp0)
 resolution_x:
Resolution for the field at the focus, in the x and y coordinates
 resolution_z:
Number of pixels for the field at the focus, in the z coordinate
 divisions_theta,divisions_phi:
Resolution for the 2D calculus (must be the same as the size of ex_lens and ey_lens)
 wavelength:
Wavelength (nm) in the medium (equals wavelength in vacuum/n)
 x0, y0:
Used for centering the XY field at an x0, y0 position (nm)
 plot_plane (string):
Available values: ‘X’ or ‘Y’, select to plot the ZX or the ZY plane respectivelly
 field_is_already_rotated (bool):
Internal variable used to tell if the field is already evaluated at phi180º and phi270º
The rest of the parameters are specified in sim
 Returns:
 arrays:
Ex_XZ,Ey_XZ,Ez_XZ,Ex_XY,Ey_XY,Ez_XY, each one is a matrix with the amplitude of each cartesian component on the XZ plane (Ex,Ey,Ez) or on the XY plane (Ex2,Ey2,Ez2)
Each index of the matrixes corresponds to a different pair of coordinates, for example:
Ex_XZ[z,x] with z each index of the coordinates np.linspace(z_range/2,z_range/2,resolution_z) and x each index for np.linspace(x_range/2**0.5,x_range/2**0.5,resolution_x) in which the field is calculated, the x range is sqrt(2) times longer for consistency with the VP and no_mask functions
Ex_XY[y,x2] with y each index of the coordinates np.linspace(x_range/2,x_range/2,resolution_x) and x each index for np.linspace(x_range/2,x_range/2,resolution_x) in which the field is calculated
The XZ plane is given by y=0 and the XZ plane by z=zp0
Functions for simulation of an interface¶

interface.
interface_custom_mask_focus_field_XY
(n_list, d_list, ex_lens, ey_lens, alpha, h, wavelength, z_int, zp0, resolution_x, divisions_theta, divisions_phi, x_range, countdown=True, x0=0, y0=0, field_is_already_rotated=False)¶ 2D integration to calculate the field focused by a high aperture lens on the XY plane with an interface
 Args:
 n_list (array):
array with the refraction index of each medium of the multilayer system
 d_list (array):
Thickness of each interface in the multilayer system (nm), must be given as a numpy array with the first and last values a np.inf. Only used if interface=True
 ex_lens,ey_lens (array):
X and Y component of the incident field. Each position must be given in the coordinates
This arrays must have the value of the amplitude of the incident field for each value of theta and phi. Example: ex_lens[phi_position,theta_position] for phi_position a value in np.linspace(0,2*np.pi,divisions_phi) and theta_position a value in np.linspace(0,alpha,divisions_theta)
 z_int:
Axial position in wich the interface is located
 zp0:
Axial position for the XY plane (given by z=zp0)
 resolution_x:
Number of pixels for the field at the focus, in the x and y coordinates
 divisions_theta,divisions_phi:
Resolution for the 2D calculus (must be the same as the size of ex_lens and ey_lens)
 wavelength:
Wavelength (nm) in the vacuum
 x0, y0:
Used for centering the XY field at an x0, y0 position
 countdown (bool):
If True, a progress bar is shown with the use of the tqdm package
The rest of the parameters are specified in sim
 Returns:
 arrays:
ex,ey,ez: Cartesian components of the focused field on the XY plane, given by z=zp0

interface.
interface_custom_mask_focus_field_XZ_XY
(n_list, d_list, ex_lens, ey_lens, alpha, h, wavelength, z_int, z_range, resolution_z, zp0, resolution_x, divisions_theta, divisions_phi, x_range, x0=0, y0=0, z0=0, plot_plane='X', field_is_already_rotated=False)¶ 2D integration to calculate the field focused by a high aperture lens on the XY and XZ planes with an interface
 Args:
 n_list (array):
array with the refraction index of each medium of the multilayer system
 d_list (array):
Thickness of each interface in the multilayer system (nm), must be given as a numpy array with the first and last values a np.inf. Only used if interface=True
 ex_lens,ey_lens (array):
X and Y component of the incident field. Each position must be given in the coordinates
This arrays must have the value of the amplitude of the incident field for each value of theta and phi. Example: ex_lens[phi_position,theta_position] for phi_position a value in np.linspace(0,2*np.pi,divisions_phi) and theta_position a value in np.linspace(0,alpha,divisions_theta)
 z_int:
Axial position in wich the interface is located
 zp0:
Axial position for the XY plane (given by z=zp0)
 resolution_x:
Number of pixels for the field at the focus, in the x and y coordinates
 resolution_z:
Number of pixels for the field at the focus, in the z coordinate
 divisions_theta,divisions_phi:
Resolution for the 2D calculus (must be the same as the size of ex_lens and ey_lens)
 wavelength:
Wavelength (nm) in the vacuum
 x0, y0:
Used for centering the XY plane at an x0, y0 position (nm)
 z0:
Used for centering the XZ plane at an z0 position (nm)
 plot_plane (string):
Available values: ‘X’ or ‘Y’, select to plot the ZX or the ZY plane respectivelly
The rest of the parameters are specified in sim
 Returns:
 arrays:
Ex_XZ,Ey_XZ,Ez_XZ,Ex_XY,Ey_XY,Ez_XY, each one is a matrix with the amplitude of each cartesian component on the XZ plane (ex_XZ,ey_XZ,ez_XZ) or on the XY plane (ex_XY,ey_XY,ez_XY)
Each index of the matrixes corresponds to a different pair of coordinates, for example:
Ex_XZ[z,x] with z each index of the coordinates np.linspace(z_range/2,z_range/2,resolution_z) and x each index for np.linspace(x_range/2**0.5,x_range/2**0.5,resolution_x) in which the field is calculated, the x range is sqrt(2) times longer for consistency with the VP and no_mask functions
Ex_XY[y,x2] with y each index of the coordinates np.linspace(x_range/2,x_range/2,resolution_x) and x each index for np.linspace(x_range/2,x_range/2,resolution_x) in which the field is calculated
The XZ plane is given by y=0 and the XZ plane by z=zp0