Morphological Granulometric Analysis of Sediment Images

Sediments are routinely analyzed in terms of the sizing characteristics of the grains of which they are composed. Via sieving methods, the grains are separated and a weight-based size distribution constructed. Various moment parameters are computed from the size distribution and these serve as sediment characteristics. This paper examines the feasibility of a fully electronic granularity analysis using digital image processing. The study uses a random model of three-dimensional grains in conjunction with the morphological method of granulometric size distributions. The random model is constructed to simulate sand, silt, and clay particle distributions. Owing to the impossibility of perfectly sifting small grains so that they do not touch, the model is used in both disjoint and non-disjoint modes, and watershed segmentation is applied in the non-disjoint model. The image-based granulometric size distributions are transformed so that they take into account the necessity to view sediment fractions at different magnifications and in different frames. Gray-scale granulometric moments are then computed using both ordinary and reconstructive granulometries. The resulting moments are then compared to moments found from real grains in seven different sediments using standard weight-based size distributions.


INTRODUCTION
Grain size characteristics and surface area measurements are of great importance to sedimentologists, whose interest lies in soil characterization and trace elemental chemistry (Horowitz and Elrik, 1987).Size distributions provide good quantification for soil studies and carry information concerning the weathering phenomena about the transport, sorting, and sediment source (Kranck and Milligan, 1985;Kranck et al., 1996a, b;Lång and Stevens, 1999).Weathering phenomena include chemical decomposition and mechanical disintergration.
The present paper examines the feasibility of an electronic granularity analysis based on digital image processing.In practice, this would be accomplished by digital image capture and the consequent computation of granulometric size distributions based on the images.It would require fairly sensitive imaging equipment and a stable lighting environment to insure consistency between samples.To emulate these conditions, the present study uses a random model of three-dimensional grains formulated to simulate natural sediment mixtures under the assumption of weightvolume proportionality.Morphological granulometric size distributions are computed digitally from the grain-model images and their moments are then compared to the moments arising from a conventional statistical analysis of the real sediment samples.
The natural sediment grains under comparison have been obtained from the Öre estuary (North Sweden).Four spots along the estuary have been selected that are prediminantly silicious type.The samples are predominantly quartz and plagioclase feldspar groups.Samples 1a and 1b sand fractions contain significant amounts of the mica group (biotite) and iron oxides (hematite or magnetite).Their silt and clay fractions also contain mica and chlorite.Sample 2 has more chlorite.Sample 3 has more montmorillonite and mica.Samples 3b, 4a, and 4b have montmorillonite in mid to small size fractions (silt -clay groups).Fig. 1 shows the structure of various sediment groups.Grain shapes are mostly uniform in silicious sediment groups.
At each of the selected spots (sample numbers 1 to 4), near surface samples (1-5 cm deep) are given a sufix 'a'.Samples of the same core taken at deeper pits (30-35 cm) are given a sufix 'b.'The silicious sediment spots have been selected considering results of (Forsgren et al., 1993).The selected samples fall under the following categories: sample 1(a) -fine and very fine sand, and coarse silt; sample 1(b) -coarse silt and very fine sand; sample 2 -coarse, medium and fine silt; samples 3(a), 3(b), 4(a), 4(b) -predominantly fine and very fine silt.
Conventional sieves have been used to obtain the size distribution for grains sizes > 32 µm (i.e) sand and coarse silt.Grains smaller than this range have been measured in the residual suspension using a coulter counter.Grain size analysis on these selected natural samples have been performed and the results compared to the image-based granulometric sizing results.
After postulating a random model and applying image-based granulometric analysis to the simulated grains, the concordance between the statistical analysis applied to both physical and simulated grains is checked.Comparisons are made on both disjoint and non-disjoint grain models.The latter case is likely more practical because it is difficult to physically separate grains.In the non-disjoint model, grains are segmented using watershed transform prior to granulometric analysis.Granulometric sizing is adapted to the sediment sieving process and both ordinary and reconstructive granulometries are applied, the latter being more akin to sediment sieving.As will be seen, the first and second moments based on image processing are very close to those derived from weight measurements of real sediment grains.The order of magnitude between samples is maintained for the third and fourth moments, there being a shift between the two methods.
The results obtained indicate the potential for image-based sediment granularity analysis.The next step would be to build a prototype image capture system and then apply granulometric analysis directly to the images.

