Benchmarking of different reconstruction algorithms for industrial cone-beam CT

In the field of industrial computed tomography, there have been a lot of interests in different reconstruction algorithms which could perform better than the standard FDK algorithm in non-ideal conditions leading to faster CT acquisitions. However, the computation efforts required by iterative algorithms are quite high, which makes them less practical in industrial environments. In this work, we propose a benchmark of the performances of various reconstruction algorithms using Tomographic Iterative GPU-based Reconstruction (TIGRE) toolbox; the benchmarking criteria used include both image quality metrics and task-oriented post-reconstruction results (measurements, porosity, etc.). Furthermore, we have also tested their performances with limited projections and limited angular range to facilitate faster CT acquisitions. The preliminary results indicate several improvements in the reconstructions (better feature definitions) with one of the iterative algorithms (CGLS: conjugate gradient least squares) as compared to FDK.


Introduction
Cone-beam computed tomography has become the default approach with lab-based CT systems, especially in industrial environments due to optimal acquisition times.In the cone-beam CT, the most widely used reconstruction algorithm is FDK (Feldkamp, Devis, and Kress), which relies on the filtered back projection (FBP) method [1].FDK being an analytical algorithm requires less intensive computation as compared to the iterative algorithms; however, the quality of the reconstruction is strongly affected by the reduction of the number of projections.Therefore, the iterative reconstruction algorithms have gained interest in the context of fast CT, especially due to their better results with limited data (sparse-view CT and limited-angle tomography) and better robustness to measurement noise [2].The basic principle of the iterative algorithm is that for every angle, a forward projection is calculated by summing the intensities of all pixels along the particular ray-path, the forward projection is then compared to the recorded projection, and the process continues iteratively until convergence to the desired result is achieved [3].This makes it computationally demanding, however, a comparison of their performances against FDK in industrial case scenarios has not been reported in the literature.In this work, we utilized the TIGRE toolbox, which is a GPU-based CT reconstruction software repository that contains a wide variety of algorithms such as the standard FDK, iterative algorithms: Gradient-based algorithms (SART, OS-SART, SIRT), and Krylov subspace algorithms (CGLS).It offers both MATLAB and python implementation.For the comparisons, several image quality matrices are provided and the most common is the L2 norm (of the difference between the forward projection of the current image estimate and the measured projection) and RMSE between iterations [4,5].In this work, two different case studies have been performed considering different objects conducted by CEA and PTB.The goal of this work is to test the performance of the different algorithms especially in the case of limited data (limited angle and/or limited projections), therefore the limited data is obtained by subsampling of the full data; generalized methodology is explained in Figure 1.

First case study
A real CT acquisition usually suffers some misalignments, eventually leading to image artifacts; hence, first it would be interesting to have a simulated dataset, which could also serve as the reference.We have obtained a simulated CT acquisition for the nominal model of hole plate (courtesy: PTB) with 900 projections (over a complete rotation cycle of 360°) with a detector size of 1024 × 1024 (pixel 0.2 mm).In addition, further undersampling is carried out to test the algorithm with reduced projections and/or limited angle trajectory (270°).The simulation was performed in CIVA, which is a multi-technique Non-destructive testing simulation software that integrates, in particular, a radiographic and computed tomographic (RT/CT) module [6].The post-processing and measurements were carried out in VGStudioMax 3.4.

Configuration
The dataset utilized for the current study is a simulated CT scan of the PTB hole plate obtained from CIVA with the following specifications (with noise and without noise):  900 projections  Image size: 1024 × 1024 (pixel 0.2 mm)  Dsd (source-to-detector distance) = 500 mm  Dso (source-to-object distance) = 200 mm Further subsampling with respect to the number of projections and the angular range was done to have multiple datasets and analyses the effect of limited projection and limited angle on the quality of reconstruction.Three algorithms (FDK, CGLS, SIRT) were implemented on the above-mentioned eight datasets and the performance results are provided here for comparison.Figure 2 shows the execution time in each run, which includes the reconstruction and saving of the file.The RMSE and L2 errors are also plotted for each run as shown in Figure 3.Although the time taken by iterative algorithms (CGLS and SIRT 20 iterations each) is considerably higher, it will be interesting to compare the measurement results for concluding the overall performance of algorithms.It is also interesting to note that the CLGS algorithms exhibited the least RMSE and L2 error.The next section reports the comparison of post-reconstruction analyses.

