1. Introduction

X-ray micro-computed tomography (μCT) is a very useful tool for imaging millimeter- to centimeter-sized objects with high spatial resolution. To provide good image quality, several parameters must be fulfilled. First, the imaging system must have the necessary high spatial resolution, which requires either an x-ray source with small emission spot or a high-resolution x-ray camera, and preferably both. Secondly, the object must provide sufficient contrast, which can be achieved either by the natural composition of the object or by staining techniques. Phase contrast is an attractive alternative, especially for soft biological tissues due to its improved contrast for elements with low atomic numbers. Thirdly, the signal-to-noise ratio must be high enough, which is inevitably limited by the total number of photons registered on the detector. Finally, the imaging system must not be disturbed by artifacts, i.e., errors degrading the final image, which arise in various ways.

One significant and classical image error in μCT is ring artifacts. They appear as concentric rings or half-rings around the center of rotation and are typically assumed to be due to uncorrected variations in the detector and the detector sensitivity. A proper flat-field correction should, thus, remove these artifacts but due to imperfections or non-linear behavior, such errors may still persist even after the flat-field correction. A vast range of correction methods to remove ring artifacts have been published, both by post-processing the sinograms [1–4], and by post-processing the reconstructed image [5, 6]. There are also some ideas that are based on moving the detector, in order to average out the measurement errors over a larger detector area [7, 8] although this leads to a loss in spatial resolution if the detector movements are not exactly known.

We propose that the μCT ring artifacts remaining after flat-field correction arise from the combined effect of beam hardening and high-spatial-frequency variations in the thickness of the scintillator detector, which is typically used in these systems. These high-resolution x-ray cameras normally consist of a thin scintillator screen, which is optically coupled to a CCD or CMOS detector, by either a fiber-optic plate or a system of lenses. Because the point-spread-function of the detector is largely determined by the thickness of the scintillator, the scintillator must be kept thin. Often gadolinium oxysulfide (Gd₂O₂S) is used as scintillator because of its high x-ray absorption and high yield of visible photons. However, Gd₂O₂S is a powder, making the thickness of the scintillator vary with grain size and distribution.

In the present paper, we first carefully determine the spatial variations in scintillator thickness by measuring the detector response for several different x-ray spectra. We then show that beam hardening induces local variations in the detector response which depend on the local scintillator thickness, leading to ring artifacts upon tomographic reconstruction. Finally, we present a correction method, which removes the induced ring artifacts by processing the projections.

The mechanism for ring-artifact formation proposed here is particularly important in μCT where the relative variations in scintillator thickness are larger than in the thicker detectors for higher-energy, lower-resolution CT. However, the basic mechanism of non-uniformity in the detector response due to beam hardening has been discussed before, albeit not in connection to variations in the scintillator thickness. Davidson et al. [9] characterized the energy-dependent response of a Medipix detector and flat-panel detectors have been corrected for large-scale errors such as bands and lines [10] and for isolated pixel clusters [11].

2. Measuring the scintillator thickness

To eliminate ring artifacts in μCT, the spatial variations in scintillator thickness must first be measured. We utilize the fact that the scintillator conversion and transfer yield (G), i.e. the number of visible photons transferred to the fiber-optic plate per incident x-ray, will depend on both scintillator thickness (T_s) and incident x-ray energy (E). This factor includes the x-ray absorption efficiency, the conversion efficiency into visible light photons and scattering of the visible light inside the scintillator. For this, the analytic expression presented by Swank [12] is used. Using a polychromatic source with spectral power S(E), filtered with aluminum with attenuation coefficient μ_f (E) and thickness T_f, the camera response will change as the spectrum hardens through the filter, and this change will depend on the scintillator thickness. Material coefficients for filters, object and scintillator were obtained from NIST [13]. We can expect a camera response R as

The images were acquired using a high-resolution x-ray CCD camera (VHR, Photonic Science Ltd, UK) with 9 μm pixels fiber-optically coupled to a Gd₂O₂S screen, specified by the manufacturer to have density 3333 mg/cm³ for the packed powder, and an area density of 5 mg/cm². The nominal thickness is thus 15 μm. The camera was 80 cm away from a liquid-metal-jet microfocus source (R6 [14], Excillum AB, Sweden), operated at 60 kV, 80 W. The spectrum was measured with an energy resolving silicon-drift detector (X123, Amptek Inc., MA). The spectrum was filtered with different thicknesses of aluminum foil, and the exposure time was adapted to achieve approximately the same number of counts in the camera, to minimize the influence of any non-linearities in the response in the CCD. Read-out noise and thermal noise were always negligible compared to photon noise. 80 images were averaged. The exposure times and filter thicknesses are shown in Table 1.

Table 1. Filter thicknesses and exposure times for the scintillator measurements. The total exposure time is 5.6 hours. This was repeated four times to assess the measurement precision.

View Table