MORPHOLOGICAL GRANULOMETRIES
Morphological granulometries model sieving processes (Matheron, 1975).The essential idea is to operate on an image in such a way that fine structure is progressively eliminated.The area of the remaining image is continuously diminished, and this decreasing area is considered as a size distribution.For a granular image, the image is sieved by successively trying to fit structuring elements into the grains and removing grains into which the probes do not fit.For compact grains, eventually all are sieved from the image and the size distribution reaches the value zero.
We begin by defining a binary granulometry and its size distribution.Consider a fixed compact, convex set B. For any positive real number t, the opening of a set S by structuring element tB is denoted by γ tB (S) and defined as the union of all translates of the structuring element that are subsets of S. As t increases, γ tB (S) diminishes, which means that for t > r, γ tB (S) ⊂ γ rB (S).The parameterized mapping γ tB is called a granulometry.For each set S, a size distribution is defined by letting Ω(t) be the area of γ tB (S).Ω(0) is the area of S. The pattern spectrum Φ of S is defined by normalizing the size distribution so that it is increasing and goes from 0 to 1, namely, Φ(t) = 1 − Ω(t)/Ω(0).The pattern spectrum is a probability distribution.Often the derivative of Φ is used.The moments of Φ, called "granulometric moments," are powerful pattern and texture descriptors (Dougherty and Pelz, 1991;Dougherty et al., 1992;Batman and Dougherty, 1997;Sand and Dougherty, 1998;Theera-Umpon and Gader, 2000).In pattern recognition, S is a random set, the pattern spectrum is a random function, and granulometric moments are random variables.
For a set S composed of a union of disjoint compact grains, a reconstructive granulometry is defined by passing in full any grain not completely eliminated by the opening.Whereas an ordinary granulometry generally diminishes each grain progressively until its elimination, a reconstructive granulometry represents a true sieve: for each grain there exists a value t 0 such that the grain is unchanged for t ≤ t 0 and is eliminated for t > t 0 .The cut-off value t 0 is called the granulometric size of the grain.A grain is passed by a reconstructive granulometry if and only if its granulometric size exceeds the parametric multiple of the primitive B. The pattern spectrum of a reconstructive granulometry is defined in the same manner as for an ordinary granulometry.Reconstructive granulometries are used both for pattern classification and image filtering (Chen and Dougherty, 1999).Pedagogical details and applications can be found (Dougherty and Astola, 1999).
A single-opening gray-scale granulometry involves opening a gray-scale image f (x, y) by a structuring element that is also a function.The most common approach involves the special case of opening a function by a set.The opening, γ B (f), of f by a compact set B at the point (x, y) is the supremum, at the point (x, y), of all translates of B in threedimensional space that fit beneath the graph of f.A granulometry generated by B takes the form γ tB (f), and the size distribution Ω(t) is defined by the volume beneath the graph of γ B (f).We assume that f has bounded support so that Ω(t) = 0 for sufficiently large t.Gray-scale granulometric moments have been used successfully for texture-based classification in medical applications (Chen et al., 1993;Chen and Dougherty, 1994;Baeg et al., 1999).For two reasons, we use a square structuring element: (1) the sieve used for the real grains has a square mesh, and (2) there are fast implementations for square structuring elements.
We consider random grain models for both disjoint and non-disjoint grains.Each grain is modeled as nonnegative function on a compact domain.The domains are disjoint in the disjoint model and are allowed to intersect according to an intersection model in the non-disjoint model.Since weight-based size distributions currently used in sedimentology are based on individual grains, we apply the watershed transformation to segment intersecting grains (Meyer andBeucher, 1990, Beucher andMeyer, 1993).Finding good markers is essential for successful watershed application.Automatic methods devised to find markers are often specific to a given set of images.To avoid over segmentation, in our simulation the modeled grains have been hand marked.

