holopy.scattering.theory package¶
Theories to compute scattering from objects.
All theories have a common interface defined by
holopy.scattering.theory.scatteringtheory.ScatteringTheory
.
Subpackages¶
Submodules¶
Compute holograms using the discrete dipole approximation (DDA). Currently uses ADDA (https://github.com/adda-team/adda) to do DDA calculations. .. moduleauthor:: Thomas G. Dimiduk <tdimiduk@physics.harvard.edu>
-
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:
holopy.scattering.theory.scatteringtheory.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
-
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)¶
-
-
class
Lens
(lens_angle, theory, quad_npts_theta=100, quad_npts_phi=100, use_numexpr=True)¶ Bases:
holopy.scattering.theory.scatteringtheory.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.
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:
holopy.scattering.theory.scatteringtheory.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
-
-
class
AberratedMieLens
(spherical_aberration=0.0, lens_angle=1.0, calculator_accuracy_kwargs={})¶ Bases:
holopy.scattering.theory.mielens.MieLens
-
parameter_names
= ('lens_angle', 'spherical_aberration')¶
-
-
class
MieLens
(lens_angle=1.0, calculator_accuracy_kwargs={})¶ Bases:
holopy.scattering.theory.scatteringtheory.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.
-
-
class
AberratedMieLensCalculator
(spherical_aberration=None, **kwargs)¶ Bases:
holopy.scattering.theory.mielensfunctions.MieLensCalculator
-
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
..math:
a_l =
- rac{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)},
b_l =
- rac{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)},
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 [1]_, which differ from those used in Bohren and Huffman [2]_.
[1] H. C. van de Hulst, “Light Scattering by Small Particles”, Dover (1981), pg 123. [2] 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: 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.
\[\xi_n(z) = z\,h^{(1)}_n(z),\]where \(h^{(1)}_n(z)\) is the first spherical Hankel function.
Parameters: 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
-
calculate_incident_field
()¶ This is here so (i) Any corrections in the theory to the scattered field
have an easy place to enter, and- Other modules can consistently use the same scattered field as this module.
-
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.
\[\]- ec{E}_{sc} = A left[ I_{12} sin(2phi) hat{y} +
- -I_{10} hat{x} +
- I_{12} cos(2phi) hat{x} +
-I_{20} hat{x} + -I_{22} cos(2phi) hat{x} + -I_{22} sin(2phi) hat{y}
ight]
- krho, 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.
- field_xcomp, field_ycomp : numpy.ndarray
- 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
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
-
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 :math`pi( heta)` and :math` au( heta) are defined as:
..math:
\pi_n( heta) =
rac{1}{sin heta} P_n^1(cos heta)
au_n( heta) =rac{mathrm{d}}{mathrm{d} heta} P_n^1(cos heta)
where \(P_n^m\) is the associated Legendre function. The functions are computed by upward recurrence using the relations
..math:
\pi_n =
rac{2n-1}{n-1}cos heta , pi_{n-1} - rac{n}{n-1}pi_{n-2}
au_n = n , cos heta , pi_n - (n+1)pi_{n-1}beginning with \(pi_0 = 0\) and \(pi_1 = 1\)
- 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
- pi, tau : ndarray
- 2D arrays with shape (len(theta), max_order) containing the values of the angular functions evaluated at theta up to order max_order
-
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
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:
holopy.scattering.theory.scatteringtheory.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)
-
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 forReturns: 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)¶
-
class
ScatteringTheory
¶ Bases:
holopy.core.holopy_object.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 Returns: Return type: ScatteringTheory instance, of the same class as self
-
parameter_names
= ()¶
-
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
-
Compute holograms using Mishchenko’s T-matrix method for axisymmetric scatterers. Currently uses
-
class
Tmatrix
¶ Bases:
holopy.scattering.theory.scatteringtheory.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
-