Simple and Robust Boolean Operations for Triangulated Surfaces

Boolean operations of geometric models is an essential issue in computational geometry. In this paper, we develop a simple and robust approach to perform Boolean operations on closed and open triangulated surfaces. Our method mainly has two stages: (1) We firstly find out candidate intersected-triangles pairs based on Octree and then compute the inter-section lines for all pairs of triangles with parallel algorithm; (2) We form closed or open intersection-loops, sub-surfaces and sub-blocks quite robustly only according to the cleared and updated topology of meshes while without coordinate computations for geometric enti-ties. A novel technique instead of inside/outside classification is also proposed to distinguish the resulting union, subtraction and intersection. Several examples have been given to illus-trate the effectiveness of our approach.


Introduction
The calculating of Boolean operations for solids plays an important role in geometric processing, which intends to obtain union, subtraction and intersection of geometric models.The most commonly used method to display geometric models probably is the boundary representation (B-Rep) [11] which can be based on parametric surfaces such as NURBS surfaces or discrete surfaces.The discrete surface is the decomposition of a specific surface domain, of which polygonal meshes especially triangular meshes is the most popular kind of representation form.
In this paper, we focus our contribution on the implementation of Boolean operations over triangular meshes.As mentioned above, triangular meshes have become the most popular surface representation.They are normally required to be manifold in most applications, however this property may be lost in some cases such as undergoing dynamical collision or large deformation.Thus reconstruction need to be carried out over those meshes to keep them manifold [26].In this paper, we perform simple and robust Boolean operations on a pair of manifold triangulated surfaces without considering the case of self-intersecting meshes.
In our algorithm, surfaces are classified into two types: open or closed.Open surfaces have outer or inner boundary loops while closed ones without.When a surface is closed, the geometric model bounded by it can be called as a volume or a block.The geometric objects we will manipulate are a pair of manifold triangulated surfaces which can be open-and-open, open-andclosed or closed-and-closed.

Related works
Varies of implementations of Boolean operations are described in the literatures.They can be roughly classified according to three main properties: the type of input data, the type of computation and the type of output data [18].About the type of computation, Tayebi et al. [22] assort them into four different categories: exact arithmetic methods [2], approximate arithmetic methods [21], interval computation methods and volumetric methods [4,24].
Most of approaches described in the literatures specially aim at one type of the input models such as those by B-Rep based on NURBS surface [22,25], curved surface [14] and polygonal meshes [2,18,20].Besides these popular surface-based representations, several other methods like Nef [9], surfel-bounded [1], L-Rep [27] for representing geometric models are also developed to construct new models or convert existing ones to target models, and then Boolean operations are also implemented.In another aspect, the input models can be manifold or non-manifold [19], and most methods perform Booleans on manifold ones.
Among the types of computation, the exact arithmetic methods and interval computation methods need directly computing Booleans over the initial surfaces, while the rest of two do not and hence deem to be indirect.
When implementing direct Booleans for mesh models, there are two key procedures which strongly affect the effectiveness and efficiency of the whole algorithm.The first one is how to robustly obtain the intersection lines and loops of all intersected entities as faster as possible.The core of such procedure is to accurately find out all potentially intersected entities in a short time to reduce computation cost.Many techniques based on BSP(binary space partitions) [2], Octree [18], OBB trees [12], bipartite graph structure [28] and tracking neighbors [8,14] have been developed to realize this goal.The second key procedure is how to correctly assemble and distinguish the union, subtraction and intersection of intersected models.The most direct method is to make a inside/outside classification [3], which checks the location of vertices or facets with the resulting volumes.Varies of algorithms of direct Boolean operations develop their own solutions for the above two problems.
Quite recently, several novel methods of Boolean operations on polygonal meshes have been presented.Pavić et al. [18] develop a hybrid method which combines polygonal and volumetric computations and representations for performing Boolean operations over polygonal meshes; similarly García et al. [7] use Extended Simplicial Chains(ESCs) to represent both boundary and volume of free-form solids when perform Booleans.Aiming at robustness of Boolean operations for large number of triangles in industrial applications, Schifko et al. [20] adopt several suitable libraries from CGAL and GNU multi precision arithmetic library to filter exact arithmetic.An approximate method based on Layered Depth Images(LDI) for polygonal models is developed in [24], in which LDI is accepted to sample and trimmed adaptive contouring is used to rebuild intersected meshes.Campen and Kobbelt [2] adopt plane-based representations, BSP techniques and a localization scheme to obtain exact and robust (self-) intersections for polygonal meshes.

