Image gradient L<sub>0</sub>-norm based PICCS for swinging multi-source CT reconstruction

Haijun Yu; Weiwen Wu; Peijun Chen; Changcheng Gong; Junru Jiang; Shaoyu Wang; Fenglin Liu; Hengyong Yu

doi:10.1364/OE.27.005264

1. Introduction

Dynamic computed tomography (CT) has been used to image a changeable object with time in industrial applications [1–4]. To improve the temporal resolution with high image quality for the dynamic CT system, one key problem is how to obtain sufficient projections in short time. However, limited by the constraints of the mechanical gantries, the scanning speed of the conventional CT can reach no more than one rotation per 200ms [5]. To further improve the scanning, the multi-source CT system is developed, in which multiple source-detector pairs installed on a gantry simultaneously rotate around the imaging object to acquire projections. For example, a novel technique called real time tomography (RTT) system [6] was developed by adopting a circular array of X-ray source rather than the complicated mechanical scanning. The C-arm CT was developed for dynamic imaging to improve the temporal resolution by reducing the scan angle range. Those systems have achieved good performance in medical applications. However, the complicated structure and expensive cost make them difficult to be applied in industrial applications.

Recently, a swinging multi-source CT (SMCT) scheme was proposed for aperiodic dynamic imaging. This system consists of $Q$ ( $Q$ is an odd integer) X-ray source/detector pairs, and they are distributed uniformly along a circumference. Different from other multi-source CT systems in medical applications [1,7–9], the SMCT system adopts the bearing technique rather than the complicated slip ring. The bearing has advantages in terms of concise and precise structure, easy manufacture and low price, and it is suitable for large-scale serial production. In addition, compared to the medical multi-source CT system equipped a slip ring for transferring current and signals, the SMCT system loads high transmission current and complicated transmission signals by cable. In the SMCT system, each source-detector pair only swings a small angle so that projections can be gathered at one time. Since there are some local dynamic regions within the objects, it is difficult for the SMCT system to collect sufficient projections within a small time period. This may cause data inconsistency and result in motion artifacts. In the case, the conventional filtered back-projection (FBP) and iterative reconstructed (ART) methods are not powerful enough to obtain high quality images.

Guided by the compressed sensing theory, high quality images can be reconstructed from undersampling data sets. For example, the total variation [10] (TV) minimization used as a common prior knowledge was introduced into CT reconstruction. Indeed, the TV minimization based reconstruction methods can improve reconstructed image quality in terms of noise suppression and artifacts reduction. Particularly, incorporating the prior image constraint into the TV-based reconstruction model, the prior image constrained compressed sensing (PICCS) was proposed for 4DCT reconstruction [11]. Then, the PICCS was also applied to dual-source CT [12] and C-arm CT systems [13].

In aperiodic dynamic imaging, only small regions of the object change with time and other part keeps unchanged, and much information within the projections is redundant. Therefore, the PICCS method was modified for the SMCT system, and it generated swinging multi-source PICCS (SM-PICCS) in our previous work [14]. The goal of SM-PICCS method mainly focuses on minimizing image gradient L₁-norm. However, the image gradient L₁-norm based reconstruction methods may lead to staircase effect with the loss of image details and features. To address the aforementioned disadvantages, lots of non-convex penalties were proposed to constrain the final reconstruction solution. For example, a more general non-convex PICCS (NCPICCS) method was proposed in [15]. Compared with the state-of-art PICCS, the NCPICCS has a better ability to preserve image edge structures and suppress noise. Indeed, the Lp-norm (p<1) constraint based reconstruction methods were proposed to obtain a better solution than the image gradient L₁-norm based methods [16–19], especially in the case of fewer projections. A homotopic L₀-minimization was proposed, and it achieved huge success in magnetic resonance imaging (MRI) reconstruction in [17]. Recently, the image gradient L₀-norm has been widely used to exploit gradient sparsity in CT image reconstruction, such as limited-angle, low-dose spectral CT, etc [20–26]. Compared with image gradient L₁-norm, the image gradient L₀-norm counts the number of nonzero pixels, which can directly measure image gradient sparsity [18]. Since image gradient L₀-norm minimizes the number of pixels that have non-zero gradient magnitudes rather than the image gradient magnitude, it can effectively preserve the image structures.

Encouraged by the advantages of image gradient L₀-norm for image reconstruction, we incorporate it into the PICCS model and generate a new reconstruction model for SMCT system, named L₀-PICCS in this study. Compared with the previous SM-PICCS method, the advantages of the proposed L₀-PICCS are three-fold. (i) It can avoid staircase effect with noise suppression in the edge areas. (ii) It can provide clearer image structure and finer features. (iii) It can reconstruct high-quality images from fewer projections. The contributions of this work are mainly in the following two aspects: (i) Image gradient L₀ norm was incorporated into the PICCS to further improve the ability for reconstructing high-quality images from undersampling projections. (ii) Because the minimization of image gradient L₀-norm is a non-convex and NP-hard problem, it is hard to be solved. Thus, the split-Bregman method [27] is employed to optimize the L₀-PICCS model.

The rest of this paper is organized as follow. In Section II, the SMCT industrial system is briefly introduced and the L₀-PICCS model is proposed and optimized. In Section III, both numerical simulations and realistic data set experiments are performed to validate and evaluate the advantages of the proposed L₀-PICCS method. In section IV, we discuss some related issues and make a conclusion.

2. Method

2.1 SMCT system

