## Abstract

This paper presents a mathematical characterization and analysis of beam-hardening artefacts in X-ray computed tomography (CT). In the field of dental and medical radiography, metal artefact reduction in CT is becoming increasingly important as artificial prostheses and metallic implants become more widespread in ageing populations. Metal artefacts are mainly caused by the beam-hardening of polychromatic X-ray photon beams, which causes mismatch between the actual sinogram data and the data model being the Radon transform of the unknown attenuation distribution in the CT reconstruction algorithm. We investigate the beam-hardening factor through a mathematical analysis of the discrepancy between the data and the Radon transform of the attenuation distribution at a fixed energy level. Separation of cupping artefacts from beam-hardening artefacts allows causes and effects of streaking artefacts to be analysed. Various computer simulations and experiments are performed to support our mathematical analysis.

## 1. Introduction

X-ray computed tomography (CT) is an important diagnostic tool for medical and dental tomographic imaging of the human body. Filtered backprojection (FBP) [1] is the most commonly used CT reconstruction algorithm. It assigns an X-ray attenuation coefficient to each pixel using the CT sinogram data. The FBP algorithm is based on the monochromaticity assumption of X-ray source so that sinogram data are an approximation of the Radon transform of an attenuation coefficient distribution at a fixed energy level.

This monochromatic assumption is acceptable for most soft tissues in the human body, and the FBP algorithm can be applied with current clinical CT to produce satisfactory reconstructions. However, metallic materials violate this assumption because their attenuation coefficients vary greatly with the energy level, and the presence of metal objects in the scan slice therefore leads to a serious mismatch between the sinogram data and the model data (i.e. the Radon transform of an image).

Because of the mismatch caused by this large variation in attenuation coefficient, CT scans of patients with metal implants (e.g. hip replacements, dental fillings, surgical clips and pacemaker wires) show streaking artefacts that result in a loss of anatomical information. Metal artefacts are caused not only by beam-hardening but also other physical effects such as scattered radiation, nonlinear partial volume and noise [2]. To suppress streaking artefacts, various metal artefact reduction (MAR) methods have been suggested in the last few decades. Most of these methods are based on inpainting-based methods (e.g. interpolation [3–7]), normalized interpolation methods [2], Poisson inpainting [8], wavelet [9–11], tissue-class models [12] and total variation [13], Euler’s elastica [14], iterative reconstruction methods [15–18] (e.g. iterative FBP, weighted least-square methods) and hybrid methods that combine the first two methods.

Despite various efforts to reduce metal artefacts, most existing methods for single-energy CT have been insufficiently effective to achieve widespread clinical use [2,12,19]. Dual-energy CT [20] provides satisfactory reconstruction via beam-hardening correction, but it requires a longer postprocessing time and a higher dose of radiation compared with single-energy CT [21].

This paper focuses on analysing beam-hardening artefacts, and aims to identify the causes and effects of streaking artefacts. We separate a *nonlinear beam-hardening factor* from the projection data. Beam-hardening artefacts can be divided into cupping artefacts and streaking artefacts. Cupping artefacts only disturb the image in the regions of problematic objects, whereas streaking artefacts corrupt the tomographic image outside the problematic region.

The nonlinear beam-hardening factor depends not only on the individual geometry of metallic objects but also on their arrangement in the field of view of CT. Based on mathematical analysis of the beam-hardening factor, we make several interesting configurations of metallic objects which do not produce the streaking artefacts although their beam-hardening effects are significant. We performed various numerical simulations and phantom experiments to support the theoretical analysis.

## 2. Analysis on beam-hardening artefact

Considering a two-dimensional parallel beam system, the projection data are given by the Lambert–Beer law [22,23]:
2.1
where *f*_{E}(**x**) is the attenuation coefficient of material at position and at energy level *E*, and is the Radon transform given by , , and *η*(*E*) represents fractional energy at photon energy *E* in the spectrum of the X-ray source [24,25]. Considering *η* as a probability density function with a mean beam energy *E*_{0} as shown in figure 1, we try to recover the *f*_{E0} from the projection data *P*(*φ*,*s*) generated by polychromatic X-ray source.

The standard CT reconstruction method is the FBP algorithm that is based on the assumption that for each *φ*. The CT image reconstructed by FBP is given by
2.2
where is the adjoint of the Radon transform given by , and is the Riesz potential given by .

### (a) Beam-hardening factor

In order to simply explain characteristics of the beam-hardening artefacts, we make some assumptions on the projection data. We assume that *η* is supported in the interval [*E*_{0}−*ε*,*E*_{0}+*ε*], i.e. and *f*_{E} is differentiable with respect to *E*. Since the projection data *P* in (2.1) are given by , the mismatch between the data *P*(*φ*,*s*) and can be estimated with the aid of ∂*f*_{E}/∂*E*:
2.3

