Improving radiation dose efficiency of X-ray differential phase contrast imaging using an energy-resolving grating interferometer and a novel rank constraint

Yongshuai Ge; Ran Zhang; Ke Li; Guang-Hong Chen

doi:10.1364/OE.24.012955

1. Introduction

In recent years, grating interferometer based x-ray differential phase contrast (DPC) imaging systems [1–3] have been extensively investigated with the hope of finding their niche in medical applications [4–9]. Combining a Talbot-Lau interferometer and a conventional x-ray imaging setup with a medical grade tube-detector assembly may put less stringent requirements on the x-ray source and detector to facilitate the translation into medical applications. However, one challenge to be overcome is to improve the radiation dose efficiency, e.g. measured by the dose normalized contrast to noise ratio (CNRD), of the DPC imaging systems with the Talbot-Lau interferometer. For a quantum limited DPC imaging system, the noise level of a DPC image is dependent on both the detected photon flux and the fringe visibility of the Talbot-Lau interferometer. The fringe visibility in Talbot-Lau interferometers has a strong dependency on energy. If the x-ray beam energy differs greatly from the energy where fringe visibility is maximized for a given interferometer (this energy is referred to as the designed energy throughout this paper), the measured visibility of the system drops considerably. For a given interferometer system with a designed energy for maximal fringe visibility, the maximal fringe visibility cannot be achieved when information from x-rays with a wide energy range is integrated at the detector. In fact, although x-ray photons with energies deviated from the designed value do carry phase information about the image object, due to the reduction of fringe visibility at these energies, they also greatly add noise to the measured DPC signals and thus significantly degrade the radiation dose efficiency of the DPC imaging systems provided that an energy integration detector (EID) is used. In this regard, an energy-resolving photon counting detector (PCD) does have potential benefits, as it enables a Talbot-Lau interferometer to be implemented as an energy-resolving grating interferometer so that one can selectively use x-ray photons to achieve the highest fringe visibility in DPC imaging.

With an energy-resolved interferometer, one can adjust the center of the energy window to match with the designed energy of the interferometer for the best fringe visibility. The next question is how to select an optimal energy window width to obtain the best DPC imaging performance. If the energy bin width is too wide, then the overall fringe visibility will drop due to the averaging effect over different beam energies within the energy window. However, if the energy bin width is too narrow, the number of photons detected within the energy window drops, leading to an increase in DPC image noise and thus a decrease in radiation dose efficiency. Therefore, an optimal energy bin width should be selected to balance the drops in fringe visibility and the number of detected photons within the selected energy bin [10, 11]. While it is straightforward to search for an optimal energy window when the system is operated at the first Talbot order, the definition of an optimal energy window may become tricky when the system is operated at higher Talbot orders [12, 13] due to the existence of multiple local maxima and minima in visibility-energy response curve.

Alternatively, when multiple energy thresholds are available for use in PCD, one can simultaneously acquire DPC images in multiple energy windows and then combine the acquired images from different energy windows with a weight, so that the signal to noise ratio in the final combined DPC image can be maximized [14]. The optimal weight in each energy bin has been found to be inversely proportional to the square of the energy and to the noise variance of the DPC signal in the corresponding energy window. A potential challenge to this method is the need to accurately calibrate the photon energy at each energy bin. It has been demonstrated that smaller weights should be assigned to both high and low energy bins that deviate from the central energy bin with highest fringe visibility.

In this paper, a novel method was developed to reduce noise of DPC images acquired using an energy-resolving interferometer system regardless of energy window partitioning schemes. This method utilizes information from all energy bins, and the method depends on neither the image object itself, nor the number of energy bins used in experiments, nor the need of accurate energy calibration at each energy bin. Numerical simulations and experimental phantom studies have been performed to validate this method.

2. Theory and method

2.1. Fundamental physics in DPC imaging

The Talbot-Lau interferometer based x-ray DPC imaging technique [3] measures the differential phase information ϕ (x,y,E) of an object at a given x-ray beam energy E. In theory, the measured DPC signal from a Talbot-Lau grating system is defined as:

ϕ (x, y, E) = \frac{C}{E^{2}} \frac{\partial}{\partial y} \int_{L} ρ_{e} (x, y, z) d z,

