Multimodal Phase-Based X-Ray Microtomography with Nonmicrofocal Laboratory Sources

Fabio A. Vittoria, Marco Endrizzi, Gibril K. Kallon, Charlotte K. Hagen, Francesco Iacoviello, Paolo De Coppi, and Alessandro Olivo Department of Medical Physics and Biomedical Engineering, University College London, Malet Place, Gower Street, London WC1E 6BT, United Kingdom Electrochemical Innovation Lab, Department of Chemical Engineering, University College London, London WC1E 7JE, United Kingdom Institute of Child Health, University College London, London WC1N 1EH, United Kingdom (Received 30 August 2017; revised manuscript received 18 October 2017; published 8 December 2017)


I. INTRODUCTION
X-ray imaging is a widely used technique for the study and investigation of specimens' internal structure in different fields of application.X rays traveling through matter are partially attenuated, and the degree of attenuation depends on the sample's chemical composition and density.This is the only contrast mechanism in standard, attenuation-based x-ray imaging.An alternative approach is represented by x-ray phase-contrast imaging (XPCI).In x-ray imaging, a sample can be described in terms of its complex refractive index n ¼ 1 − δ þ iβ.The parameter β is responsible for x-ray attenuation, while δ describes phase effects.
X rays, in fact, are electromagnetic waves which, when traveling through a medium, experience a phase shift that depends on δ and on the thickness of the medium.An XPCI system converts this phase shift into an intensity modulation that can be measured with a detector.XPCI can be implemented in different modalities with different experimental setups [1].Two XPCI techniques have differentiated themselves from the others owing to their ability to provide strong phase signals with nonmicrofocal laboratory sources: grating interferometry (GI) and edge illumination (EI).GI and EI belong to a wider class of XPCI methods sometimes referred to as "differential-phase-contrast" imaging methods.They are sensitive to the refraction induced by a sample, which is proportional to the gradient of the phase shift.
When the gradient of the phase shift varies on a scale smaller than the system resolution, the refraction signal cannot be "spatially resolved," and the radiation coming from the sample appears to be scattered at different angles.This effect has been referred to as ultra-small-angle x-ray scattering (USAXS) [2], and it can reveal small structures or inhomogeneities in a sample, on a scale smaller than the system resolution.In GI and EI, refraction sensitivity is achieved by using two or three optical elements, mutually aligned with micrometric or submicrometric precision.In a single raw image acquired with a GI or EI system, attenuation, refraction, and USAXS all contribute to the contrast between different materials.In order to separate these three contributions, at least three different acquisitions are required, with different relative positions of the optical elements [3,4].
The retrieval of attenuation, refraction, and USAXS is particularly important for the implementation of XPCI in computed tomography (CT) [5,6].Attenuation, refraction, and USAXS, in fact, can be expressed as line integrals along the x-ray propagation direction of fundamental properties of the sample [7].If multiple projections of each signal are acquired while rotating the sample with respect to the incoming radiation, the line integrals can be inverted and it is possible to reconstruct the threedimensional distributions of the linear attenuation coefficient (responsible for attenuation) of the real part of the refractive index (responsible for refraction) and of the width of the local scattering distribution (responsible for USAXS) [7].
While GI and EI are highly sensitive to small angular refractions, the need to perform multiple acquisitions while varying the relative positions of the optical elements imposes significant limitations: the total exposure time is increased proportionally with the number of images required for the retrieval procedure (at least three); any error in the position of the optical elements translates into an error in the retrieved signals; and the relative alignment between the optical elements should, ideally, remain stable for the entire duration of a measurement, which can be challenging for long acquisitions, such as those required in CT.
For the above reasons, alternative, "single-shot" XPCI methods have been proposed recently in which attenuation and refraction can be retrieved from a single acquisition, without the need to move any optical elements [8][9][10][11][12][13][14].Among those methods, speckle tracking [13] and beam tracking [14,15] can also retrieve the USAXS signal, and they have been successfully implemented with laboratory sources.The common idea behind those single-shot methods is to use a phase or an intensity modulator to create a known reference intensity pattern on an x-ray detector.When a sample is introduced into the beam, the reference pattern is locally modified.The variations induced by the sample to the intensity pattern can be modeled in terms of sample attenuation, refraction, and USAXS.
Different strategies can then be applied to invert this model and retrieve the three quantities of interest from the acquired images [15,16].In speckle tracking, the reference pattern is obtained as the near-field interference pattern generated from a random phase modulator.The amplitude of the intensity fluctuations of the speckle pattern is strongly affected by the source size, requiring the use of microfocal sources when implemented in the laboratory.Important progress recently achieved with the speckle-tracking method was its implementation in computed tomography with a polychromatic laboratory source [17].In that case, however, only attenuation and refraction were reconstructed.It has been shown, in fact, that when speckle tracking is implemented with polychromatic sources, a spurious USAXS signal can be detected due to sample attenuation [18,19], rather than its microstructures, and that the USAXS signal retrieved from microstructured samples does not increase linearly with the sample thickness [19].
A similar spurious USAXS signal was also reported for other phase-contrast imaging methods [20,21].In the beam-tracking approach, the reference pattern is created as the projection image of an absorbing mask featuring periodically repeated apertures [8,9,15].A consequence of using an intensity modulator rather than a phase modulator to create the reference pattern is that the method is less sensitive to spatial coherence and can tolerate larger source sizes.Additionally, the strong intensity modulation of the mask allows beam tracking to be used with samples simultaneously generating strong refraction and USAXS signals [22].Most importantly for CT applications, even when implemented with polychromatic sources, the USAXS signal retrieved with the beam-tracking method is proportional to the sample thickness [19].
In this work, we present a tomographic implementation of beam tracking with laboratory sources.The method is tested with a nonmicrofocal source and a high-resolution detector on a large variety of samples featuring high and low attenuation, refraction, and USAXS signals.

II. THEORETICAL BACKGROUND
A scheme of the experimental setup is shown in Fig. 1.The primary x-ray beam generated by the source is divided by the mask into several secondary beamlets.The mask apertures extend along the y direction, perpendicular to the plane in the drawing in Fig. 1, creating laminar beamlets.After propagation through the sample and in free space, the beamlets reach the detector, where their intensity profile is recorded.The measurement is repeated while rotating the sample around the y axis.
Let us consider a single row of detector pixels in the x direction and a specific angular position θ of the sample with respect to the incoming radiation.We previously showed that the intensity profile of each beamlet can be fitted with the following model function [15]: USAXS increases its standard deviation [15].Attenuation (L), refraction (R), and USAXS (S) are retrieved from the parameters A, B, and C, measured with and without the sample, using the following expressions: where z 3 is the sample-to-detector propagation distance and we introduce the transmission T ¼ A o =A i .The subscripts o and i refer to the quantities retrieved from the profiles acquired with and without the sample, respectively.
In the case of monochromatic radiation [23], the attenuation, refraction, and USAXS retrieved from a particular aperture of the mask can be written as the following line integrals: where λ is the x-ray wavelength, l is the sample thickness, x a is the position of the mask aperture along the x direction, μ ¼ ð4π=λÞβ is the linear attenuation coefficient, β is the imaginary part of the complex refractive index, δ is the difference from unity of the real part of the complex refractive index, and σ 2 ϕ is the variance of the local USAXS distribution function ϕ [23].
When polychromatic radiation is used, the following alternative expressions can be found [15]: where PðλÞ depends on the source energy spectrum, detector energy response, and transmission through the mask [15].In the last term of Eqs. ( 9) and ( 10), we make the assumption that the attenuation through the sample is weak, so that the weighted average over PðλÞTðx a ; θ; λÞ of Rðx a ; θ; λÞ and Sðx a ; θ; λÞ can be approximated with their weighted average over PðλÞ.
For weakly attenuating materials, Tðx a ; θ; λÞ ≈ 1 − Lðx a ; θ; λÞ, and it is possible to derive the following expressions: where hfi P indicates the weighted average of a function f, with weights given by P. Equations ( 11)-( 13) are the analogs of Eqs. ( 5)- (7), in the case of polychromatic radiation, where μ, δ, and σ 2 ϕ are replaced by the average quantities hμi P , hδi P , and hσ 2 ϕ i P .It is important to stress that Eqs. ( 11)-( 13) are derived in the assumption of low attenuation from the sample.For materials exhibiting strong x-ray attenuation, tomographic reconstructions performed with polychromatic radiation can result in the presence of artifacts in the reconstructed slices.A typical example is the "cupping" artifact encountered in standard attenuation-based CT performed with polychromatic radiation [24].The study of the artifacts encountered when the above assumption is not satisfied goes beyond the scope of this paper and will be the subject of future investigation.Equations ( 11)-( 13) can be inverted by using a tomographic reconstruction algorithm to obtain the three-dimensional maps of hμi P , hδi P , and hσ 2 ϕ i P from a series of projections acquired at different angles.

A. Experimental apparatus
The implemented source, a rotating tungsten tube, is operated at 46 kV(peak) and 26 mA and is characterized by a focal spot of about 70 μm full width at half maximum.The detector is a CCD camera coupled through a fiber-optic plate to a Gadox scintillator with an effective pixel size equal to 4.54 μm (Photonic Science).The mask is made of a 225-μm gold layer on a 525-μm silicon substrate.The aperture size and period of the mask are 9 and 93 μm, respectively.The source-to-mask, mask-to-sample, and sample-to-detector distances are z 1 ¼ 168 cm, z 2 ¼ 4 cm, and z 3 ¼ 30 cm, respectively.
Similar to the EI method [25], in the beam-tracking approach the intrinsic resolution in the direction perpendicular to the mask apertures is, to good approximation, determined by the size of the mask aperture, and it is therefore higher than the sampling rate of an image reconstructed from a single acquisition with the sample, given by the mask period.For this reason, we perform a ten-step scan of the sample along the x direction across one period of the mask, with the scanning step equal to 9.5 μm.
The signals extracted from each step are then recombined in a single "oversampled" image, with the effective pixel size equal to 9.5 × 3.87 μm 2 .In the direction parallel to the mask apertures, in fact, the effective pixel size is equal to the demagnified detector pixel.This process is repeated for different angular positions of the sample in the angular range (0°and 180°) to obtain a complete tomographic scan.The total number of projections acquired in the angular range above is 721 for the samples in Fig. 2, and 1441 for the samples in Figs. 3 and 4. Each frame is acquired with a 5-s exposure time.
To measure the detector dark current, ten frames are acquired before and after each scan with the source turned off, averaged, and subtracted from all of the other measurements.During the scan, after every 20 angular positions, a reference image is obtained by acquiring and averaging ten reference frames without the sample in the beam.

B. Scanning electron microscopy
A Zeiss EVO 10 system (Zeiss, Oberkochen, Germany) is utilized in backscattered electron mode.Before imaging, the sample is gold coated with an SC7620 mini sputter coater (Quorum Technologies, Laughton, England) for 30 s.

C. Tomographic reconstruction
CT reconstructions are performed with a customdesigned iterative reconstruction algorithm based on a gradient-descent minimization approach with Tikhonov-Miller regularization.For our analysis, we use the ASTRA toolbox [26], which contains an implementation of the radon transform and of its dual operator (also called back projection), together with several reconstruction algorithms.Note that the refraction signal Rðx a ; θÞ is equal to the line integral of the spatial derivative along the x direction of hδi P ðx a ; zÞ.Therefore, before performing the CT reconstruction, the retrieved refraction signal is numerically integrated along the x direction.
The reconstruction algorithm needs an initial guess of the reconstructed slice Rec 0 and the experimental sinogram Sino as input, for each of the reconstructed signals.We use a zero matrix as the starting guess.At each iteration, the reconstructed slice is updated with the following formula: RT and RT* are the Radon transform and its dual operator, respectively, ∇ 2 indicates the Laplacian operator, and a 1 and a 2 are two constant factors.The second term on the righthand side of Eq. ( 14) updates the reconstruction so that its projected sinogram matches the experimental one; a 1 determines the weight of this update.The last term, instead, represents the Tikhonov-Miller regularization and penalizes strong fluctuations in the reconstructed slice, reducing not only the noise level but also the final resolution; a 2 determines the weight of the regularization term.
The value of a 1 and a 2 are chosen ad hoc in order to obtain the best balance between the convergence speed of the algorithm, the image resolution, and the noise level.Specifically, the reconstruction of δ is performed without regularization (a 2 ¼ 0), due to the relatively low level of noise of the refraction signal.

IV. RESULTS AND DISCUSSION
We first image simple custom-made samples consisting of three wires and one rod of different known materials and diameters: maxima fishing line with a 300-μm diameter, polybutylene terephthalate (PBT) wire with a 180-μm diameter [chemical composition ðC 12 H 12 O 4 Þ n , density of 1.31 g=cm 3 ], nylon wire with a 150-μm diameter [chemical composition ðC 6 H 11 NOÞ n , density of 1.13 g=cm 3 ], and aluminum rod with a 2-mm diameter (density of 2.70 g=cm 3 ).
Figure 2 shows slices of the reconstructed attenuation (μ), phase (δ), and USAXS (σ 2 ϕ ) volumes of the different samples.The scattering signal varies from zero only at the interface between two different materials, as expected from homogeneous samples made of a single material without internal microstructures.The average values of μ and δ measured in a circular area in the middle of each wire of ten-pixel diameter are shown in Table I, together with the theoretical values calculated from Eqs. ( 11) and ( 12), and the relative errors ϵ μ and ϵ δ between the two.The energydependent quantities μðλÞ and δðλÞ are obtained from the XRAYLIB library [27] for each material using the chemical compositions and densities previously listed.PðλÞ was obtained from a simulated source spectrum [28] and detector energy response [29].
By comparing the values of μ exp with the tabulated values of μðλÞ for each material, it is possible to define an effective wavelength λ eff for which μ exp ¼ μðλ eff Þ.The same process can be repeated for hμi P , δ exp , and hδi P .The values of the corresponding effective energies E eff ¼ hc=λ eff , where h is the Planck constant and c is the speed of light in vacuum, are shown in Table I.The effective energies E eff are the energies that would be required to obtain the same values of μ or δ when performing the same experiment with monochromatic radiation [30].
The results in Table I show how, for low-attenuating materials like PBT and nylon, the retrieved values of μ and δ are in good agreement with the theoretical values predicted by Eqs. ( 11) and ( 12), and the relative errors are compatible with the experimental uncertainty due to noise, which is higher for the attenuation signal.For materials like aluminum showing higher attenuation, on the other hand, the approximations needed to express the retrieved signals as line integrals along the photon path of fundamental properties of the sample no longer hold.In this case, the retrieved μ and δ values are lower than their theoretical values, and the corresponding effective energy is higher.This is a common effect in CT performed with polychromatic radiation: as radiation propagates through the sample, the low-energy part of the spectrum is attenuated more than the high-energy part, thus increasing the mean energy of the spectrum and, therefore, decreasing the retrieved μ and δ values.
We then image another custom-made sample of three different powders: ground coffee, all-purpose flour, and icing sugar.Figure 3 shows a reconstructed slice of the attenuation, phase, and USAXS signals.At this resolution level, the coffee sample does not appear as a powder: individual bean fragments and their internal structure can be seen in the images.This is, however, an interesting sample, as it shows how the high resolution and low noise achieved in the phase signal result in a clear visualization of the small structures inside the beans.Figure 3(d) shows an TABLE I. Average values of μ exp and δ exp measured experimentally from the different materials in Fig. 2; their corresponding theoretical values hμi P and hδi P , calculated from Eqs. ( 11) and (12); and the relative errors ϵ μ ¼ jðμ exp − hμi P Þ=hμi P j × 100% and ϵ δ ¼ jðδ exp − hδi P Þ=hδi P j × 100% between the two.The theoretical values for the maxima line are not indicated due to its unknown composition.image obtained from the same coffee powder using a scanning electron (backscattered electron signal), where the same structures observed in Fig. 3(b) are visible at higher resolution.Sugar and flour, however, appear to be homogeneously distributed, with local areas of higher concentration visible in the phase image.The average values of μ, δ, and σ 2 ϕ measured in a circular area in the middle of each sample of 60-pixel diameter are shown in Table II.

Material
In the case of the coffee sample, the selected area is completely contained within a single coffee grain.For this sample, it is not possible to compare the experimental values of μ and δ with theoretical predictions from Eqs. (11) and ( 12) due to the unknown chemical compositions and densities of the powders.While the attenuation and phase signals from sugar are higher than the one from flour, the scattering contrast is reversed, confirming the complementarity of the three signals and their potential for material discrimination.
The last sample we image is a decellularized rabbit esophagus.The sample is derived using detergent enzymatic treatment [31], an innovative decellularization method which aims at eliminating the cell population while preserving tissue microarchitecture and extracellular matrix.The aim of the method is to create a cell-free scaffold from animal or human organs which can then be populated with stem cells of the organ recipient for organ transplantation.All surgical procedures and animal husbandry were carried out in accordance with the United Kingdom Home Office guidelines under the Animals (Scientific Procedures) Act 1986 and were approved by the local ethics committee (University College London, England; Project License No. 70/2719).This sample was previously imaged at a synchrotron and in the laboratory using different phase-contrast imaging techniques [32] (free-space propagation at the synchrotron and edge illumination in the laboratory).
A complete description of the sample preparation and the results obtained can be found in Ref. [32].The signals reconstructed from the esophagus are shown in Fig. 4.This sample is characterized by a complex internal structure which is resolved with a high level of detail in the phase image [Fig.4(b)].This image reveals an intact scaffold microarchitecture, and it is possible to identify all of the layers of the native tissue (mucosa, submucosa, muscularis propria, and adventitia).In the attenuation image [Fig.4(a)], the contrast is too low to differentiate among the different structures of the sample, and the noise level is relatively high.Similar considerations apply to the scattering image in Fig. 4(c).

V. CONCLUSIONS
In this paper, we present an implementation of beamtracking CT with a nonmicrofocal laboratory setup.The method can reconstruct the three-dimensional maps of attenuation, phase, and USAXS of a given sample.Starting from our previous work on beam-tracking CT with monochromatic radiation [23] and on planar imaging with polychromatic radiation [15], we build a theoretical framework which describes the signal reconstructed in beamtracking CT with polychromatic radiation.The theoretical model is tested on samples of known composition and is shown to provide quantitatively correct values for the retrieved linear attenuation coefficient μ and for the difference from unity of the real part of the refractive index δ for low-attenuating materials.
The model, however, does not consider the variations of the spectrum induced by samples featuring higher attenuation.In this case, the effective energy of the spectrum is increased by a quantity which depends on the particular sample and its attenuation properties.While this is a common problem in CT performed with polychromatic radiation, it is possible to calibrate the system and its effective energy for specific materials and applications.Alternatively, for highly attenuating samples, a source with a higher energy spectrum can be used, together with additional filters to reduce the energy spread of the spectrum.It is important to note, however, that, in this case, a thicker mask would be preferable to reduce radiation transmitted through the absorbing septa.
The method is also tested on a variety of complex samples.The results in Fig. 3 show how the phase signal and its sensitivity are not negatively affected by the presence of a strong USAXS signal, and the results obtained from the esophagus scaffold confirm the possibility of using laboratory phase-contrast CT in tissue engineering and regenerative medicine, as has already been proposed by Hagen et al. [32].
The comparison between the three signals reconstructed with beam tracking shows a clear difference in overall image quality between the phase and the attenuation and USAXS signals.Reconstructed slices of δ consistently show higher levels of detail and signal-to-noise ratio with respect to μ and σ 2 ϕ .While this behavior is promising in terms of phasecontrast imaging, future work should explore possible causes of this behavior and, specifically, whether it is an intrinsic limitation of the beam-tracking approach or if the higher level of noise of the attenuation and USAXS signals might depend on the specific retrieval method used here.
The reconstructed USAXS slices, however, are in agreement with the expectations: the USAXS signal is zero for homogeneous samples made of a single material with no internal microstructures, and it varies from zero at the interface between two different materials or for a microstructured sample.Additionally, no artifacts can be seen in the reconstructed slices, confirming the linearity of the measured signal with the sample thickness.
In regards to the attenuation signal, it is important to note that the silicon substrate of the mask filters the source spectrum, resulting in a severe reduction of low-energy photons (below approximately 10 keV).While this filtration reduces the dose delivered to the sample, it also drastically reduces the contrast in all of the reconstructed signals.This might be one of the causes for the reduced signal-to-noise ratio obtained in the attenuation and USAXS images.Future developments of the technique might include the use of masks without substrate.These masks can be obtained by laser ablation in thin metal foils, an extremely cost-effective manufacturing process, and they have been recently implemented with the EI method [33].
FIG. 1. Scheme of the experimental setup.

FIG. 3 .
FIG. 3. Reconstructed slices of a test sample made of three different powders.(a) μ with a displayed dynamic range of [0, 160].(b) δ with a displayed dynamic range of [0, 7 × 10 −7 ].(c) σ 2 ϕ with a displayed dynamic range of [0, 4.5 × 10 −7 ].The images are obtained as the average of ten detector rows, resulting in a voxel size of 9.5 × 9.5 μm 2 in the displayed plane, and 38.7 μm in the perpendicular direction.(d) Backscattered electron image of the coffee sample, showing the internal structures visible in (b) at higher resolution.Each powder is contained in a different plastic straw.

FIG. 4 .
FIG. 4. Reconstructed slices of a decellularized rabbit esophagus.(a) μ with a displayed dynamic range of [0, 180].(b) δ with a displayed dynamic range of [0, 5 × 10 −7 ].(c) σ 2 ϕ with a displayed dynamic range of [0, 2.8 × 10 −7 ].The images are obtained as the average of ten detector rows, resulting in a voxel size of 9.5 × 9.5 μm 2 in the displayed plane, and 38.7 μm in the perpendicular direction.The esophagus is contained in a plastic straw.

TABLE II .
Average values of μ, δ, and σ 2 ϕ reconstructed from the different powders in Fig.3.