Our contribution
When designing algorithms to solve a specific problem, normally a fast algorithm with better efficiency also has higher complexity than a slow one.This is certainly true for performing Boolean operations on triangular meshes.In this paper, we will try to keep a relative balance between the efficiency and complexity of such algorithms.And our objective is to develop a simple and robust approach of Boolean operations which also have satisfied efficiency.Several highlights can be drawn in our method: (1) A simple and relative fast method to calculate the intersection lines of triangles.In order to reduce computation cost, we firstly use the spatial decomposition data structure Octree to locate and search candidate intersected triangle-pairs, and then compute the intersection line of each pair of triangles parallel.
(2) A fast approach only based on operations of the indexes of entities and without computation for coordinates to obtain intersection loops, sub-surfaces and sub-blocks.
(3) A novel technique to distinguish the union, intersection and subtraction volume of two intersected closed surfaces(solid models).This simple classification method is very fast and robust since it is also only based on operations of the indexes of entities.
The rest of this paper is organized as follows.In Sect.2, we give a brief view to our algorithm; in Sect.3 we design some data structure and define several notations for our implementation; and then, all details of the approach will be described in Sect.4;Then, we give several examples to illustrate the effectiveness of the approach in Sect.5, and finally conclude and discuss our work in Sect.6.

Overview of the approach
Our method can be mainly summarized into 6 steps (Fig. 1), and can be divided into two states: the first, which is based on computation for coordinates of entities, includes searching intersected triangle-pairs, calculating intersection line of each pair of triangles, re-triangulation and updating the resulting surface meshes; the second procedure is only based on operations for indexes of entities including forming intersection loops, creating sub-surfaces, assembling and distinguishing sub-blocks.Step 1: Search the candidate intersected triangle-pairs.In order to reduce the cost on computing, a robust and fast search algorithm is needed to obtain potentially intersected triangles.In this paper, we use Octree [6] to locate and find out all candidate intersected triangle-pairs.
Step 2: Compute the intersection line for each pair of triangles.Möller [16] developed a robust and efficient algorithm to do this and we adopt his work.When there are a huge number of triangle-pairs needed to calculate their intersections, parallel computation based on OpenMP [17] is accepted.
Step 3: Merge and renumber all vertices, update and clear meshes.New vertices will be produced after intersecting of triangle-pairs, and the intersected triangles are replaced by retriangulation.In order to keep a valid topology, all vertices are merged and renumbered and all triangles are updated.
Step 4: Connect intersection lines into closed or open loops.After computing the intersection of each pairs of triangles, a set of discrete edges can be obtained, and they need to be connected into closed or open orientated loops.If there not exists at least one closed loop on an intersected surface, then no closed block bounded by triangular facets will be formed.
Step 5: Obtain sub-surfaces based on closed loops.A sub-surface includes the closed loop and all of its incident triangles.The edges of a closed loop are set as the advancing front; "grow" a new surface based on the topology until the number of faces in the sub-surface not increases (Fig. 4(b)).
Step 6: Assemble and distinguish sub-blocks.Sub-blocks can be easily created by assembling related sub-surfaces (Fig. 4(c)).In further, the boundary closed loops generated in Step 4 of sub-surface can be used to represent the sub-surface.Hence, assembling and distinguishing can be done based on the boundary closed loops.