where C is a constant, E is the x-ray beam energy, and ρ_e(x,y,z) represents the electron density distribution of the object. Herein, the derivative is taken along the y-axis, which is assumed to be perpendicular to the grooves in the gratings. In reality, the measurement is performed over a finite energy window. In this case, the beam energy E in Eq. (1) should be interpreted as the effective energy E_eff of the used energy window.

To facilitate the discussion, it is assumed that an ideal photon counting detector is used in the DPC data acquisitions, and a phase-stepping method [2] is used to acquire the data. The detected signal from a given detector pixel (p,q) (where p,q are the pixel indices corresponding to the horizontal and vertical directions on the detector plane) at the k-th phase stepping position can be modeled as [2,3,15]

N^{(k)} (p, q, E) = N_{0} (p, q, E) + N_{1} (p, q, E) \cos [2 π \frac{k}{M} + ϕ (p, q, E)],

where k = 1,2,…,M. The total number of phase steps is M, and N^k(p,q,E) is the measured photon number at energy E for the k-th phase stepping position. The three energy-specific parameters, N₀, N₁, and ϕ, which correspond to three contrast mechanisms, i.e., attenuation contrast, dark-field contrast, and differential phase contrast, respectively, can be estimated from the M measurements. For example, one of the most commonly used approach to retrieve phase signal from detector counts is given by

ϕ (p, q, E) = \tan^{- 1} [\frac{\sum_{k} N^{(k)} (p, q, E) \sin (2 π \frac{k}{M})}{\sum_{k} N^{(k)} (p, q, E) \cos (2 π \frac{k}{M})}] .

Therefore, phase retrieval is generally not a linear combination of detector counts.

The intrinsic quantum fluctuations of the measured photon number N⁽^k⁾ will be propagated into the three estimated quantities of interest (N₀, N₁, and ϕ). In particular, theoretical studies [16–22] have demonstrated that the noise variance of the extracted differential phase signal, i.e., $σ_{ϕ (E)}^{2}$ , can be written as

σ_{ϕ (E)}^{2} = \frac{1}{M N_{0} (E)} \times \frac{2}{ε^{2} (E)},

where ε(E) = N₁(E)/N₀(E) is defined as the fringe visibility of the interferometer.

2.2. Spatial-spectral DPC image matrix and the nature of its low dimensionality

If a DPC image dataset is acquired using an energy-resolving photon counting detector with Λ energy bins, the DPC image in the i-th energy window can be indexed as ϕ(p,q,E_i),1 ≤ i ≤ Λ. To simplify the presentation, each two-dimensional DPC image ϕ (p,q,E_i) is vectorized using the lexicographic order as follows:

\vec{ϕ} (p, q, E_{i}) = [\begin{matrix} \begin{array}{l} ϕ (1, 1, E_{i}) \\ ϕ (1, 2, E_{i}) \end{array} \\ ⋮ \\ \begin{array}{l} ϕ (2, 1, E_{i}) \\ ϕ (2, 2, E_{i}) \end{array} \\ ⋮ \\ ϕ (P, Q, E_{i}) \end{matrix}],

where P and Q are the largest index values along the horizontal and vertical directions of the detector plane, respectively. Using this vector notation, DPC images acquired from all energy windows can form a spatial-spectral image matrix J as follows:

J = [\vec{ϕ} (p, q, E_{1}) \dots \vec{ϕ} (p, q, E_{Λ})] = [\begin{matrix} ϕ (1, 1, E_{1}) & \dots & ϕ (1, 1, E_{Λ}) \\ ϕ (1, 2, E_{1}) & \dots & ϕ (1, 2, E_{Λ}) \\ ⋮ & ⋮ & ⋮ \\ ϕ (P, Q, E_{1}) & \dots & ϕ (P, Q, E_{Λ}) \end{matrix}] .

Note that this spatial-spectral image matrix J has P × Q rows and Λ columns. Each row of J can be regarded as the change in DPC signal for each detector pixel from one energy window to the other.

A key observation in this paper is that, the i-th and j-th columns in the spatial-spectral DPC image matrix J only differ from each other by a global scaling factor of $E_{j}^{2} / E_{i}^{2}$ [see Eq. (1)], provided that there is no noise or other non-idealities in DPC signal. Subsequently, each column is linearly depended on all other columns, and thus the rank of the matrix J is 1. This observation is independent of the number of energy windows used in DPC imaging. When noise and other non-ideal factors are present in the measured data, the rank of matrix J will be elevated. Therefore, the intrinsic rank-one nature of the spatial-spectral image matrix can be utilized to reject the contributions from noise and other non-ideal confounding factors.

