adding slicing function for volume and conversion hexahedrongrid to tetrahedron mesh

Apps/Examples/simpleGMap3.cpp

@@ -34,6 +34,8 @@ SimpleGMap3::SimpleGMap3()
 	normal = myMap.addAttribute<PFP::VEC3>(VERTEX, "normal");
 	volume = myMap.addAttribute<PFP::VEC3>(VOLUME, "volume");
 
+	CellMarker mE(myMap,EDGE);
+
 	Algo::Modelisation::Primitive3D<PFP> primCat(myMap,position);
 	Dart d = primCat.hexaGrid_topo(3,1,1);
 	primCat.embedHexaGrid(2,1,1);
@@ -56,6 +58,9 @@ SimpleGMap3::SimpleGMap3()
 	Geom::Plane3D<PFP::REAL> pl(VEC3(-1,-0.5,-0.5),VEC3(-1,-0.5,0.5),VEC3(1,0.5,0.5));
 	Algo::Modelisation::sliceConvexVolume<PFP>(myMap, position, d, pl);
 
+//	for(Dart d = myMap.begin() ; d != myMap.end() ; myMap.next(d))
+//		mE.mark(d);
+
 	for(unsigned int i = position.begin() ; i != position.end() ; position.next(i))
 		position[i] += VEC3(2,0,0);
 
@@ -77,6 +82,12 @@ SimpleGMap3::SimpleGMap3()
 	myMap.cutEdge(d);
 	position[myMap.phi1(d)] = mid;
 
+	myMap.check();
+//
+	myMap.splitFace(d,myMap.phi1(myMap.phi1(myMap.phi1(d))));
+
+	myMap.check();
+
 	for(unsigned int i = position.begin() ; i != position.end() ; position.next(i))
 		position[i] += VEC3(0,2,0);
 

Apps/Examples/volumeExplorer.h

@@ -27,7 +27,7 @@
 
 #include <iostream>
 
-#define WITH_GMAP 1
+//#define WITH_GMAP 1
 
 #include "Topology/generic/parameters.h"
 #ifdef WITH_GMAP

include/Algo/Modelisation/polyhedron.h

@@ -115,7 +115,7 @@ template <typename PFP>
 Dart createTriangularPrism(typename PFP::MAP& map);
 
 /**
- * create a 3-sided pyramid
+ * create a 4-sided pyramid
  */
 template <typename PFP>
 Dart createQuadrangularPyramid(typename PFP::MAP& map);

include/Algo/Modelisation/subdivision3.h

@@ -49,11 +49,25 @@ Dart cut3Ear(typename PFP::MAP& map, Dart d);
 /**
 * Cut a volume considering a plane
 * @param d dart of the volume
-* @return a dart from the face in the midle
+* @return a dart from the created face
+* * TODO (optimization) change to build path while splitting faces
 */
 template <typename PFP>
 Dart sliceConvexVolume(typename PFP::MAP& map, typename PFP::TVEC3& position, Dart d, Geom::Plane3D<typename PFP::REAL > pl);
 