Qualitative comparison
The reconstructed datasets were imported in VGStudioMax for the post analyses (surface determination and measurements).
Figure 4 shows a comparison of the reconstructions (center slice) for the D1 dataset, which is a baseline result of each algorithm to be compared against the subsampled datasets.The SIRT algorithm is worst performing as the contrast is extremely poor which could be problematic to define the holes followed by measurements.
Figure 6: The measurement strategy for the hole plate Figure 7 shows the deviations in the measurements of the center-to-center distances against the nominal value.As a common trend, the results obtained from the SIRT algorithms exhibit maximum deviations; however, there is no clear difference in FDK and CGLS results.The goal of the study was to investigate the performance of the several iterative reconstructions algorithm in comparison to the standard FDK algorithm.Two iterative algorithms CGLS and SIRT were tested.In terms of time, the iterative algorithms were taking considerably higher time, however, the CGLS algorithm was performing better when looking at the RMSE and L2 error results.On further comparison of the post-reconstruction analyses, it can be concluded that the SIRT resulted in a very bad reconstruction with poor contrast, which eventually leads to poor definition of the features (holes).FDK and CGLS produced very comparable results, however, FDK-based results are affected by limited projection and limited angular range approach, unlike CGLS.Therefore, CGLS iterative algorithm can be preferred over FDK when there is a small number of projections and/or partial angle datasets are available, which could then justify its longer execution time.

Measurement
Another investigation was performed using an additive manufacturing standard developed by the National Physical Laboratory (NPL).The object resembles a cuboid with another smaller cuboid and a cylinder carved out from the top and the same feature of a small cuboid and cylinder stacked on top.Therefore, the object features mostly bidirectional measurands such as diameters.
The measurement was performed on a YXLON FF35 CT system using a transmission X-ray tube in microfocus mode and with full circle trajectory, i.e., 360° sample rotation.1440 projections with an integration time of 2 s were recorded.Using a continuous scan mode results in a total scanning time of 48 min.The details about the imaging conditions are summarized in Table 2.This represents a high-quality scan rather than a fast CT scan with low noise and enough projections to eliminate sampling artifacts.To generate data equivalent to a fast CT scan the originally obtained set of projections was reduced by selecting projections at regular intervals, e.g., use only every 10 th projection.This approach realizes the reduced scanning time for fast CT by fewer projections of the same quality and angle range as opposed to other approaches like reduced integration time per projection or limited angle scans.This way four reduced sets of projections were generated with 360, 240, 180, and 144 projections.The respective scan duration for these reduced sets of projections would be 12 min, 8 min, 6 min, and 5 min.

Reconstruction algorithms and parametrization
For reconstruction of the original and reduced sets of projection data the Python implementation of the TIGRE Toolbox was used.Three different iterative algorithms were chosen in addition to the FDK algorithm, namely, OS-SART-TV, CGLS, and FISTA.These algorithms used a block size as parameter which represents the number of projections used per image update.The  black size was always chosen to be one-sixth of the respective reduced set of projections.Of the different parameters used by the iterative algorithms, only the number of iterations and the initialization were varied to optimize the quality of the reconstruction and minimize the reconstruction time.The influence of other parameters like the order strategy or the hyperparameters for total variation was not investigated.For these parameters, the default values are used.While these other parameters might be optimized to improve results, their adjustment requires more detailed knowledge of algorithms.Therefore, this investigation focuses on the variation of those parameters suitable for optimization by the end-user.The convergence behavior of the algorithms is shown in Figure 9 where the quality of the reconstruction is potted against the number of iterations exemplarily for the reduced set of 180 projections.As a quality measure, the L2 norm between the original projections and the forward projection of the reconstruction result is used.Because the value of the L2 norm is proportional to the number of projections it is reported as per projection to allow comparisons between the different reduced sets of projections.Further two different initializations are compared, i.e., no initialization where the iterative algorithm starts with an empty volume and an FDK initialization where the iterative algorithm starts from the result of a regular FDK reconstruction.The plot shows that the CGLS and OS-SART-TV algorithms converge faster reaching a stable L2 error after less than 10 iterations.In contrast, the FISTA algorithm requires about 50 iterations to reach a comparable L2 error.Initialization with an FDK volume generally leads to faster convergence of the L2 error reducing the required iterations for OS-SART-TV and CGLS to about 3 to reach a stable L2 error.For FISTA the FDK initialization also reduces the required iterations significantly, however, it needs still about 10 iterations to reach a comparable L2 error.CGLS produces the lowest L2 error of all investigated algorithms given enough iterations.This leads to the expectation that CGLS yields the best results.However, as shown in Figure 10, CGLS does not suppress streak artifacts caused by the under-sampling of the projection angles as well as OS-SART-TV or FISTA.Also, the background value is not reduced to zero which is the case for OS-SART-TV and FISTA.OS-SART-TV produces the smoothest image of the investigated algorithms.For the subsequent analysis, an optimized number of iterations and initialization was chosen for each algorithm to achieve a stable reconstruction result in the least amount of time.As the FDK initialization reduces the required number of iterations significantly it was chosen for all algorithms and the chosen number of iterations were 3 for CGLS and OS-SART-TV and 10 for FISTA.With these parameters, CGLS takes between 1.3 min and 1.7 min for the reconstruction depending on the number of projections used.This is comparable to the time needed for an FDK reconstruction.OS-SART-TV takes about three times longer with 4.4 min to 4.9 min and FISTA is the slowest algorithm taking between 7.9 min and 10. and the metrology tools to fit cylinders to the surface.As ground truth for the analysis, an FDK reconstruction of the full data set was used, and deviations of radii and form deviation are given with respect to this volume.The form deviation was evaluated as the maximum spread of all surface points used for fitting.Dimensional measurements were only investigated for volumes obtained with optimized parameters and results are analyzed in dependence of total measurement time, i.e., scan time plus reconstruction time.Results are presented in Figure 11.The smallest deviations of the radii are observed with a FISTA reconstruction using 360 projections.However, the required total time of this combination is the largest, and CGLS or FDK may perform similarly with more projections used to match the total time required for scan and reconstruction.The deviations observed with CGLS are almost the same as the ones observed with FDK.This might be due to the initialization with an FDK result, however, without FDK initialization deviations for the radii of more than 400 µm are observed.The deviation of the radii increases more and more rapidly for less than 180 projections used in the reconstruction with FDK and CGLS.The trend of the inner and outer cylinder suggests that the surface is systematically shifted outwards-For the FISTA algorithm, the radii also deteriorate for fewer projections, however, the outer cylinder decreases as well as the inner cylinder which cannot be explained by a systematic shift of the surface.For OS-SART-TV the radii stay approximately constant for the investigated numbers of projections, however, at an offset of about 6 µm.The deviation of the form error shows similar results with an increase of the deviations below 180 projections for FDK, CGLS, and FISTA while OS-SART-TV performs best here most likely due to the suppressed noise.In general, the inner cylinder shows a larger form deviation of the form error which is also due to noise which is increased at inner features.