2.3. Augmented spatial-spectral matrix J_P

When a total of Λ energy bins are used, additional DPC images can be formed by combining the detected photons from two or more energy bins. The DPC images generated in this way are referred to as prior images in this paper, which can be used to augment the spatial-spectral image matrix J to form a new spatial-spectral image matrix J_P. If a total of λ prior DPC images are formed, an augmented matrix J_P with Λ+λ columns can be written as follows:

J_{P} = [\begin{matrix} ϕ (1, 1, E_{1}) & \dots & ϕ (1, 1, E_{Λ}) & ϕ (1, 1, E_{1}^{p}) & \dots & ϕ (1, 1, E_{λ}^{p}) \\ ϕ (1, 2, E_{1}) & \dots & ϕ (1, 2, E_{Λ}) & ϕ (1, 2, E_{1}^{p}) & \dots & ϕ (1, 2, E_{λ}^{p}) \\ ⋮ & ⋮ & ⋮ & ⋮ & ⋮ & ⋮ \\ ϕ (P, Q, E_{1}) & \dots & ϕ (P, Q, E_{Λ}) & ϕ (P, Q, E_{1}^{p}) & \dots & ϕ (P, Q, E_{λ}^{p}) \end{matrix}] .

All the prior DPC images $ϕ (p, q, E_{l}^{p})$ added in this paper were generated as follows: a direct summation of the detected photon numbers of different energy windows was performed before Eq. (3) was used to extract the corresponding DPC image. Due to the nonlinear nature of phase retrieval as shown in Eq. (3), an augmented prior DPC image vector $ϕ (p, q, E_{l}^{p})$ is not equivalent to a linear combination of the first Λ columns in J_P. Clearly, the rank-one nature of the spatial-spectral image matrix J is always inherited by the augmented spatial-spectral matrix J_P provided that there is no noise or other confounding factors. However, due to the presence of noise and other confounding factors in the measured data, the augmentation of these prior images to the original spatial-spectral image matrix J can be found to be beneficial in noise reduction [23].

2.4. Rank-one approximation of the augmented spatial-spectral matrix J_P

Due to the presence of noise and other non-ideal confounding factors in DPC imaging, the rank of an experimentally measured spatial-spectral DPC image matrix J_P may deviate from the ideal value of 1 (just as it has for matrix J). Therefore, the intrinsic rank-one nature can be utilized to generate a low rank approximated spatial-spectral matrix $J_{P}^{*}$ with reduced image noise levels. The rank of a matrix can be conveniently determined using the well-known singular value decomposition (SVD) method [24], which is a generalized spectral decomposition of the matrix:

J_{P} = U \sum V^{T} = \sum_{i = 1}^{Λ + λ} σ_{i} {\vec{u}}_{i} \otimes {\vec{v}}_{i}^{T} .

In Eq. (8)U and V are two orthonormal matrices, and Σ is a diagonal matrix with entries σ₁ ≥ σ₂ ≥ ⋯ ≥ σ_Λ+_λ. After decomposition, each column in U only contains the spatial information of the object, while energy dependent information appears only in V^T. The singular values σ_i provide the weighting factors for the corresponding spatial-spectral basis ${\vec{u}}_{i} \otimes {\vec{v}}_{i}^{T}$ .

Following singular value decomposition, the rank-one nature of the augmented spatial-spectral image matrix J_P can be exploited to generate a noise reduced spatial-spectral image matrix $J_{P}^{*}$ as follows:

J_{P}^{*} \approx σ_{1} {\vec{u}}_{1} \otimes {\vec{v}}_{1}^{T} .

In total, matrix $J_{P}^{*}$ contains Λ+λ columns. Each column restores the denoised DPC information for each corresponding energy window. It will be demonstrated that, with this rank-one approximation, noise is reduced and the accuracy of the DPC signal for each energy window is maintained.

2.5. Methods to generate prior DPC image vectors

