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 even ell 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 from scipy.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 when ell_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 when ell_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 as ell_indices, and if indices is None or greater than one, the first axis should have the same length as indices
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. If None, 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. If indices contains exactly one index, the shape is (coefficients.shape[:-1], J) where J is the number of detector segments. If indices is None or contains several indices, the shape is (N, coefficients.shape[1:-1], J) where N 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 to gradient().

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 length len(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 of ell_indices.
v (ndarray[Any, dtype[TypeVar(_ScalarType_co, bound= generic, covariant=True)]]) – The second coefficient array, of the same shape as u.
resolve_spectrum (bool) – Optional. Whether to resolve the product according to each frequency band, given by the coefficients of each ell in ell_indices. Default value is False.

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 and v. Has the shape (u.shape[:-1]) if resolve_spectrum is False, and (u.shape[:-1] + (ell_max // 2 + 1,)) if resolve_spectrum is True, where ell_max is SphericalHarmonics.ell_max.

Notes

Calling this function is equivalent to calling get_inner_product() with spectral_moments=np.unique(ell_indices[ell_indices > 0]) where ell_indices is SphericalHarmonics.ell_indices. See the note to get_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 of ell_indices.
v (ndarray[Any, dtype[TypeVar(_ScalarType_co, bound= generic, covariant=True)]]) – The second coefficient array, of the same shape as u.
resolve_spectrum (bool) – Optional. Whether to resolve the product according to each frequency band, given by the coefficients of each ell in ell_indices. Defaults to False, which means that the sum of every component of the spectrum is returned. If True, components are returned in order of ascending ell. The ell included in the spectrum depends on spectral_moments.
spectral_moments (Optional[List[int]]) – Optional. List of particular values of ell to calculate the inner product for. Defaults to None, which is identical to including all values of ell in the calculation. If spectral_moments contains all nonzero values of ell and resolve_spectrum is False, the covariance of v and u will be calculated (the sum of the inner product over all non-zero ell If resolve_spectrum is True, the covariance per ell in spectral_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 and v. Has the shape (u.shape[:-1]) if resolve_spectrum is False, (u.shape[:-1] + (ell_max // 2 + 1,)) if resolve_spectrum is True and spectral_moments is None, and finally the shape (u.shape[:-1] + (np.unique(spectral_moments).size,)) if resolve_spectrum is True and spectral_moments is a list of integers found in ell_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 as ell_indices. Computations only operate over the last axis of coefficients, so derived properties in the output will have the shape (*coefficients.shape[:-1], ...).
Return type: Dict[str, Any]
Returns: A dictionary containing two sub-dictionaries, basis_set and spherical_functions. basis_set contains information particular to SphericalHarmonics, whereas spherical_functions contains information about the spherical functions represented by the coefficients which are not specific to the chosen representation.

Notes

In detail, the two sub-dictionaries basis_set and spherical_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. If ell_max is 0, 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]]. If ell_max is 0, 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 is 0, 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, and eigenvectors[..., 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. If ell_max is 0, the eigenvectors will be the Cartesian basis vectors.
main_orientations: The estimated main orientations from the largest absolute eigenvalues. If ell_max is 0, 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 is 1 (orientation represents a maximum). If ell_max is 0, 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. If ell_max is 0, the normalized standard deviations will always be zero.
power_spectra: The spectral powers of each ell in ell_indices, for each spherical function, sorted in ascending ell. If ell_max is 0, 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 in power_spectra. Equal to np.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.

Parameters

coefficients (ndarray[Any, dtype[float]]) – An array of coefficients of arbitrary shape, provided that the last dimension contains the coefficients for one function.
ell_max (Optional[int]) – The band limit of the spherical harmonic expansion.

Return type

ndarray[Any, dtype[float]]

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. If None, 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) where J is the number of detector segments. If indices is None or contains several indices, the shape is (N, coefficients.shape[1:-1], J) where N 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 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.

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, and M is the projection_matrix, then M @ 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 projection i.

If r is a residual between a projection from spherical to detector space and data from projection i, then M[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 is 1. For scalar data, the default value of 1 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 channels: int¶: The number of channels this basis supports.

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 and kwargs 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 length len(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.

Parameters

u (ndarray[Any, dtype[TypeVar(_ScalarType_co, bound= generic, covariant=True)]]) – The first coefficient array, of arbitrary shape and dimension.
v (ndarray[Any, dtype[TypeVar(_ScalarType_co, bound= generic, covariant=True)]]) – The second coefficient array, of the same shape as u.

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 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 of coefficents, so derived properties in the output will have the shape (*coefficients.shape[:-1], ...).
Return type: Dict[str, Any]
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.

Parameters

coefficients (ndarray[Any, dtype[float]]) – An array of coefficients of arbitrary shape, provided that the last dimension contains the coefficients for one function.
ell_max (Optional[int]) – The band limit of the spherical harmonic expansion.

Return type

ndarray[Any, dtype[float]]

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 and kwargs 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 is_dirty: bool¶

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 and 2 * 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 from scipy.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 is False.

sparsity_count
Number of basis set elements that can map to each segment, if enforce_sparsity is set to True. Default is 3.

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 as kernel_scale_parameter, and if indices is None or greater than one, the first axis should have the same length as indices
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. If None, 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. If indices contains exactly one index, the shape is (coefficients.shape[:-1], J) where J is the number of detector segments. If indices is None or contains several indices, the shape is (N, coefficients.shape[1:-1], J) where N 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 to gradient().

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 length len(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 the probed_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 of mumott.core.ProbedCoordinates with its vector attribute indicating the points of the sphere for which to evaluate the amplitudes.

Return type

ndarray[Any, dtype[float]]

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 this GaussianKernels instance.
v (ndarray[Any, dtype[TypeVar(_ScalarType_co, bound= generic, covariant=True)]]) – The second coefficient array, of the same shape as u.

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 the len of this instance. Computations only operate over the last axis of coefficients, so derived properties in the output will have the shape (*coefficients.shape[:-1], ...).
Return type: Dict[str, Any]
Returns: A dictionary containing two sub-dictionaries, basis_set and spherical_harmonic_analysis. basis_set contains information particular to GaussianKernels, whereas spherical_harmonic_analysis contains an analysis of the spherical function using a spherical harmonic transform.

Notes

In detail, the two sub-dictionaries basis_set and spherical_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 is 2 * 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. If None, 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) where J is the number of detector segments. If indices is None or contains several indices, the shape is (N, coefficients.shape[1:-1], J) where N 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 and GaussianKernel 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 the
polar angle is defined as \(\arccos(z)\).

property grid_hash: str¶: Returns a hash of grid.

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 is_dirty: bool¶

property kernel_scale_parameter: float¶: The scale parameter for each kernel.

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, and M is the projection_matrix, then M @ 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 projection i.

If r is a residual between a projection from Gaussian kernel to detector space and data from projection i, then M[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 from scipy.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 is False.

sparsity_count
If enforce_sparsity is set to True, the number of basis set elements that can map to each detector segment. Default value is 1.

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.

Parameters: probed_directions (ndarray[Any, dtype[float]]) – Array with length 3 along its last axis
Return type: ndarray[Any, dtype[int]]
Returns: Array with same shape as the input except for the last dimension, which contains the index of the nearest-neighbor sensitivity direction.

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. If None, the forward computation will be performed for all projections.

Return type

ndarray[Any, dtype[float]]

Returns

An array of values on the detector corresponding to the coefficients given. If indices contains exactly one index, the shape is (coefficients.shape[:-1], J) where J is the number of detector segments. If indices is None or contains several indices, the shape is (N, coefficients.shape[1:-1], J) where N 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 length len(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

coefficients (ndarray[Any, dtype[float]]) – Array of coefficients with coefficient number along its last index.
probed_directions (ndarray[Any, dtype[float]]) – Array with length 3 along its last axis.

Return type

ndarray[Any, dtype[float]]

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: ndarray[Any, dtype[float]]
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 the len of this instance. Computations only operate over the last axis of coefficients, so derived properties in the output will have the shape (*coefficients.shape[:-1], ...).
Return type: Dict[str, Any]
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: ndarray[Any, dtype[float]]
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.

Parameters

coefficients (ndarray[Any, dtype[float]]) – 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.

Return type

ndarray[Any, dtype[float]]

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 a mumott.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].

Parameters

direction_index (int) – Index of the sensitivity direction.
geometry (Geometry) – mumott.Geometry object of the full problem.
(optional) (data_container) – mumott.DataContainer compatible with Geometry from which a scalar dataset will be constructed.

Return type

tuple[Geometry, tuple[ndarray[Any, dtype[float]], ndarray[Any, dtype[float]]]]

Returns

sub_geometry – Geometry of the scalar problem.
data_tuple – Tuple containing two numpy arrays. data_tuple[0] is the data of the scalar problem. data_tuple[1] are the weights.

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. If None, the gradient computation will be performed for all projections.

Return type

ndarray[Any, dtype[float]]

Returns

An array of gradient values based on the coefficients given. If indices contains exactly one index, the shape is (coefficients.shape[:-1], J) where J is the number of detector segments. If indices is None or contains several indices, the shape is (N, coefficients.shape[1:-1], J) where N 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 the
polar angle is defined as \(\arccos(z)\).

property grid_hash: str¶: Returns a hash of grid.

property integration_mode: str¶: Mode of integration for calculating projection matrix. Accepted values are 'simpson', 'romberg', 'trapezoid', and 'midpoint'.

property is_dirty: bool¶

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, and M is the projection_matrix, then M @ 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 projection i.

property projection_matrix_hash: str¶: Returns a hash of projection_matrix.