# An Analytical Method for Reducing Metal Artifacts in X-Ray CT Images

Authors: Ming Chen; Dimeng Xia; Dan Wang; Jingqi Han; Zhu Liu
Journal: Mathematical Problems in Engineering
 (2019)

---

## Abstract

Medical CT imaging often encounters metallic implants or some metal interventional therapy apparatus. These metallic objects can produce metal artifacts in reconstruction images, which severely degrade image quality. In this paper, we analyze the difference between polychromatic projection data and Radon transform data and develop an analytical method to reduce metal artifacts. Approximate features of metal artifacts can be obtained by a simplified energy spectrum function of x-ray beam. The developed method can reduce most artifacts, and preserve more original details. It does not require prior knowledge of x-ray energy spectrum and original projection data, avoiding iterative calculation and saving reconstruction time. Simulation experimental results show that the method can greatly remove metal artifacts.


---

## Body

## 1. Introduction

X-ray CT imaging plays an indispensable role in clinical diagnosis and interventional therapy. However, due to beam-hardening caused by some metallic implants in patients or some metallic interventional therapy apparatus, such as dental fillings, orthopaedic implants, or microwave ablation needle, there are strong streak or star-shape artifacts in reconstruction images [1–3], which are usually called metal artifacts. They seriously degrade image quality and bring trouble to clinical applications. Although many metal artifact reduction (MAR) methods were developed, their applications in clinical settings are not totally successful because of the complexity of their forming cause and characteristics. Currently there is no standard solution [4, 5]. Therefore, how to reduce metal artifacts still remain a challenging problem in x-ray medical CT imaging.The effects of x-ray beam-hardening, the photon starvation, and the partial volume can all result in metal artifacts [6]. In the past few decades, a large number of MAR algorithms have been proposed to correct or reduce metal artifacts. Interpolation was widely used for data completion [7–10], where missing projection data are approximated by an interpolation technique. But due to the inaccuracy of data interpolation, additional fringe artifacts and other deformations were introduced in a new reconstruction image [10, 11]. The missing data became more accurate by use of forward projections of a prior image [12, 13]. A combination method of normalization and interpolation was proposed to remove most of the artifacts [11]. However, since artifacts are very strong for some cases, some pixels are often classified into wrong types, leading to unsatisfactory results. In order to correct beam-hardening, some iterative algorithms were proposed, which reconstructed images from some processed projections [14–16]. They can suppress some artifacts, but there is still no satisfactory result for all images. Recently, there were some researches about the deep learning strategy to reduce metal artifacts [17, 18]. One of their drawbacks is that there is no common CT image database for model training, and another is that some mild artifacts typically still remain.The forming cause of metal artifacts is mainly the high attenuation of metallic objects, which leads to x-ray beam-hardening and aggravates the scattering phenomenon for the polychromatic x-ray beam spectrum. For low-atom number metals, satisfactory results were achieved by correcting beam-hardening [19–25]. Some dual-energy correction [21] and statistical iterative correction [23, 24] were proposed to reduce beam-hardening effects. The former required longer postprocessing and higher doses of radiation, while the latter needed more prior information about the energy spectrum of the incident x-ray and the energy-dependent attenuation coefficient of the materials. Park et al. [25] put forward a MAR method by giving the approximate expression of metal artifacts, where choosing the approximate alternative of energy spectrum is very critical to obtain more accurate geometrical characterizations of artifacts.For a polychromatic x-ray spectrum CT imaging, this paper analyzes the geometrical features of metal artifacts in a two-dimensional fan-beam system, constructs an approximate energy spectrum function, and gives an approximate expression of metal artifacts by using an excess photon energy relation. The artifacts can be analytically expressed as the approximate mathematical equation in this method.The rest of this paper is organized as follows. In Section2, we summarize some background knowledge about Radon transform and the projection expression for a polychromatic x-ray spectrum. In Section 3, an analytical MAR method is developed by looking for geometrical characterizations of metal artifacts. In Section 4, we give numerical simulations to verify the effectiveness of the developed method. Finally, conclusions will be made and some related issues will be discussed.


## 2. Background Knowledge