For each pixel, Eq. (1) poses an inverse problem for the unknowns T_s and g. The inverse problem was solved by finding the values that gave the best least-square fit between calculated and experimental data for R(T_f, T_s, t)−R_d. To give an indication of the scintillator measurement precision, the measurement was repeated four times with independent sets of images.

Figure 1 shows the measured scintillator thickness, averaged over the four datasets, with an average thickness 15.2 μm and standard deviation 1.2 μm (7.6 % of average thickness) across the whole camera field of view. The appearance of the result could easily be mistaken for noise. We note that the method is very precise at measuring variations in scintillator thickness, but does not give very accurate results for the average thickness of the scintillator. For each pixel, the standard deviation was calculated in between the four measurements. The average standard deviations for all pixels was 0.7 μm (5 %). However, most of the error was due to erroneous measurement in the average thickness, not in the variations. By first normalizing the four measurements to the same average thickness, the same average standard deviation was 0.08 μm (0.5 %). The deviations between each measurement and the mean of all four measurements is demonstrated in Fig. 1(c). The distribution position corresponds to the measurement deviation in average thickness, whereas the distribution spread corresponds to the measurement deviation in thickness variations.

 

Fig. 1 The scintillator thickness measured with the proposed method. The scintillator thickness measurements with the proposed method. (a) Thickness map of the mean of all four measurements. For visibility of the small structures, we only show a smaller region (500×500 pixels).(b) Histograms of measured thickness, for the four different measurements. (c) Histograms of the difference between each measurement and the mean of all four measurements.

Download Full Size | PDF

We believe that the reason for the large variations in overall thickness are due to small variations in the flux from the source, of less than 1 %. Furthermore, we note that the scintillator response primarily relies on the absorption efficiency, which relies on the area density (mg/cm²) and not the thickness. Thus the value of the density will affect the measured absolute thickness. However, ring artifacts are introduced by the scintillator thickness variations, and not the absolute thickness.

3. Correction method for projection images

To correct tomography projections for the errors induced by the scintillator variations, we first introduce the detector response in a projection image as

The formation of a flat-field image is described similarly as

Equations (2) and (3) are combined to eliminate g and t, to give

4. Demonstrating the method in tomography

Tomography of a 20 mm cylindrical polyethylene terephthalate (PET) phantom with air filled holes was performed. Four reconstructions were made using the same tomography dataset. Firstly, using a flat-field correction without absorber. Secondly, using a 1800 μm Al absorber in the flat-field, which adapts the spectrum of the flat-field to better match the transmitted spectrum through the object. Thirdly and fourthly, by preprocessing the tomography projections with the proposed correction method considering the variations in scintillator thickness, assuming μ_obj to be aluminum and PET, respectively. The tomography was acquired as 2401 projections over 180°, with source-object distance of 100 cm, and source-detector distance of 120 cm, exposure time 15 s per projection, and using the same equipment and settings as above. After flat-fielding or making scintillator variation correction, the projections were phase retrieved using Paganin’s method [15], to remove slight phase contrast edge enhancement. The tomographic reconstruction was done using Octopus Reconstruction (Inside Matters NV, Belgium).

The results are depicted in Fig. 2. It is clear that the reconstruction from conventionally flat-fielded images in Fig. 2(a) suffers from severe ring artifacts. These ring artifacts are considerably reduced in Fig. 2(b) by using a filtered spectrum in the flat-fields, to better match the spectrum transmitted through the object. However, the thickness of the object varies across the images, and therefore the correction is not perfect. We also see that ring artifacts are introduced outside the object, where the spectrum matching is worse with the absorber in the flat-fields than without. To assess the importance of material selection in the proposed method, as explained below, we processed with materials having both erroneous and correct material properties, aluminum (Fig. 2(c)) and PET (Fig. 2(d)), respectively. It is immediately obvious that the method is capable of removing ring artifacts also for the wrong choice of material.

 

Fig. 2 A tomography of a PET phantom with air-filled holes, reconstructed with three different data processing methods. (a) Reconstructed after flat-field correction without filter. (b) Reconstructed after flat-field correction with 1800 μm Al absorber in the flat-fields. (c) Reconstructed with the proposed correction method and wrong material (Al). (d) Reconstructed with the proposed correction method and correct material (PET). Profiles are taken along the black lines, averaged over 30 pixel lines to reduce noise. The numbers are the standard deviations within the gray areas, relative to the contrast between air and PET. Scale bar is 5 mm. I–IV indicate the four domains, as described in the text.

Download Full Size | PDF

