Quantitative estimation of closed cell porosity in low density ceramic composites using X-ray microtomography

X-ray Microtomography is a proven tool for phase fraction analysis of multi-phase systems, provided that each phase is adequately partitioned by some means of data processing. For porosity in materials containing low-density ceramic phases, differentiation between pores and the low-density phase(s) can be intractable due to low scattering in the low-density phase, particularly if small pores necessitate low binning. We present a novel, combined methodology for accurate porosity analysis—despite these shortcomings. A 3-stage process is proposed, consisting of (1) Signal/noise enhancement using non-local means denoising, (2) Phase segmentation using a convolutional neural network, and (3) Quantitative analysis of the resulting 3D pore metrics. This particular combination of denoising and segmentation is robust against the fragmentation of common segmentation algorithms, while avoiding the volitional aspects of model selection associated with histogram fitting. We discuss the procedure applied to ternary phase SiC–TiC-diamond composites produced by reactive spark plasma sintering with porosity spanning 2–9 vol%.


Scientific Reports
| (2023) 13:127 | https://doi.org/10.1038/s41598-022-27114-w www.nature.com/scientificreports/ Non-local means (NLM) is one such filter effective in both denoising and edge preservation, already finding widespread use among the greater tomography community [15][16][17][18] . Following noise reduction, segmentation is the step where pores are computationally separated from the image. Segmentation is the general process whereby objective features in an image are isolated, or otherwise mapped for processing. Conventionally, this is achieved through a basic threshold of the scalar distribution, or by application of the watershed. The threshold approach is popularized by its simplicity, requiring only a single parameter for tuning. However, thresholding is susceptible to over-segmentation and is naturally semantic, unable to individualize particles from the resultant mask (without successive processing). The watershed is innately instanced and similarly prone to over-segmentation but is advantaged by the realization of boundaries between contiguous, identical particles 19 .
The advent of the convolutional neural network (CNN) as a segmentation tool has expanded the breadth of classifiable features to those most nebulous in form. Indeed, CNN-based segmentation has seen considerable adoption in the medical community for numerous biological phenomena, for which the prototypical forms of thresholding and watershed are often ill-poised 20 . Image segmentation using CNNs for microscopy purposes has greatly matured since the inception of the CNN, and there are now multiple cloud-based tools for easily leveraging these algorithms, such as APEER (https:// www. apeer. com) or CDeep3M 21 . These cloud-based tools reduce the barrier to entry for research groups not directly in the machine learning space, in terms of required technical expertise, computing hardware, and the time commitment for a minimum viable result. Here we use the Zeiss APEER tool to demonstrate that manually constructing CNN architectures for the problem is not strictly required.
In this study we explore the combination of edge-preserving NLM denoising with tailorable convolutional neural network segmentation to analyze porosity in noisy, low contrast XRM experiments. This process is discussed in the context of porous SiC-TiC-diamond composites. Different noise reduction routines are compared to NLM, then the effect of CNN segmentation is contrasted to thresholding for these denoised tomograms. Finally, statistical accuracy is discussed in terms of over/under-segmentation along with error analysis of the porosity volume fraction.
Lastly, this combined methodology is used to investigate the pore evolution in reactively spark plasma sintered SiC-TiC-diamond composites. Resulting from the reactive (liquid-phase) process, XRM reveals a pore morphology representative of the bulk microstructure. Application of the NLM-CNN processing methodology herein allows us to more confidently explicate the competing sintering mechanisms affecting the initial ceramic powder blends. We have recently reported in more detail on these sintering mechanisms in 22 . Furthermore, the combined process presented in this paper has implications for other low-density ceramic composites, such as SiC-B 4 C-diamond composites.

Materials and methodology
Materials. Spark plasma sintering (SPS) was used to fabricate intentionally porous SiC and denser SiC-TiCdiamond composites. The SPS system was manufactured by Thermal Technology LLC (Santa Rosa, California). The starting SiC powder (Thermo Fisher Scientific) with D 50 = 36 µm was sintered at a hold temperature, pressure, time, and heating rate of 1600 °C, 5 MPa, 10 min, and 100 °C/min, respectively. For the composites, powder blends of elemental Si, Ti, and TiC-coated or uncoated diamond were densified using SPS at hold temperatures of 1600, 1625, or 1650 °C at a 10 min hold time with heating rate of 100 °C/min at 50 MPa pressure. Elemental Si and Ti melt during sintering at these temperatures, reacting with diamond to form SiC/TiC. The composite microstructures are primarily SiC-diamond matrix with an intergranular TiC structure (see Fig. 2). Diamond is metastable and tends towards graphitization under these SPS conditions. See 22 for further composite synthesis and microstructural details. X-ray Microtomography. Tomographic images were acquired on a Zeiss XRADIA Versa 520 X-Ray Microscope. The source was configured for 80 kV voltage, 7 W power, and a LE6 filter. Datasets were then acquired at 20 × magnification, 2 × binning (~ 0.797 μm effective pixel size), 7 s exposure time, and − 180°-180° angular stage range. Data processing. All XRM experiments are exported from the Zeiss Reconstructor application to uncompressed *.tiff image series in 32-bit floating-point precision, valued between 0 and 1. Image series are reconciled into stacked 3D tiffs (maintaining 32-bit precision) with dimensions 956 × 992 × 964 px using ImageJ 23 . Prior to data processing, the 3D datasets are cropped to 512 × 512 × 512 px cubes positioned at the geometric center of the reconstruction volume (physical dimensions are ~ 408.1 × 408.1 × 408.1 μm or ~ 0.107 mm 3 ). Prior to training/segmentation, datasets are downsampled to 8-bit integral precision. Annotation, training, and segmentation are done using APEER (https:// www. apeer. com). 3D particle statistics are generated from segmentation results using BoneJ 24 .

Results and discussion
The proposed process is summarized in Fig. 1, which is composed of 3 main stages. In the first step, 3D nonlocal means is applied to the raw XRM datasets to improve the SNR. The resulting denoised datasets are used for labeling of pores for training the CNN model in succeeding step. The SNR improvement is intended to have the effect of improving label accuracy by more clearly revealing pores/pore boundaries that are partially obscured by noise, while also reducing the prevalence of noise as a factor which the CNN must "learn". In the second stage, the labels resulting from annotation of selected 2D slices are used to train a CNN using APEER. More specifically, APEER is a Unet-based convolutional neural network architecture using EfficientNet 25  Once training is completed, the CNN is fed whole XRM datasets for segmentation. The resulting segmentation is two-dimensional in nature, as the CNN is trained on cross-sections of pores. In the third step, the layered 2D segments of each pore are "solidified" back into 3D volumes based on their connected-ness. With the availability of each 3D pore structure, the porosity and other pore-based statistics can be easily measured.
Denoising. Improving the signal/noise ratio of an XRM dataset can have a substantial positive effect on the segmentation quality. Unfortunately, denoising algorithms tend to-at least partially-behave like low-pass filters. The suppression of higher spatial frequency features obfuscates vital details like phase boundaries and small particles/pores. As a result, a balance between meaningful noise reduction and the preservation of fine microstructural features is not strictly guaranteed. Figure 2a shows a slice through a porous SiC-TiC-diamond composite, along with various noise reduction filters ( Fig. 2b-e). The median/gaussian filters are computed using 3D kernels, and the bilateral/NLM filters using 3D search/patch distances. The gaussian filter is computed with a 3 px standard deviation and the median filter a 5 × 5 × 5 px kernel. Both the median and gaussian filters are  www.nature.com/scientificreports/ simple and fast, requiring only a single parameter for tuning (excepting truncation of the gaussian kernel). The gaussian filter is effective in denoising, at the cost of blurred boundaries and other small features. The histogram peaks belonging to diamond and SiC have partially separated due to the reduction in noise; however, visual clarity in the image has been severely reduced by erasure of higher spatial frequencies. The median filter is effective in the preservation of boundaries and smaller features; however, the achievable level of noise reduction leaves more to be desired. Figure 2d-e shows equivalent results using a bilateral and non-local means filter, respectively. The gaussian, bilateral, and NLM filters are each intimately related, in that each are based on spatial and/or scalar (gaussian) weighted averaging of the neighborhood surrounding each pixel. Given some pixel coordinate p , the (gaussian) filtered pixel I G p is the average of each other pixel coordinate q in a neighborhood around p, weighed by a gaussian function f G with standard deviation σ: where �p − q� is the spatial distance between pixels p and q . The size of is usually determined by the number of standard deviations by which the kernel is truncated. As mentioned, the gaussian filter does not have edgepreservation properties. The bilateral filter compensates for this by introducing a second gaussian weighting function r G that weighs against the difference in scalars I p and I q : where 1/C p is a normalization factor. This extra weighting function (against scalar distance) allows the bilateral filter to reduce the weight contribution from increasingly dissimilar pixels, giving the bilateral filter moderate edge preservation. This can be seen in Fig. 2d where noise reduction is achieved while maintaining phase boundary clarity. Instead of strictly weighing against spatial or scalar distance, NLM weighs against the similarity of the local neighborhood between p and q 28,29 . This is achieved using a new weighting function w , parameterized by the Euclidean distance between each neighborhood ( p and q ), and a standard deviation h , controlling the permissiveness in weighing contributions from dissimilar neighborhoods: The σ 2 value is an optional noise variance estimation, ensuring that any two neighborhoods within the noise variance are equally weighed (see 29 ). The noise variance could be estimated by the median absolute deviation (MAD), wavelet shrinkage 30 , or any other number of estimators. The effect of each denoising filter on the scalar histogram, and an individual slice is shown in Fig. 3a and Fig. 2, respectively. In Fig. 3a, the visual quality of denoising is most evident by the distinctness of the SiC and diamond peaks. I.e., no material components of the (1) www.nature.com/scientificreports/ unfiltered histogram are distinct (and the pores are entirely obfuscated). NLM gives the most visually definite peaks, due to both noise reduction and edge preservation. While NLM is quite effective, its drawbacks come in the form of parameter estimation and computation time. While not exceedingly difficult, the patch size (size of p and q ), search neighborhood (size of ), and smoothing parameter h (standard deviation of s G ) must be manually estimated. The computational complexity of NLM is O NP d R d , where N is the total number of voxels, P the patch size, R the size of the search neighborhood, and d the dimensionality. For our needs, we developed a fast multithreaded NLM implementation (C++ with a graphical user interface), with compiled binaries available 31 . However, those with access to Nvidia (or AMD ROCm-supported) graphics processing units (GPU) may find it prudent to use implementations eliciting GPU acceleration 18 . We created this software for quick, interactive parameter tuning and to lower the knowledge barrier for non-programmers.
Not all software supporting 3D NLM implement 3D neighborhoods/searching, many implementations accumulate 2D NLM filters in a paginated fashion over some axis of the 3D image. Paginated filtering will significantly reduce the computation time; however, the level of noise reduction will likewise diminish. For non-quantitative applications (i.e., XRM for purely imaging purposes), the speed gained by paginated filtering likely outweighs the boost in visual acuity. Figure 3b compares the effect of 3D NLM to 2D NLM for a SiC-TiC-diamond XRM composite. The distinctness of the SiC/Diamond peaks is best using 3D NLM; however, the other 2D NLM results are still quite favorable. For highly anisotropic microstructures, the denoising quality may be partially sensitive to the axis chosen for pagination. We find that the processing time for 3D NLM is approximately 10 × slower for the datasets and parameters used here. However, for quantitative applications-like porosity estimation-it may be prudent to favor denoising quality to minimize extraneous false positives/negatives generated during segmentation (i.e., see Fig. 7). The best choice likely depends on the data quality on a case-by-case basis, for this study we use 3D NLM for best quality denoising. There are further adaptations to the NLM filter to reduce execution time, based on Fast Fourier Transforms (FFT) 32 , omission of dissimilar neighborhoods 33 , probabilistic early termination 34 , or 35 ; all with varying levels of denoising quality.

Segmentation challenges.
Mapping the voxel-based representation of an XRM dataset into specific pore and non-pore (or pore and phasic) regions is the main pivotal step in porosity analysis. A robust method to classify regions of the dataset is critical, otherwise, phase fraction results become error-prone and reliant on manual intervention or correction by the individual. For well-behaved datasets, the segmentation step can be especially trivial if the density of each phase is well distributed. In these ideal cases, simple thresholding may be sufficient to wholly classify each phase. In this case, thresholding refers to segmentation by partitioning the dataset's scalar distribution into fixed bins characterizing each phase. The "thresholds" define the boundaries of each bin.
The shaded regions in Fig. 4a show the result of thresholding on an intentionally porous SiC specimen (4 × binning), with respect to the dataset histogram and the inset image slice. There is a sufficient difference in density and recurrence among pores and SiC, such that peaks arising from pores are isolated from non-pores. Due to the peak isolation, a well-estimated threshold value could be placed around the valley at ~ 0.52, resulting in a porosity fraction of ~ 40.36 vol%. Non-linear least squares (NLLS) fitting of the histogram components is sometimes effective in estimating phase fractions, given a suitable model. Using a split-Voigt model for each peak results in a mediocre fit (orange for pores and blue for non-pores) and the porosity fraction is computed as the relative area of the pore model: ~ 42.92 vol% porosity. Still, there is a ~ 6% relative difference between the two methods arising from uncertainty in threshold placement and least squares parameter error (imperfect model). With such a high porosity content, this error is not too significant. Figure 4b shows a similar histogram from an XRM experiment on a diamond-SiC-TiC composite (2 × binning). Unlike the porous SiC experiment, the diamond (~ 3.5 g/cm 3 ), SiC (~ 3.16 g/cm 3 ), and pore features are materially inseparable in their given state. The increase in noise is partially due to lower binning (needed to resolve smaller pores), causing the noise to increase by a factor of 8. Additionally, cutting diamond-based samples is non-trivial and can result in thicker than ideal XRM specimens, further increasing noise. Finally, the amount of porosity is much smaller (mostly < 6 vol%) resulting in a weaker, buried peak. Attempts at thresholding must be "eyeballed" to balance out the effects of erroneous false positives and negatives (see inset in Fig. 4b). Still, the result is inexact and difficult to visually confirm, depending on the accuracy needed. Further, we could not learn more about the pore distribution due to the heavily fragmented nature of the segmentation (i.e., what is the average pore size?). Also, without a physically meaningful model, NLLS deconvolution of the pore peak would likewise be hopeless.
Thresholding and NLLS fitting are not the only methods of isolating porosity information in tomographic datasets. Indeed, there are a multitude of algorithms for both instanced and semantic segmentation. These two methods, however, capture the essence of the ineffective separation of pores from non-pores in low-density ceramic materials (or otherwise noisy and low-contrast phases). The watershed algorithm-and its variationsmay be the most traditional method used for segmentation in microscopy. If the dataset is thought of as a topographical map, the watershed treats the valleys as "basins" that fill with water. As the basins "fill" and begin to touch, the boundaries between basins can are retained. Resolving the boundaries (watersheds) between connected particles (basins) is the strength of watershed segmentation: improving particle analysis by the separation of contiguous particles that would otherwise be considered singular. For pore "particles", however, there is not so much physical basis for the subdivision of contiguous pores.
Application to porous SiC-TiC-diamond composites. Twenty slices from the 3 most porous composites were chosen for annotation and training. Chosen slices were well separated within each dataset, such that no pore was resampled by any other annotated slice. Pore boundaries among the 20 slices were carefully   www.nature.com/scientificreports/ drawn for all sizes (small to large pores), then non-pore regions surrounding each pore were annotated (it is important to include pore boundaries in the background). An example slice is shown in Fig. 5a with annotated pore and background regions. Small, medium, and large pores are all represented in the annotation, and each pore boundary is included in the background regions. Every pore in a slice does not need to be annotated, we annotated ~ 80-90% of pores per slice. The composites are highly isotropic, and the axis chosen for annotation should not have an effect. The intersection-over-union (IoU) evolution of the model is shown in Fig. 5b with nominal IoU in the range 0.88-0.89. Figure 6a shows the porosity fraction for the six composites resulting from the described analysis process, parameterized by the diamond powder (coated versus uncoated) and SPS hold temperature. The porosity of the uncoated diamond is constant with temperature when accounting for the estimated error. For TiC-coated diamond, the porosity is higher with a minimum at 1625 °C, where the porosity is in the vicinity of uncoated diamond. The higher porosity in TiC-coated diamond is attributed to delamination of the TiC coating as the underlying diamond surface graphitizes, resulting in a lower bulk density. Note that the size scale of graphite produced in this process (confirmed by transmission electron microscopy in 22 ) is well below the detectable size for XRM.
For statistics of the pore distribution to be viable, there must be a minimal amount of falsely classified pores (over/under-segmentation). Much of the potential for over/under-segmentation is mitigated by the denoising step (see the effect of noise on thresholding in Fig. 4). Elevated over-segmentation will result in elevated pore number density. Figure 6b conveys the difference in number density between thresholding and CNN segmentation, for each SiC-TiC-diamond composite. To make the comparison most sensible, the threshold value was chosen so the ensuing porosity fraction was equivalent to that measured by the CNN. Specifically, a scalar histogram of 1024 bins was computed for each dataset, and the bin value closest to the CNN porosity is chosen; the largest relative error was − 0.0083% in the 1650 °C TiC-coated diamond composite. In all cases, Fig. 6b consistently shows that the number density produced by CNN segmentation is consistently 20-40% lower than by thresholding, cursorily suggesting that over-segmentation is reduced using the CNN. To further explore whether this could be a result of under-segmentation of small pores by the CNN, Fig. 7 examines the visual effects of segmentation by both methods. For clarity, Fig. 7a shows the un-segmented regions. Figures 6c and 7b show the segmentation in equivalent regions (of the same slice) by thresholding and CNN segmentation, respectively. In the upper image, thresholding results in incomplete segmentation (false negatives) of the large pore, and over-segmentation (false positives) of diamond/SiC in the lower image. However, the CNN segmentation correctly classifies each image with neither over-nor under-segmentation.

Estimation of uncertainty. Ideally, error estimations should accompany point estimates of the porosity
fraction. An error estimate establishes a frame of reference for comparison among similar point estimates of different specimens or regions. The systematic error introduced by postprocessing XRM datasets through means of denoising or segmentation is not so straight forward to estimate. However, we can estimate the statistical error associated with the measured porosity distribution. I.e., the number, size, and spread of the pore distribution should affect our confidence in the point estimate being genuinely representative of the whole. Consider two specimens of equal porosity X p = 50 vol%: one containing a single enormous pore and the other thousands of finely dispersed and disconnected pores. In this example, we should have more confidence in the estimate containing numerous small pores, as opposed to that of a single large pore. In the contrived example of a single massive pore, we have little footing to expect whether a follow-up experiment would contain another pore, or two, or any at all. Moreover, we are confident that a follow-up experiment on the specimen with numerous finely dispersed pores should likely exhibit a similar distribution if the specimen is homogenous.
If the pore distribution is known, the porosity fraction can be expressed in terms of the mean pore volume V , number of pores n , and the volume spanning the reconstruction :  where the δ operator indicates the standard error of the mean. The reconstruction volume is exact, and the number of pores should be close to exact-provided the segmentation process does not generate needless false positives. If these assumptions hold, δX p becomes, from which δX p and subsequent confidence intervals can be bootstrapped without further assumption. More concisely, if δV is taken as δV = s V / √ n , with s V being the standard deviation in pore volume, the standard error of the porosity fraction can be directly expressed as: This procedure for estimating error should only apply to closed-cell type porosity. For open-cell porosity where a single pore could feasibly span the entire field of view, we cannot effectively sample the pore metrics, invalidating Eq. (6). In such circumstances, it may be prudent to pursue gaseous adsorption-based characterization vectors.
Limitations. Compared to conventional denoising filters, non-local means requires a much longer execution time. The execution time can be considerably reduced by forgoing 3D neighborhoods in favor of 2D, at the slight expense of noise reduction (Fig. 3b). The reduced level of noise reduction may-or may not-be substantive when contrast is low. Parameter estimation may also require time if one is unfamiliar with NLM; this was partial motivation for the interactive NLM software we have provided.
Similarly, annotation and training of the neural network is a non-trivial time commitment, the number of annotations required is all but certain at the outset. With improvements like transfer learning 36 , or data augmentation such as elastic deformation fields 37 , the annotation commitment is being constantly reduced. Alternatively, annotation can be considered a tuning-knob of sort, as the segmentation quality may always be improved by simply more annotations-conventional segmentation methods rarely expose a vector for incremental improvement. However, the segmentation quality is not solely determined by the quantity of annotations, but by their quality too. Unfortunately, annotation is a human-driven process and is therefore captive to human bias.

Conclusions
A novel combined post-processing methodology for quantitative porosity analysis is described within the context of X-Ray Microtomography. Comprised of three parts: (1) Signal/noise ratio enhancement using non-local means, (2) Pore segmentation using a convolutional neural network, and (3) Statistical analysis of the pore (4) X p = Vn www.nature.com/scientificreports/ distribution, porosity metrics are robustly computed-despite elevated noise and low contrast exemplified by low-density ceramics. Porosity volume fraction estimates were computed for reactive spark plasma sintered SiC-TiC-diamond composites from 1600 to 1650 °C with, and without TiC coatings on the diamond particles (other porous low-density ceramics should work similarly well). The resulting TiC-coated diamond porosity was measured 4.0-9.0 vol%, and uncoated from 2.7 to 3.4 vol%, based on sintering temperature. Both over-and under-segmentation was significantly reduced by segmentation using a convolution neural network on nonlocal means pre-processed XRM datasets. The NLM-CNN method resulted in a 20-40% (based on temperature) reduction in estimated pores when compared to thresholding.

Data availability
All data generated or analyzed during this study are included in this published article.