In this section, we give Radon transform of a two-dimensional (2D) function and projection expressions for a polychromatic x-ray spectrum. This section will give theoretical basis for finding the relations between single-energy projection data and polychromatic projection data in next section.Mathematically, 2D Radon transforms can be seen as a line integral process from one 2D function to another. LetR denote Radon transform. For fX and a line L, R can be expressed as follows:(1)RfL=∫LfXdl,or(2)Rfr,φ=∫Lr,φfx,ydl,where X=x,y is a reconstructed point, L=Lr,φ:X·∅=r,∅=cos⁡φ,sin⁡φ,dl is a arc length on L, φ is a projection angle, and r is a sampling variable along a detector direction. In fact, (2) is also called ray sum, line integral, or the projection under single energy [26].Whenφ is fixed and r is (-∞,+∞) in (2), Rf(r,φ) is a set of parallel projections shown in Figure 1. For all φ∈[0,π), Rf(r,φ) is usually called parallel-beam projections in CT imaging. The expression of a fan-beam projection for single energy is similar to (2) by converting the parameters of parallel-beam scan mode to those of fan-beam scan mode.Figure 1
Radon transform off(x,y) under φ.For CT system with a polychromatic x-ray spectrum, the projection expression above will become difficult. LetI0,E and IE denote the incident x-ray intensity and the transmitted intensity at an energy E, respectively. Let fX,E denote the linear attenuation coefficient of the detected object at E. To emphasize an energy-dependent nature of material attenuation, Lambert-Beer law [27, 28] can be written as follows:(3)IE=I0,Eexp⁡-∫LfX,Edl.X-ray beam produced by x-ray tube operating at 120 KV generally covers a broad spectrum, which is with two sharp peaks in an example shown in Figure2 by a solid line. The spectrum distribution shows that the range of the output x-ray photon energy is [20 keV, 120 keV], where energy photons less than 20 keV are absorbed by some materials in x-ray tube [29].Figure 2
An example of x-ray energy spectrum when x-ray tube operates at 120 KV, where a horizontal axis is x-ray energy, and a longitudinal axis is a normalized output.Projection data for a polychromatic x-ray source can be expressed as follows:(4)pr,φ=-ln∫IEdE∫I0,EdE=-ln⁡∫ηEexp⁡-∫LfX,EdldE,where ηE=I0,E/∫I0,E1dE1 represents a normalized x-ray beam energy spectrum [25, 30]. When metallic implants in patients or some metallic interventional therapy apparatus exist in some CT scanned slices, pr,φ will include beam-hardening effects. Thus, projection data are not accurate for the scanned slices and not equal to the ideal projections Rfr,φ above.


## 3. An Analytical MAR method

In this section, we first analyze the forming cause of metal artifacts based on aforementioned knowledge, then construct the approximate function of an energy spectrum distribution, and finally give geometrical features of metal artifacts in 2D fan-beam system.In CT scan system with a polychromatic x-ray source, beam-hardening leads to the significant difference between raw projectionpr,φ and Radon transform Rf(r,φ). If we only consider metal artifacts caused by beam-hardening, the difference expression pr,φ-Rf(r,φ) corresponds to projections of metal artifacts.We remark the aforementioned functionfX as fE′(X) at fixed energy E′ for a polychromatic x-ray beam. Then, fE′(X) is the corrected objective function. Let R-1 denote some exact CT reconstruction algorithm, such as filtering-back-projection (FBP) algorithm or backprojection-filtering (BPF) algorithm. So, R-1pX is a reconstructed image with metal artifacts, and R-1[p-RfE′](X) seen as metal artifacts is a reconstruction result from a difference part between projection data and Radon transform data. Obviously R-1pX can be expressed as follows: (5)R-1pX=fE′X+R-1p-RfE′X.It is assumed that metallic objects have the same attenuation coefficient in a reconstruction image, which is reasonable for most cases in clinic applications. In this case, metal artifacts mainly depend on some geometrical features of metal regions, a normalized energy spectrum function, and metal attenuations about x-ray. LetfEX denote a reconstruction image with metal artifacts R-1pX; that is, pr,φ=RfEX, where E is a variable. According to (4) and (5), R-1p-RfE′X can be written as follows: (6)R-1p-RfE′X=R-1-ln⁡∫ηEexp⁡-RfE-fE′r,φdE.The water and bone beam-hardening corrections can usually be used in current medical CT imaging systems [30]. In this paper, we only consider beam-hardening effects of metallic objects, since their attenuation coefficients have more dependency on x-ray energy than soft tissues and bones. Let Ω denote a metal domain, which can include some small metal regions. Then, fE-fE′=0 outside of Ω. Since the probability of photoelectric interactions is roughly inversely proportional to the cubic of excess photon energy [31], fE-fE′ in Ω can be approximated as(7)fE-fE′≈wEE-E′,where w(E) depends on E and attenuation coefficients of metal regions.So, according to (7), (6) can be rewritten as(8)R-1p-RfE′X≈R-1-ln⁡∫ηEexp⁡-wEE-E′RλΩr,φdE,where λΩ is a characteristic function for the metal domain Ω.Because the normalized energy spectrum functionηE is unknown, R-1[p-RfE′](X) can not be accurately obtained. We choose a following alternative function in order to give a specific and relatively accurate expression of ηE:(9)η′E=Asin⁡ωE-E′+θifE-E′≤h0,otherwisewhere h is a variable. It is worth noting that the selections of an amplitude A, a frequency ω, and a phase θ are generally related to h. An example η′E is shown in Figure 2 by a dotted line, where it is the smoother one. In fact, by choosing h, we can make an approximate spectrum distribution η′E satisfy the following relationship:(10)ln⁡∫ηE-η′Eexp⁡-wEE-E′RλΩr,φdE≈1.This means that we can use ln(∫η′Eexp{-w(E)(E-E′)RλΩr,φ}dE) to approximately replace ln⁡(∫ηEexp{-w(E)(E-E′)RλΩr,φ}dE). So, (8) can be further approximated as(11)R-1p-RfE′X≈R-1-ln∫η′Eexp⁡-wEE-E′RλΩr,φdE.Through a series of mathematical derivations on the right side in (11) and supposing that w(E) is approximately independent of E, R-1[p-RfE′](X) is finally deduced as follows: (12)R-1-ln⁡Aexp-whR′wR′2+ω2sin⁡γ-ωh+ωφ-Aexp-whR′wR′2+ω2sin⁡γ--ωh+ωφ,where R′=RλΩr,φ,cos⁡γ=wR′/wR′2+ω2. According to (5), the corrected objective function fE′(X) can be obtained by the following equation: (13)fE′X=R-1pX-R-1p-RfE′X.Equation (13) is the analytical MAR algorithm expression given by approximate geometrical features of metal artifacts in (12).


