Cup products on polyhedral approximations of 3D digital images

Cup Products on Polyhedral Approximations of 3D Digital Images Rocio Gonzalez-Diaz1, Javier Lamar2, and Ronald Umble3 1Dept. of Applied Math (I), School of Computer Engineering, University of Seville, Campus Reina Mercedes, C.P. 41012, Seville, Spain [email protected] 2Pattern Recognition Department, Advanced Technologies Application Center, 7th Avenue #21812 218 and 222, Siboney, Playa, C.P. 12200, Havana City, Cuba [email protected] 3Department of Mathematics, Millersville University of Pennsylvania, P.O. Box 1002 Millersville, PA 17551-0302, Pennsylvania, USA [email protected] Abstract. Let Ibe a 3D digital image, and let Q(I) be the associated cubical complex. In this paper we show how to simplify the combinatorial structure of Q(I) and obtain a homeomorphic cellular complex P(I) with fewer cells. We introduce formulas for a diagonal approximation on a general polygon and use it to compute cup products on the cohomology H∗(P(I)). The cup product encodes important geometrical information not captured by the cohomology groups. Consequently, the ring structure of H∗(P(I)) is a finer topological invariant. The algorithm proposed here can be applied to compute cup products on any polyhedral approximation of an object embedded in 3-space. Keywords: Cellular complex, cohomology, cup product, diagonal approximation, digital image, polyhedron. 1 Introduction Throughout this paper, coefficients lie in the field Z2.LetXbe a cellular complex embedded in 3-dimensional space and constructed by gluing 3-dimensional polyhedra together along common faces (see [4]). At a most basic level, the connected components, homotopy classes of non-contractible loops, and boundaries of tunnels in Xgenerate the cellular cohomology H∗(X). At the next level, certain relationships among the generators are encoded by the cup product, which endows H∗(X) with a graded commutative ring structure. Indeed, the discriminating information encoded by the cup product improves our capability to distinguish between 3D images. For example, H∗(S1∨S1∨S2)andH∗(S1×S1) are isomorphic as vector spaces but not as rings since cup products vanish in the wedge but not in the product. Thus S1∨S1∨S2and S1×S1have quite different topological properties. To date, the cup product has seen limited application to problems in 3D image processing. In [10,11], Gonzalez-Diaz and Real used their 14-adjacency algorithm J.K. Aggarwal et al. (Eds.): IWCIA 2011, LNCS 6636, pp. 107–119, 2011. c Springer-Verlag Berlin Heidelberg 2011 108 R. Gonzalez-Diaz, J. Lamar, and R. Umble Fig. 1. Left: A digital image I=(Z3,26,6,B); the set Bconsists of 8 unit cubes (voxels). Right: The quadrangles of ∂Q(I). and the standard formulation in [17] to compute cup products on the simplicial complex K(I) associated with a given digital image I. More recently, GonzalezDiaz, Jimenez and Medrano introduced a method for computing cup products on cubical approximations Q(I). Their cup products are computed directly from the cubical complex, and no additional subdivisions are necessary [8,9]. For a geometrical interpretation of cohomology in the context of digital images, we refer the reader to [5,6,15]. In [14], Kravatz computed cup products on a general 2-dimensional polygon in terms of a combinatorial diagonal approximation, which assumes a particular ordering of the vertices. In this paper, we introduce a more general formula for computing cup products, which is independent of the ordering of vertices and computationally effective. A problem that frequently arises in 3D image processing is to efficiently encode the boundary surface of a given digital object as a set of voxels. The most popular approach to this problem uses a triangulation. While triangles are combinatorially simple, and visualization of triangulated surfaces is supported by existing hardware and software, the number of triangles required is often large and the computational analysis correspondingly slow. It is desirable, therefore, to seek more computationally economical combinatorial approximations. Adjacent coplanar triangles in a triangulation, for example, can be merged into more general polygons and become faces of more general but combinatorially simpler polyhedra. The payoff from combinatorial simplicity is improved computational efficiency. Approximating 3D objects with polyhedral complexes is a well-studied problem in the field of Computational Geometry (for example, see [1,2,3]). An algorithm for constructing polyhedral approximations in certain special cases was given by Kovalevsky and Schulz in [13,19]. Their algorithm generates the convex hull of a given object then modifies the convex hull by recursively generating convex hulls of either subsets of the given voxel set or subsets of the background voxels. The result of this method is a polyhedron that separates object voxels from background voxels. The computational methods introduced in this paper can be effectively applied to any polyhedral approximation of a 3D object. Indeed, one maximizes computational efficiency by approximating a given 3D object with a polyhedral Cup Products on Polyhedral Approximations of 3D Digital Images 109 Fig. 2. Left: Quadrangles in the boundary of cubes cand σsharing a square σ(in bold). Right: Quadrangles in ∂(c):=∂(c+σ), the boundary of the cell cafter removing σ. complex containing a minimal number of cells. To the extent that this is our long-term objective, we take a first step in this direction here. The paper is organized as follows: In Section 2 we introduce a simplification procedure, which produces a cellular complex P(I) homeomorphic to Q(I)with significantly fewer cells. In Section 3 we define a diagonal approximation on a general polygon and use it to compute the cohomology ring of P(I). Conclusions and some ideas for future work are discussed in Section 4. 2 3D Digital Pictures and Cellular Complexes Let Ibe a 3D digital image and let Q(I) be an associated cubical complex. In this section we introduce a simplification procedure, which produces a cellular complex P(I) homeomorphic to Q(I) with significantly fewer cells. Intuitively, a cellular decomposition of a 3D space Xembedded in R3is a representation of Xas a finite union of vertices (0-cells), edges (1-cells), polygons (2-cells), and polyhedra (3-cells), which have been glued together in such a way that the non-empty intersection of two cells is a cell. A k-cell is also referred to as a k-face.Acellular complex is a 3D space Xembedded in R3together with a cellular decomposition. For a precise definition of a cellular complex, which is more subtle than one might expect, see [4]. Acubical complex Qis a cellular complex whose 2-cells are squares (or quadrangles) and whose 3-cells are cubes. Note that if a cube is in Q, its bounding quadrangles are in Q; if a quadrangle is in Q, its bounding edges are in Q;and if an edge is in Q, its endpoints are in Q. Consider a 3D binary digital picture I=(Z3,26,6,B), where Z3is the underlying grid and B(the foreground) is a finite set of points of the grid fixing the 26-adjacency for the points of Band the 6-adjacency for the points of Z3\B (the background). The cells of Q(I) are unit cubes centered at the points of B with faces parallel to the coordinate planes (called the voxels of I), together with their quadrangles, edges, and vertices. Let Kbe a cellular complex. An i-cell σ∈Kis a facet of a cell σ∈Kif σis an (i+ 1)-cell and σis a face of σ.Amaximal cell of Kis not a facet of any cell of K. The boundary of K, denoted by ∂K, is the subcomplex of Kconsisting of all cells that are facets of exactly one (maximal) cell, and their faces. Note that 110 R. Gonzalez-Diaz, J. Lamar, and R. Umble Fig. 3. Critical configurations (i), (ii) and (iii) (modulo reflections and rotations) the maximal cells of ∂Q(I) are all the quadrangles of Q(I)sharedbyavoxelof Band a voxel of Z3\B(see Figure 1). Following the exposition in [8,9], given a digital image Iand its associated cubical complex Q(I), we apply a face-reduction technique to reduce the number of cells in Q(I)\∂Q(I) and obtain a cellular complex K(I) homeomorphic to Q(I) whose maximal cells are the quadrangles of ∂Q(I)(seeFigure2and Algorithm 1). Then, ∂K(I)=∂Q(I). Input: A cubical complex Q(I)associated to a 3D digital image I. Initially, K(I):=Q(I). While there exists a cell σ∈Q(I)\∂Q(I)do If σis a facet of exactly two cells c, σ ∈Q(I)do remove σand σfrom the current K(I); redefine cas c∪σ. If σis a facet of exactly one cell σ∈Q(I)do remove σand σfrom the current K(I). Output: the cellular complex K(I). Algorithm 1. Face-Reduction Process Next, we preform a simplification process in ∂K(I) to produce a cellular complex P(I) homeomorphic to K(I) such that the maximal cells of ∂P(I)are polygons. But first, we need a definition. Definition 1. A vertex v∈∂K(I)is critical if one of the following situations occurs: (i) vis a face of some edge eshared by four cubes, exactly two of which intersect along eand lie in Q(I)(see cubes w1and w2in Figure 3). (ii) vis shared by eight cubes, exactly two of which are are corner-adjacent and contained in Q(I)(see cubes s1and s2in Figure 3). (iii) vis shared by eight cubes, exactly two of which are corner-adjacent and not contained in Q(I)(cubes t1and t2in Figure 3). Cup Products on Polyhedral Approximations of 3D Digital Images 111 Fig. 4. (a) Nv←{q1,q2,q3,q4},(b)facetsofp←{e1,e2,e3,e4,e5,e6,e7,e8} It is a well-known fact that a non-critical vertex of ∂K(I) lies in a neighborhood of ∂K(I) homeomorphic to R2(see [16]). Algorithm 2 processes the non-critical vertices of ∂K(I) to obtain the cellular complex P(I). Initially, P(I)=K(I). For a vertex v∈∂P(I), let Nvbe the set of 2-cells q∈P(I)incidenttothevertexv.IfNvdefines a region Rvhomeomorphic to a disc, then Nvis replaced by a new 2-cell pin P(I), which is the union of the cells of Nv.Theedgesof∂P(I)incidenttovand the vertex vare removed from P(I) (see Figure 4). Observe that the maximal cells of the final cellular complex ∂P(I) are polygons and ∂P(I) has fewer cells than ∂K(I). We can set some terminating conditions. For example: (1) terminate when the number of edges of the polygons in ∂P(I) reach some specified maximum; or (2) terminate after merging the set Nvof coplanar 2-cells of ∂P(I) (this preserves the geometry but removes fewer cells). An example of the differences that arise from these different terminating conditions is demonstrated in Example 1. Observe that Algorithm 2 uses the ordering on the set of non-critical vertices V⊂∂K(I) to select the next non-critical vertex. To the best of our knowledge, this is the first algorithm to appear that produces a cellular complex with polygonal maximal cells by removing non-critical vertices. Example 1. Let μI be a μMRI of a trabecular bone of size: 85 ×85 ×10 voxels (see Figure 5 in which μI is given by a sequence of 10 2D digital images of size Input: The output of Algorithm 1: the cellular complex K(I). Initially, P(I):=K(I); V:= ordered set of non-critical vertices of ∂K(I). While ∃v∈Vsuch that Rvis homeomorphic to a disc do remove vfrom P(I)and V; remove the edges incident to vfrom P(I); remove the 2-cells of Nvfrom P(I); add a new 2-cell pto P(I)which is the union of the cells of Nv. Output: The cellular complex P(I). Algorithm 2. Algorithm to obtain the cellular complex P(I) 112 R. Gonzalez-Diaz, J. Lamar, and R. Umble Fig. 5. AμMRI of a trabecular bone 85 ×85). The number of quadrangles in ∂Q(μI) is 20956 (see Figure 6). After applying Algorithm 1 to obtain the cellular complex K(μI), we apply Algorithm 2toK(μI) with the terminating condition 1 (the number of the edges of the polygons in ∂P(μI) is smaller or equal to 10). Then, the number of polygons of ∂P(μI) is 1567 (see Figure 7.a). If we only consider the set of coplanar polygons Nv, then the number of polygons of ∂P(μI) after applying Algorithm 2 is 9321 (see Figure 7.b). 3 Computing the Cohomology Ring of P(I) Traditionally, one computes cup products in simplicial or cubical complex using the standard formulas in [17,20]). In this section, we give a procedure for computing cup products on P(I) (the output of Algorithm 2), which avoids triangulation by defining explicit formulas for diagonal approximations on polygons. We begin with a review of some standard definitions from Algebraic Topology (for details see [17]). Given a graded set S={Sq}q,theq-chains of S,whichare finite formal sums of elements of Sq, define an additive abelian group structure on Sq. These groups, called q-chain groups, are denoted by Cq(S). The collection of all chain groups associated with Sis denoted by C∗(S)={Cq(S)}q Fig. 6. The cubical complex ∂Q(μI) Cup Products on Polyhedral Approximations of 3D Digital Images 113 andisreferredtoasthechain group of S.Achain complex (C∗(S),∂)is a chain group C∗(S) together with a square zero homomorphism ∂={∂q: Cq(S)→Cq−1(S)}q,called the boundary operator. For example, consider a triangle vi,v j,v kwith vertices vi<v j<v k. The boundary of the triangle is the formal sum of its edges, that is, ∂2(vi,v j,v k)=vi,v j+vj,v k+vi,v k. Note that any chain group C∗(S) together with the zero boundary map ∂≡0 is a chain complex. The chain complex associated with P(I) (the output of Algorithm 2) is the collection (C∗(P(I)),∂)={Cq(P(I)),∂ q}qwhere: •each Cq(P(I)) is the chain group generated by the q-cells of P(I), •the boundary ∂q:Cq(P(I)) →Cq−1(P(I)) evaluated on a q-cell of P(I)is the formal sum of its facets, and •the boundary of a general q-chain is defined by linearly extending ∂. Givenachaincomplex(C∗(S),∂), a q-chain σ∈Cq(S) is called a q-cycle if ∂q(σ)=0.Ifσ=∂q+1(μ)forsome(q+1)-chainμthen σis called a q-boundary. Referring to the triangle vi,v j,v kabove, σ=vi,v j+vj,v k+vi,v kis both a 1-cycle and a 1-boundary since ∂1(σ)=0andσ=∂2(vi,v j,v k). Two q-cycles aand aare homologous if there exists a q-boundary bsuch that a=a+b. Denote the groups of q-cycles and q-boundaries by Zq(S)and Bq(S), respectively. All q-boundaries are q-cycles (Bq(S)⊆Zq(S)). Define the qth homology group to be the quotient group Hq(S)=Zq(S)/Bq(S), for all q. Each element of Hq(S)isaclass[a]=a+Bq(S)andais a representative q-cycle. The homology of Sis the collection of all the homology groups associated with S, i.e., H∗(S)={Hq(S)}q. Let (C∗(S),∂)and(C∗(S),∂) be chain complexes. A homomorphism f= {fq:Cq(S)→Cq(S)}qsuch that fq∂q=∂ qfqfor all qis a chain map. Note that the identity idC∗(S)={idCq(S):Cq(S)→Cq(S)}qis a chain map. Let f={fq:Cq(S)→Cq(S)}qand g={gq:Cq(S)→Cq(S)}qbe chain maps. A chain homotopy from fto gis a homomorphism φ={φq: Cq(S)→Cq+1(S)}qsuch that φq−1∂q+∂ q+1φq=fq+gqfor all q. Achain contraction of (C∗(S),∂)to (C∗(S),∂) is a triple (f={fq:Cq(S)→Cq(S)}q, g={gq:Cq(S)→Cq(S)}q,φ={φq:Cq(S)→Cq +1(S)}q) such that (i) fand gare chain maps; (ii) φis a chain homotopy from idC∗(S)to gf ={gqfq:Cq(S)→Cq(S)}q; (iii) fg ={fqgq:Cq(S)→Cq(S)}q=idC∗(S). Cochain groups are the linear duals of chain groups. Given a chain complex (C∗(S),∂), a q-cochain c∈Hom(Cq(S),Z/2). If we index the q-cells in a cellular complex from 1 to nq, their corresponding duals generate C∗(S). Thus a cochain c∈C∗(S)isaZ2-linear combination of the nqelements in the dual basis, and as such can be thought of as a bit string of length nq. The set Cq(S)ofallq-cochains is a group, and the direct sum of all cochain groups associated with Sis the graded group C∗(S)={Cq(S)}q.Thecoboundary operator δ={δq:Cq(S)→Cq+1(S)}qis defined on a q-cochain cby δq(c)= c∂q+1.Notethatδ◦δ= 0. The associated cochain complex is the pair (C∗(S),δ). 114 R. Gonzalez-Diaz, J. Lamar, and R. Umble Fig. 7. (a) The cellular complex ∂P(μI) for 10 edges as upper bound on p∈∂P(μI). (b) The cellular complex ∂P(μI) preserving geometry. Aq-cochain cis a q-cocycle if δq(c)=0.Aq-cochain bis a q-coboundary if there exists a (q−1)-cochain csuch that b=δq−1(c). Two q-cocycles cand care cohomologous if there exists a q-coboundary bsuch that c=c+b(see Figure 8). We denote the subgroup of q-cocycles by Zq(S), and the subgroup of qcoboundaries by Bq(S). The qth cohomology group is defined to be the quotient Hq(S)=Zq(S)/Bq(S). Each element of Hq(S)isaclass[c]=c+Bq(S). The element cis a representative q-cocycle of the cohomology class [c]. The cohomology of Sis the graded Z/2-vector space H∗(S)={Hq(S)}q. Since P(I) (the output of Algorithm 2) is embedded in R3,homology and cohomology of P(I) are isomorphic and torsion free. An AT-model [10,11] for a chain complex (C∗(S),∂), denoted by ((S, ∂),H, f,g,φ), consists of a chain complex (C∗(H),∂≡0) together with a chain contraction (f,g,φ)of(C∗(S),∂)to(C∗(H),∂). The following properties hold: •If σ∈Hq,then gq(σ)∈Cq(S) is a representative cycle of a class of Hq(C(S)). •The cochain ∂σf:Cq(S)→Z/2 defined by ∂σf(μ):=1,if σappears in the expression of f(μ), 0,otherwise; is a representative cocycle of a class of Hq(S). •The map Hq→Hq(S)givenbyσ→ [g(σ)] linearly extends to an isomorphism Cq(H)∼ =Hq(S). •The map Hq→Hq(S)givenbyσ→ [∂σf] linearly extends to an isomorphism Cq(H)∼ =Hq(S). Example 2. Consider the cellular complex ∂P(μI) shown in Figure 7.b obtained after applying Algorithm 2 to a μMRI of a trabecular bone of size 85 ×85 ×10 Cup Products on Polyhedral Approximations of 3D Digital Images 115 0-cochain {v1} 1-cochain {e1,e 4} 1-coboundary δ{v1}={e3,e 4} 1-cocycle c={e1,e 2} 1-cocycle d={e1,e 2,e 3,e 4} homologous cocycles cand d;sinced=c+δ{v1} Fig. 8. Example of cochain, cocycle and coboundary voxels. Table 1 shows the results of the homology computation, i.e, the number of connected components, holes and cavities obtained after computing an ATmodel for ∂P(μI). Representative 1-cycles are shown in Figure 9. Table 1. Results of the homology groups computation for the cellular complex ∂P(μI) shown in Figure 7.b. (see Figure 9). Cellular Complex H0H1H2 ∂P(μI)526 An AT-model for (C∗(S),∂) always exists and can be computed in O(m3), where mis the number of elements of S(see [10,11]). GivenanAT-model((P(I),∂),H,f,g,φ)forP(I) (the output of Algorithm 2), we have that C∗(H)∼ =H∗(P(I)) ∼ =H∗(P(I)) ∼ =Hom(H∗(P(I)),Z/2) ∼ =C∗(H). Let α∈Hn. Consider the dual elementary n-cocycle in Cn(H), α∗:Cn(H)→Z/2 such that for μ∈Hn,α ∗(μ):=1ifμ=α, 0otherwise. Proposition 1. Given an ordering {v1<··· <v n}of the vertices of P(I), each polygon p∈P(I)can be expressed as an ordered list of vertices {vi1< ···<v ik}⊆{v1<···<v n}with edges ej:= vij,v ij+1 ,if j<k, vik,v i1,if j=k. The following two theorems formulate a diagonal approximation ∇on a polygon and the cup product on H∗(P(I)) in terms of ∇. All non-trivial cup products in H∗(P(I)) are products of 1-cocycles for dimensional reasons. Theorem 1. Consider a polygon p=v1,...,v nwith edges ei=vi,v i+1, i<n,anden=vn,v 1. Then a diagonal approximation on pis given by ∇(p):=1<i<n, vi<vi+1 (e1+···+ei−1)⊗ei +1<i<n, vi>vi+1 (ei+1)+···+en)⊗ei.