## research papers

## Suppression of ring artefacts when tomographing anisotropically attenuating samples

**Sofya Titarenko,**

^{a}Valeriy Titarenko,^{b}^{*}Albrecht Kyrieleis,^{b}Philip J. Withers^{b}and Francesco De Carlo^{c}^{a}Department of Mathematics, Faculty of Physics, Moscow State University, Leninskie Gory, Moscow 119991, Russia, ^{b}Henry Moseley X-ray Imaging Facility, School of Materials, University of Manchester, Grosvenor Street, Manchester M1 7HS, England, and ^{c}X-ray Science Division, Argonne National Laboratory, 9700 South Cass Avenue, Bldg 401, Argonne, IL 60439-4837, USA^{*}Correspondence e-mail: valeriy.titarenko@manchester.ac.uk

There are many objects for which the attenuation varies significantly as they are rotated during computerized X-ray tomography, for example plate samples. This can lead to significant ring artefacts in the subsequent tomographic reconstructions. In this paper a new method is presented that can successfully suppress such ring artefacts and is applicable to both parallel and cone-beam geometries. Rapid correction is achieved *via* an analytical formula which involves only a matrix-vector multiplication, for which the matrix is known and depends on a regularization parameter. The efficacy of the method is demonstrated for a paleontological sample (calcified shark cartilage) and a carbon–carbon composite/Ti–SiC metal matrix composite test sample.

Keywords: X-ray absorption computerized tomography; laterally extended objects; ring artefacts; filtered back-projection; inverse problem; laminography.

### 1. Introduction

There are many objects which attenuate X-rays very differently in different directions. For example, this may be because geometrically they have a large aspect ratio, such as a plate, and so attenuate significantly more when viewed along the major axis, or because they comprise differently attenuating phases that are not homogeneously distributed. As a result when acquiring a series of radiographs (projections) as the sample is rotated in order to make a tomograph, the details in some projections, or regions within certain projections, will be only faintly represented, if at all. This can lead to significant ring artefacts in the subsequent reconstructions which vary in strength around the rings. In this paper a method is presented that can successfully suppress such ring artefacts.

Let us consider the following typical experimental set-up for X-ray tomography on a synchrotron X-ray beamline or laboratory X-ray imager. There is a source of X-rays which provides a beam of given geometry. The geometry of the beam may be arbitrary and is not restricted to classical parallel or cone-beam geometries. The only assumption is that the intensity of the beam is unchanged during the scan. The beam passes through a sample that is rotated about a vertical axis before falling onto a scintillator which absorbs the photon energy and re-emits it as visible light in turn. The light then passes through an optical system and is recorded by a CCD (charged couple device) camera. As a result a series of two-dimensional radiographs (projections) is obtained.

Given a set of projections one can apply various reconstruction algorithms in order to find the structure of the sample. For example, an FBP (filtered back-projection) algorithm for parallel geometry of the beam and an FDK (abbreviated from Feldkamp, Davis, Kress) algorithm for cone-beam geometry are currently the most popular algorithms (Natterer & Wübbeling, 2007; Kak & Slaney, 1988; Buzug, 2008). Let there be a monochromatic X-ray beam, having a beam path, *L*, through the sample, and *I*_{0} be the intensity of the X-ray photons before the sample. Then the intensity of the X-rays after the sample can be written as *I* = *I*_{0} ·*A* = *I*_{0}exp(-*p*), where the attenuation factor *p* is the line integral along *L* and is the at a point . This means that *p* is the path through the sample along the beam trajectory. Let *x* and *y* denote the horizontal and vertical coordinates, which are parallel to the rows and columns of the pixels on the CCD sensor, respectively, and *t* be a moment of time. Then the reconstruction algorithms use the attenuation factor *p*(*x*,*y*,*t*) in order to find the at a point .

In addition to electrons generated by the absorption of visible photons, the CCD sensor is also sensitive to electrons (so-called dark-noise electrons) generated by thermal processes within it. Therefore, to a first approximation the standard flat-field correction procedure can be defined (see, for example, Stock, 2008),

where *W* is a white (or flat) field (recorded with the sample removed from the beam) and *D* is a dark field (recorded with the beam switched off). Unfortunately, for various reasons mentioned below the measured attenuation factor differs from the real one *p*, so we can always write