To image the industrial aperiodic process, a swinging multi-source CT (SMCT) system was proposed to improve the temporal resolution [14]. Figures 1 and 2 show a typical SMCT system and the corresponding fan-beam geometry. The system consists of $Q$ ( $Q$ is an odd integer) pairs of X-ray source/detector and they are uniformly distributed on a rotating table with a concentric bearing structure. The imaging object is placed on the origin of the rotating table. The SMCT system is driven by a servo motor to swing a certain angle. For such a SMCT system, it can obtain $Q$ uniform circular projections simultaneously. Assuming the object is not changeable with time, the system can obtain full projections when each X-ray source/detector pair rotates an angle of $2 π / Q$ in one direction. However, since the object is always changeable with time, it is difficult to obtain the uniform full projections, and it results in poor quality of reconstructed images. To improve the reconstructed image quality, the SM-PICCS reconstruction method was proposed.

Fig. 1 A representative SMCT system ( $Q = 5$ ).

Download Full Size | PDF

Fig. 2 The fan-beam geometry of the SMCT system in Fig. 1.

Download Full Size | PDF

The CT imaging model can be expressed as a linear system. For noise-free projections, it can be approximated as:

A f = P,

where

A = (a_{m n})

is the system matrix with size of

M \times N

,

M

represents the total number of projections, and

N

denotes the number of image pixels.

P = {[p_{1}, p_{2}, ..., p_{M}]}^{T}

represents the measured projections vector and

f = {[f_{1}, f_{2}, \dots, f_{N}]}^{T}

is the reconstructed image, which can be obtained by solving Eq. (1). Because

M

is usually smaller than

N

in our study, Eq. (1) is ill-posed and difficult to be solved.

The regularization method has been developed to deal with the ill-posed inverse problem of CT reconstruction, which can be formulated as the following minimization problem:

\underset{f}{\arg \min} | | A f - P | |_{2}^{2} + γ R (f),

where

| | A f - P | |_{2}^{2}

denotes the data fidelity and

R (f)

represents the regularization term. The regularization term can be formulated by incorporating prior information which plays an important role in obtaining high quality image, especially for undersampling data set, such as limited-angle, sparse-view, etc.

γ

is the regularization parameter to balance the data fidelity term and the regularization term.

The SM-PICCS was proposed for SMCT in our previous work and it can be expressed as following.

\underset{f}{a r g \min} κ T V (f) + (1 - κ) T V (f - f^{p}), s . t . A f = P,

where

κ

is an empirical parameter and

f^{p}

is the prior image. TV represents the L₁-norm of the image gradient, which can be denoted as

T V (f) = \sum_{i = 2}^{I} \sum_{j = 2}^{J} \sqrt{{‖ f_{i, j} - f_{i - 1, j} ‖}^{2} + {‖ f_{i, j} - f_{i, j - 1} ‖}^{2}},

where I and J are the height and width of reconstructed image and their product equals to

N

.

f_{i, j}

represents the pixel value located at the image grid (i, j). Here, the image gradient values of image boundary are set to 0.

2.2 L₀-PICCS model

The image gradient L₀-norm was proposed to calculate the number of non-zeros in image gradient image. It can be defined as

|| \nabla f | |_{0} = \sum_{i = 2}^{I} \sum_{j = 2}^{J} φ (| f_{i, j} - f_{i - 1, j} | + | f_{i, j} - f_{i, j - 1} |),

where

φ (\cdot)

is a function and it can be read as

