Simulation description

structures

structures contains the main struct object, structure generators and helper functions related to the generation of structure geometries.

structure definition class struct

class pyGDM2.structures.struct(step, geometry, material, n1, n2, normalization, n3=None, spacing=5000.0, with_radiation_correction=True)

structure container

Defines a nanostructure to perform a GDM simulation on. This includes the geometry, the material specific dielectric function(s), as well as the environment in which the structure located.

Parameters:
  • step (float) – stepsize used for discretization
  • geometry (list of tuples) – list of coordinates of the meshpoints of the geometry. Several structure generators are provided in the module structures
  • material (instance of materials class or list of such) – use material classes provided by the materials module. If a single instance is provided, the according material will be used for the entire structure
  • n2 (n1,) – complex refractive index of environmental layers “1” (substrate), “2” (environment containing structure)
  • normalization (float) – can be obtained using structures.get_normalization(). =1 for cubic mesh, =sqrt(2) for hexagonal compact mesh.
  • n3 (complex, default: value of n2) – complex refractive index of (optional) top cladding layer (layer “3”). By default, use ref.index of structure environment (layer “2”, ref.index n2)
  • spacing (float, default = 5000) – thickness of layer “2” in nm. Layers “1” and “2” are infinitely thick (“1” from -inf to 0; “3” from spacing to +inf)
  • with_radiation_correction (bool, default: True) – Adds an optional radiative correction to the self-term, like described by B. Draine, in: Astrophys. J. 333, 848–872 (1988). Using the correction usually leads to better convergence for very large discretization stepsizes (>=20nm), and has only a weak influence on the results from fine meshes.
__init__(step, geometry, material, n1, n2, normalization, n3=None, spacing=5000.0, with_radiation_correction=True)

Initialization

getNormalization(wavelength)

Returns the normalization term at “wavelength”

Calling this function will set the internal variables struct.wavelength and struct.cnorm, self.norm_nonrad and self.norm_rad

Parameters:wavelength (float) – wavelength (in nm) at which to calculate polarizabilites
Returns:cnorm – normalization term for used mesh and wavelength
Return type:complex

Notes

For details on the normalization derivation, see e.g. Girard, C., Dujardin, E., Baffou, G. & Quidant, R. Shaping and manipulation of light fields with bottom-up plasmonic structures. New J. Phys. 10, 105016 (2008).

getPolarizability(wavelength)

Returns the polarizabilities of each meshpoint at “wavelength”

Calculate the “polarizabilites” for each meshpoint using its volume (defined by step) and its dielectric function (“material”) as well as the environment refractive index (“n2”).

For the moment, only isotropic materials (scalar espilon) are supported. (anistropic tensor support is already implemented in fortran routines but not exposed to python)

Calling this function will set the internal variables struct.wavelength, struct.alpha, struct.areal (real part of alpha) and struct.aimag (imag part of alpha)

Parameters:wavelength (float) – wavelength (in nm) at which to calculate polarizabilites
Returns:alpha – complex “polarizabilities” for each meshpoint
Return type:list of complex

Notes

These are not actual polarizabilities and in particular not identical to the ones used e.g. in the coupled dipole approximation software DDSCAT. For details, see e.g.: Girard, C. Near fields in nanostructures. Reports on Progress in Physics 68, 1883-1933 (2005).

setDtype(dtypef, dtypec)

set dtype of arrays

Structure generators

rect_wire

pyGDM2.structures.rect_wire(step, L, H, W, mesh='cube', ORIENTATION=1)

Generate rectangular wire with long axis along X on X/Y plane

Parameters:
  • step (float) – stepsize in nw
  • H, W (L,) – Length, Height and Width as Nr. of steps
  • mesh (string, default: 'cube') – meshing type. ‘cube’ or ‘hex’
  • ORIENTATION (int, default: 1) – For hex. meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

list of 3d coordinates in tuples (x,y,z)

Return type:

list of tuples