GRANULOMETRIC ANALYSIS OF SEDIMENTS
To adapt granulometric sizing to sediment samples, we assume sediments are classified into three size classes, all being disjoint with the others: sand, silt, and clay (conventional size ranges were used).For image processing, the grains cannot be analyzed without magnification because they are too small for good pixel resolution.Moreover, the different sizes cannot be viewed under a single magnification.Classes must be observed with different levels of magnification.Granulometric sizing is carried out individually on these discrete classes and later combined to form the sample sizing distribution.
Proceeding generally, if we have n nonoverlapping decreasing size classes, C 1 , C 2 ,…, C n viewed at decreasing magnifications M 1 , M 2 ,…, M n , and all grains for all classes are employed for granulometric calculation, then the overall size distribution for the magnified grains can be computed from the class size distributions.Let Ω k (t) be the size distribution for class C k at magnification M k .A change of scale is required to convert the size distribution back to the original pre-magnified scaling.This change means that t is replaced by M k t.Moreover, the total volume at magnification M k is M k 3 times the original volume.This means that, in terms of Ω k (t), the size distribution for C k at the premagnified scale is given by The overall size distribution, Ω(t), at the original scale is the sum of the size distributions: Ω k has been determined at the magnified setting, whereas Ω is at the original scale.
In practice, two adjustments have to be made.First, we do not use all grains we are given because this would require an enormous image grid at the levels of magnification necessary to accurately reveal grain structure.Hence, for each class C k we use only a fraction a k of the volume we are given.This means that the formula for Ω(t) must be adapted by multiplying by 1/a k in the C k term.In addition, we are only supplied with a sampling of each full class.If the grains we are given for class C k compose the fraction p k of the total volume of all grains (not just those given to us for granulometric analysis), then for it to be relative to the total sediment population, the size distribution must be normalized by a factor of p k for class C k .It thus takes the form The corresponding volume is given by Hence, the pattern spectrum is Differentiation yields the pattern-spectrum density ) 0 ( The pattern spectrum Φ k for Ω k satisfies the identity, Given the pattern spectrum, we use standard probabilistic definitions to compute its moments, in particular, the mean, standard deviation, skewness, and kurtosis.Fig. 2 shows the process pictorially.

SIMULATION AND DISCUSSION
In our simulations we consider the magnifications M 3 < M 2 < M 1 for sand, silt, and clay, respectively, and for these the size of the digital image is set at 1024 × 1024 pixels, 2048 × 2048 pixels, and 4096 × 4096 pixels, respectively.In our experiment the magnification factors are set at M 3 = 232, M 2 = 2250, and M 1 = 8000.
To model sand, silt, and clay, we use different parameterized three-dimensional random shaped primitives.Multisided polyhedrons are used to model sand and silt type sediments, as most of these grains closely resemble these shapes.To model mid-size grains -silt fraction, two types of polyhedrons are used.Clay-type grains were modeled as random ellipses, Chemical and physical weathering phenomena's makes them close to circular shape.We refer to Fig. 14 in Sondi et al. (1995) for photomicrograph of clay in suspension.In general, an n-sided (total sides) polyhedron has n − 1 random sides, n − 2 random angles, and n − 2 intensity functions.In this paper we prefer to address such shapes by its random sides rather than the total number of edges.
The first grain fraction, sand, is modeled by a polyhedron with six random sides and an additional fixed closing side makes it seven total sides.Six sides, five random angles and five intensity functions are randomly picked from a uniform distribution.Intensities have been chosen to provide linearly increasing and linearly decreasing slopes.The second fraction, silt is modeled by an equal proportion of five and four-sided polyhedrons.Since this fraction had a wide mixture of grain shapes and sizes, two shapes were used to model this fraction in order to better describe the randomness.The edge lengths, intensity functions and angles are randomly picked from a uniform distribution.A six-sided polyhedron has five randomly selected sides, four random angles and an additional side forms the closing edge.Similarly, a five-sided polyhedron has four randomly chosen sides and three random angles and an additional closing side.Depending on the edge lengths and angles, the polyhedron takes a variety of shapes.In most cases, polyhedrons generated are non-regularly shaped, and therefore referring to them by their random sides is appropriate.The third fraction, clay, is modeled by random ellipses.The major, minor axis, and inclination are randomly picked from a uniform range.Depending on the axis dimensions, the shape will vary from a circle to an ellipse and placed (rotated) at an arbitrary direction.The size ranges in pixels units for all the fractions (at respective magnification) are shown in the Appendix.Fig. 3 shows the parameterized grain model.
Model sizing has been chosen so that the volume ratios of simulated grains (in an imaging sense) are equal to the weight ratios of the natural sediments (in a geological sense).The assumption is that weight information is carried by grain intensity.Image Anal Stereol 2001;20:87-99 Shapes also have intensity parameters.For polyhedrons, the shape can be approximated by a finite number of triangles, with additional intensity dimensions.Each base triangle is referred to as the face of a polyhedron.Intensity levels on its faces linearly change (increase or decrease) from randomly picked start and end levels, whose ranges are shown in the appendix.The intensity ranges and the factor of increase on each face have been chosen such that the intensity profile is peaked towards the center of the grain.Intensity at any point t will be a linear function I k (t) on the k th face of the polyhedron.In this simulation, intensity ranges are preset to have a smooth increase and a decrease on each of its faces.An elliptical shape has maximum intensity level at the center and the intensity linearly slopes down to the minimum level towards the edges.The size range and intensity ranges have been chosen to closely model natural sediments with convexity in mind.In the simulation, grains were not rejected if, otherwise, certain combinations of the randomly chosen lengths and angles produce non-convex shapes (which occur in real sediments).Overall, the random grain model closely models both the randomness and structural shapes of the natural sediments.Fig. 4 illustrates fractional samples of natural grains.The simulation is finalized by performing a morphological opening with a disc-shaped structuring element to smooth the grains.Fig. 5 shows the simulated grain models.The image volumes were matched to the given prior knowledge of the sediment weight ratios of the individual fractions and granulometric sizing was computed later, following the procedure described in the previous section.
While the overall model might appear complex, it is greatly simplified in comparison to real grains.Simpler models have been tried, including binary models, but these did not provide the results comparable to sizing distributions for real sediments.Inclusion of intensity as a grain parameter played a crucial role in modeling the fractions, as weight ratios play a major role in quantifying the grain sizes in real sieving.The final model chosen yields grains that appear visually akin to the real grains, as well as comparable granulometric results.
The various parts of Table 1 compare the moment results for the seven samples listed in the leftmost column.The natural sediment results are from weight-based moment calculations using real sediments.There are four columns for image-based granulometries applied to the simulated model.Two are for reconstructive granulometries, applied to both disjoint grains and segmented non-disjoint grains.Two are for ordinary opening granulometries applied to the same two models.The results are very close for the first two moments using all four approaches.For the third and fourth moments, the imaging-based moments preserve the sizing among the moments, but are systematically lower.Owing to this systematic lowering, further study could perhaps arrive at a linear correction factor; however, this is outside our purpose of showing feasibility using a random model.The parts of Fig. 6

CONCLUSION
The purpose of this study has been to investigate the possibility of an imaging-based system to provide quantitative granulometric sediment analysis.The closeness of the image-based granulometric features to the conventional weight-based statistical moments of the seven real sediment samples shows the potential of both reconstructive and opening granulometries.Discrepancies in the weight-based and image-based methods likely result from using a random grain model, together with using the third dimension (intensity) to capture volume (weight), especially in the more sensitive higher moments.The good point is that, according to the data of this study, order is preserved in the higher order moments and the difference can be well corrected by a linear correction factor.More importantly perhaps, practical implementation of an image-based system would involve high-quality image acquisition, together with the digital tools employed here, and any calibration would be relative to moment relations between image-based granulometries and weight-based size distributions on real grains, not on modeled grains.The success of the watershed segmentation in the nondisjoint model, together with its historical success in many other contexts, indicates that sediments need not be sifted into disjoint grains for an imaging-based approach.

APPENDIX
Following parameter settings was used in the simulation: length (pixel unit at respective magnifications), angels in degrees.

Fig. 2 .
Fig. 2. Pictorial view of the steps involved to obtain granulometric sizing of real sediment grains.