Conclusions
From these results it is not clear how much benefit the iterative algorithms have over a standard FDK reconstruction in the investigated fast CT scenario.While the reconstruction results are, in general, visually more appealing the results of dimensional measurements are not significantly, if at all, better than those of an FDK reconstruction given the same amount of time.The results from the FDK reconstruction only deteriorate rapidly for a total time of less than 10 min.From this point on OS-SART-TV may prove a better option than FDK, however, this corresponds to even fewer projections (<100) to further reduce scanning time which was not investigated in this study.

2 MethodologyFigure 1 :
Figure 1: The flow diagram of the generalized methodology

Figure 2 :Figure 3 :
Figure 2: The implantation time of each of the algorithm

Figure 4 :
Figure 4: Slice view of the reconstructions for D1 dataset

Figure 5 :
Figure 5: Nominal-actual comparison with respect to FDK and CGLS

Figure 7 :
Figure 7: Deviation in the center-to-center distances measurementsFigure8apresents the deviation (from nominal value) in the diameter measurements.Overall, the FDK and CGLS algorithms resulted in smaller deviations, which is of the order of around 10µm.Whereas SIRT algorithm resulted in a very high deviation ranging from 60 to 100 µm, which could be attributed to the poor contrast and bad definition of the hole as clearly seen in the cross-sections presented in Figure8b.Moreover, the contrast is better with CGLS as compared to FDK, which results in better hole definition (better surface determination).Both the measurement results and the cross-section views indicate very similar results for CGLS algorithm even for limited projection and limited angle datasets.In contrast, FDK results in a slight degradation in the hole definition with data deficiency approaches.

Figure 8 :
Figure 8: Deviations in hole diameter measurements (a) and comparison of surface determination and hole definition (Hole 1: corresponding to the measurement shown in dotted area) (b) 3.3 Conclusions

Figure 9 :
Figure 9: Plot of L2 norm per projection between original projections and forward projection of reconstruction result against number of iterations for the reduced set of 180 projections.

Figure 10 :
Figure 10: Slice views of volumes reconstructed with the respective algorithms using 180 projections and the optimized set of parameters.
Figure11: Plot of the deviations between FDK of full data and respective reconstruction of reduced data sets for cylinder radii (left) and the cylinder's form deviation (right).Form deviation is given as maximum spread of probed surface points.

Table 1 :
Multiple datasets used for the study

Table 2 :
Detailed CT parameters for original scan of the additive manufacturing standard