Basis sets¶
- class mumott.methods.basis_sets.SphericalHarmonics(probed_coordinates=None, ell_max=0, enforce_friedel_symmetry=True, **kwargs)[source]¶
Basis set class for spherical harmonics, the canonical representation of polynomials on the unit sphere and a simple way of representing band-limited spherical functions which allows for easy computations of statistics and is suitable for analyzing certain symmetries.
- Parameters
ell_max (int) – The bandlimit of the spherical functions that you want to be able to represent. A good rule of thumb is that
ell_max
should not exceed the number of detector segments minus 1.probed_coordinates (ProbedCoordinates) – Optional. A container with the coordinates on the sphere probed at each detector segment by the experimental method. Its construction from the system geometry is method-dependent. By default, an empty instance of
ProbedCoordinates
is created.enforce_friedel_symmetry (bool) – If set to
True
, Friedel symmetry will be enforced, using the assumption that points on opposite sides of the sphere are equivalent. This results in only evenell
being used.kwargs –
Miscellaneous arguments which relate to segment integrations can be passed as keyword arguments:
- integration_mode
Mode to integrate line segments on the reciprocal space sphere. Possible options are
'simpson'
,'midpoint'
,'romberg'
,'trapezoid'
.'simpson'
,'trapezoid'
, and'romberg'
use adaptive integration with the respective quadrature rule fromscipy.integrate
.'midpoint'
uses a single mid-point approximation of the integral. Default value is'simpson'
.- n_integration_starting_points
Number of points used in the first iteration of the adaptive integration. The number increases by the rule
N
←2 * N - 1
for each iteration. Default value is 3.- integration_tolerance
Tolerance for the maximum relative error between iterations before the integral is considered converged. Default is
1e-5
.- integration_maxiter
Maximum number of iterations. Default is
10
.
- property csr_representation: tuple¶
The projection matrix as a stack of sparse matrices in CSR representation as a tuple. The information in the tuple consists of the 3 dense matrices making up the representation, in the order
(pointers, indices, data)
.
- property ell_indices: ndarray[Any, dtype[int]]¶
The
ell
associated with each coefficient and its corresponding spherical harmonic. Updated whenell_max
changes.Notes
The word
ell
is used to represent the cursive small L, also written \(\ell\), often used as an index for the degree of the Legendre polynomial in the definition of the spherical harmonics.
- property ell_max: int¶
The maximum
ell
used to represent spherical functions.Notes
The word
ell
is used to represent the cursive small L, also written \(\ell\), often used as an index for the degree of the Legendre polynomial in the definition of the spherical harmonics.
- property emm_indices: ndarray[Any, dtype[int]]¶
The
emm
associated with each coefficient and its corresponding spherical harmonic. Updated whenell_max
changes.Notes
For consistency with
ell_indices
, and to avoid visual confusion with other letters,emm
is used to represent the index commonly written \(m\) in mathematical notation, the frequency of the sine-cosine parts of the spherical harmonics, often called the spherical harmonic order.
- forward(coefficients, indices=None)[source]¶
Carries out a forward computation of projections from spherical harmonic space into detector space, for one or several tomographic projections.
- Parameters
coefficients (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array of coefficients, of arbitrary shape so long as the last axis has the same size asell_indices
, and ifindices
is None or greater than one, the first axis should have the same length asindices
indices (
Optional
[ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]]) – Optional. Indices of the tomographic projections for which the forward computation is to be performed. IfNone
, the forward computation will be performed for all projections.
- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]- Returns
An array of values on the detector corresponding to the
coefficients
given. Ifindices
contains exactly one index, the shape is(coefficients.shape[:-1], J)
whereJ
is the number of detector segments. Ifindices
isNone
or contains several indices, the shape is(N, coefficients.shape[1:-1], J)
whereN
is the number of tomographic projections for which the computation is performed.
Notes
The assumption is made in this implementation that computations over several indices act on sets of images from different projections. For special usage where multiple projections of entire fields is desired, it may be better to use
projection_matrix
directly. This also applies togradient()
.
- generate_map(coefficients, resolution_in_degrees=5, map_half_sphere=True)¶
Generate a (theta, phi) map of the function modeled by the input coefficients. If
map_half_sphere=True
(default) a map of only the z>0 half sphere is returned.- Parameters
coefficients (
ndarray
[Any
,dtype
[float
]]) – One dimensional numpy array with lengthlen(self)
containing the coefficients of the function to be plotted.resolution_in_degrees (
int
) – The resoution of the map in degrees. The map uses eqidistant lines in longitude and latitude.map_half_sphere (
bool
) – If True returns a map of the z>0 half sphere.
- Return type
tuple
[ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]]- Returns
map_intensity – Intensity values of the map.
map_theta – Polar cooridnates of the map.
map_phi – Azimuthal coordinates of the map.
- get_covariances(u, v, resolve_spectrum=False)[source]¶
Returns the covariances of the spherical functions represented by two coefficient arrays.
- Parameters
u (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – The first coefficient array, of arbitrary shape except its last dimension must be the same length as the length ofell_indices
.v (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – The second coefficient array, of the same shape asu
.resolve_spectrum (
bool
) – Optional. Whether to resolve the product according to each frequency band, given by the coefficients of eachell
inell_indices
. Default value isFalse
.
- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]- Returns
An array of the covariances of the spherical functions represented by
u
andv
. Has the shape(u.shape[:-1])
if resolve_spectrum isFalse
, and(u.shape[:-1] + (ell_max // 2 + 1,))
if resolve_spectrum isTrue
, whereell_max
isSphericalHarmonics.ell_max
.
Notes
Calling this function is equivalent to calling
get_inner_product()
withspectral_moments=np.unique(ell_indices[ell_indices > 0])
whereell_indices
isSphericalHarmonics.ell_indices
. See the note toget_inner_product()
for mathematical details.
- get_inner_product(u, v, resolve_spectrum=False, spectral_moments=None)[source]¶
Retrieves the inner product of two coefficient arrays.
Notes
The canonical inner product in a spherical harmonic representation is \(\sum_\ell N(\ell) \sum_m u_m^\ell v_m^\ell\), where \(N(\ell)\) is a normalization constant (which is unity for the \(4\pi\) normalization). This inner product is a rotational invariant. The rotational invariance also holds for any partial sums over \(\ell\). One can define a function of \(\ell\) that returns such products, namely \(S(\ell, u, v) = N(\ell)\sum_m u_m^\ell v_m^\ell\), called the spectral power function. The sum \(\sum_{\ell = 1}S(\ell)\) is equal to the covariance of the band-limited spherical functions represented by \(u\) and \(v\), and each \(S(\ell, u, v)\) is the contribution to the covariance of the band \(\ell\). See also the SHTOOLS documentation for an excellent overview of this.
- Parameters
u (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – The first coefficient array, of arbitrary shape and dimension, except the last dimension must be the same as the length ofell_indices
.v (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – The second coefficient array, of the same shape asu
.resolve_spectrum (
bool
) – Optional. Whether to resolve the product according to each frequency band, given by the coefficients of eachell
inell_indices
. Defaults toFalse
, which means that the sum of every component of the spectrum is returned. IfTrue
, components are returned in order of ascendingell
. Theell
included in the spectrum depends onspectral_moments
.spectral_moments (
Optional
[List
[int
]]) – Optional. List of particular values ofell
to calculate the inner product for. Defaults toNone
, which is identical to including all values ofell
in the calculation. Ifspectral_moments
contains all nonzero values ofell
andresolve_spectrum
isFalse
, the covariance ofv
andu
will be calculated (the sum of the inner product over all non-zeroell
Ifresolve_spectrum
isTrue
, the covariance per ell inspectral_moments
, will be calculated, i.e., the inner products will not be summed over.
- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]- Returns
An array of the inner products of the spherical functions represented by
u
andv
. Has the shape(u.shape[:-1])
ifresolve_spectrum
isFalse
,(u.shape[:-1] + (ell_max // 2 + 1,))
ifresolve_spectrum
isTrue
andspectral_moments
isNone
, and finally the shape(u.shape[:-1] + (np.unique(spectral_moments).size,))
ifresolve_spectrum
isTrue
andspectral_moments
is a list of integers found inell_indices
- get_output(coefficients)[source]¶
Returns a dictionary of output data for a given array of spherical harmonic coefficients.
- Parameters
coefficients (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array of coefficients of arbitrary shape and dimensions, except its last dimension must be the same length asell_indices
. Computations only operate over the last axis ofcoefficients
, so derived properties in the output will have the shape(*coefficients.shape[:-1], ...)
.- Return type
- Returns
A dictionary containing two sub-dictionaries,
basis_set
andspherical_functions
.basis_set
contains information particular toSphericalHarmonics
, whereasspherical_functions
contains information about the spherical functions represented by thecoefficients
which are not specific to the chosen representation.
Notes
In detail, the two sub-dictionaries
basis_set
andspherical_functions
have the following members:- basis_set
- name
The name of the basis set, i.e.,
'SphericalHarmonicParameters'
- coefficients
A copy of
coefficients
.- ell_max
A copy of
ell_max
.- ell_indices
A copy of
ell_indices
.- emm_indices
A copy of
emm_indices
.- projection_matrix
A copy of
projection_matrix
.
- spherical_functions
- means
The spherical means of each function represented by
coefficients
.- variances
The spherical variances of each function represented by
coefficients
. Ifell_max
is0
, all variances will equal zero.- r2_tensors
The traceless symmetric rank-2 tensor component of each function represented by
coefficients
, in 6-element form, in the order[xx, yy, zz, yz, xz, xy]
, i.e., by the Voigt convention. The matrix form can be recovered as r2_tensors[…, tensor_to_matrix_indices], yielding matrix elements[[xx, xy, xz], [xy, yy, yz], [xz, yz, zz]]
. Ifell_max
is0
, all tensors have elements [1, 0, -1, 0, 0, 0].- tensor_to_matrix_indices
A list of indices to help recover the matrix from the 2-element form of the rank-2 tensors, equalling precisely
[[0, 5, 4], [5, 1, 3], [4, 3, 2]]
- eigenvalues
The eigenvalues of the rank-2 tensors, sorted in ascending order in the last index. If
ell_max
is0
, the eigenvalues will always be (1, 0, -1)- eigenvectors
The eigenvectors of the rank-2 tensors, sorted with their corresponding eigenvectors in the last index. Thus,
eigenvectors[..., 0]
gives the eigenvector corresponding to the smallest eigenvalue, andeigenvectors[..., 2]
gives the eigenvector corresponding to the largest eigenvalue. Generally, one of these two eigenvectors is used to define the orientation of a function, depending on whether it is characterized by a minimum (0
) or a maximum (2
). The middle eigenvector (1
) is typically only used for visualizations. Ifell_max
is0
, the eigenvectors will be the Cartesian basis vectors.- main_orientations
The estimated main orientations from the largest absolute eigenvalues. If
ell_max
is0
, the main orientation will be the x-axis.- main_orientation_symmetries
The strength or definiteness of the main orientation, calculated from the quotient of the absolute middle and signed largest eigenvalues of the rank-2 tensor. If
0
, the orientation is totally ambiguous. The orientation is completely transversal if the value is-1
(orientation represents a minimum), and completely longitudinal if the value is1
(orientation represents a maximum). Ifell_max
is0
, the main orientations are all totally ambiguous.- normalized_standard_deviations
A relative measure of the overall anisotropy of the spherical functions. Equals \(\sqrt{\sigma^2 / \mu}\), where \(\sigma^2\) is the variance and \(\mu\) is the mean. The places where \(\mu=0\) have been set to
0
. Ifell_max
is0
, the normalized standard deviations will always be zero.- power_spectra
The spectral powers of each
ell
inell_indices
, for each spherical function, sorted in ascendingell
. Ifell_max
is0
, each function will have only one element, equal to the mean squared.- power_spectra_ell
An array containing the corresponding
ell
to each of the last indices inpower_spectra
. Equal tonp.unique(ell_indices)
.
- get_spherical_harmonic_coefficients(coefficients, ell_max=None)[source]¶
Convert a set of spherical harmonics coefficients to a different
ell_max
by either zero-padding or truncation and return the result.
- gradient(coefficients, indices=None)[source]¶
Carries out a gradient computation of projections from spherical harmonic space into detector space, for one or several tomographic projections.
- Parameters
coefficients (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array of coefficients (or residuals) of arbitrary shape so long as the last axis has the same size as the number of detector segments.indices (
Optional
[ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]]) – Optional. Indices of the tomographic projections for which the gradient computation is to be performed. IfNone
, the gradient computation will be performed for all projections.
- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]- Returns
An array of gradient values based on the coefficients given. If
indices
contains exactly one index, the shape is(coefficients.shape[:-1], J)
whereJ
is the number of detector segments. If indices isNone
or contains several indices, the shape is(N, coefficients.shape[1:-1], J)
whereN
is the number of tomographic projections for which the computation is performed.
Notes
When solving an inverse problem, one should not to attempt to optimize the coefficients directly using the
gradient
one obtains by applying this method to the data. Instead, one must either take the gradient of the residual between theforward()
computation of the coefficients and the data. Alternatively one can apply both the forward and the gradient computation to the coefficients to be optimized, and the gradient computation to the data, and treat the residual of the two as the gradient of the optimization coefficients. The approaches are algebraically equivalent, but one may be more efficient than the other in some circumstances.
- property integration_mode: str¶
Mode of integration for calculating projection matrix. Accepted values are
'simpson'
,'romberg'
,'trapezoid'
, and'midpoint'
.
- property projection_matrix: ndarray[Any, dtype[_ScalarType_co]]¶
The matrix used to project spherical functions from the unit sphere onto the detector. If
v
is a vector of spherical harmonic coefficients, andM
is theprojection_matrix
, thenM @ v
gives the corresponding values on the detector segments associated with each projection.M[i] @ v
gives the values on the detector segments associated with projectioni
.If
r
is a residual between a projection from spherical to detector space and data from projectioni
, thenM[i].T @ r
gives the associated gradient in spherical harmonic space.
- class mumott.methods.basis_sets.TrivialBasis(channels=1)[source]¶
Basis set class for the trivial basis, i.e., the identity basis. This can be used as a scaffolding class when implementing, e.g., scalar tomography, as it implements all the necessary functionality to qualify as a
BasisSet
.- Parameters
channels (
int
) – Number of channels in the last index. Default is1
. For scalar data, the default value of1
is appropriate. For any other use-case, where the representation on the sphere and the representation in detector space are equivalent, such as reconstructing scalars of multiple q-ranges at once, a different number of channels can be set.
- property csr_representation: tuple¶
The projection matrix as a stack of sparse matrices in CSR representation as a tuple. The information in the tuple consists of the 3 dense matrices making up the representation, in the order
(pointers, indices, data)
.
- forward(coefficients, *args, **kwargs)[source]¶
Returns the provided coefficients with no modification.
- Parameters
coefficients (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array of coefficients, of arbitrary shape, except the last index must specify the same number of channels as was specified for this basis.- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]- Returns
The provided :attr`coefficients` with no modification.
Notes
The
args
andkwargs
are ignored, but included for compatibility with methods that input other arguments.
- generate_map(coefficients, resolution_in_degrees=5, map_half_sphere=True)¶
Generate a (theta, phi) map of the function modeled by the input coefficients. If
map_half_sphere=True
(default) a map of only the z>0 half sphere is returned.- Parameters
coefficients (
ndarray
[Any
,dtype
[float
]]) – One dimensional numpy array with lengthlen(self)
containing the coefficients of the function to be plotted.resolution_in_degrees (
int
) – The resoution of the map in degrees. The map uses eqidistant lines in longitude and latitude.map_half_sphere (
bool
) – If True returns a map of the z>0 half sphere.
- Return type
tuple
[ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]]- Returns
map_intensity – Intensity values of the map.
map_theta – Polar cooridnates of the map.
map_phi – Azimuthal coordinates of the map.
- get_inner_product(u, v)[source]¶
Retrieves the inner product of two coefficient arrays, that is to say, the sum-product over the last axis.
- get_output(coefficients)[source]¶
Returns a dictionary of output data for a given array of coefficients.
- Parameters
coefficients (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array of coefficients of arbitrary shape and dimension. Computations only operate over the last axis ofcoefficents
, so derived properties in the output will have the shape(*coefficients.shape[:-1], ...)
.- Return type
- Returns
A dictionary containing a dictionary with the field
basis_set
.
Notes
In detail, the dictionary under the key
basis_set
contains:- basis_set
- name
The name of the basis set, i.e.,
'TrivialBasis'
- coefficients
A copy of
coefficients
.- projection_matrix
The identity matrix of the same size as the number of chanenls.
- get_spherical_harmonic_coefficients(coefficients, ell_max=None)[source]¶
Convert a set of spherical harmonics coefficients to a different
ell_max
by either zero-padding or truncation and return the result.
- gradient(coefficients, *args, **kwargs)[source]¶
Returns the provided coefficients with no modification.
- Parameters
coefficients (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array of coefficients of arbitrary shape except the last index must specify the same number of channels as was specified for this basis.- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]- Returns
The provided :attr`coefficients` with no modification.
Notes
The
args
andkwargs
are ignored, but included for compatibility with methods that input other argumetns.
- property integration_mode: str¶
Mode of integration for calculating projection matrix. Accepted values are
'simpson'
,'romberg'
,'trapezoid'
, and'midpoint'
.
- property probed_coordinates: ProbedCoordinates¶
- property projection_matrix¶
The identity matrix of the same rank as the number of channels specified.
- class mumott.methods.basis_sets.GaussianKernels(probed_coordinates=None, grid_scale=4, kernel_scale_parameter=1.0, enforce_friedel_symmetry=True, **kwargs)[source]¶
Basis set class for gaussian kernels, a simple local representation on the sphere. The kernels follow a pseudo-even distribution similar to that described by Y. Kurihara in 1965, except with offsets added at the poles.
Notes
The Gaussian kernel at location \(\rho_i\) is given by
\[N_i \exp\left[ -\frac{1}{2} \left(\frac{d(\rho_i, r)}{\sigma}\right)^2 \right]\]\[\sigma = \frac{\nu \pi}{2 (g + 1)}\]where \(\nu\) is the kernel scale parameter and \(g\) is the grid scale, and
\[d(\rho, r) = \arctan_2(\Vert \rho \times r \Vert, \rho \cdot r),\]that is, the great circle distance from the kernel location \(\rho\) to the probed location \(r\). If Friedel symmetry is assumed, the expression is instead
\[d(\rho, r) = \arctan_2(\Vert \rho \times r \Vert, \vert \rho \cdot r \vert)\]The normalization factor \(\rho_i\) is given by
\[N_i = \sum_j \exp\left[ -\frac{1}{2} \left( \frac{d(\rho_i, \rho_j)}{\sigma} \right)^2 \right]\]where the sum goes over the coordinates of all grid points. This leads to an approximately even spherical function, such that a set of coefficients which are all equal is approximately isotropic, to the extent possible with respect to restrictions imposed by grid resolution and scale parameter.
- Parameters
probed_coordinates (ProbedCoordinates) – Optional. A container with the coordinates on the sphere probed at each detector segment by the experimental method. Its construction from the system geometry is method-dependent. By default, an empty instance of
mumott.ProbedCoordinates
is created.grid_scale (int) – The size of the coordinate grid on the sphere. Denotes the number of azimuthal rings between the pole and the equator, where each ring has between
2
and2 * grid_scale
points along the azimuth.kernel_scale_parameter (float) – The scale parameter of the kernel in units of \(\frac{\pi}{2 (g + 1)}\), where \(g\) is
grid_scale
.enforce_friedel_symmetry (bool) – If set to
True
, Friedel symmetry will be enforced, using the assumption that points on opposite sides of the sphere are equivalent.kwargs –
Miscellaneous arguments which relate to segment integrations can be passed as keyword arguments:
- integration_mode
Mode to integrate line segments on the reciprocal space sphere. Possible options are
'simpson'
,'midpoint'
,'romberg'
,'trapezoid'
.'simpson'
,'trapezoid'
, and'romberg'
use adaptive integration with the respective quadrature rule fromscipy.integrate
.'midpoint'
uses a single mid-point approximation of the integral. Default value is'simpson'
.- n_integration_starting_points
Number of points used in the first iteration of the adaptive integration. The number increases by the rule
N
←2 * N - 1
for each iteration. Default value is 3.- integration_tolerance
Tolerance for the maximum relative error between iterations before the integral is considered converged. Default is
1e-5
.- integration_maxiter
Maximum number of iterations. Default is
10
.- enforce_sparsity
If
True
, makes matrix sparse by limiting the number of basis set elements that can map to each segment. Default isFalse
.- sparsity_count
Number of basis set elements that can map to each segment, if
enforce_sparsity
is set toTrue
. Default is3
.
- property csr_representation: tuple¶
The projection matrix as a stack of sparse matrices in CSR representation as a tuple. The information in the tuple consists of the 3 dense matrices making up the representation, in the order
(pointers, indices, data)
.
- property enforce_friedel_symmetry: bool¶
If
True
, Friedel symmetry is enforced, i.e., the point \(-r\) is treated as equivalent to \(r\).
- forward(coefficients, indices=None)[source]¶
Carries out a forward computation of projections from Gaussian kernel space into detector space, for one or several tomographic projections.
- Parameters
coefficients (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array of coefficients, of arbitrary shape so long as the last axis has the same size askernel_scale_parameter
, and ifindices
isNone
or greater than one, the first axis should have the same length asindices
indices (
Optional
[ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]]) – Optional. Indices of the tomographic projections for which the forward computation is to be performed. IfNone
, the forward computation will be performed for all projections.
- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]- Returns
An array of values on the detector corresponding to the
coefficients
given. Ifindices
contains exactly one index, the shape is(coefficients.shape[:-1], J)
whereJ
is the number of detector segments. Ifindices
isNone
or contains several indices, the shape is(N, coefficients.shape[1:-1], J)
whereN
is the number of tomographic projections for which the computation is performed.
Notes
The assumption is made in this implementation that computations over several indices act on sets of images from different projections. For special usage where multiple projections of entire fields are desired, it may be better to use
projection_matrix
directly. This also applies togradient()
.
- generate_map(coefficients, resolution_in_degrees=5, map_half_sphere=True)¶
Generate a (theta, phi) map of the function modeled by the input coefficients. If
map_half_sphere=True
(default) a map of only the z>0 half sphere is returned.- Parameters
coefficients (
ndarray
[Any
,dtype
[float
]]) – One dimensional numpy array with lengthlen(self)
containing the coefficients of the function to be plotted.resolution_in_degrees (
int
) – The resoution of the map in degrees. The map uses eqidistant lines in longitude and latitude.map_half_sphere (
bool
) – If True returns a map of the z>0 half sphere.
- Return type
tuple
[ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]]- Returns
map_intensity – Intensity values of the map.
map_theta – Polar cooridnates of the map.
map_phi – Azimuthal coordinates of the map.
- get_amplitudes(coefficients, probed_coordinates=None)[source]¶
Computes the amplitudes of the spherical function represented by the provided
coefficients
at theprobed_coordinates
.- Parameters
coefficients (
ndarray
[Any
,dtype
[float
]]) – An array of coefficients of arbitrary shape, provided that the last dimension contains the coefficients for one spherical function.probed_coordinates (
Optional
[ProbedCoordinates
]) – An instance ofmumott.core.ProbedCoordinates
with itsvector
attribute indicating the points of the sphere for which to evaluate the amplitudes.
- Return type
- get_inner_product(u, v)[source]¶
Retrieves the inner product of two coefficient arrays, that is to say, the sum-product over the last axis.
- Parameters
u (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – The first coefficient array, of arbitrary shape and dimension, so long as the number of coefficients equals the length of thisGaussianKernels
instance.v (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – The second coefficient array, of the same shape asu
.
- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]
- get_output(coefficients)[source]¶
Returns a dictionary of output data for a given array of basis set coefficients.
- Parameters
coefficients (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array of coefficients of arbitrary shape and dimensions, except its last dimension must be the same length as thelen
of this instance. Computations only operate over the last axis ofcoefficients
, so derived properties in the output will have the shape(*coefficients.shape[:-1], ...)
.- Return type
- Returns
A dictionary containing two sub-dictionaries,
basis_set
andspherical_harmonic_analysis
.basis_set
contains information particular toGaussianKernels
, whereasspherical_harmonic_analysis
contains an analysis of the spherical function using a spherical harmonic transform.
Notes
In detail, the two sub-dictionaries
basis_set
andspherical_harmonic_analysis
have the following members:- basis_set
- name
The name of the basis set, i.e.,
'GaussianKernels'
- coefficients
A copy of
coefficients
.- grid_scale
A copy of
grid_scale
.- kernel_Scale_paramter
A copy of
kernel_scale_parameter
.- enforce_friedel_symmetry
A copy of
enforce_friedel_symmetry
.- projection_matrix
A copy of
projection_matrix
.
- spherical_harmonic_analysis
An analysis of the spherical function in terms of spherical harmonics. See
SphericalHarmonics.get_output
for details.
- get_spherical_harmonic_coefficients(coefficients, ell_max=None)[source]¶
Computes the spherical harmonic coefficients of the spherical function represented by the provided
coefficients
using a Driscoll-Healy grid.For details on the Driscoll-Healy grid, see the SHTools page for a comprehensive overview.
- Parameters
coefficients (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array of coefficients of arbitrary shape, provided that the last dimension contains the coefficients for one function.ell_max (
Optional
[int
]) – The bandlimit of the spherical harmonic expansion. By default, it is2 * grid_scale
.
- gradient(coefficients, indices=None)[source]¶
Carries out a gradient computation of projections from Gaussian kernel space into detector space for one or several tomographic projections.
- Parameters
coefficients (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array of coefficients (or residuals) of arbitrary shape so long as the last axis has the same size as the number of detector segments.indices (
Optional
[ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]]) – Optional. Indices of the tomographic projections for which the gradient computation is to be performed. IfNone
, the gradient computation will be performed for all projections.
- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]- Returns
An array of gradient values based on the
coefficients
given. Ifindices
contains exactly one index, the shape is(coefficients.shape[:-1], J)
whereJ
is the number of detector segments. If indices isNone
or contains several indices, the shape is(N, coefficients.shape[1:-1], J)
whereN
is the number of tomographic projections for which the computation is performed.
Notes
When solving an inverse problem, one should not attempt to optimize the coefficients directly using the gradient one obtains by applying this method to the data. Instead, one must either take the gradient of the residual between the
forward()
computation of the coefficients and the data. Alternatively one can apply both the forward and the gradient computation to the coefficients to be optimized, and the gradient computation to the data, and treat the residual of the two as the gradient of the optimization coefficients. The approaches are algebraically equivalent, but one may be more efficient than the other in some circumstances. However, normally, the projection between detector andGaussianKernel
space is only a small part of the overall computation, so there is typically not much to be gained from optimizing it.
- property grid: Tuple[ndarray[Any, dtype[float]], ndarray[Any, dtype[float]]]¶
Returns the polar and azimuthal angles of the grid used by the basis.
- Returns
A
Tuple
with contents(polar_angle, azimuthal_angle)
, where thepolar angle is defined as \(\arccos(z)\).
- property grid_scale: int¶
The number of azimuthal rings from each pole to the equator in the spherical grid.
- property integration_mode: str¶
Mode of integration for calculating projection matrix. Accepted values are
'simpson'
,'romberg'
,'trapezoid'
, and'midpoint'
.
- property probed_coordinates: ProbedCoordinates¶
- property projection_matrix: ndarray[Any, dtype[_ScalarType_co]]¶
The matrix used to project spherical functions from the unit sphere onto the detector. If
v
is a vector of gaussian kernel coefficients, andM
is theprojection_matrix
, thenM @ v
gives the corresponding values on the detector segments associated with each projection.M[i] @ v
gives the values on the detector segments associated with projectioni
.If
r
is a residual between a projection from Gaussian kernel to detector space and data from projectioni
, thenM[i].T @ r
gives the associated gradient in Gaussian kernel space.
- property projection_matrix_hash: str¶
Returns a hash of
projection_matrix
.
- class mumott.methods.basis_sets.NearestNeighbor(directions, probed_coordinates=None, enforce_friedel_symmetry=True, **kwargs)[source]¶
Basis set class for nearest-neighbor interpolation. Used to construct methods similar to that presented in Schaff et al. (2015). By default this representation is sparse and maps only a single direction on the sphere to each detector segment. This can be changed; see
kwargs
.- Parameters
directions (NDArray[float]) – Two-dimensional Array containing the
N
sensitivity directions with shape(N, 3)
.probed_coordinates (ProbedCoordinates) – Optional. Coordinates on the sphere probed at each detector segment by the experimental method. Its construction from the system geometry is method-dependent. By default, an empty instance of
mumott.ProbedCoordinates
is created.enforce_friedel_symmetry (bool) – If set to
True
, Friedel symmetry will be enforced, using the assumption that points on opposite sides of the sphere are equivalent.kwargs –
Miscellaneous arguments which relate to segment integrations can be passed as keyword arguments:
- integration_mode
Mode to integrate line segments on the reciprocal space sphere. Possible options are
'simpson'
,'midpoint'
,'romberg'
,'trapezoid'
.'simpson'
,'trapezoid'
, and'romberg'
use adaptive integration with the respective quadrature rule fromscipy.integrate
.'midpoint'
uses a single mid-point approximation of the integral. Default value is'simpson'
.- n_integration_starting_points
Number of points used in the first iteration of the adaptive integration. The number increases by the rule
N
←2 * N - 1
for each iteration. Default value is 3.- integration_tolerance
Tolerance for the maximum relative error between iterations before the integral is considered converged. Default is
1e-3
.- integration_maxiter
Maximum number of iterations. Default is
10
.- enforce_sparsity
If
True
, limites the number of basis set elements that can map to each detector segemnt. Default isFalse
.- sparsity_count
If
enforce_sparsity
is set toTrue
, the number of basis set elements that can map to each detector segment. Default value is1
.
- property csr_representation: tuple¶
The projection matrix as a stack of sparse matrices in CSR representation as a tuple. The information in the tuple consists of the 3 dense matrices making up the representation, in the order
(pointers, indices, data)
.
- property enforce_friedel_symmetry: bool¶
If
True
, Friedel symmetry is enforced, i.e., the point \(-r\) is treated as equivalent to \(r\).
- find_nearest_neighbor_index(probed_directions)[source]¶
Caluculate the nearest neighbor sensitivity directions for an array of x-y-z vectors.
- forward(coefficients, indices=None)[source]¶
Carries out a forward computation of projections from reciprocal space modes to detector channels, for one or several tomographic projections.
- Parameters
coefficients (
ndarray
[Any
,dtype
[float
]]) – An array of coefficients, of arbitrary shape so long as the last axis has the same size as this basis set.indices (
Optional
[ndarray
[Any
,dtype
[int
]]]) – Optional. Indices of the tomographic projections for which the forward computation is to be performed. IfNone
, the forward computation will be performed for all projections.
- Return type
- Returns
An array of values on the detector corresponding to the
coefficients
given. Ifindices
contains exactly one index, the shape is(coefficients.shape[:-1], J)
whereJ
is the number of detector segments. Ifindices
isNone
or contains several indices, the shape is(N, coefficients.shape[1:-1], J)
whereN
is the number of tomographic projections for which the computation is performed.
- generate_map(coefficients, resolution_in_degrees=5, map_half_sphere=True)¶
Generate a (theta, phi) map of the function modeled by the input coefficients. If
map_half_sphere=True
(default) a map of only the z>0 half sphere is returned.- Parameters
coefficients (
ndarray
[Any
,dtype
[float
]]) – One dimensional numpy array with lengthlen(self)
containing the coefficients of the function to be plotted.resolution_in_degrees (
int
) – The resoution of the map in degrees. The map uses eqidistant lines in longitude and latitude.map_half_sphere (
bool
) – If True returns a map of the z>0 half sphere.
- Return type
tuple
[ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]]- Returns
map_intensity – Intensity values of the map.
map_theta – Polar cooridnates of the map.
map_phi – Azimuthal coordinates of the map.
- get_amplitudes(coefficients, probed_directions)[source]¶
Calculate function values of an array of coefficients.
- Parameters
- Return type
- Returns
Array with function values. The shape of the array is
(*coefficients.shape[:-1], *probed_directions.shape[:-1])
.
- get_function_values(probed_directions)[source]¶
Calculate the value of the basis functions from an array of x-y-z vectors.
- Parameters
probed_directions (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – Array with length 3 along its last axis- Return type
- Returns
Array with same shape as input array except for the last axis, which now has length
N
, i.e., the number of sensitivity directions.
- get_output(coefficients)[source]¶
Returns a dictionary of output data for a given array of basis set coefficients.
- Parameters
coefficients (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array of coefficients of arbitrary shape and dimensions, except its last dimension must be the same length as thelen
of this instance. Computations only operate over the last axis ofcoefficients
, so derived properties in the output will have the shape(*coefficients.shape[:-1], ...)
.- Return type
- Returns
A dictionary containing information about the optimized function.
- get_second_moments(coefficients)[source]¶
Calculate the second moments of the functions described by
coefficients
.- Parameters
coefficients (
ndarray
[Any
,dtype
[float
]]) – An array of coefficients (or residuals) of arbitrary shape so long as the last axis has the same size as the number of detector channels.- Return type
- Returns
Array containing the second moments of the functions described by coefficients, formatted as rank-two tensors with tensor indices in the last 2 dimensions.
- get_spherical_harmonic_coefficients(coefficients, ell_max=None)[source]¶
Computes and rturns the spherical harmonic coefficients of the spherical function represented by the provided
coefficients
using a Driscoll-Healy grid.For details on the Driscoll-Healy grid, see the SHTools page for a comprehensive overview.
- get_sub_geometry(direction_index, geometry, data_container=None)[source]¶
Create and return a geometry object corresponding to a scalar tomography problem for scattering along the sensitivity direction with index
direction_index
. If optionally amumott.DataContainer
is provided, the sinograms and weights for this scalar tomography problem will alse be returned.Used for an implementation of the algorithm descibed in [Schaff2015].
- gradient(coefficients, indices=None)[source]¶
Carries out a gradient computation of projections of projections from reciprocal space modes to detector channels, for one or several tomographic projections.
- Parameters
coefficients (
ndarray
[Any
,dtype
[float
]]) – An array of coefficients (or residuals) of arbitrary shape so long as the last axis has the same size as the number of detector channels.indices (
Optional
[ndarray
[Any
,dtype
[int
]]]) – Optional. Indices of the tomographic projections for which the gradient computation is to be performed. IfNone
, the gradient computation will be performed for all projections.
- Return type
- Returns
An array of gradient values based on the
coefficients
given. Ifindices
contains exactly one index, the shape is(coefficients.shape[:-1], J)
whereJ
is the number of detector segments. If indices isNone
or contains several indices, the shape is(N, coefficients.shape[1:-1], J)
whereN
is the number of tomographic projections for which the computation is performed.
- property grid: Tuple[ndarray[Any, dtype[float]], ndarray[Any, dtype[float]]]¶
Returns the polar and azimuthal angles of the grid used by the basis.
- Returns
A
Tuple
with contents(polar_angle, azimuthal_angle)
, where thepolar angle is defined as \(\arccos(z)\).
- property integration_mode: str¶
Mode of integration for calculating projection matrix. Accepted values are
'simpson'
,'romberg'
,'trapezoid'
, and'midpoint'
.
- property probed_coordinates: ProbedCoordinates¶
- property projection_matrix: ndarray[Any, dtype[_ScalarType_co]]¶
The matrix used to project spherical functions from the unit sphere onto the detector. If
v
is a vector of gaussian kernel coefficients, andM
is theprojection_matrix
, thenM @ v
gives the corresponding values on the detector segments associated with each projection.M[i] @ v
gives the values on the detector segments associated with projectioni
.
- property projection_matrix_hash: str¶
Returns a hash of
projection_matrix
.