Data structure and notation
In this section, some notations for geometric entities are defined.Additional properties such as direction are added onto some of them.
Definition 1 Directed edge is a directed segment, in which the first vertex is named as Head, while the other is Tail.
Definition 2 Orientated loop is a set of connected directed edges, which can be closed or open.It can be also represented by a set of ordered vertices.A pair of orientated loops are defined as twins if they have same vertices and opposite order.Definition 3 Normalized face is a coplanar polygon with its normal, which is triangle or polygon in this paper.
To describe our algorithm in following sections, we notate several common geometric objects and allocate arrays to store geometric entities.m aVerts: an array of vertex to store all vertices after intersection of all triangle-pairs and merging; all vertices must be checked and renumbered.m aEdges: an array of edge to store all intersection lines after intersecting of all trianglepairs; the Head and Tail of each edge in this array are updated after merging and renumbering all vertices.
m aLoops: an array of loop to store all closed/open intersection loops.The set of ordered vertices in each loop must have the newly updated IDs.m aPolys: an array of polygon to store all resulting polygons produced from the intersecting of triangle-pairs.Each polygon in this array will be decomposed into new triangles via polygon triangulation.
m aTrgls: an array of triangle to store all updated triangles after intersecting and merging; the triangles come from (1) the original triangles that do not intersect with others, or (2) the new triangles generated by re-triangulating the resulting polygons from the intersecting of trianglepairs.m aSurfs: an array of surface to store all created sub-surfaces.Each sub-surface has its boundary sub-loop(s) preparing for assembling sub-blocks.m aBlocks: an array of block to store all assembled sub-blocks.

The Methodology
In this section, we will present the detailed description to all six steps in our approach separately.Several procedures such as clearing triangular mesh, creating sub-surfaces and sub-blocks will be explained with pseudo codes.

Searching intersected triangle-pairs
Before calculating the intersection of any pair of triangles, we should know which pair of triangles potentially intersect.The direct and inefficient method is to make a bounding box intersection test for each triangle in a surface with that of another one.In order to improve the efficiency, we adopt Octree to locate and then find out candidate intersected triangle-pairs.Given two surface meshes, S A and S B , compute their smallest AABBs denoted as Box A and Box B , respectively; and then calculate the intersection Box AB of Box A and Box B (Box AB = Box A ∩Box B ); check each triangle of S A and S B whether it is outside of the volume Box AB , and divide S A and S B into two sub-arrays where S Aout + S Ain = S A , S Bout + S Bin = S B ; and then extend the volume Box AB into a cube to be an AABB for the triangles from both S Ain and S Bin (S Ain ∪S Bin ).
Let the bounding cube as the root node of an Octree, and create eight octants for each node recursively.Let the N a and N b denote the number of triangles that intersect each interior node from S Ain and S Bin respectively, the recursion of creating 8 octants for each node terminates and then the node becomes a leaf when (1)The depth of the node reaches a user-specified maximum depth; (2)Both N a and N b less than a given allowable number; (3)N a or N b is zero.
When to check whether a triangle from either S Ain or S Bin is inside a node of the Octree, a simple method is to make an intersection test between the bounding box of the triangle and the node of Octree.if they intersect, the triangle can be considered being inside the node.Noticeably, a same triangle can intersect with several nodes.To reduce the cost on computing bounding box for each triangle, we previously calculate the boxes for all triangles from both S Ain and S Bin once, and then adopt the boxes when need.

Intersecting of triangles and re-triangulating
The intersecting of a pair of triangles is not complex, Möller [16] has developed a robust and efficient algorithm and its code to do this.In this paper, we directly adopt the work from Möller.
Our interest is that, when there are still a huge number of pairs of triangles needed to calculate their intersections, parallel algorithm is meaningful and necessary here.Computing the intersection of one pair of triangles does not affect the procedure of that to another pair of triangles.The intersecting for any two pairs of triangles can be carried out at the same time.This parallelization can be realized easily by calling the OpenMP API [17] and adding several very simple routines into the common serial codes.We use the following code structure to compute the intersections of triangle-pairs parallel.
#pragma omp parallel for for each pair of triangles p i { Calculate the intersection of p i ; Save the intersection edge of p i if it exists; } After calculating the intersection for all triangle-pairs, for an individual intersected triangle, usually there are several intersection edges inside the triangle since it perhaps intersects with several other triangles, and it will be divided into several sub polygons by the intersection edges.To generate these sub polygonal faces for a triangle T i , we firstly find all of its intersection edges, and then connect the edges into one or more open or closed loops , and finally divide the newly produced polygons by a unused intersection loop iteratively until all loops of the triangle T i are adopt, as shown in Fig. 2.
(d) As described above, intersected triangles will be divided into polygonal faces.In order to manipulate the surface mesh more easily in next steps, all resulting polygonal face should be decomposed into triangles via ear-clipping algorithm [10].Since that the ear-clipping is valid for planar and count-clockwise polygons, pre-process must be carried out before re-triangulating a polygon in 3D.For a polygon P , its partition can be obtained in 3 steps: Step 1: Calculate the local coordinate system of P , and then transfer P into its local version P ; Step 2: Check whether the vertices of P orders in count-clockwise (CCW) or clockwise (CW).If CW, reverse P to make it be CCW; Step 3: Generate the polygon triangulation T of P via ear-clipping, and then directly obtain T for P according to P for that the topologies of T and T are exactly the same while the coordinates of vertices differ(Fig.2(d)).