rect_dimer

pyGDM2.structures.rect_dimer(step, L, H, W, G, mesh='cube', ORIENTATION=1)

Generate rectangular wire dimer with long axis along X on X/Y plane

Using rectwire model for generation of the two dimer particles.

Note : Using an even number gap will result in non-symmetric position with respect to (X=0,Y=0)

Parameters:
  • step (float) – stepsize in nw
  • H, W (L,) – Length, Height and Width as Nr. of steps
  • G (int (zero or positive)) – Gap between dimers, in Nr. of steps
  • mesh (string, default: 'cube') – meshing type. ‘cube’ or ‘hex’
  • ORIENTATION (int, default: 1) – For hex. meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

list of 3d coordinates in tuples (x,y,z)

Return type:

list of tuples

prism

pyGDM2.structures.prism(step, NSIDE, H, mesh='cube', ORIENTATION=1)

Generate prism on X/Y plane

Parameters:
  • step (float) – stepsize in nw
  • NSIDE (int) – sidelength of regular prism (in nrs of dipoles)
  • H (int) – height of prism (in nrs of dipoles)
  • mesh (string, default: 'cube') – meshing type. ‘cube’ or ‘hex’
  • ORIENTATION (int, default: 1) – For hex. meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

list of 3d coordinates in tuples (x,y,z)

Return type:

list of tuples

rhombus

pyGDM2.structures.rhombus(step, L, H, alpha=90.0, mesh='cube', ORIENTATION=1)

Generate planar rhombus structure on X/Y plane

Parameters:
  • step (float) – stepsize in nw
  • L (int) – sidelength of regular prism (in nrs of dipoles)
  • H (int) – height of prism (in nrs of dipoles)
  • alpha (float, default: 90) – top/bottom angle of rhombus (in degrees)
  • mesh (string, default: 'cube') – meshing type. ‘cube’ or ‘hex’
  • ORIENTATION (int, default: 1) – For hex. meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

list of 3d coordinates in tuples (x,y,z)

Return type:

list of tuples

hexagon

pyGDM2.structures.hexagon(step, NSIDE, H, rotate=0, mesh='cube', ORIENTATION=1)

Generate regular hexagon on X/Y plane

Parameters:
  • step (float) – stepsize in nw
  • NSIDE (int) – sidelength of regular hexagon (in nrs of dipoles)
  • H (int) – height of hexagon (in nrs of dipoles)
  • rotate (float, default: 0) – rotation angle (in degrees)
  • mesh (string, default: 'cube') – meshing type. ‘cube’ or ‘hex’
  • ORIENTATION (int, default: 1) – For hex. meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

list of 3d coordinates in tuples (x,y,z)

Return type:

list of tuples

double_hexagon

pyGDM2.structures.double_hexagon(step, NSIDE1, NSIDE2, DX, DY, H, rotate=0, mesh='cube', ORIENTATION=1)

Generate regular hexagon on X/Y plane

Parameters:
  • step (float) – stepsize in nw
  • NSIDE1 (int) – sidelength of hexagon1 (in nrs of dipoles)
  • NSIDE2 (int) – sidelength of hexagon2 (in nrs of dipoles)
  • DX (int) – X-offset of hexagon2 (in nrs of dipoles)
  • DY (int) – Y-offset of hexagon2 (in nrs of dipoles)
  • H (int) – height of hexagon (in nrs of dipoles)
  • rotate (float, default: 0) – rotation angle (in degrees)
  • mesh (string, default: 'cube') – meshing type. ‘cube’ or ‘hex’
  • ORIENTATION (int, default: 1) – For hex. meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

list of 3d coordinates in tuples (x,y,z)

Return type:

list of tuples

diabolo

pyGDM2.structures.diabolo(step, nside_upper, nside_lower, length, width, H, dir_prisms='in', mesh='cube', ORIENTATION=1)

generate diabolo structure on X/Y plane

