Detector-unit-dependent calibration for polychromatic projections of rock core CT

Abstract

Computed tomography (CT) plays an important role in digital rock analysis, which is a new prospective technique for oil and gas industry. But the artifacts in CT images will influence the accuracy of the digital rock model. In this study, we proposed and demonstrated a novel method to restore detector-unit-dependent functions for polychromatic projection calibration by scanning some simple shaped reference samples. As long as the attenuation coefficients of the reference samples are similar to the scanned object, the size or position is not needed to be exactly known. Both simulated and real data were used to verify the proposed method. The results showed that the new method reduced both beam hardening artifacts and ring artifacts effectively. Moreover, the method appeared to be quite robust.

Keywords

Computed tomography beam hardening artifacts ring artifacts polychromatic projection calibration

1 Introduction

Digital rock analysis is a new prospective technique for the oil and gas industry, where the creation of 3D digital rock model plays a foundation role [1]. CT, which can acquire 3D image of a rock non-destructively and effectively, becomes an important tool to establish digital rock model. However, the artifacts in CT images may affect the accuracy of the digital rock model. Therefore, artifact reduction of CT image is a cornerstone for digital rock analysis [2].

The artifacts frequently encountered in rock core CT images are ring artifacts and beam hardening artifacts, which are mainly caused by the inconsistent responses of detector units and the polychromatic feature of the X-ray. Besides, rock core CT images may also have so-called photon starvation artifacts [3]. Actually, a rock core typically is a cylindrical object that is thicker in the middle and thinner in the periphery. Limited by the dynamic range, the detector units corresponding to the middle part of the rock core have to suffer from photon starvation, so as to avoid count saturation for other detector units. To avoid photon starvation, attenuators are used to reduce the variation of beam intensity reached to the detector. However, such attenuators may cause the distortion of CT values [4].

Several methods are proposed to correct the artifacts mentioned above. The typical methods for ring artifact correction include flat field calibration methods [5, 6], detector shake methods [7,8, 7,8], filtering methods [9 –11], and etc. Generally, they can well correct the high frequency ring artifacts, but because of the local feature of the methods, they fail to deal with low frequency distortions and banding artifacts. Typical methods for beam hardening correction include linearization methods [12 –14], dual energy methods [15 –20], and iterative methods [21 –25]. The main idea of linearization methods and dual energy methods is to estimate a calibration function that converts the polychromatic projections into equivalent monochromatic projections. In iterative methods, the artifact correction is modeled as an optimization problem, and solved iteratively.

In this study, we propose a novel artifact correction method, which simultaneously eliminates both high and low frequency ring artifacts, beam hardening artifacts, as well as CT value distortion caused by attenuator. The idea of our method is similar to that of the linearization methods and the EDEC method [18] in some degree. But a striking difference is that we restore a set of detector-unit-dependent calibration functions that have taken both the detector unit response and the detector-unit-dependent attenuator thickness into consideration. The reference samples used to restore the calibration function are just simple objects, such as cylindrical samples. In addition, we do not need to know the sample’s shape, size or position, as long as their attenuation coefficients are similar to the scanned object. We verify our method using simulated data as well as real data of artificial and natural rock cores. Experimental results show that the artifacts are significantly reduced by our method.

2 Problem and model of rock core CT

In this section, we first depict the rock core CT imaging problem and build its mathematical model. A rock core CT scanner is composed of X-ray source, detector, mechanical system, attenuator, controller, and computer system. Figure 1 shows the schematic diagram of a rock core CT scanner. The CT scanning process can be modeled as

$\begin{matrix} I (β, u) = \int_{E_{min}}^{E_{max}} γ (E, u) \cdot I_{0} (E, u) \cdot exp (- μ_{0} (E) r (u)) \\ \cdot exp (- \int_{x \in L (u)} μ (x 1 pt R (β), E) dL) dE + σ (u) \end{matrix},$ (1)

Fig.1

Schematic diagram of the rock core CT.

where the symbols used in the formula are described in Table 1.

Table 1

Symbols used in formula (1)

