Portal de Periódicos da CAPES

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

2017; Soil Science Society of America; Volume: 17; Issue: 1 Linguagem: Inglês

10.2136/vzj2016.06.0049

1539-1663

Sarah De Smet, Erwan Plougonven, Angélique Leonard, Aurore Degré, Éléonore Beckers,

Hydrocarbon exploration and reservoir analysis

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 generating fractal images to be used as the background of the final constructed partial volume effect image and represented by white pixels in Fig. 3c. These pixels could not be further modified during the rest of the process. Step B involved allocating smaller pixels to the solid and pore phases; the black pixels of Fig. 3c were the pixels subjected to further fractal divisions. Step A in the generation of decreasing resolution fractal images for Image no. 1, sorted by iteration from left to right: construction of the partial volume effect background. The white pixels represent the soil matrix. The fractal dimension in the Step A iterations was set as follows: As noted above, Step B involved allocating smaller pixels to the solid and pore phases, and the black pixels in Fig. 3c are the pixels subjected to further fractal divisions. The three iterations of Step B (see Fig. 4) were produced with z = 2, n = 2, and D = 1: Step B in the generation of decreasing resolution fractal images for Image no. 1, sorted by iteration from left to right. The black pixels represent pores. The partial volume effect construction following the method of 46 was applied to the generated fractal images (from Step B), and Fig. 4c shows the first image to be processed. From a size of 1728 by 1728, the image was scaled down into an 864 by 864 image by calculating the average of 2 by 2 squares. The down-scaled image (Fig. 4c) was then overlaid on Fig. 4b, which had also been down-scaled from 1728 by 1728 to 864 by 864, by adding the corresponding pixel color to fully represent the effect of all sizes of pores. This newly created image (not shown) was scaled down to 432 by 432 using the same averaging process, overlaid on Fig. 4a, and the resulting image was then scaled down to 216 by 216 by averaging. The result is shown in Fig. 5. Figure 6 illustrates the impact of the averaging process between Fig. 4c and Fig. 5. Final partial volume effect constructed by the fractal generator for Image no. 1. The left-hand image is an enlargement of Fig. 4c (black rounded square in the upper right corner). The right-hand image is an enlargement of a portion of Fig. 5 (black rounded square in the upper right corner). The next step involved the overlaying of Fig. 5 on the corresponding RS image derived from 3. We then added random normal noise to the pure white and pure black pixels, and variance and means were calculated from our scanner noise by scanning the empty chamber (within [0; 255]: mean = 222 and variance = 15.9). The final step was to add Gaussian noise (mean = 0; variance = 0.01) to the whole image to represent high-frequency noise (Fig. 7). Final simulated grayscale image for Image no. 1. Adhering to an algebraic comparison in its strictest sense, the effect of a pre-segmentation median filter (PRE0, none; PRE1, radius one pixel; PRE2, radius two pixels) was tested on the segmentation quality of the simulated images. Median filters assign the median value of the neighboring pixels to the center pixel. These filters are less sensitive to extreme values and no grayscale value is created near the object boundary, resulting in the object edges being better preserved (44). The use of a median filter before segmentation seems to be a common step in the field of soil X-ray CT image processing. The complexity of segmentation is linked to the noise, artifacts, and partial volume effects in the grayscale images. Other sources of image degradation include ring artifacts, streak artifacts, high-frequency noise, scattered photons, and distortions (1). Therefore, besides enhancing image quality, choosing the right segmentation method is crucial. The global thresholding method described by Otsu (30) was tested because it provides acceptable results (15) and can be used in preference to the IK method where there are poorly distinguishable histogram peaks (46). Despite its non-reliability and the existence of more recent and more efficient methods, it is still a widely used method for soil images, probably because it is rapid and easy to use. It automatically chooses a TH based on the minimization of the intraclass variance between two intensity classes of pixels. In our study, it was performed with MATLAB R2015a (The MathWorks). As we had ground-truth information available, the TH that should be applied could be estimated. Through an iterative loop, the TH that minimized the difference between calculated porosity and ground-truth porosity was selected, and this value served as a benchmark. This procedure was based on the method described by 4. The MATLAB R2015a code was provided by the authors. Hereafter, we refer to the method as the PBA method. The IK method (29) has provided good results in various studies (13, 14; 36; 15; 46). Its variation, the adaptive-window indicator kriging method (13), was chosen because 13 concluded that the adaptive method required fewer computational resources than the fixed one while providing very similar results. The IK concept relies on the selection of a TH interval, T1 to T2. All grayscale values below T1 are set at 0 and all values above T2 are set at 1. The values between T1 and T2 are assigned to a specific color, namely a phase, depending on their grayscale value and their already classified neighboring pixels. The adaptive-window IK method modifies this neighboring area based on locally available information to reduce the computational costs when possible. The method was applied using the authors' software, AWIK. The choice of T1 and T2 was based on edge detection using the gradient masks (GM) method (40), an option available within AWIK software. Hereafter, we refer to the method as IK/GM. The PBA method was shown to be satisfactory, although its performance was poorer than that of IK/GM (4). The weakness in IK is the choice of the T1 to T2 interval. 40 proposed an improved automatic TH interval selection method, although it remained sensitive to noise. We therefore sought to combine the physical robustness of the PBA method with the edge-preserving IK method. The aim was to select a TH interval based on the global PBA threshold and then compute the IK method. The TH intervals tested were ±10, ±20, ±30, ±40, and ±50% of the global TH value. For example, if the PBA TH was 94 (on a 0–255 grayscale), the ±10% IK–PBA interval would be 85 to 103, the ±20% interval would be 75 to 113, and so on. Commonly, histogram percentiles would have been tested. We have, however, chosen this approach because the final objective would be to apply the IK–PBA to real soil images from which the porosity would be estimated with laboratory measurements. Because this porosity value would be uncertain, the global TH obtained with the original PBA method would also be uncertain. Therefore, the priority was that the interval include the supposedly "true" global TH value, corresponding to the "true" soil sample porosity. Hereafter, we will refer to the method as IK–PBA. For a functional comparison, a post-segmentation median filter (POST0, none; POST2, radius two pixels) was also tested on the simulated images. A post-processing cleanup was also applied by removing the pores smaller than five pixels in area. The pore characterization was performed using the Analyze Particles tool that is available in the public domain image processing ImageJ (Version 1.47c, National Institutes of Health, http://rsb.info.nih.gov/ij). where P0 is the number of pore pixels in the ground-truth image, PT is the number of pore pixels in the tested image, S0 is the number of solid pixels in the ground-truth image, and ST is the number of solid pixels in the tested image. We chose this simple indicator for its clear interpretation and because it offered the possibility of comparison with other studies (46; 39). Similarly, we used the relative error in the calculated porosity (RE_P) as a performance indicator. Calculated porosity is the ratio of black pixels (pores) to the total number of pixels. where P is the number of pore pixels, T is the total number of pixels, σp2 is the grayscale value variance in the pore pixels in the original grayscale simulated image, and σ2 is the total grayscale value variance in the original grayscale simulated image. Non-uniformity is a natural choice given the uniformity that a pore space should have, although it gives a poorer performance than the ground-truth information based indicator (48). The physical evaluation of the segmentation methods was based on the pore network modeling concept, which is effective, for example, in computing soil-water retention curves (45). Because we were dealing with two-dimensional images, we could focus only on the effects of segmentation on the two-dimensional pore network characteristics, such as pore geometry. More specifically, we focused on the irregular pore shapes, which, despite their name, tend to be the norm rather than the exception in real soils. In addition to the empty–filled dynamic within pores, the wetting film plays an important role in fluid displacement (6). Irregular pore shapes have corners where there might be an accumulation of wetting fluid. The volumetric flow rate through one pore was obtained by multiplying the conductance (g) by the fluid displacement driving force. As in an electric circuit, where resistances are summed in series, conductance values were summed in parallel. We therefore multiplied each pore's dimensionless conductance () by its squared surface area (A2) to sum all the conductance values (g) for each image, which resulted in a global conductance value. The relative error in conductance (RE_K) was calculated for each image. We also calculated the dimensionless conductance () relative error of each pore (RE_g). For the physical analysis, we then had two types of parameters: RE_K, reflecting the global conductance of the image, and RE_g, describing the pore shape accuracy. The RE_K and RE_g indicators were studied as absolute values. To assess whether or not the quality of the segmentation methods was altered by noise reduction, a three-way ANOVA was conducted to test for significant differences in the ME, NU, RE_P, and RE_K indicators for the various levels of noise reduction and the three segmentation methods. A randomized complete block design was applied, the simulated images being the random blocks. For the significant fixed interaction, three two-way ANOVAs were conducted (one per segmentation method) to test for a significant impact of noise reduction on the segmentation results. Tukey's post-hoc test was performed in cases of a significant effect (p < 0.05). To determine which segmentation method and noise reduction combination performed most accurately and to see if the IK–PBA method brought improvement, four two-way ANOVAs were conducted to test for significant differences in the ME, NU, RE_P, and RE_K values between the segmentation method and noise reduction combinations (10 levels). In cases of significance (p < 0.05), a post-hoc Dunnett test was conducted, with IK–PBA as the control. In each case, similar analyses of RE_g were conducted. Because each pore had its own shape factor, each one of the 229 pores (for all 15 images combined) was considered as a random block. Figure 8 shows the ME, NU, and RE_P averaged for the 15 simulated images. With OTSU and IK/GM, PRE1 noise reduction filtering improved the segmentation accuracy because a decrease in indicator value meant an increase in segmentation accuracy. Compared with the results obtained by 11, 46, and 39, the ME and NU values for PRE1 and PRE2 were satisfactory. Without noise reduction (PRE0), the ME value for OTSU was 60% lower than that obtained by 39 but 34% higher for IK/GM. These high differences probably arose because we considered averaged ME values with high standard errors. The OTSU_PRE0 and IK/GM_PRE0 performed very well for some images but poorly for others, with some porosity relative error values about 315% for OTSU and 500% for IK/GM. The PBA method performed well, with ME indicators below 0.01 and NU indicators below 0.05. In this case, where the exact porosity to reach was known, noise reduction did not improve the PBA method. This is consistent with the working principle of PBA and with our experimental conditions. The segmentation was not perfect, however, highlighting the main disadvantage of the selection of one threshold value for an entire two- (or three-) dimensional image. Averaged misclassification error (ME), region non-uniformity (NU), and porosity relative error (RE_P) for all segmentation methods (OTSU, porosity-based [PBA], and adaptive indicator kriging + gradient masks [IK/GM]) and for all pre-segmentation noise reductions (PRE0, PRE1, and PRE2). Table 4 presents the relative variations in performance indicator values between two noise reduction levels for each segmentation method. First, the variations among the indicators were not comparable. Therefore, to characterize the effect of one noise reduction, it is advisable to use multiple and various performance indicators. Then, the variations from a PRE1 to a PRE2 noise reduction were greater for the global segmentation methods (OTSU and PBA) than for the local IK/GM method albeit the histogram bimodality was sharpened. This is, however, consistent with the study of 14, who found that OTSU gave a greater mean difference between two noise reduction levels than IK. Figure 9 shows the resulting images after OTSU and IK/GM segmentation for both noise reduction levels. Black pixels represent the pores that match the ground-truth information, the blue pixels represent pixels that are allocated to the soil matrix but should have been allocated to pores, and the red pixels are those allocated to pores but shouldn't have been. With OTSU, from PRE1 to PRE2, small features are removed (blue pixels) and bigger pores have growing edges (red pixels). The differences between PRE1 and PRE2 for IK/GM are less striking. Resulting Image no. 10 after the OTSU and the adaptive indicator kriging + gradient masks (IK/GM) segmentation methods for two level of pre-segmentation noise reduction (PRE1, PRE2). Black pixels represent the pores that match the ground-truth information, the blue pixels represent pixels that are allocated to soil matrix but should have been allocated to pore and the red pixels are the one allocated to pore but shouldn't have been. With regard to the occurrence of the lowest ME, PBA_PRE0 provided the lowest indicator 12 times; together, PBA_PRE1 and OTSU_PRE1 provided the lowest one three times and almost always had the same ME value (13 times out of 15). In terms of NU, PBA_PRE0 always provided the lowest indicator value and, if its performance had not been taken into account, OTSU_PRE1 would always have provided the lowest ME and NU values. With regard to RE_P, OTSU_PRE1 and PBA_PRE1 provided the identical closest value to real porosity for 13 images. The IK/GM_PRE1 twice provided the lowest RE_P. This is not consistent with the findings reported by 46, who concluded that IK performed better than OTSU in the case of clear bimodal histograms, which was the case with the simulated images. The difference between the two studies was probably due to the TH interval choice when using IK. The gradient masks method (40) for calculating the TH interval was developed for unimodal images. We discuss this point further below. At no point did OTSU_PRE2 or IK/GM_PRE2 give the best performance. Because the added noise on our simulated images is uncorrelated, PRE2 noise reduction seems to be disproportionate and destroys true information. Statistical analyses confirmed that the PRE1 filter significantly improved segmentation accuracy with OTSU and IK/GM in terms of ME. With regard to NU, there was a significant difference between the three OTSU values (PRE0–PRE1–PRE2), but post-Tukey's test was not able to determine the source of the difference. Similarly, PRE0 to PRE1 and PRE0 to PRE2 were significantly different for IK/GM, but in contrast RE_P significantly differentiated PRE0 to PRE1 and PRE0 to PRE2 for OTSU but not for IK/GM. These contrasting results illustrate the variability in indicator definitions and reflect the working principles of the global and local methods. The OTSU method gives different porosities by identifying porosity within the soil matrix where grayscale values are low (porosity is represented by black pixels). This leads to porosity without physical meaning, as noted by 11. The IK/GM method identifies the right pore region, but the limits might not be accurate. Therefore, despite a high grayscale value, some pixels were taken into account, which increased the grayscale value variance and subsequently the NU. 46 concluded that the use of NU is not enough for characterizing segmentation quality but provided acceptable results in the absence of ground-truth information. They observed that selecting the best segmentation method based on both ME and NU agreed for most images. Again, the use of multiple indicators allowed us to better characterize the accuracy of the segmentation methods. From ME, NU, and RE_P analyses and according to the experimental conditions, we showed that OTSU and IK/GM were more accurate with a PRE1 noise reduction and that OTSU_PRE1 and IK/GM_PRE1 were not statistically different. With regard to the global conductance results (RE_K), Fig. 10 depicts the averaged relative error for RE_K for all 15 images. Post-segmentation noise reduction always provided higher averaged RE_K with high standard errors. Indeed, post-segmentation noise reduction alters the pores edges and has influenced the pore conductance values (Fig. 11). However, post-segmentation has the advantage of removing small features wrongly assigned to porosity, albeit a post-segmentation cleanup could also do the job if those features are small enough. Averaged conductance relative error (RE_K) for all segmentation methods (30 [OTSU], porosity-based [PBA], and adaptive indicator kriging + gradient masks [IK/GM]), for all pre-segmentation noise reductions (PRE0, PRE1, and PRE2), and for both post-segmentation noise reductions (POST0 and POST2). Resulting Image no. 6 after the OTSU (30) segmentation method. The left-hand image was obtained without a pre- or post-segmentation noise reduction. The right-hand image was obtained without a pre-segmentation noise reduction and with a two-level post-segmentation noise reduction. Black pixels represent the pores that match the ground-truth information, the blue pixels represent pixels that are allocated to the soil matrix but should have been allocated to pores, and the red pixels are those allocated to pores but shouldn't have been. Table 5 presents the variations in RE_K between OTSU or IK/GM and PBA. As noted above, PBA is here based on the images' ground-truth information and should therefore perform well. The two global methods performed in a similar way when a PRE1 noise reduction was applied (low relative variations), which is consistent with previously discussed results. In this case, where we compared one segmentation method to a supposedly accurate segmentation method (PBA), the relative variations of RE_K were fairly similar to NU values. This is an interesting point because the NU indicator is calculated without ground-truth information, while RE_K is. The statistical analyses concluded that the segmentation method and combination of noise reduction factors separately had a significant impact on RE_K. Figure 12 shows the main factor effect plot for these factors. It shows that applying a post-segmentation noise reduction without any pre-segmentation noise reduction (combination PRE0–POST2) led to a significantly higher averaged relative error when compared with any combination of pre-segmentation noise reduction (PRE0–PRE1–PRE2) without post-segmentation noise reduction (POST0). In particular, a Tukey post-hoc test concluded that the comparison of PRE0–POST0 (or PRE1–POST0) and PRE0–POST2 was highly significant (p < 0.01), and the comparison between PRE2–POST0 and PRE0–POST2 was significant (p < 0.05). Post-hoc tests also concluded that there was a significant difference between the PBA and IK/GM methods but none between OTSU and IK/GM or between OTSU and PBA, as also illustrated in Fig. 13. Main effect plots for the conductance relative error (RE_K). The upper graph displays the pre-segmentation (PRE0, PRE1, and PRE2) and post-segmentation (POST0 and POST2) noise reduction combinations as variables. The lower graph displays the segmentation methods (porosity-based [PBA]; 30 [OTSU]; and adaptive indicator kriging + gradient masks [IK/GM]) as variables. Image no. 10 at various steps: (a) simulated image; (b) image after a PRE0 pre-segmentation noise reduction and OTSU segmentation (30); (c) image after a PRE0 pre-segmentation noise reduction and the porosity-based segmentation (PBA); and (d) image after a PRE0 pre-segmentation noise reduction and the adaptive indicator kriging + gradient masks segmentation (IK/GM). Figure 14 shows the tendency between noise reduction levels and segmentation methods with regard to RE_g (the relative error across the shape factor). We did not compute the PRE0 results because some pores had merged with OTSU and IK/GM and relative errors would increase without meaning. The OTSU and PBA methods had almost the same mean. This is consistent with the fact that global TH variation would lead to porosity variation within the soil matrix but less so around the pore region edges (also illustrated in Fig. 13). Post-segmentation noise reduction led to higher relative errors, which was consistent with the RE_K results. Some small pores were removed and pore edges were smoothed when POST2 was applied. There was no significant difference, however, among segmentation methods or the noise reduction levels, which was in contrast with the statistical results of RE_K. Main effect plots for the shape factor relative error (RE_g). The upper graph displays the pre-segmentation (PRE0, PRE1, and PRE2) and post-segmentation (POST0 and POST2) noise reduction combinations as variables. The lower graph displays the segmentation methods (porosity-based [PBA]; 30 [OTSU]; and adaptive indicator kriging + gradient masks [IK/GM]) as variables. Table 6 shows the TH median values for the segmentation methods and associated noise reduction. For IK/GM, the TH interval boundaries tended to decrease from PRE0 to PRE1 and PRE2. For OTSU, this was also the case from PRE0 to PRE1. The TH then increased, however, from PRE1 to PRE2. There was indeed a right-hand shift in the lower part of the soil matrix peak (see Fig. 15). At noise reduction PRE1, OTSU and PBA even had an identical TH. This could therefore be seen as a satisfactory noise reduction for the global method. As noted above, IK/GM performed better with pre-segmentation noise reduction. In those cases, the TH interval included the global TH from OTSU and PBA. With increasing noise reduction, the IK/GM TH interval increased. The TH2 was selected as the pore–solid boundary intensity level, and this one moved to a lower value with noise reduction (see Fig. 15). The TH1 selection was based on the TH2 value and the mode value. On the basis of these findings, we investigated the choice of a TH interval around the global TH computed by PBA. Image no. 14 histograms with different pre-segmentation noise reductions. Top to bottom: PRE0, PRE1, and PRE2. The dotted line represents the global threshold obtained with the porosity-based (PBA) segmentation method; the plain bold line represents the upper threshold obtained with adaptive indicator kriging + gradient masks (IK/GM) segmentation. The IK–PBA method tended to combine the precision of the local IK method with the robustness of the global PBA method by selecting the initial local TH interval around the global TH calculated by PBA. To perform a sensitivity analysis of the TH interval impact on the IK method, the TH interval around the global TH selected with PBA ranged from ±10 to ±50%. We found that the ME indicator remained unchanged for the intervals ±10, ±20, and ±30%. After that, ME increased constantly, reaching about 25% of the initial ME value at the ±50% interval. For the following operations, we present only the segmentation results with a ±10% TH interval. First, the averaged ME (0.0023), NU (0.0217), RE_P (0.0096), and RE_K (0.064) values for IK–PBA were all in the same range as those from PBA or OTSU_PRE1. For ME, IK–PBA had the best performance twice, but for NU, RE_P, or RE_K it never had the best performance. The statistical analyses confirmed this trend by showing only OTSU_PRE0 and IK/GM_PRE0 as significantly different from IK–PBA in terms of ME, NU, and RE_P. When including post-segmentation noise reduction, RE_K analyses showed that post-segmentation noise reduction did not produce a significantly different result with IK–PBA. The IK–PBA method provided a significant improvement, however, compared with IK/GM_PRE0_POST2, which was consistent with previous findings. The RE_g again gave contrasting results. According to this indicator, IK–PBA gave the best results and differed significantly from any other combination of method and noise reduction. The IK–PBA method therefore produced the correct binary images without the use of a noise reduction process (as opposed to OTSU_PRE1) and without knowing the real characteristics of the RS image used to construct the simulated images (as opposed to PBA). This is consistent with the recommendation made by 15 and 16 that a local method could be used as an alternative to pre-segmentation processing. As noted above, the choice of the TH interval is of prime importance when using IK. With the two-peak histogram simulated image, the interval around the global PBA TH produced far better results than the interval calculated by the gradient masks method (40), and this made the original idea of IK–PBA attractive. X-ray computed tomography is widely used in soil and hydrological sciences. To be able to apply it to many situations, the prior concern is to have an accurate and correct pore space description. This comes with suitable choices of image processing that will modify the information initially available on grayscale images. The conscious and relevant processing decisions are therefore of great importance. Within this context, various noise reduction and segmentation method combinations were tested on multiple simulated grayscale images to perform statistical analyses on five types of indicators. It was shown that pre-segmentation noise reduction through a median filter led to a significant improvement in segmentation accuracy for the global segmentation method introduced by 30 and for the local adaptive-window indicator kriging (13) segmentation method whose threshold interval was calculated with the gradient masks method (40). Moreover, the statistical analyses did not significantly differentiate those two methods when a pre-segmentation median filter was applied. The PBA method calculated the global threshold value that should be applied; however, the global segmentation wasn't perfect. Therefore, a pre-segmentation noise reduction filter seems to be a necessity with global thresholding. Post-segmentation noise reduction was shown to be detrimental to segmentation quality by altering the pore shapes. The threshold interval choice with the IK method is of major importance. Adapting the calculating method to the type of image histogram is advised. Our approach (IK–PBA) based on the global threshold value calculated with the PBA method performed well by providing indicator values that were similar to those generated by PBA, but IK–PBA had the advantage of using neither ground-truth information nor noise reduction filters. We acknowledge the support of the National Fund for Scientific Research (Brussels, Belgium). We also thank Dr. Alasdair Houston for providing the AWIK software and Professor Yves Brostaux for advising us on the statistical analyses. We thank the AgriGES project, its coordinator Professor Bernard Heinesh, and Donat Regaert for supplying the fractal generator code. Finally, we wish to thank the anonymous reviewers and the editor for improving the quality of this paper.