Parameters:
  • step (float) – stepsize in nw
  • nside_upper (int) – sidelength of upper regular prism (in nrs of dipoles)
  • nside_lower (int) – sidelength of upper regular prism (in nrs of dipoles)
  • length (int) – length of bridge between prisms (in nrs of dipoles)
  • width (int) – width of bridge between prisms (in nrs of dipoles)
  • H (int) – height of prism (in nrs of dipoles)
  • dir_prisms (str, default="in") – “in” or “out”. Direction of prism tips, inwards or outwards.
  • mesh (string, default: 'cube') – meshing type. ‘cube’ or ‘hex’
  • ORIENTATION (int, default: 1) – For hex. meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

list of 3d coordinates in tuples (x,y,z)

Return type:

list of tuples

sphere

pyGDM2.structures.sphere(step, R, mesh='cube', ORIENTATION=1)

Generate sphere with radius R on X/Y plane

Parameters:
  • step (float) – stepsize in nw
  • R (float) – radius of sphere (in nrs of dipoles). Not limited to integer numbers.
  • mesh (string, default: 'cube') – meshing type. ‘cube’ or ‘hex’
  • ORIENTATION (int, default: 1) – For hex. meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

list of 3d coordinates in tuples (x,y,z)

Return type:

list of tuples

nanodisc

pyGDM2.structures.nanodisc(step, R, H, ELONGATED=0, mesh='cube', ORIENTATION=1)

Generate round nanodisc in X/Y plane. Height H, Radius R.

Parameters:
  • step (float) – stepsize in nw
  • R (float) – radius of circular crosssection (in nrs of dipoles). Not limited to integer numbers.
  • H (int) – Height of Structure (in nrs of dipoles).
  • ELONGATED (int, default: 0) – add optional elongation “bridge” between half-circles (in nrs of steps)
  • mesh (string, default: 'cube') – meshing type. ‘cube’ or ‘hex’
  • ORIENTATION (int, default: 1) – For hex. meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

list of 3d coordinates in tuples (x,y,z)

Return type:

list of tuples

nanorod

pyGDM2.structures.nanorod(step, L, R, caps='flat', mesh='cube', ORIENTATION=1)

Generate round nanorod with axis along X on X/Y plane

Parameters:
  • step (float) – stepsize in nw
  • L (int) – length of rod (in nrs of dipoles).
  • R (float) – radius of circular crosssection (in nrs of dipoles). Not limited to integer numbers.
  • caps (str, default: 'flat') – ‘flat’: flat caps, ‘round’: semispherical caps
  • mesh (string, default: 'cube') – meshing type. ‘cube’ or ‘hex’
  • ORIENTATION (int, default: 1) – For hex. meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

list of 3d coordinates in tuples (x,y,z)

Return type:

list of tuples

lshape_rect

pyGDM2.structures.lshape_rect(step, L, W, H, DELTA, mesh='cube', ORIENTATION=1)

Generate symmetric L-shaped nanoantenna in X/Y plane formed by two rectangular elements.

Parameters:
  • step (float) – stepsize in nw
  • L,W,H (int) – length,width,height of each rod (in nrs of dipoles)
  • DELTA (int (odd number)) – Gap width - only odd numbers work for symmetry reasons! (in nr of dipoles)
  • mesh (string, default: 'cube') – meshing type. ‘cube’ or ‘hex’
  • ORIENTATION (int, default: 1) – For hex. meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

list of 3d coordinates in tuples (x,y,z)

Return type:

list of tuples

lshape_rect_nonsym

pyGDM2.structures.lshape_rect_nonsym(step, L1, W1, L2, W2, H, DELTA, mesh='cube', ORIENTATION=1)

Generate symmetric L-shaped nanoantenna in X/Y plane formed by two rectangular elements.