As an example, Fig. 1(a) shows a PCD with two energy thresholds: a low energy threshold $E_{T H}^{1}$ and a high energy threshold $E_{T H}^{2} (E_{T H}^{1} < E_{T H}^{2})$ . Let $N^{(k)} (E_{T H}^{1})$ and $N^{(k)} (E_{T H}^{2})$ represent the number of photons detected above the corresponding energy threshold at k-th phase stepping position. Hence, two energy bins can be generated: Bin 1 and Bin 2. Assuming the effective mean energy of Bin 1 is E₁, then the number of detected photons in Bin 1 is given by $N^{(k)} (E_{1}) = N^{(k)} (E_{T H}^{1}) - N^{(k)} (E_{T H}^{2})$ . Similarly, the detected photon number in Bin 2 is $N^{(k)} (E_{2}) = N^{(k)} (E_{T H}^{2})$ . The DPC signals in each energy bin, ϕ (p,q,E₁) and ϕ(p,q,E₂), can be estimated using Eq. (3). For a PCD equipped with three or more energy thresholds, the detected photons in each energy bin can be measured using the same method.

Fig. 1 The full x-ray spectrum is equally divided into either two consecutive energy bins (a) or three consecutive energy bins (b). Energy threshold is denoted by $E_{T H}^{i}$ .

Download Full Size | PDF

In this paper, only one prior DPC image $ϕ (p, q, E_{1}^{p})$ has been added and examined for different energy threshold settings. It is retrieved directly from the detected photons above the lowest energy threshold $E_{T H}^{1}$ .

3. Materials and methods

3.1. Numerical simulation studies

To motivate and validate the use of the rank-one assumption in this paper, wave-optical simulations were performed in the first place. The technical details of the numerical simulation framework was presented in [25]. In this work, a 40 kVp spectrum was simulated with an energy step size of 0.1 keV. The simulated Talbot-Lau interferometer system consists of three gratings. The parameters for each grating were given as follows: The source grating G0 with 40.0 μm gold bars and a period of 36.0 μm; The G1 grating with 32.0 μm silicon bars and a period of 8.0 μm (which corresponds to a π phase shift for 25 keV x-rays); and the G2 grating with 40.0 μm gold bars and with a period of 4.5 μm. The simulated system was operated at the first order Talbot distance. An ideal photon counting detector with a pixel size of 100 μm was simulated. The phantom used in the numerical simulation consisted of three cylinders with different materials (PMMA, breast tissue and adipose) and different diameters (4.0 mm, 3.0 mm and 5.0 mm). Poisson noise was simulated for each detector pixel.

3.2. Experimental studies

To validate the proposed rank-one approximation method, the spatial-spectral DPC image datasets were acquired on a benchtop system equipped with a grating interferometer and a photon counting detector. The interferometer was operated at the first order Talbot distance, and it consists of a silicon based phase grating G1 and an absorption grating G2. Both gratings have a duty cycle of 0.5. By design, the G1 grating has a period of 8.0 μm and induces a π-phase shift for mono-energetic x-ray photons of 28 keV. The G2 grating has a period of 4.8 μm and the x-ray absorber is made from 60.0 μm thick gold. In addition, the distance between the x-ray focal spot and the G1 grating is 108.4 cm, and the G2 grating is placed at 21.7 cm downstream the G1 grating. The x-ray source used in this work is a micro-focus x-ray tube (model: L10321, Hamamatsu Photonics K.K., Japan) with a tungsten target. It was operated at 40 kVp (mean energy is about 28 keV) and 150 μA, and the nominal focal spot size is 8 μm. In the experimental data acquisitions, the full x-ray spectrum was equally divided into different energy bins: from 2 bins up to 7 bins. To study the denoising performance of the proposed rank-one approximation method at different dose levels, 30 repeated frames were acquired at each phase stepping position during a total exposure time of 20 seconds. The interferometer system is equipped with an energy-resolved single photon counting detector (model: XC-FLITE X1, XCounter AB, Sweden). The x-ray sensor is made of CdTe and the detector has two adjustable energy thresholds. Its native pixel size is 100 μm, and a 2 × 2 binning mode was used in the data acquisitions.

Two different physical phantoms were used to evaluate the performance of the proposed rank-one method. The first phantom is a polymethylmethacrylate (PMMA) tube with an outer diameter of 9.6 mm and inner diameter of 8.1 mm, see Fig. 2(a). The phantom contains five groups of spheres submerged in a vegetable oil bath. They are ordered from top to bottom as follows: Acrylic ball with a diameter of 6.3 mm; Polyoxymethylene (POM) balls with a diameter of 6.3 mm; and four POM balls with a diameter of 3.2 mm. The other phantom is a PMMA wedge phantom with a slope of 10 degrees, see Fig. 2(b). The constant refraction angle generated by this phantom enabled the contrast-to-noise ratio (CNR) to be quantified for DPC imaging.

