Cone Beam Model

Cone Beam Model#

The ConeBeamModel class implements a geometry and reconstruction model for cone beam computed tomography. This class inherits all behaviors and attributes of the Tomography Model. It also implements some cone-beam specific functions such as FDK (Feldkamp-Davis-Kress) reconstruction.

For cone-beam geometry, the default detector channel spacing is delta_det_channel is 1 ALU, and the voxels are 3D cubes with spacing delta_voxel.

The default voxel spacing is set to delta_voxel = delta_det_channel / magnification where magnification = source_detector_dist / source_iso_dist. This implies that as the magnification increases, the default voxel spacing decreases. However, these parameters can be changed by the user using the TomographyModel.set_params() method.

See the API docs for the TomographyModel class for details on a wide range of functions that can be implemented using the ConeBeamModel.

Constructor#

class mbirjax.ConeBeamModel(sinogram_shape, angles, source_detector_dist, source_iso_dist, helical_z_shifts=None, use_curved_detector=False)[source]#

Bases: TomographyModel

A class designed for handling forward and backward projections in a cone beam geometry, extending the Tomography Model. This class offers specialized methods and parameters tailored for cone beam setups.

This class inherits all methods and properties from the Tomography Model and may override some to suit parallel beam geometrical requirements. See the documentation of the parent class for standard methods like setting parameters and performing projections and reconstructions.

Parameters not included in the constructor can be set using the set_params method of Tomography Model. Refer to Tomography Model documentation for a detailed list of possible parameters.

Parameters:

  • sinogram_shape (tuple) – Shape of the sinogram as a tuple in the form (views, rows, channels), where ‘views’ is the number of different projection angles, ‘rows’ correspond to the number of detector rows, and ‘channels’ index columns of the detector that are assumed to be aligned with the rotation axis.

  • angles (ndarray or jax array) – A 1D array of projection angles, in radians, specifying the angle of each projection relative to the origin.

  • source_detector_dist (float) – Distance between the X-ray source and the detector in units of ALU.

  • source_iso_dist (float) – Distance between the X-ray source and the center of rotation in units of ALU.

  • helical_z_shifts (ndarray or jax array, optional) – Per-view axial shifts (ALU; same length as angles) for helical mode.

  • use_curved_detector (bool, optional) – Detector geometry type. Either False (default) or True implies each detector row has constant distance from source.

Note

Additional parameter:

recon_slice_offset (float, default=0) - This parameter controls the vertical offset of the reconstruction in ALU. If recon_slice_offset is positive, the region below iso is reconstructed.

See also

TomographyModel

The base class from which this class inherits.

Alternative Reconstruction#

ConeBeamModel.fdk_recon(sinogram, filter_name='ramp', view_batch_size=100)[source]#

Perform FDK reconstruction on the given sinogram.

Our implementation uses standard filtering of the sinogram, then uses the adjoint of the forward projector to perform the backprojection. This is different from many implementations, in which the backprojection is not exactly the adjoint of the forward projection. For a detailed theoretical derivation of this implementation, see the zip file linked at this page: https://mbirjax.readthedocs.io/en/latest/theory.html

Parameters:

  • sinogram (jax array) – The input sinogram with shape (num_views, num_rows, num_channels).

  • filter_name (string, optional) – Name of the filter to be used. Defaults to “ramp”

  • view_batch_size (int, optional) – Size of view batches (used to limit memory use)

Returns:

recon (jax array) – The reconstructed volume after FDK reconstruction.

ConeBeamModel.split_sino_recon(sino, weights=None, half_overlap=5, init_recon=None, max_iterations=15, stop_threshold_change_pct=0.2, first_iteration=0, compute_prior_loss=False, logfile_path='./logs/recon.log', print_logs=True)[source]#

This function reduces memory usage for cone beam MBIR reconstruction by approximately a factor of 2 by splitting the detector rows into two overlapping halves, reconstructing each half separately, and stitching the reconstructions together.

The function can be called with the same arguments as TomographyModel.recon(), and it should return a reconstruction which is approximately equal to the reconstruction returned by TomographyModel.recon().

Parameters:

  • sino (jnp.ndarray | np.ndarray) – Full sinogram of shape (num_views, num_rows, num_cols).

  • weights (jnp.ndarray | np.ndarray, optional) – Optional sinogram weights with the same shape as sino.

  • half_overlap (int) – Number of overlapping detector rows and recon slices per half. (total overlap = 2 * half_overlap)

  • init_recon (optional) – Same as in the recon method.

  • max_iterations (int, optional) – Same as in the recon method.

  • stop_threshold_change_pct (float, optional) – Same as in the recon method.

  • first_iteration (int, optional) – Same as in the TomographyModel.recon() method.

  • compute_prior_loss (bool, optional) – Same as in the TomographyModel.recon() method.

  • logfile_path (str, optional) – Same as in the TomographyModel.recon() method.

  • print_logs (bool, optional) – Same as in the TomographyModel.recon() method.

Returns:

Tuple[jnp.ndarray, dict] – - Reconstructed volume (jax array). - Dictionary of metadata containing recon and model parameters for each half.

Raises:

  • ValueError – If inputs are missing or shapes are inconsistent.

  • AssertionError – If array dimensions are invalid.

  • TypeError – If half_overlap is not an integer.

Example

>>> import jax.numpy as jnp
>>> import mbirjax as mj
>>> sino = jnp.ones((180, 64, 64))  # (views, rows, cols)
>>> model = mj.ConeBeamModel(sinogram_shape=sino.shape,
...                          angles=jnp.linspace(0, jnp.pi, 180),
...                          source_detector_dist=1000.0,
...                          source_iso_dist=500.0)
>>> recon, recon_info = model.split_sino_recon(sino, half_overlap=4)
>>> recon.shape  # Quilted reconstruction volume
(64, 64, 64)