As a first approximation the term *q* often does not depend on time, *i.e.* . This case corresponds, for example, to the case where there are some particles or pieces of dust/dirt on the surface of the scintillator which absorb X-rays. However, some other factors can lead to *q* being dependent on *p* or *t*. For example, the beam may not be monochromatic. Let the beam have an intensity dependent on the wavelength , and *p* and *q* be also functions of . The measured attenuation can be written as

and in the general case it is a function of *p*^{*}, *q*^{*} and , where

This means that , *i.e.* *q*^{*} is an implicit function of , *p*^{*} and , and is only the first term in the Taylor series of *q*^{*} for a given . Among the other causes of *q* being a function of *p* and *t* there are the following: the scintillator may have various impurities inside (the attenuation properties of these impurities/defects may be temperature dependent, so that intensity variation on passing through the sample as it is rotated may lead to temperature variations on the scintillator); scintillator performance may be non-linear; some pixels of the CCD sensor may have a non-linear response (see Bardelli *et al.*, 2006; Moszyński, 2006; Krumm *et al.*, 2008). In such cases the intensity of ring artefacts will depend on the shape of the sample: a ring artefact is more intense if the attenuation factor *p* has a greater value because the intensity of X-rays that pass through the sample is an exponential function of the attenuation factor *p* and the same errors in the measured intensity lead to greater errors in *p* if the intensity is small. Therefore for samples showing very different levels of attenuation in different directions, *i.e.* where *p* varies markedly with projection angle, the ring artefacts may have irregular shapes. See, for example, Fig. 1(*b*) which relates to a plate sample; the artefacts are most intense for arcs almost parallel to the longest axis of the sample.

In this paper we try to find *a priori* information about the attenuation factor *p*(*x*,*y*,*t*) that will help us to suppress these irregular ring artefacts. *I*(*x*,*y*) is the intensity of the X-ray beam just before the scintillator after the beam has passed through the sample. Let *I*_{v}(*x*,*y*) be the intensity of the visible light just after the scintillator. Often it can be assumed that the transform is linear and may be written as a convolution,

where *K*(*x*,*y*) is the point-spread function (PSF). In many applications *K*(*x*,*y*) has an almost Gaussian shape (Banhart, 2008; Mahajan, 2001), *i.e.*

where . Since

and , then, owing to the differentiation property of the convolution, we obtain

So it can be assumed that *I*_{v}(*x*,*y*) is a smooth function. Therefore, we can suppose that *p*(*x*,*y*,*t*) is also a smooth function of *x* and *y*.

Various assumptions on *p* and *q* have been made in a series of papers submitted by the authors. In the case of regular ring artefacts two methods have been proposed in Titarenko *et al.* (2009, 2010*a*). In the former, the correction of ring artefacts is based on a knowledge of the attenuation coefficients in some areas of a reconstructed slice. In this case the difference between the exact image and the image obtained from a measured sinogram is almost a constant over each half-ring (all full-rings are concentric with respect to the centre of rotation of the sample). In the latter, the smoothness of the sinogram along the horizontal coordinate, as well as equivalence between the first and the final rows of the sinogram (they should be the same in the case of 360° rotation or flipped with respect to the axis of rotation in the case of 180° rotation) are used to obtain a quadratic programming problem whose solution is a time-independent function, *q*(*x*). In contrast to the former, explicit values of the attenuation coefficients are not used. For some simple, but often encountered, cases an analytical formula has been proposed in Titarenko *et al.* (2010*b*); this formula is also used in the method proposed in this paper. A case of irregular ring artefacts caused by small vibrations in experimental set-up has been considered in Titarenko *et al.* (2010*c*). In the current paper we consider more general irregular ring artefacts.

