Introduction

There are various things, and I wanted to reconstruct a 3D CBCT (Cone-Beam Computed Tomography) from a projected image. However, even if I search for "CBCT reconstruction python" in Japanese, the theory of reconstruction is a hit, but there is not much information on how to do it with the code. In particular, we can find some 2D reconstructions, but few 3D reconstructions. So, this time, I reconstructed 3D CBCT with python using a library called TIGRE (https://github.com/CERN/TIGRE), so I will share how to do it.

※Notes

This article only explains how to reconstruct 3D CBCT using TIGRE, and does not create a reconstruction program from scratch. There is no theoretical explanation for CBCT reconstruction. I don't care about the theory or the details of the program, but anyway, if I can reconstruct the 3D CBCT, that's fine! !! It is for people who say. Please note. By the way, for those who want to know more about CBCT, Basics of cone beam CT image reconstruction [^ 1] is recommended.

PC environment

OS : Windows 10 Professional --Language: Python 3.7.3
GPU : Geforce RTX 2080 Ti
Python : 3.7.3
CUDA : 10.1
MSVS : 2017
TIGRE : 2.0

Difference between CT and CBCT reconstruction

To explain it roughly, CT and CBCT differ greatly in whether the X-rays entering the detector are fan beams or cone beams. Therefore, the reconstruction method also changes. In CT, X-rays are incident almost perpendicular to the detector (= fan beam), so it can be reconstructed by a method called Filtered Back Projection (FBP method). On the other hand, in the case of CBCT, X-rays enter the detector at an angle (cone angle) (= cone beam), so instead of the FBP method, the FDK method [^ 2] proposed by Feldkamp et al. Is used again. Must be configured. For a detailed explanation, see the basics of cone-beam CT image reconstruction [^ 1], but for the time being, just remember that "when reconstructing CBCT, use the FDK method instead of the FBP method."

About TIGRE

This time, we will reconstruct 3D CBCT using a library developed by CERN (European Organization for Nuclear Research) called TIGRE (https://github.com/CERN/TIGRE). It's amazing that CERN can do anything. If you search for Python's CBCT reconstruction library, you will find several others, but I personally think that TIGRE is the easiest to install and is also GPU compatible, so it is easy to use.

Install TIGRE

Follow the steps in Simple Instructions in Installation Instructions for Python on Github. First install MSVC, CUDA, python, then do the following:

git clone https://github.com/CERN/TIGRE.git
#TIGRE obtained with git/Python/After moving to
python setup.py install --user

When setup.py finishes without error, try running example.py in the TIGRE/Python folder. If you see the following image, it is installed correctly.

Execution result of example.py

Installation error

By the way, in my environment, when I ran setup.py, I got the following error.

nvcc fatal : Unsupported gpu architecture 'compute_86'

Is this because my RTX 2080 Ti doesn't support the compute_86 architecture? It seems that I solved it by commenting out a part of setup.py.

`setup.py`


#Around the 18th line
COMPUTE_CAPABILITY_ARGS = [  # '-gencode=arch=compute_20,code=sm_20', #deprecated
    #'-gencode=arch=compute_30,code=sm_30',#deprecated
    '-gencode=arch=compute_37,code=sm_37',
    '-gencode=arch=compute_52,code=sm_52',
    '-gencode=arch=compute_60,code=sm_60',
    '-gencode=arch=compute_61,code=sm_61',
    '-gencode=arch=compute_70,code=sm_70',
    '-gencode=arch=compute_75,code=sm_75',
    #'-gencode=arch=compute_86,code=sm_86',#← This
    '--ptxas-options=-v', '-c',
    '--default-stream=per-thread',
    ]

CBCT projection and reconstruction

Once you have TIGRE installed, it's easy. This is because example.py has already acquired the projected image and reconstructed it by the FDK method. However, in example.py, Geometry of CBCT is the default state, so we will allow you to change Geometry.

`FDK.py`


from __future__ import print_function
from __future__ import division
import numpy as np
from matplotlib import pyplot as plt
import tigre
import tigre.algorithms as algs
from tigre.demos.Test_data import data_loader
from tigre.utilities.Measure_Quality import Measure_Quality

#Geometry settings
geo = tigre.geometry(mode='cone',  default=False)
geo.DSD = 1500                                     # Distance Source Detector      (mm)
geo.DSO = 1000                                     # Distance Source Origin        (mm)
# Detector parameters
geo.nDetector = np.array((300, 300))               # number of pixels              (px)
geo.dDetector = np.array((0.8, 0.8))               # size of each pixel            (mm)
geo.sDetector = geo.nDetector * geo.dDetector    # total size of the detector    (mm)
# Image parameters
geo.nVoxel = np.array((256, 256, 256))             # number of voxels              (vx)
geo.sVoxel = np.array((128, 128, 128))             # total size of the image       (mm)
geo.dVoxel = geo.sVoxel/geo.nVoxel               # size of each voxel            (mm)
# Offsets
geo.offOrigin = np.array((0, 0, 0))                # Offset of image from origin   (mm)
geo.offDetector = np.array((0, 0))                 # Offset of Detector            (mm)

# Auxiliary
geo.accuracy = 0.5                                 # Accuracy of FWD proj          (vx/sample)
# Mode
geo.mode = 'cone'                                  # parallel, cone                ...
geo.filter = 'hann'                       #  None, shepp_logan, cosine, hamming, hann

#Shooting angle, 0~Shooting at 360 degrees
nangles = 360
angles = np.linspace(0, 2 * np.pi, nangles, endpoint=False, dtype=np.float32)

#Get the head phantom provided by default
head = data_loader.load_head_phantom(geo.nVoxel)

#Get the projected image of the head phantom
proj = tigre.Ax(head,geo,angles)

#Reconstruction with FDK
fdkout = algs.fdk(proj,geo,angles)

#drawing
plt.subplot(1,3,1)
plt.imshow(fdkout[geo.nVoxel[0]//2])
plt.subplot(1,3,2)
plt.imshow(fdkout[:, geo.nVoxel[1]//2, :])
plt.subplot(1,3,3)
plt.imshow(fdkout[:, :, geo.nVoxel[2]//2])
plt.show()

#Drawing method provided by tigre
tigre.plotProj(proj)
tigre.plotImg(fdkout)

Execution result

This is a reconstructed image by the FDK method.

You can also check the projected image for each angle in the tigre.plotProj (proj) part.

Change input image

The input CT image is in numpy.ndarray format, so if you want to calculate with another CT image, just change the input image to numpy.ndarray format and then replace the following part.

head = data_loader.load_head_phantom(geo.nVoxel)

Also, if you want to reconstruct using the projected image prepared by yourself, just change the projected image to numpy.ndarray format and then replace the following part.

proj = tigre.Ax(head,geo,angles)

Summary

I explained how to use a library called TIGRE to reconstruct 3D CBCT with python. I hope it helps people who want to do CBCT reconstruction but don't understand.

Please point out any mistakes.

References

[^ 1]: Basics of cone beam reconstruction (http://www.iryokagaku.co.jp/frame/03-honwosagasu/433/433.html)

Reconstruction of cone beam CT (CBCT) using python and TIGRE