## 4. Implementation of MAR Method and Numerical Simulation Experiments

Based on the aforementioned derivation of MAR method, implementations of the developed MAR technique mainly includes three steps. Firstly, raw CT images are reconstructed using FBP algorithm from the polychromatic projection data. Then, metallic objects are segmented from the images to obtain a metal domain and its characteristic function, and an approximate artifact image can be reconstructed using FBP algorithm. Finally, the corrected objective image is given by subtracting the artifact image from the raw image. Its specific implementation steps are shown in Figure3.Figure 3
Implementation flow chart of the developed MAR method.To verify the effectiveness of the developed MAR method, numerical simulations are performed on VC++ 6.0 platform, where all steps in the method are coded in C language. Fan-beam CT system parameters are as follows: the distance between x-ray source and a center of rotation is 1100mm, the distance between x-ray source and a linear-array detector is 1600mm, and the length of a detector cell is 0.87mm. Projection data are generated by using the circular trajectory equidistant fan-beam CT scan mode, where the number of detector cells is 512 and the number of 360 degree full-scan projection angles is 660. The size of image matrix is 512 × 512.We select a modified digital mandibular phantom, as shown in Figure4. Its structure definitions are shown in Table 1, which are firstly given by Lemmens et al. [32]. Three circles with maximum gray values mean the implanted mercury objects. Bone components and soft tissues in the phantom are set according to the parameters reported by ICRU44 [33]. Their attenuation coefficients are obtained by XCOM software [34]. A polychromatic x-ray beam is generated by simulations, where x-ray tube operates at 120 kVp. The polychromatic projection data are obtained shown in Figure 5.Table 1
Definitions of the digital mandibular phantom structures.
Structure sequence Number Material Dense(g/cm3) 1 Mercury 13.6 2 Bone 2.9 3 Bone 2.7 4 Soft tissue 0.8 5 Soft tissue 0.6 6 Soft tissue 0.2Figure 4
The digital mandibular phantom, where the display window is[0.0,4.5].Figure 5
The projection data with beam-hardening effects.According to the implementation flowchart in Figure3, we select the threshold value 3000HU [35, 36] to segment the metallic parts and then calculate the projection data of the characteristic function. According to (10), the parameters in the approximate energy spectrum function η′E are designed as A=π/4h,ω=π/2h,E′=65Kev,h=55Kev, θ=h and the optimized parameter w=0.042. According to (12), an artifact image is reconstructed as shown in Figure 6, where geometrical features of metal artifacts are clearly observed. The objective image is obtained by subtracting the artifact image, as shown in Figure 7(a), where CT image reconstruction takes about 85 seconds using the proposed MAR method on a computer with CPU 3.40GHz. For comparison, we give a reconstruction result using the linear interpolation MAR in Figure 7(b). We can find that the developed methods can reduce most artifacts, and preserve more original details. In facts, if the materials with lower density than the mercury are implanted in the phantom, the better corrected results can be obtained.Figure 6
The reconstructed artifact image using (12).Figure 7
The corrected images: (a) is a reconstruction result using the developed analytical MAR method in (13), and (b) is a reconstruction result using linear interpolation MAR method, where the display window is [0.0,0.6] for displaying the details clearly.
(a) (b)


## 5. Conclusion and Discussion

In this paper, for CT scan system with a polychromatic x-ray source, we propose the analytical method for reducing metal artifacts. This method can directly obtain the corrected image by reconstructing features of metal artifacts without any prior information with ray energy spectrum, largely remove structural features of the artifacts, and better retain the original information. So, the developed method can obtain high-resolution CT images and provide more accurate information for clinical applications.A good alternative expression of the probability of the photoelectric interaction is very important in the developed method. However, due to random noises from x-ray beam generation process and detection process of linear-array detectors in an actual CT system, photoelectric interactions are very difficult. If more knowledge in random system analyzation could be applied in the photoelectric interaction, the better expression may be designed. How to express the photoelectric interactions using the random process will be a challenging and meaningful work.



---

*Source: 2351878-2019-01-13.xml*