Fig. 2 (a) PMMA tube phantom containing a vegetable bath and spheres of different materials and diameters. (b) PMMA wedge phantom with a slope of 10 degrees. Three regions correspond to different PMMA thicknesses were selected, A, B, and C, to evaluate the DPC image CNRs.

Download Full Size | PDF

4. Results

4.1. Numerical simulation results

Numerical simulation results are presented in this section to demonstrate the rank-one nature of the noiseless DPC signal and to validate the proposed low rank approximation method to reduce noise and thus improve dose efficiency.

4.1.1. Demonstration of rank-one nature of the DPC signal

Figure 3 shows the singular value distribution of the spatial-spectral DPC image matrix J at several different number of energy bins when the quantum noise is absent and present. As can be observed, when noise is absent, only the first singular value is dominant, which demonstrated the rank-one nature of the DPC signal. When noise is present, smaller singular values may be elevated and the nature of the proposed rank-one approximation is to discard the contributions from these smaller singular values since the corresponding basis image contain primarily noise.

Fig. 3 The singular value distributions of the spatial-spectral DPC image matrix J with and without quantum noise for energy bin number Λ = 3,5,7,9.

Download Full Size | PDF

4.1.2. Validation of rank-one approximation in noise reduction

When quantum noise was introduced, the proposed rank-one approximation was able to generate noise reduced DPC images while well preserve the image information content, as shown in Fig. 4 for three energy bins. The results for more energy bins showed the same trend and thus were not presented here.

Fig. 4 DPC images generated by numerical simulations. Denoising using rank-one approximation. All images are displayed with a range of [−0.9 0.9].

Download Full Size | PDF

4.2. Experimental results

4.2.1. Measured beam spectrum and fringe visibility of the grating interferometer

By adjusting the energy thresholds of the PCD with a step size of 3 keV, both the x-ray transmission spectrum at the detector and the fringe visibility ε(E) of the used interferometer were experimentally measured. Results shown in Fig. 5 verify that both the spectrum of x-rays incident on the detector plane and the fringe visibility peaked near the designed operation energy (28 keV) of the grating interferometer.

Fig. 5 (a) Number of detected x-ray photons at different energies. The peak corresponds to 27 keV. (b) Measured fringe visibility at different energies. The peak corresponds to 28 keV.

Download Full Size | PDF

4.2.2. Singular values and spatial basis of the augmented spatial-spectral image matrix

To validate the assumption that the DPC information is primarily contained in the first singular value and the corresponding basis, while higher order singular values and bases contribute primarily to noise, the singular values and the corresponding spatial bases for the PMMA tube phantom are presented in Fig. 6. As clearly shown in the results, the first singular value σ₁ is much larger than the second and third two singular values. Similarly, only the first spatial basis ${\vec{u}}_{1}$ contains meaningful object information. In contrast, the other two spatial bases, ${\vec{u}}_{2}$ and ${\vec{u}}_{3}$ , contain negligible information of the image object. These experimental results validated the rank-one approximation and the feasibility of noise reduction by neglecting higher order singular values and spatial bases.

Fig. 6 (a) Magnitude of singular values σ_i of the spatial-spectral matrix J_P. (b)–(d) correspond to the three columns ( ${\vec{u}}_{1}$ , ${\vec{u}}_{2}$ , and ${\vec{u}}_{3}$ ) of the decomposed spatial basis matrix U. Clearly, the object information is predominantly contained in the first basis ${\vec{u}}_{1}$ .

Download Full Size | PDF

4.2.3. Experimental DPC images

DPC images generated without the proposed rank-one approximation for both phantoms were shown in the top row of Fig. 7. Clearly, the image noise levels varied from one energy bin to another. Images generated using the proposed rank-one approximation method are presented in the second row. Noise in these images was significantly reduced. Difference images between the “Original” images and the “SVD:J_P” images are also presented in the bottom row, and there is no noticeable loss of structural information associated with the proposed rank-one approximation method. The line profiles of the tube phantom in Fig. 8 demonstrate that the spatial edge profile of the phantom is well preserved by the proposed rank-one approximation method.

Fig. 7 Experimental imaging results for both the tube and the wedge phantom. The first image row show the original images of different energy bins. Images processed by the proposed rank-one approximation method are listed in the second row. The difference between the first two rows are listed in the bottom row. Images of the wedge phantom, the upper dark part corresponds to air, and the lower bright part corresponds to the PMMA wedge. The three ROIs used for quantitative measurements are marked as: A, B, and C. All images are displayed with a range of [−0.5 1.5].