+/**
+* Cut a volume considering a set of marked edges and vertices
+* marked edges and vertices must form a simple path
+* @param d dart of the volume
+* @param edgesToCut marker to identify edges along the slice
+* @param verticesToSplit marker to identify edges on the slice
+* @return a dart from the created face
+* TODO (optimization) change to build path while splitting faces
+*/
+template <typename PFP>
+Dart sliceConvexVolume(typename PFP::MAP& map, typename PFP::TVEC3& position, Dart d, CellMarker& edgesToCut, CellMarker& verticesToSplit);
+
+
 
 /**
 * catmull clark volumic : do not move the original vertices

include/Algo/Modelisation/subdivision3.hpp

@@ -60,6 +60,8 @@ Dart cut3Ear(typename PFP::MAP& map, Dart d)
 	}
 	else
 	{
+		std::vector<Dart> vPath;
+
 		//triangulate around the vertex
 		do
 		{
@@ -67,17 +69,14 @@ Dart cut3Ear(typename PFP::MAP& map, Dart d)
 			if(map.template phi<111>(e) != e)
 				map.splitFace(map.phi_1(e), map.phi1(e));
 
-			dRing = map.phi1(e);
-			dRing2 = map.phi2(dRing);
+			dRing = map.phi2(map.phi1(e));
 
-			map.unsewFaces(dRing);
+			vPath.push_back(dRing); //remember all darts from the ring
 
 			e = dN;
 		} while (e != d);
 
-		map.closeHole(dRing);
-		map.closeHole(dRing2);
-		map.sewVolumes(map.phi2(dRing), map.phi2(dRing2));
+		map.splitVolume(vPath);
 	}
 
 	return map.phi2(dRing);
@@ -156,13 +155,15 @@ Dart sliceConvexVolume(typename PFP::MAP& map, typename PFP::TVEC3& position, Da
 				Dart dS = dd;
 				bool split=false;
 
-				do {
+				do
+				{
 					//find the new vertex
 					if(vs.isMarked(dS))
 					{
 						Dart dSS = map.phi1(dS);
 						//search an other new vertex (or an existing vertex intersected with the plane) in order to split the face
-						do {
+						do
+						{
 							if(vs.isMarked(dSS))
 							{
 								nbSplit++;
@@ -205,6 +206,223 @@ Dart sliceConvexVolume(typename PFP::MAP& map, typename PFP::TVEC3& position, Da
 	return dRes;
 }
 
+template <typename PFP>
+Dart sliceConvexVolume(typename PFP::MAP& map, typename PFP::TVEC3& position, Dart d, CellMarker& edgesToCut, CellMarker& verticesToSplit)
+{
+	typedef typename PFP::VEC3 VEC3;
+
+	Dart dRes;
+	unsigned int nbInter = 0;
+	unsigned int nbVertices = 0;
+	CellMarkerStore vs(map, VERTEX);			//marker for new vertices from edge cut
+	CellMarkerStore cf(map, FACE);
+	Dart dPath;
+
+	MarkerForTraversor<typename PFP::MAP::ParentMap > mte(map,EDGE);
+	MarkerForTraversor<typename PFP::MAP::ParentMap > mtf(map,FACE);
+
+	//search edges and vertices crossing the plane
+	Traversor3WE<typename PFP::MAP::ParentMap > te(map,d);
+	for(Dart dd = te.begin() ;dd != te.end() ; dd = te.next())
+	{
+		if(!mte.isMarked(dd) && edgesToCut.isMarked(dd))
+		{
+			nbInter++;
+			VEC3 p = (position[dd]+position[map.phi1(dd)])*0.5f;
+			cf.mark(dd);			//mark face and opposite face to split
+			cf.mark(map.phi2(dd));
+
+			map.cutEdge(dd);
+			Dart dN = map.phi1(dd);
+
+			mte.mark(dN);
+
+			vs.mark(dN);		//mark vertex for split
+			position[dN] = p;
+		}
+	}
+
+//	std::cout << "edges cut: " << nbInter << std::endl;
+	unsigned int nbSplit=0;
+
+	//at least two edges are concerned
+	assert(nbInter>1);
+
+	Traversor3WF<typename PFP::MAP::ParentMap > tf(map,d);
+	for(Dart dd = tf.begin() ; dd != tf.end() ; dd = tf.next())
+	{
+		//for faces with a new vertex
+		if(cf.isMarked(dd))
+		{
+			cf.unmark(dd);
+
+			Dart dS = dd;
+			bool split=false;
+
+			do {
+				//find the new vertex
+				if(vs.isMarked(dS) || verticesToSplit.isMarked(dS))
+				{
+					Dart dSS = map.phi1(dS);
+					//search an other new vertex (or an existing vertex intersected with the plane) in order to split the face
+					do {
+						if(vs.isMarked(dSS) || verticesToSplit.isMarked(dSS))
+						{
+							nbSplit++;
+							map.splitFace(dS,dSS);
+							dPath=map.phi_1(dS);
+							split=true;
+						}
+						dSS = map.phi1(dSS);
+					} while(!split && dSS!=dS);
+				}
+				dS = map.phi1(dS);
+			} while(!split && dS!=dd);
+		}
+
+		//define the path to split
+		std::vector<Dart> vPath;
+		vPath.reserve((nbSplit+nbVertices)+1);
+		vPath.push_back(dPath);
+		for(std::vector<Dart>::iterator it = vPath.begin() ;it != vPath.end() ; ++it)
+		{
+			Dart dd = map.phi1(*it);
+
+			Dart ddd = map.phi1(map.phi2(dd));
+
+			while(!vs.isMarked(map.phi1(ddd)) && ddd!=dd)
+				ddd = map.phi1(map.phi2(ddd));
+
+			if(vs.isMarked(map.phi1(ddd)) && !map.sameVertex(ddd,*vPath.begin()))
+				vPath.push_back(ddd);
+		}
+
+		assert(vPath.size()>2);
+		map.splitVolume(vPath);
+		dRes = map.phi2(*vPath.begin());
+	}
+
+	return dRes;
+}
+
+template <typename PFP>
+Dart sliceConvexVolumes(typename PFP::MAP& map, typename PFP::TVEC3& position,CellMarker& volumesToCut, CellMarker& edgesToCut, CellMarker& verticesToSplit)
+{
+	Dart dRes=NIL;
+
+    typedef typename PFP::VEC3 VEC3;
+    CellMarkerStore localVerticesToSplit(map, VERTEX); //marker for new vertices from edge cut
+
+    //Step 1: Cut the edges and mark the resulting vertices as vertices to be face-split
+    TraversorE<typename PFP::MAP> te(map);
+    CellMarkerStore cf(map, FACE);
+
+    for(Dart d = te.begin(); d != te.end(); d=te.next()) //cut all edges
+    {
+        if(edgesToCut.isMarked(d))
+        {
+            VEC3 p = (position[d]+position[map.phi1(d)])*0.5f;
+            cf.mark(d);			//mark face and opposite face to split
+            cf.mark(map.phi2(d));
+
+            map.cutEdge(d);
+            Dart dN = map.phi1(d);
+
+            localVerticesToSplit.mark(dN);		//mark vertex for split
+            position[dN] = p;
+        }
+    }
+
+    //Step 2: Split faces with cut edges
+    TraversorF<typename PFP::MAP> tf(map);
+    for(Dart d = tf.begin(); d != tf.end(); d=tf.next())
+    {
+        if(cf.isMarked(d))
+        {
+            cf.unmark(d);
+            Dart dS = d;
+            bool split=false;
+            do
+            {
+                //find the new vertex
+                if(localVerticesToSplit.isMarked(dS))
+                {
+                    Dart dSS = map.phi1(dS);
+                    //search an other new vertex (or an existing vertex to split) in order to split the face
+                    do
+                    {
+                        if(localVerticesToSplit.isMarked(dSS) || verticesToSplit.isMarked(dSS))
+                        {
+                            map.splitFace(dS,dSS);
+                            split=true;
+                        }
+                        dSS = map.phi1(dSS);
+                    } while(!split && dSS!=dS);
+                }
+                dS = map.phi1(dS);
+            } while(!split && dS!=d);
+        }
+    }
+
+    //Step 3 : Find path and split volumes
+    TraversorW<typename PFP::MAP> tw(map);
+    for(Dart d = tw.begin(); d != tw.end(); d=tw.next()) //Parcours des volumes
+    {
+        if(volumesToCut.isMarked(d)){
+
+            Traversor3WV<typename PFP::MAP> t3wv(map,d);
+            int nbVInPath=1;
+            Dart dPath;
+            bool found=false;
+
+            //find a vertex of the volume to start the path to split
+            for(Dart dd = t3wv.begin(); dd != t3wv.end() && !found; dd=t3wv.next())
+            {
+                if(localVerticesToSplit.isMarked(dd) || verticesToSplit.isMarked(dd))
+                {
+                    Dart ddd = map.phi1(map.phi2(dd));
+                    while(!localVerticesToSplit.isMarked(map.phi1(ddd))
+                    		&& !verticesToSplit.isMarked(map.phi1(ddd))
+                    		&& ddd!=dd)
+                        ddd = map.phi1(map.phi2(ddd));
+                    found=true;
+                    nbVInPath++;
+                    dPath=ddd;
+                }
+            }
+            //define the path to split
+            std::vector<Dart> vPath;
+            vPath.reserve(nbVInPath);
+            vPath.push_back(dPath);
+            CellMarker cmf(map,FACE);
+            int nbFacesVisited=0;
+
+            //define the path to split for the whole volume
+            for(std::vector<Dart>::iterator it = vPath.begin() ;it != vPath.end(); ++it)
+            {
+                Dart dd = map.phi1(*it);
+                Dart ddd = map.phi1(map.phi2(dd));
+                while(!localVerticesToSplit.isMarked(map.phi1(ddd))
+                		&& !verticesToSplit.isMarked(map.phi1(ddd))
+                		&& ddd!=dd)
+                    ddd = map.phi1(map.phi2(ddd));
+
+                if( !map.sameVertex(ddd,*vPath.begin()) )
+                {
+                    vPath.push_back(map.phi1(ddd));
+                    localVerticesToSplit.unmark(map.phi1(ddd));
+                }
+            }
+
+            assert(vPath.size()>2);
+            map.splitVolume(vPath);
+            dRes = map.phi2(*vPath.begin());
+        }
+    }
+
+    return dRes;
+}
+
 template <typename PFP, typename EMBV, typename EMB>
 void catmullClarkVol(typename PFP::MAP& map, EMBV& attributs, const FunctorSelect& selected)
 {

include/Algo/Modelisation/tetrahedralization.h

@@ -39,11 +39,17 @@ namespace Tetrahedralization
 
 /**
 * subdivide a hexahedron into 5 tetrahedron
-* @param d dart of the hexahedron
 */
 template <typename PFP>
 void hexahedronToTetrahedron(typename PFP::MAP& map, Dart d);
 
