Do Not Ignore the ASTRA Toolbox!

ASTRA logo.

If you are into computed tomography (CT) from the perspective of algorithm development, or if you want to do the reconstruction yourself instead of using a standard software package (e.g., the one that was included with your scanner), you cannot ignore the ASTRA Tomography Toolbox. Well, I am confident you can, if you put your mind to it, but you shouldn’t. This toolbox is developed by the ASTRA research group (“All Scale Tomographic Reconstruction Antwerp”) of the Vision Lab of the University of Antwerp, together with the Computational Imaging group of the Centrum Wiskunde & Informatica (CWI). Full disclosure: I work at the Vision Lab as a postdoc. However, this article is not sponsored or endorsed by the Vision Lab, it is purely my own opinion.

Quoting from the official description:

It supports 2D parallel and fan beam geometries, and 3D parallel and cone beam. All of them have highly flexible source/detector positioning. A large number of 2D and 3D algorithms are available, including FBP, SIRT, SART, and CGLS. The basic forward and backward projection operations are GPU-accelerated, and directly callable from MATLAB to enable building new algorithms.

During my own PhD research, I have used the ASTRA Toolbox for algorithm development. Because the Toolbox provides GPU implementations of basic components such as projectors, I could use GPU processing already during the experimental stage, allowing me to use real-world (i.e., large) datasets during testing. This is not obvious, since it is often only when an algorithm is complete that a GPU version is written, restricting initial testing to the aptly named academic examples.

MATLAB Wrapper

First of all, to put your mind at rest if you are not a hardcore programmer, the ASTRA Toolbox comes with a user-friendly MATLAB wrapper. All functionality of the Toolbox is available from MATLAB. Since MATLAB is the main way to use the Toolbox, I’ll leave you in the good hands of the official documentation for that, and add an example of my own below that uses the contributed Python wrapper.

Python Wrapper

If you are more of a Python person, there is also a Python wrapper. As an example, I’ve written a Python script that recreates the example with the two white squares on the Wikipedia page on the Radon transform (and additionally also reconstructs the sinogram). The example shows that you can do a lot with only a few lines of code. To keep things brief, there are not a lot of comments, but there are more detailed examples for each part of this script in the samples directory of the Toolbox.

import numpy as np
from scipy import misc
 
import astra
 
# Create phantom.
phantom = np.zeros((128, 128))
phantom[32 : 64, 32 : 64] = np.ones((32, 32))
phantom[64 : 96, 64 : 96] = np.ones((32, 32))
misc.imsave('phantom.png', phantom)
 
# Create geometries and projector.
vol_geom = astra.create_vol_geom(128, 128)
angles = np.linspace(0, np.pi, 180, endpoint=False)
proj_geom = astra.create_proj_geom('parallel', 1., 128, angles)
projector_id = astra.create_projector('linear', proj_geom, vol_geom)
 
# Create sinogram.
sinogram_id, sinogram = astra.create_sino(phantom, projector_id)
misc.imsave('sinogram.png', sinogram)
 
# Create reconstruction.
reconstruction_id = astra.data2d.create('-vol', vol_geom)
cfg = astra.astra_dict('SIRT')
cfg['ReconstructionDataId'] = reconstruction_id
cfg['ProjectionDataId'] = sinogram_id
cfg['ProjectorId'] = projector_id
cfg['option'] = {}
cfg['option']['MinConstraint'] = 0.  # Force solution to be nonnegative.
algorithm_id = astra.algorithm.create(cfg)
astra.algorithm.run(algorithm_id, 100)  # 100 iterations.
reconstruction = astra.data2d.get(reconstruction_id)
misc.imsave('reconstruction.png', reconstruction)
 
# Cleanup.
astra.algorithm.delete(algorithm_id)
astra.data2d.delete(reconstruction_id)
astra.data2d.delete(sinogram_id)
astra.projector.delete(projector_id)

The output of this script is shown below. On the left is the phantom with the two white squares. In the middle is the sinogram, which contains a (2D) projection in each horizontal line. From this sinogram, the reconstruction on the right was created through the algebraic algorithm SIRT.

Phantom (left), sinogram (middle), and reconstruction (right)

For more information on tomography itself, have a look at my series of articles on that subject. There is also a tutorial on creating a 3D reconstruction from 2D cone-beam projections that you might find interesting if you are working with your own datasets.

[update] This post was selected for Writing Clinic at DailyBlogTips! As a result, I was able to make several improvements to the text. Thanks again, Ali!

[update 2] Updated the homepage of the Toolbox; added CWI as a developer; reformulated the text to reflect that the Python wrapper is now an integral part of the Toolbox, instead of a contributed module.