Download Full Size | PDF

Fig. 8 Line profiles of the tube phantom (central image region) in each energy bin.

Download Full Size | PDF

4.2.4. Quantitative assessment of signal value and contrast-to-noise ratio

Quantitative measurements were performed with images of the wedge phantom to evaluate the performance of the proposed rank-one approximation method. Mean DPC signal values and standard deviations were measured in three rectangular regions (A, B, and C) in Fig. 7, and the results are tabulated in Table 1. The CNRs of DPC images are presented in Fig. 9.

Table 1. List of mean±standard deviation values of DPC images generated with conventional method and with the proposed rank-one approximation methods.

View Table

Fig. 9 CNR of the three energy bins measured at ROIs A-C that correspond to regions in the wedge phantom with different thicknesses.

Download Full Size | PDF

There are several observations based on these results: (1) the proposed rank-one approximation method is able to maintain the signal accuracy (less than ±5% on average); (2) the proposed rank-one approximation method reduces noises for all energy bins, and thus is able to improve CNR; (3) the image CNRs cross all energy bins become identical after the rank-one approximation; (4) the augmented rank-one approximation method using matrix J_P can further reduce image noise and thus improve image CNRs.

4.3. Dependence on radiation exposure levels

Figure 10 shows CNR values of DPC images (region B) acquired at different radiation exposure levels, and three processing methods were compared: The first method, denoted as “Original-Max”, used only a single energy bin that achieved the maximum CNR across different energy bins. The second method, denoted as “SVD:J”, used the rank-one approximation and matrix J without the augmented image vectors. The third method, referred to as “SVD:J_P”, used the augmented matrix J_P and the proposed rank-one approximation.

Fig. 10 Comparisons of the measured CNRs from three different methods for different exposure levels. The linear relationship between the square of CNR and the exposure level is well maintained after the rank-one approximations.

Download Full Size | PDF

As shown in Fig. 10, all three methods demonstrated linear relationships between the measured squared CNR and radiation exposure level, and the proposed rank-one approximation method improved CNR at all exposure levels involved in this study. Note that the CNR improvements of the DPC images processed by the proposed rank-one approximation methods might be limited at very low exposure levels.

4.4. Dependence on total number of energy bins

Figure 11 plots the measured CNR value as a function of total number of energy bins. For a given radiation dose level, the CNR of DPC images decreased with increasing total number of energy bins for all three methods, as the number of photons in each bin decreased. For a given total number of energy bins and given exposure level, the proposed rank-one approximation method with the augmented matrix always generated the highest CNR (and therefore the highest radiation dose efficiency).

Fig. 11 Comparisons of the measured CNRs from three different methods for different energy window numbers.

Download Full Size | PDF

5. Discussion and summary

With an energy-resolving photon counting detector, a novel rank-constrained method was developed to reduce DPC image noise without sacrificing the DPC signal accuracy. Through both the numerical simulations and experimental phantom studies, it was demonstrated that the CNR of DPC images can be improved by the proposed method, therefore the radiation dose efficiency of the DPC imaging systems can be potentially improved.

It was demonstrated that the rank-one approximation method always leads to matched CNR across all energy windows. This is a direct result of the rank-one constraint. Because of the rank-one nature of the spatial-spectral image matrix J and J_P, a common spatial basis, u₁ is shared by images in different energy bins after the rank-one approximation, and this spatial basis is scaled by a numerical scaling factor determined by the effective energy of each energy bin to generate the corresponding DPC image for that energy bin.

The interferometer system used in this study was operated at the first fractional Talbot order, at which the fringe visibility of the diffraction fringe pattern drops gradually when the beam energy deviates from the designed operation energy of the interferometer. At higher fractional Talbot orders, the fringe visibility may reach zero at multiple energy levels [12,13]. In practice, however, the width of an energy bin width is always finite, and the energy resolution of the PCD is always limited, As a result, these theoretical zero-crossings of the visibility-energy curve are usually smeared out with polychromatic x-ray radiations. Due to this reason, we expect the proposed method is still applicable to high order Talbot interferometer setup, although it needs to be experimentally validated in future work.