We think it is impossible in principle to develop an algorithm that suppresses ring artefacts for any dataset. The reason lies in the nature of tomography. Image reconstruction from a given set of projections is an ill-posed problem. This means that at least one of the following three requirements is not satisfied: the solution exists, is unique and stable. Image reconstruction is always an unstable problem, *i.e.* small errors in input data lead to significant variations in the solution, see examples by Lavrent'ev *et al.* (2001) and Natterer & Wübbeling (2007). Also from the theory of ill-posed problems it is known that without *a priori* information about properties of a sample and noise levels in input data you may achieve a solution which varies significantly from the exact one. Note that many well known methods of image reconstruction assume some properties of the solution, which are often not met in real problems. In the case of filtered back-projection the solution is to be an infinite number multiplied by the differentiable function, *i.e.* it cannot be used to reconstruct a sample with any edges/boundaries, even a sphere or cylinder. In addition, the methods have some internal parameters allowing a user to `regularize' a solution in order to avoid high noise level in a solution. This `regularization' is often implemented as cutting high frequencies in the input data.

Similar ideas of `filtering' are used in all the methods we know to suppress ring artefacts, which have been developed over the last 20 years. Of course, this `filtering' is implemented in different ways. There seem to be three main approaches. For the first, a data acquisition is modified in order to `blur' ring artefacts. In Davis & Elliott (1997) time-delay integration is used, so a detector is moved laterally during acquisition. This approach allows users to avoid strong rings on a given slice; however, it may also increase the noise level on several neighbouring slices. The second approach is based on post-processing of the reconstructed images (see, for example, Sijbers & Postnov, 2004; Walls *et al.*, 2005; Yang *et al.*, 2008; Titarenko *et al.*, 2009). A polar-to-Cartesian coordinates transform is often used to find the rings. The third approach is based on applying various filtering techniques to the input data (often sinograms) (see Antoine *et al.*, 2002; Raven, 1998; Boin & Haibel, 2006; Titarenko *et al.*, 2010*b*; Titarenko & Yagola, 2010). Of course, some methods may combine several approaches. The method proposed in this paper belongs to the last group.

In our opinion, owing to the ill-posedness of image reconstruction, data processing does not allow ring artefacts to be suppressed without real features being suppressed at the same time. A good example is a sample with a small dense particle near the centre of rotation. In this case a sinogram will have a column which is brighter than the neighbouring ones. The same sinogram may be formed in the case of a ring artefact and no dense particle. Other examples may involve samples with concentric features. It seems the best approach is to combine data processing and acquisition, so that if there are any doubts as to whether we see an artefact it is possible to acquire more data (moving/rotating/inclining a sample) to check this. Of course, it may be possible to develop a ring-suppression method specifically for a particular class of sample such that real features may be easily distinguished. However, in the case of artefacts developed during acquisition or from some non-linear effects, *e.g.* beam-hardening, it can be very hard to suppress ring artefacts.

### 2. Ring artefact suppression

Let us take a horizontal slice of data, *i.e.* the data recorded by a single row of pixels over all projections. This can be represented by the matrix , where *t*_{i}, , and *x*_{j}, , are grids (*m* and *n* are numbers), *i* is the projection number and *j* the pixel in the row. In many cases, *e.g.* when the sample rotates at constant speed, the variable *t* (time) may be replaced by the angle of rotation of the sample . However, in the case of variations in beam intensity, the temperature of the scintillator or the camera or the strength of the artefacts depends not only on so it is better to use *t*. The grid {*x*_{j}}_{1}^{n} is uniform (*x*_{j}-*x*_{j-1} is often the pixel size), {*t*_{i}}_{1}^{m} may be arbitrary, since we have not assumed that *p*(*x*,*y*,*t*) is a smooth function of *t*. Therefore

where *p*_{ij} must be found. Finding *p*_{ij} is equivalent to finding *q*_{ij}, which we will use thereafter. *q* denotes the matrix *q*_{ij}. Since *p*(*x*,*y*,*t*) is a smooth function of *x*, we require that

attains its minimal value for a given error in measurement. Another assumption is that the error in measurements is relatively small, *i.e.*

should tend to zero.

Let us introduce a smoothing (Tikhonov's) functional

where is a regularization parameter. Then the minimizer of this equation tends to the exact solution as all errors in the measurements tend to zero. Note that the regularization parameter depends on errors and cannot be chosen arbitrarily. However, for given non-zero errors one could vary in order to satisfy some restrictions on the image structure. More information about methods of solving ill-posed problems can be found in Titarenko *et al.* (1995, 1998), Engl *et al.* (1996), Bakushinsky & Kokurin (2004) and Ivanov *et al.* (2002).

Now we assume that *q*_{ij} can be written as a final sum

where *f*_{is} are orthonormal vectors, *i.e.*

The vector is defined for . Let and . Then the vectors *c*_{w} and *g*_{w} are introduced such that

with elements defined as

Now our aim is to find a global minimizer *c*_{sj} of the smoothing functional (12). As shown in *Appendix A* the property (14) allows us to obtain *S* systems of linear equations (*n* equations with *n* variables)

where *T* is a transposed vector, *i.e.* a column. All these systems have the same matrix *A*.

To find the solution of (18) we use the analytical formula found in Titarenko *et al.* (2010*b*). Then

where

and . To avoid large numbers the last formula can be rewritten as

for . The matrix *H*_{jk} is symmetric, so is trivial.

For the practical applications described below we use the following set of orthonormal vectors (see Hsu, 1995),

These vectors may be abbreviated to

Note that there are only *m* independent vectors, so .

### 3. Practical applications

#### 3.1. Test object 1: a high-aspect-ratio (laterally extended) sample

Test object 1 is calcified cartilage from the pectoral girdle of a fossilized symmoriid shark. A series of projections was acquired at the microtomography beamline 2-BM at the Advanced Photon Source (APS) at Argonne National Laboratory, USA. A double multilayer monochromator was used. The X-ray beam had a parallel geometry so a standard filtered back-projection (FBP) algorithm (see *e.g.* Natterer & Wübbeling, 2007; Kak & Slaney, 1988) was applied to the sinograms preprocessed by the algorithm proposed above. Although the sample comprises phases of similar attenuation coefficients it has a large aspect ratio () in the projections.

From Figs. 1(*a*) and 1(*b*) it can be seen that the standard flat-field correction technique is not effective owing to the large aspect ratio. As a result the strength of ring artefacts seen in the reconstructed horizontal slice in Fig. 1(*b*) varies with angle about the centre of the slice, being largest normal to the long axis of the sample, *i.e.* the arcs are most significant essentially where they are parallel to the long axis. After applying the proposed ring-suppression method [taking , *S* = 21 and the functions *f*_{is} defined in (23)] the artefacts have been suppressed, as demonstrated in Fig. 1(*c*). Note that weak ring artefacts still persist near the edge of the reconstructed slice. However, they are not intense and can be suppressed further if a larger value of *S* is used. The value *S* = 21 means that the minimal period of the orthogonal functions defined in (23) is one-tenth of the height of the sinogram (the height is measured along variable *t*). Roughly speaking one can expect that artefacts whose period is greater than the minimal period will be suppressed, whereas artefacts with a smaller period may still persist but with decreased strength. The parameter *S* controls how ring artefacts are suppressed along a radial coordinate (in a polar frame) and along the angular coordinate. *S* = 21 is taken, since the minimal period of artefacts seems to be about one-tenth of the sinogram's height, and seems to be optimal for the suppression along the angular coordinate (for a larger some rings remain while for a smaller some genuine features become suppressed). Note that there is also a `spider' artefact near the centre of Fig. 1(*c*) which is discussed below.