Let *D*_{1},…,*D*_{j} be segments of high-attenuation materials with their derivatives of attenuation coefficients being on *D*_{j}. Assuming that ∂*f*_{E}/∂*E*|_{E=E0}=0 outside of ∪*D*_{j}, the distribution of the attenuation coefficient, *f*_{E}, can be expressed as
2.4
where *χ*_{Dj} denotes the characteristic function of *D*_{j}; *χ*_{Dj}=1 in *D*_{j} and 0 otherwise. With this expression, we have the following observation on a nonlinear beam-hardening factor.

### Observation 2.1

Assuming that *f*_{E} is expressed as (2.4), the projection data *P*(*φ*,*s*) can be expressed as
2.5
where represents a nonlinear beam-hardening factor depending on with and such that
2.6
and *δ*_{D,αE0}=0 outside the set .

### Proof.

For simplicity, we abbreviate by *α*_{j} and *δ*_{D,αE0} by *δ*_{D,α}. Under the assumption of (2.4), projection data *P*(*φ*,*s*) can be approximated by
2.7

Since , it follows from the generalized mean-value theorem for integration that for each , there exists *δ*_{D,α}(*φ*,*s*) in the interval [−*ε*,*ε*] such that
2.8
Using (2.8), the projection data *P*(*φ*,*s*) in (2.7) can be expressed as
2.9
▪

Let us briefly discuss the effects of *α*_{j} and *δ*_{D,α} in (2.5). In a monochromatic X-ray model with *δ*≈0, it follows from the estimate in (2.6) that *δ*_{D,α}≈0, and therefore . In high-energy CT, the attenuation coefficients of most materials do not vary much with *E* after 100 keV electron energy, and therefore we may assume *α*_{j}≈0. In these cases, we have and FBP method based on the analytic formula of provides artefact-free CT images with *f*_{CT}≈*f*_{E0}. Figure 2 illustrates numerical simulations of beam-hardening artefacts in CT images at different energy levels. The *f*_{CT} is computed when sinogram data are obtained at the three different energy ranges 20–40 keV (figure 2*b*), 40–60 keV (figure 2*c*) and 80–150 keV (figure 2*d*), respectively. Streaking artefacts in *f*_{CT} are significantly reduced in the energy range between 80 and 150 keV, because *α*_{j} in energy range of 80–150 keV is negligibly small. In this simulation, the projection data *P* in (2.1) were simulated based on the experimental findings of attenuation coefficient in [26] and source spectra described in [27].

Figure 3 shows nonlinear beam-hardening factor *δ*_{D,α} changing nonlinearly with (*φ*,*s*). We consider *δ*_{D,α} for the simplest case where *D* is a single disc. From (2.8),
2.10
Taking advantage of the radial symmetry of the disc *D*, we have
2.11
where is the radial function given by
2.12
We should note that *δ*_{D,α}(*φ*,**x**⋅** θ**)=

*δ*

_{D,α}(0,

*x*

_{1}) for all

*φ*∈(−

*π*,

*π*] and

**x**∈

*D*.

The following observation provides structure of the artefacts.

### Observation 2.2

Assume that *f*_{E} is expressed as (2.4) and *D*_{1},…,*D*_{N} are discs. Let *Φ* be the set of the projection angles such that
2.13
Here, *Φ* represents the set of projection angles such that the projection data over the set *Φ* are expressed by a Radon transform of a function. (Note that the streaking artefacts occur when *P* over the entire projection angles (−*π*,*π*] is not in the range space of the Radon transform.) We have the following:

Figure 4 shows the structure of in the case when *D* consists of two discs. To validate the observations (ii)-(a) and (ii)-(b) in observation 2.2, we compute the ratio of them with respect to or . In the case of two discs being iron in figure 4, we have
and
where .

Observation 2.2 can be applied to a more general geometry of *D*. For an example, figure 5 shows streaking artefacts caused by a square-shaped metal box. At angles of *φ*=0,*π*/2, the discrepancy between the projection data *P* and the Radon transform of the monochromatic image is large, and this mismatch generates the streaking artefacts *Γ* in *f*_{CT} along the lines in the directions *φ*=0,*π*/2. In this case, *Φ*^{c} in (2.13) could be a small neighbourhood of {0,*π*/2}.

### (b) Streaking artefact due to beam-hardening factor

Viewing the streaking artefacts as the propagation of singularities, microlocal analysis [28] can be used to describe the locations and orientations of singularities, simultaneously. Recently, Park *et al.* [29] have characterized the streaking artefacts owing to beam-hardening effects using the wavefront set [28]. The crucial observation in [29] is that metal streaking artefacts are produced only when the wavefront set of the Radon transform of does not contain the wavefront set of the square of the Radon transform. The metal streaking artefacts are generated by the severe nonlinear mismatch of the projection data *P* that depends on the geometry of the high attenuation materials .

