Helical sample-stepping for faster speckle-based multi-modal tomography with the Unified Modulated Pattern Analysis (UMPA) model

Speckle-based imaging (SBI) is a multi-modal X-ray imaging technique that gives access to absorption, phase-contrast, and dark-field signals from a single dataset. However, it is often difficult to disentangle the different signals from a single measurement. Having complementary data obtained by repeating the scan under slightly varied conditions (multiframe approach) can significantly enhance the accuracy of signal extraction and, consequently, improve the overall quality of the final reconstruction. In order to retrieve the different channels, SBI relies on a reference pattern, generated by the addition of a wavefront marker in the beam (i.e., a sandpaper or gratings). Here, we show how a continuous helical acquisition can extend the field of view (FOV) and speed up the acquisition while maintaining a multiframe approach for the signal retrieval of a test object.


Introduction
Fast multimodal X-ray imaging is valuable for expanding the capabilities of X-ray tomography beyond attenuation-based imaging.By combining different imaging modalities, it allows to obtain a more comprehensive view of internal structures, materials, and physical properties of objects.Several multimodal X-ray imaging techniques have been developed over the past decades: analyzer-based imaging (ABI) [1,2], grating-based imaging (GBI) [3,4], edge illumination (EI) [5], and more recently, speckle-based imaging (SBI) [6][7][8].These techniques have led to advancements in a wide range of applications (i.e, materials science, optics characterisation, and biomedical imaging).
SBI is a multi-modal X-ray technique that gives access to attenuation, phase contrast, and dark-field signals.Signal retrieval is based on the use of a reference pattern, generated by a wavefront marker (e.g., sandpaper [9,10] or 2D Talbot array illuminators [11][12][13]).Among others [14], Unified Modulated Pattern Analysis (UMPA) [15,16] is an algorithm capable of modelling local distortions of the reference pattern, which occur when a sample is inserted in the beam.Although the method allows signal retrieval from a single-shot acquisition, higher spatial resolution and better image quality can be achieved by acquiring multiple sets of reference and sample images with relative transverse locations of the sample or modulator [15].
UMPA can process datasets with two stepping modalities: diffuser stepping, and the more recently added sample stepping modality.The latter allows extending the reconstructed field of view (FOV) since it combines frames with different lateral sample positions.We here present "helical stepping", a new SBI acquisition scheme which uses helical tomography to generate multiframe tomographic data compatible with UMPA.Helical stepping allows expanding the FOV while cutting down on acquisition times.

Experimental setup
The measurements were performed at the synchrotron PETRA III in Hamburg (DESY) [17].The imaging beamline P05 uses X-rays from an undulator located d 0 = 85 m from the microtomography endstation.A double-crystal silicon monochromator delivered a highly monochromatic beam (Δ/ ≈ 10 −4 ) with an energy of 20 keV.A CMOSIS CMV 20000 camera from Karlsruhe Institute of Technology (KIT) with a pixel size of 6.4 µm was used, connected to an optical microscope with a 100 µm CdWO 4 scintillator.Using 5× magnification, we obtained a useful area of 2.50 mm × 6.55 mm (1951 × 5120 pixels), with an effective pixel size of 1.28 µm in the sample plane.Spatial resolution was estimated from a slanted-edge MTF measurement, resulting in a value of 2.14 µm (1.67 pixels).The setup's main components are shown in figure 1(a).The only change to the conventional setup was the addition of a modulator on piezoelectric translation stages.
We used 6 layers of 1000-grit silicon carbide sandpaper (mean particle size 5.8 µm) as a wavefront marker for the near-field speckle tomography scans.The sandpaper was positioned d 1 = 115 mm upstream of the sample and the sample-detector distance was d 2 = 175 mm.The main experimental parameters are reported in table 1.
It also included a toothpick and a pipette tip containing polystyrene microspheres (300 nm).