φ (| f_{i, j} - f_{i - 1, j} | + | f_{i, j} - f_{i, j - 1} |) = {\begin{matrix} 1 & | f_{i, j} - f_{i - 1, j} | + | f_{i, j} - f_{i, j - 1} | \neq 0 \\ 0 & otherwise \end{matrix} .

Compared with the TV, the image gradient L₀-norm counts the number of non-zeros for gradient image. Consequently, the image edge information and sharp structure can be effectively retained. Thus, the L₀-norm of image gradient is employed to replace the L₁-norm in the SM-PICCS algorithm. The reconstruction model of L₀-PICCS can be formulated as following

\begin{matrix} \begin{matrix} \underset{f}{arg \min} & κ | | \nabla f | |_{0} \end{matrix} + (1 - κ) | | \nabla (f - f^{p}) | |_{0} & s \end{matrix} . t . A f = P .

Equation (7) is an equality constraint problem, which can be converted into an unconstrained optimization problem.

\begin{matrix} \underset{f}{\arg \min} & \frac{1}{2} | | A f - P | |_{2}^{2} + λ (κ | | \nabla f | |_{0} + (1 - κ) | | \nabla (f - f^{p}) | |_{0}) \end{matrix} .

To solve Eq. (8), two auxiliary variables

u_{1}

and

u_{2}

are introduced to replace

f

and

f - f^{p}

respectively. Then, the Eq. (8) can be rewritten as

\begin{matrix} \begin{matrix} \underset{f}{\arg min} & \frac{1}{2} \end{matrix} || A f - P | |_{2}^{2} + λ (κ | | \nabla u_{1} | |_{0} + (1 - κ) | | \nabla u_{2} | |_{0}) & s . t . u_{1} \end{matrix} = f, u_{2} = f - f^{p} .

Equation (9) can be further evolved into an unconstrained minimization problem

\begin{array}{l} \begin{matrix} \underset{f, u_{1}, u_{2}, t_{1} . t_{2}}{\arg min} & {\frac{1}{2} || A f - P | |_{2}^{2} + λ (κ | | \nabla u_{1} | |_{0} + (1 - κ) | | \nabla u_{2} | |_{0}) \end{matrix} \\ + \frac{δ_{1}}{2} | | f - u_{1} - t_{1} | |_{2}^{2} + \frac{δ_{2}}{2} | | f - f^{p} - u_{2} - t_{2} | |_{2}^{2}}, \end{array}

where

t_{1}

and

t_{2}

represent two error feedback variables respectively. The problem Eq. (10) can be divided into three sub-problems by employing the split-Bregman strategy [29]:

Sub-problem 1:

f^{(k + 1)} = \underset{f}{\arg \min} \frac{1}{2} | | A f - P | |_{2}^{2} + \frac{δ_{1}}{2} | | f - u_{1}^{(k)} - t_{1}^{(k)} | |_{2}^{2} + \frac{δ_{2}}{2} | | f - f^{p} - u_{2}^{(k)} - t_{2}^{(k)} | |_{2}^{2} .

Sub-problem 2:

u_{1}^{(k + 1)} = \underset{u_{1}}{\arg \min} \frac{δ_{1}}{2} | | f^{(k + 1)} - u_{1} - t_{1}^{(k)} | |_{2}^{2} + λ κ | | \nabla u_{1} | |_{0},

u_{2}^{(k + 1)} = \underset{u_{2}}{\arg \min} \frac{δ_{2}}{2} | | f^{(k + 1)} - f^{p} - u_{2} - t_{2}^{(k)} | |_{2}^{2} + λ (1 - κ) | | \nabla u_{2} | |_{0} .

Sub-problem 3:

t_{1}^{(k + 1)} = t_{1}^{(k)} - (f^{(k + 1)} - u_{1}^{(k + 1)}),

t_{2}^{(k + 1)} = t_{2}^{(k)} - (f^{(k + 1)} - f^{p} - u_{2}^{(k + 1)}) .

Because the mathematical model of Eq. (11a) is a strictly convex problem, the sub-problem 1 can be easily solved. To determine the minimization point of Eq. (11a), the derivative of Eq. (11a) should be equivalent to zero. We have,

A^{T} (A f - P) + δ_{1} (f - u_{1}^{(k)} - t_{1}^{(k)}) + δ_{2} (f - f^{p} - u_{2}^{(k)} - t_{2}^{(k)}) = 0,

which can be simplified as

(A^{T} A + δ_{1} I + δ_{2} I) f = A^{T} P + δ_{1} (u_{1}^{(k)} + t_{1}^{(k)}) + δ_{2} (f^{p} + u_{2}^{(k)} + t_{2}^{(k)}) .

Equation (13) is equivalent to

\begin{array}{l} (A^{T} A + δ_{1} I + δ_{2} I) f = (A^{T} A + δ_{1} I + δ_{2} I) f^{(k)} - (A^{T} A + δ_{1} I + δ_{2} I) f^{(k)} \\ + A^{T} P + δ_{1} (u_{1}^{(k)} + t_{1}^{(k)}) + δ_{2} (f^{p} + u_{2}^{(k)} + t_{2}^{(k)}), \end{array}

Therefore, the

f^{(k + 1)}

can be updated by

\begin{array}{l} f^{(k +1)} = f^{(k)} - {(A^{T} A + δ_{1} I + δ_{2} I)}^{- 1} \\ \times [A^{T} (A f^{(k)} - P) + δ_{1} (f^{(k)} - u_{_{1}}^{(k)} - t_{_{1}}^{(k)}) + δ_{2} (f^{(k)} - u_{_{2}}^{(k)} - t_{_{2}}^{(k)} - f^{p})] . \end{array}

The sub-problem 2 includes an L₀-norm minimization and it results in non-convex and NP-hard problem. Fortunately, this problem can be solved by utilizing an approximate method [28]. Here, since Eq. (11c) is similar to Eq. (11b), only Eq. (11b) is discussed. Equation (11b) is equivalent to the following problem:

u_{1}^{(k + 1)} = \underset{u_{1}}{\arg \min} | | f^{(k + 1)} - u_{1} - t_{1}^{(k)} | |_{2}^{2} + \frac{2 λ κ}{δ_{1}} | | \nabla u_{1} | |_{0} .

First, let

{(\partial_{x} u_{1})}_{n}

and

{(\partial_{y} u_{1})}_{n}

be defined as

f_{i, j} - f_{i - 1, j}

and

f_{i, j} - f_{i, j - 1}

, where

n = (i - 1) \times J + j

. Then, two auxiliary variables

(h_{n}, v_{n})

are introduced to replace

({(\partial_{x} u_{1})}_{n}, {(\partial_{y} u_{1})}_{n})

. Equation (16) can be further transformed into the following problem:

\underset{u_{1}, {h_{n}, v_{n}}}{\arg \min} | | f^{(k + 1)} - u_{1} - t_{1}^{(k)} | |_{2}^{2} + λ_{1}^{*} | | \nabla u_{1} | |_{0} + β_{1} ({({(\partial_{x} u_{1})}_{n} - h_{n})}^{2} + {({(\partial_{y} u_{1})}_{n} - v_{n})}^{2}),

where the

λ_{1}^{*} = 2 λ κ / δ_{1}

is a smoothing parameter, which is designed to control the image edge and finer structures.

β_{1}

represents an automatically adaptive parameter to control the similarity between

(h_{n}, v_{n})

and

({(\partial_{x} u_{1})}_{n}, {(\partial_{y} u_{1})}_{n})

.

According to [22,28], Eq. (17) can be updated by the following two alternative steps:

\begin{array}{l} {{(h_{n})}^{(i + 1)}, {(v_{n})}^{(i + 1)}} = \\ {\begin{matrix} (0, 0) & {({({(\partial_{x} (u_{1}))}_{n})}^{(i)})}^{2} + {({({(\partial_{y} (u_{1}))}_{n})}^{(i)})}^{2} \leq λ_{1}^{*} / β_{1} \\ {{({(\partial_{x} (u_{1}))}_{n})}^{(i)}, {({(\partial_{y} (u_{1}))}_{n})}^{(i)}} & otherwise \end{matrix}, \end{array}

u_{1}^{(i + 1)} = F^{- 1} {\frac{F (f^{(i + 1)} - d^{(i)}) + β_{1} (F^{*} (\partial x) F (h^{(i + 1)}) + F^{*} (\partial y) F (v^{(i + 1)}))}{F (1) + β_{1} (F^{*} (\partial x) F (\partial x) + F^{*} (\partial y) F (\partial y))}},

where

F (\cdot)

is Fast Fourier Transform (FFT), and

F^{*} (\cdot)

is the complex conjugate transform of Fourier Transform. The problem Eq. (11c) can be updated with the same steps. The main steps of L₀-PICCS are summarized in Table 1.

Algorithm 1: L0-PICCS

View Table | View all tables in this article

3. Experiments

The more X-ray source/detector pairs mounted on SMCT system, the better reconstructed image quality and the higher the temporal resolution that can be obtained [14]. However, a large number of X-ray source/detector pairs can reduce radius of the field of view (FOV). In addition, the large number of source/detector pairs can increase the manufacturing costs of the system. Meanwhile, the scattering effect of X-ray source will be aggravated from different directions, which will degrade the reconstructed image quality. Therefore, the number of source/detector pairs are set as 7 and 5 and the corresponding sampling views are 700 and 500 in the simulation. In the realistic data set experiment, the sampling views are 450 and 448 respectively.

There are four parameters in the L₀-PICCS need to be selected, i.e., the Lagrangian multiplier $δ_{1}$ and $δ_{2}$ , smoothing parameters $λ_{1}^{*}$ and $λ_{2}^{*}$ . The optimized parameters of all reconstruction methods are picked out by comparing the RMSE and visual evaluation. Tables 2 and 3 list the optimal parameters of L₀-PICCS for different cases. The iteration number for all the algorithms is set as 300. The sampling time is set as 0.1s and other numerical simulations parameters are listed in Table 4.

Table 2. The parameters of L₀-PICCS ( $δ_{1}$ , $δ_{2}$ , $λ_{1}^{*}$ , $λ_{2}^{*}$ ) for realistic simulation.

View Table | View all tables in this article

Table 3. Numerical simulation parameters

View Table | View all tables in this article

Table 4. Quantitative results in terms of RMSEs

View Table | View all tables in this article

3.1 Numerical simulation

To evaluate the performance of the proposed L₀-PICCS, three different undersampling factors, i.e., $w = 10, 25, 50$ are set. The undersampling factor 10 means 70 projections are collected at one time frame for the SMCT system equipped with 7 source/detector pairs. In that case, the temporal resolution of the SMCT system is 1.0s. Similarly, the temporal resolution of the SMCT system can be 0.4s and 0.2s for the undersampling factors 25 and 50.

To demonstrate the advantages of the SMCT system, a phantom (see Fig. 3) was designed to simulate the generation and disappear of bubbles within 5.0s for casting cooling in our previous study [14]. More details about the phantom can refer to [14]. In this work, the range of greyscale value is set as [0 1] rather than [0 1000] in [14]. In this study, to evaluate the denoising performance of the L₀-PICCS, a uniformly distributed Gaussian noise, where the standard variance is set as 10% of the maximum value of noise-free projection, is added to the projection data.

Fig. 3 From left to right columns represent the statues of the casting cooling process at 0.1s, 0.6s and 5.0s, respectively.

Download Full Size | PDF

The high quality prior image plays an important role for both SM-PICCS and L₀-PICCS reconstruction techniques. For each time frame, the SMCT system can only collect few projections. The previous proposed SM-PICCS algorithm achieved a great success by incorporating the prior image. The feasibility strategy is to utilize the traditional FBP method to reconstruct the object image from complete data set. However, it is impossible to obtain such a complete data set for a dynamic changeable object even with the SMCT. The projections from different time frames usually can be treated as a “quasi-complete” data set within a short time window. Finally, the prior image can be reconstructed using FBP from such “quasi-complete” data set. However, the used “quasi-complete” data set is corrupted by data inconsistency, which further implies the prior image with severe artifacts from the data inconsistence.

The reconstructed images consist of 256 $\times$ 256 pixels and each of them covers 1 $\times$ 1 mm². For simplicity, we only analyze the reconstructed images at 0.6s in the case of 5 X-ray sources from the noisy projections and the results are displayed in Fig. 4. From Fig. 4, it can be seen that the results reconstructed by ART [29] and total variation minimization with the steepest descent search (TVM-SD) [30] are seriously destroyed by shadow and motion artifacts. Compared with the ART and TVM-SD results, both SM-PICCS and L₀-PICCS can effectively restore the image information. Furthermore, compared with the SM-PICCS, the L₀-PICCS can more effectively reduce shadow artifacts and preserve the image edge. These advantages have been confirmed by the region indicated by red arrows in Fig. 4. It can be observed that the SM-PICCS cannot be able to separate the adjacent small bubbles with the increase of the sampling factors. However, these blurred bubbles in the SM-PICCS results can be easily distinguished by the L₀-PICCS. Indeed, L₀-PICCS can provide better spatial resolution than the SM-PICCS. To further demonstrate the out-performance of the L₀-PICCS, the differences between reconstructed and true images are also illustrated in Fig. 5. From Fig. 5, we can see that the images from the L₀-PICCS are closer to true image than those reconstructed by other competitors.

Fig. 4 Reconstructed results of the time frame 6 from the noisy projections acquired with 5 X-ray sources. The 1st −3rd rows represents the reconstructed results of the undersampling factors 10, 25, 50, respectively. 1st-3rd columns are the results reconstructed by ART, TVM-SD, SM-PICCS, L₀-PICCS methods, respectively. The display window is [0.04,0.1].

Download Full Size | PDF

Fig. 5 The difference images between reconstructed images in Fig. 4 and true image. The display window is [0.04,0.1].

Download Full Size | PDF

In order to quantitatively compare the reconstructed results, the index values in terms of RMSE, PSNR and SSIM of all reconstruction algorithms with 5 and 7 X-ray sources are listed in Tables 5, 6 and 7, respectively. From Tables 5 and 6, we can see that the L₀-PICCS always can obtain the smallest RMSE and the largest PSNR in the cases of noise-free and noisy projections. Meanwhile, the proposed L₀-PICCS can achieved the highest SSIM values in all cases. This means that the images reconstructed by L₀-PICCS are closer to the ideal phantom than other comparisons. From Tables 5-7, it can be seen that the results reconstructed by ART and TVM-SD from 7 X-ray sources are better than those obtained from 5 X-ray sources. This is because the projections from 7 X-ray sources are more than that achieved from 5 X-ray sources. However, regarding the L₀-PICCS or SM-PICCS results, we can observe that the RMSEs and SSIMs are very close in the cases of 5 and 7 sources. Unexpectedly, the results using the proposed L₀-PICCS in 5 X-ray sources are even better than those reconstructed by SM-PICCS with 7 X-ray sources. In other words, the L₀-PICCS may be has a potential to reduce the number of X-ray source with the same image quality. The conclusion can be also confirmed by the index of PSNRs.

Table 5. Quantitative results in terms of PSNRs

View Table | View all tables in this article

Table 6. Quantitative results in terms of SSIMs

View Table | View all tables in this article

3.2 Realistic data set

To demonstrate the usefulness of the L₀-PICCS algorithm in practical applications, a specimen is scanned by a micro cone-beam CT system with one X-ray source and one flat panel detector (FPD) in Chongqing University. In this study, the tube voltage and current are set as at 60 kV and 200 mAs. The FPD consists of 1024 $\times$ 1024 units and each of them covers an area of 0.127 $\times$ 0.127 mm². The distance starting from the X-ray source to FPD and object are 502.0 mm and 132.3 mm respectively. 450 projections are uniformly acquired over a full scan range. Figure 6(a) shows the reconstructed 3D structure of the specimen using FDK algorithm and Fig. 6(b) shows the 50 slices around the middle plane. From Fig. 6, it can be seen that the structure of the specimen is similar along rotation axis direction, i.e., most parts of structure have not changed except for some small regions. Since the distance from the source to object far outweigh that between source and object, the circular cone-beam geometry can be treated as a good approximation of 3D parallel-beam geometry for such micro-CT system. In fact, this is true for most of the nano-CT and micro-CT scanners. Hence, the extracted 50 slices of the FPD along the z-axis can be treated as the 2D dynamic object. Some representative slices reconstructed from different FPD rows are shown in Fig. 7.

Fig. 6 The 3D structure of the specimen reconstructed by FDK algorithm (left) and extracted 50 reconstructed slices (right).

Download Full Size | PDF

Fig. 7 50 reconstructed slices from 50 FPD slices by using FBP and each reconstructed slice consists of 256 $\times$ 256 pixels. The display window is [0 0.05].

Download Full Size | PDF

To imaging the aperiodic dynamic process in industrial applications, we extract 450 and 448 projections for SMCT systems equipped with 5 and 7 source/detector pairs respectively. The undersampling factors are set as 10, 18, 30 for SMCT system with 5 source/detector pairs and result in the number of projections being 45, 25, 15. Similarly, the undersampling factors are set as 8, 16, 32 for SMCT system with 7 source/detector pairs and lead to the number of projections being 56, 28, 14. To clarify the projection extraction, here the undersampling factor is assumed as 10 for SMCT with 5 source/detector pairs. In the case, the interlaced strategy is adopted to extract 45 projections for one row of FPD.

For simplicity, we only analyze the reconstructed results at time frame 6. The reconstructed images from 6th row of the FPD by using FBP is served as a benchmark (see Fig. 8). Figures 8 and 9 show the reconstructed results for different system configuration with different undersampling factors. From Figs. 8 and 9, we can obviously observe that the results of ART and TVM-SD contain too severe artifacts to the image information is completely destroyed. Compared with ART and TVM-SD results, the SM-PICCS technique can recover high quality image to some extent. However, the image edges and structures indicated with the red arrows in Figs. 8 and 9 are blurred so that some small gaps are disappeared. In contrast, the results are constructed by L₀-PICCS have clearer image edges and much finer structures. Especially, the gaps between small structures can be distinguished clearly. This is because that L₀-PICCS minimizes the L₀-norm of image gradient. It penalizes the non-zero number rather than the amplitude, which is good for image edge protection and finer features recovering.

Fig. 8 Reconstructed results of the time frame 6 from realistic data set acquired with 5 X-ray sources. The 1st-3rd rows represents the reconstructed results of the undersampling factors 10, 18 and 30, respectively. 1st-3rd columns are the results reconstructed by ART, TVM-SD, SM-PICCS, L0-PICCS methods, respectively. The display window is [0,0.05].

Download Full Size | PDF

Fig. 9 Same as Fig. 8 but for 7 X-ray sources.

Download Full Size | PDF

4. Discussion and conclusion

To study the locally aperiodic dynamic imaging in industrial applications, we previously designed a SMCT system. For the SMCT system, the conventional compressed sensing based reconstruction method may be failed to provide high-quality image. To address this problem, we first applied the PICCS method to the SMCT system. In this study, to further improve the temporal resolution of the SMCT system, the L₀-norm of image gradient is incorporated into the PICCS method to obtain an L₀-PICCS algorithm. Both numerical and realistic data set experiments demonstrate the L₀-PICCS method can provide a higher image quality with sharper image edge and finer structures than the SM-PICCS. This is because image gradient L₁-norm minimization can lead to a global smoothing by punishing image gradients coefficients. In contrast, the image gradient L₀-norm minimization penalizes the non-zero number of the small image gradient rather than the amplitude, which can be beneficial for image edge information protection and finer structure recovering.

From the results with different undersampling factors in numerical simulations, it can be seen that the reconstructed image quality is decreased with the undersampling factor increased. Especially, the undersampling factor is degraded to 25, the reconstructed images from L₀-PICCS and SM-PICCS are degraded than those obtained from undersampling factor 10. However, the images reconstructed by the L₀-PICCS provide a higher image quality than those obtained by the SM-PICCS. As the aforementioned, increasing the undersampling factor will improve the temporal resolution. Thus, L₀-PICCS is able to acquire better images quality than SM-PICCS when the temporal resolution increased. Again, compared with the SM-PICCS, the proposed L₀-PICCS method can reconstructed a satisfied image with fewer projections. It indicates the developed method has a potential to reduce the number of the X-ray sources with the same temporal resolution. When the undersampling factor is fixed, the results of L₀-PICCS from 5 X-ray sources are more accurate than that obtained by SM-PICCS from 7 X-ray sources. The realistic data set experiments demonstrate that L₀-PICCS can acquire a better image quality than SM-PICCS. The realistic data set experiments demonstrate that if there are changing region rapidly within the object, both L₀-PICCS and SM-PICCS may fail to reconstruct high quality images. For example, the undersampling factor 32 for SMCT system equipped with 7 X-ray source-detector pairs, the projections are contaminated by data inconsistency and result in image quality are further comprised. Fortunately, the proposed L₀-PICCS can still achieved better reconstructed image with clear image edge and finer structures.

Although the L₀-PICCS algorithm has obtained excellent performance, there are still some issues. First, the L₀-PICCS method contains four parameters, which may become the biggest problem in practical applications. In this study, those parameters are optimized based on RMSE and visual evaluation. Figure 10 further demonstrates the results in terms of the RMSE v.s. four parameters in the case of undersampling factor 10 for SMCT with 5 source/detector pairs. The theoretical analyses and furtherly optimizations are still open problems that need to be fully addressed in our future work [31]. Second, the realistic data set is collected from the micro-CT system equipped with one X-ray source/detector pair. The system ignores the forward and backward X-ray scattering effects. Thus, how to obtain a high quality reconstructed image is a challenge considering scattering effect in practical. Third, the L₀-PICCS is first proposed for SMCT system in this study. In fact, it can be also applied to cardiac CT (4DCT) [1,32,33], spectral CT [34,35]. Particularly, when the proposed L₀-PICCS is applied to spectral CT reconstruction, the averaged image can be considered as the prior image.

Fig. 10 RMSEs v.s. different parameters settings. (a)-(d) represent $δ_{1}, λ_{1}^{*}, δ_{2}, λ_{2}^{*}$ respectively.

Download Full Size | PDF

In summary, based on the advantages of the image gradient L₀-norm minimization in image edge protection and finer features recovering, it was employed to improve SMCT system temporal resolution by combining the prior image. Considering the image gradient L₀-norm optimization is a non-convex and NP-hard problem, the split-Bregman strategy was used to solve the L₀-PICCS model. The experiment results demonstrate the proposed L₀-PICCS can obtain better reconstructed image quality than other comparisons. This will be extremely meaningful for dynamic CT reconstruction.

Funding

National Natural Science Foundation of China (61471070), the NIH/NIBIB U01 (EB017140), and the National Instrumentation Program of China (2013YQ030629).

References

1. M. Holbrook, D. P. Clark, and C. T. Badea, “Low-dose 4D cardiac imaging in small animals using dual source micro-CT,” Phys. Med. Biol. 63(2), 025009 (2018). [PubMed]

2. M. Nakano, A. Haga, J. Kotoku, T. Magome, Y. Masutani, S. Hanaoka, S. Kida, and K. Nakagawa, “Cone-beam CT reconstruction for non-periodic organ motion using time-ordered chain graph model,” Radiat. Oncol. 12(1), 145 (2017). [CrossRef] [PubMed]

3. A. Bauereiß, T. Scharowsky, and C. Körner, “Defect generation and propagation mechanism during additive manufacturing by selective beam melting,” J. Mater. Process Tech. 214(11), 2522–2528 (2014). [CrossRef]

4. R. Schwarze and J. Klostermann, “Computational Fluid Dynamic (CFD) Simulations of Liquid Steel Infiltration in Porous Ceramic Structures: Dynamics of the Penetrating Melt Surface,” Steel Res. Int. 87(4), 465–471 (2016). [CrossRef]

5. P. FitzGerald, P. Edic, H. Gao, Y. Jin, J. Wang, G. Wang, and B. Man, “Quest for the ultimate cardiac CT scanner,” Med. Phys. 44(9), 4506–4524 (2017). [CrossRef] [PubMed]

6. W. M. Thompson, W. R. B. Lionheart, E. J. Morton, M. Cunningham, and R. D. Luggar, “High speed imaging of dynamic processes with a switched source x-ray CT system,” Meas. Sci. Technol. 26(5), 55401 (2015). [CrossRef]

7. Y. Liu, H. Liu, Y. Wang, and G. Wang, “Half-scan cone-beam CT fluoroscopy with multiple x-ray sources,” Med. Phys. 28(7), 1466–1471 (2001). [CrossRef] [PubMed]

8. R. A. Robb, E. A. Hoffman, L. J. Sinak, L. D. Harris, and E. L. Ritman, “High-speed three-dimensional X-ray computed tomography: The dynamic spatial reconstructor,” Proc. IEEE 71(3), 308–319 (1983). [CrossRef]

9. B. Liu, G. Wang, E. L. Ritman, G. Cao, J. Lu, O. Zhou, L. Zeng, and H. Yu, “Image reconstruction from limited angle projections collected by multisource interior x-ray imaging systems,” Phys. Med. Biol. 56(19), 6337–6357 (2011). [CrossRef] [PubMed]

10. E. Y. Sidky and X. Pan, “Image reconstruction in circular cone-beam computed tomography by constrained, total-variation minimization,” Phys. Med. Biol. 53(17), 4777–4807 (2008). [CrossRef] [PubMed]

11. G. H. Chen, J. Tang, and S. Leng, “Prior image constrained compressed sensing (PICCS): a method to accurately reconstruct dynamic CT images from highly undersampled projection data sets,” Med. Phys. 35(2), 660–663 (2008). [CrossRef] [PubMed]

12. G. H. Chen, J. Tang, and J. Hsieh, “Temporal resolution improvement using PICCS in MDCT cardiac imaging,” Med. Phys. 36(6), 2130–2135 (2009). [CrossRef] [PubMed]

13. G. H. Chen, P. Theriault-Lauzier, J. Tang, B. Nett, S. Leng, J. Zambelli, Z. Qi, N. Bevins, A. Raval, S. Reeder, and H. Rowley, “Time-resolved interventional cardiac C-arm cone-beam CT: an application of the PICCS algorithm,” IEEE Trans. Med. Imaging 31(4), 907–923 (2012). [CrossRef] [PubMed]

14. W. Wu, H. Yu, C. Gong, and F. Liu, “Swinging multi-source industrial CT systems for aperiodic dynamic imaging,” Opt. Express 25(20), 24215–24235 (2017). [CrossRef] [PubMed]

15. J. C. Ramirez-Giraldo, J. Trzasko, S. Leng, L. Yu, A. Manduca, and C. H. McCollough, “Nonconvex prior image constrained compressed sensing (NCPICCS): theory and simulations on perfusion CT,” Med. Phys. 38(4), 2157–2167 (2011). [CrossRef] [PubMed]

16. R. Chartrand, “Exact Reconstruction of Sparse Signals via Nonconvex Minimization,” IEEE Signal Proc. Lett. 14(10), 707–710 (2007). [CrossRef]

17. J. Trzasko and A. Manduca, “Highly undersampled magnetic resonance image reconstruction via homotopic l(0) -minimization,” IEEE Trans. Med. Imaging 28(1), 106–121 (2009). [CrossRef] [PubMed]

18. E. Y. Sidky, R. Chartrand, J. M. Boone, and X. Pan, “Constrained TpV Minimization for Enhanced Exploitation of Gradient Sparsity: Application to CT Image Reconstruction,” IEEE J. Transl. Eng. Health Med. 2, 1–18 (2014). [CrossRef] [PubMed]

19. E. Y. Sidky, R. Chartrand, and X. Pan, “Image reconstruction from few views by non-convex optimization,” in IEEE Nuclear Science Symposium Conference Record (2007), p. 3526. [CrossRef]

20. W. Yu, C. Wang, and M. Huang, “Edge-preserving reconstruction from sparse projections of limited-angle computed tomography using ℓ₀-regularized gradient prior,” Rev. Sci. Instrum. 88(4), 043703 (2017). [CrossRef] [PubMed]

21. L. Zhang, L. Zeng, and Y. Guo, “l0 regularization based on a prior image incorporated non-local means for limited-angle X-ray CT reconstruction,” J. XRay Sci. Technol. 26(3), 481–498 (2018). [CrossRef] [PubMed]

22. W. Wu, Y. Zhang, Q. Wang, F. Liu, P. Chen, and H. Yu, “Low-dose spectral CT reconstruction using ℓ₀ image gradient and tensor dictionary,” Appl. Math. Model. 63, 538–557 (2018). [CrossRef]

23. C. Gong, L. Zeng, C. Wang, and L. Ran, “Design and Simulation Study of a CNT-Based Multisource Cubical CT System for Dynamic Objects,” Scanning 2018, 6985698 (2018). [CrossRef] [PubMed]

24. Y. Sun and J. Tao, “Image reconstruction from few views by l(0)-norm optimization,” Chinese Phys. B 23(7), 078703 (2014). [CrossRef]

25. M. Li, C. Zhang, C. Peng, Y. Guan, P. Xu, M. Sun, and J. Zheng, “Smoothed l₀ Norm Regularization for Sparse-View X-Ray CT Reconstruction,” BioMed Res. Int. 2016, 2180457 (2016). [CrossRef] [PubMed]

26. L. Zeng and C. Wang, “Error bounds and stability in the l0 regularized for CT reconstruction from small projections,” Inverse Probl. Imaging (Springfield) 10(3), 829–853 (2016). [CrossRef]

27. T. Goldstein and S. Osher, “The Split Bregman Method for L1-Regularized Problems,” SIAM J. Imaging Sci. 2(2), 323–343 (2009). [CrossRef]

28. L. Xu, C. Lu, Y. Xu, and J. Jia, “Image smoothing via L0 gradient minimization,” in SIGGRAPH Asia Conference (2011), p. 174.

29. R. Gordon, R. Bender, and G. T. Herman, “Algebraic reconstruction techniques (ART) for three-dimensional electron microscopy and x-ray photography,” J. Theor. Biol. 29(3), 471–481 (1970). [PubMed]

30. E. Y. Sidky and X. Pan, “Image reconstruction in circular cone-beam computed tomography by constrained, total-variation minimization,” Phys. Med. Biol. 53(17), 4777–4807 (2008). [CrossRef] [PubMed]

31. P. R. Johnston and R. M. Gulrajani, “Selecting the corner in the L-curve approach to Tikhonov regularization,” IEEE Trans. Biomed. Eng. 47(9), 1293–1296 (2000). [CrossRef] [PubMed]

32. D. Richter, H. I. Lehmann, A. Eichhorn, A. M. Constantinescu, R. Kaderka, M. Prall, P. Lugenbiel, M. Takami, D. Thomas, C. Bert, M. Durante, D. L. Packer, and C. Graeff, “ECG-based 4D-dose reconstruction of cardiac arrhythmia ablation with carbon ion beams: application in a porcine model,” Phys. Med. Biol. 62(17), 6869–6883 (2017). [CrossRef] [PubMed]

33. S. Kim, Y. Chang, and J. B. Ra, “Cardiac motion correction based on partial angle reconstructed images in x-ray CT,” Med. Phys. 42(5), 2560–2571 (2015). [PubMed]

34. W. Wu, F. Liu, Y. Zhang, Q. Wang, and H. Yu, “Non-local Low-rank Cube-based Tensor Factorization for Spectral CT Reconstruction,” IEEE Trans. Med. Imaging1 (2018).

35. W. Wu, Y. Zhang, Q. Wang, F. Liu, F. Luo, and H. Yu, “Spatial-spectral cube matching frame for spectral CT reconstruction,” Inverse Probl. 34(10), 104003 (2018). [CrossRef]

	Sources	$w$	$δ_{1}$ (10⁻²)	$δ_{2}$ (10⁻²)	$λ_{1}^{*}$ (10⁻²)	$λ_{2}^{*}$ (10⁻⁸)
Noise-free	5	10	5	2	2.00	1.5
		25	5	2	2.00	1.5
		50	3	500	2.00	4
	7	10	14.3	333	2.00	5000
		25	5	2	2.00	5000
		50	2	2	2.00	5000
Noise	5	10	6	200	2.50	500
		25	5	450	2.50	1.5
		50	5	450	2.90	1.5
	7	10	27.3	85.9	1.58	10000
		25	8.2	85.9	1.58	500
		50	5	2	2.00	1.5

Sources	$w$	$δ_{1}$	$δ_{2}$ (10⁻²)	$λ_{1}^{*}$ (10⁻⁵)	$λ_{2}^{*}$ (10⁻⁷)
5	10	0.2	2.4	0.5	0.11
	18	0.2	2.4	0.6	0.14
	30	0.16	2.24	0.6	0.11
7	8	0.2	3.5	1.3	743
	16	0.2	2.75	1.4	238
	32	0.2	2.75	1.4	238

Parameters	Q = 5	Q = 7
Distance from source to rotation center (mm)	577	640
Distance from source to detector (mm)	977	1080
Diameter of field of view (mm)	60	60
Efficient detector array length (mm)	204.8	204.8
Detector pixel size (mm)	0.8	0.8
Interval angle between two adjacent projection views	$π / 250$	$π / 350$
Reconstructed image size (pixel)	256 $\times$ 256	256 $\times$ 256

		Q = 5			Q = 7
Undersampling factors		10	25	50	10	25	50
Noise-free	ART	0.255	0.289	0.322	0.239	0.240	0.262
	TVM-SD	0.232	0.281	0.321	0.173	0.217	0.253
	SM-PICCS	0.056	0.058	0.068	0.056	0.059	0.069
	L₀-PICCS	0.043	0.052	0.064	0.041	0.046	0.059
Noisy	ART	0.259	0.290	0.324	0.254	0.244	0.264
	TVM-SD	0.232	0.282	0.322	0.175	0.217	0.254
	SM-PICCS	0.066	0.068	0.076	0.069	0.067	0.075
	L₀-PICCS	0.055	0.058	0.072	0.059	0.062	0.070

		Q = 5			Q = 7
Undersampling factors		10	25	50	10	25	50
Noise-free	ART	12.368	10.908	9.850	14.317	12.790	11.709
	TVM-SD	12.711	11.025	9.862	15.258	13.293	11.939
	SM-PICCS	25.056	24.694	23.353	25.068	24.743	23.533
	L₀-PICCS	27.217	25.729	23.810	27.788	26.720	24.508
Noisy	ART	12.242	10.869	9.817	13.964	12.660	11.647
	TVM-SD	12.681	11.011	9.843	15.129	13.260	11.915
	SM-PICCS	23.627	23.396	22.407	23.651	23.592	22.705
	L₀-PICCS	25.112	24.583	22.809	24.549	24.179	23.035

Image gradient L₀-norm based PICCS for swinging multi-source CT reconstruction

Abstract

1. Introduction

2. Method

2.1 SMCT system

2.2 L₀-PICCS model

3. Experiments

3.1 Numerical simulation

3.2 Realistic data set

4. Discussion and conclusion

Funding

References

Cited By

Figures (10)

Tables (7)

Equations (23)

Optics Express

		Q = 5			Q = 7
Undersampling factors		10	25	50	10	25	50
Noise-free	ART	0.8448	0.7568	0.6586	0.9100	0.8606	0.8044
	TVM-SD	0.8555	0.7615	0.6567	0.9272	0.8747	0.8125
	SM-PICCS	0.9935	0.9930	0.9904	0.9936	0.9931	0.9908
	L₀-PICCS	0.9960	0.9943	0.9913	0.9965	0.9955	0.9925
Noisy	ART	0.8409	0.7552	0.6561	0.9032	0.8569	0.8023
	TVM-SD	0.8545	0.7608	0.6549	0.9250	0.8738	0.8117
	SM-PICCS	0.9909	0.9904	0.9878	0.9909	0.9908	0.9887
	L₀-PICCS	0.9935	0.9926	0.9887	0.9924	0.9918	0.9893

Abstract

1. Introduction

2. Method

2.1 SMCT system

2.2 L0-PICCS model

3. Experiments

3.1 Numerical simulation

3.2 Realistic data set

4. Discussion and conclusion

Funding

References

Cited By

Figures (10)

Tables (7)

Equations (23)

Optics Express

2.2 L₀-PICCS model