In summary, a method that exploits the intrinsic rank-one nature of x-ray differential phase contrast images generated by different energy bins of an energy-resolving PCD was developed. Both numerical simulations and experimental validation studies were performed, and it was demonstrated that the proposed method can improve the contrast-to-noise ratio and radiation dose efficiency of differential phase contrast imaging systems.

Acknowledgments

The authors would like to thank John Hayes for his editorial assistance and stimulating discussions. This work was partially supported by the National Institute of Biomedical Imaging and Bioengineering of the National Institutes of Health under Award Number R01EB020521, and by the Office of the Assistant Secretary of Defense for Health Affairs, through the Breast Cancer Research Program, under Award No. W81XWH-16-1-0031. Opinions, interpretations, conclusions and recommendations are those of the author; they does not necessarily represent the official views of the National Institutes of Health, and are not necessarily endorsed by the Department of Defense.

References and links

1. A. Momose, S. Kawamoto, I. Koyama, Y. Hamaishi, H. Takai, and Y. Suzuki, “Demonstration of x-ray Talbot interferometry,” Jpn. J. Appl. Phys. Part 2 42(7B), 866–868 (2003). [CrossRef]

2. T. Weitkamp, A. Diaz, C. David, F. Pfeiffer, M. Stampanoni, E. Zeigler, and P. Cloetens, “X-ray phase imaging with a grating interferometer,” Opt. Express 12(16), 6296–6304 (2005). [CrossRef]

3. F. Pfeiffer, T. Weitkamp, O. Bunk, and C. David, “Phase retrieval and differential phase-contrast imaging with low-brilliance X-ray sources,” Nat. Phys. 2(4), 258–261 (2006). [CrossRef]

4. D. Stutman, T. Beck, J. Carrino, and C. Bingham, “Talbot phase-contrast x-ray imaging for the small joints of the hand,” Phys. Med. Bio. 56(17), 5697–5720 (2011). [CrossRef]

5. J. Herzen, M. Willner, A. Fingerle, P. Noël, T. Köhler, E. Drecoll, E. Rummeny, and F. Pfeiffer, “Imaging liver lesions using grating-based phase-contrast computed tomography with bi-lateral filter post-processing,” PloS one 9(1), e83369 (2014). [CrossRef] [PubMed]

6. K. Li, Y. Ge, J. Garrett, N. Bevins, J. Zambelli, and G-H. Chen, “Grating-based phase contrast tomosynthesis imaging: Proof-of-concept experimental studies,” Med. Phys. 41(1), 011903 (2014). [CrossRef] [PubMed]

7. N. Hauser, Z. Wang, R. Kubik-Huch, M. Trippel, G. Singer, M. Hohl, E. Roessl, T. Köhler, U. van Stevendaal, N. Wieberneit, Nataly, and M. Stampanoni, , “A study on mastectomy samples to evaluate breast imaging quality and potential clinical relevance of differential phase contrast mammography,” Invest. Radiol. 49(3), 131–137 (2014). [CrossRef]

8. Y. Ge, K. Li, J. Garrett, and G-H. Chen, “Grating based x-ray differential phase contrast imaging without mechanical phase stepping,” Opt. Express 22(12), 14246–14252 (2014). [CrossRef] [PubMed]

9. J. Garrett, Y. Ge, K. Li, and G-H. Chen, “Anatomical background noise power spectrum in differential phase contrast and dark field contrast mammograms,” Med. Phys. 41(12), 120701 (2014). [CrossRef] [PubMed]

10. T. Thuering, W. Barber, Y. Seo, F. Alhassen, J. Iwanczyk, and M. Stampanoni, “Energy resolved X-ray grating interferometry,” Appl. Phys. Lett. 102(19), 191113 (2013). [CrossRef]

11. Y. Ge, K. Li, and G-H. Chen, “Cramér-Rao lower bound in differential phase contrast imaging and its application in the optimization of data acquisition systems,” Proc. SPIE 9033, 90330F (2014). [CrossRef]

12. M. Engelhardt, C. Kottler, O. Bunk, C. David, C. Schroer, Joachim Baumann, M. Schuster, and F. Pfeiffer, “The fractional Talbot effect in differential x-ray phase-contrast imaging for extended and polychromatic x-ray sources,” J. Microsc. 232(1), 145–157 (2008). [CrossRef] [PubMed]

13. T. Thuering and M. Stampanoni, “Performance and optimization of X-ray grating interferometry,” Phil. Trans. R. Soc. A 372(2010), 20130027 (2014). [CrossRef] [PubMed]