#### 3.2. Test object 2: compound sample

In order to more seriously challenge the correction algorithm a sample was devised showing more extreme heterogeneity of the ^{-3}; the material has a flat two-dimensional woven architecture of high modulus fibres) sandwiched between two (2.5 ×1.0 and 6.5 ×1.2 mm) blocks of Ti matrix/unidirectional 140 µm-diameter SiC monofilament composite all wrapped up in (30 µm) aluminium foil. The dataset was acquired at the B-16 beamline of the Diamond Light Source, UK. The X-ray beam had a white spectrum and parallel geometry. FBP was used to reconstruct slices. While the aspect ratio of the sample is relatively small (1.5 : 1.0) the is highly anisotropic because the carbon composite at the centre attenuates only slightly whereas the metal composite attenuates the X-ray beam significantly. In addition, the presence of the SiC fibres in the metal composite makes this region locally heterogeneous. Beam hardening arising from the white beam is also apparent as black lines parallel to the metal composite blocks and diagonal lines in Fig. 2, making the treatment of the ring artefacts more complicated.

### 4. Discussion

In applying our approach a key issue is the choice of and *S*. First let us discuss how should depend on *s*. One may use the same for all *f*_{is}. However, the function *q*(*x*_{i},*y*,*t*) is assumed to be smooth along *t*. Note that Fourier coefficients of a smooth function decrease as 1/*s* (see Beerends *et al.*, 2003). Since the functions *f*_{is} defined in (23) are discrete Fourier coefficients, the use of = would seem to be more reasonable than = . This strategy can be assessed by comparing images from the right- and left-hand columns of Fig. 3. For small values of *S* the quality of the reconstruction of the titanium composite blocks is almost the same for both cases, while for large *S* the case of being constant for all *s* provides more wave artefacts. For example, the carbon cores of the 140 µm SiC fibres are clearer in many of the fibres in Fig. 3(*b*), while there are additional wave features around the blocks in Fig. 3(*c*) which decrease the quality of the image.

Another important question is the most appropriate value of *S*. From the arguments of the previous section one might suppose that increasing *S* indefinitely might give better and better results. However, increasing *S* can have important implications with regard to periodic structures. Compare Figs. 3(*b*) (*S* = 5) and (*d*) (*S* = 51) where there are new dark arcs near the C cores of the SiC fibres. In Figs. 4(*a*)–4(*c*) it can be seen how genuine features (thin curves) in the sinogram are suppressed and how the sinogram flattens as *S* increases. Note the horizontal lines in Figs. 4(*a*)–4(*c*) are due to intensity profile variations (see Titarenko *et al.*, 2010*c*). Therefore, if the value of *S* is increased then the genuine features in the reconstructed image may be suppressed.

An interesting feature is the `spider' feature that is absent for *S* = 5, but becomes increasingly prominent as *S* increases (see Fig. 4). In order to consider its origin let us undertake a gedanken experiment. Consider imaging a simple ball for which we vary the distance between the ball and the centre of rotation, and for which there are no errors in the simulations (see Fig. 5). The proposed correction is applied to both sinograms using the same parameters (*S* = 100, = 0.0001*s*^{2}). Comparing Figs. 5(*a*) and 5(*b*) one can see that the footprint of the ball over the sinogram is affected more drastically when it is closer to the axis of rotation. It would seem that features which are closer to the centre of rotation are suppressed earlier when *S* increases and/or decreases. Clearly this is unsatisfactory and leads to flattening of the detail very near to the centre of rotation, giving rise to the characteristic `spider' noted above [see Figs. 4(*d*)–4(*f*)].

In order to circumvent this problem it is preferable to increase *S* with the distance from the centre of the sinogram. For example, see Fig. 6 where *S* = 5 inside the dashed circle and *S* = 30 outside. Of course, one can also vary *S* continuously from the centre to the edges; however, even the stepwise approach gives good results. In the general case when varying *S* a compromise between suppression of the ring artefacts towards the periphery of the slice and the suppression of genuine features towards the centre of the slice should be found.

Now we provide some estimates for the times required to preprocess a sinogram with the method proposed in this paper. We have used a desktop computer with a dual core processor (E7200, Intel Core 2 Duo, 2.53 GHz) and Intel performance libraries. To form the *n*×*n* matrix *H*, see (21), it takes 0.11 s for *n* = 2048 (the case of the shark cartilage) and 0.44 s for *n* = 3880 (the case of the compound sample). Once the matrix is found the solution can be obtained in less than 0.01 s in both cases. Therefore, in order to process a three-dimensional volume it is better to find *S* *n*×*n* matrices first and store them to memory.

### 5. Conclusions

We have developed a new method for the suppression of ring artefacts suited to anisotropically attenuating objects which can be applied to both parallel and cone-beam geometries. This anisotropy may be because the object is laterally extended such as a plate, or because the attenuation across the sample varies heterogeneously. The series of projections can be processed as a set of parallel two-dimensional slices. The slices can be chosen by fixing either rows or columns of the CCD sensor. A solution can be quickly realised because an analytical solution exists that is a product of a known matrix and a vector. Critical to the method are appropriate choices of and *S*. We have found that works well. As to the number of terms, *S*, we found that increasing *S* ends to preferentially suppress genuine features towards the centre of the slice, giving rise to a `spider' artefact. As a result a better solution is to increase *S* with distance from the centre of the sinogram.

Two datasets preprocessed by the proposed method were acquired at different synchrotron sources, *i.e.* using a parallel geometry of the beam, and reconstructed by FBP. Of course, the method can also be applied to data obtained by laboratory tomography, *i.e.* in cone geometry of the beam, if one preprocesses the data taken by a given row of the CCD sensor. The method is based on preprocessing two-dimensional data, *i.e.* data taken by one row of the sensor. This means that neighbouring slices/sinograms are preprocessed independently. However, in forthcoming papers the authors intend to extend the proposed method to three-dimensional data, *i.e.* use the assumption that two neighbouring slices should vary slightly. The authors believe that the studies presented here suggest that the method will be widely applicable to biological and medical specimens taken using white X-ray beams.

### APPENDIX A

### Formula derivation

Let there be a smoothing functional

with the regularization parameter . We know that an element *q*_{ij} can be written as a final sum

where vectors *f*_{is} are orthonormal, *i.e.*

Function (24) is continuously differentiable at any point *c* with elements *c*_{sj}, , . A point *c* can be a local minimizer of the smoothing functional only if

see, for example, Dennis & Schnabel (1996).

Taking (26) into account we obtain

Since we have defined *p*_{i0} = *p*_{i1} and *p*_{i,n+1} = *p*_{in}, then

Therefore

Using vectors *c*_{w}, *g*_{w} defined in (15) and (16) the local minimizer of the smoothing functional (24) should satisfy systems

where the matrix is the same for all *w*,

Note that the Hessian matrix of with respect to variables *c*_{wu} is a diagonal matrix with diagonal elements , *i.e.* it is positive definite. Therefore, the solution found is a local minimizer (see Dennis & Schnabel, 1996). Since the systems have only unique solutions, the element found, *c*, is the global minimum of Tikhonov's functional .

### Acknowledgements

The authors would like to thank the Diamond Light Source for use of the X-ray light source. The use of the Advanced Photon Source was supported by the US Department of Energy, Office of Science, Office of Basic Energy Sciences, under Contract No. DE-AC02-06CH11357. The piece of shark cartilage was provided by Professor Michael Coates (University of Chicago) and scanned by Mason Dean (Max Planck Institute). ST is grateful to EPSRC for funds through the `Collaborating for Success' grant. VT and AK are grateful to STFC for funding the project.

### References

Antoine, C., Nygåard, P., Gregersen, Ø. W., Holmstad, R., Weitkamp, T. & Rau, C. (2002). *Nucl. Instrum. Methods Phys. Res. A*, **490**, 392–402. CrossRef CAS Google Scholar

Bakushinsky, A. B. & Kokurin, M. Y. (2004). *Iterative Methods for Approximate Solution of Inverse Problems*, Vol. 577, *Mathematics and its Applications*. Dordrecht: Springer. Google Scholar

Banhart, J. (2008). Editor. *Advanced Tomographic Methods in Materials Research and Engineering.* No. 66, *Monographs on the Physics and Chemistry of Materials*. Oxford University Press. Google Scholar

Bardelli, L. *et al.* (2006). *Nucl. Instrum. Methods Phys. Res. A*, **569**, 743–753. CrossRef CAS Google Scholar

Beerends, R. J., ter Morsche, H. G., van den Berg, J. C. & van de Vrie, E. M. (2003). *Fourier and Laplace Transforms*. Cambridge University Press. Google Scholar

Boin, M. & Haibel, A. (2006). *Opt. Express*, **14**, 12071–12075. Web of Science CrossRef PubMed Google Scholar

Buzug, T. M. (2008). *Computed Tomography: From Photon Statistics to Modern Cone-Beam CT*. Berlin: Springer. Google Scholar

Davis, G. R. & Elliott, J. S. (1997). *Nucl. Instrum. Methods Phys. Res. A*, **394**, 157–162. CrossRef CAS Web of Science Google Scholar

Dennis, J. E. & Schnabel, R. B. (1996). *Numerical Methods for Unconstrained Optimization and Nonlinear Equations*. Philadelphia: SIAM. Google Scholar

Engl, H. W., Hanke, M. & Neubauer, A. (1996). *Regularization of Inverse Problems*, Vol. 375, *Mathematics and its Applications*. Dordrecht: Kluwer Academic Publishers. Google Scholar

Hsu, H. P. (1995). *Schaum's Outline of Signals and Systems*. New York: McGraw-Hill. Google Scholar

Ivanov, V. K., Vasin, V. V. & Tanana, V. P. (2002). *Theory of Linear Ill-posed Problems and its Applications*. Utrecht: VSP. Google Scholar

Kak, A. C. & Slaney, M. (1988). *Principles of Computerized Tomographics Imaging*, Vol. 33, *Classics in Applied Mathematics*. Philadelphia: SIAM. Google Scholar

Krumm, M., Kasperl, S. & Franz, M. (2008). *NDT E Int.* **41**, 242–251. CrossRef Google Scholar

Lavrent'ev, M. M., Zerkal, S. M. & Trofimov, O. E. (2001). *Computer Modelling in Tomography and Ill-Posed Problems*. Utrecht: VSP. Google Scholar

Mahajan, V. N. (2001). *Optical Imaging and Aberrations*, Part II, *Wave Diffraction Optics*. Bellingham: SPIE Press. Google Scholar

Moszyński, M. (2006). *Radiation Detectors for Medical Applications*, pp. 293–315. *NATO Security through Science Series.* Berlin: Springer. Google Scholar

Natterer, F. & Wübbeling, F. (2007). *Mathematical Methods in Image Reconstruction*, Vol. 5, *SIAM Monographs on Mathematical Modeling and Computation.* Philadelphia: SIAM. Google Scholar

Raven, C. (1998). *Rev. Sci. Instrum.* **69**, 2978–2980. Web of Science CrossRef CAS Google Scholar

Sijbers, J. & Postnov, A. (2004). *Phys. Med. Biol.* **49**, N247–N253. Web of Science CrossRef PubMed Google Scholar

Stock, S. R. (2008). *MicroComputed Tomography: Methodology and Applications*. Boca Raton: CRC Press. Google Scholar

Tikhonov, A. N., Goncharsky, A. V., Stepanov, V. V. & Yagola, A. G. (1995). *Numerical Methods for the Solution of Ill-posed Problems*, Vol. 328, *Mathematics and its Applications*. Dordrecht: Kluwer Academic Publishers. Google Scholar

Tikhonov, A. N., Leonov, A. S. & Yagola, A. G. (1998). *Nonlinear Ill-Posed Problems*, Vols. 1, 2 and 14, *Applied Mathematics and Mathematical Computation*. London: Chapman and Hall. Google Scholar

Titarenko, S., Titarenko, V., Kyrieleis, A. & Withers, P. J. (2009). *Appl. Phys. Lett.* **95**, 071113. Web of Science CrossRef Google Scholar

Titarenko, S., Titarenko, V., Kyrieleis, A. & Withers, P. J. (2010*a*). *J. Synchrotron Rad.* **17**, 540–549. Web of Science CrossRef IUCr Journals Google Scholar

Titarenko, S., Withers, P. J. & Yagola, A. (2010*b*). *Appl. Math. Lett.* **23**, 1489–1495. Web of Science CrossRef Google Scholar

Titarenko, S. S. & Yagola, A. G. (2010). *Moscow University Phys. Bull.* **65**, 65–67. Web of Science CrossRef Google Scholar

Titarenko, V., Titarenko, S., Withers, P. J., De Carlo, F. & Xiao, X. (2010*c*). *J. Synchrotron Rad.* **17**, 689–699. Web of Science CrossRef CAS IUCr Journals Google Scholar

Walls, J. R., Sled, J. G., Sharpe, J. & Henkelman, R. M. (2005). *Phys. Med. Biol.* **50**, 4645–4665. Web of Science CrossRef PubMed Google Scholar

Yang, J., Zhen, X., Zhou, L., Zhang, S., Wang, Z., Zhu, L. & Lu, W. (2008). *The 2nd International Conference on Bioinformatics and Biomedical Engineering*, pp. 2386–2389. Piscataway: IEEE. Google Scholar

© International Union of Crystallography. Prior permission is not required to reproduce short quotations, tables and figures from this article, provided the original authors and source are cited. For more information, click here.