Residual calculators

class mumott.methods.residual_calculators.GradientResidualCalculator(data_container, basis_set, projector, use_scalar_projections=False, scalar_projections=None)

Class that implements the GradientResidualCalculator method. This residual calculator is an appropriate choice for SAXS tensor tomography, as it relies on the small-angle approximation. It relies on inverting the John transform (also known as the X-ray transform) of a tensor field (where each tensor is a representation of a spherical function) by comparing it to scattering data which has been corrected for transmission.

Parameters

• data_container (DataContainer) – Container for the data which is to be reconstructed.

• basis_set (BasisSet) – The basis set used for representing spherical functions.

• projector (Projector) – The type of projector used together with this method.

• use_scalar_projections (bool) – Whether to use a set of scalar projections, rather than the data in data_container.

• scalar_projections (NDArray[float]) – If use_scalar_projections is true, the set of scalar projections to use. Should have the same shape as data_container.data, except with only one channel in the last index.

property coefficients: ndarray[Any, dtype[_ScalarType_co]]

Optimization coefficients for this method.

property dtype: dtype

dtype used by the projector object attached to this instance.

get_estimate_of_matrix_norm()

Calculate the matrix norm of the matrix implied by the basis set and projector associated with this residual calculator.

Return type

An estimate of the matrix norm (largest singular value)

get_gradient_from_residual_gradient(residual_gradient)

Projects a residual gradient into coefficient and volume space. Used to get gradients from more complicated residuals, e.g., the Huber loss. Assumes that any weighting to the residual gradient has already been applied.

Parameters

residual_gradient (ndarray[Any, dtype[float]]) – The residual gradient, from which to calculate the gradient.

Return type

Dict

Returns

An NDArray containing the gradient.

get_residuals(get_gradient=False, get_weights=False, gradient_part=None)

Calculates a residuals and possibly a gradient between coefficients projected using the BasisSet and Projector attached to this instance.

Parameters

• get_gradient (bool) – Whether to return the gradient. Default is False.

• get_weights (bool) – Whether to return weights. Default is False. If True along with get_gradient, the gradient will be computed with weights.

• gradient_part (Optional[str]) – Used for the zonal harmonics resonstructions to determine what part of the gradient is being calculated. Default is None. Raises a NotImplementedError for any other value.

Return type

dict[str, ndarray[Any, dtype[float]]]

Returns

A dictionary containing the residuals, and possibly the gradient and/or weights. If gradient and/or weights are not returned, their value will be None.

property is_dirty: bool

True if stored hashes of geometry or basis set objects do not match their current hashes. Used to trigger updates

property probed_coordinates: ndarray[Any, dtype[_ScalarType_co]]

An array of 3-vectors with the (x, y, z)-coordinates on the reciprocal space map probed by the method. Structured as (N, K, I, 3), where N is the number of projections, K is the number of detector segments, I is the number of points to be integrated over, and the last axis contains the (x, y, z)-coordinates.

Notes

The region of the reciprocal space map spanned by each detector segment is represented as a parametric curve between each segment. This is intended to simulate the effect of summing up pixels on a detector screen. For other methods of generating data (e.g., by fitting measurements to a curve), it may be more appropriate to only include a single point, which will then have the same coordinates as the center of a detector segments. This can be achieved by setting the property integration_samples.

The number of detector segments is len(geometry.detecor_angles)*len(geometry.two_theta) i.e. the product of the number of two_theta bins times the number of azimuthal bins. As a default, only on two theta bin is used. When several two_theta bins are used, the second index corresponds to a raveled array, where the azimuthal is the fast index and two theta is the slow index.

class mumott.methods.residual_calculators.ZHTTResidualCalculator(data_container, basis_set, projector)

Class that implements the gradient calculations for a model that uses a SphericalHarmonics basis set restricted to zonal harmonics parametrized by a primary axis with polar coordinates \(\theta_0\) and \(\phi_0\) ,defined as:

\[\begin{split}\begin{pmatrix} x_0\\ y_0\\ z_0\end{pmatrix} = \begin{pmatrix} \sin(\theta_0) \sin(\phi_0) \\ \sin(\theta_0) \cos(\phi_0) \\ \cos(\theta_0) \end{pmatrix}\end{split}\]

This model is equivalent to the one used in [Liebi2015], but uses a different approach to computation.

This implementation avoids doing some of the expensive calculations of trigonometric functions and Legendre polynomials by doing the rotation in the space of the spherical harmonics using Wigner (small) d-matrices. The forward model only involves a small number of trigonometric functions to evaluate the \(d_z(\text{angle})\) matrices for the \(\theta\) and \(\phi\) rotations. Everything else is expressed as matrix products with precomputed matrices.

The full forward model may be written as:

\[\boldsymbol{I} = \boldsymbol{W} \boldsymbol{P} \boldsymbol{d}_z(\phi_0) \boldsymbol{d}_y(\frac{\pi}{4})^T \boldsymbol{d}_z(-\theta_0) \boldsymbol{d}_y(\frac{\pi}{4}) \boldsymbol{a}'_{l0},\]

where \(\boldsymbol{W}\) is the mapping from spherical harmonic modes to detector segments, which can be precomputed. \(\boldsymbol{P}\) is the typical projector from normal 3D tomography and \(\boldsymbol{d}_i(\text{angle})\) with \(i = x,y,z\) are Wigner (small) d matrices for real spherical harmonics. \(\theta_0\), \(\phi_0\), and \(\boldsymbol{a}_{l0}\) are the model parameters for each voxel.

Derivatives are easy to evaluate because the angles only appear in the \(\boldsymbol{d}_z(\text{angle})\)-matrices. All the expensive trigonometric and spherical harmonics calculations have been put into the precomputation of \(\boldsymbol{W}\) and \(\boldsymbol{d}_y(\frac{\pi}{4})\).

Parameters

• data_container (DataContainer) – Container holding the data to be reconstructed.

• basis_set (SphericalHarmonics) – The basis set used for representing spherical functions.

• projector (Projector) – The type of projector used together with this method.

property coefficients: ndarray[Any, dtype[_ScalarType_co]]

Optimization coefficients for this method. Contains both the zonal coefficients and the angles. The first N-2 elements are zonal coefficients. The N-1th element is the polar angle and the last element is the azimuthal angle.

property directions

Returns the direction of symmetry as a unit vector in in xyz coordinates. The vector index is the last index of the output.

property dtype: dtype

dtype used by the projector object attached to this instance.

property ell_max: int

l max

get_estimate_of_matrix_norm()

Calculate the matrix norm of the matrix implied by the basis set and projector associated with this residual calculator.

Return type

An estimate of the matrix norm (largest singular value)

get_gradient(get_weights=False, gradient_part='full')

Calculates the gradient of half the sum of residuals squared.

Parameters

• get_gradient – Whether to return the gradient. Default is False.

• gradient_part (str) – If gradient_part is 'full' (Default) the gradient is computed with respect to all parameters; if gradient_part is 'angles' only the gradient with respect to the angles is computed; if gradient_part is 'coefficients' only the gradient with respect to the zonal spherical harmonics coefficients is computed.

Return type

dict[str, ndarray[Any, dtype[float]]]

Returns

A dictionary containing the residuals of the gradient. If only a part of the gradient is computed, the rest of the elements will be filled with zeros.

get_residuals(get_gradient=False, get_weights=False, gradient_part='full')

Calculates the residuals and possibly the gradient of the residual square sum (without the factor of -2!) with respect to the parameters. The coefficients are projected using the SphericalHarmonics and Projector attached to this instance.

Parameters

• get_gradient (bool) – Whether to return the gradient. Default is False.

• get_weights (bool) – Whether to return weights. Default is False. If True along with get_gradient, the gradient will be computed with weights.

• gradient_part (str) – If gradient_part is 'full' (Default) the gradient is computed with respect to all parameters; if gradient_part is 'angles' only the gradient with respect to the angles is computed; if gradient_part is 'coefficients' only the gradient with respect to the zonal spherical harmonics coefficients is computed.

Return type

dict[str, ndarray[Any, dtype[float]]]

Returns

A dictionary containing the residuals, and possibly the gradient and/or weights. If gradient and/or weights are not returned, their value will be None.

property is_dirty: bool

True if stored hashes of geometry or basis set objects do not match their current hashes. Used to trigger updates

property probed_coordinates: ndarray[Any, dtype[_ScalarType_co]]

An array of 3-vectors with the (x, y, z)-coordinates on the reciprocal space map probed by the method. Structured as (N, K, I, 3), where N is the number of projections, K is the number of detector segments, I is the number of points to be integrated over, and the last axis contains the (x, y, z)-coordinates.

Notes

The region of the reciprocal space map spanned by each detector segment is represented as a parametric curve between each segment. This is intended to simulate the effect of summing up pixels on a detector screen. For other methods of generating data (e.g., by fitting measurements to a curve), it may be more appropriate to only include a single point, which will then have the same coordinates as the center of a detector segments. This can be achieved by setting the property integration_samples.

The number of detector segments is len(geometry.detecor_angles)*len(geometry.two_theta) i.e. the product of the number of two_theta bins times the number of azimuthal bins. As a default, only on two theta bin is used. When several two_theta bins are used, the second index corresponds to a raveled array, where the azimuthal is the fast index and two theta is the slow index.

property rotated_coefficients

Returns the real spherical harmonics coefficients.

property volume_shape: int

Shape of voxel volume