Six dimensional X-ray Tensor Tomography with a compact laboratory setup

Y. Sharma, M. Wieczorek, F. Schaff, S. Seyyedi, F. Prade, F. Pfeiffer, and T. Lasser Lehrstuhl f€ ur Biomedizinische Physik, Physik-Department and Institut f€ ur Medizintechnik, Technische Universit€ at M€ unchen, 85748 Garching, Germany Computer Aided Medical Procedures, Technische Universit€ at M€ unchen, 85748 Garching, Germany Institut f€ ur diagnostische und interventionelle Radiologie, Klinikum rechts der Isar, Technische Universit€ at M€ unchen, 81675 M€ unchen, Germany

Conventional X-ray imaging is based on the attenuation of X-rays when they pass through a material. Recently, grating interferometer based approaches [1][2][3] are being used to exploit the scattering and refraction of X-rays by the material in order to obtain additional contrast modalities, namely, phase contrast and dark-field contrast. Dark-field contrast originates from the ultra-small angle scattering of X-rays by sub-micron and micron sized structures. [4][5][6][7][8] Unlike attenuation, the dark-field signal obtained in a grating interferometer is anisotropic in nature, [9][10][11] which means that the measured signal depends on the orientation of the scattering structure with respect to the grating interferometer. Thus, by varying the relative orientation of a scattering structure within the interferometer, it is possible to deduce information about its orientation from the variations in the dark-field signal.
Two dimensional structure orientations can be recovered from several projections acquired by rotating the sample around the beam propagation direction. 12,13 However, in order to measure the structure orientations in three dimensions, Malecki et al. 14 introduced several non-standard acquisition poses in a technique called X-ray Tensor Tomography (XTT). This was achieved by placing the sample on an Eulerian cradle (Figure 1), which allows rotation of the sample with three degrees of freedom, and not just around one axis as done in conventional tomography. The three degrees of freedom are achieved by rotating the cradle around the y axis (angle w), rotating the sample stage around the z 0 axis by moving it on the cradle (angle h), and rotating the sample around the y 0 axis (angle /); here, the coordinate system (x, y, z) denotes world coordinates, and ðx 0 ; y 0 ; z 0 Þ is the sample coordinate frame ( Figure 1). A standard XTT dense sampling scheme S comprises several acquisition poses s i S ¼ fs i :¼ ðw; h; /Þ; w 2 f0 ; 20 ; 40 g; h 2 H; / 2 Ug; where H ¼ f0 ; 30 ; 60 ; 90 g; U ¼ f0 ; 2:01 ; …; 360 g. Such an acquisition geometry enables the measurement of the anisotropic dark-field signal from several poses spread over the unit sphere. These dark-field projections are then fed into an iterative reconstruction algorithm, 15 which splits the combined signal measured in all the projections into several auxiliary scattering components, each of which determines the strength of scattering along a unique orientation. The scattering orientations are then scaled with their respective scattering strengths, and a Principal Component Analysis is performed to fit a scattering tensor in each voxel. At last, the structure orientation is approximated as the smallest half axis of the tensor. An illustration of the XTT acquisition and reconstruction method is shown in Figure 2.
XTT is a promising technique with a number of applications such as visualization of the microstructure of fibrous materials and dentinal tubules of human teeth. 14,15 It is, however, limited by long acquisition times owing to the large number of measurements and the time required for the movement of several motors. Moreover, the Eulerian cradle is a bulky component and poses challenges for a compact setup design. Hence, there is a need to optimize the XTT acquisition geometry. In order to do so, it is essential to establish a relationship between the acquisition geometry and the expected result.
In order to visualize the acquisition geometries, we formulate an example tomographic trajectory rotating about y 0 as a function of w and h Aðw; hÞ :¼ fx :¼ ðw; h; /Þ; / 2 ½0 ; 18:95 ; …; 360 g: A common sparse acquisition scheme for conventional computed tomography (CT) using 20 equally spaced projection angles is then expressed as A(0 , 0 ). In other words, the sample is rotated around the fixed y axis and line integrals through the sample along z, also known as the optical axis, are recorded. Let us definetðxÞ as the unit vector along the optical axis for a given acquisition pose x :¼ ðw; h; /Þ. For visualization, we plot 6tðxÞ on a unit sphere, assuming that the sample is fixed at the center of the sphere. The points 6tðxÞ for A(0 , 0 ) are shown in blue in Figure 3(a). Another important variable for XTT is the normalized sensitivity axiŝ sðxÞ. The sensitivity axis of a grating interferometer is the direction in which the phase shift of X-rays is measured. It lies in the plane of the gratings, orthogonal to the grating bars ( Figure 1). We haveŝðxÞ ¼ y when the grating bars are horizontal. Assuming again that the sample is fixed at the center of a sphere, we also plot 6ŝðxÞ on the unit sphere. The points 6ŝðxÞ for A(0 , 0 ) are shown in red in Figure 3(a). It can be seen that only one direction of scattering can be measured with A(0 , 0 ). Therefore, additional tomographic trajectories for XTT (as shown in Figure 3) are obtained by rotating y 0 using w and h.
The visualization of the trajectories of the optical axiŝ tðxÞ and sensitivity axisŝðxÞ shown in Figure 3 gives a qualitative understanding of the orientations that are measured with a given acquisition scheme. Below, we present a method for the quantification of such acquisition schemes.
Let R and T be the two sets of evenly distributed points chosen from a hemisphere where N ¼ 4843. The mean vector of R is given by Let X be an acquisition scheme using n different acquisition poses x i The objective is to define a quantity, Coverage(X, k), that represents how well the orientation r k 2 R is measured by X.
Our proposed procedure for computing Coverage(X, k) is outlined in the following: (1) Define X k & X as the set of poses that measure the orientation r k 2 R X k :¼ fx 2 X; wðr k ; xÞ > T c g; where w(r k , x) is a weight factor that specifies how well the orientation r k is measured by the acquisition pose x, and T c ¼ 0.7 is an arbitrarily chosen threshold. w(r k , x) is given by where j Á Â Á j denotes the magnitude of the cross product and hÁ; Ái denotes the standard scalar product. (2) R k is the mean resulting length of the vectorstðxÞ for all where k Á k denotes the Euclidean length of a vector. (3) T k & T is the set of points that are measured by X k , assuming that each vectortðxÞ contributes to the measurement of points lying within a cone of opening angle 5 around it T k :¼ fs 2 T ; hs;tðxÞi > cosð2:5 Þ 8x 2 X k g: (4) The Coverage(X, k) is then defined as The quantity 2 Â ð1 À R k Þ is the spherical variance, 16 while jT k j=N is related to the cumulative solid angle spanned by thetðxÞ; x 2 X k . Thus, a higher value of Coverage(X, k) is achieved when thetðxÞ is distributed widely over the unit sphere.
Next, we scale the unit vectors r k 2 R with the values of Coverage(X, k) for all values of k resulting in CðXÞ :¼ fc k :¼ r k Â CoverageðX; kÞ; k ¼ 1; …; Ng; and its mean The angular deviation of the normalized mean vectorm CðXÞ from the normalized mean vectorm R is a measure of the non-uniformity of the Coverage Sphere. Finally, we can now introduce the proposed Coverage Metric CM(X) as Coverage X; k ð Þ: CMðXÞ 2 ½0; 1 is a measure of the efficiency of the acquisition protocol X. Higher values of CM imply a more comprehensive and uniform measurement of all orientations on the unit sphere.
In order to correlate the concept of coverage with experimental observations and to study the effect of sparse acquisition schemes, a Short Fiber Reinforced Polymer (SFRP) sample made of glass fibers (18 lm in diameter) is analyzed. SFRPs are widely used in the automotive industry for their improved mechanical properties, which depend heavily on the orientation and length distribution of the reinforcing carbon or glass fibers. A very high resolution X-ray micro computed tomography (XCT) technique combined with fiber tracking techniques is commonly used to analyze the fiber orientation distribution of such materials. 17 XTT is a very useful technique for these materials as it directly resolves the orientations of the fibers within a much larger sample size. We measure a SFRP sample (dimensions 10 Â 10 Â 2 mm) with the scheme S employing 3 axes of rotation, 7 phase steps, and 3 s exposure per phase step resulting in a total acquisition time of approximately 12.6 h. The setup 7 used for the XTT measurement comprises of a micro-focus X-ray tube (operated at voltage 45 kVp and power 25 W) and a Varian PaxScan 2520DX detector (pixel size 127 lm). The three gratings with periods of 10 lm, 5 lm, and 10 lm, respectively, were arranged in the first fractional Talbot configuration.
Scattering tensors T S ðr i Þ 2 R 3Â3 þ for voxels r i , i ¼ 1,…, I, discretizing the volume of interest (isotropic voxel size 64 lm 3 ), are calculated using the method proposed by Vogel et al. 15 and demonstrated in Figure 2; the subscript (here S) indicates the acquisition protocol used to compute the scattering tensors. Figure 4(a) shows the scattering tensors T S (r i ) along with the corresponding structure orientation overlaid on the XCT volume of the sample described above. The tensors are color coded with their orientation. It is evident from Figure 4(a) that XTT reveals the three dimensional orientations of the fibers in this sample. The resulting fiber orientations are comparable qualitatively to the fiber tracking results presented by Hannesschl€ ager et al. 18 using a voxel size of 6.5 lm 3 .
In the following, we will study sparse sampling schemes for XTT. S represents a densely sampled standard scheme within the practical limits of the hardware (jwj > 40 is not achievable since the cradle would block the X-ray beam). XTT reconstructions of the SFRP sample are obtained using different acquisition schemes X and the resulting tensors T X (r i ) are compared to the base T S (r i ). For comparison, we compute a distance metric between corresponding tensors in each voxel. We use the Riemannian manifold Sym þ 3 , i.e., the manifold of positive-definite symmetric matrices, equipped with a Riemannian metric. 20 The distance d : Sym þ 3 Â Sym þ 3 ! R þ between two tensors T 1 ; T 2 2 Sym þ 3 on this manifold is then computed as where r i denotes the i-th eigenvalue of the matrix T À 1 2 1 T 2 T À 1 2 1 . Two corresponding tensors T X (r i ) and T S (r i ) are considered to be correlated if dðT X ðr i Þ; T S ðr i ÞÞ is less than a certain threshold T p . Based on this, we introduce a Performance Metric PM(X), which determines how well the acquisition scheme X performs with respect to S   FIG. 3. Spherical representation of acquisition scheme A(w, h) for different combinations of w and h. On each sphere, blue points represent the trajectories of the optical axis 6tðxÞ; and red points represent the trajectories of the sensitivity axis 6ŝðxÞ assuming that the sample is fixed at the center of the sphere.

PM X
where T p ¼ 0.4 and j Á j is the number of elements in a set. We now use the presented numerical (CM) and experimental (PM) approach for the relative quantification of acquisition schemes for two studies: (b) Effect of Hardware Complexity: Hardware complexity is induced by the Eulerian cradle, which is required to provide the two additional axes of rotation (w and h). Therefore, in this study, we begin with the acquisition scheme B(4.2) as representative of a complex setup employing 3 rotation axes and rename it S 3 , the subscript here indicating the number of rotation axes. Next, we remove rotation axes in steps to obtain the schemes S 2 and S 1 . Scheme S 2 represents a medium complexity setup employing two rotation axes with fixed w (optimally chosen using CM), while scheme S 1 represents a simple setup employing only one rotation axis (fixed w and h) The Coverage Spheres for S 2 and S 1 are shown in Figure  4(b), and the reconstructed tensors are shown in Figure 4(a). The plot of CM and PM versus the number of rotation axes is shown in Figure 4(d). It can be seen that the quality of the reconstruction deteriorates slightly but not enough to affect their interpretation in an application when only a single value of w is used (i.e., a medium complexity setup using two rotation axes). However, using a simple setup with just one rotation axis and fixed w, h (scheme S 1 ) markedly deteriorates the quality of the reconstructions. This indicates that a compact setup of medium complexity, employing only two rotation axes, can be sufficient for XTT. The CM calculation can be extended in future work to incorporate prior information about the sample and to determine the exact value(s) of h required to measure all the structure orientations within the sample.
We conclude that it is possible to perform a full sixdimensional X-ray Tensor Tomography technique with at most two axes of sample rotation. Therefore, by adding only one additional rotation axis to a conventional X-ray Computed Tomography device, it is possible to obtain compact XTT setups for industrial and potential medical applications.