UMPA acquisition schemes
As reported in [16], it is possible to model reference pattern distortions induced by the sample for each pixel r = (, ) to obtain different signals.The model accounts for absorption, small-angle scattering and refraction of the incident X-rays.It is defined as where  (r),  (r) and u(r) = (  (r),   (r)) correspond to the transmittance, dark-field, horizontal and vertical differential-phase signals. 0 (r) is the reference image and ⟨ 0 (r)⟩ is the local mean of  0 in the vicinity of r.The method's performance relies on reference pattern quality and the pattern's sampling over the object features.Hence, by acquiring multiple sample-reference image sets {  }, { 0, } (1 ≤  ≤ ) with different lateral translations between the sample and modulator, we can obtain more accurate reconstructions and improve image quality [9,14].UMPA fits the model to the sample images by minimizing the sum of squared differences in a window Γ(w): where  pairs of images contribute to the final cost function (r).If the sample is moved instead of the diffuser, the difference in sample position s  between the different steps has to be considered, so that identical portions of the sample add to the local minimization and signal retrieval in UMPA.
The set of sample and reference frames in this study were processed using UMPA with a 3 × 3 analysis window after drift correction for each tomography angle.We obtained volumetric reconstructions using conventional filtered back-projection tomography with ASTRA-toolbox [18,19].Differential-phase signals were integrated [20] for each projection before tomographic reconstruction.UMPA can process datasets with two stepping modalities described in the following sections.
Diffuser stepping.Diffuser stepping consists in moving the diffuser transversely to the beam between tomographic scans, while keeping the sample fixed.In practical terms, this means setting s  = 0 in eq.(2.2).A way to combine this with tomographic acquisition, in a stable setup, is to acquire a whole tomographic dataset at a fixed diffuser position, move the diffuser, perform the next scan, and so on.
For the diffuser-stepping acquisition, we collected 3001 angular views of the phantom in a continuous 180 • tomographic scan for each of 20 diffuser positions.20 dark images and 70 reference images were taken before and after each scan.The raw images were dark current and beam profile corrected.Bad pixel outliers were replaced by the median of their closest neighbors, and the beam profile was estimated from low-pass filtered (50 × 50 px kernel size) reference images.
-3 -Sample stepping.In sample stepping, the modulator remains fixed while the sample is moved.Thus, different parts of the sample are superimposed with the same reference pattern.If the sample positions s  (in pixels) are known, it is possible to process the dataset in a way analogous to diffuser stepping.The minimization of  (r) in eq.(2.2) requires the analysis window to move so that the correct portion of the sample is accounted for in each   .This modality allows to extend the reconstructed FOV by intrinsically stitching the displaced sample images during the signal extraction process.However, in the final reconstruction, some pixels receive contributions by less than  frames, and thus have different noise levels.This can be visualized with a frame coverage map [see figure 2(c)].Helical stepping.Helical stepping is a method to reorder and transform data from a single helical tomography scan into a multi-frame sample-stepping dataset compatible with UMPA.This can accelerate the acquisition, but requires determining the exact sample position for each projection.
In the helical-stepping acquisition, the sample was rotated while continuously moving the tomographic stage at a constant speed in the vertical direction by 200 µm per turn.In total, 48000 projections were acquired for 10 full turns of the sample.When considering the subset of projections with an interval of 180 • between them (flipping every other frame horizontally to reverse the effect of the half-rotation), the change between frames is a net vertical movement of the sample without rotation, see figure 1(b).We thus obtained, for each of 2400 different tomographic angles, 20 frames differing only in the vertical location of the sample.

Drift corrections
Motor accuracy and repeatability, along with instabilities in the experimental setup, can introduce vibrations that significantly impact image quality in high-resolution speckle-based scans.Since SBI relies on matching speckle patterns with and without the sample [see figure 3(c) and (d)], it is essential that the pattern remains stable over time in both images to ensure accurate reconstructions.
-4 -Sample motion might occur due to setup stability issues.Usually, multiple reference images would be acquired between every tomographic acquisition (i.e., every diffuser position).The reference pattern can change if the source properties change, since speckles form from free-space propagation interference effects of the sandpaper grains.By using reference images taken closest in time to the actual sample image, we attempt to compensate for any beam instability or temporal variations.However, this requires moving the sample out of the beam after every tomographic acquisition, potentially leading to increased positioning inaccuracies.We here present a method to correct random fluctuations of both diffuser and sample position between individual frames.
We begin with the correction of sample position fluctuations.Eq. ( 2.3) gives an overview of the procedure.Here, {•  } signifies the set of all values with subscripts  = 1, . . ., , and {• , ′ } represents the set of all values with subscripts (,  ′ ) as defined in eq.(2.4).
We assume a complete set of sample images {  } and reference images { 0, }, where either the position of the diffuser (diffuser stepping) or the sample (sample stepping) varies.Note that, in sample stepping,  0,1 =  0,2 = . . .=  0,  .UMPA allows multimodal signal retrieval from a single-shot acquisition.This means that only one pair of images is needed, one with the sandpaper in the beam  0, and another with the sandpaper and the sample in the beam   .The first step is to apply UMPA to all pairs (  ,  0, ) of single frames, yielding a series of transmission images {  }, which should not contain any visible speckles.Differences between the   should therefore only be due to variations in sample position, and differ only in a translation of the image content.We can determine these translation vectors {s , ′ } between any two (  ,   ′ ), e.g., by identifying the peak of the cross-correlation of the two.For diffuser stepping, the resulting translation vector is entirely due to positioning errors.In a sample-stepping dataset however, the sample is intentionally stepped by large distances (see s  in eq.(2.2)), which means that only part of the sample data is shared between any pair of frames.To account for this, we identify the overlap area for every (  ,   ′ ) image pair based on the known {s  }, and use only this area for the calculation of the translation vector.Furthermore, we only calculate the translation vector for the image pairs whose overlap is sufficiently large (here we require at least 50% of the full FOV).This defines a subset  of frame pairs: 1 ≤  <  ′ ≤ , and   and   ′ have more than 50% overlap (2.4) Next, we retrieve an absolute sample position ŝ for every frame  from the relative translation vectors {s , ′ }.To do so, we define a cost function and minimize: The final step consists in shifting both   and  0, by ŝ via sub-pixel interpolation, yielding the images  ′  and  ′ 0, (for all ).Since the impact of sample drifts has been eliminated from these images, they can be used as input for the diffuser drift correction.This procedure is shown schematically in eq.(2.6).
Similarly to the sample drift correction, we apply UMPA to pairs of single frames ( ′  ,  ′ 0, ), but we now focus on the retrieved differential-phase channels u  = ( , ,  , ).Any drift of the diffuser between  ′  and  ′ 0, should result in a uniform offset of u  across the entire image.Consider any pair of frames (,  ′ ) from .The image content of u  and u  ′ in the overlap area should be identical, except for an offset.We thus define u , ′ as the mean of u  − u  ′ across the overlap area.It is equal to the difference in sample-to-reference diffuser drifts between  and  ′ .
Using an approach equivalent to eq. (2.5), we can then calculate the absolute (per-frame) diffuser drifts { û } from the relative drifts {u , ′ }.For the minimization, we use a starting estimate û0, , which we derive via a cross-correlation-based estimate of shift between  0, and  0,0 [see expected positions in figure 2(a)].Finally, we apply the diffuser drift correction by translating the reference frames { ′ 0, } by { û } via sub-pixel interpolation, yielding the frames { ′′ 0, }.A fully corrected dataset, ready to be processed with UMPA, thus consists of the sample-reference sets { ′  }, { ′′ 0, }.

