- Email: [email protected]
Image analysis and fractal geometry to characterize soil desiccation cracks
Geoderma 154 (2009) 153–163
Contents lists available at ScienceDirect
Geoderma j o u r n a l h o m e p a g e : w w w. e l s ev i e r. c o m / l o c a t e / g e o d e r m a
Image analysis and fractal geometry to characterize soil desiccation cracks J.U. Baer a, T.F. Kent b, S.H. Anderson c,⁎ a b c
Yoder Baer & Associates, Villa Ridge, Missouri, United States Department of Mathematics, Marywood University, Scranton, PA, United States Department of Soil, Environmental, and Atmospheric Sciences, 302 ABNR Building, Univ. of Missouri, Columbia, MO 65211-7250, United States
a r t i c l e
i n f o
Article history: Received 10 September 2008 Received in revised form 9 October 2009 Accepted 15 October 2009 Available online 8 November 2009 Keywords: Fractal dimension Fractal lacunarity Claypan soils Dessication cracking
a b s t r a c t The rapid transport of water and solute through desiccation soil cracks can lead to crop water and nutrient stress as well as ground and surface water contamination. The objectives of this research were to determine (1) soil surface crack area, (2) fractal dimension of crack edges (D1), and (3) mass fractal dimension (D2) and lacunarity of soil surface cracks. Photographs were used to record in situ soil cracking in a Mexico silt loam (fine, smectitic, mesic Vertic Epiaqualf) at two locations for three sampling dates. Photographs of the cracks were digitized and an algorithm was developed to determine crack area, and the fractal dimension and lacunarity of the digitized images. Location 1 produced 7% crack area, a mass fractal dimension (D2) of 1.44, and a crack edge fractal dimension (D1) of 1.08. Location 2 produced 25% crack area, a mass fractal dimension (D2) of 1.64, and a crack edge fractal dimension of 1.08 (D1). Lacunarity curves as a function of box size for crack boundary fractal dimension were more variable at location 2. Soil profile moisture content and texture analysis indicated the presence of an argillic horizon at the 20-cm depth at Location 1 while at Location 2 the argillic horizon was exposed to the soil surface. Photographic image analysis techniques and fractal analysis appear to distinguish differences in crack area and patterns which may be useful in characterizing soil cracking. © 2009 Elsevier B.V. All rights reserved.
1. Introduction Desiccation cracking is an important phenomena in soils. Bronswijk and Evers-Vermeer (1990) cite several positive and negative effects of cracking on agricultural soils. The positive effects include improved drainage (Bouma et al., 1979), improved infiltration (Swartz, 1966), and improved soil structure (Wilding and Hallmark, 1984; Bronswijk, 1989). Soil cracks have also been utilized as a factor to alter pore size through intermittent wetting (Jalai-Farahani et al., 1993), and during reclamation of salt-affected high shrink-swell soils as pathways for water and salts (Malik et al., 1991). One negative effect of cracking on agricultural soils is the rapid transport of water through soil cracks, leading to crop water stress (Thomas and Phillips, 1979). Another negative effect is the rapid transport of solute through cracks, leading to crop nutrient stress (Coles and Trudgill, 1985), ground water contamination (Kazemi et al., 2008), and surface water contamination if transport is to a drainage system (Harris et al., 1994; Brown et al., 1995). In the Midwestern U.S., there are an estimated 4 million ha of soils with an argillic (Bt) horizon high in expanding clay occurring abruptly below the soil surface (Blanco-Canqui et al., 2002). Desiccation soil cracking is commonly observed in expanding clay soils. As these soils dry, there is a soil matrix volume decrease caused by the collapsing layers of the clay micelles. Cracking occurs as shrinkage produces tensile forces ⁎ Corresponding author. E-mail address: [email protected] (S.H. Anderson). 0016-7061/$ – see front matter © 2009 Elsevier B.V. All rights reserved. doi:10.1016/j.geoderma.2009.10.008
acting in the horizontal plane. Limited tension develops in the vertical direction because gravity causes the surface to lower as the volume of the soil matrix decreases (Billings, 1954). The development and extent of desiccation cracking are determined by the interaction of numerous factors, including clay content and mineralogy (Aitchinson and Holmes, 1953), soil thickness, surface configuration, rate of drying, total drying time, depth of groundwater table, and direction of surface drainage (Plummer and Gostin, 1981). Because soil structure is such an important property affecting water storage and movement, modeling water retention and flow requires measurement of crack size and pattern. Examples of basic crack measurements are reported (Zein el Abedine and Robinson, 1971; Lima and Grismer, 1992). Zein El Abedine and Robinson (1971) devised a bilateral device to aid in measuring crack width and used a wire probe to measure crack depth along transects in vertisols of several regions of the Sudan. Logsdon et al. (1990) characterized soil crack distribution in situ by exposing selected horizontal planes and tracing cracks on clear polyethylene sheets. They found this method superior to photographic slides because film noise was eliminated and analysis time reduced. Lima and Grismer (1992) directly measured crack island width, crack width, and crack depth along a transect in a field irrigated with three different salinity waters to determine the effects of soil salinity on soil crack morphology. Imaging techniques have also been developed to make basic crack measurements. Guidi et al. (1978) developed a laboratory method based on electro-optical determinations to estimate size distribution of cracks in artificially dried soil samples. They found that the analysis could be
154
J.U. Baer et al. / Geoderma 154 (2009) 153–163
carried out on photographs without loss of precision, making the method applicable for field measurements of undisturbed soil surfaces. Ringrose-Voase and Bullock (1984) used a commercial image analyzing computer to make basic pore space measurements of resin-impregnated, undisturbed soil samples. They also translated the two-dimensional measurements into three-dimensional volumes using the mathematical tools of stereology. In order to study the role soil pores play in rapid infiltration in a long-term, continuous, no-till corn watershed, Edwards et al. (1988) developed a technique to count and size macropores. They photographed cleaned, horizontal surfaces at several depths in the watershed and scanned the images with a commercial image analyzer. Lima and Grismer (1992) also used photographic image analysis techniques to determine soil surface cracked. Descriptions of soil surface crack patterns are reported in the literature (Bronswijk, 1991; Chertkov, 2001). Bronswijk (1991) discusses the use of a polygonal (hexagon) crack pattern in the case of bare or pastured soils in the absence of tillage, and a cubical pattern in the case of row-cropped soils to simulate water movement in cracking soils. Lachenbruch (1962) detailed the mechanics of contraction cracking of a two-dimensional system and classified crack networks in terms of their dominant intersection angle. The intuitive view that contraction cracks should form regular hexagons was reinforced by showing that such a system releases the maximum strain energy per unit crack area. However, the variation of polygonal shapes formed by actual contraction cracks made it clear that energy considerations are only useful in explaining the fact that the polygons tend to be equidimensional. He concluded that since a crack tends to propagate in the direction perpendicular to maximum tension, a crack advancing in an oblique direction will intersect a pre-existing crack perpendicularly after entering the stress relief area caused by the pre-existing crack. Research has also been conducted to evaluate the crack network of swelling clay soils (Chertkov and Ravina, 1998; Chertkov, 2000, 2001). Because of the complexities involved in desiccation soil crack formation, no single classification scheme is able to adequately describe all the crack patterns found in nature. Fractal geometry (Mandelbrot, 1982) provides a mathematical framework to quantify the irregular shapes found in nature, such as soil crack patterns. The introduction of fractal geometry has had an enormous impact on modeling and analysis in the natural sciences and there are numerous examples in the literature, such as work done by Young and Crawford (1992), to characterize soil dessication using fractal dimension. It has been used in models to characterize solute transport (Ross, 1986; Rieu and Sposito, 1991a,b). Baveye et al. (1998) observed that stain pattern and measured fractal dimension were highly correlated when a dye was used to study preferential flow. Some research has been conducted to evaluate the fractal properties of crack systems (Yu et al., 2000; Alves et al., 2001; Wang, 2004). Velde (1999) evaluated surface crack patterns in soils and found that cultivated soils showed the greatest irregularity in fractal dimension for a given porosity. Soils with a high degree of cracking may provide more rapid transport of water and pollutants to the subsoil. With appropriate soil properties and crack descriptions, simulated crack networks may be generated to model water and pollutant transport in cracking soils. Evaluation of soil crack patterns relative to their fractal properties would be useful as a characterization tool. The objectives of this research were to use photographic image analysis techniques to record soil surface cracks in order to determine (1) surface crack areas, (2) fractal dimension of crack edges (D1), and (3) mass fractal dimension (D2) and lacunarity of soil surface cracks. 2. Materials and methods
Environmental Quality research site near Centralia in central Missouri (39° 13′ N, 92° 07′ W). Soils at the site are mapped as Mexico silt loam (fine, smectitic, mesic Vertic Epiaqualf). Location 1 was on a summit (1–2% slopes) landscape position and Location 2 was on a backslope (1–3% slopes, eroded) landscape position. These soils are typical claypan soils characterized by an abrupt occurrence of a clay-rich layer at varying depths from 10 cm to 30 cm below the soil surface. A rainout shelter (Fig. 1) was constructed at each location to isolate a 2.44 × 2.44 m soil surface area from precipitation and surface run on. The perimeter of the area was also lined with plastic to a depth of 1.22 m to minimize subsurface flow into the area. Vegetation inside the rainout shelters was removed and the top 2.5 cm of the soil was homogenized and leveled. A grid was laid over the soil surface to define a 1.88 × 1.16 m plot area, leaving an adequate buffer zone around the perimeter. An access tube was installed to a depth of 1.50 m in the buffer zone of each rainout shelter to monitor soil water content by neutron attenuation. Before the experiment began, the profile under each rainout shelter was saturated by sprinkling water on the soil surface for several days. When saturated, the soil profiles were swollen with no cracks visible. Soil particle-size determination by the pipette method (Gee and Bauder, 1986) was performed on soil samples taken at 0.15-m depth increments from the soil surface to the 0.75 m depth in the buffer zone at each location. 2.2. Photography and image files Photographs were taken of the 1.88 × 1.16 m plot surfaces as the soil profiles dried and cracks formed. A Pentax K-1000 camera with 28–70 mm lens was mounted from a fixed position 1.52 m above and parallel to the soil surface. Images were exposed onto ISO 400 color slide film using early morning ambient lighting. Soil water content readings were taken at each rainout shelter at 0.15-m depth increments to 1.50 m after photographing the soil surfaces. Data were collected in this manner at three times: 19 (14 June), 54 (19 July), and 82 (16 August) days after experiment initiation. The photographic images were developed to 20.3 × 30.5 cm (8 × 12 in) glossy color prints and scanned at 300 × 300 dpi, 256 grey scale using a color scanner. The digitized images were cropped to 2000 × 2000 pixels (113.1 × 113.1 cm) with a resolution of 0.57 mm pixel− 1. The 2000 × 2000 pixel images were split into 16 equal sub-areas giving more manageable (500 × 500 pixels) sub-images with which to work. Individual sub-images were edited by hand so that crack pixels could be rendered black and noncrack pixels rendered white. The subimages were reassembled into full 2000 × 2000 pixel binary (crack– noncrack) images, and the number of crack pixels counted. To isolate the crack boundaries, a deletion process was used to eliminate all crack pixels which were not on a crack boundary, leaving only crack boundary pixels. A crack pixel was considered to be a boundary pixel if at least one of its eight immediate neighbors was a noncrack pixel. If a crack pixel failed to be a boundary pixel, then it was considered to be a nonboundary pixel and eliminated. 2.3. Fractal dimension and lacunarity An algorithm developed by Chen (1987) for image texture description and segmentation, based on the box counting probability method (Voss, 1986), was used to determine the fractal dimension, D, and fractal lacunarity, C(L), of the crack images. To calculate these parameters, an image with M points is overlaid with boxes of side L. The number of boxes needed to cover the entire image can be calculated using the following equation: M Pðm; LÞ m=1 m K
2.1. Field experimental details
NðLÞ = ∑
Two locations on a continuous cool-season grass and legume plot were chosen for study at the USDA-ARS Agricultural Systems for
where K is the maximum number of possible points inside a box, M is the total number of points, and P(m,L) is the probability that m points
ð1Þ
J.U. Baer et al. / Geoderma 154 (2009) 153–163
155
Fig. 1. Schematic diagram of rainout shelter with camera access hole and access tube for water content readings.
fall within a box of side L. P(m,L) is determined by centering a box of size L over each point of the image. The number of image points that fall within the box is then counted. P(m,L) is the number of boxes of side L that contain m points divided by the total number of boxes of side L that were considered. The fractal dimension, a measure of the space filling nature of an image, is calculated as the negative slope of the linear regression of log(N) vs. log(L). Fractal lacunarity, a parameter describing the heterogeneity of a fractal image (greater lacunarity implies greater heterogeneity), is calculated as: CðLÞ =
M2⁎ðLÞ−½M1⁎ðLÞ2 ½M1⁎ðLÞ2
ð2Þ
where M1⁎ and M2⁎ are defined as: K
M1⁎ðLÞ = ∑ mPðm; LÞ m=1 K
M2⁎ðLÞ = ∑ m Pðm; LÞ 2
m=1
ð3Þ
ð4Þ
images of known fractal dimension that had the same image dimension as an average crack island from the plot data (approximately 300 × 300 pixels). These images included the following: a hexagon (D = 1.000; Fig. 2), a curve with known low fractal dimension (LDC1, D = 1.086; Fig. 2), and another curve with known low fractal dimension (LDC2, D = 1.161; Fig. 2). These images were chosen since a preliminary analysis of fractal dimension for both plots gave estimated values of less than 1.15. Using a minimum box size of 5 and a maximum box size of 45, it was feasible to both minimize the sum of the squares of the differences between the theoretical dimension and the calculated dimension, as well as maximize the r2 value of the linear regression for the experimental data while giving a wide range of box sizes (relative to the average island size). Using these parameters, the algorithm was tested against images of known fractal dimension greater than 1.161, namely the Triadic Koch, the Quadratic Koch, and the Sierpinski Gasket (Fig. 2). These image curves were selected due to their theoretical fractal dimensions spanning the expected range of surface soil crack boundary images (1 < D < 2). These test images are available at the following website: http://ac.marywood.edu/ tfkent/www/research/fractal_soil/. 2.5. Crack area mass fractal dimension (D2) and lacunarity
This parameter is similar in some ways to an estimate of the coefficient of variation (deviations between the second moment and the square of the first moment divided by the square of the first moment). 2.4. Crack edge fractal dimension (D1) The first method used to estimate the fractal dimension evaluated crack edge boundaries in the dessication crack images, the crack edge fractal dimension (referred to as D1). To minimize interference from other crack areas, the algorithm was modified to count a crack boundary pixel only if there was a path that was completely contained within the maximum box size from it to the center boundary pixel. Further details on this algorithm are located at the following website: http://ac.marywood. edu/tfkent/www/research/fractal_soil/. In order to choose the minimum and maximum box sizes to use as parameters for the program, the modified code was calibrated against
The second method used to estimate the fractal dimension evaluated crack areas in the dessication crack images, the mass fractal dimension (referred to as D2). This analysis was similar to that used for the mass fractal dimension. The D2 code was developed first in order to count all crack boundaries within the box. Then the D2 code was modified to estimate D1 by counting only those crack boundaries where there is a path which can be traced out by moving between adjacent points that lie totally within the maximum box size. Further details on this algorithm are located at the following website: http:// ac.marywood.edu/tfkent/www/research/fractal_soil/. Preliminary analysis showed that the crack area mass fractal dimensions of each image were distributed between 1.30 and 1.70. A series of test images were generated with known fractal dimension based on a two dimensional analog of the Cantor Dust (Fig. 3). These images were each 2000 × 2000 pixels and represented both complete
156
J.U. Baer et al. / Geoderma 154 (2009) 153–163
Fig. 2. Boundaries of sample curves used to test fractal program (D1): (a) hexagon, (b) low dimension curve one (LDC1), (c) low dimension curve two (LDC2), (d) Triadic Koch curve, (e) Quadratic Koch curve, and (f) Sierpinski gasket.
cracking (D ranging from 1.261 to 1.712, Table 1) and incomplete cracking (D ranging from 1.280 to 1.549, Table 1). These test images are available at the following website: http://ac.marywood.edu/ tfkent/www/research/fractal_soil/. In order to generate each of the individual images, a scaling value r was fixed between 0.333 and 0.499. The scaling between cracks at each succeeding level of the image is r, while the largest cracks have width 2000 × (1 − 2r). For the complete cracking, this gives a theoretical dimension of D = −log(4)/log(r), while for the incomplete cracking, a theoretical dimension of D = −log(3)/log(r) is given. These images were generated as partial fractals that only have cracking three levels deep since the experimental images of the soil show only two to three distinct levels of cracking.
Similar to the analysis for the crack edge fractal dimension, these images were used to find that a box size range of 5 to 401 (about 20% of the width of the individual images) minimized the sum of squares of the differences between the theoretical dimension and the calculated dimension, as well as maximized the coefficient of determination of the linear regression of experimental data while providing a wide range of box sizes. 3. Results and discussion 3.1. Theoretical curves The comparison of the calculated fractal dimensions with the theoretical fractal dimensions obtained from testing the crack edge
J.U. Baer et al. / Geoderma 154 (2009) 153–163
157
Fig. 3. Patterns developed to test mass fractal program (D2): (a) random with D = 1.28, (b) boundaries of random with D = 1.28, (c) square with D = 1.62, and (d) boundaries of square with D = 1.62.
(D1) algorithm on selected images is illustrated in Table 2. Fractal dimensions of these curves ranged from 1.000 to 1.585 which covered the range of expected values for the dessication crack images. Coefficients of determination for the log (number of boxes) vs. log (box size) plots estimating fractal dimension were all higher than 0.999. The calculated fractal dimension was slightly higher than the theoretical dimension for four of the curves and slightly lower for two of the curves; however the differences were small. The percent error between the calculated and theoretical values ranged from −0.34% (Sierpinski Gasket) to 1.72% (Quadratic Koch). These values are well within the expected range for estimating fractal dimensions. Examples of the comparison of the calculated fractal dimension with the theoretical fractal dimension obtained from testing the crack area mass fractal dimension (D2) are given in Table 1. For the full cracking test images, the percent error ranged between −5.04% to
Table 1 Theoretical and calculated mass fractal dimensions (D2) with coefficients of determination. Fractal image
Full cracking Full cracking Full cracking Partial cracking Partial cracking Partial cracking
1.77% with an average of.1.22% and a standard deviation of 1.80%. For the partial cracking test images, the percent error ranged between −2.93% to 9.67% with an average error of 5.04% and a standard deviation of 3.35%. These error ranges are higher than those found on the crack edge test; however, they are still within an expected range for estimating fractal dimension, especially from incomplete fractal images. Additionally, both sets of test images gave coefficients of determination ranging from 0.9886 to 0.9934. 3.2. Dessication crack area Fig. 4 shows the digitized, binary crack area images from both Locations for all sampling times. It is apparent from these images that Location 2 had more crack area for all sampling dates compared to
Table 2 Theoretical and calculated fractal dimensions (D1) for the boundary of hexagon, low dimension curve one (LDC1), low dimension curve two (LDC2), Triadic Koch Curve, Quadratic Koch Curve, and Sierpinski Gasket. Fractal image
Fractal dimension Theoretical
Calculated
r2
Difference
− log(4)/log(0.333) = 1.261 − log(4)/log(0.387) = 1.460 − log(4)/log(0.445) = 1.712 − log(3)/log(0.424) = 1.280 − log(3)/log(0.462) = 1.423 − log(3)/log(0.492) = 1.549
1.283 1.445 1.626 1.404 1.501 1.504
.9926 .9925 .9886 .9934 .9901 .9919
0.022 − 0.015 − 0.086 0.124 0.078 − 0.045
Also included is the difference (calculated − theoretical).
Hexagon LDC1 LDC2 Triadic Koch curve Quadratic Koch curve Sierpinski gasket
Fractal dimension Theoretical
Calculated
r2
Difference
log(2)/log(2) = 1.0000 log(7)/log(6) = 1.0860 log(5)/log(4) = 1.1610 log(4)/log(3) = 1.2619 log(8)/log(4) = 1.5000 log(3)/log(2) = 1.5850
1.0092 1.0924 1.1577 1.2699 1.5258 1.5796
.9999 .9999 .9999 .9999 .9999 .9996
0.0092 0.0064 − 0.0033 0.0080 0.0258 − 0.0054
Also included is the difference (calculated − theoretical).
158
J.U. Baer et al. / Geoderma 154 (2009) 153–163
Fig. 4. Digitized surface soil crack area images for two locations and three dates after experiment initiation: (a) Location 1, Day 19; (b) Location 1, Day 54; (c) Location 1, Day 82; (d) Location 2, Day 19; (e) Location 2, Day 54; and (f) Location 2, Day 82. Black regions are cracks and white regions are non-crack areas.
Location 1. This was attributed to higher clay concentrations occurring closer to the surface for Location 2. In addition, crack area seems to increase over time for both locations which was attributed to additional drying of the soil surface with time. In Fig. 5, the digitized, binary crack boundary images from both locations are shown for all three sampling dates. Similar patterns described relative to the binary crack images (Fig. 4) were found for the two locations and sampling dates. Crack areas as a percent of the total image area are given in Table 3 (data obtained from Fig. 4) from both locations for all three dates. At both locations, an increase in crack area with time occurs with linear correlation coefficients of 0.95 (Location 1) and 0.82 (Location 2). It is noted these correlations indicate trends with time since probabilities
are low due to three data points used for each correlation. Most striking is the large difference in crack area between Locations 1 and 2. Location 1 developed nearly 7% crack area while Location 2 developed over 25% crack area over the 82 day desiccation period. Soil profile water content and texture data from both locations provide some explanation for this large difference in crack area between locations. Fig. 6 shows soil profile water content as a function of time at both locations. The large difference between the 7.5-cm depth and the 22.5-cm depth water content values at Location 1, as compared to Location 2, indicates an abrupt textural change between these depths at Location 1. Clay content analysis (Fig. 7) revealed the presence of an argillic horizon beginning approximately at the 20-cm depth at Location 1, while at Location 2 the argillic horizon is exposed to the
J.U. Baer et al. / Geoderma 154 (2009) 153–163
159
Fig. 5. Digitized surface soil crack boundary images for two locations and three dates after experiment initiation: (a) Location 1, Day 19; (b) Location 1, Day 54; (c) Location 1, Day 82; (d) Location 2, Day 19; (e) Location 2, Day 54; and (f) Location 2, Day 82. Images were created from Fig. 4.
soil surface due to erosion. These data indicate that, when not exposed to the surface, the higher water holding capacity of the argillic horizon reduced soil crack area even when surface soil water contents were low. Conversely, if exposed, the argillic horizon produced significantly more cracking. Quantifying crack area alone is not sufficient for modeling water retention, water flow, and solute transport in cracking soil. Consider the case of Location 1 which produced nearly 900 cm2 crack area (7% of the 113 × 113 cm image area) by Day 82. If this crack area were consolidated into one rectilinearly shaped feature, 30 cm wide × 30 cm long × 10 cm deep, the wall area would be 1200 cm2. If, on the other hand, the same crack area was consolidated into one recti-
linearly shaped feature, 1 cm wide × 900 cm long × 10 cm deep, the wall area would be 18,020 cm2, nearly 15 times greater. These extremes illustrate the importance of quantifying crack pattern as well as area. Many soil processes are influenced by the amount of crack wall area. 3.3. Crack edge fractal dimension (D1) Table 3 contains the crack edge fractal dimensions (D1) for the crack boundary images (entire file) from both locations for the sampling dates. Fig. 8 shows the linear regression of log(N) vs. log(L) used to determine the fractal dimension of the crack boundary
160
J.U. Baer et al. / Geoderma 154 (2009) 153–163
Table 3 Surface soil crack area and fractal dimension of crack edge boundaries (D1) estimated using the entire file at the two locations for three dates after experiment initiation. Fractal dimension (D1) Crack area
Entire file
r2
% of entire area Location 1 Day 19 Day 54 Day 82
3.59 4.54 6.96
1.0701 1.0698 1.0801
0.9983 0.9999 0.9999
Location 2 Day 19 Day 54 Day 82
18.13 26.42 25.11
1.0459 1.0569 1.0827
0.9999 0.9999 0.9999
images; all linear regressions show good fit to the data (r2 = 0.998 or higher). The measured fractal dimension for Location 1 was higher compared to Location 2 for all dates except day 82. The fractal
Fig. 7. Mass fraction of clay-size particles as a function of soil depth at Locations 1 and 2.
dimension increased as a linear function of time for Location 2 (r = 0.96). For Location 1, the dimension was constant between day 19 and day 54 and then increased between day 54 and day 82. The fractal dimension (to two significant figures) for boundary images is 1.08 for both Location 1 and Location 2 at day 82. The crack edge fractal dimension (D1) is an indication of the fractal nature of the edges of the cracks as compared to the mass fractal dimension of cracks. 3.4. Crack area mass fractal dimension (D2) and lacunarity
Fig. 6. Soil water content as a function of days after experiment initiation and soil depth for two locations: (a) Location 1 and (b) Location 2.
Mass fractal dimension (D2) values for the crack boundary images (entire file) from both locations for all sampling dates are shown in Table 4. Linear regressions of log(N) vs. log(L) showed good fit to the data (r2 = 0.99 or higher, Fig. 9). The measured fractal dimension for Location 2 was higher compared to Location 1 for all dates. The fractal dimension increased as a linear function of time for Location 1 (r2 = 0.977), but did not follow this pattern for Location 2. For Location 2, the dimension increased between day 19 and day 54 and then decreased slightly between day 54 and day 82. The fractal dimension (to two significant figures) for crack area is 1.44 at Location 1 and 1.64 at Location 2 on day 82. It is interesting to note the differences between the mass (D2) and edge fractal dimension (D1) values; mass values are higher compared to edge values. As expected, the mass fractal dimension is correlated with crack area (r = 0.978, P < 0.01)) while there is no correlation between edge fractal dimension and crack area (r = −0.26, P = 0.62). The fact that the fractal dimension of crack boundary images from Location 2 are greater than Location 1 seems reasonable from viewing the images, since the surface of Location 2 will fragment more under drying than the surface of Location 1 due to the exposed argillic horizon at Location 2. The crack area mass fractal dimension (D2) is an indication of the fractal nature of the crack mass as compared to the fractal dimension (D1) of the crack edges. Rachman et al. (2005) found mass fractal dimensions for computed tomography-measured macropores for three different management systems. They found fractal dimension values that were either comparable or slightly higher than those of the current study with 1.70 under grass hedges and 1.49 under row crop management. A slightly smaller mass fractal dimension to values in this study of 1.16 was determined for a deposition zone. The value reported by Gantzer
J.U. Baer et al. / Geoderma 154 (2009) 153–163
161
Fig. 8. Crack edge boundary fractal dimension (D1) estimation curves of log(number of boxes) versus log(box size) for two locations and three dates after experiment initiation: (a) Location 1, Day 19; (b) Location 1, Day 54; (c) Location 1, Day 82; (d) Location 2, Day 19; (e) Location 2, Day 54; and (f) Location 2, Day 82.
and Anderson (2002) for a no-till treatment was 1.26, which is slightly less than the values found in this study. Fractal lacunarity relationships estimated for crack areas (same images as for D2) for sampling dates at both locations are given in Fig. 10. Both Locations 1 and 2 had an increase in lacunarity over the test period (for box size < 100). For all dates, Location 1 had smaller variation than Location 2. This seems reasonable from viewing the images since Location 1 produced a more heterogeneous crack pattern compared to Location 2. The lacunarity curves at Location 1 are fairly consistent across box size (40 < box size < 125) for the first two dates, while at Location 2, the curves increase nearly proportional to increases in box size initially (box size < 35 to 40) and then decrease with increasing box size (35 < box size < 125). This may indicate that
Location 1 for the first two dates more closely resembles a true fractal pattern compared to Location 2. Other researchers have shown the benefits of assessing fractal lacunarity for describing the fractal patterns in soil properties (Zeng et al., 1996). These researchers indicated that fractal lacunarity (aggregate and pore uniformity) of X-ray computed tomography images was useful in discriminating soils of similar but slightly different structure. Zeng et al. (1996) recommended that use of both lacunarity with the fractal dimension was necessary to discriminate some soil structures. Fractal lacunarity is an assessment of the space filling pattern of the crack area; high lacunarity indicates more space between crack areas and low lacunarity indicates a more dense pattern of crack areas. 3.5. Potential applications
Table 4 Mass fractal dimension (D2) estimated using the entire file at the two locations for three dates after experiment initiation. Fractal dimension (D2) Entire file
r2
Location 1 Day 19 Day 54 Day 82
1.3135 1.3979 1.4358
0.9930 0.9906 0.9911
Location 2 Day 19 Day 54 Day 82
1.5456 1.6458 1.6370
0.9945 0.9926 0.9916
The two different fractal dimensions determined in this study have useful applications for quantifying water and solute interactions with soil crack features. Some simulation models such as the Soil Water Assessment Tool (SWAT, Neitsch et al., 2005) consider soil cracking. Some other models developed to directly include soil cracks in simulating solute transport include CRACK-NP (Armstrong et al., 2000b) and MACRO (Armstrong et al., 2000a). These models either simulate the amount of soil cracking per area or volume of soil using an input parameter indicating the degree of soil cracking (SWAT) or include an estimate of crack macroporosity along with associated hydraulic parameters (CRACK-NP, MACRO). However, most models do not simulate how cracks are distributed or the roughness of crack edges. A recent model was developed to simulate transport from predicted
162
J.U. Baer et al. / Geoderma 154 (2009) 153–163
Fig. 9. Mass fractal dimension (D2) estimation curves of log(number of boxes) versus log(box size) for two locations and three dates after experiment initiation: (a) Location 1, Day 19; (b) Location 1, Day 54; (c) Location 1, Day 82; (d) Location 2, Day 19; (e) Location 2, Day 54; and (f) Location 2, Day 82.
crack distributions (Vogel et al., 2005a). This model was developed using measured crack distribution data (Vogel et al., 2005b). This latter model may include data from the current study to do future predictions of crack distribution patterns as well as crack edge roughness. The crack edge fractal dimension quantified in the current study estimates a parameter related to the degree of roughness along the crack surface. This dimension can be further evaluated and potentially used in models as a parameter on how much area (roughness along the crack edge) is exposed for interfacing with water and solutes which is important for transport predictions with these types of soils (Kazemi et al., 2008). The crack area mass fractal dimension quantified in the current study is related to the amount of crack area or volume. These crack features are currently used in some models (Neitsch et al., 2005). Further evaluation of the application of these fractal dimension estimates relative to their use in models is needed. 4. Conclusions
Fig. 10. Fractal lacunarity calculated from mass fractal dimension (D2) data as a function of box size for two locations and three dates after experiment initiation: (a) Location 1 for three dates and (b) Location 2 for three dates.
This research was conducted to evaluate soil crack patterns evolving over time for exposed surfaces with claypan soils. Surface soil crack area and pattern were quantified at two locations for three sampling dates. Location 2 produced more crack area compared to Location 1 because of an exposed argillic horizon high in smectitic clay. Crack edge fractal dimension (D1) values were similar (1.08) for the two locations indicating the crack edges had similar variations even though the crack areas were quite different between locations.
J.U. Baer et al. / Geoderma 154 (2009) 153–163
Location 2 also produced a more space-filling, heterogeneous crack pattern compared to Location 1 as determined by the mass fractal dimension (D2) of dessication crack areas (1.44 vs. 1.64, respectively). Areas with a high degree of soil cracking, such as Location 2, may provide more rapid transport of water and solute to the subsoil. With appropriate soil properties and crack descriptions, artificial crack networks could be generated to more successfully model water and solute transport in cracking soils. Image analysis techniques and fractal geometry provide a way to quantify soil cracks and their patterns. Acknowledgements The authors wish to express their appreciation to Dr. Yinghui Zeng for modifying the code used to estimate fractal dimension and lacunarity, and to the following individuals for assisting with both field and laboratory work: Kim McGinty, Peter Baer, and David Keller. References Aitchinson, G.D., Holmes, J.W., 1953. Aspects of swelling in the soil profile. Aust. J. Appl. Sci. 4, 244–259. Alves, L.M., da Silva, R.V., Mokross, B.J., 2001. Influence of crack fractal geometry on elastic-plastic fracture mechanics. Physica. A 295, 144–148. Armstrong, A., Aden, K., Amraoui, N., Diekkruger, B., Jarvis, N., Mouvet, C., Nicholls, P., Wittwer, C., 2000a. Comparison of the performance of pesticide-leaching models on a cracking clay soil: results using the Brimstone Farm dataset. Agric. Water Manag. 44, 85–104. Armstrong, A.C., Matthews, A.M., Portwood, A.M., Leeds-Harrison, P.B., Jarvis, N.J., 2000b. CRACK-NP: a pesticide leaching model for cracking clay soils. Agric. Water Manag. 44, 183–199. Baveye, P., Boast, C.W., Ogawa, S., Parlange, J., Steenhuis, T., 1998. Influence of image resolution and thresholding on the apparent mass fractal characteristics of preferential flow patterns in field soils. Water Resour. Res. 34, 2783–2796. Billings, M.P., 1954. Structural geology, 2nd ed. Prentice-Hall, Englewood Cliffs, NJ. Blanco-Canqui, H., Gantzer, C.J., Anderson, S.H., Alberts, E.E., Ghidey, F., 2002. Saturated hydraulic conductivity and its impact on simulated runoff for claypan soils. Soil Sci. Soc. Am. J. 66, 1596–1602. Bouma, J., Dekker, L.W., Haans, J.C.F.M., 1979. Drainability of some Dutch clay soils: a case study of soil survey interpretation. Geoderma 22, 193–203. Bronswijk, J.J.B., 1989. Prediction of actual cracking and subsidence in clay soils. Soil Sci. 148, 87–93. Bronswijk, J.J.B., 1991. Magnitude, modeling and significance of swelling and shrinkage processes in clay soils. Ph.D. Diss. Wageningen Agricultural University, Wageningen, The Netherlands. Bronswijk, J.J.B., Evers-Vermeer, J.J., 1990. Shrinkage of Dutch clay soil aggregates. Neth. J. Agric. Sci. 38, 175–194. Brown, C.D., Hodgkinson, R.A., Rose, D.A., Syers, J.K., Wilcockson, S.J., 1995. Movement of pesticides to surface waters from a heavy clay soil. Pestic. Sci. 43, 131–140. Chen, S.S., 1987. Fractal geometry in image understanding. Ph.D. Diss. Univ. of MissouriColumbia. Chertkov, V.Y., 2000. Using surface crack spacing to predict crack network geometry in swelling soils. Soil Sci. Soc. Am. J. 64, 1918–1921. Chertkov, V.Y., 2001. Reply to “Comments on ‘Using surface crack spacing to predict crack network geometry in swelling soils’”. Soil Sci. Soc. Am. J. 65, 1574–1575. Chertkov, V.Y., Ravina, I., 1998. Modeling the crack network of swelling clay soils. Soil Sci. Soc. Am. J. 62, 1162–1171. Coles, N., Trudgill, S., 1985. The movement of nitrogen fertilizer from the soil surface to drainage waters by preferential flow in weakly structured soils Slapton South Devon UK. Agric. Ecosyst. Environ. 13, 241–259. Edwards, W.M., Norton, L.D., Redmond, C.E., 1988. Characterizing macropores that affect infiltration into nontilled soil. Soil Sci. Soc. Am. J. 52, 483–487. Gantzer, C.J., Anderson, S.H., 2002. Computed tomographic measurement of macroporosity in chisel-disk and no-tillage seedbeds. Soil Tillage Res. 64, 101–111.
163
Gee, G.W., Bauder, J.W., 1986. Particle-size analysis, In: Klute, A. (Ed.), Methods of soil analysis. Part 1, 2nd ed. : Agronomy, vol. 9, pp. 383–411. Guidi, G., Pagliai, M., Petruzzelli, G., 1978. Quantitative size evaluation of cracks and clods in artificially dried soil samples. Geoderma 19, 105–113. Harris, G.L., Nicholls, P.H., Bailey, S.W., Howse, K.R., Mason, D.J., 1994. Factors influencing the loss of pesticides in drainage from cracking clay soils. J. Hydrol. (Amst.) 159, 235–253. Jalai-Farahani, H.R., Heermann, D.F., Duke, H.R., 1993. Physics of surge irrigation. II. Relationships between soil physical and hydraulic parameters. Trans. ASAE 36, 45–50. Kazemi, H.V., Anderson, S.H., Goyne, K.W., Gantzer, C.J., 2008. Atrazine and alachlor transport in claypan soils as influenced by differential antecedent soil water content. J. Environ. Qual. 37, 1599–1607. Lachenbruch, A.H., 1962. Mechanics of thermal contraction cracks and ice-wedge polygons in permafrost. Spec. Pap. - Geol. Soc. Am. 70. Lima, L.A., Grismer, M.E., 1992. Soil crack morphology and soil salinity. Soil Sci. 153, 149–153. Logsdon, S.D., Allmaras, R.R., Wu, L., Swan, J.B., Randall, G.W., 1990. Macroporosity and its relation to saturated hydraulic conductivity under different tillage practices. Soil Sci. Soc. Am. J. 54, 1096–1101. Malik, M., Amrheim, C., Letey, J., 1991. Polyacrylamide to improve water flow and salt removal in a high shrink-swell soil. Soil Sci. Soc. Am. J. 55, 1664–1667. Mandelbrot, B.B., 1982. The fractal geometry of nature. Freeman, San Francisco. Neitsch, S.L., Arnold, J.G., Kiniry, J.R., Williams, J.R., 2005. Soil and water assessment tool theoretical documentation, 2005 version. USDA-ARS, Grassland, Soil and Water Research Laboratory, Temple, Texas. 566 pp. Plummer, P.S., Gostin, V.A., 1981. Shrinkage cracks: desiccation or synaeresis. J. Sediment. Petrol. 54, 1147–1156. Rachman, A., Anderson, S.H., Gantzer, C.J., 2005. Computed tomographic measurement of soil macroporosity parameters as affected by stiff-stemmed grass hedges. Soil Sci. Soc. Am. J. 69, 1609–1616. Rieu, M., Sposito, G., 1991a. Fractal fragmentation, soil porosity, and soil water properties; I, theory. Soil Sci. Soc. Am. J. 55, 1231–1238. Rieu, M., Sposito, G., 1991b. Fractal fragmentation, soil porosity, and soil water properties; II, applications. Soil Sci. Soc. Am. J. 55, 1239–1244. Ringrose-Voase, A., Bullock, P., 1984. The measurement of soil structural parameters by image analysis. In: Bouma, J., Raats, P.A.C. (Eds.), International symposium on water and solute movement in heavy clay soils. Proc. ISSS, Wageningen, The Netherlands. 27–30 Aug. 1984. ILRI, Wageningen, The Netherlands, pp. 350–354. Ross, B., 1986. Dispersion in fractal fracture networks. Water Resour. Res. 22, 823–827. Swartz, G.L., 1966. Water entry into a black earth under flooding. Queensl. J. Agric. Anim. Sci. 23, 407–422. Thomas, G.W., Phillips, R.E., 1979. Consequences of water movement in macropores. J. Environ. Qual. 8, 149–152. Velde, B., 1999. Structure of surface cracks in soil and muds. Geoderma 93, 101–124. Vogel, H.J., Hoffmann, H., Leopold, A., Roth, K., 2005a. Studies of crack dynamics in clay soil: II. A physically based model for crack formation. Geoderma 125, 213–223. Vogel, H.J., Hoffmann, H., Roth, K., 2005b. Studies of crack dynamics in clay soil: I. Experimental methods, results and morphological quantification. Geoderma 125, 203–211. Voss, R., 1986. Random fractals: characterization and measurement. In: Pynn, R., Skeltorp, A. (Eds.), Phenomena in disordered systems. Plenum, New York. Wang, S.G., 2004. The effect of quenched disorder on the fractal dimension of crack for brittle fracture in two dimensions. Phys. A: Stat. Mech. Appl. 335, 1–8. Wilding, L.P., Hallmark, C.T., 1984. Development of structural microfabric properties in shrinking and swelling clays. In: Bouma, J., Raats, P.A.C. (Eds.), Water and solute movement in heavy clay soils. Proc. Symposium ISSS. International Institute for Land Reclamation and Improvement, Wageningen, Netherlands, pp. 1–22. Young, I.M., Crawford, J.W., 1992. The analysis of fracture profiles using fractal geometry. Aust. J. Soil Res. 30, 291–295. Yu, G., Xie, H., Zhao, J., Yang, L., 2000. Fractal evolution of a crack network in overburden rock stratum. Discrete Dyn. Nat. Soc. 5, 47–52. Zein el Abedine, A., Robinson, G.H., 1971. A study on cracking in some vertisols of the Sudan. Geoderma 5, 229–241. Zeng, Y., Gantzer, C.J., Peyton, R.L., Anderson, S.H., 1996. Fractal dimension and lacunarity of bulk density determined with X-ray computed tomography. Soil Sci. Soc. Am. J. 60, 1718–1724.
Contents lists available at ScienceDirect
Geoderma j o u r n a l h o m e p a g e : w w w. e l s ev i e r. c o m / l o c a t e / g e o d e r m a
Image analysis and fractal geometry to characterize soil desiccation cracks J.U. Baer a, T.F. Kent b, S.H. Anderson c,⁎ a b c
Yoder Baer & Associates, Villa Ridge, Missouri, United States Department of Mathematics, Marywood University, Scranton, PA, United States Department of Soil, Environmental, and Atmospheric Sciences, 302 ABNR Building, Univ. of Missouri, Columbia, MO 65211-7250, United States
a r t i c l e
i n f o
Article history: Received 10 September 2008 Received in revised form 9 October 2009 Accepted 15 October 2009 Available online 8 November 2009 Keywords: Fractal dimension Fractal lacunarity Claypan soils Dessication cracking
a b s t r a c t The rapid transport of water and solute through desiccation soil cracks can lead to crop water and nutrient stress as well as ground and surface water contamination. The objectives of this research were to determine (1) soil surface crack area, (2) fractal dimension of crack edges (D1), and (3) mass fractal dimension (D2) and lacunarity of soil surface cracks. Photographs were used to record in situ soil cracking in a Mexico silt loam (fine, smectitic, mesic Vertic Epiaqualf) at two locations for three sampling dates. Photographs of the cracks were digitized and an algorithm was developed to determine crack area, and the fractal dimension and lacunarity of the digitized images. Location 1 produced 7% crack area, a mass fractal dimension (D2) of 1.44, and a crack edge fractal dimension (D1) of 1.08. Location 2 produced 25% crack area, a mass fractal dimension (D2) of 1.64, and a crack edge fractal dimension of 1.08 (D1). Lacunarity curves as a function of box size for crack boundary fractal dimension were more variable at location 2. Soil profile moisture content and texture analysis indicated the presence of an argillic horizon at the 20-cm depth at Location 1 while at Location 2 the argillic horizon was exposed to the soil surface. Photographic image analysis techniques and fractal analysis appear to distinguish differences in crack area and patterns which may be useful in characterizing soil cracking. © 2009 Elsevier B.V. All rights reserved.
1. Introduction Desiccation cracking is an important phenomena in soils. Bronswijk and Evers-Vermeer (1990) cite several positive and negative effects of cracking on agricultural soils. The positive effects include improved drainage (Bouma et al., 1979), improved infiltration (Swartz, 1966), and improved soil structure (Wilding and Hallmark, 1984; Bronswijk, 1989). Soil cracks have also been utilized as a factor to alter pore size through intermittent wetting (Jalai-Farahani et al., 1993), and during reclamation of salt-affected high shrink-swell soils as pathways for water and salts (Malik et al., 1991). One negative effect of cracking on agricultural soils is the rapid transport of water through soil cracks, leading to crop water stress (Thomas and Phillips, 1979). Another negative effect is the rapid transport of solute through cracks, leading to crop nutrient stress (Coles and Trudgill, 1985), ground water contamination (Kazemi et al., 2008), and surface water contamination if transport is to a drainage system (Harris et al., 1994; Brown et al., 1995). In the Midwestern U.S., there are an estimated 4 million ha of soils with an argillic (Bt) horizon high in expanding clay occurring abruptly below the soil surface (Blanco-Canqui et al., 2002). Desiccation soil cracking is commonly observed in expanding clay soils. As these soils dry, there is a soil matrix volume decrease caused by the collapsing layers of the clay micelles. Cracking occurs as shrinkage produces tensile forces ⁎ Corresponding author. E-mail address: [email protected] (S.H. Anderson). 0016-7061/$ – see front matter © 2009 Elsevier B.V. All rights reserved. doi:10.1016/j.geoderma.2009.10.008
acting in the horizontal plane. Limited tension develops in the vertical direction because gravity causes the surface to lower as the volume of the soil matrix decreases (Billings, 1954). The development and extent of desiccation cracking are determined by the interaction of numerous factors, including clay content and mineralogy (Aitchinson and Holmes, 1953), soil thickness, surface configuration, rate of drying, total drying time, depth of groundwater table, and direction of surface drainage (Plummer and Gostin, 1981). Because soil structure is such an important property affecting water storage and movement, modeling water retention and flow requires measurement of crack size and pattern. Examples of basic crack measurements are reported (Zein el Abedine and Robinson, 1971; Lima and Grismer, 1992). Zein El Abedine and Robinson (1971) devised a bilateral device to aid in measuring crack width and used a wire probe to measure crack depth along transects in vertisols of several regions of the Sudan. Logsdon et al. (1990) characterized soil crack distribution in situ by exposing selected horizontal planes and tracing cracks on clear polyethylene sheets. They found this method superior to photographic slides because film noise was eliminated and analysis time reduced. Lima and Grismer (1992) directly measured crack island width, crack width, and crack depth along a transect in a field irrigated with three different salinity waters to determine the effects of soil salinity on soil crack morphology. Imaging techniques have also been developed to make basic crack measurements. Guidi et al. (1978) developed a laboratory method based on electro-optical determinations to estimate size distribution of cracks in artificially dried soil samples. They found that the analysis could be
154
J.U. Baer et al. / Geoderma 154 (2009) 153–163
carried out on photographs without loss of precision, making the method applicable for field measurements of undisturbed soil surfaces. Ringrose-Voase and Bullock (1984) used a commercial image analyzing computer to make basic pore space measurements of resin-impregnated, undisturbed soil samples. They also translated the two-dimensional measurements into three-dimensional volumes using the mathematical tools of stereology. In order to study the role soil pores play in rapid infiltration in a long-term, continuous, no-till corn watershed, Edwards et al. (1988) developed a technique to count and size macropores. They photographed cleaned, horizontal surfaces at several depths in the watershed and scanned the images with a commercial image analyzer. Lima and Grismer (1992) also used photographic image analysis techniques to determine soil surface cracked. Descriptions of soil surface crack patterns are reported in the literature (Bronswijk, 1991; Chertkov, 2001). Bronswijk (1991) discusses the use of a polygonal (hexagon) crack pattern in the case of bare or pastured soils in the absence of tillage, and a cubical pattern in the case of row-cropped soils to simulate water movement in cracking soils. Lachenbruch (1962) detailed the mechanics of contraction cracking of a two-dimensional system and classified crack networks in terms of their dominant intersection angle. The intuitive view that contraction cracks should form regular hexagons was reinforced by showing that such a system releases the maximum strain energy per unit crack area. However, the variation of polygonal shapes formed by actual contraction cracks made it clear that energy considerations are only useful in explaining the fact that the polygons tend to be equidimensional. He concluded that since a crack tends to propagate in the direction perpendicular to maximum tension, a crack advancing in an oblique direction will intersect a pre-existing crack perpendicularly after entering the stress relief area caused by the pre-existing crack. Research has also been conducted to evaluate the crack network of swelling clay soils (Chertkov and Ravina, 1998; Chertkov, 2000, 2001). Because of the complexities involved in desiccation soil crack formation, no single classification scheme is able to adequately describe all the crack patterns found in nature. Fractal geometry (Mandelbrot, 1982) provides a mathematical framework to quantify the irregular shapes found in nature, such as soil crack patterns. The introduction of fractal geometry has had an enormous impact on modeling and analysis in the natural sciences and there are numerous examples in the literature, such as work done by Young and Crawford (1992), to characterize soil dessication using fractal dimension. It has been used in models to characterize solute transport (Ross, 1986; Rieu and Sposito, 1991a,b). Baveye et al. (1998) observed that stain pattern and measured fractal dimension were highly correlated when a dye was used to study preferential flow. Some research has been conducted to evaluate the fractal properties of crack systems (Yu et al., 2000; Alves et al., 2001; Wang, 2004). Velde (1999) evaluated surface crack patterns in soils and found that cultivated soils showed the greatest irregularity in fractal dimension for a given porosity. Soils with a high degree of cracking may provide more rapid transport of water and pollutants to the subsoil. With appropriate soil properties and crack descriptions, simulated crack networks may be generated to model water and pollutant transport in cracking soils. Evaluation of soil crack patterns relative to their fractal properties would be useful as a characterization tool. The objectives of this research were to use photographic image analysis techniques to record soil surface cracks in order to determine (1) surface crack areas, (2) fractal dimension of crack edges (D1), and (3) mass fractal dimension (D2) and lacunarity of soil surface cracks. 2. Materials and methods
Environmental Quality research site near Centralia in central Missouri (39° 13′ N, 92° 07′ W). Soils at the site are mapped as Mexico silt loam (fine, smectitic, mesic Vertic Epiaqualf). Location 1 was on a summit (1–2% slopes) landscape position and Location 2 was on a backslope (1–3% slopes, eroded) landscape position. These soils are typical claypan soils characterized by an abrupt occurrence of a clay-rich layer at varying depths from 10 cm to 30 cm below the soil surface. A rainout shelter (Fig. 1) was constructed at each location to isolate a 2.44 × 2.44 m soil surface area from precipitation and surface run on. The perimeter of the area was also lined with plastic to a depth of 1.22 m to minimize subsurface flow into the area. Vegetation inside the rainout shelters was removed and the top 2.5 cm of the soil was homogenized and leveled. A grid was laid over the soil surface to define a 1.88 × 1.16 m plot area, leaving an adequate buffer zone around the perimeter. An access tube was installed to a depth of 1.50 m in the buffer zone of each rainout shelter to monitor soil water content by neutron attenuation. Before the experiment began, the profile under each rainout shelter was saturated by sprinkling water on the soil surface for several days. When saturated, the soil profiles were swollen with no cracks visible. Soil particle-size determination by the pipette method (Gee and Bauder, 1986) was performed on soil samples taken at 0.15-m depth increments from the soil surface to the 0.75 m depth in the buffer zone at each location. 2.2. Photography and image files Photographs were taken of the 1.88 × 1.16 m plot surfaces as the soil profiles dried and cracks formed. A Pentax K-1000 camera with 28–70 mm lens was mounted from a fixed position 1.52 m above and parallel to the soil surface. Images were exposed onto ISO 400 color slide film using early morning ambient lighting. Soil water content readings were taken at each rainout shelter at 0.15-m depth increments to 1.50 m after photographing the soil surfaces. Data were collected in this manner at three times: 19 (14 June), 54 (19 July), and 82 (16 August) days after experiment initiation. The photographic images were developed to 20.3 × 30.5 cm (8 × 12 in) glossy color prints and scanned at 300 × 300 dpi, 256 grey scale using a color scanner. The digitized images were cropped to 2000 × 2000 pixels (113.1 × 113.1 cm) with a resolution of 0.57 mm pixel− 1. The 2000 × 2000 pixel images were split into 16 equal sub-areas giving more manageable (500 × 500 pixels) sub-images with which to work. Individual sub-images were edited by hand so that crack pixels could be rendered black and noncrack pixels rendered white. The subimages were reassembled into full 2000 × 2000 pixel binary (crack– noncrack) images, and the number of crack pixels counted. To isolate the crack boundaries, a deletion process was used to eliminate all crack pixels which were not on a crack boundary, leaving only crack boundary pixels. A crack pixel was considered to be a boundary pixel if at least one of its eight immediate neighbors was a noncrack pixel. If a crack pixel failed to be a boundary pixel, then it was considered to be a nonboundary pixel and eliminated. 2.3. Fractal dimension and lacunarity An algorithm developed by Chen (1987) for image texture description and segmentation, based on the box counting probability method (Voss, 1986), was used to determine the fractal dimension, D, and fractal lacunarity, C(L), of the crack images. To calculate these parameters, an image with M points is overlaid with boxes of side L. The number of boxes needed to cover the entire image can be calculated using the following equation: M Pðm; LÞ m=1 m K
2.1. Field experimental details
NðLÞ = ∑
Two locations on a continuous cool-season grass and legume plot were chosen for study at the USDA-ARS Agricultural Systems for
where K is the maximum number of possible points inside a box, M is the total number of points, and P(m,L) is the probability that m points
ð1Þ
J.U. Baer et al. / Geoderma 154 (2009) 153–163
155
Fig. 1. Schematic diagram of rainout shelter with camera access hole and access tube for water content readings.
fall within a box of side L. P(m,L) is determined by centering a box of size L over each point of the image. The number of image points that fall within the box is then counted. P(m,L) is the number of boxes of side L that contain m points divided by the total number of boxes of side L that were considered. The fractal dimension, a measure of the space filling nature of an image, is calculated as the negative slope of the linear regression of log(N) vs. log(L). Fractal lacunarity, a parameter describing the heterogeneity of a fractal image (greater lacunarity implies greater heterogeneity), is calculated as: CðLÞ =
M2⁎ðLÞ−½M1⁎ðLÞ2 ½M1⁎ðLÞ2
ð2Þ
where M1⁎ and M2⁎ are defined as: K
M1⁎ðLÞ = ∑ mPðm; LÞ m=1 K
M2⁎ðLÞ = ∑ m Pðm; LÞ 2
m=1
ð3Þ
ð4Þ
images of known fractal dimension that had the same image dimension as an average crack island from the plot data (approximately 300 × 300 pixels). These images included the following: a hexagon (D = 1.000; Fig. 2), a curve with known low fractal dimension (LDC1, D = 1.086; Fig. 2), and another curve with known low fractal dimension (LDC2, D = 1.161; Fig. 2). These images were chosen since a preliminary analysis of fractal dimension for both plots gave estimated values of less than 1.15. Using a minimum box size of 5 and a maximum box size of 45, it was feasible to both minimize the sum of the squares of the differences between the theoretical dimension and the calculated dimension, as well as maximize the r2 value of the linear regression for the experimental data while giving a wide range of box sizes (relative to the average island size). Using these parameters, the algorithm was tested against images of known fractal dimension greater than 1.161, namely the Triadic Koch, the Quadratic Koch, and the Sierpinski Gasket (Fig. 2). These image curves were selected due to their theoretical fractal dimensions spanning the expected range of surface soil crack boundary images (1 < D < 2). These test images are available at the following website: http://ac.marywood.edu/ tfkent/www/research/fractal_soil/. 2.5. Crack area mass fractal dimension (D2) and lacunarity
This parameter is similar in some ways to an estimate of the coefficient of variation (deviations between the second moment and the square of the first moment divided by the square of the first moment). 2.4. Crack edge fractal dimension (D1) The first method used to estimate the fractal dimension evaluated crack edge boundaries in the dessication crack images, the crack edge fractal dimension (referred to as D1). To minimize interference from other crack areas, the algorithm was modified to count a crack boundary pixel only if there was a path that was completely contained within the maximum box size from it to the center boundary pixel. Further details on this algorithm are located at the following website: http://ac.marywood. edu/tfkent/www/research/fractal_soil/. In order to choose the minimum and maximum box sizes to use as parameters for the program, the modified code was calibrated against
The second method used to estimate the fractal dimension evaluated crack areas in the dessication crack images, the mass fractal dimension (referred to as D2). This analysis was similar to that used for the mass fractal dimension. The D2 code was developed first in order to count all crack boundaries within the box. Then the D2 code was modified to estimate D1 by counting only those crack boundaries where there is a path which can be traced out by moving between adjacent points that lie totally within the maximum box size. Further details on this algorithm are located at the following website: http:// ac.marywood.edu/tfkent/www/research/fractal_soil/. Preliminary analysis showed that the crack area mass fractal dimensions of each image were distributed between 1.30 and 1.70. A series of test images were generated with known fractal dimension based on a two dimensional analog of the Cantor Dust (Fig. 3). These images were each 2000 × 2000 pixels and represented both complete
156
J.U. Baer et al. / Geoderma 154 (2009) 153–163
Fig. 2. Boundaries of sample curves used to test fractal program (D1): (a) hexagon, (b) low dimension curve one (LDC1), (c) low dimension curve two (LDC2), (d) Triadic Koch curve, (e) Quadratic Koch curve, and (f) Sierpinski gasket.
cracking (D ranging from 1.261 to 1.712, Table 1) and incomplete cracking (D ranging from 1.280 to 1.549, Table 1). These test images are available at the following website: http://ac.marywood.edu/ tfkent/www/research/fractal_soil/. In order to generate each of the individual images, a scaling value r was fixed between 0.333 and 0.499. The scaling between cracks at each succeeding level of the image is r, while the largest cracks have width 2000 × (1 − 2r). For the complete cracking, this gives a theoretical dimension of D = −log(4)/log(r), while for the incomplete cracking, a theoretical dimension of D = −log(3)/log(r) is given. These images were generated as partial fractals that only have cracking three levels deep since the experimental images of the soil show only two to three distinct levels of cracking.
Similar to the analysis for the crack edge fractal dimension, these images were used to find that a box size range of 5 to 401 (about 20% of the width of the individual images) minimized the sum of squares of the differences between the theoretical dimension and the calculated dimension, as well as maximized the coefficient of determination of the linear regression of experimental data while providing a wide range of box sizes. 3. Results and discussion 3.1. Theoretical curves The comparison of the calculated fractal dimensions with the theoretical fractal dimensions obtained from testing the crack edge
J.U. Baer et al. / Geoderma 154 (2009) 153–163
157
Fig. 3. Patterns developed to test mass fractal program (D2): (a) random with D = 1.28, (b) boundaries of random with D = 1.28, (c) square with D = 1.62, and (d) boundaries of square with D = 1.62.
(D1) algorithm on selected images is illustrated in Table 2. Fractal dimensions of these curves ranged from 1.000 to 1.585 which covered the range of expected values for the dessication crack images. Coefficients of determination for the log (number of boxes) vs. log (box size) plots estimating fractal dimension were all higher than 0.999. The calculated fractal dimension was slightly higher than the theoretical dimension for four of the curves and slightly lower for two of the curves; however the differences were small. The percent error between the calculated and theoretical values ranged from −0.34% (Sierpinski Gasket) to 1.72% (Quadratic Koch). These values are well within the expected range for estimating fractal dimensions. Examples of the comparison of the calculated fractal dimension with the theoretical fractal dimension obtained from testing the crack area mass fractal dimension (D2) are given in Table 1. For the full cracking test images, the percent error ranged between −5.04% to
Table 1 Theoretical and calculated mass fractal dimensions (D2) with coefficients of determination. Fractal image
Full cracking Full cracking Full cracking Partial cracking Partial cracking Partial cracking
1.77% with an average of.1.22% and a standard deviation of 1.80%. For the partial cracking test images, the percent error ranged between −2.93% to 9.67% with an average error of 5.04% and a standard deviation of 3.35%. These error ranges are higher than those found on the crack edge test; however, they are still within an expected range for estimating fractal dimension, especially from incomplete fractal images. Additionally, both sets of test images gave coefficients of determination ranging from 0.9886 to 0.9934. 3.2. Dessication crack area Fig. 4 shows the digitized, binary crack area images from both Locations for all sampling times. It is apparent from these images that Location 2 had more crack area for all sampling dates compared to
Table 2 Theoretical and calculated fractal dimensions (D1) for the boundary of hexagon, low dimension curve one (LDC1), low dimension curve two (LDC2), Triadic Koch Curve, Quadratic Koch Curve, and Sierpinski Gasket. Fractal image
Fractal dimension Theoretical
Calculated
r2
Difference
− log(4)/log(0.333) = 1.261 − log(4)/log(0.387) = 1.460 − log(4)/log(0.445) = 1.712 − log(3)/log(0.424) = 1.280 − log(3)/log(0.462) = 1.423 − log(3)/log(0.492) = 1.549
1.283 1.445 1.626 1.404 1.501 1.504
.9926 .9925 .9886 .9934 .9901 .9919
0.022 − 0.015 − 0.086 0.124 0.078 − 0.045
Also included is the difference (calculated − theoretical).
Hexagon LDC1 LDC2 Triadic Koch curve Quadratic Koch curve Sierpinski gasket
Fractal dimension Theoretical
Calculated
r2
Difference
log(2)/log(2) = 1.0000 log(7)/log(6) = 1.0860 log(5)/log(4) = 1.1610 log(4)/log(3) = 1.2619 log(8)/log(4) = 1.5000 log(3)/log(2) = 1.5850
1.0092 1.0924 1.1577 1.2699 1.5258 1.5796
.9999 .9999 .9999 .9999 .9999 .9996
0.0092 0.0064 − 0.0033 0.0080 0.0258 − 0.0054
Also included is the difference (calculated − theoretical).
158
J.U. Baer et al. / Geoderma 154 (2009) 153–163
Fig. 4. Digitized surface soil crack area images for two locations and three dates after experiment initiation: (a) Location 1, Day 19; (b) Location 1, Day 54; (c) Location 1, Day 82; (d) Location 2, Day 19; (e) Location 2, Day 54; and (f) Location 2, Day 82. Black regions are cracks and white regions are non-crack areas.
Location 1. This was attributed to higher clay concentrations occurring closer to the surface for Location 2. In addition, crack area seems to increase over time for both locations which was attributed to additional drying of the soil surface with time. In Fig. 5, the digitized, binary crack boundary images from both locations are shown for all three sampling dates. Similar patterns described relative to the binary crack images (Fig. 4) were found for the two locations and sampling dates. Crack areas as a percent of the total image area are given in Table 3 (data obtained from Fig. 4) from both locations for all three dates. At both locations, an increase in crack area with time occurs with linear correlation coefficients of 0.95 (Location 1) and 0.82 (Location 2). It is noted these correlations indicate trends with time since probabilities
are low due to three data points used for each correlation. Most striking is the large difference in crack area between Locations 1 and 2. Location 1 developed nearly 7% crack area while Location 2 developed over 25% crack area over the 82 day desiccation period. Soil profile water content and texture data from both locations provide some explanation for this large difference in crack area between locations. Fig. 6 shows soil profile water content as a function of time at both locations. The large difference between the 7.5-cm depth and the 22.5-cm depth water content values at Location 1, as compared to Location 2, indicates an abrupt textural change between these depths at Location 1. Clay content analysis (Fig. 7) revealed the presence of an argillic horizon beginning approximately at the 20-cm depth at Location 1, while at Location 2 the argillic horizon is exposed to the
J.U. Baer et al. / Geoderma 154 (2009) 153–163
159
Fig. 5. Digitized surface soil crack boundary images for two locations and three dates after experiment initiation: (a) Location 1, Day 19; (b) Location 1, Day 54; (c) Location 1, Day 82; (d) Location 2, Day 19; (e) Location 2, Day 54; and (f) Location 2, Day 82. Images were created from Fig. 4.
soil surface due to erosion. These data indicate that, when not exposed to the surface, the higher water holding capacity of the argillic horizon reduced soil crack area even when surface soil water contents were low. Conversely, if exposed, the argillic horizon produced significantly more cracking. Quantifying crack area alone is not sufficient for modeling water retention, water flow, and solute transport in cracking soil. Consider the case of Location 1 which produced nearly 900 cm2 crack area (7% of the 113 × 113 cm image area) by Day 82. If this crack area were consolidated into one rectilinearly shaped feature, 30 cm wide × 30 cm long × 10 cm deep, the wall area would be 1200 cm2. If, on the other hand, the same crack area was consolidated into one recti-
linearly shaped feature, 1 cm wide × 900 cm long × 10 cm deep, the wall area would be 18,020 cm2, nearly 15 times greater. These extremes illustrate the importance of quantifying crack pattern as well as area. Many soil processes are influenced by the amount of crack wall area. 3.3. Crack edge fractal dimension (D1) Table 3 contains the crack edge fractal dimensions (D1) for the crack boundary images (entire file) from both locations for the sampling dates. Fig. 8 shows the linear regression of log(N) vs. log(L) used to determine the fractal dimension of the crack boundary
160
J.U. Baer et al. / Geoderma 154 (2009) 153–163
Table 3 Surface soil crack area and fractal dimension of crack edge boundaries (D1) estimated using the entire file at the two locations for three dates after experiment initiation. Fractal dimension (D1) Crack area
Entire file
r2
% of entire area Location 1 Day 19 Day 54 Day 82
3.59 4.54 6.96
1.0701 1.0698 1.0801
0.9983 0.9999 0.9999
Location 2 Day 19 Day 54 Day 82
18.13 26.42 25.11
1.0459 1.0569 1.0827
0.9999 0.9999 0.9999
images; all linear regressions show good fit to the data (r2 = 0.998 or higher). The measured fractal dimension for Location 1 was higher compared to Location 2 for all dates except day 82. The fractal
Fig. 7. Mass fraction of clay-size particles as a function of soil depth at Locations 1 and 2.
dimension increased as a linear function of time for Location 2 (r = 0.96). For Location 1, the dimension was constant between day 19 and day 54 and then increased between day 54 and day 82. The fractal dimension (to two significant figures) for boundary images is 1.08 for both Location 1 and Location 2 at day 82. The crack edge fractal dimension (D1) is an indication of the fractal nature of the edges of the cracks as compared to the mass fractal dimension of cracks. 3.4. Crack area mass fractal dimension (D2) and lacunarity
Fig. 6. Soil water content as a function of days after experiment initiation and soil depth for two locations: (a) Location 1 and (b) Location 2.
Mass fractal dimension (D2) values for the crack boundary images (entire file) from both locations for all sampling dates are shown in Table 4. Linear regressions of log(N) vs. log(L) showed good fit to the data (r2 = 0.99 or higher, Fig. 9). The measured fractal dimension for Location 2 was higher compared to Location 1 for all dates. The fractal dimension increased as a linear function of time for Location 1 (r2 = 0.977), but did not follow this pattern for Location 2. For Location 2, the dimension increased between day 19 and day 54 and then decreased slightly between day 54 and day 82. The fractal dimension (to two significant figures) for crack area is 1.44 at Location 1 and 1.64 at Location 2 on day 82. It is interesting to note the differences between the mass (D2) and edge fractal dimension (D1) values; mass values are higher compared to edge values. As expected, the mass fractal dimension is correlated with crack area (r = 0.978, P < 0.01)) while there is no correlation between edge fractal dimension and crack area (r = −0.26, P = 0.62). The fact that the fractal dimension of crack boundary images from Location 2 are greater than Location 1 seems reasonable from viewing the images, since the surface of Location 2 will fragment more under drying than the surface of Location 1 due to the exposed argillic horizon at Location 2. The crack area mass fractal dimension (D2) is an indication of the fractal nature of the crack mass as compared to the fractal dimension (D1) of the crack edges. Rachman et al. (2005) found mass fractal dimensions for computed tomography-measured macropores for three different management systems. They found fractal dimension values that were either comparable or slightly higher than those of the current study with 1.70 under grass hedges and 1.49 under row crop management. A slightly smaller mass fractal dimension to values in this study of 1.16 was determined for a deposition zone. The value reported by Gantzer
J.U. Baer et al. / Geoderma 154 (2009) 153–163
161
Fig. 8. Crack edge boundary fractal dimension (D1) estimation curves of log(number of boxes) versus log(box size) for two locations and three dates after experiment initiation: (a) Location 1, Day 19; (b) Location 1, Day 54; (c) Location 1, Day 82; (d) Location 2, Day 19; (e) Location 2, Day 54; and (f) Location 2, Day 82.
and Anderson (2002) for a no-till treatment was 1.26, which is slightly less than the values found in this study. Fractal lacunarity relationships estimated for crack areas (same images as for D2) for sampling dates at both locations are given in Fig. 10. Both Locations 1 and 2 had an increase in lacunarity over the test period (for box size < 100). For all dates, Location 1 had smaller variation than Location 2. This seems reasonable from viewing the images since Location 1 produced a more heterogeneous crack pattern compared to Location 2. The lacunarity curves at Location 1 are fairly consistent across box size (40 < box size < 125) for the first two dates, while at Location 2, the curves increase nearly proportional to increases in box size initially (box size < 35 to 40) and then decrease with increasing box size (35 < box size < 125). This may indicate that
Location 1 for the first two dates more closely resembles a true fractal pattern compared to Location 2. Other researchers have shown the benefits of assessing fractal lacunarity for describing the fractal patterns in soil properties (Zeng et al., 1996). These researchers indicated that fractal lacunarity (aggregate and pore uniformity) of X-ray computed tomography images was useful in discriminating soils of similar but slightly different structure. Zeng et al. (1996) recommended that use of both lacunarity with the fractal dimension was necessary to discriminate some soil structures. Fractal lacunarity is an assessment of the space filling pattern of the crack area; high lacunarity indicates more space between crack areas and low lacunarity indicates a more dense pattern of crack areas. 3.5. Potential applications
Table 4 Mass fractal dimension (D2) estimated using the entire file at the two locations for three dates after experiment initiation. Fractal dimension (D2) Entire file
r2
Location 1 Day 19 Day 54 Day 82
1.3135 1.3979 1.4358
0.9930 0.9906 0.9911
Location 2 Day 19 Day 54 Day 82
1.5456 1.6458 1.6370
0.9945 0.9926 0.9916
The two different fractal dimensions determined in this study have useful applications for quantifying water and solute interactions with soil crack features. Some simulation models such as the Soil Water Assessment Tool (SWAT, Neitsch et al., 2005) consider soil cracking. Some other models developed to directly include soil cracks in simulating solute transport include CRACK-NP (Armstrong et al., 2000b) and MACRO (Armstrong et al., 2000a). These models either simulate the amount of soil cracking per area or volume of soil using an input parameter indicating the degree of soil cracking (SWAT) or include an estimate of crack macroporosity along with associated hydraulic parameters (CRACK-NP, MACRO). However, most models do not simulate how cracks are distributed or the roughness of crack edges. A recent model was developed to simulate transport from predicted
162
J.U. Baer et al. / Geoderma 154 (2009) 153–163
Fig. 9. Mass fractal dimension (D2) estimation curves of log(number of boxes) versus log(box size) for two locations and three dates after experiment initiation: (a) Location 1, Day 19; (b) Location 1, Day 54; (c) Location 1, Day 82; (d) Location 2, Day 19; (e) Location 2, Day 54; and (f) Location 2, Day 82.
crack distributions (Vogel et al., 2005a). This model was developed using measured crack distribution data (Vogel et al., 2005b). This latter model may include data from the current study to do future predictions of crack distribution patterns as well as crack edge roughness. The crack edge fractal dimension quantified in the current study estimates a parameter related to the degree of roughness along the crack surface. This dimension can be further evaluated and potentially used in models as a parameter on how much area (roughness along the crack edge) is exposed for interfacing with water and solutes which is important for transport predictions with these types of soils (Kazemi et al., 2008). The crack area mass fractal dimension quantified in the current study is related to the amount of crack area or volume. These crack features are currently used in some models (Neitsch et al., 2005). Further evaluation of the application of these fractal dimension estimates relative to their use in models is needed. 4. Conclusions
Fig. 10. Fractal lacunarity calculated from mass fractal dimension (D2) data as a function of box size for two locations and three dates after experiment initiation: (a) Location 1 for three dates and (b) Location 2 for three dates.
This research was conducted to evaluate soil crack patterns evolving over time for exposed surfaces with claypan soils. Surface soil crack area and pattern were quantified at two locations for three sampling dates. Location 2 produced more crack area compared to Location 1 because of an exposed argillic horizon high in smectitic clay. Crack edge fractal dimension (D1) values were similar (1.08) for the two locations indicating the crack edges had similar variations even though the crack areas were quite different between locations.
J.U. Baer et al. / Geoderma 154 (2009) 153–163
Location 2 also produced a more space-filling, heterogeneous crack pattern compared to Location 1 as determined by the mass fractal dimension (D2) of dessication crack areas (1.44 vs. 1.64, respectively). Areas with a high degree of soil cracking, such as Location 2, may provide more rapid transport of water and solute to the subsoil. With appropriate soil properties and crack descriptions, artificial crack networks could be generated to more successfully model water and solute transport in cracking soils. Image analysis techniques and fractal geometry provide a way to quantify soil cracks and their patterns. Acknowledgements The authors wish to express their appreciation to Dr. Yinghui Zeng for modifying the code used to estimate fractal dimension and lacunarity, and to the following individuals for assisting with both field and laboratory work: Kim McGinty, Peter Baer, and David Keller. References Aitchinson, G.D., Holmes, J.W., 1953. Aspects of swelling in the soil profile. Aust. J. Appl. Sci. 4, 244–259. Alves, L.M., da Silva, R.V., Mokross, B.J., 2001. Influence of crack fractal geometry on elastic-plastic fracture mechanics. Physica. A 295, 144–148. Armstrong, A., Aden, K., Amraoui, N., Diekkruger, B., Jarvis, N., Mouvet, C., Nicholls, P., Wittwer, C., 2000a. Comparison of the performance of pesticide-leaching models on a cracking clay soil: results using the Brimstone Farm dataset. Agric. Water Manag. 44, 85–104. Armstrong, A.C., Matthews, A.M., Portwood, A.M., Leeds-Harrison, P.B., Jarvis, N.J., 2000b. CRACK-NP: a pesticide leaching model for cracking clay soils. Agric. Water Manag. 44, 183–199. Baveye, P., Boast, C.W., Ogawa, S., Parlange, J., Steenhuis, T., 1998. Influence of image resolution and thresholding on the apparent mass fractal characteristics of preferential flow patterns in field soils. Water Resour. Res. 34, 2783–2796. Billings, M.P., 1954. Structural geology, 2nd ed. Prentice-Hall, Englewood Cliffs, NJ. Blanco-Canqui, H., Gantzer, C.J., Anderson, S.H., Alberts, E.E., Ghidey, F., 2002. Saturated hydraulic conductivity and its impact on simulated runoff for claypan soils. Soil Sci. Soc. Am. J. 66, 1596–1602. Bouma, J., Dekker, L.W., Haans, J.C.F.M., 1979. Drainability of some Dutch clay soils: a case study of soil survey interpretation. Geoderma 22, 193–203. Bronswijk, J.J.B., 1989. Prediction of actual cracking and subsidence in clay soils. Soil Sci. 148, 87–93. Bronswijk, J.J.B., 1991. Magnitude, modeling and significance of swelling and shrinkage processes in clay soils. Ph.D. Diss. Wageningen Agricultural University, Wageningen, The Netherlands. Bronswijk, J.J.B., Evers-Vermeer, J.J., 1990. Shrinkage of Dutch clay soil aggregates. Neth. J. Agric. Sci. 38, 175–194. Brown, C.D., Hodgkinson, R.A., Rose, D.A., Syers, J.K., Wilcockson, S.J., 1995. Movement of pesticides to surface waters from a heavy clay soil. Pestic. Sci. 43, 131–140. Chen, S.S., 1987. Fractal geometry in image understanding. Ph.D. Diss. Univ. of MissouriColumbia. Chertkov, V.Y., 2000. Using surface crack spacing to predict crack network geometry in swelling soils. Soil Sci. Soc. Am. J. 64, 1918–1921. Chertkov, V.Y., 2001. Reply to “Comments on ‘Using surface crack spacing to predict crack network geometry in swelling soils’”. Soil Sci. Soc. Am. J. 65, 1574–1575. Chertkov, V.Y., Ravina, I., 1998. Modeling the crack network of swelling clay soils. Soil Sci. Soc. Am. J. 62, 1162–1171. Coles, N., Trudgill, S., 1985. The movement of nitrogen fertilizer from the soil surface to drainage waters by preferential flow in weakly structured soils Slapton South Devon UK. Agric. Ecosyst. Environ. 13, 241–259. Edwards, W.M., Norton, L.D., Redmond, C.E., 1988. Characterizing macropores that affect infiltration into nontilled soil. Soil Sci. Soc. Am. J. 52, 483–487. Gantzer, C.J., Anderson, S.H., 2002. Computed tomographic measurement of macroporosity in chisel-disk and no-tillage seedbeds. Soil Tillage Res. 64, 101–111.
163
Gee, G.W., Bauder, J.W., 1986. Particle-size analysis, In: Klute, A. (Ed.), Methods of soil analysis. Part 1, 2nd ed. : Agronomy, vol. 9, pp. 383–411. Guidi, G., Pagliai, M., Petruzzelli, G., 1978. Quantitative size evaluation of cracks and clods in artificially dried soil samples. Geoderma 19, 105–113. Harris, G.L., Nicholls, P.H., Bailey, S.W., Howse, K.R., Mason, D.J., 1994. Factors influencing the loss of pesticides in drainage from cracking clay soils. J. Hydrol. (Amst.) 159, 235–253. Jalai-Farahani, H.R., Heermann, D.F., Duke, H.R., 1993. Physics of surge irrigation. II. Relationships between soil physical and hydraulic parameters. Trans. ASAE 36, 45–50. Kazemi, H.V., Anderson, S.H., Goyne, K.W., Gantzer, C.J., 2008. Atrazine and alachlor transport in claypan soils as influenced by differential antecedent soil water content. J. Environ. Qual. 37, 1599–1607. Lachenbruch, A.H., 1962. Mechanics of thermal contraction cracks and ice-wedge polygons in permafrost. Spec. Pap. - Geol. Soc. Am. 70. Lima, L.A., Grismer, M.E., 1992. Soil crack morphology and soil salinity. Soil Sci. 153, 149–153. Logsdon, S.D., Allmaras, R.R., Wu, L., Swan, J.B., Randall, G.W., 1990. Macroporosity and its relation to saturated hydraulic conductivity under different tillage practices. Soil Sci. Soc. Am. J. 54, 1096–1101. Malik, M., Amrheim, C., Letey, J., 1991. Polyacrylamide to improve water flow and salt removal in a high shrink-swell soil. Soil Sci. Soc. Am. J. 55, 1664–1667. Mandelbrot, B.B., 1982. The fractal geometry of nature. Freeman, San Francisco. Neitsch, S.L., Arnold, J.G., Kiniry, J.R., Williams, J.R., 2005. Soil and water assessment tool theoretical documentation, 2005 version. USDA-ARS, Grassland, Soil and Water Research Laboratory, Temple, Texas. 566 pp. Plummer, P.S., Gostin, V.A., 1981. Shrinkage cracks: desiccation or synaeresis. J. Sediment. Petrol. 54, 1147–1156. Rachman, A., Anderson, S.H., Gantzer, C.J., 2005. Computed tomographic measurement of soil macroporosity parameters as affected by stiff-stemmed grass hedges. Soil Sci. Soc. Am. J. 69, 1609–1616. Rieu, M., Sposito, G., 1991a. Fractal fragmentation, soil porosity, and soil water properties; I, theory. Soil Sci. Soc. Am. J. 55, 1231–1238. Rieu, M., Sposito, G., 1991b. Fractal fragmentation, soil porosity, and soil water properties; II, applications. Soil Sci. Soc. Am. J. 55, 1239–1244. Ringrose-Voase, A., Bullock, P., 1984. The measurement of soil structural parameters by image analysis. In: Bouma, J., Raats, P.A.C. (Eds.), International symposium on water and solute movement in heavy clay soils. Proc. ISSS, Wageningen, The Netherlands. 27–30 Aug. 1984. ILRI, Wageningen, The Netherlands, pp. 350–354. Ross, B., 1986. Dispersion in fractal fracture networks. Water Resour. Res. 22, 823–827. Swartz, G.L., 1966. Water entry into a black earth under flooding. Queensl. J. Agric. Anim. Sci. 23, 407–422. Thomas, G.W., Phillips, R.E., 1979. Consequences of water movement in macropores. J. Environ. Qual. 8, 149–152. Velde, B., 1999. Structure of surface cracks in soil and muds. Geoderma 93, 101–124. Vogel, H.J., Hoffmann, H., Leopold, A., Roth, K., 2005a. Studies of crack dynamics in clay soil: II. A physically based model for crack formation. Geoderma 125, 213–223. Vogel, H.J., Hoffmann, H., Roth, K., 2005b. Studies of crack dynamics in clay soil: I. Experimental methods, results and morphological quantification. Geoderma 125, 203–211. Voss, R., 1986. Random fractals: characterization and measurement. In: Pynn, R., Skeltorp, A. (Eds.), Phenomena in disordered systems. Plenum, New York. Wang, S.G., 2004. The effect of quenched disorder on the fractal dimension of crack for brittle fracture in two dimensions. Phys. A: Stat. Mech. Appl. 335, 1–8. Wilding, L.P., Hallmark, C.T., 1984. Development of structural microfabric properties in shrinking and swelling clays. In: Bouma, J., Raats, P.A.C. (Eds.), Water and solute movement in heavy clay soils. Proc. Symposium ISSS. International Institute for Land Reclamation and Improvement, Wageningen, Netherlands, pp. 1–22. Young, I.M., Crawford, J.W., 1992. The analysis of fracture profiles using fractal geometry. Aust. J. Soil Res. 30, 291–295. Yu, G., Xie, H., Zhao, J., Yang, L., 2000. Fractal evolution of a crack network in overburden rock stratum. Discrete Dyn. Nat. Soc. 5, 47–52. Zein el Abedine, A., Robinson, G.H., 1971. A study on cracking in some vertisols of the Sudan. Geoderma 5, 229–241. Zeng, Y., Gantzer, C.J., Peyton, R.L., Anderson, S.H., 1996. Fractal dimension and lacunarity of bulk density determined with X-ray computed tomography. Soil Sci. Soc. Am. J. 60, 1718–1724.