14. G. Pelzer, T. Weber, G. Anton, R. Ballabriga, F. Bayer, M. Campbell, T. Gabor, W. Haas, F. Horn, X. Llopart, and N. Michel, “Grating-based x-ray phase-contrast imaging with a multi energy-channel photon-counting pixel detector,” Opt. Express 21(22), 25677–25684 (2013). [CrossRef] [PubMed]

15. F. Pfeiffer, M. Bech, O. Bunk, P. Kraft, E.F. Eikenberry, CH. Brönnimann, C. Grünzweig, and C. David, “Hard-X-ray dark-field imaging using a grating interferometer,” Nat. Mater. 7(2), 134–137 (2008). [CrossRef] [PubMed]

16. V. Revol, C. Kottler, R. Kaufmann, U. Straumann, and C. Urban, “Noise analysis of grating-based x-ray differential phase contrast imaging,” Rev. Sci. Instrum. 81, 073709 (2010). [CrossRef] [PubMed]

17. T. Köhler, K. Engel, and E. Roessl, “Noise properties of grating-based x-ray phase contrast computed tomography,” Med. Phys. 38, S106–S116 (2011). [CrossRef] [PubMed]

18. R. Raupach and T.G. Flohr, “Analytical evaluation of the signal and noise propagation in x-ray differential phase-contrast computed tomography,” Phys. Med. Bio. 56, 2219–2244 (2011). [CrossRef]

19. G-H. Chen, J. Zambelli, K. Li, N. Bevins, and Z. Qi, “Scaling law for noise variance and spatial resolution in differential phase contrast computed tomography,” Med. Phys. 38, 584–588 (2011). [CrossRef] [PubMed]

20. T. Weber, P. Bartl, F. Bayer, J. Durst, W. Haas, T. Michel, A. Ritter, and G. Anton, “Noise in x-ray grating-based phase-contrast imaging,” Med. Phys. 38(7), 4133–4140 (2011). [CrossRef] [PubMed]

21. K. Engel, D. Geller, T. Köhler, G. Martens, S. Schusser, G. Vogtmeier, and E. Rössl, “Contrast-to-noise in X-ray differential phase contrast imaging,” Nucl. Instrum. Meth. A 648, S202–S207 (2011). [CrossRef]

22. R. Raupach and T.G. Flohr, “Performance evaluation of x-ray differential phase contrast computed tomography (PCT) with respect to medical imaging,” Med. Phys. 39(8), 4761–4774 (2012). [CrossRef] [PubMed]

23. G-H. Chen and Y. Li, “Synchronized multiartifact reduction with tomographic reconstruction (SMART-RECON): A statistical model based iterative image reconstruction method to eliminate limited-view artifacts and to mitigate the temporal-average artifacts in time-resolved CT,” Med. Phys. 42(8), 4698–4707 (2015). [CrossRef] [PubMed]

24. G. Strang, Introduction to linear algebra, (Wellesley-Cambridge Press, 1993).

25. R. Zhang, L. Zhang, Z. Chen, W. Peng, and R. Li, “Sensitivity of a non-interferometric grating-based x-ray imaging system,” Phys. Med. Bio. 59(7), 1573–1588 (2014). [CrossRef]

PMMA thickness	Method	Bin 1	Bin 2	Bin 3
1.5 cm	Original	0.59 ± 0.17	0.48 ± 0.18	0.21 ± 0.31
	SVD:J	0.57 ± 0.15	0.48 ± 0.12	0.25 ± 0.06
	SVD:J_P	0.57 ± 0.12	0.48 ±0.10	0.24 ± 0.05

3.5 cm	Original	0.48 ± 0.21	0.35 ± 0.21	0.17 ± 0.31
	SVD:J	0.45 ± 0.15	0.38 ± 0.13	0.20± 0.07
	SVD:J_P	0.46 ± 0.13	0.39 ± 0.11	0.20 ± 0.05

5.5 cm	Original	0.34 ± 0.27	0.31 ± 0.26	0.14 ± 0.33
	SVD:J	0.35 ± 0.20	0.29 ± 0.16	0.15 ± 0.09
	SVD:J_P	0.38 ± 0.16	0.32 ± 0.14	0.16 ± 0.07

Improving radiation dose efficiency of X-ray differential phase contrast imaging using an energy-resolving grating interferometer and a novel rank constraint

Abstract

1. Introduction