Parameters:
  • step (float) – stepsize in nw
  • L,W,H (int) – length,width,height of each rod (in nrs of dipoles)
  • DELTA (int (odd number)) – Gap width - only odd numbers work for symmetry reasons! (in nr of dipoles)
  • mesh (string, default: 'cube') – meshing type. ‘cube’ or ‘hex’
  • ORIENTATION (int, default: 1) – For hex. meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

list of 3d coordinates in tuples (x,y,z)

Return type:

list of tuples

lshape_round

pyGDM2.structures.lshape_round(step, L, W, H, DELTA, RAD, mesh='cube', ORIENTATION=1)

Generate symmetric L-shaped nanoantenna in X/Y plane formed by two rectangular elements.

Parameters:
  • step (float) – stepsize in nm
  • L,W,H (int) – length,width,height of each rod (in nrs of dipoles)
  • DELTA (int (odd number)) – Gap width - only odd numbers work for symmetry reasons! (in nr of dipoles)
  • RAD (float) – radius of curvature in steps (nr of dipoles)
  • mesh (string, default: 'cube') – meshing type. ‘cube’ or ‘hex’
  • ORIENTATION (int, default: 1) – For hex. meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

list of 3d coordinates in tuples (x,y,z)

Return type:

list of tuples

split_ring

pyGDM2.structures.split_ring(step, R, W, H, G=-1, alphaG=0, mesh='cube', ORIENTATION=1)

Generate splitring structure on X/Y plane

If G (gap) is set: close the resonator up to “gap” (G=)

Parameters:
  • step (float) – stepsize in nw
  • R,W,H (int,int,int) – Dimensions: Radius, linewidth and height of structure (in nr. of dipoles)
  • G ((optional) int, default: -1) – Gap width in numbers of steps. If G == -1, alphaG is taken instead
  • alphaG ((optional) float, default: 0) – Gap width as angle in degrees. Allowed angles: Between 0 and 180. ‘0’: entirely closed
  • mesh (string, default: 'cube') – meshing type. ‘cube’ or ‘hex’
  • ORIENTATION (int, default: 1) – For hex. meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

list of 3d coordinates in tuples (x,y,z)

Return type:

list of tuples

rect_split_ring

pyGDM2.structures.rect_split_ring(step, L1, L2, H, W, G=False, mesh='cube', ORIENTATION=1)

Generate rectangular splitring structure on X/Y plane

If G (gap) is set: close the resonator up to “gap” (G=)

Parameters:
  • step (float) – stepsize in nw
  • H, W (L1,L2,) – Length (L1:X,L2:Y), Height and Width in Nr. of steps
  • G ((optional) int, default: False) –
    Gap of a closed splitring.
    If False: totally open splitring, If ‘0’: entirely closed
  • mesh (string, default: 'cube') – meshing type. ‘cube’ or ‘hex’
  • ORIENTATION (int, default: 1) – For hex. meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

list of 3d coordinates in tuples (x,y,z)

Return type:

list of tuples

Mesher

_meshCubic

pyGDM2.structures._meshCubic(XRANGE, YRANGE, ZRANGE, condition)

mesh on a cubic grid

Mesh a 3d structure on a cubic grid within given spacial limits using a boolean selection function

Parameters:
  • XRANGE (tuple) – tuple of 2 ints (x0, x1): Indices of x-coordinates on mesh
  • YRANGE (tuple) – tuple of 2 ints (y0, y1): Indices of y-coordinates on mesh
  • ZRANGE (tuple) – tuple of 2 ints (z0, z1): Indices of z-coordinates on mesh
  • condition (function) – function of mesh-coord.-indices, whether to place a meshpoint or not: ‘func(xi,yi,zi)’, must returns a Boolean
Returns:

list of 3d coordinate-index tuples

Return type:

list of tuples

_meshHexagonalCompact

pyGDM2.structures._meshHexagonalCompact(XRANGE, YRANGE, ZRANGE, condition, ORIENTATION=1)

mesh on a hexagonal compact grid