Merging and updating
After intersecting of triangle-pairs, new vertices will be produced and the original intersected triangles are replaced by re-triangulated triangles.In order to get ready for next operations that based on valid topology, the surfaces must be merged and updated.It needs to carry out the following clearings: (1) merge all vertices and renumber them; (2) update the indexes of the vertices for each triangle and loop; (3) check each triangle and loop whether there are same indexed vertices; (4) reverse some newly produced triangles (Fig. 3).
Several requirements must be conformed after clearing: no same vertices, no degenerated triangles, no same vertices in a loop and no same edges.For that the next procedures such as connecting loops and creating sub-surfaces are all based on a valid topology of meshes; thus it is necessary to clear the topology besides merging and renumbering vertices.As shown in Fig. 3, (a) is the original triangular meshes with 3 unallowable faces for that they have same edges with their adjacent triangles and (b) is the cleared version.

Forming intersection loops
As defined in Sect.3,Orientated loop is a set of connected directed edges, which can be closed as a cycle or open as a polyline.Open loop is the intersection loop in which either the first or the last vertex is only shared by one intersection edge while each of the rest is shared by two edges.Close loop is the intersection loop in which each vertex is shared by at least two intersection edges.we classify the closed intersection loops into: hard closed and soft closed: Hard closed loop is the one in which each vertex is shared by only two intersection edges; for example, there are six hard closed loops in Fig. 4.
Soft closed loop is the one in which the first and last vertices are shared by more than two intersection edges while each of the rest is shared by two edges, see 4 soft closed loop in Fig. 9.
Only closed intersection loops(hard closed and soft closed ones) can be used to create subsurfaces (Sect.4.5).The closed or open loops can be assembled by strictly connecting the Head of an edge with the Tail of another edge, or the Tail of an edge with the Head with another edge.The iteration of assembling will stop when no more connected edges can be found.Besides being displayed by a set of edges, the resulting closed/open loops can be also represented by a set of ordered vertices.For an original surface, there may be several loops on it; and we should know whether a loop is closed or open simply by comparing the Head of the first edge with the Tail of the last edge.

Creating sub-surfaces
Assuming there is at least one closed loops in an original surface, then sub-surfaces can be created beginning from the closed loop.A sub-surface includes the closed loop and its incident triangles.Set the edges of a closed loop as the advancing edge front, "growing" until the number of faces in the sub-surface not increases.The growing can be done according the topology of the updated surface.Finding the next face is to search the incident faces of an advancing front edge.After a growing step, all advancing front must be updated to prepare for next growing until a sub-surface is formed; when all closed loops are used, the growing can terminate(Algorithm 1).while the number of triangles of newSf increases do

5:
for edge e j of loop L i do 6: let egHead and egT ail be the first and second vertices of e j ; 7: for triangle T k in m aTrgls do  add newSf into m aSurfs; 18: end for Generally, there are more than one sub-surfaces can be formed in a original surface, and we classify them into public and private.Definition 4 Private sub-surface: a sub-surface is private when there is only one closed loops in it.Private means it is owned by a loop privately.
Definition 5 Public sub-surface: a sub-surface is public when there is more than one closed loops in it."Public" means several closed loops share this sub-surface together.Noticeably, if this public sub-surface is generated from an original open surface, then it also must have a closed boundary loop besides the intersection closed loop(s).
Definition 6 Sub-surface owner: the only closed loop which owns the private sub-surface or several closed loops which share the public sub-surface.
Remark 1 There is no more than one public sub-surface.
Remark 2 There is at least one private sub-surface.