Symbol	Definition
u	the coordinate of detector unit
L (u)	the X-ray path from the source spot to detector unit u
r (u)	the path length that X-ray passes through the attenuator corresponding to detector unit u
β	the rotation angle
R (β)	rotation matrix
E	_{the energy of the X - ray photon}
E_min/ - E_max	the minimum/maximum energy of the X-ray photons
x	the coordinate of sample scanned
μ (x, E)	the linear attenuation coefficient of the sample at the point x related to energy E
μ₀ (E)	the linear attenuation coefficient of the attenuator
γ (E, u)	the efficiency of detector unit u
I₀ (E, u)	the incident photons intensity emitted from X-ray source
I (β, u)	the photon intensity measured by detector unit u when the sample in rotation angle β
σ (u)	the intensity of the scattering photons

Let ${\bar{I}}_{0} (u)$ represent the photon intensity measured by detector unit u without the sample, that is

${\bar{I}}_{0} (u) = \int_{E_{min}}^{E_{max}} γ (E, u) \cdot I_{0} (E, u) \cdot exp (- μ_{0} (E) r (u)) dE + σ_{0} (u),$ (2)

where σ₀ (u) is the intensity of scattered photons. Then, the polychromatic projections can be written as

$p (u, β) = - log \frac{I (u, β)}{{\bar{I}}_{0} (u)} .$ (3)

However, it is well-known that the commonly used CT reconstruction algorithms, such as FBP or ART are based on monochromatic assumption. Therefore, direct application of these algorithms for image reconstruction from polychromatic projections p (u, β) will face a risk of importing some artifacts into the images, such as ring artifacts and beam hardening artifacts.

3 Calibration of polychromatic data

In this section, we first analyze the causes of ring artifacts and beam hardening artifacts, and then propose our calibration method.

3.1 Analysis on artifact causes

3.1.1 Ring artifacts

If the scattered photons are ignored, then we have projections for the case of monochromatic X-ray as $\begin{matrix} p (u, β) \\ = - log \frac{γ (E_{0}, u) \cdot I_{0} (E_{0}, u) \cdot exp (- μ_{0} (E_{0}) r (u)) \cdot exp (- A (u, β, E_{0}))}{γ (E_{0}, u) \cdot I_{0} (E_{0}, u) \cdot exp (- μ_{0} (E_{0}) r (u))} \\ = A (u, β, E_{0}) \end{matrix},$ (4) where A (u, β, E₀) = ∫_x∈L(u)μ (xR (β), E₀) dL. The γ (E₀, u) in numerator and denominator are canceled out, which means that the projection p (u, β) is independent of the detection efficiency γ (E₀, u). However, for the case of polychromatic X-ray, the related projection $p (u, β) = - log \frac{\int_{E_{min}}^{E_{max}} γ (E, u) \cdot I_{0} (E, u) \cdot exp (- μ_{0} (E) r (u)) \cdot exp (- A (u, β, E)) dE}{\int_{E_{min}}^{E_{max}} γ (E, u) \cdot I_{0} (E, u) \cdot exp (- μ_{0} (E) r (u)) dE}$ (5) always depends on γ (E, u). Therefore, the projection inconsistency for adjacent detector units is caused not only by the difference of their detection efficiency, but also by the X-ray polychromatic feature. The accumulation of such projection inconsistency between adjacent detector units in the reconstructed images forms the so-called ring artifacts.

3.1.2 Beam hardening artifacts

A polychromatic X-ray beam is composed of individual photons with a range of energies. As the beam passes through an object, it becomes “harder”, that is to say its mean energy increases, because the lower-energy photons are absorbed more rapidly than the higher-energy photons. For a uniform cylindrical object, beam hardening usually results in a cupping artifact. In practice, a bowtie shaped attenuator (or filters) is often employed to correct the beam hardening artifacts, and also to reduce the dynamic range requirement on the detector [4]. It is noticed from (3) that the essential role of an attenuator is to modify the emission spectra I₀ (E, u) of X-rays for individual detector units into the following effective spectra $I_{0} (E, u) \cdot exp (- μ_{0} (E) r (u)) .$ (6)

Some kinds of image distortions, such as “capping” or “Mexican hat” artifacts [26], may be caused by the bowtie attenuator. According to the analysis above, the X-ray polychromatic feature is the substantial reason causing the ring artifacts and the beam hardening artifacts. To eliminate the artifacts, the polychromatic projections should be well calibrated for each detector unit.

3.2 The calibration method

The essence of calibration is to convert the polychromatic projections into the so called “desired data” [14], which are usually the line integral of reference material. In clinical CT, water is commonly used as such reference material, as most of body tissues have similar attenuation coefficients with water. For rock core CT, we choose aluminum serving as the reference material, since the attenuation coefficient of aluminum is similar to those of the prominent compositions of nature rock core [27], in addition, aluminum is inexpensive and easy to process.

