RTK/Scripts/3DCG
This script uses the file SheppLogan.txt as input.
# Create a simulated geometry
rtksimulatedgeometry -n 180 -o geometry.xml
# You may add "--arc 200" to make the scan short or "--proj_iso_x 200" to offset the detector
# Create projections of the phantom file
rtkprojectgeometricphantom -g geometry.xml -o projections.mha --spacing 2 --dimension 256 --phantomfile SheppLogan.txt
# Reconstruct
rtkconjugategradient -p . -r projections.mha -o 3dcg.mha -g geometry.xml --spacing 2 --dimension 256 -n 20
# Create a reference volume for comparison
rtkdrawgeometricphantom --spacing 2 --dimension 256 --phantomfile SheppLogan.txt -o ref.mha
In the presence of noise, all projection data may not be equally reliable. The conjugate gradient algorithm can be modified to take this into account, and each pixel of the projections can be associated with a weight. The higher the weight, the more reliable the pixel data. Download noisy projections and the associated weights, as well as the geometry, and run the following to compare the regular least squares reconstruction (without weights) and the weighted least squares reconstruction.
# Perform least squares reconstruction
rtkconjugategradient -p . -r noisyLineIntegrals.mha -o LeastSquares.mha -g geom.xml -n 20
# Perform weighted least squares reconstruction
rtkconjugategradient -p . -r noisyLineIntegrals.mha -o WeightedLeastSquares.mha -g geom.xml -w weightsmap.mha -n 20
Taking the weights into account slows down convergence. This can be corrected by using a preconditioner in the conjugate gradient algorithm. The preconditioner is computed automatically from the weights map, you just need to activate the flag :
# Perform preconditioned conjugate gradient reconstruction with weighted least squares cost function
rtkconjugategradient -p . -r noisyLineIntegrals.mha -o WeightedLeastSquares.mha -g geom.xml -w weightsmap.mha -n 20 --preconditioned