For mathematical characterization of streaking artefacts, let us confine the streaking artefacts to wavefront set of straight line satisfying *Σ*_{x}(*f*_{CT})≠∅, for all **x**∈*L*_{φ,s}, where *Σ*(*v*) is the smallest closed conic subset of outside of which decays rapidly, and *Σ*_{x}(*u*) is a closed conic subset in defined as
According to their analysis of streaking artefacts [29], the beam-hardening factor *δ*_{D,α} produces the streaking artefacts along the tangent line to ∂*D*, the boundary of region *D*, provided that
2.14
where WF(*u*) denotes the wavefront set of *u*, which is a closed conic subset in defined as
In the case when *η*(*E*)=1/2*ε* and , Park *et al.* [29] proved that the streaking artefacts are straight lines *L*_{φ,s} tangent to the boundaries with touching at least two different points in the boundaries. In this case, the beam-hardening factor *δ*_{D,α} is determined by
2.15

According to this observation, if the region *D* is strictly convex or it is an annulus, there are no tangent lines touching two different points on the boundary of *D*, so that the CT image *f*_{CT} does not have the streaking artefacts as shown in figures 4 (top) and 6*b*. On the other hand, in the case of figures 4 (bottom) and 6*d*, there are four tangent lines touching two different points on the boundary of *D*, along which streaking artefacts are produced in the reconstructed CT image by FBP.

According to the observations 2.1 and 2.2, the streaking artefacts occur by way of minimizing the mismatch between projection and a Radon transform of monochromatic image. With this wavefront set point of view, the streaking artefacts are produced to compensate for the abrupt mismatch between projection data and Radon transform of *f*_{E0}.

## 3. Experiments and suggestions

Being used for imaging subjects containing high attenuation materials, a dual detector CT (Discovery CT750 HD, GE Healthcare, USA) is used to see the influence of beam-hardening artefacts. For a given polychromatic X-ray source, the absolute value of ∂*f*_{E}/∂*E*|_{E0} at an effective energy *E*_{0} is the major factor of beam hardening. The main difficulty in reducing beam-hardening artefacts comes from the fact that the discrepancy between the projection *P* and *f*_{E} is related to the geometry of *D*. It is because the intensities of the projection data do not change linearly with the thickness of *D*. Figure 7 illustrates virtual monochromatic images using dual energy CT. The beam-hardening artefacts are seen in figure 7*b*, although the monochromatic image, in principle, should not have the beam-hardening artefacts [30]. According to the wavefront set-based characterizations of streaking artefacts in §2b, singular lines representing streaking artefacts are included in (all intersections of tangent spaces for two different points in ∂*D*). These descriptions are observed from various experiments, even though they may contain other effects such as scattering, noise or photon starvation effects.

We can reduce metal artefacts significantly by acquisition of tilted CT reconstructions from different angles of the same object. The use of titled CT may provide a transverse CT slice containing less metallic objects, and therefore streaking artefacts are significantly reduced. This is because metal streaking artefacts are included in (all intersections of tangent spaces for two different points in ∂*D*). As shown in figure 8*c*, streaking artefacts occur along the tangent line to boundary of two metallic objects. By contrast, in the case of single metallic object, the streaking artefacts are rarely seen in the reconstructed image as shown in figure 8*b*. Indeed, Ballhausen *et al.* [31] have proposed a metal artefact reduction method based on similar observations about streaking artefacts obtained through the experiments. Since the streaking artefacts are highly dependent on the geometry of the objects, they can be reduced by combining the CT images of the same objects which are obtained from different angles.

## 4. Discussion and conclusion

High-attenuation object-related artefacts seriously degrade CT images, resulting in a loss of anatomical information for medical diagnosis. Therefore, patients with metal implants (e.g. hip replacements, dental fillings, surgical clips, and pacemaker wires) may not benefit from CT scanning. A major cause of metal artefacts is beam hardening, which causes the projection data *P*(*φ*,*s*) to deviate from the range of the Radon transform, rendering the discrepancy between *P* and *f*_{E0}. Other causes of metal artefacts include scatter, nonlinear partial volume effects, Poisson noise, motion and photon starvation.

This paper addresses mainly mathematical theories on beam-hardening artefacts so that the theories can guide metal artefact reduction. Mathematical characterization of metal streaking artefacts and the geometries of the high attenuation objects that cause them can be used to develop an efficient MAR method with reasonable computational cost. Although there are various iterative MAR methods using data fitting algorithms, massive computational loads that they require make them unacceptable for clinical applications at present.

The development of artefact reduction methods continues to pose difficulties because the inverse problem is highly nonlinear, depending on the geometry of high attenuation objects. However, knowledge of the locations and shapes of implants within the patient can assist in reducing metal streak artefacts. Incorporating prior knowledge of metallic implant shapes could aid in the development of a novel artefact correction method based on the mathematical analysis reported here. Our future research will aim to develop a prior-knowledge-based MAR method for use in the field of diagnosis, pre-operative and presurgical assessment and surgical navigation.

## Funding statement

This work was supported by Samsung Science and Technology Foundation (no. SSTF-BA1402-01). Y.E.C. and J.K.S. were supported (in part) by the Yonsei University Future-leading Research Initiative of 2014.

## Footnotes

One contribution of 11 to a theme issue ‘X-ray tomographic reconstruction for materials science’.

- Accepted January 27, 2015.

- © 2015 The Author(s) Published by the Royal Society. All rights reserved.