A zoneplate is a lens based on diffraction of light. We will use PyPropagate to calculate the stationary field in a binary zoneplate and its focal point in paraxial approximation.
In [1]:
from pypropagate import *
%matplotlib inline
In [2]:
zoneplate_settings = presets.settings.create_paraxial_wave_equation_settings()
s = zoneplate_settings.symbols
In [3]:
zp = zoneplate_settings.create_category("zoneplate",info="Zoneplate parameters")
zp.create_symbol("F",info="Focal length")
zp.create_symbol("z_F",zp.F+s.sz/2,info="Focal spot")
zp.create_symbol("n1",info="Index of refraction of the zoneplate")
zp.create_symbol("n2",info="Index of refraction")
zp.create_symbol("N",info="Number of zones",type=pc.Types.Integer)
zp.create_symbol("i",info="Zone index",type=pc.Types.Integer)
zp.create_function("r",(zp.i,),info="radius of zone i")
zp.create_symbol("R",zp.r.function(2*zp.N+1),info="radius of the zoneplate");
We define a formula for the zone radii.
In [4]:
zp.r = pc.sqrt(zp.i*s.wavelength*zp.F + zp.i**2*s.wavelength**2 / 4)
zoneplate_settings.get_definition(zp.r)
Out[4]:
We set the refractive indices of the zones to Tantalum and Vacuum. There will be $50$ Zones with a focal point at $1\text{mm}$.
In [5]:
zp.F = 1*units.mm
zp.n1 = presets.medium.create_material('Ta',zoneplate_settings)
zp.n2 = 1
zp.N = 50
The simulation box is defined large enough to contain all zones with some extra space we call padding. The simulation box size in $z$ direction corresponds to the thickness of the zoneplate.
In [6]:
padding = 0.2
s.zmin = 0
s.zmax = 20*units.um
s.Nz = 10000
s.xmax = (1+padding)*zp.R
s.xmin = -s.xmax
s.Nx = 20000
zoneplate_settings.wave_equation.set_energy(12*units.keV)
presets.boundaries.set_plane_wave_initial_conditions(zoneplate_settings)
We save the zones information in a one-dimensional binary array where the index corresponds to the distance from the center. Since we want to be able to change zoneplate parameters at any time the zoneplate array should be created as a last step before the simulation. Therefore we define the function that creates the zoneplate array as an initializer function.
In [7]:
def create_zoneplate_array(settings):
s = settings.symbols
key = settings.get_unitless(tuple(zp.keys()+[s.dx]))
# Avoid unnecessary re-calculation of the zones
cached = settings.get_cache().get('zoneplate_array')
if cached is None or cached != key:
size = settings.get_as(zp.R/s.dx+0.5,int)
ro = 0
zf = pc.numpyfy(settings.get_unitless(pc.floor(zp.r/s.dx)))
zones = zf(i_zoneplate = range(settings.get_as(2*(zp.N+1),int)))
zoneplate_array = np.ones(size,dtype=np.bool)
for ri,ro in zip(zones[::2],zones[1::2]):
zoneplate_array[ri:ro] = np.zeros(ro-ri,dtype=np.bool)
zoneplate_function = pc.array("zones",zoneplate_array)
inside_zoneplate_r = pc.equal(zoneplate_function(pc.floor(pc.sqrt(s.x**2+s.y**2)/s.dx + 0.5)),0)
s.n = pc.piecewise((zp.n1,inside_zoneplate_r),(zp.n2,True))
settings.get_cache()['zoneplate_array'] = key
zoneplate_settings.initializers['create_zoneplate_array'] = create_zoneplate_array
To verify the settings we plot the refraction index at $y=0$.
In [8]:
n = s.n.function
plot(expression_to_array(pc.unequal(n(s.x,0,s.zmin),1),zoneplate_settings),figsize=(15,4));
In [9]:
zppropagator = propagators.FiniteDifferencesCS(zoneplate_settings)
zpfield = zppropagator.run_slice()[0::,::s.sz/1000]
In [10]:
plot(zpfield,figsize=(12,8));
plot(zpfield[3.1*units.um:zp.R],figsize=(12,4))
Out[10]:
We will use the exit field of the finite difference propagator as the initial field for propagation in vacuum. Since we are not interested in the pre-focal field we will do a single step up to 80% of the focal length. To avoid reflections from the boundary we will add further padding to the simulation box.
In [11]:
def add_absorbing_padding_x(settings,size,n = None):
from pypropagate import pc,presets
if n is None:
n = presets.medium.create_material('Pb',settings)
s = settings.symbols
r = pc.Max(pc.Max(0,(s.x - s.xmax + size)/size),(-s.x + s.xmin + size)/size)
s.n = settings.get_numeric(s.n) * (1-r) + r * n
In [12]:
vacuum_settings = presets.settings.create_next_settings(zoneplate_settings)
vacuum_settings.wave_equation.n = 1
add_absorbing_padding_x(vacuum_settings,s.sx/2 - zp.R)
vacuum_settings.simulation_box.sz = 1.3 * zoneplate_settings.get_numeric(zp.z_F - s.sz/2)
vacuum_settings.simulation_box.Nz = 100000
initial = presets.boundaries.add_padding(zppropagator.get_field(),0.5)
presets.boundaries.set_initial(vacuum_settings,initial)
In [13]:
vapropagator = propagators.FiniteDifferencesCS(vacuum_settings)
vafield = vapropagator.run_slice()[::,::s.sz/1000]
In [14]:
from matplotlib.colors import LogNorm
plot(vafield[:],figsize=(10,10),norm=LogNorm(10**-10))
Out[14]:
In [15]:
from pypropagate.colormaps import fire_colormap
import matplotlib.pyplot as plt
fig,(ax1,ax2) = plt.subplots(2,figsize=(15,10),sharex=True)
plot(vafield[-0.5*units.um:0.5*units.um],ax=ax1,cmap=fire_colormap(),vmax=100)
plot(vafield[0],ax=ax2)
plt.xlim(0.1,1.3);
In [16]:
zF = zoneplate_settings.get_numeric(zp.F + s.sz/2)
plot(vafield[-0.25*units.um:0.25*units.um,zF - 0.25*units.mm:zF+0.25*units.mm],cmap=fire_colormap(),figsize=(15,5));