Image Anal Stereol 2011;30:31-38 Original Research Paper ON THE DILATED FACETS OF A POISSON-VORONOI TESSELLATION Claudia Redenbach University of Kaiserslautern, Mathematics Department, Erwin-Schrödinger-Straße, 67653 Kaiserslautern, Germany e-mail: redenbach@mathematik.uni-kl.de (Accepted February 15, 2011) ABSTRACT In this paper, the parallel set Er of the facets ({d— 1)-faces) of a stationary Poisson-Voronoi tessellation in R^ and M? is investigated. An analytical formula for the spherical contact distribution function of the tessellation allows for the derivation of formulae for the volume density and the specific surface area of Er. The densities of the remaining intrinsic volumes are studied by simulation. The results are used for fitting a dilated Poisson-Voronoi tessellation to the microstructure of a closed-cell foam. Keywords: contact distribution function, intrinsic volumes, parallel set, surface density, stochastic geometry. INTRODUCTION Random tessellations are widely used to model cellular or polycrystalline materials such as foams or sintered ceramics (see, e.g., Telley e/a/., 1996 Coster e/a/., 2005; Ribeiro-Ayeh, 2005; Redenbach, 2009). For instance, the system of facets of a random tessellation can model the wall system in closed-cell foams or the grain boundaries in polycrystalline materials (see Fig. 1 for examples of such materials). Random tessellations are space-filling cell systems by definition, i.e., their cell facets have zero thickness. In contrast, real materials typically consist of two components such that the facets have a certain thickness. Models for such microstructures can be obtained using the parallel set (or dilation) of the facets of the tessellation's cells. Using these models, relations between the microstructure of a material and its macroscopic properties can be investigated (Redenbach, 2009). Fig. 1. Sectional microscopic image of an AhOi, sinter material (left) and visualisation of a tomographic image of a closed polymer foam (right). A well-known tessellation model is the Voronoi tessellation generated by a stationary Poisson point process. A huge number of analytic results for geometric characteristics of this model is available (see Okabe etal, 2000 for an overview). Here, we will use explicit formulae for the spherical contact distribution function to derive formulae for the volume and surface area density of the dilated system of cell facets. These characteristics, which are easy to estimate and to interpret, are sufificient for the estimation of the model parameters, the intensity A of the Poisson process and the dilation radius R. The densities of the remaining intrinsic volumes can then be used to validate a fitted model. So far, however, they can only be studied by simulation. The paper is organised as follows: First, we give an introduction into the concepts and the notation used in this paper. Then we obtain some results for the contact distribution function of the random closed set Er formed by the dilated facet system of a stationary Poisson-Voronoi tessellation. In the following section, the densities of the intrinsic volumes of Er are given as functions of the dilation radius R. Section "Application" is concerned with the application of the obtained results to the modelling of cellular materials. We conclude with a discussion of the results and an outlook on possible future work. DILATED POISSON-VORONOI TESSELLATIONS Throughout this paper we are working in d-dimensional Euclidean space R^ equipped with the Euclidean norm || • ||. For x e and r > 0 let b{x, r) denote the closed (i-dimensional ball of radius r centred in x The unit sphere in R^ is denoted by 5^-1. Write ^ for the Borel sets in R^, X^ for the č/-dimensional Lebesgue measure, and Od-\ for the surface measure on Let ^ and be the system of closed and compact convex sets in R^, respectively. Elements of JT are called convex bodies. Equip ^ with the topology of closed convergence (Schneider and Weil, 2008, p. 18) and denote the corresponding system of Borel sets on ^ by A random closed set S is a measurable mapping from a probability space to It is called stationary if its distribution is invariant under translations of R^. For sets K ^ JT, the intrinsic volumes k = 0,..., č/, can be defined using the Steiner formula k=0 Here, 0 denotes Minkowski addition, i.e., A(BB = {a+ b : ae A,b e B} for A,B e and Kj, is the volume of the A-dimensional unit ball. For č/ = 2 the intrinsic volumes are - up to constant factors - the area ^ = l^, the boundary length I = 214, and the Euler characteristic x = for d = 3 they are the volume l^, the surface area 5 = 21^, the integral of mean curvature M = kV\, and the Euler characteristic = Ko, see e. g. (Schneider, 1993, p. 210). Important characteristics for stationary random closed sets S are the densities of the intrinsic volumes. They are defined as the limits IV,^(S) = lim (1) where H/ G JT is a convex body with K/(H/) > 0 and E denotes expectation. A sufficient condition for the existence of the limit is that S is almost surely locally polyconvex, i.e., fox K ^ the set S H iT can almost surely be written as a finite union of convex bodies, and satisfies ^ < oo (Schneider and Weil, 2008, Theorem 9.2.1). Here N{X) denotes the smallest number m such that Z = iCi U ... U i^j with iCi,..., iiffl G JT . Analogously to the intrinsic volumes, their densities have the following meaning: For d=2, Aa = Vy^2 is the area density. La = 2Vy^\ is the density of the boundary length, and Xa = Vv,o is the density of the Euler characteristic, for rf = 3, Vv= VV,3 is the volume density, Sy = 2 Vy^2 is the surface density (or specific surface area). My = nVy^i is the density of the integral of mean curvature, and Xv = is the density of the Euler characteristic. Denote by Č the interior of the bounded set C C R^. Let r be a set of bounded convex d-dimensional subsets ofR^ with (i) if Ci,C2 G Tand Q ^ C2 then A HQ = 0, (ii) (iii) T is locally finite, i.e., for every bounded B G ^ it holds that #{C G T: Cn 5 0} < Then T is called a tessellation of R^. The cells are d-polytopes due to their convexity and the fact that the tessellation is space-filling (Schneider and Weil, 2008, Lemma 10.1.1). Denote by ^^{C) the set of all k-faces, A = 0,..., č/, of a č/-polytope C. The space of all tessellations of R^ can also be equipped with the topology of closed convergence. Hence, we can define a random tessellation as a random variable X with values in Then the union of all facets, i.e., {d- l)-faces, of the tessellation forms a random closed set u where = The Poisson-Voronoi tessellation is maybe the most well-known random tessellation model (Schneider and Weil, 2008; Stoyan et al, 1995). It is generated by the points of a stationary Poisson process O by assigning to a point G O the cell C{x, O) = {7 G R^ : 117- I < 11/- z| | for all z G O}. In the following, let O be a stationary Poisson process with intensity A, Vo(0) the induced Voronoi tessellation, and S the stationary random closed set formed by the facets of Vo(0). For i? > 0 denote by S7e = Seö(0,i?)= U {F(Bb{Q,R)) (2) the dilated facet system of Vo(0). Finally, denote by the almost surely unique cell of Vo(0) which contains the origin. In the following, we will investigate geometric characteristics of the random closed set S^e- Realisations of this set in both R^ and R^ are shown in Fig. 2. Fig. 2. Realisations of dilated Poisson-Voronoi tessellations in R^ (left) and R^ (right). CONTACT DISTRIBUTION FUNCTIONS Let 5 G Jf^ be a convex body containing the origin and let 0 be a random closed set. Then the contact distribution function Hb of 0 is defined via HB{r) = p{@r\rB^%\Qi@\ r>0. V / Important special cases are the spherical contact distribution function Hg, where B = b{0,1) is the unit ball centred in the origin, and the linear contact distribution function Hj, where B = l{u) is a line segment of unit length in direction u G Explicit formulae for contact and chord length distribution Sanctions of the Poisson-Voronoi tessellation are given in Heinrich (1998) and Muche and Stoyan (1992). Theorem 1. (Heinrich, 1998, Theorem 1) LetBd'^'' be a compact, star-shaped set containing the origin. For r>0 the contact distribution function HB{r) ofE is given by l-HB{r) 0 5£/-i Remark 1. For d =2 and d =3 the spherical contact distribution function ofE has been computed explicitly in Muche and Stoyan (1992). For d =2, it is given by -2nX / pe and for d = 3 we have -X (^(4i^+2p2)(7i-arccos ^)+6ri/p^-i^j dp, (4) / Hlr) = 1-AkI j 0 Remark 2. The volume density Vy^fj of Er equals the spherical contact distribution function ofE, i.e.. For any random closed set 0 and a compact set 5 G the capacity fianctional is defined as T(q{B) = P{Q n 5 / 0). It uniquely determines the distribution of 0 (Schneider and Weil, 2008, Theorem 2.1.3). If 5 G Jf^ is a convex body containing the origin, we have T^{rB) = P{Ef^rB^%) = \-P{rBc Frf(O)) = //^(r) for any r > 0. The capacity fianctional of the dilated tessellation is related to the contact distribution fianction of the original tessellation via T^i^irB) = \-P[rB®ö(0,R) c Frf(O) \ (6) Theorem 2. The contact distribution function Hb^r of Er is given by l-HB,R{r) _ osd-i_ 0 Sd-i Proof The application of the law of total probability to the definitions of Hb^r and Tb,j^{B) yields Now use Eqs. 6, 5, and Theorem 1. (7) □ (3) Vv,ö{ER) = HS{R). (5) Example 1 (Spherical contact distribution fianction). For B= b{Q,\), we have Hence, Hs,R{r) = 1 - ^^^ explicit formulae for Hg in d = 2 and d =3 are given in Remark 1. DENSITIES OF THE INTRINSIC VOLUMES In this section we study the densities of the intrinsic volumes of the random closed set in R^ and R^. These exist since Er fulfills the above-mentioned conditions for the existence of the limit (Eq. 1) which follows from the representation (Eq. 2). We derive integral formulae for the densities of area and boundary length in R^ and for the volume density and the specific surface area in R^. The interpretation of these two densities is straightforward and they are sufficient for the estimation of the parameters A and R. Theorem 3. For d = 2, the area density Aa^r is given by Hs{ R) in Eg. 3. The density of the boundary length La,R of Er is given by La,r = 47rÄ - + 4A J pe -I ((4i?2+2p2)(7i-arccos S)+eRy/p2^R2"j I-\ R[n- arccos - j + ^ p^ - j dp ^ Proof By Stoyan et af (1995, Equation (6.2.4)), La,r can be obtained as La,r = (1 — Aa,r)H'^ . □ The integrals can be evaluated numerically. Plots for the densities of area and boundary length for a tessellation with unit intensity as functions of R are shown in Fig. 3. Theorem 4. For d =3, the volume density Vv,r of Er is given by Hs{R) in Eg. 4. The surface area density Sv,R ofE R is given by Sv,r = / ^ \ 0 \ pe {R+p)'dp Proof By Stoyan et al (1995, Equation (6.2.4)), Sv,r can be obtained as Sv,r = (1 - Vv,r). □ As in the planar case, the integrals can be evaluated numerically. Plots for the densities of volume and surface area for a tessellation with unit intensity as functions of Rare shown in Fig. 4. For the interpretation of the remaining intrinsic volumes we notice that the complement of the random closed set 5r forms a particle process in R'^. It consists of non-overlapping convex polytopes which are eroded versions of the cells of the tessellation (see Fig. 2). We denote its intensity by Nv,r. Alternatively, we can write Nv,r = Pr^ , where pR is the probability that the typical cell 'survives' an erosion with a ball of radius R. For particle processes X of non-overlapping objects with intensity Ny, the additivity of the intrinsic volumes yields Vv,,{X) = NvEV,{Xo), where X) is the typical grain of X (Schneider and Weil, 2008). Let now A' C M'^ be a compact, polyconvex and topologically regular set. Then the Euler characteristics of X and the topological closure X'^ of the complement of X are related via d+l, see (Ohserand Schladitz, 2009, Formula (3.11)). Since ErDK with K e Jf fulfils the requirements listed above, the density of the Euler characteristic can be interpreted as For the interpretation of Mv,r, we remark that the ^nsities of the integral of mean curvature of X and X'= just differ by a factor -1 (Ohser and Schladitz, 2009, p. 159). Furthermore, the mean width b{K) of a convex body K e Jf is related to_ its integral of mean curvature M{K} via M{K} = 2nb{K) (Ohser and Schladitz, 2009, p. 21). Hence, we get MV,R = -2nNv,RbR, where bR is the expected mean width of the eroded cells of the tessellation. In applications, the values of Mv,r and Xv,r (or Xa,r) can be used to validate a Poisson-Voronoi model which has been fitted using Vv,r and Sv,r (or Aa,r and La,r). This is of particular interest since the estimators for the intrinsic volume densities used in this paper estimate all densities simultaneously (see also Ohser and Schladitz (2009)). Hence, characteristics for parameter estimation and model validation can be obtained in one single estimation step. For these reasons, we believe that also the intrinsic volume densities Xa,r for c/ = 2 and Mv,r and Xv,r for d =3 deserve an investigation. Since we do not see a way to obtain analytical formulae for them, they are studied by simulation. For c/ = 2, the density of the Euler characteristic was estimated from 50 realisations of a stationary Poisson-Voronoi tessellation with intensity A = 1 in a square with edge length 40. The dilated edge systems of the tessellations were discretized in binary images with 2000 X 2000 pixels. The density of the Euler characteristic was estimated from the image data using the estimator proposed in Ohser and Schladitz (2009) which is implemented in the MAVI software package (Fraunhofer ITWM, 2005); see Fig. 3 and Table 1. In the case c/ = 3, the densities of the integral of mean curvature and the Euler characteristic were estimated from 30 realisations of Poisson-Voronoi tessellations with unit intensity in a cube with edge length 10, i.e., containing 1000 cells on average. The dilated facet systems of the tessellation were discretized in binary images with 10003 voxels. Again, the densities of the intrinsic volumes were estimated from the image data using the estimators proposed in Ohser and Schladitz (2009); see Fig. 4 and Table 2. Fig. 3. Densities of area, boundary length, and Euler characteristic of a dilated Poisson-Voronoi tessellation in R2 w^^it^h int^ensity 1 = 1 as a ^unct^ion of R. Cii^cles correspond t^o the mean v^alues estimat^ed from 50 realisat^ions of t^he model. T^e envelopes plot^ted in g^r^ey a^r^e t^he minimal ^nd max^imal lvalues l^r^om t^he 50 i^ealisat^ions. Fig. 4. Densities of v^olume, sur^^ace ar^ea, int^eg^r^al of mean curvature, and Euler characteristic of a dilated Poi^so^-Vor^o^oi t^essellat^ion inR3 w^it^h int^ensity 1 = 1 as a function of R. Circles correspond to the mean lvalues est^imat^ed^r^om 30 r^ealisat^ions of t^he model. T^e envelopes plot^t^ed in g^r^ey ar^e t^he minimal and max^imal lvalues fr^om t^he 30 r^ealisat^ions. Remark 3. The limits of the intrinsic volume densities for SR for R ^ 0 can be computed from geometric characteristics of the Poisson-Voronoi tessellation. For d =2, we have Aa,O = 0, = = and where L denotes the expected total edge length of the Poisson Voronoi tessellation per unit volume. For d =3, we get = 0, 5^,0 = 25^ 5.820^5, M^o = -2nb^ -9.16U"^ and where S denotes the expected total surface are^of the facets of the tessellation per unit volume, and b is the expected mean width of its typical cell. APPLICATION In this section we apply the formulae for Vy^R and Sy^R derived above to fit a Poisson-Voronoi model to the microstructure of a closed-cell polymer foam which is used for the insulation of buildings. We analysed a CT image of the material consisting of 480 X 480 X 360 voxels with a voxel edge length of 10.21 /im. The densities of volume and surface area were estimated from a binarised version of the image using the estimators of Ohser and Schladitz (2009). The estimated values are Vy = 7.91% and Sy = 6.983 mm"^. Based on the formulae given in Theorem 4 we fitted the parameters of the dilated Poisson- Voronoi tessellation to the estimated values using a Nelder-Mead simplex optimisation procedure (Neider and Mead, 1965). The distance function which was minimised was Vv-VP{X,R) % + Sv-SP{X,R) ŠV where Vf (A, R) and , R) are the characteristics for a dilated Poisson-Voronoi tessellation with intensity X and dilation radius R. The results are X = 1.989 mm"^ and R = 9.904 /im which corresponds to = 7.08% and = 6.984 mm-^. The fit of these values is relatively close. However, the visualisations of the original sample and the model shown in Fig. 5 show a certain deviation in the cell shape. This indicates that a tessellation with more regular cell shapes might be a preferable model for this material. Unfortunately, the validation method suggested above is not applicable for this data set: it is well-known that the estimators for My and Xv are sensitive to the resolution of the image (Ohser and Schladitz, 2009). In the current image, the resolution is not high enough to fully resolve the cell walls which are slightly perforated in the binary image. While the effect on the estimates of the volume and surface density (loss of several pixels) is negligible, the estimates My = 10.61 mm"^ and Xv = —401.2 mm"-^ are highly affected. In particular, the discussion in the previous section shows that a negative value of My should be expected for the dilated facet system of a three-dimensional tessellation while Xv should be positive. In the foam example, the signs of both characteristics are reversed. Consequently, their interpretation is questionable and they cannot be used for model validation. Fig. 5. Visualisations and sectional images of the original foam (left) and the fitted Poisson-Voronoi model (right). Sample and CT imaging: R. Schlimper, Fraunhofer IMW Halle. Remark 4. //S is a stationary and isotropic random closed set in R^ then its volume and surface area density are related to the area density A^j^ and the boundary length L^ of the two-dimensional section Sn(R^x{0}) via Vy = A^ and Hence, a three-dimensional model ca^n also be fitted based on measurements in 2D images. Whenever possible, three-dimensional data should be preferred for the model fit. In the data presented above, the area fr^actions estimated from sectional images vary betw^een 6.74% and 9.41%. Hence, the choice of a single slice might lead to a significantly different estimation of the model parameters. Nevertheless, the possibility to fit a model using 2D observations is valuable for materials which are not tr^actable to 3D imaging techniques, e.g., the sinter mate^al shown in Fig. 1. DISCUSSION We have studied the densities of the intrinsic volumes of the dilated facet systems of stationary Poisson-Voronoi tessellations in R2 and R3. We obtained analytical formulae for the densities of the area and boundary length in R2 and the volume and surface area in R3. The densities of the remaining intrinsic volumes could only be studied by simulation. The results presented in this paper allow for a quick and easy fitting of stationary Poisson-Voronoi tessellations to cellular materials. Since two of the intrinsic volume densities are sufficient for the parameter estimation, the remaining ones can be used for model validation. In many applications, tessellation models generated by regular point processes might provide better model fitting results. Nevertheless, the Poisson-Voronoi tessellation could be sufficient in cases where the cell shape is of minor importance. Furthermore, thanks to its analytical tractability it plays an important role as reference model. Due to the independent placement of the cells' generators in the Poisson-Voronoi tessellation this model yields a very disordered structure. Hence, it could also be seen as an extreme case which can possibly provide bounds for characteristics of more regular tessellation models. Formulae for the spherical contact distribution function of more general Voronoi tessellations and of Poisson Laguerre tessellations were given in Heinrich (1998) and Lautensack (2007), respectively. These formulae could be used to extend the results from this paper to further tessellation models. So far, however, these formulae are given in less explicit form than the formulae for the Poisson-Voronoi tessellation which hinders their applicability. APPENDIX Tables 1 and 2 contain the means and standard deviations of Xa,r, Mvr, and Xv,R estimated by simulation. Table 1. Values of Xa,r in R2 estimated by simulation. R Xa,r Xa,r mean sd 0.025 -0.9993 0.0224 0.050 -0.9995 0.0226 0.100 -0.9976 0.0224 0.150 -0.9853 0.0210 0.200 -0.9484 0.0194 0.250 -0.8784 0.0159 0.300 -0.7682 0.0114 0.350 -0.6308 0.0094 0.400 -0.4829 0.0105 0.450 -0.3435 0.0108 0.500 -0.2257 0.0116 0.550 -0.1359 0.0098 0.600 -0.0771 0.0072 0.650 -0.0406 0.0044 Table 2. Values of Mvr and CVR in R3 estimated by simulation. R Mvrr Mvrr XV,R XV,R mean sd mean sd 0.025 -8.9695 0.2039 1.0025 0.0335 0.050 -8.5529 0.1901 1.0022 0.0347 0.100 -7.7027 0.1598 1.0026 0.0341 0.150 -6.8224 0.1268 1.0014 0.0348 0.200 -5.8910 0.0928 0.9989 0.0329 0.250 -4.8739 0.0550 0.9758 0.0327 0.300 -3.7635 0.0250 0.9218 0.0283 0.350 -2.6042 0.0390 0.7827 0.0193 0.400 -1.5528 0.0516 0.5707 0.0141 0.450 -0.7594 0.0484 0.3366 0.0161 0.500 -0.2938 0.0311 0.1542 0.0138 0.550 -0.0872 0.0160 0.0532 0.0069 0.600 -0.0199 0.0082 0.0134 0.0043 ACKNOWLEDGEMENTS This work was supported by the German Research Foundation (DFG) under grant RE 3002/1-1. REFERENCES Coster M, Arnould X, Chermant JL, Moataz AE, Chartier T (2005). A microstructural model by space tessellation for a sintered ceramic: cerine. Image Anal Stereol 24:105-16. Fraunhofer ITWM (2005). MAVl - Modular Algorithms for Volume Images. Heinrich L (1998). Contact and chord length distribution of a stationary Voronoi tessellation. Adv Appl Prob 30:603-18. Lautensack C (2007). Random Laguerre tessellations. Ph.D. thesis, Universitat Karlsruhe. Weiler bei Bingen: Verlag Lautensack. Muche L, Stoyan D (1992). Contact and chord length distributions of the Poisson Voronoi tessellation. J Appl Probab 29:467-71. Nelder J, Mead R (1965). A simplex method for function minimization. Comput J 7:308-13. Ohser J, Schladitz K (2009). 3D images of materials structures - Processing and analysis. Weinheim: Wiley. Okabe A, Boots B, Sugihara K, Chiu SN (2000). Spatial tessellations - Concepts and applications of Voronoi diagrams, 2nd ed. Chichester: Wiley. Redenbach C (2009). Modelling foam structures using random tessellations. In: Capasso V, Aletti G, Micheletti A, eds. Stereology and Image Analysis. Proc 10th Eur Conf ISS (ECS10), vol. 4 of The MlRlAM Project Series. Bologna: Esculapio. Ribeiro-Ayeh S (2005). Finite Element Modelling of the Mechanics of Solid Foam Materials. Ph.D. thesis, Kungliga Tekniska Hoogskolan, Stockholm. Schneider R (1993). Convex bodies: the Brunn-Minkowski theory. Cambridge: Cambridge University Press. Schneider R, Weil W (2008). Stochastic and integral geometry. Berlin: Springer. Stoyan D, Kendall WS, Mecke J (1995). Stochastic geometry and its applications, 2nd ed. Chichester: Wiley. Telley H, Liebling TM, Mocellin A (1996). The Laguerre model of grain growth in two dimensions, l. Cellular structures viewed as dynamical laguerre tessellations. Philos Mag B 73:395-408.