Parameters:
  • XRANGE (tuple) – tuple of 2 ints (x0, x1): Indices of x-coordinates on mesh
  • YRANGE (tuple) – tuple of 2 ints (y0, y1): Indices of y-coordinates on mesh
  • ZRANGE (tuple) – tuple of 2 ints (z0, z1): Indices of z-coordinates on mesh
  • condition (function) – function of mesh-coord.-indices, whether to place a meshpoint or not: ‘func(xi,yi,zi)’, must returns a Boolean
  • ORIENTATION (int, default: 1) – For hex. meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

list of 3d coordinate-index tuples

Return type:

list of tuples

Other functions

image_to_struct

pyGDM2.structures.image_to_struct(img_name, nm_per_pixel, stepsize, H, threshold=100, useDarkPixel=True, returnImage=False, mesh='cube', ORIENTATION=1)

Convert an image to a planar structure

Might be useful with SEM images

Parameters:
  • img_name (string) – path to image-file
  • nm_per_pixel (float) – number of nanometers corresponding to one pixel in image
  • stepsize (float) – target stepsize
  • H (int) – Height of final structure in numbers of stepsize
  • threshold (float, (default: 100)) – threshold value between [0, 255] to declare a pixel as structure all brighter pixels will be used
  • useDarkPixel (bool, default: True) – if False, use bright pixels as structure (below threshold), instead of the dark ones
  • returnImage (bool, default: False) – if True, returns a 2D numpy array corresponding to the image AND the structure
  • mesh (string, default: 'cube') – meshing type. ‘cube’ or ‘hex’
  • ORIENTATION (int, default: 1) – For hexagonal meshing only. ‘1’ or ‘2’. Orientation of Hexagons.
Returns:

  • struct (list of tuples) – list of 3d coordinates in tuples (x,y,z)
  • if returnImage==True – returns tuple: (img_array, struct)

struct_to_image

pyGDM2.structures.struct_to_image(STRUCT, RGB=True, minimum_size=-1, projection='XY')

generate image-compatible array from structure geometry

Parameters:
  • struct (list of tuple or core.simulation) – list of coordinate tuples or instance of core.simulation
  • RGB (bool, default: True) – If True, generates an RGB array. If False, generate a 2D array with simply 0 and 1 as entries (0: no material, 1: material)
  • minimum_size (int, default: -1) – minimum canvas size in pixels. -1: use structure extensions. If >0: enlarge canvas and center structure if structure extensions smaller than image.
  • projection (str, default: "XY") – 2D plane for projection. One of [‘XY’, ‘XZ’, ‘YZ’]
Returns:

  • np.array for image conisting of rgb tuples for each pixel.
  • Plot e.g. using matplotlib’s imshow.
  • Save to bitmap or png file e.g. using scipy.misc.imsave.

get_normalization

pyGDM2.structures.get_normalization(mesh)

Provide normalization factor for mesh type

Parameters:mesh (string) – mesh used for structure generator
Returns:normalization – normalization factor corresponding to mesh
Return type:float

center_struct

pyGDM2.structures.center_struct(STRUCT, returnOffsets=False)

Center a structure on the x/y plane around (0,0)

Parameters:
  • STRUCT (list of tuples, or struct-object) – list of (x,y,z) coordinate tuples, as generated by the functions of the structure module; or struct or simulation object (which includes the geometry information)
  • returnOffsets (bool, default: False) – if True, return tuple of applied (X,Y)-offset together with structure
Returns:

    • list of 3d coordinates as tuples (x,y,z); or struct or – simulation object with adapted geometry.
  • - optionally (Offset-tuple as second return value)

rotate_XY

pyGDM2.structures.rotate_XY(STRUCT, ALPHA)

Rotate a structure around the Z-axis on the X/Y plane

Parameters:
  • STRUCT (list of tuples) – list of (x,y,z) coordinate tuples, as generated by the functions of the structure module
  • ALPHA (float) – rotation angle in degrees
Returns:

list of 3d coordinates in tuples (x,y,z)

Return type:

list of tuples

materials

materials is a collection of classes defining the dispersion of materials.

fromFile class

class pyGDM2.materials.fromFile(refindex_file, unit_wl='micron', interpolate_order=1)

tabulated dispersion

Use tabulated data provided from textfile for the complex material refractive index

Parameters:
  • refindex_file (str) –

    path to text-file with the tabulated refractive index (3 whitespace separated columns: #1 wavelength, #2 real(n), #3 imag(n))

    Data can be obtained e.g. from https://refractiveindex.info/ using the “Full database record” export function.

  • unit_wl (str, default: 'micron') – Units of the wavelength in file, one of [‘micron’, ‘nm’]
  • interpolate_order (int, default: 1) – interpolation order for data (1: linear, 2: square, 3: cubic) “1” uses numpy.interp, “2” and “3” require scipy (using scipy.interpolate.interp1d)
__init__(refindex_file, unit_wl='micron', interpolate_order=1)

Use tabulated dispersion

epsilon(wavelength)

Tabulated interpolated dielectric function

constant dielectric function material

Parameters:wavelength (real) – wavelength at which to evaluate dielectric function (in nm)

Predefined materials

class pyGDM2.materials.dummy(n=(2+0j))

constant index

Material with spectrally constant refractive index

Parameters:n (complex, default: (2.0 + 0.0j)) – complex refractive index of constant material (returned dielectric function will be n**2)
class pyGDM2.materials.gold(interpolate_order=1)

gold index

Complex dielectric function of gold from: P. B. Johnson and R. W. Christy. Optical Constants of the Noble Metals, Phys. Rev. B 6, 4370-4379 (1972)

Parameters:interpolate_order (int, default: 1) – interpolation order for data (1: linear, 2: square, 3: cubic) “1” uses numpy, “2” and “3” require scipy (scipy.interpolate.interp1d)
class pyGDM2.materials.alu(interpolate_order=1)

alu index

Complex dielectric function of aluminium from: A. D. Rakić, A. B. Djurišic, J. M. Elazar, and M. L. Majewski. Optical properties of metallic films for vertical-cavity optoelectronic devices, Appl. Opt. 37, 5271-5283 (1998)

Parameters:interpolate_order (int, default: 1) – interpolation order for data (1: linear, 2: square, 3: cubic) “1” uses numpy, “2” and “3” require scipy (scipy.interpolate.interp1d)
class pyGDM2.materials.silicon(interpolate_order=1)

silicon index

Complex dielectric function of silicon from: Edwards, D. F. in Handbook of Optical Constants of Solids (ed. Palik, E. D.) 547–569 (Academic Press, 1997).

Parameters:interpolate_order (int, default: 1) – interpolation order for data (1: linear, 2: square, 3: cubic) “1” uses numpy, “2” and “3” require scipy (scipy.interpolate.interp1d)

fields

fields contains the main efield object and field-generator functions for frequently used fundamental fields.

Incident electro-magnetic field class efield

class pyGDM2.fields.efield(field_generator, wavelengths, kwargs)

incident electromagnetic field container class

Defines an incident electric field including information about wavelengths, polarizations, focal spot or whatever optional parameter is supported by the used field generator.

Parameters:
  • field_generator (callable) –
    field generator function. Mandatory arguments are
  • wavelengths (list) – list of wavelengths (wavelengths in nm)
  • kwargs (list of dict or dict) –

    possible additional keyword arguments, passed to field_generator. Either dict or list of dicts.

    • If list of dicts, each entry must correspond exactly to one parameters-set for field-generator.
    • If dict, maybe contain lists for configurations of the parameters. In that case, all possible parameter-permutations will be generated.

Examples

>>> kwargs = dict(theta = [0.0,45,90])
[{'theta': 0.0}, {'theta': 45.0}, {'theta': 90.0}]

is equivalent to:

>>> kwargs = [dict(theta=0.0), dict(theta=45.0), dict(theta=90.0)]
[{'theta': 0.0}, {'theta': 45.0}, {'theta': 90.0}]
__init__(field_generator, wavelengths, kwargs)

initialize field container

setDtype(dtypef, dtypec)

set dtype of arrays

field generators

plane wave

pyGDM2.fields.planewave(struct, wavelength, theta, kSign=-1, consider_substrate_reflection=False, returnField='E')

Normally incident (along Z) planewave

Amplitude = 1 for both, E and B

Parameters:
  • struct (structures.struct) – structure definition, including environment setup
  • wavelength (float) – Wavelength in nm
  • theta (float) – polarization angle in degrees, 0deg = 0X direction
  • kSign (int, default: -1) – sign of wavenumber. +1: propagation from bottom to top (towards increasing z) -1: propagation from top to bottom (towards smaller z, default)
  • consider_substrate_reflection (bool, default: False) – Whether to consider the reflection / transmission coefficient at the substrate for adjusting the field amplitude
  • returnField (str, default: 'E') – if ‘E’: returns electric field; if ‘B’: magnetic field
Returns:

E0 (B0) – list of (complex) 3-tuples: [(Ex1, Ey1, Ez1), …]

Return type:

Complex E-(B-)Field at each dipole position as

oblique plane wave for total internal reflection sim.

pyGDM2.fields.evanescent_planewave(struct, wavelength, theta_inc=0, polar='s', returnField='E')

oblique incident planewave

Oblique incidence (from bottom to top) through n1/n2/n3 layer interfaces. May be used to simulate evanescent fields in the total internal reflection configuration. Amplitude = 1 for both, E and B.

Parameters:
  • struct (structures.struct) – structure definition, including environment setup
  • wavelength (float) – Wavelength in nm
  • theta_inc (float, default: 0) –
    incident angle in the XZ plane with respect to e_z, in degrees.
    • 0deg = along Z (form neg to pos Z)
    • 90deg = along X (form pos to neg X)
  • polar (str, default: 's') – incident polarization. Either ‘s’ or ‘p’. At 0 degree incident angle, ‘s’ is polarized along x, ‘p’ along y.
  • returnField (str, default: 'E') – if ‘E’: returns electric field; if ‘B’: magnetic field
Returns:

E0 (B0) – list of (complex) 3-tuples: [(Ex1, Ey1, Ez1), …]

Return type:

Complex E-(B-)Field at each dipole position as

focused plane wave

pyGDM2.fields.focused_planewave(struct, wavelength, theta, xSpot, ySpot, NA=-1, spotsize=-1, kSign=-1, phase=0, consider_substrate_reflection=False, returnField='E')

Normally incident (along Z) planewave with gaussian intensity profile

Maximum amplitude = 1 for both, E and B, focused at (x0,y0)

Parameters:
  • struct (structures.struct) – structure definition, including environment setup
  • wavelength (float) – Wavelength in nm
  • theta (float) – polarization angle in degrees, 0deg = 0X direction
  • kSign (int, default: -1) – sign of wavenumber. +1: propagation from bottom to top (towards increasing z) -1: propagation from top to bottom (towards smaller z, default)
  • phase (float, default: 0) – additional phase of the beam, in degrees
  • consider_substrate_reflection (bool, default: False) – Whether to consider the reflection / transmission coefficient at the substrate for adjusting the field amplitude
  • returnField (str, default: 'E') – if ‘E’: returns electric field; if ‘B’: magnetic field
Returns:

E0 (B0) – list of (complex) 3-tuples: [(Ex1, Ey1, Ez1), …]

Return type:

Complex E-(B-)Field at each dipole position as

gaussian

pyGDM2.fields.gaussian(struct, wavelength, theta, xSpot, ySpot, zSpot=0, NA=-1, spotsize=-1, kSign=-1, paraxial=False, phase=0, returnField='E')

Normal incident (along Z) Gaussian Beam Field

obligatory “einKwargs” are: ‘xSpot’, ‘ySpot’ and one of ‘NA’ or ‘spotsize’

Parameters:
  • struct (structures.struct) – structure object including environment
  • wavelength (float) – Wavelength in nm
  • theta (float) – polarization angle in degrees, 0deg = 0X direction
  • xSpot,ySpot,zSpot (float,float) – x/y/z coordinates of focal point
  • NA (float) – Numerical aperture to calculate beamwaist
  • spotsize (float (optional)) – Gaussian beamwaist (overrides “NA”)
  • kSign (float, default: -1) – Direction of Beam. -1: top to Bottom, 1 Bottom to top
  • paraxial (bool, default: False) – Use paraxial Gaussian beam: No longitudinal fields. If “False”, longitudinal components are obtained using Maxwell equation div(E)=0 as condition
  • phase (float, default: 0) – additional phase of the beam, in degrees
  • returnField (str, default: 'E') – if ‘E’: returns electric field; if ‘B’: magnetic field
Returns:

E0

Return type:

Complex E-Field at each dipole position

double gaussian

pyGDM2.fields.double_gaussian(struct, wavelength, theta1, theta2, xSpot1, ySpot1, zSpot1, xSpot2, ySpot2, zSpot2, phase1=0, phase2=0, beam1_amplitude=1.0, beam2_amplitude=1.0, paraxial=True, spotsize=-1, kSign=-1, returnField='E')

Two focused beams, based on gaussian()

Parameters:
  • theta2 (theta1,) – polarization angle of beam 1 and 2
  • ySpot1, zSpot1 (xSpot1,) – beam 1 focal spot
  • ySpot2, zSpot2 (xSpot2,) – beam 2 focal spot
  • phase2 (phase1,) – relative phase, in degrees
  • beam2_amplitude (beam1_amplitude,) – max. amplitude of beam 1, beam 2

:param For the other parameters, see gaussian():

Returns:E0
Return type:Complex E-Field at each dipole position ( (Ex,Ey,Ez)-tuples )

electric dipolar emitter

pyGDM2.fields.dipole_electric(struct, wavelength, x0, y0, z0, mx, my, mz, returnField='E')

field emitted by an electric dipole at (x0,y0,z0) with complex amplitude (mx,my,mz)

obligatory kwargs are: wavelength, x0,y0,z0, mx,my,mz

From: G. S. Agarwal, Phys. Rev. A, 11, 230 (1975). (Eqs. 4.5/4.6)

Parameters:
  • struct (structures.struct) – structure object including environment
  • wavelength (float) – Wavelength in nm
  • x0,y0,z0 (float) – x/y/z coordinates of magnetic dipole position
  • mx,my,mz (float) – x/y/z amplitude of mag. dipole vector
  • returnField (str, default: 'E') – if ‘E’: returns electric field; if ‘B’: magnetic field
Returns:

E0

Return type:

Complex E-Field at each dipole position ( (Ex,Ey,Ez)-tuples )

magnetic dipolar emitter

pyGDM2.fields.dipole_magnetic(struct, wavelength, x0, y0, z0, mx, my, mz, returnField='E')

field emitted by a magnetic dipole at (x0,y0,z0) with complex amplitude (mx,my,mz)

obligatory “einKwargs” are: wavelength, x0,y0,z0, mx,my,mz

From: G. S. Agarwal, Phys. Rev. A, 11, 230 (1975). (Eqs. 4.5/4.6)

Parameters:
  • struct (structures.struct) – structure object including environment
  • wavelength (float) – Wavelength in nm
  • x0,y0,z0 (float) – x/y/z coordinates of magnetic dipole position
  • mx,my,mz (float) – x/y/z amplitude of mag. dipole vector
  • returnField (str, default: 'E') – if ‘E’: returns electric field; if ‘B’: magnetic field
Returns:

E0

Return type:

Complex E-Field at each dipole position ( (Ex,Ey,Ez)-tuples )