+/**
+* WARNING : assume all volumes to be hexahedrons
+* subdivide a hexahedron mesh into a tetrahedron mesh
+*/
+template <typename PFP>
+void hexahedronsToTetrahedrons(typename PFP::MAP& map);
+
 /************************************************************************************************
  * 																Tetrahedron functions																	   *
  ************************************************************************************************/

include/Algo/Modelisation/tetrahedralization.hpp

@@ -22,6 +22,8 @@
  *                                                                              *
  *******************************************************************************/
 
+#include "Algo/Modelisation/subdivision3.h"
+
 namespace CGoGN
 {
 
@@ -37,18 +39,54 @@ namespace Tetrahedralization
 template <typename PFP>
 void hexahedronToTetrahedron(typename PFP::MAP& map, Dart d)
 {
-	//Splitting Path
-	std::vector<Dart> sp;
-	sp.reserve(32);
+	Dart d1 = d;
+	Dart d2 = map.phi1(map.phi1(d));
+	Dart d3 = map.phi_1(map.phi2(d));
+	Dart d4 = map.phi1(map.phi1(map.phi2(map.phi_1(d3))));
+
+	cut3Ear<PFP>(map,d1);
+	cut3Ear<PFP>(map,d2);
+	cut3Ear<PFP>(map,d3);
+	cut3Ear<PFP>(map,d4);
+}
 
-	Traversor3VE<typename PFP::MAP> tra(map, d);
-	for (Dart d = tra.begin() ; d != tra.end() ; d = tra.next())
+template <typename PFP>
+void hexahedronsToTetrahedrons(typename PFP::MAP& map)
+{
+	TraversorV<typename PFP::MAP> tv(map);
+	for(Dart s=tv.begin() ; s!=tv.end() ; s=tv.next() )
 	{
-		map.splitFace(map.phi1(d), map.phi_1(d));
-		sp.push_back(map.phi1(d));
-	}
+		bool allDegreesAroundEqual3=true;
+
+		std::vector<Dart> dov;
+
+		FunctorStore fs(dov);
+		map.foreach_dart_of_vertex(s,fs);
+
+		CellMarkerStore cmv(map,VOLUME);
+
+		for(std::vector<Dart>::iterator it=dov.begin();allDegreesAroundEqual3 && it!=dov.end();++it)
+		{
+			if(!cmv.isMarked(*it) && !map.isBoundaryMarked(*it))
+			{
+				cmv.mark(*it);
+				allDegreesAroundEqual3 = (map.phi1(map.phi2(map.phi1(map.phi2(map.phi1(map.phi2(*it))))))==*it);
+			}
+		}
 
-	map.splitVolume(sp);
+		if( allDegreesAroundEqual3)
+		{
+			CellMarkerStore cmv(map,VOLUME);
+			for(std::vector<Dart>::iterator it=dov.begin();it!=dov.end();++it)
+			{
+				if(!cmv.isMarked(*it) && !map.isBoundaryMarked(*it))
+				{
+					cmv.mark(*it);
+					cut3Ear<PFP>(map,*it);
+				}
+			}
+		}
+	}
 }