Stereology in Measurement and Modeling: a case study

•Stereology is the part of microstructology that quantifies microstructural geometry

Stereology provides a way of measuring geometric properties from sampled data such as micrographs and tomograms. Stereology is a tool to link material properties to microstructural geometry.

•Classical stereology is based on the assumption of a continuum; digital image processing and mathematical morphology adapt it to modern datasets

Stereology deals with the statistics of volume, surface area, mean curvature and Gaussian curvature in three dimensions. In two dimensions it deals correspondingly with area, line lingth and curvature. These geometric quantities correspond to (after applying suitable constants) the Minkowski functionals for three and two dimensions respectively. Thus stereology has roots in integral geometry. Digital images are made up of pixels; they are not continuous. Mathematical morphology provides a formalism for estimating on digital images the values of the measures associated with the Minkowski functionals.

•This case study illustrates (1) how stereology, mathematical morphology and image processing work together to estimate the length per area of microstructural-size cracks using digital images and (2) how stereology provides a way of estimating the probability of cracks growing into one another as an example of its use in microstructural modeling.



 
presented at the
3rd International
Conference on
Microstructology
held at the
University of Alabama
Birmingham, AL
May 15-20, 2005

NEXT SLIDE

What subjects are covered in classical stereology and its contemporary extensions?

•The book, Stereology and Stochastic Geometry, by John Hilliard and Lawrence Lawson summarizes in simple terms, most of what classical stereology is about.

•Stereology got its name from its ability to infer three dimensional quantities from measurements on two dimensional sections.

This is what most people think of when they hear the term "stereology". Generalized stereology formulas apply to spaces with more than three dimensions. The derivations of these stereology formulas are based on integrating over some Grassmann manifold of possible sections or projections to infer such results as the one that says the volume fraction of some phase in a microstructure is equal to the expected areal fraction. Practical stereologists cannot actually averageover all possible sections. Hence they often assume an ergodic property that says that averaging over a single large slice will give the same result as an ensemble average over many slices. If this is true, the microstructure is said to be "isotropic" in that sense.

•Stereology is the application side of stochastic geometry.

Tesselation properties, estimating statistical distributions of bodies of known shape and estimating measures of properties of pattern not tied to the geometric solids or the convex ring are also examples of stereology.


NEXT SLIDE

The Case Study: microcracks and their coalescence.

•Microcracks are microstructural features.

Microcracks can be viewed as elastically soft planes in common with sliding grain or interphase boundaries and slip bands. They often form along the surface interfaces of elastically hard phases. Observation of elastic-plastic materials confirms the prediction of elasticity theory that spontaneous fracture will occur at certain geometric corners due to their singular stress distributions. Once started, these cracks propagate in fatigue according to the Paris-Erdogan relationship like cracks of macroscopic size.

•The linking of microcracks is an important growth mechanism.

Microcracks are apt to encounter obstacles to their growth such as local stress variations. These obstacles cause arrest. However, a growing microcrack may grow into an arrested crack and join to it. The resulting sudden length increase causes a jump in the stress intensity. This is called crack coalescence. The enlarged crack grows more quickly and is less likely to become arrested.


•Predicting the liability to fracture and the constitutive properties of material having microcracks

Knowing the state of microcracking is important to understanding a material's future. There are a few materials that can disintegrate from microcrack coalescence alone. For many others, coalescence is one of the major events in (macro-) crack initiation. Stereology becomes important in measuring the state of microcracking and in modeling its evolution in time.


NEXT SLIDE

Stereology to Measure and Model Microcracks: Crofton's First Theorem

A Short History:

"...the mean breadth of any convex body is equal to the diameter of a circle whose corcumference equals the length of the boundary."

Crofton (1868)



The theorem is expressed in terms of the intersection of a probe set E with a set A. The probe set could be a line.

As seen by Hadwigger (1957)



The "W"s are the Minkowski functionals where the superscript indicates the dimensionality of the space and the subscript is the index.

As seen by Stoyan et al. (1995)



The set being probed is called X and the probe set is a line e. The l.h.s. is the perimeter divided by pi. The r.h.s. is the integral of the Euler characteristic of the intersection.

As seen by Ohser et al. (2000)



 
Hilliard uses words to express the interpretation of this equation for cases where the set has holes or other nonconformities in it.

As seen by Hilliard and Lawson (2003)