To quantify the errors, we generated a line profile radially from the center of rotation. To reduce the influence of noise and other artifacts than ring artifacts, the line profiles were averaged over 30 pixel lines, corresponding to the width of the black line in the figure. Ring artifacts are perpendicular to the extracted profile, and are therefore unaffected by the averaging. In the profile, standard deviations in four different domains were calculated: outside the phantom (I), in the outer part of the phantom (II), in a air-filled hole (III), and in the plastic close to the center of rotation (IV). The standard deviations were normalized to the contrast difference between the air and plastic. In domain I, the proposed method is equivalent to flat-fielding without filter, whereas flat-fielding with a filter introduces ring artifacts. In all the other three domains, the proposed method (Fig. 2(d)) clearly exhibits the least amount of ring artifacts. In domains III and IV, the proposed method reduced the standard deviations from 4.2 % and 9.1 % to 1.9 % and 3.5 % respectively. We note that the standard deviations also include noise, other types of artifacts and ring artifacts of other origin than scintillator variations. Thus, the actual reduction of ring artifacts from scintillator variations are consequently larger than indicated by the numbers. Furthermore, detection errors close to the center of rotation will generate rings with smaller circumference, making the total error spread over a smaller domain in the reconstructed image [16]. Rings closer to the center of rotation are thus intrinsically more difficult to correct for.

Beam hardening can introduce so called cupping artifacts [17], making a uniform material appear darker in the center. This is due to the fact that the total material attenuation coefficient is reduced as the beam hardens. In Fig. 2(a) and Fig. 2(b), the built-in beam hardening correction in the reconstruction software was applied, but still leaves slight cupping artifacts, visible as a brighter area at the outermost part of the object. The proposed correction method already incorporates beam hardening correction. In Fig. 2(d), beam hardening correction was not applied, without leaving any residual artifacts. Still, in Fig. 2(c), the wrong material selection generated inverse cupping artifacts, which were successfully removed by applying the beam hardening correction. Regarding the ring artifacts, comparing Fig. 2(c) and Fig. 2(d) shows that the method is surprisingly insensitive to the choice of material.

For the correction to work appropriately, the local spectrum at every point in the projection image must be estimated. The estimate must not be exact, but any estimate that is more accurate than just assuming the empty-beam spectrum will reduce ring artifacts. In Eq. (2), S(E) exp(−μ_obj(E)T_obj) is the spectrum attenuated through the object. Here, we assume that μ_obj (E) is known, i.e. the material composition of the object is known. To demonstrate the effect on ring artifacts from the selection of wrong material, we compare Fig. 2(c) and Fig. 2(d) to find that selecting wrong material gives only slightly more ring artifacts. Strictly, the method is only valid for a single-material object where the material is known. To test the validity also for multi-material objects and erroneous material selection, we simulated projection images where the full spectrum was registered in every pixel. The simulation results are shown in Fig. 3. The simulated object consists of mostly water, with aluminum inserts. The energy estimate using water shows a small underestimation of the mean energy where aluminum is present. Using the heavier elements aluminum and iron slightly overestimates the beam mean energy. Using beryllium underestimates the change in mean energy from beam hardening by around 30 %.

 

Fig. 3 A simulation to calculate the deviation between true spectrum and the estimated spectrum from camera intensity. (a) The 5 mm water phantom diameter water cylinder used in the simulation, with 0.2 mm and 1 mm diameter aluminum and air inserts. (b) The mean energies of the true spectrum in the simulation and the local spectrum estimation in the simulation using four different materials. The estimate is based on water absorption. The method gives slight errors from wrong material or where there is phase contrast. The deviations between the estimated mean energies and the true mean energy is shown in (c).

Download Full Size | PDF

5. Conclusions

The proposed method to remove ring artifacts is based on a true physical cause that generates ring artifacts. Many published methods for removing or reducing ring artifacts rely only on image post-processing. In such methods, there is always a risk that true features in the object are identified as artifacts and removed, or that false features are introduced. By intervening directly on the true physical cause of ring artifacts, this risk should be considerably smaller. Naturally, the method is limited to removing rings that originate from such scintillator variations. However, the method removes most of the ring artifacts in the phantom tomography example, so we conclude that scintillator variations are the dominant cause of ring artifacts in this imaging system. We note that applying the method does not hinder application of ring artifact removal methods based on image post-processing.

In order to implement the correction method, the scintillator variations must be measured. This measurement requires stable and uniform illumination from the source as well as qualitative knowledge about the spectrum. In our setup, this was easily obtained by setting the detector at large enough distance from the source. We note that the method provided very good precision in measuring variations in the scintillator area density, but not as precise for the average thickness.

When the scintillator variations have been measured, the only extra information required compared to a normal tomography, is the empty-beam spectrum and the object material. However, the spectrum does not have to be quantitative, nor has to be the same as the spectrum used for measuring the scintillator area density. The simulated study on different materials showed that also material choices very different from the true sample can provide energy estimates well within 1 keV, compared to the maximum beam hardening of around 10 keV. The exception is materials with very low atomic numbers. However, in most cases the object material is known much better than that (e.g. absorption in biomedical samples is mostly from oxygen and carbon, and in some parts calcium).

Reducing artifacts generally improves the image quality and perception. This can specifically be of importance when imaging slight variations in object density at high resolution, e.g. to differentiate soft tissue types [18] or structures within a soft tissue [19] in a biomedical sample.