To convert the polychromatic projections to the line integrals, a map from polychromatic projection to line integral is usually built by scanning some reference samples with different thickness. Existing methods ignore the inconsistency of the detector response, and build one unified map, which usually works fine on beam hardening artifact correction for CT systems with no attenuator, but fails to correct ring artifacts and CT value distortion. In order to correct simultaneously all the artifacts mentioned above, in this paper, we propose a method to build detector-unit-dependent maps, i.e., build one map for each individual detector unit.

Suppose that the reference sample consists of homogeneous materials with attenuation coefficient μ (E). Then its polychromatic projection is $\bar{p} (t, u) \overset{def}{=} - log \frac{\int_{E_{min}}^{E_{max}} γ (E, u) \cdot I_{0} (E, u) \cdot exp (- μ_{0} (E) r (u)) \cdot exp (- μ (E) t) dE}{\int_{E_{min}}^{E_{max}} γ (E, u) \cdot I_{0} (E, u) \cdot exp (- μ_{0} (E) r (u)) dE},$ (7) where t = l (u, β) is the intersection length of the reference sample and the X-ray path L (u) at rotation angle β. In order to calibrate the polychromatic projections for each detector unit u, we need to estimate the inverse function $t = t (\bar{p}, u)$ , that converts the polychromatic projection $\bar{p}$ of reference sample to its line integral t. Generally, the inverse function is approximated by a linear combination of some basis functions ${t_{n} (\bar{p}) | n = 0, 1, \dots, N}$ as $t (\bar{p}, u) \overset{def}{=} \sum_{n = 0}^{N} c_{n} (u) 1 pt 1 pt \cdot t_{n} (\bar{p}),$ (8) where c_n (u), 1ptn = 0, 1, …, N are the detector-unit-dependent coefficients. The basis functions are usually chosen to be the monomial $t_{n} (\bar{p}) = {\bar{p}}^{n}$ , 1ptn = 0, 1, …, N. In the existing related papers [14 , 28], each coefficient c_n (u) is taken as a constant c_n, which implies that the inconsistences of detection efficiencies on detector units are ignored. Then, from the relation of (8), it is followed that the image f (x) of the reference sample can be expressed by a linear combination of series images f_n (x) reconstructed from $t_{n} (\bar{p})$ as $f (x) = R^{- 1} (t (\bar{p}, u)) = \sum_{n = 0}^{N} c_{n} 1 pt \cdot f_{n} (x),$ (9) where R^-1 denotes the FBP reconstruction operator, $f_{n} = R^{- 1} (1 pt {\bar{p}}^{n})$ . The set of coefficients {c_n|n = 0, 1, …, N } is determined by solving an optimization problem that minimizes the deviation between f (x) and a given reference image. Those methods can significantly reduce the beam hardening artifacts in the reconstructed images. However, the ignorance of the inconsistences of detection efficiencies on detector units may cause some kind artifacts, such as ring artifacts and streak artifacts. These image-domain-based methods as describe in [14] are difficult to estimate the detector-unit-dependent coefficient c_n (u).

In the following, we introduce our projection-domain-based method of estimating the detector-unit-dependent coefficients {c_n (u) |n = 0, 1, …, N } in formula (8). For a given u, if a set of sampling data $Ω_{u} = {(\bar{p} (t, u), t) 1 pt 1 pt | 1 pt 1 pt t = l (u, β), 1 pt 1 pt β = 1, 2, \dots, B}$ is obtained in some way, then these coefficients can be solved by the least squares method ${c_{n} (u) | 1 pt n = 0, 1, \dots, N} = arg min \sum_{β = 1}^{B} {(\sum_{n = 0}^{N} c_{n} (u) 1 pt 1 pt \cdot {(\bar{p} (l (u, β), u))}^{n} 1 pt - l (u, β) 1 pt 1 pt)}^{2} .$ (10)