(Added note (2007): In the Babylonian Talmud, M'nachot 29b, Moses is confounded by the complex words of Akiva. He asks where these complicated and confusing things came from and is told that they were given to Moses at Sinai. Moses finds reassuring.)

NEXT SLIDE

Interpreting Crofton's First Theorem (Ohser)

The Euler characteristic equals unity if there is an intersection (with a convex figure) and zero otherwise.

The inner integral simply reduces to the projection of the figure onto the radius line at the specified angle. The outer integral with its prefactor averages these projections to obtain the mean breadth. The result is just that given by Crofton in 1868. These more generalized notations for stereological formulas have become common. More precise than previous ones, they are better suited to advancing the mathematical infrastructure of stereology, but they are often obscure and cumbersome in applications.

•Thus, the notation of Ohser is demystified. Similar constructions can be applied to the other formulations.

•Buffon's equation can be derived from Crofton's First Theorem.

NEXT SLIDE

Mathematical Morphology Implements Stereological Measurement

•Digital images result from hit-or-miss transformations in which a structuring element of dots is ANDed with a scene to produce an image divided into pixels.
•Structuring elements consisting of matrices are then convolved with the image to accomplish image processing and measurement.
•This example shows the processing step called "dilation", Minkowski addition of a ball. The matrix consists of ones. The origin, at (2,2), is translated to each pixel. The products of the matrix elements and pixels are boolean added and placed in the pixel at the origin of the structuring element (matrix).
•At the left, the structuring element is shown slightly displaced. The value of the pixel at its origin is zero, but there is a one at (1,3). Consequently a one is placed at the origin. The ones resulting from dilation are shown in blue.

The morphological approach is not limited to boolean operations. For example, edge finding by convolution with the Laplacian operator is implemented in this way.

NEXT SLIDE

Streamlining Mixed Boolean Operations

The length of a lineal array or boundary can be measured using a morphological form of Crofton's First Theorem.

The morphological transformation given below of the image assists making the measurement.

•Each neighborhood of a black and white image pixel can be stored in itself to reduce the number of pixels accessed in a convolution.
•A matrix convolution becomes a bitwise operation that requires access to each pixel only once.
•To measure a boundary only 2x2 neighborhoods are required because all line segments must pass through only adjacent pixels.
•For example, a 2x2 structuring element with (1,1) as its origin (blue dot) has its origin translated to each pixel. The structuring element contains powers of 2 as shown. These powers of 2 are multiplied by their associated (one or zero) picture elements and the products are added. The result is stored in the picture element corresponding to the origin. This transforms a b/w image to a gray scale image.
NEXT SLIDE

Implementing Boundary Length Measurement

Crofton's First Theorem, in a variation ascribed to Cauchy by Serra, is converted into a morphological operator. Recognize that the formula simply repeats that circumference is pi times diameter. Curly P denotes "projection".
Expressed in terms of the transformed image described on the previous slide,with weights c-sub-l and projections given by the formula below, the formula applies to non-convex figures but has a small bias almost always of the sign implied by the greater than or equals sign.
The bitwise operations use delta functions and weights as indicated. The essential idea is to count a projection if it has zeros and ones in key positions while ignoring others.
For example, a 45 degree projection of a boundary on the bottom of a figure would correspond with the structure element on the left. The dots indicate don't-care conditions. This structure element corresponds in detail to four others on the right. Notice that their indices correspond to their numerical value as found in the previous slide.
NEXT SLIDE

Measuring a Lineal Array's Length Using Buffon's Equation

Buffon's estimator uses the number of interceptions per unit length of a linear probe to find the length per unit area of a lineal array. It is an unbiased estimator for an isotropic random array. For a digitized anisotropic array, bias apart from digitization error is correctable.

•Buffon's estimator is suited to lineal features of variable thickness.
•Four scan lines go through each pixel. The number of interceptions is counted.
• Click for Details
 The chart below shows the accuracy in measuring the length of oriented line segments with a four-orientation scan.
The process of applying Buffon's estimator differs from the standard practice of mathematical morphology in that it involves correlation in a linear space and the use of conditional matrices, shown at left. The first structuring element slides along the scan line until it records a hit. It is imediately transformed into the second structuring element and continues along the scan line. When the second records a hit, its form reverts to the initial one and a complete interception is counted.
NEXT SLIDE

The Microcrack Coalescence Mechanism
increases crack length through a jump process

•Microcracks encounter others at random and join together (coalescence).
•The sum of all microcracks on a surface can be regarded as a lineal array.
•Both theoretically and experimentally the probability of coalescence is proportional to the probability of the array overlapping itself and hence to its length per unit area.
•Elastic interactions favor coalescence. Anisotropy of the crack array discourages coalescence.

NEXT SLIDE

A Typical Quadrat of Microcracks in an Al/SiC(w) Composite

•Microcracks typically align themselves into a brick wall pattern on account of mutual elastic sheilding.
•Microcracks often change direction to move toward nearby microcracks. Notice how cracks 1 and 2 have already merged. Given their brick wall pattern, coalescence often takes place diagonally. This is one reason why "newly initiated" cracks often are observed to have "a shear orientation".
•Microcracking can be approximated as a marked point process with the tensile lengths being the marks. The distribution of marks often controls the probability of macrocrack initiation.

•Unlike most other microstructural features, microcracks are both rare and clustered. Stating their average density on a surface tells little.
•Because coalescence is a binary process, two-point correlation functions should be useful predictors. The variogram and kriging methods of geostatistics should be applicable.
•Dividing a surface into quadrats where in each the point process intensity is uniform, is probably the best way of characterizing microcracking.

NEXT SLIDE

Steps in Measuring Length Per Unit Area

a general recipe

1. Removal of artifacts.

In addition to removing scratches, correcting for the point spread function of the imaging system (MTF correction) is an example of artifact removal. Noise is anything that gets in the way of the picture. Hence noise reduction is part of this step.
 

2. Enhancement of contrast between "features" and background.

The separation of the features of interest from everything else tends to be a gradual process. This step may include opening, closing and filtering as well as transformations of tone and color.
 

3. Isolation of features; removal of background. (Segmentation)

Segmentation is often accomplished by filtering or thresholding but depends on the nature of the feature to be isolated.
 

4. Idealizaton of features when necessary to fit the measurement technique.

In this example the features to be measured are lines. Lines are a geometric idea and do not really exist as features in images. Their best approximation is a sequence of adjacent pixels. Lineal features often have considerable thickness. Some form of thinning may be needed. Boundary extraction may offer another approach.
 

5. Making the measurement.

Measuring areal fraction can be performed by pixel counting. Measuring lineal features requires some way of feeling their length. The boundary of a feature can be measured using mathematical morphology as discussed. Finding the length of a lineal feature can be performed by the intersection counting method discussed, since this method tends to reduce the liklihood of thickness contributing a false length.
NEXT SLIDE

Removal of Artifacts: MTF Correction

replica image of Al/SiC(p) composite

In this enlarged portion of the circled image at left, the cracks have diffuse edges and are indistinct from the background.
Matrix convolution with a gradient-enhancing matrix or "mask" (shown) makes the cracks distinct but enhances the noise.
Adaptive unsharp masking gives a better compromise but still leaves some background clutter that would appear as short cracks unless removed by opening.
NEXT SLIDE

Segmentation of the Image

Thresholding and opening to clean the previous image give the image above. The intention is to have a binary image composed only of black features on a white background. However, perhaps not all the black objects in the image are cracks. There will always be errors of false inclusion and false omission. Experience with human versus machine grading of features in images has shown that machine grading gives more consistant results.
There remains the issue of finite thickness. None of the cracks in this image particularly resemble lines. Some are so rounded that they have no apparent orientation. In order to be measured, cracks must have a defined direction at every pixel that consititutes the length direction. 		One approach is to identify the boundary of a feature. The image above shows the boundary of each feature as white on a black background. This image was obtained by matrix convolution of the image with a mask containing a numerical representation of the Laplacian operator, the one that sums the second derivitives in x and y. The length of this boundary can be measured using the Buffon estimator perviously discussed. The boundary could also be measured (on the image at the left) directly using the Crofton morphological estimator. The results would be essentially the same. Were the cracks very line-like, half their boundary would be a good estimate of their length. But, they are not. Another approach is to skeletonize the image. Above, the skeletons appear in black while the (far left) unprocessed image is superimposed in green. Skeletonization is essentially a curve fitting process in which midpoints are selected at intervals and then connected with lines. Such a fit is not always perfect. As you can see above, some excursions can be missed while some probably false branches may be added. Small dot-like features are ignored. One might wish to argue each fitting, but this would be pointless. At this level of resolution we have no way of knowing that the machine fit is not optimal. We can be certain that it is not far from optimal and at least repeatable and rapidly obtained. The resulting skeleton image lends itself to measurement using the Buffon estimator.
NEXT SLIDE

Results of Length Measurements Using the Buffon estimator

Test Pattern

True Length = 51 pixel widths.
Measurement (uncorrected) = 49 pixel widths.   Error:3.7%
Corrected measurement = 51.6 pixel widths.   Error: 1.3% <1 pixel width.

The test pattern is the one shown previously.

Skeletonized Full Quadrat, Al/SiC(p) composite --as shown

Measured Length (uncorrected) = 10331 pixel widths.
Corrected Length = 10899 pixel widths   Measured Area = 3.22 megapixels.

 

Non-Skeletonized Full Quadrat

Measured Length (uncorrected) 13379 pixel widths.
Half Measured Boundary Length (uncorrected) = 15565 pixel widths.

The effect of finite thickness
can be noted by comparing the three
measurements.
NEXT SLIDE

Modeling Microcrack Coalescence

The probability of microcrack coalescence is closely related to the probability of self overlap of the microcrack lineal array.

In a truly isotropically random array, the marginal probability of a self overlap occurring when the array is extended by an increment delta-el is given by Buffon's equation and depends only on length per unit area. Surprisingly this equation is a good predictor of crack coalescence even though the crack array is anisotropic and stress field interactions are present.
Even with anisotropy and interactions the expectation of the probability of coalescence can be shown to still be proportional to the length per unit area of the crack array. However the proportionality constant is no longer 2/pi.
The effect of anisotropy can be accounted for by carrying out the derivation of Buffon's equation through numerically integrating over the rose of directions. The integrals for doing this (at left) are given in Hilliard and Lawson as Eqn.5.33.
NEXT SLIDE

Modeling the Crack Interaction Component as a Geometric Feature

To account for elastic interactions, filamants are added to the end of each crack to represent the salient edge of Yokobori et.al's zone of no return for coalescence.

These extend the length of the cracks in the perpendicular direction for the purpose of calculating of the probability of coalescence.

With both corrections,
k=0.53 which is close to k=2/pi=0.637

The interaction correction largely offsets that due to anisotropy thus explaining why Buffon's equation was a good predictor.
NEXT SLIDE

Monte Carlo Markov Chain Model of Crack Length

A model using the results just obtained predicts the evolution of the crack length histogram for a given quadrat having a relatively uniform point-process intensity allowing calculation of the rate of formation of significant macro-cracks

The matching of the two histograms of crack length, round dots for the actual material, diamond dots for the model, shows that the model reproduces the material behavior. Component life prediction using replicas for inspections and models of this type is realistic and would be cost-efficient for many high-performance systems.

NEXT SLIDE

Conclusions and Acknowledgements

•Classical stereology is easily adapted to the task of measuring geometrical properties as represented through digital images.

•Stereology is useful for modeling those microstructural properties that are capable of geometric interpretation.


I especially wish to thank: Burton R. (Pat) Patterson and the sponsors of the conference for inviting me and for their support, Edward Chen and Masahiro Meshii for the micrographs used, Rick CulpepperIII of Microsoft for donating VisualC++.NET that works so seamlessly with CImg, David Tschumperle of European Project PrestoSpace for help with setting up CImg, Maryellen Brooks and Marietta Frank of the University of Pittsburgh Library System, and especially, John Cahn for his continued encouragement.

NEXT SLIDE

References

M. Crofton, "On the Theory of Local Probability, Applied to Straight Lines Drawn at Random on a Plane, The Methods Used Being Also Extended to the Proof of Certain New Theorems in the Integral Calculus", Phil. Trans. Roy. Soc. London, 158 (1868) 181-199.
R. Gonzalez and R. Woods, Digital Image Processing (Addison-Wesley, Reading, MA 2002).
H. Hadwiger, Vorlesungen ueber Inhalt Oberflaeche und Isoperimetrie (Springer, Berlin 1957).
J. Hilliard and L. Lawson, Stereology and Stochastic Geometry (Springer, Dordrecht 2003).
L.R. Lawson, "Talk Notes: Stereology in Measurement and Modeling: A Case Study", handout at the Third International Conference on Microstructology, Birmingham, AL, May 15-20, 2005.
L.R. Lawson,"FOURWAY.CPP", C++ code for implementing the Buffon estimator on 2-D images, available from the author, free, at larrylawson@bethelcrofter.com.
L.R. Lawson and E.Y. Chen, "Fatigue Crack Coalescence in Discontinuously Reinforced Metal Matrix Composites: Implications for Reliability Prediction" Journal of Composites Technology and Research (ASTM), 21(3) (1999) 147-152.
J. Ohser and F. Muecklich, Statistical Analysis of Microstructures in Materials Science (John Wiley & Sons, Chichester 2000).
C. Picu and V. Gupta, "Stress Singularities at Triple Junctions with Freely Sliding Grains" Int. J. Solids and Structures, 33 (1996) 1535-1541.
J . Serra, Image Analysis and Mathematical Morphology, (Academic Press, London 1982).
D. Stoyan, W. Kendall and J. Mecke, Stochastic Geometry and its Applications (John Wiley, Chichester 1995).
T. Yokobori, M. Uuozumi and M. Ichikawa, "Interaction Between Non-coplanar Parallel Staggered Elastic Cracks", Reports of the Research Institute for Strength and Fracture of Materials, Tohku University, 7 (1971) 25-47.