Assembling and distinguishing sub-blocks 4.6.1 Assembling all possible sub-blocks
Given two original surface S A and S B , a set of sub-surfaces will be obtained after intersecting.For a sub-surface SS from the original surface S A , assuming it has n closed loops, noted as L i (i = 0, . . ., n − 1).It also means that there are n owners share the sub-surface SS.If we can find n private surfaces from the original surface S B , and each private sub-surface is only owned by the closed loop L i , then a sub-block can be easily assembled by the sub-surface SS from the surface S A and n private surfaces from S B .Details of assembling sub-blocks is listed in Algorithm 2. And we can draw some conclusions for assembling sub-blocks.
Remark 3 A private sub-surface can be used once or twice.
Remark 4 A public sub-surface without boundary loop can be used once or twice.
Remark 5 A public sub-surface with a boundary loop, which is generated from an open original surface, cannot be adopted to assemble sub-blocks.For that it is unable to find a subsurface owned also by the boundary closed loop in its intersected original surface.

Distinguishing
According to Algorithm 2, we can obtain all possibly produced sub-blocks.Now, there is a question that given two blocks B A and B B (Fig. 4(a)), and after creating all sub-blocks, how we distinguish the union, intersection and subtraction.We develop a novel and simple method to solve this problem.
Algorithm 2 Create Sub-blocks From Sub-surfaces Input: A set of sub-surfaces stored in m aSurfs Output: A set of sub-blocks stored in m aBlocks 1: while all sub-surfaces in m aSurfs are not tested completely do 2: find a untested sub-surface startSS as the starting one; 3: set startSS as being tested; 4: if startSS is a public sub-surface with boundary loop(s) then create a new empty sub-block newBlk; 8: add startSS into newBlk; 9: for the i th closed loop L i in startSS do 10: for each sub-surface SS i in m aSurfs do if (SS i is untested and also owned only by L i ) and (SS i and startSS not come from same surface) then 12: set SS i as being tested; add newBlk into m aBlocks; 18: end while Step 1: Obtain the non-subtraction(union and intersection) according to the orientation of closed intersection loops when assemble sub-blocks In this step, we can only determine whether the sub-blocks are subtraction or not, and cannot distinguish they are the union or intersection clearly.
For a sub-block SB, assuming it totally has n closed loops, noted as L i (i = 0, . . ., n−1).For each closed loop L i , it has two opposite versions in which L + i represents this loop with its original orientation and L − i is the loop with same vertices as L + i but opposite orientation (Fig. 5).These two loops are defined as twins in Sect. 3  Step 2: Distinguish the union firstly and then the intersections As mentioned above, we let (B A ∪B B ), (B A ∩B B ), (B A −B B ) and (B B −B A ) represent the union, intersection, and subtractions of a pair of blocks B A and B B .Obviously, we can receive the following relationship: (B A ∩B B ) ≤ (B A ∪B B ).The equivalent case holds only when B A and B B coincides, which has to be recognized and handled previously and thus is unnecessary to be considered at this time.Hence, we can in further has (B A ∩B B ) < (B A ∪B B ).This means the intersection is always smaller than its corresponding union.This also means the minimum coordinates of all vertices of the intersection never smaller than those of the union while the maximum ones of the intersection never larger than those of the union.Therefore, the minimum and maximum coordinates of the vertices from both the union and intersection is exactly equal those only from the union.Thus, we can distinguish the union and the intersection via following three sub-procedure.
Step 2-1: Obtain the maximum and minimum coordinates of both intersection and union; this can be done very easily by sorting all the vertices.In fact, this has been realized when merge and update the vertices.We can store their indexes there and call them here.
Step 2-2: Check each candidate union sub-block by comparing its maximum and minimum coordinates.If all of them are equivalent to their corresponding values obtained in above Step 2.1, then it must be the only union volume (Fig. 4(c)).
The above procedures are invalid only in the case that mesh B A is in B B or B B is in B A .Hence, these special cases must be handled in pre-process.
Step 3: Determine all subtractions After classifying the union and intersection, the subtractions (B A −B B ) and (B B −B A ) can be easily determined: define all sub-surfaces which compose the only union volume as outer ones, while those from the intersection volume(s) as inner, for each undistinguished sub-block, (1) If one of the sub-surfaces from the sub-block is outer, and if the sub-surface belongs the original block B A , then the sub-block is part of (B A −B B ); else the sub-block is part of (B B −B A ).
(2) If one of the sub-surfaces from the sub-block is inner, and if the sub-surface belongs the original block B A , then the sub-block is part of (B B −B A ); else the sub-block is part of (B A −B B ).
After distinguishing, all triangles in inner sub-surfaces of subtractions have to reverse to obtain a valid topology of whole mesh model and make all the normal of facets outwards.