Results and discussion
Reconstructed projections from the diffuser-stepping and helical-stepping datasets are shown in figure 3(a) and (b), highlighting the larger FOV in the latter.Figure 3(e)-(l) show the effect of various corrections on the transmission and differential-phase signals.Figure 3(g) and (h) show improvements due to sample drift correction (drifts of ≈ 18.5 px) compared to figure 3(e) and (f) (no correction).The toothpick exhibits sharper and better delimited features.Figure 3(j) shows additional improvements due to diffuser drift correction.Moreover, figure 3(l) shows that the standard deviation further decreases in an empty region after UMPA bias correction [16].Tomograms for both stepping approaches are shown in figure 4. It is important to mention that, between each of the reconstructed 2400 projections over 180 • in the helical-stepping dataset, there is an offset in height since the first sample position moves vertically as the sample rotates.This offset has to be estimated via cross-correlation and the projections have to be shifted accordingly before tomographic reconstruction in order to yield well-aligned sinograms.Artifacts arising from propagation fringes at material interfaces are evident in all the volumes.The UMPA model accounts only for the first derivative of the phase, thus, second-order phase effects such as fringes are not taken into account by the reconstruction and end up mainly in the attenuation and dark-field tomograms.To mitigate these artifacts, the sample could have been positioned closer to the detector, which would however have reduced phase sensitivity.
Figure 4(b) and (e) show similar  histogram distributions, demonstrating that both approaches yield comparable information.Beam profile and scintillator imperfection corrections are more effective in a sample-stepping acquisition.This is because a consistent portion of the detector's FOV moves during the multi-frame reconstruction process, allowing to average out such issues.This approach also helps in reducing ring artifacts [compare figure 4

Conclusion
In order to develop speckle-based tomography into a routine imaging technique, high acquisition speeds are required.Furthermore, short acquisition times allow capturing dynamic processes in real-time.A continuous helical scan approach is a straightforward choice for extending the detector FOV while also maintaining short acquisition times [21].Spiral computed tomography is a routine procedure in hospitals, where time is crucial [22].However, conventional CT only provides attenuation information with a spatial resolution around 1 mm.Higher-resolution scans are more susceptible to instabilities, and have a smaller FOV because of limitations in beam or detector size.
Here we show that high-resolution, multimodal tomography with an extended FOV is possible but challenging.A similar measurement has previously been performed at a lab-based GBI setup with lower spatial resolution (100 µm) [23].In our case, the reconstruction requires an elaborate analysis of big datasets (> 10 million pixels per frame), making it computationally demanding.Moreover, at high resolution, beam instabilities and/or vibrations can greatly influence the quality of the results.Recent methods make use of deep learning algorithms to improve image quality of smaller datasets from fast scans [24].This could be considered for improving speed and image -7 -  quality in the future.By using UMPA with the proposed drift correction method, we show that a continuous acquisition scheme can already cut acquisition time in half, while maintaining image quality comparable to a standard SBI scan.

Figure 1 .
Figure 1.(a) Experimental setup.(b) Helical acquisition, the phantom was rotated by 10 full rotations to obtain 20 overlapping frames of each tomographic angle with different vertical positions of the sample.

Figure 2 .
Figure 2. Relative diffuser and sample movements for one projection.(a) Diffuser positions for the diffuser-stepping scan.Refined positions indicate displacement after alignment corrections.(b) Sample positions in helical acquisition.We expected net vertical movement, but estimated positions indicate a small drift.(c) Frame coverage map for the helical acquisition.The FOV is extended (axes in pixels), but the number of contributing images to the signal retrieval is not constant across the FOV.More frames are used in the central region (up to 20 images).The black dots mark the movement trajectory of the sample during the scan.

Table 1 .
Experimental parameters for the standard SBI acquisition and the continuous helical acquisitions.