holopy.scattering.theory package

Theories to compute scattering from objects.

All theories have a common interface defined by holopy.scattering.theory.scatteringtheory.ScatteringTheory.

Subpackages

• holopy.scattering.theory.mie_f package
  • Submodules
  • holopy.scattering.theory.mie_f.mie_specfuncs module
  • holopy.scattering.theory.mie_f.miescatlib module
  • holopy.scattering.theory.mie_f.multilayer_sphere_lib module
• holopy.scattering.theory.tmatrix_f package

Submodules

holopy.scattering.theory.dda module

Compute holograms using the discrete dipole approximation (DDA). Currently uses ADDA (https://github.com/adda-team/adda) to do DDA calculations.

class DDA(n_cpu=1, use_gpu=False, gpu_id=None, max_dpl_size=None, use_indicators=True, keep_raw_calculations=False, addacmd=[], suppress_C_output=True)

Bases: ScatteringTheory

Computes scattering using the the Discrete Dipole Approximation (DDA). It can (in principle) calculate scattering from any arbitrary scatterer. The DDA uses a numerical method that represents arbitrary scatterers as an array of point dipoles and then self-consistently solves Maxwell’s equations to determine the scattered field. In practice, this model can be extremely computationally intensive, particularly if the size of the scatterer is larger than the wavelength of light. This model requires an external scattering code: a-dda

n_cpu

Number of threads to use for the DDA calculation

Type:

int (optional)

max_dpl_size

Force a maximum dipole size. This is useful for forcing extra dipoles if necessary to resolve features in an object. This may make dda calculations take much longer.

Type:

float (optional)

use_indicators

If true, a scatterer’s indicators method will be used instead of its built-in adda definition

Type:

bool

keep_raw_calculations

If true, do not delete the temporary file we run ADDA in, instead print its path so you can inspect its raw results

Type:

bool

Notes

Does not handle near fields. This introduces ~5% error at 10 microns. This can in principle handle any scatterer, but in practice it will need excessive memory or computation time for particularly large scatterers.

classmethod can_handle(scatterer)

Given a scatterer, returns a bool

raw_scat_matrs(scatterer, pos, medium_wavevec, medium_index)

Given a (3, N) array pos etc, returns an (N, 2, 2) array

required_spacing(bounds, medium_wavelen, medium_index, n)

holopy.scattering.theory.lens module

class Lens(lens_angle, theory, quad_npts_theta=100, quad_npts_phi=100, use_numexpr=True)

Bases: ScatteringTheory

Wraps a ScatteringTheory and overrides the raw_fields to include the effect of an objective lens.

can_handle(scatterer)

Given a scatterer, returns a bool

desired_coordinate_system = 'cylindrical'

numexpr_integrand_prefactor1 = 'exp(1j * krho_p * sintheta * cos(phi_relative))'

numexpr_integrand_prefactor2 = 'exp(1j * kz_p * (1 - costheta))'

numexpr_integrand_prefactor3 = 'sqrt(costheta) * sintheta * phi_wts * theta_wts'

numexpr_integrandl = 'prefactor * (cosphi * (cosphi * S2 + sinphi * S3) + sinphi * (cosphi * S4 + sinphi * S1))'

numexpr_integrandr = 'prefactor * (sinphi * (cosphi * S2 + sinphi * S3) - cosphi * (cosphi * S4 + sinphi * S1))'

parameter_names = ('lens_angle',)

raw_fields(positions, scatterer, medium_wavevec, medium_index, illum_polarization)

Given a (3, N) array pos, etc, returns a (3, N) array

gauss_legendre_pts_wts(a, b, npts=100)

Quadrature points for integration on interval [a, b]

pts_wts_for_phi_integrals(npts)

Quadrature points for integration on the periodic interval [0, pi]

Since this interval is periodic, we use equally-spaced points with equal weights.

holopy.scattering.theory.mie module

Calculates holograms of spheres using Fortran implementation of Mie theory. Uses superposition to calculate scattering from multiple spheres. Uses full radial dependence of spherical Hankel functions for scattered field.

class Mie(compute_escat_radial=True, full_radial_dependence=True, eps1=0.01, eps2=1e-16)

Bases: ScatteringTheory

Compute scattering using the Lorenz-Mie solution.

This theory calculates exact scattering for single spheres and approximate results for groups of spheres. It does not account for multiple scattering, hence the approximation in the case of multiple spheres. Neglecting multiple scattering is a good approximation if the particles are sufficiently separated.

This model can also calculate the exact scattered field from a spherically symmetric particle with an arbitrary number of layers with differing refractive indices, using Yang’s recursive algorithm ([Yang2003]).

By default, calculates radial component of scattered electric fields, which is nonradiative.

Currently, in calculating the Lorenz-Mie scattering coefficients, the maximum size parameter x = ka is limited to 1000.

can_handle(scatterer)

Given a scatterer, returns a bool

raw_cross_sections(scatterer, medium_wavevec, medium_index, illum_polarization)

Calculate scattering, absorption, and extinction cross sections, and asymmetry parameter for spherically symmetric scatterers.

Parameters:

scatterer (scatterpy.scatterer object) – spherically symmetric scatterer to compute for (Calculation would need to be implemented in a radically different way, via numerical quadrature, for sphere clusters)

Returns:

cross_sections – Dimensional scattering, absorption, and extinction cross sections, and <cos heta>

Return type:

array (4)

Notes

The radiation pressure cross section C_pr is given by C_pr = C_ext - <cos heta> C_sca.

The radiation pressure force on a sphere is

F = (n_med I_0 C_pr) / c

where I_0 is the incident intensity. See van de Hulst, p. 14.

raw_fields(positions, scatterer, medium_wavevec, medium_index, illum_polarization)

Given a (3, N) array pos, etc, returns a (3, N) array

raw_scat_matrs(scatterer, pos, medium_wavevec, medium_index)

Returns far-field amplitude scattering matrices (with theta and phi dependence only) – assume spherical wave asymptotic r dependence

holopy.scattering.theory.mielens module

class AberratedMieLens(spherical_aberration=0.0, lens_angle=1.0, calculator_accuracy_kwargs={})

Bases: MieLens

parameter_names = ('lens_angle', 'spherical_aberration')

class MieLens(lens_angle=1.0, calculator_accuracy_kwargs={})

Bases: ScatteringTheory

Exact scattering from a sphere imaged through a perfect lens.

Calculates holograms of spheres using an analytical solution of the Mie scattered field imaged by a perfect lens (see [Leahy2020]). Can use superposition to calculate scattering from multiple spheres.

See also

mielensfunctions.MieLensCalculator

can_handle(scatterer)

Given a scatterer, returns a bool

desired_coordinate_system = 'cylindrical'

parameter_names = ('lens_angle',)

raw_fields(positions, scatterer, medium_wavevec, medium_index, illum_polarization)

Parameters:

• positions ((3, N) numpy.ndarray) – The (k * rho, phi, z) coordinates, relative to the sphere, of the points to calculate the fields. Note that the radial coordinate is rescaled by the wavevector.

• scatterer (scatterer.Sphere object)

• medium_wavevec (float)

• medium_index (float)

• illum_polarization (2-element tuple) – The (x, y) field polarizations.

holopy.scattering.theory.mielensfunctions module

class AberratedMieLensCalculator(spherical_aberration=None, **kwargs)

Bases: MieLensCalculator

This class handles lenses with spherical aberration. It is based on MieLensCalculator, which documents the methods and parameters.

Parameters:

• spherical_aberration (float or array-like of floats) – The spherical aberration, up to arbitrary order. If a float, just the coefficient of the 3rd-order aberration (4th-order in wavefront). When an array, the coefficients of aberrations in ascending order (3rd, 5th, 7th, etc), where the wavefront distortion for the nth-order aberration is of the form (cos(theta) - 1)^(n+1), where n = 3, 5, 7, etc Default is None, which raises an error.

• MieLensCalculator (See)

must_be_specified = ['particle_kz', 'index_ratio', 'size_parameter', 'lens_angle', 'spherical_aberration']

class AlBlFunctions

Bases: object

Group of functions for calculating the Mie scattering coefficients, used for expressing the scattered field in terms of vector spherical harmonics.

The coefficients a_l, b_l are defined as

\[ \begin{align}\begin{aligned}\begin{split}a_l = \frac{\psi_l(x) \psi_l'(nx) - n \psi_l(nx) \psi_l'(x)} {\\xi_l(x) \psi_l'(nx) - n \psi_l(nx) \\xi_l'(x)},\end{split}\\\begin{split}b_l = \frac{\psi_l(nx) \psi_l'(x) - n \psi_l(x) \psi_l'(nx)} {\psi_l(nx) \\xi_l'(x) - n \\xi_l(x) \psi_l'(nx)},\end{split}\end{aligned}\end{align} \]

where \(\psi_l\) and \(\\xi_l\) are the Riccati-Bessel functions of the first and third kinds, respectively. The definitions used here follow those of van de Hulst [vdh], which differ from those used in Bohren and Huffman [bandh].

References

[vdh]

H. C. van de Hulst, “Light Scattering by Small Particles”, Dover (1981), pg 123.

[bandh]

C. F. Bohren and Donald R. Huffman, “Absorption and Scattering of Light by Small Particles”, Wiley (2004), pg 101.

static calculate_al_bl(index_ratio, size_parameter, l)

Returns a_l and b_l; see class docstring.

Parameters:

• index_ratio (float) – relative index of refraction

• size_paramter (float) – Size parameter

• l (int, array-like) – Order of scattering coefficient

Returns:

a_l, b_l

Return type:

numpy.ndarray

static riccati_psin(n, z, derivative=False)

Riccati-Bessel function of the first kind or its derivative.

\[\psi_n(z) = z\,j_n(z),\]

where \(j_n(z)\) is the spherical Bessel function of the first kind.

Parameters:

• n (int, array_like) – Order of the Bessel function (n >= 0).

• z (complex or float, array_like) – Argument of the Bessel function.

• derivative (bool, optional) – If True, the value of the derivative (rather than the function itself) is returned.

Returns:

psin

Return type:

ndarray

static riccati_xin(order, z, derivative=False)

Riccati-Bessel function of the third kind or its derivative.

\[\begin{split}\\xi_n(z) = z\,h^{(1)}_n(z),\end{split}\]

where \(h^{(1)}_n(z)\) is the first spherical Hankel function.

Parameters:

• n (int, array_like) – Order of the Bessel function (n >= 0).

• z (complex or float, array_like) – Argument of the Bessel function.

• derivative (bool, optional) – If True, the value of the derivative (rather than the function itself) is returned.

Returns:

xin

Return type:

ndarray

class MieLensCalculator(particle_kz=None, index_ratio=None, size_parameter=None, lens_angle=None, quad_npts=100, interpolate_integrals='check', interpolator_window_size=30.0, interpolator_degree=32)

Bases: object

Calculates the field from a Mie scatterer imaged in a high-NA lens.

The incindent electric field is E e^{ikz}, with the particle position at z. The scattered field takes into account the varying phase of the incoming field.

Parameters:

• particle_kz (float) –

  • z is away from the lens

• index_ratio (float > 0)

• size_parameter (float > 0)

• lens_angle (float on (0, pi/2))

calculate_scattered_field(krho, phi)

tuple of 2 numpy.ndarrays, of shape krho

calculate_total_field(krho, phi)

tuple of 2 numpy.ndarrays, of shape krho

calculate_total_intensity(krho, phi)

numpy.ndarray, of shape krho

Parameters:

• quad_npts (int, optional) – The number of points for numerical quadrature of the integrals over the lens pupil.

• interpolate_integrals ({'check', True, False}) – Whether or not to interpolate the internally-evaluated integrals for speed. Default is ‘check’, which interpolates if it will be faster or does direct numerical quadrature otherwise. Interpolation is done via a piecewise Chebyshev approximant.

• interpolator_window_size (float, optional) – The spacing, in units of 1/k, for the windows of the piecewise Chebyshev approximants. A lower value gives more accurate results, although accuracy depends on the interpolator_degree parameter as well. The default is 39. which gives 5e-13 relative accuracy.

• interpolator_degree (int, optional) –

  The polynomial degree for the piecewise Chebyshev approximants. A higher value gives more accurate results, although accuracy depends on the interpolator_window_size parameter as well. The default is 32, which gives 5e-13 relative accuracy.

  It is best to leave the interpolator parameters as-is; they are only exposed for testing and advanced usage.

calculate_scattered_field(krho, phi)

Calculates the field from a Mie scatterer imaged through a high-NA lens and excited with an electric field of unit strength directed along the optical axis.

\[\vec{E}_{sc} = A \left[ I_{12} \sin(2\phi) \hat{y} + -I_{10} \hat{x} + I_{12} \cos(2\phi) \hat{x} + -I_{20} \hat{x} + -I_{22} \cos(2\phi) \hat{x} + -I_{22} \sin(2\phi) \hat{y} \right]\]

Parameters:

• krho (numpy.ndarray) – The position of the particle relative to the focal point of the lens, in (i) cylindrical coordinates and (ii) dimensionless wavevectur units. Must all be the same shape.

• phi (numpy.ndarray) – The position of the particle relative to the focal point of the lens, in (i) cylindrical coordinates and (ii) dimensionless wavevectur units. Must all be the same shape.

Returns:

field_xcomp, field_ycomp – The (x, y) components of the electric field at the detector, where the initial field is polarized in the x-direction. Same shape as krho, phi

Return type:

numpy.ndarray

Notes

This will have problems for large rho, z, because of the quadrature points. Empirically this problem happens for rho >~ 4 * quad_npts. Could be adaptive if needed….

calculate_total_field(krho, phi)

The total (incident + scattered) field at the detector

calculate_total_intensity(krho, phi)

must_be_specified = ['particle_kz', 'index_ratio', 'size_parameter', 'lens_angle']

class MieScatteringMatrix(parallel_or_perpendicular='perpendicular', index_ratio=None, size_parameter=None, max_l=None)

Bases: object

class PiecewiseChebyshevApproximant(function, degree, window_breakpoints, *args)

Bases: object

calculate_al_bl(index_ratio, size_parameter, l)

calculate_pil_taul(theta, max_order)

The 1st through Nth order angle dependent functions for Mie scattering, evaluated at theta. The functions \(\pi(\theta)\) and \(\tau(\theta)\) are defined as:

\[\begin{split}\pi_n(\theta) &= \frac{1}{\sin \theta} P_n^1(\cos\theta) \\ \tau_n(\theta) &= \frac{\mathrm{d}}{\mathrm{d}\theta} P_n^1(\cos\theta)\end{split}\]

where \(P_n^m\) is the associated Legendre function. The functions are computed by upward recurrence using the relations

\[\begin{split}\pi_n &= \frac{2n-1}{n-1}\cos\theta \, \pi_{n-1} - \frac{n}{n-1}\pi_{n-2} \\ \tau_n &= n \, \cos\theta \, \pi_n - (n+1)\pi_{n-1}\end{split}\]

beginning with \(\pi_0 = 0\) and \(\pi_1 = 1\)

Parameters:

• theta (array_like) – angles (in radians) at which to evaluate the angular functions

• max_order (int > 0) – Order at which to halt iteration. Must be > 0

Returns:

pi, tau – 2D arrays with shape (len(theta), max_order) containing the values of the angular functions evaluated at theta up to order max_order

Return type:

ndarray

gauss_legendre_pts_wts(a, b, npts=100)

Quadrature points for integration on interval [a, b]

j2(x)

A fast J_2(x) defined in terms of other special functions

spherical_h1n(n, z, derivative=False)

Spherical Hankel function H_n(z) or its derivative

spherical_h2n(n, z, derivative=False)

Spherical Hankel function H_n(z) or its derivative

holopy.scattering.theory.multisphere module

Defines Multisphere theory class, which calculates scattering for multiple spheres using the (exact) superposition method implemented in modified version of Daniel Mackowski’s SCSMFO1B.FOR. Uses full radial dependence of spherical Hankel functions for the scattered field.

class Multisphere(niter=200, eps=1e-06, meth=1, qeps1=1e-05, qeps2=1e-08, compute_escat_radial=False, suppress_fortran_output=True)

Bases: ScatteringTheory

Exact scattering from a cluster of spheres.

Calculate the scattered field of a collection of spheres through a numerical method that accounts for multiple scattering and near-field effects (see [Fung2011], [Mackowski1996]). This approach is much more accurate than Mie superposition, but it is also more computationally intensive. The Multisphere code can handle any number of spheres; see notes below for details.

niter

maximum number of iterations to use in solving the interaction equations

Type:

integer (optional)

meth

method to use to solve interaction equations. Set to 0 for biconjugate gradient; 1 for order-of-scattering

Type:

integer (optional)

eps

relative error tolerance in solution for interaction equations

Type:

float (optional)

qeps1

error tolerance used to determine at what order the single-sphere spherical harmonic expansion should be truncated

Type:

float (optional)

qeps2

error tolerance used to determine at what order the cluster spherical harmonic expansion should be truncated

Type:

float (optional)

Notes

According to Mackowski’s manual for SCSMFO1B.FOR [1] and later papers [2], the biconjugate gradient is generally the most efficient method for solving the interaction equations, especially for dense arrays of identical spheres. Order-of-scattering may converge better for non-identical spheres.

Multisphere does not check for overlaps becaue overlapping spheres can be useful for getting fits to converge. The results to be sensible for small overlaps even though mathemtically speaking they are not xstrictly valid.

Currently, Multisphere does not calculate the radial component of scattered electric fields. This is a good approximation for large kr, since the radial component falls off as 1/kr^2.

scfodim.for contains three parameters, all integers:

• nod: Maximum number of spheres

• nod: Maximum order of individual sphere expansions. Will depend on

  size of largest sphere in cluster.

• notd: Maximum order of cluster-centered expansion. Will depend on

  overall size of cluster.

Changing these values will require recompiling Fortran extensions.

The maximum size parameter of each individual sphere in a cluster is currently limited to 1000, indepdently of the above scfodim.for parameters.

References

[1]

Daniel W. Mackowski, SCSMFO.FOR: Calculation of the Scattering Properties for a Cluster of Spheres, ftp://ftp.eng.auburn.edu/pub/dmckwski/scatcodes/scsmfo.ps

[2]

D.W. Mackowski, M.I. Mishchenko, A multiple sphere T-matrix Fortran code for use on parallel computer clusters, Journal of Quantitative Spectroscopy and Radiative Transfer, In Press, Corrected Proof, Available online 11 March 2011, ISSN 0022-4073, DOI: 10.1016/j.jqsrt.2011.02.019.

can_handle(scatterer)

Given a scatterer, returns a bool

raw_cross_sections(scatterer, medium_wavevec, medium_index, illum_polarization)

Calculate scattering, absorption, and extinction cross sections, and asymmetry parameter for sphere clusters with polarized incident light.

The extinction cross section is calculated from the optical theorem. The scattering cross section is calculated by numerical quadrature of the scattered field, and the absorption cross section is calculated from the difference of the extinction cross section and the scattering cross section.

Parameters:

scatterer (scatterpy.scatterer object) – sphere cluster to compute for

Returns:

cross_sections – Dimensional scattering, absorption, and extinction cross sections, and <cos heta>

Return type:

array (4)

raw_fields(positions, scatterer, medium_wavevec, medium_index, illum_polarization)

Given a (3, N) array pos, etc, returns a (3, N) array

raw_scat_matrs(scatterer, pos, medium_wavevec, medium_index)

Calculate far-field amplitude scattering matrices at multiple positions

normalize_polarization(illum_polarization)

holopy.scattering.theory.scatteringtheory module

class ScatteringTheory

Bases: HoloPyObject

Defines common interface for all scattering theories.

Subclasses must implement: * can_handle * raw_fields or raw_scat_matrs or both. * (optionally) raw_cross_sections,

Notes

Subclasses should implement the following methods to create a scatteringtheory that works for the following user-facing methods: * calc_holo: raw_fields or raw_scat_matrs * calc_intensity: raw_fields or raw_scat_matrs * calc_field: raw_fields or raw_scat_matrs * calc_scat_matrix: raw_scat_matrs * calc_cross_sections: raw_cross_sections

By default, ScatteringTheories computer raw_fields from the raw_scat_matrs; over-ride the raw_fields method to compute the fields in a different way.

can_handle(scatterer)

Given a scatterer, returns a bool

desired_coordinate_system = 'spherical'

from_parameters(parameters)

Creates a ScatteringTheory like the current one, but with different parameters. Used for fitting

Parameters:

dict – keys should be valid self.parameter_names fields, values should be the corresponding kwargs

Return type:

ScatteringTheory instance, of the same class as self

parameter_names = ()

property parameters

raw_cross_sections(scatterer, medium_wavevec, medium_index, illum_polarization)

Returns cross-sections, as an array [cscat, cabs, cext, asym]

raw_fields(pos, scatterer, medium_wavevec, medium_index, illum_polarization)

Given a (3, N) array pos, etc, returns a (3, N) array

raw_scat_matrs(scatterer, pos, medium_wavevec, medium_index)

Given a (3, N) array pos etc, returns an (N, 2, 2) array

holopy.scattering.theory.tmatrix module

Compute holograms using Mishchenko’s T-matrix method for axisymmetric scatterers. Currently uses

class Tmatrix

Bases: ScatteringTheory

Computes scattering using the axisymmetric T-matrix solution by Mishchenko with extended precision.

It can calculate scattering from axisymmetric scatterers such as cylinders and spheroids. Calculations for particles that are very large or have high aspect ratios may not converge.

Notes

Does not handle near fields. This introduces ~5% error at 10 microns.

can_handle(scatterer)

Given a scatterer, returns a bool

raw_fields(pos, scatterer, medium_wavevec, medium_index, illum_polarization)

Given a (3, N) array pos, etc, returns a (3, N) array

raw_scat_matrs(scatterer, pos, medium_wavevec, medium_index)

Given a (3, N) array pos etc, returns an (N, 2, 2) array