Tests and results
As described in Sect.1, we aim at performing Boolean operations only over a pair of surfaces.Either of the surfaces can be open or closed.In this section, we obtain the Boolean operations of several pairs of surfaces including open-and-open (Fig. 6), open-and-closed (Fig. 7 and Fig. 8) and closed-and-closed (Fig. 9 and Fig. 10) to demonstrate our approach.
The intersection of a pair of open triangulated surfaces is tested in Fig. 6.Before considering the boundary outer loops of the original surfaces, four open intersection loops can be formed for each surface; and then the boundary loop is divided into five closed loops; hence five corresponding sub-surfaces can be created for each surface.The division of a planar meshed surface with a triangulated Chinese lion is presented to test the Booleans for a open surface with a closed surface.The model Chinese lion (Fig. 7) is obtained from http://shapes.aim-at-shape.net/.From Fig. 7, we can see that three closed loops can be formed in both the Chinese lion and the surface.Figure.9displays the dividing results that four sub-surfaces can be created for both of original surfaces, and the lion is divided into four sub-blocks.
To test Booleans on two closed surfaces, we give the examples of Boolean operations for a pair of cylinders (Fig. 9) and a pair of torus (Fig. 10).After computing the intersection edges of triangles for cylinders, four soft closed (defined in Sect.4.4) intersection loops can be formed.Both the first vertex and the last vertex of each soft closed loop is shared by four intersection edges.Based on these soft closed loops, four sub-surfaces can be created for each cylinder, and then sub-blocks such as union, intersection and subtractions can be easily assembled and distinguished according to the sub-faces.We present a simple and robust approach to perform Boolean operations on a pair of manifold triangulated surfaces.Our approach falls into the type of exact arithmetic methods in the aspect of computation.We use Octree to locate and find potentially intersected triangles and adopt parallel algorithm to calculate the intersection lines of triangles.We form and define the intersection loops into open, hard closed and soft closed ones, and then create sub-surfaces by growing while only the closed intersection loops are set as boundary loops of sub-surfaces and deemed as the initial advancing fronts.Sub-blocks are assembled very easily according to boundary loops of sub-surfaces, and then distinguished by comparing min and max coordinates of updated vertices.
There are two stages in our approach: the first, which is to calculate intersection lines of triangles and then merge and update the resulting surfaces, needs coordinates calculations; the second, which includes forming intersection loops, creating sub-surfaces, assembling and distinguishing sub-blocks, is only based on the cleared and updated topology of triangular meshes while without coordinate computations of geometric entities.Effectiveness of our approach has been illustrated by several test examples.

Figure 1 :
Figure 1: Flow of the approach

Figure 2 :
Figure 2: Intersected triangle and the re-triangulation.(a) A triangle and its edge-loops; (b) Original triangle divided into 2 polygons; (c) Original triangle divided completely; (d) Triangular partition of the sub polygons.

1 T 2 Figure 3 :
Figure 3: Clear the topology of triangular meshes.(a) Original meshes with invalid triangles T 0 , T 1 and T 2 ; (b) Cleared meshes without same edges Input: A set of closed loops stored in m aLoops and a triangular surface S Output: A set of sub-surfaces stored in m aSurfs 1: for loop L i in m aLoops do 2: create a new empty surface newSf ; 3: copy and add L i as a boundary loop into newSf ; 4:

Figure 4 :Figure 5 :
Figure 4: Boolean operations of a pair of blocks

Figure 6 :
Figure 6: Intersection of a pair of open surfaces Figure.10 shows the Boolean operations for a pair of torus which have the same inner radius.It is similar to that of cylinders but a little more complex.The above two examples only produced soft closed loops, a simple example in which hard closed loops are obtained is presented in Fig.4.

Figure 9 :
Figure 9: Boolean operations of a pair of cylinders

Figure 10 :
Figure 10: Boolean operations of a pair of torus . The loop L i owns or shares two sub-surfaces denoted as SS A L i and SS B L i which come from two different original surfaces S A and S B , respectively.we can easily determine the sub-block SB according to following conditions: