Computation of a Consistent System Matrix for Cone-beam Computed Tomography

TL;DR

This work addresses inconsistencies in cone-beam CT forward models arising from line-based discretizations by deriving explicit analytic formulas for a consistent area-/volume-based system matrix $\mathbf{W}$ in both 2D and 3D. The method decomposes cone–voxel intersections into basic subvolumes (triangles and trapezoids in 2D; tetrahedra, prismatoids, and their combinations in 3D) and computes exact weighting factors, enabling $\mathbf{p}=\mathbf{W}\mathbf{u}+\bm{\varepsilon}$ with a regularized objective $F(\mathbf{u})$ to optimize image reconstructions. Numerical experiments on synthetic and real data in 2D and 3D show that the area-/volume-based discretization improves reconstruction quality over traditional line-based approaches, particularly in underdetermined or low-regularization regimes, and a CUDA implementation is provided for practical use. The approach offers a path toward more quantitative CT and potentially better performance in low-dose or dynamic imaging, at the expense of greater computational burden that the authors plan to address in future work.

Abstract

We propose a method for the computation of a consistent system matrix for two- and three-dimensional cone-beam computed tomography (CT). The method relies on the decomposition of the cone-voxel intersection volumes into subvolumes that contribute to distinct detector elements and whose contributions to the system matrix admit exact formulae that can be evaluated without the invocation of costly iterative subroutines. We demonstrate that the reconstructions obtained when using the proposed system matrix are superior to those obtained when using common line-based integration approaches with numerical experiments on synthetic and real CT data. Moreover, we provide a CUDA implementation of the proposed method.

Figures (10)

• Figure 1: Discretization of the measurement process for the fan beam geometry.
• Figure 2: Schematic subdivision of a pixel by rays to the detector edges and/or pixel vertices. We denote $\beta_{\mathbf{m}}$ the angles to the rays to detector edges and $\gamma_{\mathbf{n},i,\varphi}$ the angles of rays to vertices of the pixel.
• Figure 3: Trapezoid shape sub-voxel delimited by two rays to subsequent detector elements.
• Figure 4: Various cases of contribution of a sub-voxel (delimited by green lines) to a detector element. Rays to corners of the detector element are shown in blue. (a) Sub-voxel VPQW with the coefficient vector $\mathbf{C_{m}} = (1, 0, 0, 0)$. (b) Sub-voxel $V_{1}V_{2}W_{2}Q_{2}Q_{1}W_{1}$ with the coefficient vector $\mathbf{C_{m}} = (1, 1, 0, 0)$. The volume is the difference between the volumes $V_{1}PQ_{1}W_{1}$ and $V_{2}PW_{2}Q_{2}$. (c) Sub-voxel $P_{1}Q_{1}Q_{0}P_{0}$$V_{1}W_{1}V_{0}W_{0}$ with the coefficient vector $\mathbf{C_{m}} = (1, 1, 1, 1)$. (d) Sub-voxel $PQ_{1}Q_{2}V_{12}W_{1}W_{2}$ with the coefficient vector $\mathbf{C_{m}} = (1, 1, 1)$.
• Figure 5: Determining the volume of sub-voxels contributing to a detector element by proper segmentation. The volume of the blue delimited sub-voxels with the coefficient vectors $\mathbf{C_{m+1}}=(1,1,0,1)$ and $\mathbf{C_{m}}=(0,0,1,0)$ can be calculated directly. The volume of the green delimited sub-voxels are calculated by subtraction. The plane in magenta color splits the voxel in an upper and lower part.