Parallel-beam tomography systems at synchrotron facilities have limited field of view (FOV) determined by the available beam size and detector system coverage. Scanning the full size of samples bigger than the FOV requires various data acquisition schemes such as grid scan, 360-degree scan with offset center-of-rotation (COR), helical scan, or combinations of these schemes. Though straightforward to implement, these scanning techniques have not often been used due to the lack of software and methods to process such types of data in an easy and automated fashion. The ease of use and automation is critical at synchrotron facilities where using visual inspection in data processing steps such as image stitching, COR determination, or helical data conversion is impractical due to the large size of datasets. Here, we provide methods and their implementations in a Python package, named Algotom, for not only processing such data types but also with the highest quality possible. The efficiency and ease of use of these tools can help to extend applications of parallel-beam tomography systems.
© 2021 Optical Society of America under the terms of the OSA Open Access Publishing Agreement
A major drawback of a parallel-beam tomography system is that the reconstructable area  is limited by the beam size of the X-ray source and the size of the detector, which together determine the field of view (FOV). In practice, it is not always possible to prepare samples to fit to the FOV such as valuable fossils [2,3] or cultural heritage specimens . In such cases, different scanning schemes to extend the FOV are needed. There are three basic data acquisition schemes which can be used: a grid scan [5,6], a 360-degree scan with offset center-of-rotation (COR), known as a half-acquisition scan [7,8], and a helical scan [9,10]. Each of these scanning schemes requires different pre-processing methods before the data can be reconstructed using the traditional approach, which reconstructs a slice from a 180-degree sinogram. In a grid scan and half-acquisition scan, image stitching methods are required. In a helical scan, data need to be converted to an equivalent circular scan dataset. With high-flux X-ray sources at synchrotron facilities, acquiring a tomographic dataset can be done in minutes, seconds or even faster  making in-situ and time-resolved experiments possible. Scientists therefore need to reconstruct a few slices quickly for assessment and decision making about the next step of their experiment. In some cases, users would like to reconstruct slices from a whole sample to locate the region-of-interest (ROI) for subsequent, higher resolution scans. This requires a processing pipeline that is easy-to-use, reliable, and automated.
Although the aforementioned data acquisition schemes are easy to implement, not many applications employing these techniques at synchrotron facilities can be found in the literature [12,13]. One of the main reasons is that there is a lack of data processing methods and open-source software to process such data. Some open-source software [14,15] can process a half-acquisition scan but users must input the COR manually based on visual inspection. Tomosaic software  can be used for a grid scan but is currently limited to processing of 180-degree scans. Additionally, metadata about spatial positions of the scans are required for auto-stitching methods to work. For helical parallel-beam scans, we found no Python open-source software to reconstruct this type of data.
The scope of this work is to provide the research community with data processing methods and open-source software for advanced scanning schemes, which are: easy-to-use, automatable, and capable of the highest possible reconstruction quality. Specifically, we introduce:
- - Methods for a half-acquisition scan which: 1) enable fully automated determination of the COR or auto-conversion of a half-acquired sinogram to a 180-degree sinogram, and 2) allow reconstruction directly from a half-acquired sinogram instead of using an extra step of converting to a 180-degree sinogram.
- - A grid scan with the rotation axis offset with respect to the FOV of the grid as a whole, to reduce the number of sample translations and make use of the methods for auto-stitching and auto-conversion of sinograms. For this type of scan, a sample detection method is proposed because a sample can be occasionally outside the FOV. This method allows automation of the processing pipeline.
- - A detailed procedure for data acquisition and calibration for a helical scan. This is important to make sure that the acquired data is reconstructable. For reconstruction, a step-by-step algorithm to convert the data to an equivalent circular acquisition data is presented. We also propose a method for removing a characteristic artifact caused by unresponsive detector pixels for this type of data.
- - A Python software package, named Algotom (Code 1, Ref. ), which implements the methods introduced above, including a full data processing pipeline of pre-processing methods, reconstruction methods, and post-processing methods. It includes implementations of our methods published before, improvements of methods published in the literature, and some developments not presented elsewhere such as a zinger removal method and a double-wedge filter. This package has been being developed in-house and on-demand at Beamline I12, Diamond Light Source (DLS) UK , to improve tomographic data quality, make use all of beamline capabilities, and process tomographic data acquired by non-standard scans.
Considering that high-speed experiments and large datasets are possible using synchrotron-based tomography systems, the proposed methods and methods implemented in Algotom have been developed for working in the sinogram space to reduce the computational cost and to parallelize the computational work. To process sinograms independently, the information recorded in each sinogram must be independent. This requires that a parallel-beam tomography system must be aligned, i.e the rotation axis of a sample stage is parallel to the plane of the detector and perpendicular to the rows of the pixel-grid detector . Alignment is a standard routine at synchrotron facilities [19,20]. If the independence condition cannot be met due to image rotation or radial lens distortion of a detecting system, Algotom provides methods for extracting tilted sinograms and unwarped sinograms with low computational cost.
Details of the methods are presented in Section 2. It is divided into two sub-sections: a section about methods used for processing a half-acquisition scan, a grid scan, and a helical scan; and a section about methods in the data processing pipeline implemented in Algotom. Section 3 show results of applying these methods to process two datasets: one acquired by using a grid scan with offset COR; and one acquired by using a helical scan with offset COR. Data used for this paper were collected at Beamline I12 and Beamline I13-2 of DLS .
2.1 Pre-processing techniques for data acquisition schemes extending the FOV
2.1.1 Methods for a 360-degree scan with offset center-of-rotation
220.127.116.11 Determination of the overlap-side and the overlap-area between sinograms
In a half-acquisition scan, a common routine is that a 360-degree sinogram is converted to an equivalent 180-degree sinogram (Fig. 1(a) and 1(c)), by stitching two halves of the 360-degree sinogram, before reconstruction. Image stitching is a well-established research area in computational photography. There are many choices of algorithms where feature-based approaches have been proven to be fast, robust, and fully automated [22–24]. However, applying these methods to sinogram images of X-ray tomography is not straightforward due to three major challenges. Firstly, the feature space of sinogram images is often not rich enough to find similar features as can be seen in Fig. 1(b) (the blue frames). Secondly, overlap areas are often chosen to be small to maximize the FOV, e.g., 10% to 20% of the area of a sinogram. This strongly affects the reliability of both feature-based and Fourier transform-based methods . Lastly, there are regions of free space beyond the object (Fig. 1(d)), the features in the beam profile and optical artifacts are often very similar and thus highly correlated (Fig. 2(b)). They can cause incorrectly determined direction of stitching of the sinograms by automatic correlation-based methods (Fig. 1(e)).
Here we use a correlation-based method with a technique proposed to solve the problem of wrong-side image-stitching caused by free space. A routine for determining the overlap-side and overlap-area between sinograms is as follows, where the first sinogram is used as the reference; the overlap-side, left or right, is referred to this sinogram:
- - Select a region at the outermost left side of the 2nd sinogram (the blue frame in the sinogram numbered 2 in Fig. 1(b)). The width of the region is insensitive, but it must be smaller than the overlap-width. The same size region in the 1st sinogram is selected and the correlation metric between two regions is calculated. Repeating the calculation for each position of the selection window across the 1st sinogram (the blue frames in the sinogram numbered 1 in Fig. 1(b)) results in a list of metric values (Fig. 2(a)).
- - Apply parabolic fit to 3 points around the minimum-metric point for each of the resulting metric lists (Fig. 2(c) and 2(d)) where the metrics are normalized by dividing to the minimum value before the fitting. The fitted result which has the higher quadratic coefficient (i.e., Fig. 2(c)) corresponds to the sharper minimum. This determines the overlap-side between two sinograms. The overlap-area is calculated from the frame position giving the minimum metric (the green dot in Fig. 2(c)).
The routine described above performs 1D searching, which has low computational cost, across the full width of a sinogram. This enables automated processing. The reason for using a quadratic coefficient to determine the overlap-side is that it represents the sharpness of the true overlap correlation (Fig. 2(c)); a small mismatch in the overlap destroys the correlation whereas accidental strong correlation due to low contrast features within space region persists over a broader range (Fig. 2(d)). The robustness of the algorithm can be improved by applying some processing steps such as denoising using a smoothing filter and normalizing intensities between two regions before calculating the correlation metric. Details of these implementations are in Code 1 (Ref. ).
For automation to work there must be a certain level of overlap between recording areas which can be coarsely adjusted, given the FOV and the sample stage position, before data are acquired. In our experience, an overlap of 10% to 20% of an image width of 2k pixels works for most cases. In the example  shown in Fig. 1 and Fig. 2, the width of the sinogram is 2.8k pixels. The calculated overlap-area (or overlap-width) is 418 pixels. The width of the search window is 100 pixels but insensitive. A larger window size increases the reliability of the result. However, it must be smaller than the overlap-area which can be estimated from data acquisition settings.
Here the correlation metric used is the inversion of the Pearson coefficient , a common metric to measure linear correlation between two images. Other metrics available in Python open-source software such as Spearman, Kendall’s tau, or point-biserial correlation coefficient , can be used. However, we found that the Pearson coefficient is reliable enough for our data. The inversion means that the minimum metric corresponds to the highest correlation. This is done to be consistent with other metrics used by our software (Ref. ) where the minimum values are the best values. The metric is calculated as
If the overlap-side is the left side of the 1st sinogram,
If the overlap-side is the right side of the 1st sinogram,
From the determined overlap-side and the overlap-area, a linear ramp to blend and stitch sinograms is used. Details of these implementations are shown in Code 1 (Ref. ). Although the above routine is demonstrated on sinograms, it certainly can be used on projection images. However, this may reduce the reliability of the method because projection images may exhibit repeated patterns, e.g., projections of a rectangular grid.
18.104.22.168 Reconstruction without converting to a 180-degree sinogram
A common approach to processing data from a half-acquisition scan is that two halves of a 360-degree sinogram are stitched to form a 180-degree sinogram. Then, it can be reconstructed using a standard reconstruction method. This approach works for most of the cases. However, stitching a sinogram suffering from intensity-fluctuation problems can result in sharp changes of intensities in the overlap-area which may be enhanced by the ramp-filter in the FBP (filtered back-projection) reconstruction method . This can give rise to ring artifacts or streak artifacts. Here we introduce a very simple approach where there is no need to use image stitching or conversion of a 360-degree sinogram to a 180-degree sinogram. It needs to be used with back-projection-based reconstruction methods such as the FBP method or the BPF (back-projection and filtering) method .
The pre-processing steps apply linear ramp weighting row-by-row across the overlap area of a 360-degree sinogram (Fig. 3(a) and 3(b)), where the negative of the logarithm of the intensities is used; then zero-padding is applied to the result to expand the width of the sinogram to size of (2×W −OA) (Fig. 3(c)). Details of these implementations are in Code 1 (Ref. ). Back-projection in the full range of 360 degrees is applied to the expanded sinogram to reconstruct a slice as can be seen in Fig. 4.
2.1.2 Methods for grid scans
22.214.171.124 Combination of grid scans and half-acquisition scans with image stitching applied in the sinogram domain
Using a grid scanning technique to scan a sample larger than the FOV of a detecting system is an easy implementation from the data acquisition perspective. However, to reconstruct the whole sample efficiently in term of storage usage and computational cost is challenging. Data processing steps commonly used in X-ray micro-tomography are flat-field correction, zinger removal, ring artefact removal, distortion correction, denoising, and tomographic reconstruction . They need to be combined properly with an image stitching step to produce a full reconstructed sample.
In most synchrotron-based tomography systems, the rotation speed of a sample stage is much higher than the translation speed, so a common data collection routine is to perform a full tomographic fly-scan at each translation position instead of doing translation scanning at each angle position. Different grid-scanning schemes [6,13,16] can be used depending on how the rotation axis is positioned with respect to the FOV and a sample. This dictates how a data processing pipeline should be constructed.
In the first approach, sub-volumes of a sample are scanned and reconstructed independently then the stitching is performed on the reconstruction space. Scanning each sub-volume can be done using a traditional 180-degree scan with the rotation axis at the middle of the FOV or with a half-acquisition scan . A disadvantage of this approach is that the intermediate reconstruction results of each scan need to be stored for later use in the stitching step. More importantly, stitching reconstructed images is not straightforward. Reconstructing from truncated projections introduces artifacts [31–33]. Additionally, the low-frequency components and the DC-component of slices between overlapped regions may vary depending on the contribution of sample parts outside the FOV of these slices. This affects the quality of the blended images.
In a synchrotron-based tomography system, the detector position is often fixed relative to the source. In the previous approach, each sub-volume of the sample is translated in 3D to the reconstructable volume, while the position of the rotation axis is unchanged with respect to the FOV. In the second approach, the rotation axis is translated horizontally in 1D relative to the FOV and the sample is translated vertically in 1D to scan at different heights (Fig. 5(a)). Thus, the translation dimension is 1D less than the first approach. To process data acquired by this type of scanning, stitching is performed in projection space or sinogram space before reconstruction. Data can be acquired in the range of [0; 180°] with the rotation axis positioned at the horizontal center of the grid as shown in Fig. 5(a) and used in . Here we propose a modification of this approach by acquiring data in the range of [0;360°] (Fig. 5(b)) to reconstruct the same volume but with a half number of translations. This scanning technique is simply derived from a half-acquisition scan where the rotation axis is positioned with respect to the FOV of the grid as a whole.
To avoid storing intermediate results, i.e., stitched data, the data processing pipeline is conducted in the sinogram domain. If there is radial lens distortion  it needs to be corrected prior to the auto-stitching method or the results will not be reliable. We developed a correction routine, accompanied by Python software [35,36], where a reference pattern is imaged to characterize the distortion produced by the detector system. Our approach is to extract information from this pattern to calculate parameters of a distortion correction model. From the parameters determined, an undistorted sinogram is generated using a few adjacent rows of projection images (Code 1, Ref. ).
A demonstration of how to reconstruct a slice of a full sample acquired using our approach is shown in Fig. 6. 360-degree sinograms of each tomograph (Fig. 6(a) and 6(b)) are stitched to form a single 360-degree sinogram. This sinogram is converted to a 180-degree sinogram (Section 126.96.36.199), then the FBP reconstruction method is applied resulting in the slice in Fig. 6(e). It can also be reconstructed directly from the sinogram in Fig. 6(c) using the method described in Section 188.8.131.52. Python implementations of these steps are in Code 1 (Ref. ). A practical advantage of the proposed scanning technique is that the pre-processing steps of COR determination and image stitching can all use the same method described in Section 184.108.40.206. This method allows search of the whole width of a sinogram with a low computational cost. It thus enables a fully automated data processing pipeline. Details of the data processing pipeline to reconstruct a whole sample using our proposed data acquisition method are presented in Section 3.
220.127.116.11 Sample detection
Depending on the shape of a sample there are parts of a grid scan where no portions of the sample are present (e.g., the 2nd row of the 1st column of the grid in Fig. 5(a)). For automating a processing pipeline, one needs to detect these areas and avoid using them for determining the overlap-side and overlap-area between cells of the grid. A method to indicate if there is a sample in an image is needed. It can happen that a part of a sample may stay inside the projection images at certain angles but outside at other angles, so using sinogram images for detection is more reliable. Here we propose a simple detection technique by exploiting an interesting feature of a sinogram is that its 2D Fourier transform exhibits a double-wedge region of significant coefficients (Fig. 7(a) and 7(c)). This feature shows not only in a 360-degree sinogram but also in a 180-degree sinogram or a sinogram with the COR out of the FOV. A mathematical derivation and thorough discussion about this feature are fully presented in Ref. .
Knowing that significant Fourier coefficients of a free space image do not exhibit the same feature as a sinogram (Fig. 7(b) and 7(d)), our technique aims to detect the boundary of the double-wedge region. If such a boundary exists, the sinogram contains a sample. Due to the symmetry of the Fourier space, using a quarter of an image is enough for detecting the boundary. The detection algorithm works as follows:
- - Apply 2D Fourier transform to a sinogram (Fig. 7(c)).
- - Extract the top-left quadrant of the resulting image excluding a few rows and columns around the middle (∼20) to remove zero-frequency components. This is to avoid the impact of these values on the step for locating the boundary using maximum values.
- - Apply a strong Gaussian filter, e.g., a sigma of 20, to the logarithmic scaling of the previous image (Fig. 7(e)). This is to remove local maxima for the next step.
- - Locate the position of the maximum value of each row at the bottom of the resulting image (∼ 50% of total rows).
- - If the distance between the maximum point of each row and the right edge of the image consistently increases when moving upwards (Fig. 7(e)), there is a sample in the sinogram.
2.1.3 Methods for a helical scan
18.104.22.168 Calibration and data acquisition
In a helical scan, a sample is rotated and translated along the rotation axis during the scanning process. To process data acquired using the helical scan, a proposed approach is to convert the data to an equivalent circular scan data before reconstruction. It is important here to determine the pixel size precisely with respect to the translation increment; taking into account that the X-ray beam is slightly divergent in a synchrotron-based tomography system, and the mechanical error of the translation stage needs to be smaller than the pixel size. To do that, we scan a dense sphere as shown in Fig. 8. The sphere is segmented and its center of mass (COM) is calculated. Using the positions of the COM of the sphere (Fig. 8(c)), referred to the coordinate of an image, in correlation to the translation positions of the stage we can determine the pixel size and evaluate the accuracy of the translation stage.
The routine for collecting helical data is demonstrated in Fig. 9 and described as follows:
Suppose that one would like to scan a sample with a length of (ye−ys) where ye and ys are determined by moving the translation stage to see the bottom and the top of the sample, respectively. To be able to reconstruct the top and the bottom of the sample, it needs to be inside the FOV at least in the angle range of [0, 180°]. This is done by shifting the translation stage to the position of12,38], Pt needs to be large enough (> 50% of the FOV).
22.214.171.124 Converting from helical acquisition data to equivalent circular acquisition data
Figure 10 demonstrates the process of generating a circular sinogram from helical scan data. As the sample is translated, to extract 1D-projections of the same slice of the sample one must extract data at different rows of the collected images. Details of the conversion routine are as follows.
Suppose that one needs to generate a circular sinogram at the height of
Knowing that at each movement the stage is translated by
In this projection, the row position of the slice yk is given by
If j0 is not an integer, the 1D-projection is extracted by interpolation between the nearest rows. The routine is repeated for all projections in the range of [i0; i0 + N180] resulting in the circular sinogram as shown in Fig. 10(c). Note that the corresponding angle range must be used for reconstruction instead of resetting to the range of [0°; 180°]. Otherwise, the reconstructed images are rotated between slices. In general, one can extract all projections of the slice yk which is in the range of [i0; i0 +FOV/ystep] and use them all for reconstruction. This helps to improve the signal-to-noise ratio of the reconstructed image. Python implementations of the mentioned steps are in Code 1 (Ref. ).
Equation (12) shows that the precise determination of the pixel size, Δx, in relation to the stage translation is crucial to the quality of the extracted sinogram. If the pixel size cannot be determined using the above calibration routine, an alternative way is to calculate the total number of vertically shifted pixels between a projection at its counterpart at 180° or 360°. This requires that the pitch needs to be smaller than the FOV, to be able to record the same part of the sample after 180° or 360° rotation. If this approach cannot yield a good result, one can generate multiple sinograms using estimated pixel sizes and measure the quality metrics of reconstructed images  to find the best value.
126.96.36.199 Ring artifact removal and blob removal
A significant advantage of using a helical scan is that it helps to clean ring artifacts caused by partially defective pixels and scintillator defects of a detecting system (Fig. 11) [12,38]. It can also remove low-frequency ring artifacts caused by the flat-field correction from an X-ray source having an unstable beam profile, or from a strongly absorbing sample recorded by a low dynamic range camera. However, in some cases where flat-field images are noisy there are still some ring artifacts around the middle of the reconstructed image (Fig. 11(c)). These ring artifacts are caused by the flat-field correction process and can be cleaned (Fig. 11(d)) using the sorting-based method, algorithm 3, in .
Another problem is that unresponsive pixels can give rise to streak artifacts in a reconstructed image as can be seen in Fig. 12. To remove these artifacts we adapt algorithm 4 and algorithm 6 as used in  for removing ring artifacts caused by unresponsive pixels in a circular scan. First, a binary blob map (Fig. 13(a)) is generated by applying algorithm 4 to each row of a flat-field image (Fig. 12(a)) where the intensity profile of each row is normalized by dividing with the median-filtered version of itself. Then, the intensities inside the blobs are interpolated from neighboring pixels using the blob mask in combination with Eq. (12). Details of these implementations are in Code 1 (Ref. ). As can be seen in Fig. 13(b) and 13(c), all the streak artifacts are removed.
2.2 Data processing pipeline and the Algotom package
In the tomographic reconstruction process, a raw dataset goes through a series of data processing steps , commonly called a data processing pipeline, before the final result can be used for post analysis. The pipeline is designed to satisfy two basic requirements: the reconstructed images show the least artifacts and the measured refractive indices are closest to those calculated. As shown in Fig. 14, the reconstruction process is divided into three stages: pre-processing, reconstruction, and post-processing. In each stage, many data processing methods can be used where each of them may work in different spaces: the projection space, the sinogram space, or the reconstruction space. In the pre-processing stage, there are different ways of building a pipeline to optimize the computational cost and the quality of the reconstructed images. It is important that pre-processing methods are chosen to work in the same space, because switching between the projection space and the sinogram space requires intermediate storage of data.
Algotom development was started at the I12-JEEP beamline in 2014 as Python codes to process data acquired by the beamline’s large FOV detector, which uses two imaging sensors to cover a rectangular FOV . Images from these cameras must be stitched, using the same approach as shown in Section 188.8.131.52, before tomographic reconstruction can take place. Data processing methods for improving the quality of tomographic data; removing artifacts caused by imperfections of hardware components; making use the beamline capabilities; processing data acquired by non-traditional scanning techniques; and automating data processing pipeline have been actively developed at I12 over the years. These methods have been used internally by I12’s users and refined further for publication and sharing with the research community through open-source software [36,41]. Specifically, methods for automated determination of the center-of-rotation , distortion correction , and removing ring artifacts  have been contributed to Tomopy  and Savu , which are software frameworks for high-throughput processing of tomographic data.
In contrast to Savu and Tomopy, which are optimized for speed, Algotom is a package of data processing algorithms and tools which are designed to be easy-to-use and easy-to-deploy prototype methods. The development of Algotom has focused on pre-processing methods which work in the sinogram space to reduce computational cost. Methods working in the projection space such as phase filter, distortion correction, or rotation correction have been adapted to work in the sinogram space. Figure 14 shows the variety of data processing methods which can be added to the pipeline depending on the set-up and the hardware issues of a tomography system. The following sections present data processing methods available in the Algotom package and commonly used at the I12 beamline. The software is the implements of the methods presented in the previous section, methods developed by us earlier [18,34,40–42], methods published in the literature [28,43], and some of our methods not published elsewhere. It provides users with a complete set of tools to process data acquired by the data acquisition schemes for extending the FOV, aiming for the highest quality possible.
2.2.1 Pre-processing methods
184.108.40.206 MTF (modulation transfer function) deconvolution
The loss of resolution caused by the scattering of scintillation photons in a scintillator-based detector system is reduced by determining the MTF of the scintillator using a calibration image and following the routine described in . Then, raw images are deconvolved by the division in the Fourier domain.
220.127.116.11 Flat-field correction
The intensity ratio of the image to the incoming beam at each pixel is found simply by the division of an image-with-sample by an image-without-sample, i.e., flat-field image, where a dark-field image (camera’s dark noise) is subtracted from each image before the division.
18.104.22.168 Zinger removal
Zingers are prominent bright dots in projection images caused by scattered X-rays hitting the detector system CCD or CMOS chip (Fig. 15(a) and 15(b)). They produce line artifacts across a reconstructed image (Fig. 15(c)). There are some ways to remove them such as normalization-based techniques [44,45] or statistical methods [15,30]. They may work better in the projection space than the sinogram space or vice versa. Here we propose a unified method which can work in both spaces. It can be used on raw images or flat-field-corrected images. Zingers in flat-field image, if not removed, cause stripe artifacts in a sinogram and ring artifacts in a reconstructed image due to the flat-field correction process. These artifacts can be removed using the corresponding removal methods .
Our approach is a normalization-based method. However, in contrast to other techniques [14,45], which use a smoothing filter for producing a smoothed version of the original image, we use the average of the images shifted around the original position but not including the original image (Fig. 16(a) and 16(b)). The shifting step is the same size as a zinger. In the normalized image, obtained by dividing the original image by the average of the shifted versions, the zingers are located using a threshold parameter (Fig. 16(c)). Then, their values in the original image are replaced by the values in the averaged image, only at these locations with size equal to or smaller than the given size (Fig. 16(d)). A major advantage of our approach is that it does not significantly alter the image, even when used to remove low-contrast zingers.
22.214.171.124 Ring artifact removal
Ring artifacts may come from the irregular response of a detecting system caused by the non-linear response of the sensor chip, dust on the optics components, and defects in the scintillator. Our work on classifying types of ring artifacts and removing them efficiently were published in  and available in open source software such as Sarepy , Tomopy , and Savu . In Algotom, including with these methods we also provide supplementary tools for customizing ring removal methods (Code 1, Ref. ).
Another source of ring artifacts is introduced by the flat-field correction process. It is caused by features in a flat-field image which are uncorrelated to projection images. In Fig. 17, there are features (arrowed) in the flat-field image caused by the X-ray optics . During tomographic scanning, these features change over time, so that they cannot be removed using the flat-field correction technique (Fig. 17(c)). As a result, there are very large ring artifacts (Fig. 17(d)) which are difficult to remove without altering the image. In such cases, using a helical scan is a good solution.
126.96.36.199 Distortion correction
The distortion correction routine is used to correct radial lens distortion of a lens-coupled detector system. Distortion parameters are calculated from a reference pattern using our Vounwarp software [35,36]. From the parameters determined, projection images are undistorted by a backward mapping process. Directly generating an undistorted sinogram using a few adjacent rows of projection images is also available in Algotom.
188.8.131.52 Tilted sinogram generation
The tilt correction routine is used to correct the problem of projection-image rotation but working in sinogram space using a few adjacent sinograms.
184.108.40.206 Center-of-rotation determination
A method developed by us and published in  is used to calculate the center-of-rotation from a 180-sinogram. A method for a half-acquisition sinogram presented in Section 2.1.1 is also implemented in Algotom.
220.127.116.11 Double-wedge filter
The previously described double-wedge property of the 2D Fourier transform of a sinogram (Section 18.104.22.168)  can used to remove sample parts outside the FOV of a sinogram. In fast tomography experiments [11,46], the number of projections collected is often insufficient for direct reconstruction methods to work well . As a result, iterative reconstruction methods, mainly linear algebra-based, are preferred. Because these methods construct and solve a set of linear equations from a sinogram, they require that the sums of intensities of each sinogram-row are the same. This condition is not satisfied in the case of a sample  larger than the FOV as can be seen in Fig. 18(a). It shows a sinogram of the sample where there are parts coming in and out the field of few of a single row of a detector. Applying iterative reconstruction methods [47–50] on this sinogram results in various types of artifacts as shown in Fig. 18(b)–18(d). The slice from the FBP method (Fig. 18(e)) is not affected by the problem. However, as the number of projections is insufficient, the slice is noisy.
We can reduce the effect of extended objects on the performance of linear algebra-based methods by removing parts of a sample larger than the width of its sinogram. This can be done by masking the 2D Fourier transform of this sinogram. The double-wedge region in the Fourier space is symmetric and well-defined if the sinogram used is 360-degree (Fig. 19(a)). Because a sinogram is often collected in a 180-degree range in parallel-beam tomography, forming a 360-degree sinogram from this sinogram can be done  by: flipping the 180-degree sinogram from left to right; shifting this image to the COR; then appending the shifted image vertically to the original sinogram (Fig. 19(a)).
The double-wedge property has been exploited to calculate the COR , sinogram restoration , or sample detection as shown in Section 22.214.171.124. A mathematical explanation of the formation of a double-wedge region is fully presented in Ref. . Here we only use a formula which defines the wedge-angle of a double-wedge region to design the mask. The angle of a wedge (Fig. 19(c)) is calculated by13). This means that only Fourier coefficients of the sample inside the FOV, i.e., the width of the sinogram, are kept. Because altering coefficients in the Fourier space results in global change of an image in the real space and there is frequency leakage, applying the double-wedge mask needs to be repeated a few times, e.g., 5-10 times.
- - Convert a 180-degree sinogram to a 360-degree sinogram.
- - Apply 2D Fourier transform to the 360-sinogram (Fig. 19(b)).
- - Multiply the result with a double-wedge binary mask (Fig. 19(c)).
- - Perform inverse 2D Fourier transform.
- - Repeat the last 3 steps a few times (e.g., 5-10 times).
- - Remove the bottom half of the resulting sinogram (Fig. 19(d)).
The result of applying the double-wedge filter to the sinogram in Fig. 18(a) is shown in Fig. 19(d). As can be seen the external parts of the sample are removed. However, there is a side effect that the low-frequency components are altered as seen by comparing the difference between the sinograms before and after the filtering (Fig. 19(e)). We can reduce this effect by applying normalization to the processed sinogram (Fig. 19(d)) before reconstruction. Figure 20 shows improvements of the iterative reconstruction methods after the double-wedge filter is applied. Strong artifacts disappear as compared to Fig. 18.
The double-wedge filter can also be used to remove horizontal stripes in a sinogram which may be caused by the fluctuation of an X-ray source (Fig. 21(a) and 21(d)) or by the extinction contrast of crystalline samples (Fig. 21(b) and 21(e)) . Note that these horizontal stripes do not have significant impact on the reconstruction quality of the FBP method but they may affect iterative reconstruction methods.
126.96.36.199 Sinogram stitching
This module of methods is used for processing data acquired using a grid scan, a half-acquisition scan, or their combination. It includes the methods of determining the overlap-side and the overlap-area between sinograms, as shown in Section 188.8.131.52, and the method for stitching them. In practice, it may happen that there are no overlaps between images. In such a case, Algotom also provides methods for joining sinograms based on interpolation.
184.108.40.206 Helical sinogram transformation
This module is used for a helical scan as presented in Section 220.127.116.11.
18.104.22.168 Denoising or contrast enhancement
In a synchrotron-based tomography system, the partial coherence of the X-ray source introduces a propagation-based phase contrast effect resulting in edge-enhanced images as shown in Fig. (1), Fig. (3), or Fig. (6). For applications which only use absorption-contrast information, this effect is unwanted as it causes streak artifacts in reconstructed images (e.g., Fig. (4)). To reduce the effect, a low-pass filter can be used to remove the high-frequency edges in the images. One of the well-known filters was published by Paganin et al. , implemented in many tomographic software packages [14,15]. As a low-pass filter it reduces noise and the dynamic range of an image. This helps to enhance the contrast between low-contrast features which can be confused if this enhancement comes from the phase effect. To properly retrieve phase shift from edge-enhanced images, a class of methods known as the propagation-based methods need at least two measurements at different sample-detector distances [54,55]. Paganin’s method is considered to be a phase retrieval method if it is applied to the image of a single homogeneous material . However, the advantage of phase-contrast imaging over absorption-contrast imaging is only shown in a multi-material sample with small-density differences [56,57]. Considering Paganin’s formula which is rewritten as [53,58]
As R, i.e., δ/β, is undetermined for a multi-material sample, users must manually select this value as a trade-off between the contrast enhancement and the sharpness of a reconstructed image. In practice, R can be arbitrarily adjusted in the range of 10 to 104. Because of that, the constant of (λz)/(4πΔx2) has no meaning and can be omitted which gives a much simpler formula
The adjustable range of the R parameter is not changed much in Eq. (18) compared to Eq. (17) because the value of (λz)/(4πΔx2) is about between 0.1 and 1 in the near-field hard X-ray regime with a micro-resolution detector. The 1D version of this filter for working in the sinogram space is given as
This 1D filter is useful to denoise or improve the contrast of reconstructed images from under-sampled sinograms in grid scans because applying Eq. (18) to all projection images before sinogram generation is very computationally expensive. Implementations of the filters in Eqs. (18) and (19) are shown in Code 1 (Ref. ).
2.2.2 Reconstruction methods
22.214.171.124 Direct reconstruction methods and edge-padding for removing cupping artifacts
Large-size images and high number of datasets limit the choices of reconstruction methods when users have to process their data in a reasonable timescale. Here, we mainly use the FBP method  or the DFI (direct Fourier inversion) method , depending on the availability of GPU resources. The under-sampling problem, caused by applying these methods to data having a low ratio between the number of projections and the image width, can be improved using denoising methods before or after reconstruction [53,59].
In FFT (fast Fourier transform)-based reconstruction methods, such as the FBP or DFI method, the 1D FFT is applied to a sinogram row-by-row before another processing step, e.g., interpolation, is used. FFT requires that an input is periodic, which is not the case for a sinogram of a sample larger than the FOV. This non-periodicity leads to artifacts in the processed sinogram. Zero-padding is often used to reduce these artifacts . However, the sharp changes in intensities at padded edges introduce errors in the Fourier domain which give rise to cupping artifacts in reconstructed images (Fig. 22). Published techniques can be used to tackle the problem, known as the truncated projection problem [31–33]. They mainly utilize extrapolation techniques to compensate for missing data which are not straightforward to implement. Here we use a simple edge-padding technique to avoid the problem caused by sharp changes in intensities at the edges. It removes the cupping artifacts in the slices of both the FBP method and the DFI method as can be seen in Fig. 23. This technique may alter the DC-component of a sinogram. However, this variance is not significant and can be corrected by normalization. The technique is already used elsewhere for the FBP method  but not for the DFI method. Here, we demonstrate this good practice because some tomographic open-source software currently uses built-in zero-padding for the FFT-based reconstruction methods [14,47,61]. To show the advantage of using edge-padding, we implement the FBP method and the DFI method ourselves in the Algotom package (Ref. ).
2.2.3 Post-processing methods
126.96.36.199 Ring artifact removal
Ring artefacts still present in reconstructed images may be reduced by transforming reconstructed images between the polar coordinates and the Cartesian coordinates and applying FFT-based filters to remove ring artifacts. Routines for these procedures are included here.
188.8.131.52 Image downsampling and image rescaling
Methods for downsampling a 2D or 3D image dataset include the mean, the median, the maximum, or the minimum value of a pixel grid or a voxel grid. Methods for rescaling 32-bit data to 16-bit data or 8-bit data are also provided. It is important to be aware that rescaling causes information loss. The global extrema or user-chosen percentile of a 3D data or a 4D data (time-series tomography) need to be used for rescaling to limit the loss.
3.1 Using a grid scan in combination with half-acquisition scans
Datasets used for demonstrating this combined scanning technique were collected at DLS, Beamline I12 with the following conditions: Monochromatic X-ray beam energy 76 keV. Imaging module 3 was used with a PCO.edge 5.5 sCMOS camera having a 2,560 × 2,160 pixels sensor. The optics module gives a pixel size of 3.2 µm resulting in a FOV of 8.2 mm wide and 6.9-mm high. To cover a whole sample, which is a pebble having a size of ∼3.6 cm × 5.1 cm × 2.2 cm (width × height × depth), a grid scan of 3 columns by 8 rows was performed giving 24 tomographic datasets in total. Each tomographic dataset contains 3601 projections evenly distributed in the range of [0°; 360°]. The rotation axis was offset to increase the width of the stitched FOV to approximately 6 times (2 × 3 image widths) of the FOV of the detector. Linear translations of the sample were selected to get the overlap between views of ∼0.8 mm, which is equivalent to ∼250 pixels. The scanning was done row-by-row to avoid the mismatch between sinograms of the same sample slice, which may be caused by the mechanical error of the vertical translation stage. Raw datasets were stored in HDF (Hierarchical Data Format) files accompanying by NeXus metadata files . A critical advantage of using the HDF file format for tomography is that subsets of a 3D dataset can be extracted quickly and efficiently. This is important for a data processing pipeline because different methods in a pipeline may work in different slice dimensions.
Figure 24 shows how the data are processed. The following methods are used for removing artifacts: flat-field correction, zinger removal, and ring artifact removal methods which are the combination of algorithm 3, 4, 5, and 6 in . Sinograms are stitched and converted to 180-degree sinograms using the method described in Section 184.108.40.206. The CORs are determined using Eq. (5). The denoising method with R = 300 is applied (Eq. (18)). Finally, the images are reconstructed using the FBP method.
In practice, synchrotron users often reconstruct a few slices of a sample for quick decision-making during an experiment, or when optimizing parameters of data processing methods prior to performing full reconstruction of the sample. As shown in Fig. 24, reconstructing a few slices from data acquired using the combination of a grid scan and half-acquisition scans is straightforward and fully automated using the methods described in Section 2.1.1. Specifically, reconstructing a slice of the sample shown in Fig. 24 from the raw data takes around 2.5 minutes on a normal workstation at Beamline I12 having 16GB RAM memory and Intel Xeon E-2124 CPU. The computational cost of the auto-stitching and auto-centering steps is most expensive taking nearly 2 minutes. The pre-processing step of flat-field correction and the gridrec reconstruction  cost only 0.5 minutes, giving the reconstructed slice with a size of ∼14k × 14k pixels. One needs to calculate the overlap-areas and overlap-sides for stitching sinograms only once. After that time for reconstructing another slice is below 1 minute.
For reconstructing the whole sample, the approach is slightly different. Firstly, overlap-areas and overlap-sides for stitching and converting sinograms are calculated column-by-column and row-by-row between cells in a grid. This is necessary to determine the global width of reconstructed images. Because the overlap-areas may change between the rows of the grid, possibly caused by the mechanical error of the translation stage, the widths of the reconstructed images may vary correspondingly. The global width chosen is the maximum value among the rows of the grid. Sinograms are then padded before reconstruction to compensate for the mismatch. For stitching in the vertical direction, the overlap-sides and overlap- areas between the grid rows are calculated using the method in Section 220.127.116.11, which is applied to the projection images at the same angle between the grid rows. Sinograms at the top and bottom of two grid cells at the same column can also be used for determining the overlaps, but this approach is more costly due to the overhead of extracting sinograms. Knowing the overlaps, sinograms of the same slice are combined using the distance-weighted average before reconstruction as shown in Fig. 24.
Although reconstructing the whole sample into a massive single HDF file is straightforward, it is not practical for post-analysis or transferring the file. Here, we divided the full sample into 12 sub-sections and reconstructed them independently using 12 GPU nodes of the DLS Hamilton cluster. It took nearly a day for processing all the data using all the processing methods described above. The shape of the 3D reconstructed image (32-bit) is of ∼15k × 14k × 14k pixels occupying 11.3 TB of the storage system. The 3D rendering of the 3D reconstructed image, from the Avizo software, is shown in Fig. 25(a) and 25(b) where the image volume was downsampled by a factor of 8 to make the visualization possible.
3.2 Using a helical scan in combination with a half-acquisition scan
A dataset used for demonstration this type of scanning was also collected at the beamline I12 using the same setup as in Section 3.1. The monochromatic X-ray beam energy was 54 keV. Parameters of the helical scan are: ystart = 3.75 mm; ystop = 19.5 mm; Pt = 6.5 mm; and N180 = 1,801. The pixel size determined using the routine in Section 18.104.22.168 is of 3.24 µm. The rotation axis was offset to approximately double the FOV. The slice thickness, Δy, is chosen to be equal the pixel size. Figure 26 shows some projections of the tomographic data where Fig. 26(c) is created by combining the last columns of all projection images. It helps to show the direction of extracting sinograms for reconstruction and the areas (the blue overlays) where the data have no use or are not reconstructable. Figure 26(c) also shows another possible way of generating sinograms which includes two steps : applying a shearing transform in the same plane as in Fig. 26(c); extracting columns of the resulting images. This approach, however, requires storing of the intermediate data and its performance relies on the IO performance of a storage system.
Data processing methods used in the pipeline (Fig. 27) for this type of scan are: flat-field correction, circular sinogram generation including the blob removal (Section 2.1.3), partial ring artifact removal (algorithm 3 in ), zinger removal, 180-degree sinogram conversion using the auto-centering and auto-stitching methods in Section 22.214.171.124, and FBP method for reconstruction. A bottleneck in the processing pipeline is the step for generating the circular sinograms where the program walks through the 3D dataset and extract subsets of the data. Here, the time cost was around 10s per sinogram for this dataset on the DLS Hamilton cluster.
4. Discussion and conclusion
Although there are many data processing methods presented in the paper and they are all important in each scanning technique and data processing pipeline, we would like to highlight a few important contributions of our work.
Firstly, there is a method for determining the overlap-side between sinograms. This method leads to the determination of the overlap-area between sinograms of a grid scan and/or the center of rotation in a half-acquisition scan. Practically, the method searches the full width of a sinogram in both left and right directions with low computational cost. This enables automated data processing. Explanations of the feasibility of the proposal method, in comparison to other approaches, for tomographic imaging are shown in Section 126.96.36.199.
Secondly, there is a combination of a grid scan with half-acquisition scans. This scanning technique allows a multiple doubling of the FOV of a parallel-beam tomography system with the least number of stage translations. Practically, it makes use of the proposed method for both sinogram stitching and 360-degree sinogram conversion. This enables to reconstruct a slice of the grid scan fully automatically with low computational cost.
Thirdly, it is a detailed procedure of calibration, data acquisition, and data conversion for helical scans. Although it is easy to perform a helical scan, processing the acquired data is very cumbersome without a proper calibration and acquisition. The pixel size of a detector in relation to the translation of a stage needs to be calibrated as X-ray beams are not perfectly parallel. The start and stop position of the translation stage need to be calculated properly to be able to reconstruct the region of interest. Details of the conversion algorithm and the methods of removing characteristic artifacts are also provided. All these contributions allow users or beamline scientists to perform helical scans with ease.
Lastly, it is software which not only implements all the proposed methods but also common methods in a data processing pipeline. It provides users with a complete set of tools to reconstruct data acquired by the scanning techniques described in this paper with ease and the highest possible quality. They all help to extend applications of parallel-beam tomography systems.
As stated about the scope of this work, we have demonstrated data processing methods and software (Ref. ) developed for reconstructing data acquired by the scanning techniques, which extend the FOV of a parallel-beam tomography system, including half-acquisition scans, grid scans, helical scans, and their combinations. Methods in a full data processing pipeline are provided in the Algotom package to give users a complete set of tools for tomographic reconstruction with the highest possible quality. Considering the importance of decision making during an in-situ experiment, a very large data can be produced in a single beamtime, and a high demand of delivering a high-throughput system of data acquisition and data processing at a synchrotron facility; our approaches were designed to be easy-to-use, automatable, and efficient; e.g., a large slice of a grid scan can be reconstructed fully automatically in a few minutes using a normal workstation computer as shown in Section 3.1.
Although the methods and the software have been developed and tested on parallel-beam X-ray tomography systems, they should be applicable to other types of parallel-beam tomography such as neutron tomography, X-ray fluorescence tomography, or electron tomography. Certainly, characteristic pre-processing methods of image alignment for electron tomography, scattering correction for neutron tomography, or self-absorption correction for fluorescence tomography need to be added to the software. To distribute the work and further development, the software is also made available at: https://github.com/algotom/algotom. Details of how to use the Algotom package to reconstruct images from the raw data presented in Section 3.1 can be found in the “examples” folder of Code 1 (Ref. ). Links to the raw data uploaded to zenodo.org are provided inside the codes.
This work was carried out with the support of the Diamond Light Source. Works related to the grid scan were performed for beamtimes: MG22109-1, “Textual contrast enhancement from multimodal tomographic morphology,” Brent Seales, 2019; SW24740-1, “Reading a Herculaneum scroll via X-ray tomography,” Brent Seales, 2019. Works related to the helical scan were performed for beamtimes: EE14033-2, “X-ray high energy phase and dark field contrast imaging using speckle based technique,” Hongchang Wang, 2016; MT19526-1, “Densitometry: a fundamental study of electron density mapping and calibration,” Rhian Jones, 2018; MG21999-1, “The effect of gravity on the differentiation dynamics of rocky planets,” Matt Pankhurst, 2019. This research used a resource of the National Synchrotron Light Source II, a U.S. Department of Energy (DOE) Office of Science User Facility operated for the DOE Office of Science by Brookhaven National Laboratory under Contract No. DE-SC0012704. The authors would like to thank Mark Basham and Nicola Wadeson for fruitful discussions.
The authors declare that there are no conflicts of interest related to this article.
1. A. C. Kak and M. Slaney, Principles of Computerized Tomographic Imaging, (IEEE, 1988).
2. R. J. Garwood, J. A. Dunlop, P. A. Selden, A. R. T. Spencer, R. C. Atwood, N. T. Vo, and M. Drakopoulos, “Almost a spider: A 305-million-year-old fossil arachnid and spider origins,” Proc. R. Soc. London, Ser. B 283(1827), 20160125 (2016). [CrossRef]
3. C. Baars, M. Ghobadi Pour, and R. C. Atwood, “The earliest rugose coral,” Geol. Mag. 150(2), 371–380 (2013). [CrossRef]
4. W. B. Seales, J. Griffioen, R. Baumann, and M. Field, “Analysis of Herculaneum papyri with X-ray computed tomography,” in International Conference on Nondestructive Investigations and Microanalysis for the Diagnostics and Conservation of Cultural and Environmental Heritage, (2011).
5. J. W. Eberhard and K. C. Tam, “Translate rotate scanning method for X-ray imaging,” US Patent 5,032,990 (July 16, 1991).
6. A. Kyrieleis, M. Ibison, V. Titarenko, and P. J. Withers, “Image stitching strategies for tomographic imaging of large objects at high resolution at synchrotron sources,” Nucl. Instrum. Methods Phys. Res., Sect. A 607(3), 677–684 (2009). [CrossRef]
7. P. S. Cho, R. H. Johnson, and T. W. Griffin, “Cone-beam CT for radiotherapy applications,” Phys. Med. Biol. 40(11), 1863–1883 (1995). [CrossRef]
8. S. R. Stock, “Recent advances in X-ray microtomography applied to materials,” Int. Mater. Rev. 53(3), 129–181 (2008). [CrossRef]
9. P. E. Slavin, “X-ray helical scanning means for displaying an image of an object within the body being scanned,” US Patent 3,432,657 (March 11, 1969).
10. Y. Bresler and C. J. Skrabacz, “Optimum interpolation in helical scan computerized tomography,” in Proceedings of IEEE International Conference on Acoustics, Speech and Signal Processing (IEEE, 1989), pp. 1472–1475.
11. M. Drakopoulos, T. Connolley, C. Reinhard, R. Atwood, O. Magdysyuk, N. Vo, M. Hart, L. Connor, B. Humphreys, G. Howell, S. Davies, T. Hill, G. Wilkin, U. Pedersen, A. Foster, N. De Maio, M. Basham, F. Yuan, and K. Wanelik, “I12: The Joint Engineering, Environment and Processing (JEEP) beamline at Diamond Light Source,” J. Synchrotron Radiat. 22(3), 828–838 (2015). [CrossRef]
12. M. J. Pankhurst, N. T. Vo, A. R. Butcher, H. Long, H. Wang, S. Nonni, J. Harvey, G. Guđfinnsson, R. Fowler, R. Atwood, R. Walshaw, and P. D. Lee, “Quantitative measurement of olivine composition in three dimensions using helical-scan X-ray micro-tomography,” Am. Min. 103(11), 1800–1811 (2018). [CrossRef]
13. E. Borisova, G. Lovric, A. Miettinen, L. Fardin, S. Bayat, A. Larsson, M. Stampanoni, J. C. Schittny, and C. M. Schlepütz, “Micrometer-resolution X-ray tomographic full-volume reconstruction of an intact post-mortem juvenile rat lung,” Histochem. Cell Biol. 155(2), 215–226 (2021). [CrossRef]
14. D. Gürsoy, F. De Carlo, X. Xiao, and C. Jacobsen, “Tomopy: a framework for the analysis of synchrotron tomographic data,” J. Synchrotron Radiat. 21(5), 1188–1193 (2014). [CrossRef]
15. N. Wadeson and M. Basham, “Savu: a Python-based, MPI framework for simultaneous processing of multiple, N-dimensional, large tomography datasets,” (2016). https://arxiv.org/abs/1610.08015.
16. R. Vescovi, M. Du, V. de Andrade, W. Scullin, D. Gursoy, and C. Jacobsen, “Tomosaic: efficient acquisition and reconstruction of teravoxel tomography data using limited-size synchrotron X-ray beams,” J. Synchrotron Radiat. 25(5), 1478–1489 (2018). [CrossRef]
17. N. T. Vo, “Algotom: Data processing algorithms for tomography,” figshare (2021https://figshare.com/s/7342960859dd73e98124.
18. N. T. Vo, M. Drakopoulos, R. C. Atwood, and C. Reinhard, “Reliable method for calculating the center of rotation in parallel-beam tomography,” Opt. Express 22(16), 19078–19086 (2014). [CrossRef]
19. N. T. Vo and R. C. Atwood, “Tilt alignment using ImageJ macro,” https://confluence.diamond.ac.uk/display/I12Tech/Tilt+alignment+using+ImageJ+macro.
20. A. Rack, “What to know/do before starting tomo experiment” http://www.esrf.eu/UsersAndScience/Experiments/Imaging/ID19/BeamlineDescription/Experimental-hutch/microtomography/what_to_know#alignxc.
21. C. Rau, U. Wagner, Z. Pešić, and A. De Fanis, “Coherent imaging at the Diamond beamline I13,” Phys. Status Solidi A 208(11), 2522–2525 (2011). [CrossRef]
22. D. Gledhill, G. Y. Tian, D. Taylor, and D. Clarke, “Panoramic imaging—A review,” Comput. Graphics 27(3), 435–445 (2003). [CrossRef]
23. R. Szeliski, “Image alignment and stitching: A tutorial,” Found. Trends Comput. Graph. Vis. 2(1), 1–104 (2007). [CrossRef]
24. H. Joshi and K. Sinha, “A survey on image mosaicing techniques,” IJARCET 2, 365–369 (2013).
25. Y. Chai, Y. Wang, Z. Yousaf, N. T. Vo, T. Lowe, P. Potluri, and P. J. Withers, “Damage evolution in braided composite tubes under torsion studied by in-situ X-ray computed tomography,” Compos. Sci. Technol. 188, 107976 (2020). [CrossRef]
26. K. Pearson, “Contributions to the mathematical theory of evolution. III. Regression, heredity, and panmixia,” Philos. Trans. R. Soc. Lond. Ser. A 187, 253–318 (1896). [CrossRef]
27. P. Virtanen, R. Gommers, T. E. Oliphant, M. Haberland, T. Reddy, D. Cournapeau, E. Burovski, P. Peterson, W. Weckesser, J. Bright, S. J. van der Walt, M. Brett, J. Wilson, K. J. Millman, N. Mayorov, A. R. J. Nelson, E. Jones, R. Kern, E. Larson, C. J. Carey, İ. Polat, Y. Feng, E. W. Moore, J. V. Plas, D. Laxalde, J. Perktold, R. Cimrman, I. Henriksen, E. A. Quintero, C. R. Harris, A. M. Archibald, A. H. Ribeiro, F. Pedregosa, P. van MulbregtSciPy 1.0 Contributors, “SciPy 1.0: Fundamental algorithms for scientific computing in Python,” Nat. Methods , 17(3), 261–272 (2020). [CrossRef]
28. G. N. Ramachandran and A. V. Lakshminarayanan, “Three dimensional reconstructions from radiographs and electron micrographs: Application of convolution instead of Fourier transforms,” Proc. Nat. Acad. Sci. 68(9), 2236–2240 (1971). [CrossRef]
29. G. T. Gullberg, “The reconstruction of fan-beam data by filtering the back-projection,” Comput. Graphics Image Process. 10(1), 30–47 (1979). [CrossRef]
30. R. C. Atwood, A. J. Bodey, S. W. T. Price, M. Basham, and M. Drakopoulos, “A high-throughput system for high-quality tomographic reconstruction of large datasets at diamond light source,” Philos. Trans. A. Math. Phys. Eng. Sci. 373(2043), 20140398 (2015). [CrossRef]
31. B. E. Oppenheim, “Reconstruction tomography from incomplete projections,” in Proceeding of Reconstruction Tomography in Diagnostic Radiology and Nuclear Medicine (University Park, 1977), pp. 155–183.
32. R. M. Lewitt, “Processing of incomplete measurement data in computed tomography,” Med. Phys. 6(5), 412–417 (1979). [CrossRef]
33. R. M. Ranggayyan, A. T. Dhawan, and R. Gordon, “Algorithms for limited-view computed tomography: An annotated bibliography and a challenge,” Appl. Opt. 24(23), 4000–4012 (1985). [CrossRef]
34. N. T. Vo, R. C. Atwood, and M. Drakopoulos, “Radial lens distortion correction with sub-pixel accuracy for X-ray micro-tomography,” Opt. Express 23(25), 32859–32868 (2015). [CrossRef]
35. N. T. Vo, “Python implementation of distortion correction methods for X-ray tomography,” Zenodo (2018) https://doi.org/10.5281/zenodo.1322720.
36. N. T. Vo, “Distortion correction for X-ray tomography,” https://github.com/nghia-vo/vounwarp.
37. P. R. Edholm, R. M. Lewitt, and B. Lindholm, “Novel properties of the Fourier decomposition of the sinogram,” Proc. SPIE 0671, 8–18 (1986). [CrossRef]
38. D. M. Pelt and D. Y. Parkinson, “Ring artifact reduction in synchrotron x-ray tomography through helical acquisition,” Meas. Sci. Technol. 29(3), 034002 (2018). [CrossRef]
39. T. Donath, F. Beckmann, and A. Schreyer, “Automated determination of the center of rotation in tomography data,” J. Opt. Soc. Am. A 23(5), 1048–1057 (2006). [CrossRef]
40. N. T. Vo, R. C. Atwood, and M. Drakopoulos, “Superior techniques for eliminating ring artifacts in X-ray micro-tomography,” Opt. Express 26(22), 28396–28412 (2018). [CrossRef]
41. N. T. Vo, “Numerical techniques for eliminating ring artifacts in X-ray micro-tomography,” https://sarepy.readthedocs.io.
42. N. T. Vo, R. C. Atwood, and M. Drakopoulos, “Preprocessing techniques for removing artifacts in synchrotron-based tomographic images,” Proc. SPIE 11113, Developments in X-Ray Tomography XII, 111131I (2019).
43. R. M. Mersereau, “Direct Fourier transform techniques in 3-D image reconstruction,” Comput. Biol. Med. 6(4), 247–IN4 (1976). [CrossRef]
44. M. Rivers, “Tutorial introduction to X-ray computed microtomography data processing,” http://www.mcs.anl.gov/research/projects/X-ray-cmt/rivers/tutorial.html.
45. J. C. E. Mertens, J. J. Williams, and N. Chawla, “A method for zinger artifact reduction in high-energy x-ray computed tomography,” Nucl. Instrum. Methods Phys. Res. A 800, 82–92 (2015). [CrossRef]
46. S. Hasan, V. Niasar, N. K. Karadimitriou, J. R. Godinho, N. T. Vo, S. An, A. Rabbani, and H. Steeb, “Direct characterization of solute transport in unsaturated porous media using fast x-ray synchrotron microtomography,” PNAS 117(38), 23443–23449 (2020). [CrossRef]
47. W. V. Aarle, W. J. Palenstijn, J. Cant, E. Janssens, F. Bleichrodt, A. Dabravolski, J. De Beenhouwer, K. J. Batenburg, and J. Sijbers, “Fast and flexible X-ray tomography using the ASTRA toolbox,” Opt. Express 24(22), 25129–25147 (2016). [CrossRef]
48. P. Gilbert, “Iterative methods for the reconstruction of three dimensional objects from their projections,” J. Theor. Biol. 36(1), 105–117 (1972). [CrossRef]
49. A. H. Andersen and A. C. Kak, “Simultaneous algebraic reconstruction technique (SART): A superior implementation of the art algorithm,” Ultrason. Imaging 6(1), 81–94 (1984). [CrossRef]
50. M. R. Hestenes and E. Stiefel, “Methods of conjugate gradients for solving linear systems,” J. Res. Natl. Bur. Stand. (U.S). 49(6), 409–436 (1952). [CrossRef]
51. J. S. Karp, G. Muellehner, and R. M. Lewitt, “Constrained Fourier space method for compensation of missing data in emission computed tomography,” IEEE Trans. Med. Imaging 7(1), 21–25 (1988). [CrossRef]
52. W. H. Zachariasen, “The secondary extinction correction,” Acta Cryst. 16(11), 1139–1144 (1963). [CrossRef]
53. D. Paganin, S. C. Mayo, T. E. Gureyev, P. R. Miller, and S. W. Wilkins, “Simultaneous phase and amplitude extraction from a single defocused image of a homogeneous object,” J. Microsc. 206(1), 33–40 (2002). [CrossRef]
54. M. R. Teague, “Irradiance moments: Their propagation and use for unique retrieval of phase,” J. Opt. Soc. Am. 72(9), 1199–1209 (1982). [CrossRef]
55. N. T. Vo, R. C. Atwood, H. O. Moser, P. D. Lee, M. B. H. Breese, and M. Drakopoulos, “A fast-converging iterative method for X-ray in-line phase contrast tomography,” Appl. Phys. Lett. 101(22), 224108 (2012). [CrossRef]
56. A. Momose, T. Takeda, Y. Itai, and K. Hirano, “Phase–contrast X-ray computed tomography for observing biological soft tissues,” Nat. Med. 2(4), 473–475 (1996). [CrossRef]
57. M.-C. Zdora, “State of the art of X-ray speckle-based phase-contrast and dark-field imaging,” J. Imaging 4(5), 60 (2018). [CrossRef]
58. T. Weitkamp, D. Haas, D. Wegrzynek, and A. Rack, “ANKAphase: Software for single-distance phase-retrieval from inline X-ray phase contrast radiographs,” J. Synchrotron Radiat. 18(4), 617–629 (2011). [CrossRef]
59. D. Kazantsev, E. Pasca, M. J. Turner, and P. J. Withers, “CCPi-Regularisation toolkit for computed tomographic image reconstruction with proximal splitting algorithms,” SoftwareX 9, 317–323 (2019). [CrossRef]
60. A. Kyrieleis, V. Titarenko, M. Ibison, T. Connolley, and P. J. Withers, “Region-of-interest tomography using filtered backprojection: assessing the practical limits,” J. Microsc. 241(1), 69–82 (2011). [CrossRef]
61. S. van der Walt, J. L. Schönberger, J. Nunez-Iglesias, F. Boulogne, J. D. Warner, N. Yager, E. Gouillart, T. Yuthe scikit-image contributors, “scikit-image: Image processing in Python,” PeerJ 2, e453 (2014). [CrossRef]
62. M. Konnecke, F. A. Akeroyd, H. J. Bernstein, A. S. Brewster, S. I. Campbell, B. Clausen, S. Cottrell, J. U. Hoffmann, P. R. Jemian, D. Mannicke, R. Osborn, P. F. Peterson, T. Richter, J. Suzuki, B. Watts, E. Wintersberger, and J. Wuttke, “The NeXus data format,” J. Appl. Cryst. 48(1), 301–305 (2015). [CrossRef]