X‐ray Micro‐CT: How Soil Pore Space Description Can Be Altered by Image Processing

Abstract

A physically accurate conversion of the X-ray tomographic reconstructions of soil into pore networks requires a certain number of image processing steps. An important and much discussed issue in this field relates to segmentation, or distinguishing the pores from the solid, but pre- and post-segmentation noise reduction also affects the pore networks that are extracted. We used 15 two-dimensional simulated grayscale images to quantify the performance of three segmentation algorithms. These simulated images made ground-truth information available and a quantitative study feasible. The analyses were based on five performance indicators: misclassification error, non-region uniformity, and relative errors in porosity, conductance, and pore shape. Three levels of pre-segmentation noise reduction were tested, as well as two levels of post-segmentation noise reduction. Three segmentation methods were tested (two global and one local). For the local method, the threshold intervals were selected from two concepts: one based on the histogram shape and the other on the image visible-porosity value. The results indicate that pre-segmentation noise reduction significantly (p < 0.05) improves segmentation quality, but post-segmentation noise reduction is detrimental. The results also suggest that global and local methods perform in a similar way when noise reduction is applied. The local method, however, depends on the choice of threshold interval. Characterizing the soil's physical properties and understanding the resulting functions of the soil is of major importance for many agricultural and environmental issues. The soil is at the interface of most physical, chemical, and biological processes. In this regard, there is increasing interest in the use of noninvasive X-ray microtomography to obtain a microscopic three-dimensional view of the inner soil pore space (for a full description of the technology, see 23). Several reviews (42; 12; 47) have discussed the use of X-ray microtomography in soil and hydrological sciences. In these fields, the technique has been used at both the core scale (e.g., 9; 17; 8; 25; 37; 5; 20; 10; 24; 19) and the aggregate scale (e.g., 28; 36; 32; 21; 49) for describing the pore space and studying the impact of land use and agricultural management on soil structure (9; 28; 17; 36; 32; 25; 5; 21; 10; 49), as well as for analyzing the relationships between soil pore networks and soil physical properties (8; 20; 24; 19). Flow simulations in observed pore networks (7) or a similar constructed pore network (45) have also been conducted. These analyses assumed that the pore space description generated from the image processing accurately represents the physical reality of the sample microstructure, but the choice of X-ray computed tomography (CT) image processing methodology has a visible impact on the resulting structure. Figure 1 shows an example of the processing steps from sample acquisition to binary image. Each step involves choosing the appropriate method and parameters, which are numerous and can have a profound effect on the resulting structure. These choices ultimately depend on the experience of the operator. Processing steps from sample to binary image. Some sources of variability are written in lowercase. What is important here is not only the diversity of these choices but also the fact that they are often inadequately described or justified. Table 1 shows an example of the diversity of methodologies used in a selection of soil science research papers (selection based on number of citations and diversity of research teams). Within Table 1, the pre-segmentation and post-segmentation steps are differentiated. Pre-segmentation steps are varied and are more efficient at handling image degradation than post-segmentation processing; a general rule (for more than just image analysis) is that the more upstream a problem is corrected, the easier is it to process the data downstream. Segmentation is the essential step when pixels are assigned to either the solid or porous phase. There are numerous segmentation methods; a review of those used in soil science was provided by 44. In this study, we differentiated global and local thresholding methods. The aim of a thresholding method is to select a grayscale value, manually or automatically, that separates the image gray levels into two groups: greater than or equal to the threshold (TH) and less than the TH. In soil science, these two groups are often defined as the solid phase (soil matrix) and the void phase (pore space). With global thresholding, a constant TH is chosen for the entire image, whereas with local thresholding, the value is computed for every pixel based on the local neighborhood (44). Segmentation precision depends on the initial quality of the grayscale images. Enhancing the projections before reconstruction and the reconstructed images before segmentation is the typical approach, but each research team has its own procedures (see Table 1). An efficient method for enhancing image quality is to apply noise reduction filters (18; 47) as mentioned in six of the 15 studies listed in Table 1. Some researchers have shown (4; 22; 36; 34; 43) that, in most practical cases, the choice of segmentation method plays a crucial role in the resulting pore structure, but no standards have yet been proposed. Several studies have sought to classify thresholding techniques based on information available from the resulting binary images (2; 14; 15; 39). So far as we know, only 46 have used synthetic soil aggregate images, from which ground-truth information was available, to compare thresholding methods. Even these studies were based on image-by-image analyses and did not provide a tool with which to properly evaluate the processing methodologies. Within this context, our study sought to provide a statistical analysis of the segmentation processing effects on the resulting data. By evaluating Otsu's global method (30), the local adaptive-window indicator kriging (IK) method (13), and the porosity-based (PBA) global method (4) on two-dimensional simulated soil images from which ground-truth information was available, we could also objectively support existing reviews. The first objective of our study was to quantify the effects of pre-segmentation noise reduction on the accuracy of the thresholding method based on the performance indicators. The second objective was to evaluate the impact of post-segmentation processing on pore functionalities. The third objective was to propose an approach for calculating the initial TH interval necessary using the local IK method based on the global TH calculated using the PBA method (4), considering that IK is sensitive to the initial choice of TH interval (15; 40; 46; 13). Here we focus initially on the construction of our simulated images. The general framework was based on the methodology described by 46. It involved superimposing a realistic binary pore image (real soil [RS] images) on an image representing partial volume effects and then adding Gaussian noise (see Fig. 2 for a detailed illustration). We created 15 simulated images from the combination of 15 selected RS binary images and 15 generated partial volume effect images using a method based on fractals and the method of 46. The thresholding methods tested should identify the pore region from the original RS image. Detailed illustration of the simulated image construction. The RS images were derived from the 3 study. We selected 15 two-dimensional images from silt loam soil. Details about the materials, sampling, and X-ray acquisition parameters can be found in 3. Reconstructions were performed using NRecon software provided free of charge by Brucker micro-CT. This software provides tomographic artifact correction methods, which were not tested in this study. Automatic misalignment compensation was used, along with a Level 7 (out of 20) ring artifact correction. The RS images were not subjected to a beam hardening correction. In X-ray microtomography, the most commonly cited artifact is beam hardening due to the polychromatic nature of the X-ray beam, implying a deviation from the Beer Lambert law. For cylindrical objects, it results in a radial grayscale intensity variation from the edges to the center. The beam hardening effect is barely distinguishable from the circular compaction that occurs when sampling soil, and removing beam hardening effects might create noise. Finally, an intensity rescaling was applied to increase contrast (44). The partial volume effect images were generated through the overlaying of decreasing resolution images, as proposed by 46. Our addition to Wang's method was to produce decreasing resolution fractal images with a fractal dimension calculated from the RS images' fractal dimension (Steps A and B). Those images were then combined to form one partial volume effect image (Step C). where N is the constant number of transformed elements at each iteration and r is the ratio between the dimension of the parent element and the dimension of the transformed element. Because power-law dependencies have been observed in soil science, researchers have applied fractal geometry to the study of soil behavior (31). For example, 38 successfully derived a soil-water retention curve from the pore-size distribution fractal dimension of a silt loam soil. Many studies have reported that this concept provides a good description of the complexity of soil microstructure (e.g., 21). where d is the Euclidian dimension. To construct partial volume effect images, we first generated decreasing resolution fractal images. This process was based on two main steps (A and B), each consisting of three fractal iterations. Step A involved fractal images to be used as the of the constructed partial volume effect image and by pixels in Fig. These pixels could not be the of the Step involved pixels to the solid and pore the pixels of Fig. were the pixels subjected to fractal Step A in the of decreasing resolution fractal images for 1, by iteration from to construction of the partial volume effect The pixels the soil The fractal dimension in the Step A was as Step involved pixels to the solid and pore and the pixels in Fig. are the pixels subjected to fractal The three of Step (see Fig. were with and Step in the of decreasing resolution fractal images for 1, by iteration from to The pixels The partial volume effect construction the method of 46 was applied to the generated fractal images Step B), and Fig. shows the first image to be a of by the image was into an by image by calculating the of 2 by 2 The image was then on Fig. which also been from by to by by adding the pixel to the effect of of This created image was to by using the on Fig. and the resulting image was then to by by The is shown in Fig. Figure the impact of the process between Fig. and Fig. partial volume effect constructed by the fractal for 1. The image is an of Fig. in the The image is an of a of Fig. in the The step involved the overlaying of Fig. on the RS image derived from 3. We then noise to the and and and were calculated from our noise by the and The step was to Gaussian noise to the image to noise simulated grayscale image for 1. to an in its the effect of a pre-segmentation one two was tested on the segmentation quality of the simulated images. filters the value of the pixels to the These filters are less sensitive to and no grayscale value is created the resulting in the edges (44). The use of a before segmentation to be a step in the field of soil X-ray image The complexity of segmentation is to the and partial volume effects in the grayscale images. sources of image degradation ring and enhancing image quality, choosing the segmentation method is The global thresholding method described by was tested it provides results and can be used in to the IK method where there are distinguishable histogram its and the of more and more efficient methods, it is a used method for soil images, it is and to It a TH based on the of the between two intensity of In our study, it was performed with we ground-truth information available, the TH that should be applied could be an the TH that the between calculated and ground-truth was and this value as a This was based on the method described by The was provided by the we to the method as the PBA The IK method has provided good results in studies 14; 36; 15; the adaptive-window indicator kriging method (13), was chosen that the method than the one similar The IK concept on the selection of a TH to grayscale are at and are at 1. The between and are assigned to a a on grayscale value and The adaptive-window IK method this based on available information to the when The method was applied using the The choice of and was based on using the method an available we to the method as The PBA method was shown to be its performance was than that of The in IK is the choice of the to interval. proposed an TH interval selection method, it sensitive to noise. We sought to the physical of the PBA method with the IK The aim was to select a TH interval based on the global PBA threshold and then the IK The TH intervals tested were and of the global TH value. For example, the PBA TH was a the interval be to the interval be to and histogram have been We however, chosen this approach the objective be to apply the to soil images from which the be with Because this value be the global TH with the original PBA method also be the was that the interval the global TH value, to the soil sample we to the method as For a a post-segmentation two was also tested on the simulated images. A was also applied by removing the pores than five pixels in The pore was performed using the tool that is available in the image processing of where is the number of pore pixels in the ground-truth image, is the number of pore pixels in the tested image, is the number of solid pixels in the ground-truth image, and is the number of solid pixels in the tested image. We this indicator for its and it the of with other studies 39). we used the relative in the calculated as a performance is the ratio of pixels to the number of where is the number of pore is the number of is the grayscale value in the pore pixels in the original grayscale simulated image, and is the grayscale value in the original grayscale simulated image. is a choice the that a pore space should it a performance than the ground-truth information based indicator The physical of the segmentation methods was based on the pore network which is for example, in soil-water retention Because we were with two-dimensional images, we could focus only on the effects of segmentation on the two-dimensional pore network as pore we on the pore to be the than the in In addition to the the plays an important role in pore have where there might be an of The through one pore was by the by the in an where are in were in We each by its to the for each image, which in a global value. The relative in was calculated for each image. We also calculated the relative of each pore For the physical we then two of the global of the image, and describing the pore shape The and were as To or not the quality of the segmentation methods was by noise a was to for in the and for the levels of noise reduction and the three segmentation methods. A was the simulated images the For the three were segmentation to for a impact of noise reduction on the segmentation was performed in of a effect (p < To which segmentation method and noise reduction combination performed most accurately and to see the method were to for in the and between the segmentation method and noise reduction In of (p < a was with as the In each similar analyses of were conducted. Because each pore its own shape each one of the pores (for 15 images was as a Figure shows the and for the 15 simulated images. With and noise reduction the segmentation accuracy a in indicator value an increase in segmentation with the results by and the and for and were noise reduction the value for was than that by but for These we with The and performed well for images but for with relative about for and for The PBA method performed with and In this where the to was noise reduction did not the PBA This is with the of PBA and with our The segmentation was not however, the main of the selection of one threshold value for an entire image. misclassification region and relative for segmentation methods porosity-based and indicator kriging and for pre-segmentation noise and Table the relative in performance indicator between two noise reduction levels for each segmentation the the were not to the effect of one noise it is to use and performance indicators. the from a to a noise reduction were greater for the global segmentation methods and than for the local method the histogram was This however, with the study of found that a greater between two noise reduction levels than Figure shows the resulting images and segmentation for both noise reduction pixels the pores that the ground-truth the pixels pixels that are to the soil but should have been to and the pixels are those to pores but have With from to are and pores have edges The between and for are less the and the indicator kriging segmentation methods for two of pre-segmentation noise reduction pixels the pores that the ground-truth the pixels pixels that are to soil but should have been to pore and the pixels are the one to pore but have With to the of the provided the indicator and provided the one three and the value of In of provided the indicator value its performance not been into have provided the and With to and provided the value to for images. The provided the This is not with the reported by that IK performed than in the of which was the with the simulated images. The between the two studies was due to the TH interval choice when using The method for calculating the TH interval was for images. We this no did or the Because the noise on our simulated images is noise reduction to be and analyses that the significantly segmentation accuracy with and in of With to there was a between the three but was not to the of the to and to were significantly for but in contrast significantly differentiated to and to for but not for These results the variability in indicator and the of the global and local methods. The method by the soil where grayscale are is by This to physical as by The method the pore but the might not be a grayscale value, pixels were into which the grayscale value and the 46 that the use of is not for segmentation quality but provided results in the of ground-truth observed that the segmentation method based on both and for most images. the use of to the accuracy of the segmentation methods. and analyses and to the we that and were more accurate with a noise reduction and that and were not With to the global results Fig. the relative for for 15 images. noise reduction provided with post-segmentation noise reduction the pores edges and has the pore post-segmentation has the of removing assigned to porosity, a post-segmentation could also the those are relative for segmentation methods porosity-based and indicator kriging for pre-segmentation noise and and for both post-segmentation noise and the segmentation The image was a pre- or post-segmentation noise reduction. The image was a pre-segmentation noise reduction and with a post-segmentation noise reduction. pixels the pores that the ground-truth the pixels pixels that are to the soil but should have been to and the pixels are those to pores but have Table the in between or and PBA is here based on the images' ground-truth information and should perform The two global methods performed in a similar way when a noise reduction was applied relative which is with discussed In this where we one segmentation method to a accurate segmentation method the relative of were similar to This is an the indicator is calculated ground-truth The statistical analyses that the segmentation method and combination of noise reduction a impact on Figure shows the main effect for these It shows that a post-segmentation noise reduction pre-segmentation noise reduction to a significantly relative when with combination of pre-segmentation noise reduction post-segmentation noise reduction In a that the of and was (p < and the between and was (p < also that there was a between the PBA and methods but between and or between and as also in Fig. effect for the relative The the pre-segmentation and and post-segmentation and noise reduction as The the segmentation methods and indicator kriging as at simulated image a pre-segmentation noise reduction and segmentation image a pre-segmentation noise reduction and the porosity-based segmentation and image a pre-segmentation noise reduction and the indicator kriging segmentation Figure shows the between noise reduction levels and segmentation methods with to relative the shape We did not the results pores with and and relative errors increase The and PBA methods the This is with the fact that global TH variation to variation the soil but less the pore region edges in Fig. 13). noise reduction to relative which was with the Some pores were and pore edges were when was applied. There was no however, segmentation methods or the noise reduction which was in contrast with the statistical results of effect for the shape relative The the pre-segmentation and and post-segmentation and noise reduction as The the segmentation methods and indicator kriging as Table shows the TH for the segmentation methods and noise reduction. For the TH interval to from to and For this was also the from to The TH then however, from to There was a in the of the soil (see Fig. noise reduction and PBA an TH. This could be as a noise reduction for the global performed with pre-segmentation noise reduction. In those cases, the TH interval the global TH from and With increasing noise the TH interval The was selected as the intensity and this one to a value with noise reduction (see Fig. The selection was based on the value and the value. the of these we the choice of a TH interval the global TH computed by with pre-segmentation noise to and The represents the global threshold with the porosity-based (PBA) segmentation the represents the threshold with indicator kriging The method to the precision of the local IK method with the of the global PBA method by the initial local TH interval the global TH calculated by To perform a analysis of the TH interval impact on the IK method, the TH interval the global TH selected with PBA from to We found that the indicator for the intervals and that, about of the initial value at the interval. For the we only the segmentation results with a TH interval. the and for were in the as those from PBA or For the performance but for or it the The statistical analyses this by only and as significantly from in of and post-segmentation noise analyses that post-segmentation noise reduction did not produce a significantly with The method provided a however, with which was with The to this the results and significantly from other combination of method and noise reduction. The method the binary images the use of a noise reduction process to and the of the RS image used to construct the simulated images to This is with the made by 15 and that a local method could be used as an to pre-segmentation the choice of the TH interval is of importance when using With the histogram simulated image, the interval the global PBA TH far results than the interval calculated by the method and this made the original of X-ray computed tomography is used in soil and hydrological sciences. To be to apply it to many the is to have an accurate and pore space This with choices of image processing that the information initially available on grayscale images. The and processing are of Within this context, noise reduction and segmentation method were tested on simulated grayscale images to perform statistical analyses on five of indicators. It was shown that pre-segmentation noise reduction through a to a in segmentation accuracy for the global segmentation method by and for the local adaptive-window indicator kriging segmentation method threshold interval was calculated with the method the statistical analyses did not significantly those two methods when a pre-segmentation was applied. The PBA method calculated the global threshold value that should be however, the global segmentation a pre-segmentation noise reduction to be a with global noise reduction was shown to be to segmentation quality by the pore The threshold interval choice with the IK method is of major the calculating method to the of image histogram is Our approach based on the global threshold value calculated with the PBA method performed well by indicator that were similar to those generated by but the of using ground-truth information noise reduction We the support of the for We also for the software and for on the statistical We the its and for the fractal Finally, we to the and the for the quality of this