Projectors¶
mumott provides four different implementations of projectors suitable for SAXS tomography. While they should yield approximately equivalent results, they differ with respect to the resources they require.
SAXSProjectorCUDA
and SAXSProjectorBilinear
implement an equivalent algorithm for GPU and CPU resources, respectively.
These two options are recommended in terms of speed and ease of use.
SAXSProjectorAstra
provides an interface to a GPU implementation in the ASTRA toolbox.
Finally, SAXSProjectorNumba
provides a CPU implementation using the numba library.
- class mumott.methods.projectors.SAXSProjectorBilinear(geometry)[source]¶
Projector for transforms of tensor fields from three-dimensional space to projection space using a bilinear interpolation algorithm that produces results similar to those of
SAXSProjectorCUDA
using CPU computation.- Parameters
geometry (Geometry) – An instance of
Geometry
containing the necessary vectors to compute forwared and adjoint projections.
- adjoint(projections, indices=None)[source]¶
Compute the adjoint of a set of projections according to the system geometry.
- Parameters
projections (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array containing coefficients in its last dimension, from e.g. the residual of measured data and forward projections. The first dimension should matchindices
in size, and the second and third dimensions should match the system projection geometry. The array must be contiguous and row-major.indices (
Optional
[ndarray
[Any
,dtype
[int
]]]) – A one-dimensional array containing one or more indices indicating from which projections the adjoint is to be computed.
- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]- Returns
The adjoint of the provided projections. An array with four dimensions
(X, Y, Z, P)
, where the first three dimensions are spatial and the last dimension runs over coefficients.
- property dtype: Union[dtype[Any], None, type[Any], _SupportsDType[dtype[Any]], str, tuple[Any, int], tuple[Any, Union[SupportsIndex, collections.abc.Sequence[SupportsIndex]]], list[Any], _DTypeDict, tuple[Any, Any]]¶
Preferred dtype of this
Projector
.
- forward(field, indices=None)[source]¶
Compute the forward projection of a tensor field.
- Parameters
field (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array containing coefficients in its fourth dimension, which are to be projected into two dimensions. The first three dimensions should match thevolume_shape
of the sample.indices (
Optional
[ndarray
[Any
,dtype
[int
]]]) – A one-dimensional array containing one or more indices indicating which projections are to be computed. IfNone
, all projections will be computed.
- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]- Returns
An array with four dimensions
(I, J, K, L)
, where the first dimension matchesindices
, such thatprojection[i]
corresponds to the geometry of projectionindices[i]
. The second and third dimension contain the pixels in theJ
andK
dimension respectively, whereas the last dimension is the coefficient dimension, matchingfield[-1]
.
- property is_dirty: bool¶
Returns
True
if the system geometry has changed without the projection geometry having been updated.
- property number_of_projections: int¶
The number of projections as defined by the length of the
Geometry
object attached to this instance.
- class mumott.methods.projectors.SAXSProjectorCUDA(geometry)[source]¶
Projector for transforms of tensor fields from three-dimensional space to projection space. Uses a projection algorithm implemented in
numba.cuda
.- Parameters
geometry (Geometry) – An instance of
Geometry
containing the necessary vectors to compute forwared and adjoint projections.
- adjoint(projections, indices=None)¶
Compute the adjoint of a set of projections according to the system geometry.
- Parameters
projections (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array containing coefficients in its last dimension, from e.g. the residual of measured data and forward projections. The first dimension should matchindices
in size, and the second and third dimensions should match the system projection geometry. The array must be contiguous and row-major.indices (
Optional
[ndarray
[Any
,dtype
[int
]]]) – A one-dimensional array containing one or more indices indicating from which projections the adjoint is to be computed.
- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]- Returns
The adjoint of the provided projections. An array with four dimensions
(X, Y, Z, P)
, where the first three dimensions are spatial and the last dimension runs over coefficients.
- property dtype: Union[dtype[Any], None, type[Any], _SupportsDType[dtype[Any]], str, tuple[Any, int], tuple[Any, Union[SupportsIndex, collections.abc.Sequence[SupportsIndex]]], list[Any], _DTypeDict, tuple[Any, Any]]¶
Preferred dtype of this
Projector
.
- forward(field, indices=None)¶
Compute the forward projection of a tensor field.
- Parameters
field (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array containing coefficients in its fourth dimension, which are to be projected into two dimensions. The first three dimensions should match thevolume_shape
of the sample.indices (
Optional
[ndarray
[Any
,dtype
[int
]]]) – A one-dimensional array containing one or more indices indicating which projections are to be computed. IfNone
, all projections will be computed.
- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]- Returns
An array with four dimensions
(I, J, K, L)
, where the first dimension matchesindices
, such thatprojection[i]
corresponds to the geometry of projectionindices[i]
. The second and third dimension contain the pixels in theJ
andK
dimension respectively, whereas the last dimension is the coefficient dimension, matchingfield[-1]
.
- property is_dirty: bool¶
Returns
True
if the system geometry has changed without the projection geometry having been updated.
- property number_of_projections: int¶
The number of projections as defined by the length of the
Geometry
object attached to this instance.
- class mumott.methods.projectors.SAXSProjectorAstra(*args, **kwargs)[source]¶
Projector for transforms of tensor fields from three-dimensional space to projection space.
This class provides wrappers for interoperability with the python bindings for ASTRA, a tomography toolbox which implements GPU-based computations for tomography.
- Parameters
geometry (Geometry) – An instance of
Geometry
containing the necessary vectors to compute forwared and adjoint projections.
- adjoint(projections, indices=None)[source]¶
Compute the adjoint of a set of projections according to the system geometry.
- Parameters
projections (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array containing coefficients in its last dimension, from e.g. the residual of measured data and forward projections. The first dimension should matchindices
in size, and the second and third dimensions should match the system projection geometry.indices (
Optional
[ndarray
[Any
,dtype
[int
]]]) – A one-dimensional array containing one or more indices indicating from which projections the adjoint is to be computed.
- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]- Returns
The adjoint of the provided projections. An array with four dimensions
(X, Y, Z, P)
, where the first three dimensions are spatial and the last dimension runs over coefficients.
- property dtype: Union[dtype[Any], None, type[Any], _SupportsDType[dtype[Any]], str, tuple[Any, int], tuple[Any, Union[SupportsIndex, collections.abc.Sequence[SupportsIndex]]], list[Any], _DTypeDict, tuple[Any, Any]]¶
The dtype input fields and projections require.
- forward(field, indices=None)[source]¶
Compute the forward projection of a tensor field.
- Parameters
field (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array containing coefficients in its fourth dimension, which are to be projected into two dimensions. The first three dimensions should match thevolume_shape
of the sample.indices (
Optional
[ndarray
[Any
,dtype
[int
]]]) – Must beNone
since the computation of individual projections is not supported for this projector.
- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]- Returns
An array with four dimensions
(I, J, K, L)
, where the first dimension corresponds to projection index. The second and third dimension contain the pixels in theJ
andK
dimension respectively, whereas the last dimension is the coefficient dimension, matchingfield[-1]
.
- property is_dirty: bool¶
Returns
True
if the system geometry has changed without the projection geometry having been updated.
- property number_of_projections: int¶
The number of projections as defined by the length of the
Geometry
object attached to this instance.
- class mumott.methods.projectors.SAXSProjectorNumba(*args, **kwargs)[source]¶
Projector for transforms of tensor fields from three-dimensional space to projection space.
- Parameters
geometry (Geometry) – An instance of
Geometry
containing the necessary vectors to compute forwared and adjoint projections.step_size (np.float64) – A number
0 < step_size <= 1
that indicates the degree of upsampling when computing projections. Default is0.5
.
- adjoint(projections, indices=None)[source]¶
Compute the adjoint of a set of projections according to the system geometry.
- Parameters
projections (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array containing coefficients in its last dimension, from e.g. the residual of measured data and forward projections. The first dimension should matchindices
in size, and the second and third dimensions should match the system projection geometry.indices (
Optional
[ndarray
[Any
,dtype
[int
]]]) – A one-dimensional array containing one or more indices indicating from which projections the adjoint is to be computed.
- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]- Returns
The adjoint of the provided projections. An array with four dimensions
(X, Y, Z, P)
, where the first three dimensions are spatial and the last dimension runs over coefficients.
- create_sampling_kernel(kernel_dimensions=(3, 3), kernel_width=(1.0, 1.0), kernel_type='bessel')[source]¶
Creates a kernel to emulate the point spread function (PSF) of the beam. This can improve the accuracy of the projection function.
- Parameters
kernel_dimensions (
Tuple
[int
,int
]) – A tuple of how many points should be sampled in each direction of the kernel. The total number of line integrals per pixel in the data iskernel_dimensions[0] * kernel_dimensions[1]
.kernel_width (
Tuple
[float
,float
]) – Width parameter for the kernel in units of pixels. Typically the full-width-half-maximum (FWHM) of the beam used for data acquisition.kernel_type (
str
) –The type of kernel to use. Acceptable values are
'bessel'
,'rectangular'
, and'gaussian'
.'bessel'
uses a sinc function multiplied by a Lanczos window with the width parameter being the first Bessel zero. This gives a sharply peaked distribution that goes to zero at twice the full width half maximum, and samples are taken up to this zero.'rectangular'
samples uniformly in a rectangle of sizekernel_width
.'gaussian'
uses a normal distribution with the FWHM given bykernel_width
sampled up to twice the FWHM.
- Return type
- property dtype: Union[dtype[Any], None, type[Any], _SupportsDType[dtype[Any]], str, tuple[Any, int], tuple[Any, Union[SupportsIndex, collections.abc.Sequence[SupportsIndex]]], list[Any], _DTypeDict, tuple[Any, Any]]¶
Preferred dtype of this
Projector
.
- forward(field, indices=None)[source]¶
Compute the forward projection of a tensor field.
- Parameters
field (
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]) – An array containing coefficients in its fourth dimension, which are to be projected into two dimensions. The first three dimensions should match thevolume_shape
of the sample.indices (
Optional
[ndarray
[Any
,dtype
[int
]]]) – A one-dimensional array containing one or more indices indicating which projections are to be computed. IfNone
, all projections will be computed.
- Return type
ndarray
[Any
,dtype
[TypeVar
(_ScalarType_co
, bound=generic
, covariant=True)]]- Returns
An array with four dimensions
(I, J, K, L)
, where the first dimension matchesindices
, such thatprojection[i]
corresponds to the geometry of projectionindices[i]
. The second and third dimension contain the pixels in theJ
andK
dimension respectively, whereas the last dimension is the coefficient dimension, matchingfield[-1]
.
- property is_dirty: bool¶
Returns
True
if the system geometry has changed without the projection geometry having been updated.
- property number_of_projections: int¶
The number of projections as defined by the length of the
Geometry
object attached to this instance.
- property projection_shape: Tuple[int]¶
The shape of each projection defined by the
Geometry
object attached to this instance, as a tuple.
- property sampling_kernel: ndarray[Any, dtype[float64]]¶
Convolution kernel for the sampling profile.
- property step_size: float64¶
One-dimensional upsampling ratio for the integration of each projection line.