Commonly the least squares method cannot get optimal result when the sampling data in Ω_u are insufficient. In order to improve the performance of the least squares method, further information is needed. According to the characteristics that $t (\bar{p}, u)$ is a monotonically increasing convex function on $\bar{p}$ , we have $\frac{\partial t (\bar{p}, u)}{\partial \bar{p}} \geq 0, 1 pt 1 pt \frac{\partial^{2} t (\bar{p}, u)}{\partial {\bar{p}}^{2}} \geq 0,$ or $\begin{matrix} \sum_{n = 0}^{N - 1} (n + 1) \cdot c_{n + 1} (u) 1 pt \cdot 1 pt {\bar{p}}^{n} \geq 0 \\ \sum_{n = 0}^{N - 2} (n + 2) \cdot (n + 1) c_{n + 2} (u) 1 pt 1 pt \cdot {\bar{p}}^{n} \geq 0, \forall 1 pt 1 pt \bar{p} \geq 0 \end{matrix},$ (11)

Adding the conditions (11) into the model, we obtain the following constrained optimization problem: $\begin{matrix} {c_{n} (u) | n = 0, 1, \dots, N} = arg min \sum_{β = 1}^{B} {(\sum_{n = 0}^{N} c_{n} (u) \cdot {(\bar{p} (l (u, β), u))}^{n} - l (u, β))}^{2} \\ s . t . \sum_{n = 0}^{N - 1} (n + 1) \cdot c_{n + 1} (u) \cdot {\bar{p}}^{n} \geq 0 \\ \sum_{n = 0}^{N - 2} (n + 2) \cdot (n + 1) c_{n + 2} (u) \cdot {\bar{p}}^{n} \geq 0, \forall 1 pt 1 pt \bar{p} \geq 0 \end{matrix}$ (12)

Which can be solved with classical optimization algorithms [29]. We will call Eq. 10 unconstrained fit model (UFM), and Eq. 12 constrained fit model (CFM) in the rest of the paper. Please note that each detector unit corresponds to an optimization problem Eq. 10 or Eq. 12, and should be solved individually.

How to acquire an efficient set of data ${(\bar{p} (t, u), t) | t = l (u, β), β = 1, 2, \dots, B}$ is very important for the optimization model. A simple method is to measure a series of stacked plate samples with known thicknesses and compositions. But in this method all the samples have to be sequentially exposed to the X-rays. Hence the measurement process is time consuming, potentially error prone, and also difficult to use in the gantry rotation CT system.

In this study, we used CT image-based method to acquire the data set as follows:

Place some simple shaped homogeneous samples, for instance, one cylindrical sample or two combined cylindrical samples, in CT platform.

Scan the samples to acquire the polychromatic projections.

Reconstruct the CT image from the polychromatic projections.

Calculate the intersection length of X-ray through the samples from their segmented CT image.

In order to acquire the data set with a range of intersection lengths as large as possible for each detector unit u, the position and diameter of samples on the CT platform should be appropriately chosen.

4 Experiment

In this section, both simulated data and real data are used to evaluate the performance of the proposed detector-unit-dependent calibration method. For simplify, all the experiments are involved in a fan beam CT. It is straightforward to extend the method to cone-beam CT.

4.1 Numerical simulations

In this subsection, we verify the accuracy of UFM and CFM with simulated data. The polychromatic projections of a reference sample (Fig. 3(a)) and a test sample (Fig. 5(a)) are simulated numerically according to Eq. (1), Eq. (2) and Eq. (3), where r (u) = 0, i.e., no attenuator is used. The X-ray spectrum (Fig. 2) is obtained by SpectrumGUI (see http://spectrumgui.sourceforge.net/) and the mass attenuation coefficient is obtained from NIST [30]. The Poisson noises corresponding to emission flux of 10⁷ photons per measurement are added to the simulated data.

Fig.2

Simulated X-ray spectrum.

Fig.3

(a) the reference sample. (b) the polychromatic projections of (a). (c) the variation range of the intersection length.

Fig.4

Comparison of the restored calibration functions corresponding to the four detector units: (a) u = 1; (b) u = 200;(c) u = 400; (d) u = 612.

Fig.5

Simulation results of the test sample: (a) The test sample; (b) The image reconstructed from the polychromatic data; (c) and (d) the images reconstructed from the data calibrated with UFM and CFM respectively. The display windows of (a), (c), and (d) are [0.9, 1.1], and that of (b) is [0.04, 0.06].

The parameters of the scanning configuration are set as follows: the focus-to-object distance (FOD) is 1000 mm; the focus-to-detector distance (FDD) is 1200 mm; the line detector is of 1224 units with pitch 0.083 mm; and the full sampling angle B is 1080. With the above parameters, the diameter of the field of view is about 130 mm.

The reference sample is a homogeneous aluminum cylinder with a diameter of 100 mm. In order to obtain a range of intersection lengths as large as possible for each u, the reference sample is positioned off-center in the FOV about 21.20 mm away, as shown in Fig. 3(a). Figure 3(b) shows the simulated polychromatic projections of the reference sample. Figure 3(c) shows the variation range of intersection length for each u (marked by dark points) where the horizontal axis denotes detector unit u, and the vertical axis denotes the intersection length corresponding to u.

As mentioned previously, the variation range of intersection length for each u can affect the accuracy of the restored calibration function. We select four typical detector units (u = 1, 200, 400, 612), whose positions are marked with the blue lines in Fig. 3(b) and (c), to illustrate the corresponding results. The variation ranges of the four detector units’ intersection lengths are [0, 39.32], [0, 86.07], [59.77, 99.54], and [93.26, 99.92], respectively. Compared to the total range [0, 100], the four sets are insufficient in a certain degree. The data of u = 1 lacks relatively big intersection length, while the data of u = 612 lacks relatively small intersection length.

Figure 4 shows the restored calibration functions (with degree N = 4) corresponding to the four selected detector units. In each subfigure, the top part displays the histograms of sampling data, the middle part displays the calibration functions restored with UFM and CFM, and the bottom part displays the error between the restored calibration functions and the theoretical ones. These results indicate that, compared with UFM, the calibration function restored with CFM have a good tolerance to the insufficient sampling data.

We designed a cylindrical test sample containing four cubes inside, as shown in Fig. 5(a), to verify the effect of the restored calibration functions. The diameter of the cylinder is 100 mm. The four cubes with density 0.00, 1.35, 2.57, and 2.67 g·cm^–3 from left to right are used to simulate the non-homogeneity of the rock core. The density of the remaining part is 2.70 g·cm^–3.

We simulate polychromatic projections of the test sample and convert them into the line integral data using the calibration functions restored with UFM and CFM respectively. Figure 5 (b)-(d) show the comparison results. The profiles of the central rows of the images are plotted inFig. 6.

Fig.6

Profiles of the central rows of the images shown in Fig. 5: (a) The profile of Fig. 5 (b); (b) The profiles of Fig. 5 (a), (c), and (d).

We can see from Fig. 5 to Fig. 6 that the calibration functions restored with CFM are much more accurate than that restored with UFM. In addition, new artifacts are shown in the CT image (Fig. 5 (c)) that is calibrated with UFM. We can see from Fig. 5(c) that the artifacts gather in two lines which are close to the center of the sample. The artifacts, we analyze, are caused by the relative big errors of the UFM along the two lines. When restoring the calibration functions, the detector units near the center of the detector are lack of relatively small intersection length, the errors are relatively big when the intersection lengths are small (see Fig. 4(d)). While for the test sample, because of the tubes with small attenuation, the intersection lengths along the two lines are just relatively small, and result in the artifacts.

4.2 Real data experiments

Several types of rock cores are scanned and calibrated with the proposed method, and achieves satisfactory results. We show here the results of three typical rock cores: a sedimentary rock core, a metamorphic rock core, and a road core made of asphalt and gravels, which are shown in Fig. 7. As a comparison, we have also implemented the wavelet-Fourier filter method [11], which is taken as a practical method for ring artifact reduction.

Fig.7

Photographs of the rock core samples: (a) A sedimentary rock; (b) A metamorphic rock; (c) A road core made of asphalt and gravels.

The rock cores are scanned by the CT system of our laboratory. The X-ray source is YXLON Y.TU-450-D05, which has a focal spot size 3.0 mm. The detector is a NTB linear detector, which has 3710 detector units with 0.083 mm unit pitch and 12-bit counter dynamic range. The X-ray tube is operated at 140 kV, current 10 mA, and the exposure time is 30 ms for each projection. The FOD and the FDD are 900 mm and 1196 mm, respectively. 1080 projections are scanned perrotation.

In order to reduce the requirement of the detector dynamic range, an attenuator, as shown in Fig. 8(a), is placed in front of the detector. Two homogeneous aluminum cylinders with diameters 80 mm and 40 mm, as shown in Fig. 8(b), are chosen as the reference samples. The bigger one is placed close to the center of the FOV, and the smaller one is placed near the bigger one to extend the variation range of intersection length.

Fig.8

The attenuator (a) and the reference samples (b) used for the experiments.

Figure 9(b) shows the sampling values and the restored calibration functions (with degree N = 4) corresponding to the three typical units marked with blue, black, and red lines in Fig. 9(a). It is seen from Fig. 9(b) that the curvatures of the three calibration functions are significantly different due to the effect of the attenuator. In fact, the mean energy of X-ray beam becomes higher as the attenuator thickness increases. Therefore, the curvature of the calibration function corresponding to the peripheral detector units is smaller than that of the central ones.

Fig.9

(a) the polychromatic projections of reference samples and (b) the sampling values and the restored calibration functions corresponding to the three typical units marked with blue, black, and red lines in (a).

The restored calibration functions are used to correct the polychromatic projections of the rock cores. Figure 10 to Fig. 12 show the CT images of the three rock cores, where subfigure (a) displays the image directly reconstructed from the polychromatic projections; subfigure (b) displays the image corrected by the wavelet-Fourier method (L = 3, wname = ‘db5’, and sigma = 2.5) [11]; subfigure (c) displays the image reconstructed from the image reconstructed from the projections calibrated by the proposed method.

Fig.10

CT images of the sedimentary rock core: (a) The image reconstructed from the polychromatic projections; (b) The image corrected by the wavelet-Fourier method (L = 3, wname = ‘db5’ and sigma = 2.5); (c) The image reconstructed from the projections calibrated by the proposed method. The display window of images (a)-(c) is [0, 0.09].

We can see from Fig. 10 to Fig. 12 that the wavelet-Fourier method reduces high frequency ring artifacts effectively, but introduces low frequency ring artifacts at the same time. While there are no obvious high or low frequency ring artifacts in the images reconstructed with the proposedmethod.

Fig.11

CT images of the metamorphic rock core: (a) The image reconstructed from the polychromatic projections; (b) The image corrected by the wavelet-Fourier method (L = 3, wname = ‘db5’ and sigma = 2.5); (c) The image reconstructed from the projections calibrated by the proposed method. The display window of images (a)-(c) is [0, 0.1].

Fig.12

CT images of the road core: (a) The image reconstructed from the polychromatic projections; (b) The image corrected by the wavelet-Fourier method (L = 3, wname = ‘db5’ and sigma = 2.5); (c) The image reconstructed from the projections calibrated by the proposed method. The display window of images (a)-(c) is [0, 0.1].

5 Conclusion

In this paper, we presented a simple and accurate method to restore the detector-unit-dependent calibration functions, which convert the polychromatic projections into the line integral data of the reference material. Experimental results show that both beam hardening artifacts and ring artifacts can be significantly reduced by our calibration method.

Footnotes

Acknowledgments

This work was supported by the Beijing Education Committee, and the Natural National Science Foundation of China [grant numbers 61401289, 61671311 and 61372179].

References

Pyrak-Nolte

L.J.

, Montemagno

C.D.

and Nolte

D.D.

, Volumetric imaging of aperture distributions in connected fracture networks, Geophysical Research Letters 24(18) (1997), 2343–2346.

Ketcham

R.A.

and Carlson

W.D.

, Acquisition, optimization and interpretation of X-ray computed tomographic imagery: Applications to the geosciences, Computers & Geosciences 27(4) (2001), 381–400.

Cnudde

and Boone

M.N.

, High-resolution X-ray computed tomography in geosciences: A review of the current technology and applications, Earth-Science Reviews 123(4) (2013), 1–17.

Kontson

and Jennings

R.J.

, Bowtie filters for dedicated breast CT: Theory and computational implementation, Medical Physics 42(3) (2015), 1453–1462.

Herman

G.T.

, Image reconstruction from projections, Real-Time Imaging 1(1) (1995), 3–18.

Chen

, Xi

, Cong

, Liu

, Wei

and Wang

, X-ray CT geometrical calibration via locally linear embedding, Journal of X-ray Science and Technology 24(2) (2016), 241–256.

Davis

G.R.

and Elliott

J.C.

, X-ray microtomography scanner using time-delay integration for elimination of ring artefacts in the reconstructed image, Nuclear Instruments and Methods in Physics Research Section A: Accelerators, spectrometers, Detectors and Associated Equipment 394(1) (1997), 157–162.

Zhu

, Zhao

, Li

and Zhang

, Micro-CT artifacts reduction based on detector random shifting and fast data inpainting, Medical Physics 40(10) (2013), 31114–31634.

Tang

, Ning

, Yu

and Conover

, Cone beam volume CT image artifacts caused by defective cells in x-ray flat panel imagers and the artifact removal using a wavelet-analysis-based algorithm, Medical Physics 28(5) (2001), 812–825.

10.

Boin

and Haibel

, Compensation of ring artefacts in synchrotron tomographic images, Optics Express 14(25) (2006), 12071–12075.

11.

Münch

, Trtik

, Marone

and Stampanoni

, Stripe and ring artifact removal with combined wavelet–Fourier filtering, Optics Express 17(10) (2009), 8567–8591.

12.

Herman

G.T.

, Correction for beam hardening in computed tomography, Physics in Medicine and Biology 24(1) (1979), 81–106.

13.

Hammersberg

and Mångård

, Correction forbeam hardening artefacts in computerised tomography, Journalof X-ray Science and Technology 8(1) (1998), 75–93.

14.

Kachelrieß

, Sourbelle

and Kalender

W.A.

, Empirical cupping correction: A first-order raw data precorrection for cone-beam computed tomography, Medical Physics 33(5) (2006), 1269–1274.

15.

De Man

, Nuyts

, Dupont

, Marchal

and Suetens

, An iterative maximum-likelihood polychromatic algorithm for CT, IEEE Transactions on Medical Imaging 20(10) (2001), 999–1008.

16.

Casteele

E.V.D.

, Dyck

D.V.

, Sijbers

and Raman

, An energy-based beam hardening model in tomography, Physics in Medicine and Biology 47(23) (2002), 4181–4190.

17.

Casteele

E.V.D.

, Dyck

D.V.

, Sijbers

and Raman

, A model-based correction method for beam hardening artefacts in X-ray microtomography, Journal of X-Ray Science and Technology 12(1) (2004), 43–57.

18.

Stenner

, Berkus

and Kachelrieß

, Empirical dual energy calibration (EDEC) for cone-beam computed tomography, Medical Physics 34(9) (2007), 3630–3641.

19.

Zhao

, Hu

, Zhao

and Zhang

, A novel iterative reconstruction method for dual-energy computed tomography based on polychromatic forward-projection calibration, Insight 56(10) (2014), 541–548.

20.

Zhao

, Zhao

and Zhang

, An extended algebraic reconstruction technique (E-ART) for dual spectral CT, IEEE Transactions on Medical Imaging 34(3) (2015), 761–768.

21.

Kyriakou

, Meyer

, Prell

and Kachelrieß

, Empirical beam hardening correction (EBHC) for CT, Medical Physics 37(10) (2010), 5179–5187.

22.

Krumm

, Kasperl

and Franz

, Reducing non-linear artifacts of multi-material objects in industrial 3D computed tomography, NDT and E International 41(4) (2008), 242–251.

23.

Van Gompel

, Van Slambrouck

, Defrise

, Batenburg

K.J.

, de Mey

, Sijbers

and Nuyts

, Iterative correction of beamhardening artifacts in CT, Medical Physics 38(S1) (2011), S36–S49.

24.

Maaß

, Meyer

and Kachelrieß

, Exact dual energy material decomposition from inconsistent rays (MDIR), Medical Physics 38(2) (2011), 691–700.

25.

Kano

and Koseki

, A new metal artifact reduction algorithm based on a deteriorated CT image, Journal of X-ray Science and Technology 24(6) (2016), 901–912.

26.

Lück

, Kolditz

, Hupfer

and Kalender

W.A.

, Effect of shaped filter design on dose and image quality in breast CT, Physics in Medicine and Biology 58(12) (2013), 4205–4223.

27.

Wedepohl

K.H.

, The composition of the continental crust, Geochimica Et Cosmochimica Acta 59(7) (1995), 1217–1232.

28.

Ritschl

, Bergner

, Fleischmann

and Kachelrieß

, Water calibration for CT scanners with tube voltage modulation, Physics in Medicine and Biology 55(14) (2010), 4107–4117.

29.

Waltz

R.A.

, Morales

J.L.

, Nocedal

and Orban

, An interior algorithm for nonlinear optimization that combines line search and trust region steps, Mathematical Programming 107(3) (2006), 391–408.

30.

Hubbell

J.H.

and Seltzer

S.M.

, NIST: X-Ray Mass Attenuation Coefficients, 1996, pp. 1–112.