From 03f47a06a3ca1cba6c4285973993f189a93284b7 Mon Sep 17 00:00:00 2001 From: untereiner Date: Mon, 12 Mar 2012 10:05:00 +0100 Subject: [PATCH] Adding 3maps MR with Regular Subdivision --- Apps/Examples/simpleGMap3.cpp | 2 +- Apps/Tuto/tuto5.cpp | 6 +- .../ImplicitHierarchicalMesh/subdivision3.hpp | 46 +- include/Algo/Import/import.h | 32 +- include/Algo/Import/import2tables.h | 5 +- include/Algo/Import/import2tablesVolume.hpp | 172 +--- include/Algo/Import/importMesh.hpp | 317 +------- include/Algo/Import/importTet.hpp | 4 +- include/Algo/Import/importTs.hpp | 114 +-- include/Algo/Modelisation/subdivision3.hpp | 463 ++++++++--- include/Algo/Render/GL2/topo3Render.h | 3 + include/Algo/Render/GL2/topo3Render.hpp | 224 ++++++ include/Topology/generic/genericmap.hpp | 1 - include/Topology/gmap/gmap2.h | 2 +- include/Topology/map/embeddedMap2.h | 5 + include/Topology/map/map3.h | 5 + include/Topology/map/map3MR/filters_Primal.h | 752 ++++++++++++++++++ .../Topology/map/map3MR/map3MR_PrimalAdapt.h | 196 +++++ .../map/map3MR/map3MR_PrimalRegular.h | 159 ++++ src/Topology/gmap/gmap2.cpp | 17 +- src/Topology/gmap/gmap3.cpp | 22 +- src/Topology/map/embeddedMap2.cpp | 47 +- src/Topology/map/embeddedMap3.cpp | 15 +- src/Topology/map/map2.cpp | 19 +- src/Topology/map/map3.cpp | 16 +- .../map/map3MR/map3MR_PrimalAdapt.cpp | 665 ++++++++++++++++ .../map/map3MR/map3MR_PrimalRegular.cpp | 671 ++++++++++++++++ 27 files changed, 3284 insertions(+), 696 deletions(-) create mode 100644 include/Topology/map/map3MR/filters_Primal.h create mode 100644 include/Topology/map/map3MR/map3MR_PrimalAdapt.h create mode 100644 include/Topology/map/map3MR/map3MR_PrimalRegular.h create mode 100644 src/Topology/map/map3MR/map3MR_PrimalAdapt.cpp create mode 100644 src/Topology/map/map3MR/map3MR_PrimalRegular.cpp diff --git a/Apps/Examples/simpleGMap3.cpp b/Apps/Examples/simpleGMap3.cpp index d992c966..f354bc13 100644 --- a/Apps/Examples/simpleGMap3.cpp +++ b/Apps/Examples/simpleGMap3.cpp @@ -54,7 +54,7 @@ SimpleGMap3::SimpleGMap3() SelectorMarked sm(markOrient); std::cout << "AAA"<< std::endl; - Algo::Modelisation::catmullClarkVol(myMap, position, sm); + //Algo::Modelisation::catmullClarkVol(myMap, position, sm); for(unsigned int i = position.begin() ; i != position.end() ; position.next(i)) position[i] += VEC3(2,0,0); diff --git a/Apps/Tuto/tuto5.cpp b/Apps/Tuto/tuto5.cpp index 260e9345..c6714269 100644 --- a/Apps/Tuto/tuto5.cpp +++ b/Apps/Tuto/tuto5.cpp @@ -263,8 +263,8 @@ void MyQT::cb_keyPress(int code) if(code == 'c') { - SelectorDartNoBoundary nb(myMap); - Algo::Modelisation::catmullClarkVol(myMap, position, nb); + //SelectorDartNoBoundary nb(myMap); + Algo::Modelisation::catmullClarkVol(myMap, position); m_positionVBO->updateData(position); m_render->initPrimitives(myMap, allDarts, Algo::Render::GL2::TRIANGLES); @@ -272,7 +272,7 @@ void MyQT::cb_keyPress(int code) m_render->initPrimitives(myMap, allDarts, Algo::Render::GL2::POINTS); - m_render_topo->updateData(myMap, position, 0.9f, 0.9f, 0.9f, nb); + m_render_topo->updateData(myMap, position, 0.9f, 0.9f, 0.9f, allDarts); } } diff --git a/include/Algo/ImplicitHierarchicalMesh/subdivision3.hpp b/include/Algo/ImplicitHierarchicalMesh/subdivision3.hpp index 4405f1c7..f05494d8 100644 --- a/include/Algo/ImplicitHierarchicalMesh/subdivision3.hpp +++ b/include/Algo/ImplicitHierarchicalMesh/subdivision3.hpp @@ -466,31 +466,31 @@ Dart subdivideVolumeClassic2(typename PFP::MAP& map, Dart d, typename PFP::TVEC3 - //Third step : 3-sew internal faces - for (std::vector >::iterator it = subdividedfaces.begin(); it != subdividedfaces.end(); ++it) - { - Dart f1 = (*it).first; - Dart f2 = (*it).second; - - - if(map.isBoundaryFace(map.phi2(f1)) && map.isBoundaryFace(map.phi2(f2))) - { - std::cout << "boundary" << std::endl; - //id pour toutes les faces interieures - map.sewVolumes(map.phi2(f1), map.phi2(f2)); - +// //Third step : 3-sew internal faces +// for (std::vector >::iterator it = subdividedfaces.begin(); it != subdividedfaces.end(); ++it) +// { +// Dart f1 = (*it).first; +// Dart f2 = (*it).second; // -// //Fais a la couture !!!!! -// unsigned int idface = map.getNewFaceId(); -// map.setFaceId(map.phi2(f1),idface, FACE); - } // -// //FAIS a la couture !!!!!!! -// //id pour toutes les aretes exterieurs des faces quadrangulees -// unsigned int idedge = map.getEdgeId(f1); -// map.setEdgeId(map.phi2(f1), idedge, DART); -// map.setEdgeId( map.phi2(f2), idedge, DART); - } +// if(map.isBoundaryFace(map.phi2(f1)) && map.isBoundaryFace(map.phi2(f2))) +// { +// std::cout << "boundary" << std::endl; +// //id pour toutes les faces interieures +// map.sewVolumes(map.phi2(f1), map.phi2(f2)); +// +//// +//// //Fais a la couture !!!!! +//// unsigned int idface = map.getNewFaceId(); +//// map.setFaceId(map.phi2(f1),idface, FACE); +// } +//// +//// //FAIS a la couture !!!!!!! +//// //id pour toutes les aretes exterieurs des faces quadrangulees +//// unsigned int idedge = map.getEdgeId(f1); +//// map.setEdgeId(map.phi2(f1), idedge, DART); +//// map.setEdgeId( map.phi2(f2), idedge, DART); +// } //LA copie de L'id est a gerer avec le sewVolumes normalement !!!!!! //id pour les aretes interieurs : (i.e. 6 pour un hexa) diff --git a/include/Algo/Import/import.h b/include/Algo/Import/import.h index 2924a97b..6e423466 100644 --- a/include/Algo/Import/import.h +++ b/include/Algo/Import/import.h @@ -44,10 +44,9 @@ namespace Import * import a mesh * @param map the map in which the function imports the mesh * @param filename (*.{trian,trianbgz,off,obj,ply}) -* @param positions table of vertices positions attribute -* @param m a marker that will be set by the function. If closeObject=false the phi2 that have fixed point are marked, else the created darts of the boundary are marked. +* @param attrNames attribute names * @param kind what kind of mesh is the file (if none (-1) determined by filename extension) (cf enum in Mesh2Tables for other kind values) -* @param closeObject a boolean indicating if the imported mesh should be closed +* @param mergeColseVertces a boolean indicating if the imported mesh should be closed * @return a boolean indicating if import was successfull */ template @@ -55,14 +54,25 @@ bool importMesh(typename PFP::MAP& map, const std::string& filename, std::vector /** * import a volumetric mesh + * @param map the map in which the function imports the mesh + * @param filename (*.{tet,off,ts}) + * @param attrNames attribute names + * @param kind what kind of mesh is the file (if none (-1) determined by filename extension) (cf enum in Mesh2Tables for other kind values) + * @param mergeColseVertces a boolean indicating if the imported mesh should be closed + * @return a boolean indicating if import was successfull */ -//template -//bool importMesh(typename PFP::MAP& map, const std::string& filename, typename PFP::TVEC3& positions, ImportVolumique::ImportType kind = ImportVolumique::UNKNOWNVOLUME); - -//template -//bool importObjWithTex(typename PFP::MAP& map, const std::string& filename); -// +template +bool importMeshV(typename PFP::MAP& map, const std::string& filename, std::vector& attrNames, ImportVolumique::ImportType kind = ImportVolumique::UNKNOWNVOLUME, bool mergeCloseVertices=false); +/** + * import a mesh and extrud it + * @param map the map in which the function imports the mesh + * @param filename (*.{trian,trianbgz,off,obj,ply}) + * @param attrNames attribute names + * @param kind what kind of mesh is the file (if none (-1) determined by filename extension) (cf enum in Mesh2Tables for other kind values) + * @param mergeColseVertces a boolean indicating if the imported mesh should be closed + * @return a boolean indicating if import was successfull + */ template bool importMeshToExtrude(typename PFP::MAP& map, const std::string& filename, std::vector& attrNames, ImportSurfacique::ImportType kind = ImportSurfacique::UNKNOWNSURFACE); @@ -80,6 +90,10 @@ template bool importTs(typename PFP::MAP& the_map, const std::string& filename, std::vector& attrNames, float scaleFactor = 1.0f); +//template +//bool importObjWithTex(typename PFP::MAP& map, const std::string& filename); +// + } // namespace Import } // namespace Algo diff --git a/include/Algo/Import/import2tables.h b/include/Algo/Import/import2tables.h index c1b66e18..20a18ec7 100644 --- a/include/Algo/Import/import2tables.h +++ b/include/Algo/Import/import2tables.h @@ -57,7 +57,7 @@ namespace Import namespace ImportVolumique { - enum ImportType { UNKNOWNVOLUME ,TET ,TRIANBGZ ,PLY ,OFF, OBJ }; + enum ImportType { UNKNOWNVOLUME , TET, ELE, TS }; } @@ -184,9 +184,6 @@ public: bool importTet(const std::string& filename, std::vector& attrNames, float scaleFactor); - bool importPly(const std::string& filename, std::vector& attrNames); - - bool importTrianBinGz(const std::string& filename, std::vector& attrNames); MeshTablesVolume(typename PFP::MAP& map): m_map(map) diff --git a/include/Algo/Import/import2tablesVolume.hpp b/include/Algo/Import/import2tablesVolume.hpp index 0560b7d7..9ffed880 100644 --- a/include/Algo/Import/import2tablesVolume.hpp +++ b/include/Algo/Import/import2tablesVolume.hpp @@ -37,11 +37,11 @@ ImportVolumique::ImportType MeshTablesVolume::getFileType(const std::string if ((filename.rfind(".tet")!=std::string::npos) || (filename.rfind(".TET")!=std::string::npos)) return ImportVolumique::TET; - if ((filename.rfind(".trianbgz")!=std::string::npos) || (filename.rfind(".TRIANBGZ")!=std::string::npos)) - return ImportVolumique::TRIANBGZ; + if ((filename.rfind(".ele")!=std::string::npos) || (filename.rfind(".ELE")!=std::string::npos)) + return ImportVolumique::ELE; - if ((filename.rfind(".ply")!=std::string::npos) || (filename.rfind(".PLY")!=std::string::npos)) - return ImportVolumique::PLY; + if ((filename.rfind(".ts")!=std::string::npos) || (filename.rfind(".TS")!=std::string::npos)) + return ImportVolumique::TS; return ImportVolumique::UNKNOWNVOLUME; } @@ -54,14 +54,12 @@ bool MeshTablesVolume::importMesh(const std::string& filename, std::vector< switch (kind) { - case ImportVolumique::PLY: - return importPly(filename, attrNames); - break; case ImportVolumique::TET: return importTet(filename, attrNames, scaleFactor); break; - case ImportVolumique::TRIANBGZ: - return importTrianBinGz(filename, attrNames); + case ImportVolumique::ELE: + break; + case ImportVolumique::TS: break; default: CGoGNerr << "Not yet supported" << CGoGNendl; @@ -73,48 +71,44 @@ bool MeshTablesVolume::importMesh(const std::string& filename, std::vector< template bool MeshTablesVolume::importTet(const std::string& filename, std::vector& attrNames, float scaleFactor=1.0f) { - AttributeContainer& container = m_map.getAttributeContainer(VERTEX) ; - -// AttributeHandler positions = m_map.template addAttribute(VERTEX, "position") ; - AttributeHandler positions = m_map.template getAttribute(VERTEX, "position") ; + AttributeHandler positions = m_map.template getAttribute(VERTEX, "position") ; if (!positions.isValid()) - positions = m_map.template addAttribute(VERTEX, "position") ; + positions = m_map.template addAttribute(VERTEX, "position") ; attrNames.push_back(positions.name()) ; - // open file + AttributeContainer& container = m_map.getAttributeContainer(VERTEX) ; + + //open file std::ifstream fp(filename.c_str(), std::ios::in); if (!fp.good()) { - CGoGNerr << "Unable to open file " << filename<< CGoGNendl; + CGoGNerr << "Unable to open file " << filename << CGoGNendl; return false; } - std::string ligne; - int nbv,nbt; - // lecture des nombres de sommets/tetra + std::string ligne; + unsigned int nbv, nbt; + // reading number of vertices std::getline (fp, ligne); std::stringstream oss(ligne); oss >> nbv; - //CGoGNout << "nbV = " << nbv << CGoGNendl; - + // reading number of tetrahedra std::getline (fp, ligne); std::stringstream oss2(ligne); oss2 >> nbt; - //CGoGNout << "nbT = " << nbt << CGoGNendl; - - //lecture sommets + //reading vertices std::vector verticesID; verticesID.reserve(nbv); - for(int i=0; i::importTet(const std::string& filename, std::vector> x; oss >> y; oss >> z; - // on peut ajouter ici la lecture de couleur si elle existe + // TODO : if required read other vertices attributes here VEC3 pos(x*scaleFactor,y*scaleFactor,z*scaleFactor); - //CGoGNout << "VEC3 = " << pos << CGoGNendl; - unsigned int id = container.insertLine(); positions[id] = pos; verticesID.push_back(id); } - m_nbVertices = verticesID.size(); - //CGoGNout << "nbVertices = " << m_nbVertices << CGoGNendl; + m_nbVertices = nbv; m_nbVolumes = nbt; - //CGoGNout << "nbVolumes = " << m_nbVolumes << CGoGNendl; + // lecture tetra // normalement m_nbVolumes*12 (car on ne charge que des tetra) @@ -154,12 +145,7 @@ bool MeshTablesVolume::importTet(const std::string& filename, std::vector> n; // nb de faces d'un volume ? -// m_nbVerticesPerFace.push_back(3); -// m_nbVerticesPerFace.push_back(3); -// m_nbVerticesPerFace.push_back(3); -// m_nbVerticesPerFace.push_back(3); -// m_nbFacesPerVolume.push_back(4); m_nbEdges.push_back(3); m_nbEdges.push_back(3); m_nbEdges.push_back(3); @@ -193,112 +179,6 @@ bool MeshTablesVolume::importTet(const std::string& filename, std::vector -bool MeshTablesVolume::importTrianBinGz(const std::string& filename, std::vector& attrNames) -{ -// // open file -// igzstream fs(filename.c_str(), std::ios::in|std::ios::binary); -// -// if (!fs.good()) -// { -// CGoGNerr << "Unable to open file " << filename << CGoGNendl; -// return false; -// } -// // read nb of points -// int nbp; -// fs.read(reinterpret_cast(&nbp), sizeof(int)); -// -// // read points -// std::vector vertices; -// {// juste pour limiter la portee des variables -// vertices.reserve(nbp); -// int labEmb=0; -// float* buffer = new float[nbp*3]; -// fs.read(reinterpret_cast(buffer),3*nbp*sizeof(float)); -// float *ptr = buffer; -// for (int i=0; isetLabel(labEmb++); -// vertices.push_back(em); -// } -// delete[] buffer; -// } -// m_nbVertices = vertices.size(); -// -// // read nb of faces -// fs.read(reinterpret_cast(&m_nbFaces), sizeof(int)); -// m_nbVolumes=1; -// m_nbFacesPerVolume.push_back(m_nbFaces); -// m_nbVerticesPerFace.reserve(m_nbFaces); -// m_emb.reserve(3*m_nbFaces); -// -// // read indices of faces -// {// juste pour limiter la portee des variables -// int* buffer = new int[m_nbFaces*6]; -// fs.read(reinterpret_cast(buffer),6*m_nbFaces*sizeof(float)); -// int *ptr = buffer; -// -// for (unsigned i=0; i -bool MeshTablesVolume::importPly(const std::string& filename, std::vector& attrNames) -{ -// PlyImportData pid; -// -// if (! pid.read_file(filename) ) -// { -// CGoGNerr << "Unable to open file " << filename<< CGoGNendl; -// return false; -// } -// -// // lecture des nombres de sommets/faces/ar�tes -// m_nbVertices = pid.nbVertices(); -// m_nbFaces = pid.nbFaces(); -// m_nbVolumes=1; -// -// -// //lecture sommets -// std::vector vertices; -// vertices.reserve(m_nbVertices); -// for (unsigned i=0; isetLabel(i); -// vertices.push_back(em); -// } -// -// m_nbVerticesPerFace.reserve(m_nbFaces); -// m_emb.reserve(m_nbVertices*8); -// m_nbFacesPerVolume.push_back(m_nbFaces); -// -// for (unsigned i=0; i& mtv) AutoAttributeHandler< NoMathIONameAttribute< std::vector > > vecDartsPerVertex(map, VERTEX, "incidents"); - //ALGO - // précondition - // indexer le tableau du nombre d'arrete par faces par le nombre de faces par volume - // - //creation d'une sous-carte (2-carte) pour les faces d'un volume en fonction du nombre d'arete grace a un meshtableSurface - //avoir un tableau_de_sommets[nbsommets] indexe par les sommets dans lequel on comptera le nombre de fois qu'on a croisé un sommet - // - // - //fusion de toutes les composantes connexes avec couture par phi3 - - - unsigned int nbf = mtv.getNbFaces(); //getNbVolume() - unsigned int indexNbFaces = 0; - - //for each volume - for(unsigned int i = 0 ; i < nbf ; ++i) - { - unsigned int indexNbFacesMin = indexNbFaces; - - indexNbFaces += mtv.getNbFacesVolume(i) - 1; - - //for each face of this volume - for(unsigned int j = indexNbFacesMin ; j<= indexNbFaces ; ++j) - { - //reconstructing the faces for this volume - - unsigned int nbEdges = mtv.getNbEdgesFace(j); - - - -// std::vector edgesBuffer; -// edgesBuffer.reserve(8); -// -// edgesBuffer.clear(); -// unsigned int prec = EMBNULL; -// for (unsigned int j = 0; j < nbe; ++j) -// { -// unsigned int em = mtv.getEmbIdx(index++); -// if (em!=prec) -// { -// prec = em; -// edgesBuffer.push_back(em); -// } -// } -// // check first/last vertices -// if (edgesBuffer.front() == edgesBuffer.back()) -// edgesBuffer.pop_back(); -// -// // create only non degenerated faces -// nbe = edgesBuffer.size(); -// if (nbe >2) -// { -// Dart d = map.newFace(nbe); -// for (unsigned int j = 0; j < nbe; ++j) -// { -// unsigned int em = edgesBuffer[j]; // get embedding -// map.setDartEmbedding(VERTEX, d, em); // associate to dart -// vecDartsPerVertex[em].push_back(d); // store incident darts for fast phi2 reconstruction -// d = map.phi1(d); -// } -// m.markOrbit(FACE, d); // mark on the fly to unmark on second loop -// } - - } - - indexNbFaces++; - } - - - - - - - - - - int index = 0; - // buffer for tempo faces (used to remove degenerated edges) - std::vector edgesBuffer; - edgesBuffer.reserve(8); - - DartMarkerNoUnmark m(map) ; - - // for each face of table - for(unsigned int i = 0; i < nbf; ++i) - { - // store face in buffer, removing degenerated edges - unsigned int nbe = mtv.getNbEdgesFace(i); - edgesBuffer.clear(); - unsigned int prec = EMBNULL; - for (unsigned int j = 0; j < nbe; ++j) - { - unsigned int em = mtv.getEmbIdx(index++); - if (em!=prec) - { - prec = em; - edgesBuffer.push_back(em); - } - } - // check first/last vertices - if (edgesBuffer.front() == edgesBuffer.back()) - edgesBuffer.pop_back(); - - // create only non degenerated faces - nbe = edgesBuffer.size(); - if (nbe >2) - { - Dart d = map.newFace(nbe); - for (unsigned int j = 0; j < nbe; ++j) - { - unsigned int em = edgesBuffer[j]; // get embedding - map.setDartEmbedding(VERTEX, d, em); // associate to dart - vecDartsPerVertex[em].push_back(d); // store incident darts for fast phi2 reconstruction - d = map.phi1(d); - } - m.markOrbit(FACE, d); // mark on the fly to unmark on second loop - } - } - //unsigned int nbnm = 0; - // reconstruct neighbourhood - - for (Dart d = map.begin(); d != map.end(); map.next(d)) - { - if (m.isMarked(d)) - { - std::vector& vec = vecDartsPerVertex[map.phi1(d)]; - - //Phi3 reconstruction - unsigned int emb1d = map.getEmbedding(VERTEX, map.phi_1(d)); - unsigned int emb2d = map.getEmbedding(VERTEX, map.phi1(map.phi1(d))); - - for (typename std::vector::iterator it = vec.begin(); it != vec.end(); ++it) - { - Dart dprim=*it; - unsigned int emb1dprim = map.getEmbedding(VERTEX, map.phi1(map.phi1(dprim))); - unsigned int emb2dprim = map.getEmbedding(VERTEX, map.phi_1(dprim)); - - if(emb1d == emb1dprim && emb2d == emb2dprim) - { - map.sewVolumes(d,dprim); - } - } - - //Phi2 reconstruction - emb1d = map.getEmbedding(VERTEX, d); - emb2d = map.getEmbedding(VERTEX, map.phi1(d)); - - for (typename std::vector::iterator it = vec.begin(); it != vec.end(); ++it) - { - Dart dprim=*it; - unsigned int emb1dprim = map.getEmbedding(VERTEX, map.phi1(dprim)); - unsigned int emb2dprim = map.getEmbedding(VERTEX, dprim); - - if(emb1d == emb1dprim && emb2d == emb2dprim) - { - map.sewFaces(d,dprim); - } - } - - m.unmarkOrbit(DART,d); - } - } - - return true ; - -// -// //creer un polyhedre de nbf faces -// //stocker pour un dart d tout les darts autour de l'arete -// // reconstruire le voisine de l'arete pour phi3 -// -// AutoAttributeHandler< NoMathIONameAttribute< std::vector > > vecDartsPerEdge(map, VERTEX, "incidents"); -// -// unsigned nbv = mtv.getNbVolumes(); -// int index = 0; -// // buffer for tempo faces (used to remove degenerated edges) -// std::vector edgesBuffer; -// edgesBuffer.reserve(8); -// -// DartMarkerNoUnmark m(map) ; -// -// //for each volume of table -// for(unsigned int i = 0; i < nbv ; ++i) -// { -// unsigned int nbf = mtv.getNbFacesPerVolume(i); -// -// //for each face of a volume -// for(unsigned int j = 0; j < nbf; ++j) -// { -// // store face in buffer, removing degenerated edges -// unsigned int nbe = mtv.getNbVerticesPerFace(j); -// edgesBuffer.clear(); -// unsigned int prec = EMBNULL; -// for (unsigned int k = 0; k < nbe; ++k) -// { -// unsigned int em = mtv.getEmbIdx(index++); -// if (em!=prec) -// { -// prec = em; -// edgesBuffer.push_back(em); -// } -// } -// -// // check first/last vertices -// if (edgesBuffer.front() == edgesBuffer.back()) -// edgesBuffer.pop_back(); -// -// // create only non degenerated faces -// nbe = edgesBuffer.size(); -// -// if (nbe >2) -// { -// Dart d = map.newFace(nbe); -// Dart d3 = map.newFace(nbe); -// map.sewVolumes(d,d3); -// -// for (unsigned int l = 0; l < nbe; ++l) -// { -// unsigned int em = edgesBuffer[l]; // get embedding -// map.setDartEmbedding(VERTEX, d, em); // associate to dart -// vecDartsPerVertex[em].push_back(d); // store incident darts for fast phi2 reconstruction -// d = map.phi1(d); -// -//// if(l == nbe-1) -//// { -//// unsigned int em2 = edgesBuffer[0]; // get embedding -//// map.setDartEmbedding(VERTEX, d3, em2); // associate to dart -//// vecDartsPerVertex[em2].push_back(d3); // store incident darts for fast phi2 reconstruction -//// } -//// else -//// { -//// unsigned int em2 = edgesBuffer[l+1]; // get embedding -//// map.setDartEmbedding(VERTEX, d3, em2); // associate to dart -//// vecDartsPerVertex[em2].push_back(d3); // store incident darts for fast phi2 reconstruction -//// } -//// -//// d3 = map.phi_1(d3); -// } -// m.markOrbit(FACE, d); // mark on the fly to unmark on second loop -// //m.markOrbit(FACE, d3); // mark on the fly to unmark on second loop -// } -// } -// } -// -// -// unsigned int nbnm = 0; -// // reconstruct neighbourhood -// for (Dart d = map.begin(); d != map.end(); map.next(d)) -// { -// if (m.isMarked(d)) -// { -// // darts incident to end vertex of edge -// std::vector& vec = vecDartsPerVertex[map.phi1(d)]; -// -// unsigned int embd = map.getEmbedding(VERTEX, d); -// unsigned int nbf = 0; -// Dart good_dart = d; -// -// CGoGNout << "vec size" << vec.size() << CGoGNendl; -// -// for(typename std::vector::iterator it = vec.begin(); it != vec.end(); ++it) -// { -// if (map.getEmbedding(VERTEX, map.phi1(*it)) == embd) -// { -// good_dart = *it; -// nbf++; -// } -// } -// -// if (nbf == 2) -// { -// if (good_dart == map.phi2(good_dart)) -// { -// map.sewFaces(d, good_dart); -// m.unmarkOrbit(EDGE, d); -// } -// } -// else -// { -// ++nbnm; -// } -// } -// } -// -// if (nbnm > 0) -// { -// CGoGNout << "Warning " << nbnm << " darts with phi2 fix points" << CGoGNendl; -// } -// -// return true ; + return false; } template @@ -549,6 +261,33 @@ bool importMesh(typename PFP::MAP& map, const std::string& filename, std::vector return importMesh(map, mts); } +template +bool importMeshV(typename PFP::MAP& map, const std::string& filename, std::vector& attrNames, ImportVolumique::ImportType kind, bool mergeCloseVertices) +{ + switch (kind) + { + case ImportVolumique::TET: + return Algo::Import::importTet(map, filename, attrNames, 1.0f); + break; + case ImportVolumique::ELE: + { + size_t pos = filename.rfind("."); + std::string fileOFF = filename; + fileOFF.erase(pos); + fileOFF.append(".off"); + return Algo::Import::importOFFWithELERegions(map, fileOFF, filename, attrNames); + break; + } + case ImportVolumique::TS: + Algo::Import::importTs(map, filename, attrNames, 1.0f); + break; + default: + CGoGNerr << "Not yet supported" << CGoGNendl; + break; + } + return false; +} + template bool importMeshToExtrude(typename PFP::MAP& map, const std::string& filename, std::vector& attrNames, ImportSurfacique::ImportType kind) diff --git a/include/Algo/Import/importTet.hpp b/include/Algo/Import/importTet.hpp index a7ce0dda..327d183e 100644 --- a/include/Algo/Import/importTet.hpp +++ b/include/Algo/Import/importTet.hpp @@ -91,8 +91,8 @@ bool importTet(typename PFP::MAP& map, const std::string& filename, std::vector< verticesID.push_back(id); } - m_nbVertices = nbv; + m_nbVertices = nbv; m_nbVolumes = nbt; CGoGNout << "nb points = " << m_nbVertices << " / nb tet = " << m_nbVolumes << CGoGNendl; @@ -129,7 +129,6 @@ bool importTet(typename PFP::MAP& map, const std::string& filename, std::vector< for(unsigned int j = 0 ; j < 3 ; ++j) { FunctorSetEmb fsetemb(map, VERTEX, verticesID[pt[2-j]]); -// foreach_dart_of_orbit_in_parent(&map, VERTEX, d, fsetemb) ; map.foreach_dart_of_orbit( PFP::MAP::ORBIT_IN_PARENT(VERTEX), d, fsetemb); @@ -149,7 +148,6 @@ bool importTet(typename PFP::MAP& map, const std::string& filename, std::vector< d = map.phi_1(map.phi2(d)); FunctorSetEmb fsetemb(map, VERTEX, verticesID[pt[3]]); -// foreach_dart_of_orbit_in_parent(&map, VERTEX, d, fsetemb) ; map.foreach_dart_of_orbit( PFP::MAP::ORBIT_IN_PARENT(VERTEX), d, fsetemb); diff --git a/include/Algo/Import/importTs.hpp b/include/Algo/Import/importTs.hpp index 0ec59125..59435bed 100644 --- a/include/Algo/Import/importTs.hpp +++ b/include/Algo/Import/importTs.hpp @@ -45,12 +45,13 @@ bool importTs(typename PFP::MAP& map, const std::string& filename, std::vector position = map.template addAttribute(VERTEX, "position") ; attrNames.push_back(position.name()) ; + AttributeHandler scalaire = map.template addAttribute(VERTEX, "scalar"); attrNames.push_back(scalaire.name()) ; AttributeContainer& container = map.getAttributeContainer(VERTEX) ; - unsigned int m_nbVertices = 0, m_nbFaces = 0, m_nbEdges = 0, m_nbVolumes = 0; + unsigned int m_nbVertices = 0, m_nbVolumes = 0; AutoAttributeHandler< NoMathIONameAttribute< std::vector > > vecDartsPerVertex(map, VERTEX, "incidents"); // open file @@ -61,21 +62,15 @@ bool importTs(typename PFP::MAP& map, const std::string& filename, std::vector> nbv; - - //CGoGNout << "nbV = " << nbv << CGoGNendl; oss >> nbt; - //CGoGNout << "nbT = " << nbt << CGoGNendl; - - //lecture sommets + //reading vertices std::vector verticesID; verticesID.reserve(nbv); for(unsigned int i = 0; i < nbv;++i) @@ -92,10 +87,8 @@ bool importTs(typename PFP::MAP& map, const std::string& filename, std::vector> x; oss >> y; oss >> z; - // on peut ajouter ici la lecture de couleur si elle existe - VEC3 pos(x*scaleFactor,y*scaleFactor,z*scaleFactor); - //CGoGNout << "VEC3 = " << pos << CGoGNendl; + VEC3 pos(x*scaleFactor,y*scaleFactor,z*scaleFactor); unsigned int id = container.insertLine(); position[id] = pos; @@ -106,22 +99,19 @@ bool importTs(typename PFP::MAP& map, const std::string& filename, std::vector> nbe; //number of vertices = 4 + assert(nbe == 4); - //Algo::Modelisation::Polyhedron::createTetra(map); Dart d = Algo::Modelisation::createTetrahedron(map); + Geom::Vec4ui pt; oss >> pt[0]; oss >> pt[1]; oss >> pt[2]; oss >> pt[3]; - //regions ? - oss >> nbe; - -// CGoGNout << "\t embedding number : " << pt[0] << " " << pt[1] << " " << pt[2] << " " << pt[3] << CGoGNendl; + //if regions are defined use this number + oss >> nbe; //ignored here - // Embed three vertices + // Embed three "base" vertices for(unsigned int j = 0 ; j < 3 ; ++j) { -// CGoGNout << "\t embedding number : " << pt[j]; + FunctorSetEmb fsetemb(map, VERTEX, verticesID[pt[2-j]]); + map.foreach_dart_of_orbit( PFP::MAP::ORBIT_IN_PARENT(VERTEX), d, fsetemb); - FunctorSetEmb femb(map, VERTEX, verticesID[pt[j]]); + //store darts per vertices to optimize reconstruction Dart dd = d; do { - femb(dd); - //vecDartPtrEmb[pt[j]].push_back(dd); - vecDartsPerVertex[pt[j]].push_back(dd); + m.mark(dd) ; + vecDartsPerVertex[pt[2-j]].push_back(dd); dd = map.phi1(map.phi2(dd)); } while(dd != d); d = map.phi1(d); - -// CGoGNout << " done" << CGoGNendl; } - //Embed the last vertex -// CGoGNout << "\t embedding number : " << pt[3] << CGoGNendl; + //Embed the last "top" vertex d = map.phi_1(map.phi2(d)); - FunctorSetEmb femb(map, VERTEX, verticesID[pt[3]]); + FunctorSetEmb fsetemb(map, VERTEX, verticesID[pt[3]]); + map.foreach_dart_of_orbit( PFP::MAP::ORBIT_IN_PARENT(VERTEX), d, fsetemb); + + //store darts per vertices to optimize reconstruction Dart dd = d; do { - femb(dd); -// CGoGNout << "embed" << CGoGNendl; - //vecDartPtrEmb[pt[3]].push_back(dd); + m.mark(dd) ; vecDartsPerVertex[pt[3]].push_back(dd); dd = map.phi1(map.phi2(dd)); } while(dd != d); -// CGoGNout << "end tetra" << CGoGNendl; + //end of tetra } -// CGoGNout << "end 1/2" << CGoGNendl; - //Association des phi3 + unsigned int nbBoundaryFaces = 0 ; for (Dart d = map.begin(); d != map.end(); map.next(d)) { - std::vector& vec = vecDartsPerVertex[d]; - - for(typename std::vector::iterator it = vec.begin(); it != vec.end(); ++it) + if (m.isMarked(d)) { - if(map.phi3(*it)==*it) + std::vector& vec = vecDartsPerVertex[map.phi1(d)]; + + Dart good_dart = NIL; + for(typename std::vector::iterator it = vec.begin(); it != vec.end() && good_dart == NIL; ++it) { - bool sewn=false; - for(typename std::vector::iterator itnext = it+1; itnext != vec.end() && !sewn; ++itnext) + if(map.getEmbedding(VERTEX, map.phi1(*it)) == map.getEmbedding(VERTEX, d) && + map.getEmbedding(VERTEX, map.phi_1(*it)) == map.getEmbedding(VERTEX, map.phi_1(d)) /*&& + map.getEmbedding(VERTEX, *it) == map.getEmbedding(VERTEX, map.phi1(d)) */) { - if(map.getEmbedding(VERTEX,map.phi1(*it))==map.getEmbedding(VERTEX,map.phi_1(*itnext)) - && map.getEmbedding(VERTEX,map.phi_1(*it))==map.getEmbedding(VERTEX,map.phi1(*itnext))) - { - map.sewVolumes(*it, map.phi_1(*itnext)); - sewn = true; - } + good_dart = *it ; } } + + if (good_dart != NIL) + { + map.sewVolumes(d, good_dart, false); + m.unmarkOrbit(FACE, d); + } + else + { + m.unmarkOrbit(PFP::MAP::ORBIT_IN_PARENT(FACE), d); + ++nbBoundaryFaces; + } } } + + if (nbBoundaryFaces > 0) + { + std::cout << "closing" << std::endl ; + map.closeMap(); + CGoGNout << "Map closed (" << nbBoundaryFaces << " boundary faces)" << CGoGNendl; + } + fp.close(); return true; } diff --git a/include/Algo/Modelisation/subdivision3.hpp b/include/Algo/Modelisation/subdivision3.hpp index 1e69a85d..d837d7ad 100644 --- a/include/Algo/Modelisation/subdivision3.hpp +++ b/include/Algo/Modelisation/subdivision3.hpp @@ -82,188 +82,395 @@ Dart cut3Ear(typename PFP::MAP& map, Dart d) return map.phi2(dRing); } +//template +//void catmullClarkVol(typename PFP::MAP& map, EMBV& attributs, const FunctorSelect& selected) +//{ +// std::vector l_centers; +// std::vector l_vertices; +// +// DartMarkerNoUnmark mv(map); +// CellMarkerNoUnmark me(map, EDGE); +// CellMarker mf(map, FACE); +// +// AutoAttributeHandler< EMB > attBary(map, VOLUME); +// CellMarker vol(map, VOLUME); +// +// //pre-computation : compute the centroid of all volume +// for (Dart d = map.begin(); d != map.end(); map.next(d)) +// { +// if(selected(d) && !vol.isMarked(d)) +// { +// vol.mark(d); +// attBary[d] = Algo::Geometry::volumeCentroidGen(map, d, attributs); +// } +// } +// +// // first pass: cut edges +// for (Dart d = map.begin(); d != map.end(); map.next(d)) +// { +// //memorize each vertices per volumes +// if(selected(d) && !mv.isMarked(d)) +// { +// l_vertices.push_back(d); +//// mv.markOrbitInParent(VERTEX,d); +// mv.markOrbit(PFP::MAP::ORBIT_IN_PARENT(VERTEX),d); +// } +// +// //cut edges +// if (selected(d) && !me.isMarked(d)) +// { +// Dart f = map.phi1(d); +// map.cutEdge(d); +// Dart e = map.phi1(d) ; +// +// attributs[e] = attributs[d]; +// attributs[e] += attributs[f]; +// attributs[e] *= 0.5; +// +// me.mark(d); +// me.mark(e); +// +// //mark new vertices +// mv.markOrbit(VERTEX, e); +// +// Dart dd = d; +// do +// { +// mf.mark(dd) ; +// mf.mark(map.phi2(dd)); +// dd = map.alpha2(dd); +// } while(dd != d); +// } +// } +// +// // second pass: quandrangule faces +// std::map toSew; +// for (Dart d = map.begin(); d != map.end(); map.next(d)) +// { +// mv.unmark(d); +// +// if (selected(d) && mf.isMarked(d)) // for each face not subdivided +// { +// mf.unmark(d); +// // compute center skip darts of new vertices non embedded +// EMB center = AttribOps::zero(); +// unsigned int count = 0 ; +// Dart it = d; +// +// do +// { +// me.unmark(it); +// me.unmark(map.phi1(it)); +// +// center += attributs[it]; +// ++count ; +// +// it = map.template phi<11>(it) ; +// } while(it != d) ; +// +// center /= double(count); +// +// Dart cf = quadranguleFace(map, d); // quadrangule the face +// attributs[cf] = center; // affect the data to the central vertex +// } +// } +// +// //third pass : create the inner faces +// for (std::vector::iterator it = l_vertices.begin(); it != l_vertices.end(); ++it) +// { +// Dart d = *it; +// //unsew all around the vertex +// //there are 2 links to unsew for each face around (-> quadrangulation) +// std::vector v; +// do +// { +// v.push_back(map.phi1(map.phi1(d))); +// v.push_back(map.phi1(d)); +// +// d = map.phi2(map.phi_1(d)); +// } +// while(d != *it); +// +//// do +//// { +//// Dart dN = map.phi1(map.phi2(d)); +//// +//// Dart dRing = map.phi1(d); +//// +//// if(map.phi2(dRing)!=dRing) +//// { +//// toSew.insert(std::pair(dRing,map.phi2(dRing))); +//// v.push_back(dRing); +//// } +//// +//// dRing = map.phi1(dRing); +//// +//// if(map.phi2(dRing)!=dRing) +//// { +//// toSew.insert(std::pair(dRing,map.phi2(dRing))); +//// v.push_back(dRing); +//// } +//// +//// d = dN; +//// } while (d != *it); +// +//// //close the generated hole and create the central vertex +//// //unsigned int degree = map.closeHole(map.phi1(d)); +// +// //TODO : pb de face en trop avec splitVolume +// //map.splitVolume(v); +// +// // +// +//// Dart e = v.front(); +//// for(std::vector::iterator it = v.begin() ; it != v.end() ; ++it) +//// if(map.Map2::isBoundaryEdge(*it)) +//// map.unsewFaces(*it); +// +//// Dart dd = map.phi1(map.phi2(map.phi1(d))); +//// map.splitFace(map.phi_1(dd),map.phi1(dd)); +//// Dart dS = map.phi1(dd); +//// map.cutEdge(dS); +// +//// attributs[map.phi1(dS)] = attBary[d]; +// +// +// //TODO : test with vertices with degree higher than 3 +//// for(unsigned int i=0; i < (degree/2)-2; ++i) +//// { +//// map.splitFace(map.phi2(dS),map.template phi<111>(map.phi2(dS))); +//// dS = map.template phi<111>(map.phi2(dS)); +//// } +// } +// +//// map.deleteVolume(map.phi3(map.phi2(map.phi1(l_vertices.front())))); +// +// map.check(); +// +// //sew all faces leading to the central vertex +// for (std::map::iterator it = toSew.begin(); it != toSew.end(); ++it) +// { +// +//// Dart f1 = map.phi2((*it).first); +//// Dart f2 = map.phi2((*it).second); +//// if(map.isBoundaryFace(f1) && map.isBoundaryFace(f2)) +//// { +//// map.sewVolumes(f1, f2); +//// } +// +// //Dart dT = map.phi2(it->first); +//// if(dT==map.phi3(dT)) +//// { +//// map.sewVolumes(dT,map.phi2(it->second)); +//// } +// } +//} + + template void catmullClarkVol(typename PFP::MAP& map, EMBV& attributs, const FunctorSelect& selected) { - std::vector l_centers; + //std::vector l_centers; std::vector l_vertices; - DartMarkerNoUnmark mv(map); - CellMarkerNoUnmark me(map, EDGE); - CellMarker mf(map, FACE); - - AutoAttributeHandler< EMB > attBary(map, VOLUME); - CellMarker vol(map, VOLUME); - //pre-computation : compute the centroid of all volume + AutoAttributeHandler< EMB > attBary(map, VOLUME); + TraversorW traW(map, selected); for (Dart d = map.begin(); d != map.end(); map.next(d)) { - if(selected(d) && !vol.isMarked(d)) - { - vol.mark(d); - attBary[d] = Algo::Geometry::volumeCentroidGen(map, d, attributs); - } + attBary[d] = Algo::Geometry::volumeCentroidGen(map, d, attributs); } - // first pass: cut edges - for (Dart d = map.begin(); d != map.end(); map.next(d)) + //subdivision + //1. cut edges + DartMarkerNoUnmark mv(map); + TraversorE travE(map); + for (Dart d = travE.begin(); d != travE.end(); d = travE.next()) { //memorize each vertices per volumes if(selected(d) && !mv.isMarked(d)) { l_vertices.push_back(d); -// mv.markOrbitInParent(VERTEX,d); mv.markOrbit(PFP::MAP::ORBIT_IN_PARENT(VERTEX),d); } - //cut edges - if (selected(d) && !me.isMarked(d)) - { - Dart f = map.phi1(d); - map.cutEdge(d); - Dart e = map.phi1(d) ; - - attributs[e] = attributs[d]; - attributs[e] += attributs[f]; - attributs[e] *= 0.5; - - me.mark(d); - me.mark(e); + Dart f = map.phi1(d); + map.cutEdge(d) ; + Dart e = map.phi1(d) ; - //mark new vertices - mv.markOrbit(VERTEX, e); + attributs[e] = attributs[d]; + attributs[e] += attributs[f]; + attributs[e] *= 0.5; - Dart dd = d; - do - { - mf.mark(dd) ; - mf.mark(map.phi2(dd)); - dd = map.alpha2(dd); - } while(dd != d); - } + travE.mark(d) ; + travE.mark(e) ; } - // second pass: quandrangule faces - std::map toSew; - for (Dart d = map.begin(); d != map.end(); map.next(d)) + //2. split faces - quadrangule faces + TraversorF travF(map) ; + for (Dart d = travF.begin(); d != travF.end(); d = travF.next()) { - mv.unmark(d); - - if (selected(d) && mf.isMarked(d)) // for each face not subdivided - { - mf.unmark(d); - // compute center skip darts of new vertices non embedded - EMB center = AttribOps::zero(); - unsigned int count = 0 ; - Dart it = d; + EMB center = Algo::Geometry::faceCentroidGen(map,d,attributs); - do - { - me.unmark(it); - me.unmark(map.phi1(it)); - - center += attributs[it]; - ++count ; + Dart dd = map.phi1(d) ; + Dart next = map.phi1(map.phi1(dd)) ; + map.splitFace(dd, next) ; - it = map.template phi<11>(it) ; - } while(it != d) ; + Dart ne = map.phi2(map.phi_1(dd)) ; + map.cutEdge(ne) ; + travF.mark(dd) ; - center /= double(count); + attributs[map.phi1(ne)] = center; - Dart cf = quadranguleFace(map, d); // quadrangule the face - attributs[cf] = center; // affect the data to the central vertex + dd = map.phi1(map.phi1(next)) ; + while(dd != ne) + { + Dart tmp = map.phi1(ne) ; + map.splitFace(tmp, dd) ; + travF.mark(tmp) ; + dd = map.phi1(map.phi1(dd)) ; } + + travF.mark(ne) ; } - //third pass : create the inner faces + //3. create inside volumes + + std::vector > subdividedFaces; + subdividedFaces.reserve(2048); for (std::vector::iterator it = l_vertices.begin(); it != l_vertices.end(); ++it) { - Dart d = *it; - //unsew all around the vertex - //there are 2 links to unsew for each face around (-> quadrangulation) - std::vector v; + Dart e = *it; + std::vector v ; + do { - v.push_back(map.phi1(map.phi1(d))); - v.push_back(map.phi1(d)); + v.push_back(map.phi1(e)); + v.push_back(map.phi1(map.phi1(e))); + + if(!map.PFP::MAP::ParentMap::isBoundaryEdge(map.phi1(e))) + subdividedFaces.push_back(std::pair(map.phi1(e),map.phi2(map.phi1(e)))); - d = map.phi2(map.phi_1(d)); + if(!map.PFP::MAP::ParentMap::isBoundaryEdge(map.phi1(map.phi1(e)))) + subdividedFaces.push_back(std::pair(map.phi1(map.phi1(e)),map.phi2(map.phi1(map.phi1(e))))); + + e = map.phi2(map.phi_1(e)); } - while(d != *it); + while(e != *it); -// do -// { -// Dart dN = map.phi1(map.phi2(d)); -// -// Dart dRing = map.phi1(d); -// -// if(map.phi2(dRing)!=dRing) -// { -// toSew.insert(std::pair(dRing,map.phi2(dRing))); -// v.push_back(dRing); -// } -// -// dRing = map.phi1(dRing); -// -// if(map.phi2(dRing)!=dRing) -// { -// toSew.insert(std::pair(dRing,map.phi2(dRing))); -// v.push_back(dRing); -// } -// -// d = dN; -// } while (d != *it); + // + // SplitSurfaceInVolume + // -// //close the generated hole and create the central vertex -// //unsigned int degree = map.closeHole(map.phi1(d)); + std::vector vd2 ; + vd2.reserve(v.size()); - //TODO : pb de face en trop avec splitVolume - //map.splitVolume(v); + // save the edge neighbors darts + for(std::vector::iterator it2 = v.begin() ; it2 != v.end() ; ++it2) + { + vd2.push_back(map.phi2(*it2)); + } - // + assert(vd2.size() == v.size()); -// Dart e = v.front(); -// for(std::vector::iterator it = v.begin() ; it != v.end() ; ++it) -// if(map.Map2::isBoundaryEdge(*it)) -// map.unsewFaces(*it); + map.PFP::MAP::ParentMap::splitSurface(v, true, false); -// Dart dd = map.phi1(map.phi2(map.phi1(d))); -// map.splitFace(map.phi_1(dd),map.phi1(dd)); -// Dart dS = map.phi1(dd); -// map.cutEdge(dS); + // follow the edge path a second time to embed the vertex, edge and volume orbits + for(unsigned int i = 0; i < v.size(); ++i) + { + Dart dit = v[i]; + Dart dit2 = vd2[i]; -// attributs[map.phi1(dS)] = attBary[d]; + // embed the vertex embedded from the origin volume to the new darts + if(map.isOrbitEmbedded(VERTEX)) + { + map.copyDartEmbedding(VERTEX, map.phi2(dit), map.phi1(dit)); + map.copyDartEmbedding(VERTEX, map.phi2(dit2), map.phi1(dit2)); + } + // embed the edge embedded from the origin volume to the new darts + if(map.isOrbitEmbedded(EDGE)) + { + unsigned int eEmb = map.getEmbedding(EDGE, dit) ; + map.setDartEmbedding(EDGE, map.phi2(dit), eEmb); + map.setDartEmbedding(EDGE, map.phi2(dit2), eEmb); + } - //TODO : test with vertices with degree higher than 3 -// for(unsigned int i=0; i < (degree/2)-2; ++i) -// { -// map.splitFace(map.phi2(dS),map.template phi<111>(map.phi2(dS))); -// dS = map.template phi<111>(map.phi2(dS)); -// } - } + // embed the volume embedded from the origin volume to the new darts + if(map.isOrbitEmbedded(VOLUME)) + { + map.copyDartEmbedding(VOLUME, map.phi2(dit), dit); + map.copyDartEmbedding(VOLUME, map.phi2(dit2), dit2); + } + } -// map.deleteVolume(map.phi3(map.phi2(map.phi1(l_vertices.front())))); + // + // + // - map.check(); + Dart dd = map.phi2(map.phi1(*it)); + Dart next = map.phi1(map.phi1(dd)) ; + map.PFP::MAP::ParentMap::splitFace(dd, next); - //sew all faces leading to the central vertex - for (std::map::iterator it = toSew.begin(); it != toSew.end(); ++it) - { + if (map.isOrbitEmbedded(VERTEX)) + { + map.copyDartEmbedding(VERTEX, map.phi_1(next), dd) ; + map.copyDartEmbedding(VERTEX, map.phi_1(dd), next) ; + } + + Dart ne = map.phi2(map.phi_1(dd)); + map.PFP::MAP::ParentMap::cutEdge(ne); -// Dart f1 = map.phi2((*it).first); -// Dart f2 = map.phi2((*it).second); -// if(map.isBoundaryFace(f1) && map.isBoundaryFace(f2)) + +// dd = map.phi1(map.phi1(next)) ; +// while(dd != ne) // { -// map.sewVolumes(f1, f2); +// Dart tmp = map.phi1(ne) ; +// map.PFP::MAP::ParentMap::splitFace(tmp, dd); +// +// if (map.isOrbitEmbedded(VERTEX)) +// { +// map.copyDartEmbedding(VERTEX, map.phi_1(dd), tmp) ; +// map.copyDartEmbedding(VERTEX, map.phi_1(tmp), dd) ; +// } +// +// dd = map.phi1(map.phi1(dd)) ; // } +// + } - //Dart dT = map.phi2(it->first); -// if(dT==map.phi3(dT)) +// setCurrentLevel(getMaxLevel()) ; +// //4 couture des relations precedemment sauvegarde +// for (std::vector >::iterator it = subdividedFaces.begin(); it != subdividedFaces.end(); ++it) // { -// map.sewVolumes(dT,map.phi2(it->second)); +// Dart f1 = phi2((*it).first); +// Dart f2 = phi2((*it).second); +// +// //if(isBoundaryFace(f1) && isBoundaryFace(f2)) +// if(phi3(f1) == f1 && phi3(f2) == f2) +// sewVolumes(f1, f2, false); // } - } +// embedOrbit(VERTEX, centralDart, getEmbedding(VERTEX, centralDart)); + //attributs[map.phi1(ne)] = attBary[*it]; +// +// setCurrentLevel(getMaxLevel() - 1) ; +// } +// +// //A optimiser +// +// TraversorE travE2(map); +// for (Dart d = travE2.begin(); d != travE2.end(); d = travE2.next()) +// { +// map.embedOrbit(VERTEX, map.phi1(d), map.getEmbedding(VERTEX, map.phi1(d))); +// } +// +// TraversorF travF2(map) ; +// for (Dart d = travF2.begin(); d != travF2.end(); d = travF2.next()) +// { +// map.embedOrbit(VERTEX, map.phi2(map.phi1(d)), map.getEmbedding(VERTEX, map.phi2(map.phi1(d)))); +// } + + } } //namespace Modelisation diff --git a/include/Algo/Render/GL2/topo3Render.h b/include/Algo/Render/GL2/topo3Render.h index 0bad6748..d3cd4c1a 100644 --- a/include/Algo/Render/GL2/topo3Render.h +++ b/include/Algo/Render/GL2/topo3Render.h @@ -290,6 +290,9 @@ public: template void computeDartMiddlePositions(typename PFP::MAP& map, typename PFP::TVEC3& posExpl, const FunctorSelect& good = allDarts); + template + void updateDataMap3OldFashioned(typename PFP::MAP& mapx, const typename PFP::TVEC3& positions, float ke, float kf, float kv, const FunctorSelect& good); + protected: /** * update all drawing buffers to render a dual map diff --git a/include/Algo/Render/GL2/topo3Render.hpp b/include/Algo/Render/GL2/topo3Render.hpp index 32044cc7..d4e83e9e 100644 --- a/include/Algo/Render/GL2/topo3Render.hpp +++ b/include/Algo/Render/GL2/topo3Render.hpp @@ -537,6 +537,230 @@ void Topo3Render::computeDartMiddlePositions(typename PFP::MAP& map, typename PF } +template +void Topo3Render::updateDataMap3OldFashioned(typename PFP::MAP& mapx, const typename PFP::TVEC3& positions, float ke, float kf, float kv, const FunctorSelect& good) +{ + Map3& map = reinterpret_cast(mapx); + + typedef typename PFP::VEC3 VEC3; + typedef typename PFP::REAL REAL; + + + if (m_attIndex.map() != &map) + { + m_attIndex = map.template getAttribute(DART, "dart_index"); + if (!m_attIndex.isValid()) + m_attIndex = map.template addAttribute(DART, "dart_index"); + } + + m_nbDarts = 0; + + // table of center of volume + std::vector vecCenters; + vecCenters.reserve(1000); + // table of nbfaces per volume + std::vector vecNbFaces; + vecNbFaces.reserve(1000); + // table of face (one dart of each) + std::vector vecDartFaces; + vecDartFaces.reserve(map.getNbDarts()/4); + + unsigned int posDBI=0; + + DartMarker mark(map); // marker for darts + for (Dart d = map.begin(); d != map.end(); map.next(d)) + { + if (good(d)) + { + CellMarkerStore markVert(map, VERTEX); //marker for vertices + VEC3 center(0, 0, 0); + unsigned int nbv = 0; + unsigned int nbf = 0; + std::list visitedFaces; // Faces that are traversed + visitedFaces.push_back(d); // Start with the face of d + + // For every face added to the list + for (std::list::iterator face = visitedFaces.begin(); face != visitedFaces.end(); ++face) + { + if (!mark.isMarked(*face)) // Face has not been visited yet + { + // store a dart of face + vecDartFaces.push_back(*face); + nbf++; + Dart dNext = *face ; + do + { + if (!markVert.isMarked(dNext)) + { + markVert.mark(dNext); + center += positions[dNext]; + nbv++; + } + mark.mark(dNext); // Mark + m_nbDarts++; + Dart adj = map.phi2(dNext); // Get adjacent face + if (adj != dNext && !mark.isMarked(adj)) + visitedFaces.push_back(adj); // Add it + dNext = map.phi1(dNext); + } while(dNext != *face); + } + } + center /= typename PFP::REAL(nbv); + vecCenters.push_back(center); + vecNbFaces.push_back(nbf); + } + } + + // debut phi1 + AutoAttributeHandler fv1(map, DART); + // fin phi1 + AutoAttributeHandler fv11(map, DART); + + // phi2 + AutoAttributeHandler fv2(map, DART); + AutoAttributeHandler fv2x(map, DART); + + m_vbo4->bind(); + glBufferData(GL_ARRAY_BUFFER, 2*m_nbDarts*sizeof(VEC3), 0, GL_STREAM_DRAW); + GLvoid* ColorDartsBuffer = glMapBuffer(GL_ARRAY_BUFFER, GL_READ_WRITE); + VEC3* colorDartBuf = reinterpret_cast(ColorDartsBuffer); + + m_vbo0->bind(); + glBufferData(GL_ARRAY_BUFFER, 2*m_nbDarts*sizeof(VEC3), 0, GL_STREAM_DRAW); + GLvoid* PositionDartsBuffer = glMapBuffer(GL_ARRAY_BUFFER, GL_READ_WRITE); + VEC3* positionDartBuf = reinterpret_cast(PositionDartsBuffer); + + + + std::vector::iterator face = vecDartFaces.begin(); + for (unsigned int iVol=0; iVol vecPos; + vecPos.reserve(16); + + // store the face & center + VEC3 center(0, 0, 0); + Dart dd = d; + do + { + const VEC3& P = positions[d]; + vecPos.push_back(P); + center += P; + d = map.phi1(d); + } while (d != dd); + center /= REAL(vecPos.size()); + + //shrink the face + unsigned int nb = vecPos.size(); + float okf = 1.0f - kf; + float okv = 1.0f - kv; + for (unsigned int i = 0; i < nb; ++i) + { + vecPos[i] = vecCenters[iVol]*okv + vecPos[i]*kv; + vecPos[i] = center*okf + vecPos[i]*kf; + } + vecPos.push_back(vecPos.front()); // copy the first for easy computation on next loop + + // compute position of points to use for drawing topo + float oke = 1.0f - ke; + for (unsigned int i = 0; i < nb; ++i) + { + VEC3 P = vecPos[i]*ke + vecPos[i+1]*oke; + VEC3 Q = vecPos[i+1]*ke + vecPos[i]*oke; + + m_attIndex[d] = posDBI; + posDBI+=2; + + *positionDartBuf++ = P; + *positionDartBuf++ = Q; + *colorDartBuf++ = VEC3(1.,1.,1.0); + *colorDartBuf++ = VEC3(1.,1.,1.0); + + fv1[d] = P*0.1f + Q*0.9f; + fv11[d] = P*0.9f + Q*0.1f; + + fv2[d] = P*0.52f + Q*0.48f; + fv2x[d] = P*0.48f + Q*0.52f; + + d = map.phi1(d); + } + + } + } + + m_vbo0->bind(); + glUnmapBuffer(GL_ARRAY_BUFFER); + + // phi1 + m_vbo1->bind(); + glBufferData(GL_ARRAY_BUFFER, 2*m_nbDarts*sizeof(typename PFP::VEC3), 0, GL_STREAM_DRAW); + GLvoid* PositionBuffer1 = glMapBufferARB(GL_ARRAY_BUFFER, GL_READ_WRITE); + + //phi2 + m_vbo2->bind(); + glBufferData(GL_ARRAY_BUFFER, 2*m_nbDarts*sizeof(typename PFP::VEC3), 0, GL_STREAM_DRAW); + GLvoid* PositionBuffer2 = glMapBufferARB(GL_ARRAY_BUFFER, GL_READ_WRITE); + + //phi3 + m_vbo3->bind(); + glBufferData(GL_ARRAY_BUFFER, 2*m_nbDarts*sizeof(typename PFP::VEC3), 0, GL_STREAM_DRAW); + GLvoid* PositionBuffer3 = glMapBufferARB(GL_ARRAY_BUFFER, GL_READ_WRITE); + + VEC3* positionF1 = reinterpret_cast(PositionBuffer1); + VEC3* positionF2 = reinterpret_cast(PositionBuffer2); + VEC3* positionF3 = reinterpret_cast(PositionBuffer3); + + m_nbRel2=0; + m_nbRel3=0; + + for(std::vector::iterator face = vecDartFaces.begin(); face != vecDartFaces.end(); ++face) + { + Dart d = *face; + do + { + Dart e = map.phi2(d); + if (d < e) + { + *positionF2++ = fv2[d]; + *positionF2++ = fv2x[e]; + *positionF2++ = fv2[e]; + *positionF2++ = fv2x[d]; + m_nbRel2++; + } + e = map.phi3(d); + if (!map.isBoundaryMarked(e) && (d < e)) + { + *positionF3++ = fv2[d]; + *positionF3++ = fv2x[e]; + *positionF3++ = fv2[e]; + *positionF3++ = fv2x[d]; + m_nbRel3++; + } + e = map.phi1(d); + *positionF1++ = fv1[d]; + *positionF1++ = fv11[e]; + + d = map.phi1(d); + } while (d != *face ); + } + + m_vbo3->bind(); + glUnmapBufferARB(GL_ARRAY_BUFFER); + + m_vbo2->bind(); + glUnmapBufferARB(GL_ARRAY_BUFFER); + + m_vbo1->bind(); + glUnmapBuffer(GL_ARRAY_BUFFER); + + m_vbo4->bind(); + glUnmapBuffer(GL_ARRAY_BUFFER); +} + }//end namespace VBO diff --git a/include/Topology/generic/genericmap.hpp b/include/Topology/generic/genericmap.hpp index 89ac687a..95b335cd 100644 --- a/include/Topology/generic/genericmap.hpp +++ b/include/Topology/generic/genericmap.hpp @@ -115,7 +115,6 @@ inline void GenericMap::deleteDart(Dart d) for(unsigned int i = m_mrCurrentLevel; i < m_mrDarts.size(); ++i) { unsigned int index = (*m_mrDarts[i])[d.index] ; - std::cout << "deleteDart index = " << index << std::endl << std::endl; if(isDartValid(index)) deleteDartLine(index) ; } diff --git a/include/Topology/gmap/gmap2.h b/include/Topology/gmap/gmap2.h index 341ca757..b347cc7c 100644 --- a/include/Topology/gmap/gmap2.h +++ b/include/Topology/gmap/gmap2.h @@ -257,7 +257,7 @@ public: /*! * */ - virtual void splitCC(std::vector& vd); + virtual void splitSurface(std::vector& vd, bool firstSideClosed = true, bool secondSideClosed = true); //@} diff --git a/include/Topology/map/embeddedMap2.h b/include/Topology/map/embeddedMap2.h index ab094815..3a931ee8 100644 --- a/include/Topology/map/embeddedMap2.h +++ b/include/Topology/map/embeddedMap2.h @@ -135,6 +135,11 @@ public: */ virtual bool mergeVolumes(Dart d, Dart e) ; + /** + * + */ + virtual void splitSurface(std::vector& vd, bool firstSideClosed = true, bool secondSideClosed = true); + /** * No attribute is attached to the new face */ diff --git a/include/Topology/map/map3.h b/include/Topology/map/map3.h index da825e0a..578143a8 100644 --- a/include/Topology/map/map3.h +++ b/include/Topology/map/map3.h @@ -233,6 +233,11 @@ public: */ unsigned int vertexDegree(Dart d) ; + //! Compute the number of edges of the vertex of d on the boundary + /*! @param d a dart + */ + unsigned int vertexDegreeOnBoundary(Dart d); + //! Tell if the vertex of d is on the boundary /*! @param d a dart */ diff --git a/include/Topology/map/map3MR/filters_Primal.h b/include/Topology/map/map3MR/filters_Primal.h new file mode 100644 index 00000000..d1d6b93f --- /dev/null +++ b/include/Topology/map/map3MR/filters_Primal.h @@ -0,0 +1,752 @@ +/******************************************************************************* +* CGoGN: Combinatorial and Geometric modeling with Generic N-dimensional Maps * +* version 0.1 * +* Copyright (C) 2009-2012, IGG Team, LSIIT, University of Strasbourg * +* * +* This library is free software; you can redistribute it and/or modify it * +* under the terms of the GNU Lesser General Public License as published by the * +* Free Software Foundation; either version 2.1 of the License, or (at your * +* option) any later version. * +* * +* This library is distributed in the hope that it will be useful, but WITHOUT * +* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or * +* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License * +* for more details. * +* * +* You should have received a copy of the GNU Lesser General Public License * +* along with this library; if not, write to the Free Software Foundation, * +* Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. * +* * +* Web site: http://cgogn.unistra.fr/ * +* Contact information: cgogn@unistra.fr * +* * +*******************************************************************************/ + +#ifndef __3MR_FILTERS_PRIMAL__ +#define __3MR_FILTERS_PRIMAL__ + +#include +#include "Algo/Geometry/centroid.h" + +namespace CGoGN +{ + +namespace Multiresolution +{ + +class MRFilter +{ +public: + MRFilter() {} + virtual ~MRFilter() {} + virtual void operator() () = 0 ; +} ; + + +/********************************************************************************* + * LOOP BASIC FUNCTIONS + *********************************************************************************/ +template +typename PFP::VEC3 loopOddVertex(typename PFP::MAP& map, const typename PFP::TVEC3& position, Dart d1) +{ + Dart d2 = map.phi2(d1) ; + Dart d3 = map.phi_1(d1) ; + Dart d4 = map.phi_1(d2) ; + + typename PFP::VEC3 p1 = position[d1] ; + typename PFP::VEC3 p2 = position[d2] ; + typename PFP::VEC3 p3 = position[d3] ; + typename PFP::VEC3 p4 = position[d4] ; + + p1 *= 3.0 / 8.0 ; + p2 *= 3.0 / 8.0 ; + p3 *= 1.0 / 8.0 ; + p4 *= 1.0 / 8.0 ; + + return p1 + p2 + p3 + p4 ; +} + +template +typename PFP::VEC3 loopEvenVertex(typename PFP::MAP& map, const typename PFP::TVEC3& position, Dart d) +{ + map.incCurrentLevel() ; + + typename PFP::VEC3 np(0) ; + unsigned int degree = 0 ; + Traversor2VVaE trav(map, d) ; + for(Dart it = trav.begin(); it != trav.end(); it = trav.next()) + { + ++degree ; + np += position[it] ; + } + + map.decCurrentLevel() ; + + float mu = 3.0/8.0 + 1.0/4.0 * cos(2.0 * M_PI / degree) ; + mu = (5.0/8.0 - (mu * mu)) / degree ; + np *= 8.0/5.0 * mu ; + + return np ; +} + +/********************************************************************************* + * SHW04 BASIC FUNCTIONS : tetrahedral/octahedral meshes + *********************************************************************************/ +template +typename PFP::VEC3 SHW04Vertex(typename PFP::MAP& map, const typename PFP::TVEC3& position, Dart d) +{ + typename PFP::VEC3 res(0); + + if(map.isTetrahedron(d)) + { + Dart d1 = map.phi1(d) ; + Dart d2 = map.phi_1(d); + Dart d3 = map.phi_1(map.phi2(d)); + + typename PFP::VEC3 p = position[d]; + typename PFP::VEC3 p1 = position[d1] ; + typename PFP::VEC3 p2 = position[d2] ; + typename PFP::VEC3 p3 = position[d3] ; + + p *= -1; + p1 *= 17.0 / 3.0; + p2 *= 17.0 / 3.0; + p3 *= 17.0 / 3.0; + + res += p + p1 + p2 + p3; + res *= 1.0 / 16.0; + } + else + { + Dart d1 = map.phi1(d); + Dart d2 = map.phi_1(d); + Dart d3 = map.phi_1(map.phi2(d)); + Dart d4 = map.phi_1(map.phi2(d3)); + Dart d5 = map.phi_1(map.phi2(map.phi_1(d))); + + typename PFP::VEC3 p = position[d]; + typename PFP::VEC3 p1 = position[d1] ; + typename PFP::VEC3 p2 = position[d2] ; + typename PFP::VEC3 p3 = position[d3] ; + typename PFP::VEC3 p4 = position[d4] ; + typename PFP::VEC3 p5 = position[d5] ; + + p *= 3.0 / 4.0; + p1 *= 1.0 / 6.0; + p2 *= 1.0 / 6.0; + p3 *= 1.0 / 6.0; + p4 *= 7.0 / 12.0; + p5 *= 1.0 / 6.0; + + res += p + p1 + p2 + p3 + p4 + p5; + res *= 1.0 / 2.0; + } + + return res; +} + +/********************************************************************************* + * BSXW02 BASIC FUNCTIONS : polyhedral meshes + *********************************************************************************/ + + + +/********************************************************************************* + * ANALYSIS FILTERS + *********************************************************************************/ + +/* Loop + *********************************************************************************/ +template +class LoopEvenAnalysisFilter : public MRFilter +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + LoopEvenAnalysisFilter(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + void operator() () + { + TraversorV trav(m_map) ; + for (Dart d = trav.begin(); d != trav.end(); d = trav.next()) + { + if(m_map.isBoundaryVertex(d)) + { + Dart db = m_map.findBoundaryFaceOfVertex(d); + typename PFP::VEC3 p = loopEvenVertex(m_map, m_position, db) ; + m_position[db] -= p ; + } + } + } +} ; + +template +class LoopNormalisationAnalysisFilter : public MRFilter +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + LoopNormalisationAnalysisFilter(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + void operator() () + { + TraversorV trav(m_map) ; + for (Dart d = trav.begin(); d != trav.end(); d = trav.next()) + { + if(m_map.isBoundaryVertex(d)) + { + Dart db = m_map.findBoundaryFaceOfVertex(d); + + unsigned int degree = m_map.vertexDegreeOnBoundary(db) ; + float n = 3.0/8.0 + 1.0/4.0 * cos(2.0 * M_PI / degree) ; + n = 8.0/5.0 * (n * n) ; + + m_position[db] /= n ; + } + } + } +} ; + +template +class LoopOddAnalysisFilter : public MRFilter +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + LoopOddAnalysisFilter(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + void operator() () + { + TraversorE trav(m_map) ; + for (Dart d = trav.begin(); d != trav.end(); d = trav.next()) + { + if(m_map.isBoundaryEdge(d)) + { + Dart db = m_map.findBoundaryFaceOfEdge(d); + typename PFP::VEC3 p = loopOddVertex(m_map, m_position, db) ; + + m_map.incCurrentLevel() ; + + Dart oddV = m_map.phi2(db) ; + m_position[oddV] -= p ; + + m_map.decCurrentLevel() ; + } + else + { + typename PFP::VEC3 p = (m_position[d] + m_position[m_map.phi2(d)]) * typename PFP::REAL(0.5); + + m_map.incCurrentLevel() ; + + Dart midV = m_map.phi2(d) ; + m_position[midV] -= p ; + + m_map.decCurrentLevel() ; + } + } + } +} ; + +/********************************************************************************* + * SYNTHESIS FILTERS + *********************************************************************************/ + +/* Linear Interpolation + *********************************************************************************/ +template +class LerpEdgeSynthesisFilter : public MRFilter +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + LerpEdgeSynthesisFilter(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + void operator() () + { + TraversorE trav(m_map) ; + for (Dart d = trav.begin(); d != trav.end(); d = trav.next()) + { + typename PFP::VEC3 p = (m_position[d] + m_position[m_map.phi2(d)]) * typename PFP::REAL(0.5); + + m_map.incCurrentLevel() ; + + Dart midV = m_map.phi2(d) ; + m_position[midV] = p ; + + m_map.decCurrentLevel() ; + } + } +} ; + +template +class LerpFaceSynthesisFilter : public MRFilter +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + LerpFaceSynthesisFilter(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + void operator() () + { + TraversorF trav(m_map) ; + for (Dart d = trav.begin(); d != trav.end(); d = trav.next()) + { + typename PFP::VEC3 p = Algo::Geometry::faceCentroid(m_map, d, m_position); + + m_map.incCurrentLevel() ; + + Dart midF = m_map.phi2(m_map.phi1(d)); + m_position[midF] = p ; + + m_map.decCurrentLevel() ; + } + } +} ; + +template +class LerpVolumeSynthesisFilter : public MRFilter +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + LerpVolumeSynthesisFilter(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + void operator() () + { + TraversorW trav(m_map) ; + for (Dart d = trav.begin(); d != trav.end(); d = trav.next()) + { + typename PFP::VEC3 p = Algo::Geometry::volumeCentroid(m_map, d, m_position); + + m_map.incCurrentLevel() ; + + Dart midV = m_map.phi_1(m_map.phi2(m_map.phi1(d))); + m_position[midV] = p ; + + m_map.decCurrentLevel() ; + } + } +} ; + +/* Loop on Boundary Vertices + *********************************************************************************/ +template +class LoopEvenSynthesisFilter : public MRFilter +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + LoopEvenSynthesisFilter(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + void operator() () + { + TraversorV trav(m_map) ; + for (Dart d = trav.begin(); d != trav.end(); d = trav.next()) + { + if(m_map.isBoundaryVertex(d)) + { + Dart db = m_map.findBoundaryFaceOfVertex(d); + typename PFP::VEC3 p = loopEvenVertex(m_map, m_position, db) ; + m_position[db] += p ; + } + } + } +} ; + +template +class LoopNormalisationSynthesisFilter : public MRFilter +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + LoopNormalisationSynthesisFilter(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + void operator() () + { + TraversorV trav(m_map) ; + for (Dart d = trav.begin(); d != trav.end(); d = trav.next()) + { + if(m_map.isBoundaryVertex(d)) + { + Dart db = m_map.findBoundaryFaceOfVertex(d); + + unsigned int degree = m_map.vertexDegreeOnBoundary(db) ; + float n = 3.0/8.0 + 1.0/4.0 * cos(2.0 * M_PI / degree) ; + n = 8.0/5.0 * (n * n) ; + + m_position[db] *= n ; + } + } + } +} ; + +template +class LoopOddSynthesisFilter : public MRFilter +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + LoopOddSynthesisFilter(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + void operator() () + { + TraversorE trav(m_map) ; + for (Dart d = trav.begin(); d != trav.end(); d = trav.next()) + { + if(m_map.isBoundaryEdge(d)) + { + Dart db = m_map.findBoundaryFaceOfEdge(d); + typename PFP::VEC3 p = loopOddVertex(m_map, m_position, db) ; + + m_map.incCurrentLevel() ; + + Dart oddV = m_map.phi2(db) ; + m_position[oddV] += p ; + + m_map.decCurrentLevel() ; + } + else + { + typename PFP::VEC3 p = (m_position[d] + m_position[m_map.phi2(d)]) * typename PFP::REAL(0.5); + + m_map.incCurrentLevel() ; + + Dart midV = m_map.phi2(d) ; + m_position[midV] += p ; + + m_map.decCurrentLevel() ; + } + } + } +} ; + +template +class LoopVolumeSynthesisFilter : public MRFilter +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + LoopVolumeSynthesisFilter(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + void operator() () + { + TraversorW trav(m_map) ; + for (Dart d = trav.begin(); d != trav.end(); d = trav.next()) + { + if(!m_map.isTetrahedron(d)) + { + typename PFP::VEC3 p = Algo::Geometry::volumeCentroid(m_map, d, m_position); + + m_map.incCurrentLevel() ; + + Dart midV = m_map.phi_1(m_map.phi2(m_map.phi1(d))); + m_position[midV] = p ; + + m_map.decCurrentLevel() ; + } + } + } +} ; + +template +class SHW04VolumeNormalisationSynthesisFilter : public MRFilter +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + SHW04VolumeNormalisationSynthesisFilter(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + void operator() () + { + m_map.incCurrentLevel() ; + TraversorV trav(m_map) ; + for (Dart d = trav.begin(); d != trav.end(); d = trav.next()) + { + if(!m_map.isBoundaryVertex(d)) + { + typename PFP::VEC3 p = typename PFP::VEC3(0); + unsigned int degree = 0; + + Traversor3VW travVW(m_map, d); + for(Dart dit = travVW.begin() ; dit != travVW.end() ; dit = travVW.next()) + { + p += SHW04Vertex(m_map, m_position, dit); + ++degree; + } + + p /= degree; + + m_position[d] = p ; + } + } + m_map.decCurrentLevel() ; + } +} ; + +/********************************************************************************* + * FUNCTORS + *********************************************************************************/ + +/* Linear Interpolation + *********************************************************************************/ +template +class LerpVertexVertexFunctor : public FunctorType +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + LerpVertexVertexFunctor(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + bool operator() (Dart d) + { + m_map.decCurrentLevel() ; + typename PFP::VEC3 p = m_position[d] ; + m_map.incCurrentLevel() ; + + m_position[d] = p ; + + return false ; + } +}; + +template +class LerpEdgeVertexFunctor : public FunctorType +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + LerpEdgeVertexFunctor(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + bool operator() (Dart d) + { + Dart d1 = m_map.phi2(d); + + m_map.decCurrentLevel(); + Dart d2 = m_map.phi2(d1); + typename PFP::VEC3 p = (m_position[d1] + m_position[d2]) * typename PFP::REAL(0.5); + m_map.incCurrentLevel(); + + m_position[d] = p; + + return false; + } +} ; + +template +class LerpFaceVertexFunctor : public FunctorType +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + LerpFaceVertexFunctor(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + bool operator() (Dart d) + { + Dart df = m_map.phi1(m_map.phi1(d)) ; + + m_map.decCurrentLevel() ; + typename PFP::VEC3 p = Algo::Geometry::faceCentroid(m_map, df, m_position); + m_map.incCurrentLevel() ; + + m_position[d] = p ; + + return false ; + } +} ; + +template +class LerpVolumeVertexFunctor : public FunctorType +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + LerpVolumeVertexFunctor(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + bool operator() (Dart d) + { + Dart df = m_map.phi_1(m_map.phi2(m_map.phi1(d))) ; + + m_map.decCurrentLevel() ; + typename PFP::VEC3 p = Algo::Geometry::volumeCentroid(m_map, df, m_position); + m_map.incCurrentLevel() ; + + m_position[d] = p ; + + return false ; + } +} ; + + +/* SHW04 basic functions : tetrahedral/octahedral meshes + *********************************************************************************/ +template +class SHW04VertexVertexFunctor : public FunctorType +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + SHW04VertexVertexFunctor(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + bool operator() (Dart d) + { + if(m_map.isBoundaryVertex(d)) + { + Dart db = m_map.findBoundaryFaceOfVertex(d); + + m_map.decCurrentLevel() ; + + typename PFP::VEC3 np(0) ; + unsigned int degree = 0 ; + Traversor2VVaE trav(m_map, db) ; + for(Dart it = trav.begin(); it != trav.end(); it = trav.next()) + { + ++degree ; + np += m_position[it] ; + } + float tmp = 3.0 + 2.0 * cos(2.0 * M_PI / degree) ; + float beta = (5.0 / 8.0) - ( tmp * tmp / 64.0 ) ; + np *= beta / degree ; + + typename PFP::VEC3 vp = m_position[db] ; + vp *= 1.0 - beta ; + + m_map.incCurrentLevel() ; + + m_position[d] = np + vp ; + } + else + { + std::cout << "not a boundary vertex vertex" << std::endl; +// typename PFP::VEC3 p = typename PFP::VEC3(0); +// unsigned int degree = 0; +// +// Traversor3VW travVW(m_map, d); +// for(Dart dit = travVW.begin() ; dit != travVW.end() ; dit = travVW.next()) +// { +// p += SHW04Vertex(m_map, m_position, dit); +// ++degree; +// } +// +// p /= degree; +// +// m_position[d] = p ; + } + return false ; + } +} ; + +template +class SHW04EdgeVertexFunctor : public FunctorType +{ +protected: + typename PFP::MAP& m_map ; + typename PFP::TVEC3& m_position ; + +public: + SHW04EdgeVertexFunctor(typename PFP::MAP& m, typename PFP::TVEC3& p) : m_map(m), m_position(p) + {} + + bool operator() (Dart d) + { + if(m_map.isBoundaryVertex(d)) + { + Dart dd = m_map.phi2(d) ; + m_map.decCurrentLevel() ; + + Dart d1 = m_map.findBoundaryFaceOfEdge(dd); + + Dart d2 = m_map.phi2(d1) ; + Dart d3 = m_map.phi_1(d1) ; + Dart d4 = m_map.phi_1(d2) ; + + typename PFP::VEC3 p1 = m_position[d1] ; + typename PFP::VEC3 p2 = m_position[d2] ; + typename PFP::VEC3 p3 = m_position[d3] ; + typename PFP::VEC3 p4 = m_position[d4] ; + + p1 *= 3.0 / 8.0 ; + p2 *= 3.0 / 8.0 ; + p3 *= 1.0 / 8.0 ; + p4 *= 1.0 / 8.0 ; + + m_map.incCurrentLevel() ; + + m_position[d] = p1 + p2 + p3 + p4 ; + } + else + { + std::cout << "not a boundary edge vertex" << std::endl; +// typename PFP::VEC3 p = typename PFP::VEC3(0); +// unsigned int degree = 0; +// +// Traversor3VW travVW(m_map, d); +// for(Dart dit = travVW.begin() ; dit != travVW.end() ; dit = travVW.next()) +// { +// p += SHW04Vertex(m_map, m_position, dit); +// ++degree; +// } +// +// p /= degree; +// +// m_position[d] = p ; + } + + return false ; + } +} ; + + +} // namespace Multiresolution + +} // namespace CGoGN + + +#endif /* __3MR_FILTERS_PRIMAL__ */ diff --git a/include/Topology/map/map3MR/map3MR_PrimalAdapt.h b/include/Topology/map/map3MR/map3MR_PrimalAdapt.h new file mode 100644 index 00000000..879b208d --- /dev/null +++ b/include/Topology/map/map3MR/map3MR_PrimalAdapt.h @@ -0,0 +1,196 @@ +/******************************************************************************* +* CGoGN: Combinatorial and Geometric modeling with Generic N-dimensional Maps * +* version 0.1 * +* Copyright (C) 2009-2012, IGG Team, LSIIT, University of Strasbourg * +* * +* This library is free software; you can redistribute it and/or modify it * +* under the terms of the GNU Lesser General Public License as published by the * +* Free Software Foundation; either version 2.1 of the License, or (at your * +* option) any later version. * +* * +* This library is distributed in the hope that it will be useful, but WITHOUT * +* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or * +* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License * +* for more details. * +* * +* You should have received a copy of the GNU Lesser General Public License * +* along with this library; if not, write to the Free Software Foundation, * +* Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. * +* * +* Web site: http://cgogn.unistra.fr/ * +* Contact information: cgogn@unistra.fr * +* * +*******************************************************************************/ + +#ifndef __MAP3MR_PRIMAL_ADAPT__ +#define __MAP3MR_PRIMAL_ADAPT__ + +#include "Topology/map/embeddedMap3.h" +#include "Topology/generic/traversorCell.h" +#include "Topology/generic/traversor3.h" +#include + +namespace CGoGN +{ + +/*! \brief The class of adaptive 3-map MR + */ + +class Map3MR_PrimalAdapt : public EmbeddedMap3 +{ +protected: + enum SubdivideType + { + S_TRI, + S_QUAD + } ; + + bool shareVertexEmbeddings; + + FunctorType* vertexVertexFunctor ; + FunctorType* edgeVertexFunctor ; + FunctorType* faceVertexFunctor ; + FunctorType* volumeVertexFunctor ; + +public: + Map3MR_PrimalAdapt(); + + std::string mapTypeName() { return "Map3MR_PrimalAdapt"; } + + /*! @name Topological helping functions + * + *************************************************************************/ + //@{ + //! + /*! + */ + void swapEdges(Dart d, Dart e); + + //! + /*! + * + */ + bool isTetrahedron(Dart d); + + //! + /*! + */ + void saveRelationsAroundVertex(Dart d, std::vector >& vd); + + //! + /*! + */ + void unsewAroundVertex(std::vector >& vd); + //@} + + /*! @name Cells informations + * + *************************************************************************/ + //@{ + //! Return the level of the edge of d in the current level map + /* @param d Dart from the edge + */ + unsigned int edgeLevel(Dart d); + + //! Return the level of the face of d in the current level map + /* @param d Dart from the face + */ + unsigned int faceLevel(Dart d); + + //! Return the level of the volume of d in the current level map + /* @param d Dart from the volume + */ + unsigned int volumeLevel(Dart d); + + //! Return the oldest dart of the face of d in the current level map + /* @param d Dart from the edge + */ + Dart faceOldestDart(Dart d); + + //! Return the level of the edge of d in the current level map + /* @param d Dart from the edge + */ + Dart volumeOldestDart(Dart d); + + //! Return true if the edge of d in the current level map + //! has already been subdivided to the next level + /*! @param d Dart from the edge + */ + bool edgeIsSubdivided(Dart d) ; + + //! Return true if the face of d in the current level map + //! has already been subdivided to the next level + /*! @param d Dart from the face + */ + bool faceIsSubdivided(Dart d) ; + + //! Return true if the volume of d in the current level map + //! has already been subdivided to the next level + /*! @param d Dart from the volume + */ + bool volumeIsSubdivided(Dart d); + //@} + + /*! @name Subdivision + * + *************************************************************************/ +protected: + //@{ + //! Add a new resolution level + /*! + */ + void addNewLevel() ; + + //! Subdivide the edge of d to the next level + /*! @param d Dart from the edge + */ + void subdivideEdge(Dart d) ; + + //! Subdivide the edge of d to the next level + /*! @param d Dart frome the face + */ + void subdivideFace(Dart d, SubdivideType sType = S_QUAD) ; + +public: + //! Subdivide the volume of d to hexahedral cells + /*! @param d Dart from the volume + */ + void subdivideVolume(Dart d) ; + + //! Subdivide the volume of d to hexahedral cells + /*! @param d Dart from the volume + */ + void subdivideVolumeTetOcta(Dart d) ; + //@} + + /*! @name Vertices Attributes management + * + *************************************************************************/ + //@{ + //! + /*! + */ + void setVertexVertexFunctor(FunctorType* f) { vertexVertexFunctor = f ; } + + //! + /* + * + */ + void setEdgeVertexFunctor(FunctorType* f) { edgeVertexFunctor = f ; } + + //! + /*! + */ + void setFaceVertexFunctor(FunctorType* f) { faceVertexFunctor = f ; } + + //! + /*! + */ + void setVolumeVertexFunctor(FunctorType* f) { volumeVertexFunctor = f ; } + //@} + +}; + +} //end namespace CGoGN + +#endif /* __MAP3MR_PRIMAL__ */ diff --git a/include/Topology/map/map3MR/map3MR_PrimalRegular.h b/include/Topology/map/map3MR/map3MR_PrimalRegular.h new file mode 100644 index 00000000..91056b64 --- /dev/null +++ b/include/Topology/map/map3MR/map3MR_PrimalRegular.h @@ -0,0 +1,159 @@ +/******************************************************************************* +* CGoGN: Combinatorial and Geometric modeling with Generic N-dimensional Maps * +* version 0.1 * +* Copyright (C) 2009-2012, IGG Team, LSIIT, University of Strasbourg * +* * +* This library is free software; you can redistribute it and/or modify it * +* under the terms of the GNU Lesser General Public License as published by the * +* Free Software Foundation; either version 2.1 of the License, or (at your * +* option) any later version. * +* * +* This library is distributed in the hope that it will be useful, but WITHOUT * +* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or * +* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License * +* for more details. * +* * +* You should have received a copy of the GNU Lesser General Public License * +* along with this library; if not, write to the Free Software Foundation, * +* Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. * +* * +* Web site: http://cgogn.unistra.fr/ * +* Contact information: cgogn@unistra.fr * +* * +*******************************************************************************/ + +#ifndef __MAP3MR_PRIMAL__ +#define __MAP3MR_PRIMAL__ + +#include "Topology/map/embeddedMap3.h" +#include "Topology/generic/traversorCell.h" +#include "Topology/generic/traversor3.h" + +#include "Topology/map/map3MR/filters_Primal.h" + +namespace CGoGN +{ + +/*! \brief The class of regular 3-map MR + */ + +class Map3MR_PrimalRegular : public EmbeddedMap3 +{ +protected: + bool shareVertexEmbeddings; + + std::vector synthesisFilters ; + std::vector analysisFilters ; + + +public: + Map3MR_PrimalRegular(); + + std::string mapTypeName() { return "Map3MR_PrimalRegular"; } + + /*! @name Topological helping functions + * + *************************************************************************/ + //@{ + //! + /* + * + */ + void swapEdges(Dart d, Dart e); + + //! + /*! + * + */ + bool isTetrahedron(Dart d); + + void splitFaceInSurface(Dart d, Dart e); + + Dart cutEdgeInSurface(Dart d); + + //! + /*! + * + */ + void saveRelationsAroundVertex(Dart d, std::vector >& vd); + + //! + /*! + * + */ + void unsewAroundVertex(std::vector >& vd); + + Dart quadranguleFace(Dart d); + + void splitSurfaceInVolume(std::vector& vd, bool firstSideClosed = true, bool secondSideClosed = false); + + //@} + + /*! @name Level creation + * + *************************************************************************/ + //@{ + //! + /* + * + */ + void addNewLevelTetraOcta(bool embedNewVertices); + + //! + /* + * + */ + void addNewLevel(bool embedNewVertices); + + //@} + + /*! @name Geometry modification + * + *************************************************************************/ + //@{ + //! + /* + * + */ + void addSynthesisFilter(Multiresolution::MRFilter* f) { synthesisFilters.push_back(f) ; } + + //! + /* + * + */ + void addAnalysisFilter(Multiresolution::MRFilter* f) { analysisFilters.push_back(f) ; } + + //! + /* + * + */ + void clearSynthesisFilters() { synthesisFilters.clear() ; } + + //! + /* + * + */ + void clearAnalysisFilters() { analysisFilters.clear() ; } + //@} + + /*! @name + * + *************************************************************************/ + //@{ + //! + /* + * + */ + void analysis() ; + + //! + /* + * + */ + void synthesis() ; + //@} +}; + +} //end namespace CGoGN + +#endif /* __MAP3MR_PRIMAL__ */ diff --git a/src/Topology/gmap/gmap2.cpp b/src/Topology/gmap/gmap2.cpp index 86b67e65..a11195a5 100644 --- a/src/Topology/gmap/gmap2.cpp +++ b/src/Topology/gmap/gmap2.cpp @@ -584,9 +584,24 @@ bool GMap2::mergeVolumes(Dart d, Dart e) return true ; } -void GMap2::splitCC(std::vector& vd) +void GMap2::splitSurface(std::vector& vd, bool firstSideClosed, bool secondSideClosed) { //assert(checkSimpleOrientedPath(vd)) ; + Dart e = vd.front() ; + Dart e2 = phi2(e) ; + + //unsew the edge path + for(std::vector::iterator it = vd.begin() ; it != vd.end() ; ++it) + { + if(!GMap2::isBoundaryEdge(*it)) + unsewFaces(*it) ; + } + + if(firstSideClosed) + GMap2::fillHole(e) ; + + if(secondSideClosed) + GMap2::fillHole(e2) ; } /*! @name Topological Queries diff --git a/src/Topology/gmap/gmap3.cpp b/src/Topology/gmap/gmap3.cpp index 0855727b..407d4ac7 100644 --- a/src/Topology/gmap/gmap3.cpp +++ b/src/Topology/gmap/gmap3.cpp @@ -410,15 +410,23 @@ void GMap3::splitVolume(std::vector& vd) Dart e = vd.front(); Dart e2 = phi2(e); - //unsew the edge path - for(std::vector::iterator it = vd.begin() ; it != vd.end() ; ++it) - GMap2::unsewFaces(*it); - - GMap2::fillHole(e) ; - GMap2::fillHole(e2) ; + GMap2::splitSurface(vd,true,true); //sew the two connected components - GMap3::sewVolumes(beta2(e), beta2(e2), false); + GMap3::sewVolumes(phi2(e), phi2(e2), false); + +// Dart e = vd.front(); +// Dart e2 = phi2(e); +// +// //unsew the edge path +// for(std::vector::iterator it = vd.begin() ; it != vd.end() ; ++it) +// GMap2::unsewFaces(*it); +// +// GMap2::fillHole(e) ; +// GMap2::fillHole(e2) ; +// +// //sew the two connected components +// GMap3::sewVolumes(beta2(e), beta2(e2), false); } /*! @name Topological Queries diff --git a/src/Topology/map/embeddedMap2.cpp b/src/Topology/map/embeddedMap2.cpp index b46b7ac5..85fda104 100644 --- a/src/Topology/map/embeddedMap2.cpp +++ b/src/Topology/map/embeddedMap2.cpp @@ -425,6 +425,51 @@ bool EmbeddedMap2::mergeVolumes(Dart d, Dart e) return false ; } +void EmbeddedMap2::splitSurface(std::vector& vd, bool firstSideClosed, bool secondSideClosed) +{ + std::vector darts ; + darts.reserve(vd.size()); + + // save the edge neighbors darts + for(std::vector::iterator it = vd.begin() ; it != vd.end() ; ++it) + { + darts.push_back(phi2(*it)); + } + + assert(darts.size() == vd.size()); + + Map2::splitSurface(vd, firstSideClosed, secondSideClosed); + + // follow the edge path a second time to embed the vertex, edge and volume orbits + for(unsigned int i = 0; i < vd.size(); ++i) + { + Dart dit = vd[i]; + Dart dit2 = darts[i]; + + // embed the vertex embedded from the origin volume to the new darts + if(isOrbitEmbedded(VERTEX)) + { + copyDartEmbedding(VERTEX, phi2(dit), phi1(dit)); + copyDartEmbedding(VERTEX, phi2(dit2), phi1(dit2)); + } + + // embed the edge embedded from the origin volume to the new darts + if(isOrbitEmbedded(EDGE)) + { + unsigned int eEmb = getEmbedding(EDGE, dit) ; + setDartEmbedding(EDGE, phi2(dit), eEmb); + setDartEmbedding(EDGE, phi2(dit2), eEmb); + } + + // embed the volume embedded from the origin volume to the new darts + if(isOrbitEmbedded(VOLUME)) + { + copyDartEmbedding(VOLUME, phi2(dit), dit); + copyDartEmbedding(VOLUME, phi2(dit2), dit2); + } + } +} + unsigned int EmbeddedMap2::closeHole(Dart d, bool forboundary) { unsigned int nbE = Map2::closeHole(d, forboundary) ; @@ -434,7 +479,7 @@ unsigned int EmbeddedMap2::closeHole(Dart d, bool forboundary) { if (isOrbitEmbedded(VERTEX)) { - copyDartEmbedding(VERTEX, f, alpha1(f)) ; + copyDartEmbedding(VERTEX, f, phi2(phi_1(f))) ; } if (isOrbitEmbedded(EDGE)) { diff --git a/src/Topology/map/embeddedMap3.cpp b/src/Topology/map/embeddedMap3.cpp index d94e527a..afef1a85 100644 --- a/src/Topology/map/embeddedMap3.cpp +++ b/src/Topology/map/embeddedMap3.cpp @@ -53,18 +53,15 @@ Dart EmbeddedMap3::cutEdge(Dart d) copyCell(EDGE, nd, d) ; } -// if(isOrbitEmbedded(ORIENTED_FACE)) if(isOrbitEmbedded(FACE2)) { Dart f = d; do { Dart f1 = phi1(f) ; -// copyDartEmbedding(ORIENTED_FACE, f1, f); + copyDartEmbedding(FACE2, f1, f); Dart e = phi3(f1); -// copyDartEmbedding(ORIENTED_FACE, phi1(e), e); -// copyDartEmbedding(FACE2, f1, f); copyDartEmbedding(FACE2, phi1(e), e); f = alpha2(f); } while(f != d); @@ -438,10 +435,14 @@ bool EmbeddedMap3::check() { if(isOrbitEmbedded(VERTEX)) { - if( getEmbedding(VERTEX, d) != getEmbedding(VERTEX, alpha1(d)) || - getEmbedding(VERTEX, d) != getEmbedding(VERTEX, alpha2(d)) ) + if( getEmbedding(VERTEX, d) != getEmbedding(VERTEX, alpha1(d))) + { + std::cout << "Embedding Check : different embeddings on vertex (alpha1(d) != d)" << std::endl ; + return false ; + } + if(getEmbedding(VERTEX, d) != getEmbedding(VERTEX, alpha2(d)) ) { - std::cout << "Embedding Check : different embeddings on vertex" << std::endl ; + std::cout << "Embedding Check : different embeddings on vertex (alpha2(d) != d)" << std::endl ; return false ; } } diff --git a/src/Topology/map/map2.cpp b/src/Topology/map/map2.cpp index 07dc57d4..11131528 100644 --- a/src/Topology/map/map2.cpp +++ b/src/Topology/map/map2.cpp @@ -329,8 +329,8 @@ void Map2::swapEdges(Dart d, Dart e) { assert(!Map2::isBoundaryEdge(d) && !Map2::isBoundaryEdge(e)); - Dart d2 = phi2(d); - Dart e2 = phi2(e); + //Dart d2 = phi2(d); + //Dart e2 = phi2(e); phi2unsew(d); phi2unsew(e) ; @@ -549,7 +549,10 @@ void Map2::splitSurface(std::vector& vd, bool firstSideClosed, bool second //unsew the edge path for(std::vector::iterator it = vd.begin() ; it != vd.end() ; ++it) - unsewFaces(*it) ; + { + if(!Map2::isBoundaryEdge(*it)) + unsewFaces(*it) ; + } if(firstSideClosed) Map2::fillHole(e) ; @@ -569,7 +572,7 @@ bool Map2::sameOrientedVertex(Dart d, Dart e) { if (it == e) // Test equality with e return true; - it = alpha1(it); + it = phi2(phi_1(it)); } while (it != d); return false; // None is equal to e => vertices are distinct } @@ -581,7 +584,7 @@ unsigned int Map2::vertexDegree(Dart d) do { ++count ; - it = alpha1(it) ; + it = phi2(phi_1(it)) ; } while (it != d) ; return count ; } @@ -593,7 +596,7 @@ bool Map2::isBoundaryVertex(Dart d) { if (isBoundaryMarked(it)) return true ; - it = alpha1(it) ; + it = phi2(phi_1(it)) ; } while (it != d) ; return false ; } @@ -605,7 +608,7 @@ Dart Map2::findBoundaryEdgeOfVertex(Dart d) { if (isBoundaryMarked(it)) return it ; - it = alpha1(it) ; + it = phi2(phi_1(it)) ; } while (it != d) ; return NIL ; } @@ -781,7 +784,7 @@ bool Map2::foreach_dart_of_vertex(Dart d, FunctorType& f, unsigned int thread) { if (f(dNext)) return true; - dNext = alpha1(dNext); + dNext = phi2(phi_1(dNext)); } while (dNext != d); return false; } diff --git a/src/Topology/map/map3.cpp b/src/Topology/map/map3.cpp index cd7f71fb..7c10f53c 100644 --- a/src/Topology/map/map3.cpp +++ b/src/Topology/map/map3.cpp @@ -511,13 +511,6 @@ void Map3::splitVolume(std::vector& vd) Dart e = vd.front(); Dart e2 = phi2(e); -// //unsew the edge path -// for(std::vector::iterator it = vd.begin() ; it != vd.end() ; ++it) -// Map2::unsewFaces(*it); -// -// Map2::fillHole(e) ; -// Map2::fillHole(e2) ; - Map2::splitSurface(vd,true,true); //sew the two connected components @@ -603,6 +596,13 @@ unsigned int Map3::vertexDegree(Dart d) return count; } +unsigned int Map3::vertexDegreeOnBoundary(Dart d) +{ + assert(Map3::isBoundaryVertex(d)); + + return Map2::vertexDegree(d); +} + bool Map3::isBoundaryVertex(Dart d) { DartMarkerStore mv(*this); // Lock a marker @@ -756,7 +756,7 @@ bool Map3::check() if(phi1(d3) != phi3(phi_1(d))) { std::cout << "Check: phi3 , faces are not entirely sewn" << std::endl; - return false; + //return false; } Dart d2 = phi2(d); diff --git a/src/Topology/map/map3MR/map3MR_PrimalAdapt.cpp b/src/Topology/map/map3MR/map3MR_PrimalAdapt.cpp new file mode 100644 index 00000000..44939860 --- /dev/null +++ b/src/Topology/map/map3MR/map3MR_PrimalAdapt.cpp @@ -0,0 +1,665 @@ +/******************************************************************************* +* CGoGN: Combinatorial and Geometric modeling with Generic N-dimensional Maps * +* version 0.1 * +* Copyright (C) 2009-2012, IGG Team, LSIIT, University of Strasbourg * +* * +* This library is free software; you can redistribute it and/or modify it * +* under the terms of the GNU Lesser General Public License as published by the * +* Free Software Foundation; either version 2.1 of the License, or (at your * +* option) any later version. * +* * +* This library is distributed in the hope that it will be useful, but WITHOUT * +* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or * +* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License * +* for more details. * +* * +* You should have received a copy of the GNU Lesser General Public License * +* along with this library; if not, write to the Free Software Foundation, * +* Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. * +* * +* Web site: http://cgogn.unistra.fr/ * +* Contact information: cgogn@unistra.fr * +* * +*******************************************************************************/ + +#include "Topology/map/map3MR/map3MR_PrimalAdapt.h" + +namespace CGoGN +{ + +Map3MR_PrimalAdapt::Map3MR_PrimalAdapt() : + shareVertexEmbeddings(false), + vertexVertexFunctor(NULL), + edgeVertexFunctor(NULL), + faceVertexFunctor(NULL), + volumeVertexFunctor(NULL) +{ + initMR(); +} + +/*! @name Topological helping functions + *************************************************************************/ +void Map3MR_PrimalAdapt::swapEdges(Dart d, Dart e) +{ + if(!Map2::isBoundaryEdge(d) && !Map2::isBoundaryEdge(e)) + { + Dart d2 = phi2(d); + Dart e2 = phi2(e); + + Map2::swapEdges(d,e); + +// Map2::unsewFaces(d); +// Map2::unsewFaces(e); +// +// Map2::sewFaces(d, e); +// Map2::sewFaces(d2, e2); + + if(isOrbitEmbedded(VERTEX)) + { + copyDartEmbedding(VERTEX, d, phi2(phi_1(d))); + copyDartEmbedding(VERTEX, e, phi2(phi_1(e))); + copyDartEmbedding(VERTEX, d2, phi2(phi_1(d2))); + copyDartEmbedding(VERTEX, e2, phi2(phi_1(e2))); + } + + if(isOrbitEmbedded(EDGE)) + { + + } + + if(isOrbitEmbedded(VOLUME)) + embedNewCell(VOLUME, d); + } +} + +bool Map3MR_PrimalAdapt::isTetrahedron(Dart d) +{ + unsigned int nbFaces = 0; + + //Test the number of faces end its valency + Traversor3WF travWF(*this, d); + for(Dart dit = travWF.begin() ; dit != travWF.end(); dit = travWF.next()) + { + //increase the number of faces + nbFaces++; + if(nbFaces > 4) //too much faces + return false; + + //test the valency of this face + if(faceDegree(dit) != 3) + return false; + } + + return true; +} + +void Map3MR_PrimalAdapt::saveRelationsAroundVertex(Dart d, std::vector >& vd) +{ + assert(getDartLevel(d) <= getCurrentLevel() || !"saveRelationsAroundVertex : called with a dart inserted after current level") ; + + //le brin est forcement du niveau cur + Dart dit = d; + + do + { + vd.push_back(std::pair(dit,phi2(dit))); + + dit = phi2(phi_1(dit)); + + }while(dit != d); +} + +void Map3MR_PrimalAdapt::unsewAroundVertex(std::vector >& vd) +{ + //unsew the edge path + for(std::vector >::iterator it = vd.begin() ; it != vd.end() ; ++it) + { + Dart dit = (*it).first; + Dart dit2 = (*it).second; + + Map2::unsewFaces(dit); + + if(isOrbitEmbedded(VERTEX)) + { + copyDartEmbedding(VERTEX, phi2(dit2), dit); + copyDartEmbedding(VERTEX, phi2(dit), dit2); + } + + if(isOrbitEmbedded(EDGE)) + { + + } + } +} + +/*! @name Cells informations + *************************************************************************/ +unsigned int Map3MR_PrimalAdapt::edgeLevel(Dart d) +{ + assert(getDartLevel(d) <= getCurrentLevel() || !"edgeLevel : called with a dart inserted after current level") ; + + unsigned int r = 0; + Dart dit = d; + do + { + unsigned int ld = getDartLevel(dit) ; + unsigned int ldd = getDartLevel(phi2(dit)) ; + unsigned int min = ld < ldd ? ldd : ld; + + r = r < min ? min : r ; + + dit = alpha2(dit); + } while(dit != d); + + return r; +} + +unsigned int Map3MR_PrimalAdapt::faceLevel(Dart d) +{ + assert(getDartLevel(d) <= getCurrentLevel() || !"faceLevel : called with a dart inserted after current level") ; + + if(getCurrentLevel() == 0) + return 0 ; + + Dart it = d ; + unsigned int min1 = getDartLevel(it) ; // the level of a face is the second minimum of the + it = phi1(it) ; + unsigned int min2 = getDartLevel(it) ; // insertion levels of its darts + + if(min2 < min1) + { + unsigned int tmp = min1 ; + min1 = min2 ; + min2 = tmp ; + } + + it = phi1(it) ; + while(it != d) + { + unsigned int dl = getDartLevel(it) ; + if(dl < min2) + { + if(dl < min1) + { + min2 = min1 ; + min1 = dl ; + } + else + min2 = dl ; + } + it = phi1(it) ; + } + + return min2 ; +} + +unsigned int Map3MR_PrimalAdapt::volumeLevel(Dart d) +{ + assert(getDartLevel(d) <= getCurrentLevel() || !"volumeLevel : called with a dart inserted after current level") ; + + if(getCurrentLevel() == 0) + return 0 ; + + Dart oldest = d ; + unsigned int vLevel = std::numeric_limits::max(); //hook sioux + + //First : the level of a volume is the minimum of the levels of its faces + Traversor3WF travF(*this, d); + for (Dart dit = travF.begin(); dit != travF.end(); dit = travF.next()) + { + // in a first time, the level of a face + //the level of the volume is the minimum of the + //levels of its faces + unsigned int fLevel = faceLevel(dit); + vLevel = fLevel < vLevel ? fLevel : vLevel ; + } + + //Second : the case of all faces regularly subdivided but not the volume itself + + return vLevel; +} + +Dart Map3MR_PrimalAdapt::faceOldestDart(Dart d) +{ + assert(getDartLevel(d) <= getCurrentLevel() || !"faceOldestDart : called with a dart inserted after current level") ; + + Dart it = d ; + Dart oldest = it ; + unsigned int l_old = getDartLevel(oldest) ; + do + { + unsigned int l = getDartLevel(it) ; + if(l < l_old || (l == l_old && it < oldest)) + { + oldest = it ; + l_old = l ; + } + it = phi1(it) ; + } while(it != d) ; + return oldest ; +} + +Dart Map3MR_PrimalAdapt::volumeOldestDart(Dart d) +{ + assert(getDartLevel(d) <= getCurrentLevel() || !"volumeOldestDart : called with a dart inserted after current level") ; + + Dart oldest = d; + Traversor3WF travF(*this, d); + for (Dart dit = travF.begin(); dit != travF.end(); dit = travF.next()) + { + //for every dart in this face + Dart old = faceOldestDart(dit); + if(getDartLevel(old) < getDartLevel(oldest)) + oldest = old; + } + + return oldest; +} + +bool Map3MR_PrimalAdapt::edgeIsSubdivided(Dart d) +{ + assert(getDartLevel(d) <= getCurrentLevel() || !"edgeIsSubdivided : called with a dart inserted after current level") ; + + if(getCurrentLevel() == getMaxLevel()) + return false ; + + Dart d2 = phi2(d) ; + incCurrentLevel() ; + Dart d2_l = phi2(d) ; + decCurrentLevel() ; ; + if(d2 != d2_l) + return true ; + else + return false ; +} + +bool Map3MR_PrimalAdapt::faceIsSubdivided(Dart d) +{ + assert(getDartLevel(d) <= getCurrentLevel() || !"faceIsSubdivided : called with a dart inserted after current level") ; + + if(getCurrentLevel() == getMaxLevel()) + return false ; + + // a face whose level in the current level map is lower than + // the current level can not be subdivided to higher levels + unsigned int fLevel = faceLevel(d) ; + if(fLevel < getCurrentLevel()) + return false ; + + bool subd = false ; + incCurrentLevel() ; + if(getDartLevel(phi1(phi1(d))) == getCurrentLevel()) + subd = true ; + decCurrentLevel() ; + + return subd ; +} + +bool Map3MR_PrimalAdapt::volumeIsSubdivided(Dart d) +{ + assert(getDartLevel(d) <= getCurrentLevel() || !"volumeIsSubdivided : called with a dart inserted after current level") ; + + unsigned int vLevel = volumeLevel(d); + if(vLevel <= getCurrentLevel()) + return false; + + //Test if all faces are subdivided + bool facesAreSubdivided = faceIsSubdivided(d) ; + Traversor3WF travF(*this, d); + for (Dart dit = travF.begin(); dit != travF.end(); dit = travF.next()) + { + facesAreSubdivided &= faceIsSubdivided(dit) ; + } + + //But not the volume itself + bool subd = false; + incCurrentLevel(); + if(facesAreSubdivided && getDartLevel(phi2(phi1(phi1(d)))) == getCurrentLevel()) + subd = true; + decCurrentLevel() ; + return subd; +} + +/*! @name Subdivision + *************************************************************************/ +void Map3MR_PrimalAdapt::addNewLevel() +{ + addLevel(); +} + + +void Map3MR_PrimalAdapt::subdivideEdge(Dart d) +{ + assert(getDartLevel(d) <= getCurrentLevel() || !"subdivideEdge : called with a dart inserted after current level") ; + assert(!edgeIsSubdivided(d) || !"Trying to subdivide an already subdivided edge") ; + + incCurrentLevel() ; + + //Duplicate all darts around the edge + Dart dit = d; + do + { + duplicateDart(dit); + + Dart dd = phi2(dit) ; + duplicateDart(dd) ; + + Dart d1 = phi1(dit) ; + duplicateDart(d1) ; + + Dart dd1 = phi1(dd) ; + duplicateDart(dd1) ; + + dit = alpha2(dit); + }while(dit != d); + + cutEdge(d) ; + (*edgeVertexFunctor)(phi1(d)); + + decCurrentLevel() ; +} + +void Map3MR_PrimalAdapt::subdivideFace(Dart d, SubdivideType sType) +{ + assert(getDartLevel(d) <= getCurrentLevel() || !"subdivideFace : called with a dart inserted after current level") ; + assert(!faceIsSubdivided(d) || !"Trying to subdivide an already subdivided face") ; + + unsigned int fLevel = faceLevel(d) ; + Dart old = faceOldestDart(d) ; + + pushLevel() ; + setCurrentLevel(fLevel) ; // go to the level of the face to subdivide its edges + + unsigned int degree = 0 ; + Dart it = old ; + do + { + ++degree ; // compute the degree of the face + + if(!edgeIsSubdivided(it)) + subdivideEdge(it) ; // and cut the edges (if they are not already) + it = phi1(it) ; + } while(it != old) ; + + setCurrentLevel(fLevel + 1) ; // go to the next level to perform face subdivision + + if(degree == 3 && sType == S_TRI) // if subdividing a triangle + { + Dart dd = phi1(old) ; + Dart e = phi1(dd) ; + (*vertexVertexFunctor)(e) ; + e = phi1(e) ; + splitFace(dd, e) ; + + dd = e ; + e = phi1(dd) ; + (*vertexVertexFunctor)(e) ; + e = phi1(e) ; + splitFace(dd, e) ; + + dd = e ; + e = phi1(dd) ; + (*vertexVertexFunctor)(e) ; + e = phi1(e) ; + splitFace(dd, e) ; + } + else // if subdividing a polygonal face + { + Dart dd = phi1(old) ; + Dart next = phi1(dd) ; + (*vertexVertexFunctor)(next) ; + next = phi1(next) ; + splitFace(dd, next) ; // insert a first edge + Dart ne = phi2(phi_1(dd)); + + cutEdge(ne) ; // cut the new edge to insert the central vertex + + dd = phi1(next) ; + (*vertexVertexFunctor)(dd) ; + dd = phi1(dd) ; + while(dd != ne) // turn around the face and insert new edges + { // linked to the central vertex + splitFace(phi1(ne), dd) ; + dd = phi1(dd) ; + (*vertexVertexFunctor)(dd) ; + dd = phi1(dd) ; + } + + (*faceVertexFunctor)(phi2(ne)) ; + } + + popLevel() ; +} + +void Map3MR_PrimalAdapt::subdivideVolume(Dart d) +{ + assert(getDartLevel(d) <= getCurrentLevel() || !"subdivideVolume : called with a dart inserted after current level") ; + assert(!volumeIsSubdivided(d) || !"Trying to subdivide an already subdivided face") ; + + unsigned int vLevel = volumeLevel(d); + Dart old = volumeOldestDart(d); + + Traversor3WV traV(*this, old); + + pushLevel() ; + setCurrentLevel(vLevel) ; // go to the level of the face to subdivide its edges + + if(getCurrentLevel() == getMaxLevel()) + addNewLevel() ; + + // + // Subdivide Faces + // + std::vector > subdividedFaces; + subdividedFaces.reserve(128); + + Traversor3WF traF(*this, old); + for(Dart dit = traF.begin(); dit != traF.end(); dit = traF.next()) + { + //if needed subdivide face + if(!faceIsSubdivided(dit)) + subdivideFace(dit); + + //save darts from the central vertex of each subdivided face + incCurrentLevel() ; + saveRelationsAroundVertex(phi2(phi1(dit)), subdividedFaces); + decCurrentLevel() ; + } + + incCurrentLevel() ; + unsewAroundVertex(subdividedFaces); + decCurrentLevel() ; + + // + // Create inside volumes + // + Dart centralDart = NIL; + setCurrentLevel(vLevel + 1) ; // go to the next level to perform volume subdivision + for(Dart dit = traV.begin(); dit != traV.end(); dit = traV.next()) + { + Map2::fillHole(phi1(dit)); + + Dart old = phi2(phi1(dit)); + Dart bc = newBoundaryCycle(faceDegree(old)); + sewVolumes(old, bc, false); + + if (isOrbitEmbedded(VERTEX)) + { + Dart it = bc; + do + { + copyDartEmbedding(VERTEX, it, phi1(phi3(it))); + it = phi1(it) ; + } while(it != bc) ; + } + + Dart dd = phi1(phi1(old)) ; + splitFace(old,dd) ; + + Dart ne = phi1(phi1(old)) ; + + cutEdge(ne); + centralDart = phi1(ne); + + Dart stop = phi2(phi1(ne)); + ne = phi2(ne); + do + { + dd = phi1(phi1(phi1(ne))); + + splitFace(ne, dd) ; + + ne = phi2(phi_1(ne)); + dd = phi1(phi1(dd)); + } + while(dd != stop); + } + + // + // Sew inside volumes + // + for (std::vector >::iterator it = subdividedFaces.begin(); it != subdividedFaces.end(); ++it) + { + Dart f1 = phi2((*it).first); + Dart f2 = phi2((*it).second); + + if(isBoundaryFace(f1) && isBoundaryFace(f2)) + { + std::cout << "plop" << std::endl; + sewVolumes(f1, f2, false); + } + } + + (*volumeVertexFunctor)(centralDart) ; + popLevel() ; +} + +void Map3MR_PrimalAdapt::subdivideVolumeTetOcta(Dart d) +{ + assert(getDartLevel(d) <= getCurrentLevel() || !"subdivideVolumeTetOcta : called with a dart inserted after current level") ; + assert(!volumeIsSubdivided(d) || !"Trying to subdivide an already subdivided face") ; + + unsigned int vLevel = volumeLevel(d); + Dart old = volumeOldestDart(d); + + pushLevel() ; + setCurrentLevel(vLevel) ; // go to the level of the face to subdivide its edges + + if(getCurrentLevel() == getMaxLevel()) + addNewLevel() ; + + // + // Subdivide Faces + // + std::vector > subdividedFaces; + subdividedFaces.reserve(128); + + Traversor3WF traF(*this, old); + for(Dart dit = traF.begin(); dit != traF.end(); dit = traF.next()) + { + //if needed subdivide face + if(!faceIsSubdivided(dit)) + subdivideFace(dit, S_TRI); + } + + // + // Create inside volumes + // + Dart centralDart = NIL; + + bool isNotTet = false; + Traversor3WV traV(*this, old); + setCurrentLevel(getMaxLevel()) ; + for(Dart dit = traV.begin(); dit != traV.end(); dit = traV.next()) + { +// //check if neighboring volumes have to be subdivided +// Traversor3VW traVN(*this, dit); +// for(Dart ditN = traVN.begin(); ditN != traVN.end(); ditN = traVN.next()) +// { +// if(volumeLevel(ditN) == vLevel -1) +// std::cout << "to subdivide" << std::endl;//subdivideVolumeTetOcta(ditN); +// } + + Dart f1 = phi1(dit); + Dart e = dit; + std::vector v ; + + do + { + v.push_back(phi1(e)); + e = phi2(phi_1(e)); + } + while(e != dit); + + splitVolume(v) ; + + //if is not a tetrahedron + unsigned int fdeg = faceDegree(phi2(f1)); + if(fdeg > 3) + { + isNotTet = true; + Dart old = phi2(phi1(dit)); + Dart dd = phi1(old) ; + splitFace(old,dd) ; + + Dart ne = phi1(old); + + cutEdge(ne); + centralDart = phi1(ne); + (*volumeVertexFunctor)(centralDart) ; + + Dart stop = phi2(phi1(ne)); + ne = phi2(ne); + do + { + dd = phi1(phi1(ne)); + + splitFace(ne, dd) ; + + ne = phi2(phi_1(ne)); + dd = phi1(dd); + } + while(dd != stop); + } + } + + //switch inner faces + if(isNotTet) + { + + for(Dart dit = traV.begin(); dit != traV.end(); dit = traV.next()) + { + Dart x = phi_1(phi2(phi1(dit))); + Dart f = x; + + DartMarkerStore me(*this); + + do + { + Dart f3 = phi3(f); + + if(!me.isMarked(f3)) + { + Dart tmp = phi_1(phi2(phi_1(phi2(phi_1(f3))))); //future voisin par phi2 + centralDart = tmp; + + Dart f32 = phi2(f3); + swapEdges(f3, tmp); + + me.markOrbit(EDGE, f3); + me.markOrbit(EDGE, f32); + } + + f = phi2(phi_1(f)); + }while(f != x); + } + + embedOrbit(VERTEX, centralDart, getEmbedding(VERTEX, centralDart)); + + //(*volumeVertexFunctor)(centralDart) ; + } + + popLevel(); +} + +} //end namespace CGoGN diff --git a/src/Topology/map/map3MR/map3MR_PrimalRegular.cpp b/src/Topology/map/map3MR/map3MR_PrimalRegular.cpp new file mode 100644 index 00000000..e999c38b --- /dev/null +++ b/src/Topology/map/map3MR/map3MR_PrimalRegular.cpp @@ -0,0 +1,671 @@ +/******************************************************************************* +* CGoGN: Combinatorial and Geometric modeling with Generic N-dimensional Maps * +* version 0.1 * +* Copyright (C) 2009-2012, IGG Team, LSIIT, University of Strasbourg * +* * +* This library is free software; you can redistribute it and/or modify it * +* under the terms of the GNU Lesser General Public License as published by the * +* Free Software Foundation; either version 2.1 of the License, or (at your * +* option) any later version. * +* * +* This library is distributed in the hope that it will be useful, but WITHOUT * +* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or * +* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License * +* for more details. * +* * +* You should have received a copy of the GNU Lesser General Public License * +* along with this library; if not, write to the Free Software Foundation, * +* Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. * +* * +* Web site: http://cgogn.unistra.fr/ * +* Contact information: cgogn@unistra.fr * +* * +*******************************************************************************/ + +#include "Topology/map/map3MR/map3MR_PrimalRegular.h" + +namespace CGoGN +{ + +Map3MR_PrimalRegular::Map3MR_PrimalRegular() : + shareVertexEmbeddings(true) +{ + initMR(); +} + +/*! @name Topological helping functions + *************************************************************************/ +void Map3MR_PrimalRegular::swapEdges(Dart d, Dart e) +{ + if(!Map2::isBoundaryEdge(d) && !Map2::isBoundaryEdge(e)) + { + Dart d2 = phi2(d); + Dart e2 = phi2(e); + + Map2::swapEdges(d,e); + + if(isOrbitEmbedded(VERTEX)) + { + copyDartEmbedding(VERTEX, d, phi2(phi_1(d))); + copyDartEmbedding(VERTEX, e, phi2(phi_1(e))); + copyDartEmbedding(VERTEX, d2, phi2(phi_1(d2))); + copyDartEmbedding(VERTEX, e2, phi2(phi_1(e2))); + } + + if(isOrbitEmbedded(EDGE)) + { + + } + + if(isOrbitEmbedded(VOLUME)) + embedNewCell(VOLUME, d); + } +} + +bool Map3MR_PrimalRegular::isTetrahedron(Dart d) +{ + unsigned int nbFaces = 0; + + //Test the number of faces end its valency + Traversor3WF travWF(*this, d); + for(Dart dit = travWF.begin() ; dit != travWF.end(); dit = travWF.next()) + { + //increase the number of faces + nbFaces++; + if(nbFaces > 4) //too much faces + return false; + + //test the valency of this face + if(faceDegree(dit) != 3) + return false; + } + + return true; +} + + + +void Map3MR_PrimalRegular::splitSurfaceInVolume(std::vector& vd, bool firstSideClosed, bool secondSideClosed) +{ + std::vector vd2 ; + vd2.reserve(vd.size()); + + // save the edge neighbors darts + for(std::vector::iterator it = vd.begin() ; it != vd.end() ; ++it) + { + vd2.push_back(phi2(*it)); + } + + assert(vd2.size() == vd.size()); + + Map2::splitSurface(vd, firstSideClosed, secondSideClosed); + + // follow the edge path a second time to embed the vertex, edge and volume orbits + for(unsigned int i = 0; i < vd.size(); ++i) + { + Dart dit = vd[i]; + Dart dit2 = vd2[i]; + + // embed the vertex embedded from the origin volume to the new darts + if(isOrbitEmbedded(VERTEX)) + { + copyDartEmbedding(VERTEX, phi2(dit), phi1(dit)); + copyDartEmbedding(VERTEX, phi2(dit2), phi1(dit2)); + } + } + +} + +void Map3MR_PrimalRegular::addNewLevelTetraOcta(bool embedNewVertices) +{ + pushLevel(); + + addLevel(); + setCurrentLevel(getMaxLevel()); + + //create the new level with the old one + for(unsigned int i = m_mrattribs.begin(); i != m_mrattribs.end(); m_mrattribs.next(i)) + { + unsigned int newindex = copyDartLine((*m_mrDarts[m_mrCurrentLevel])[i]) ; // duplicate all darts + (*m_mrDarts[m_mrCurrentLevel])[i] = newindex ; // on the new max level + if(!shareVertexEmbeddings) + (*m_embeddings[VERTEX])[newindex] = EMBNULL ; // set vertex embedding to EMBNULL if no sharing + } + + //subdivision + //1. cut edges + TraversorE travE(*this); + for (Dart d = travE.begin(); d != travE.end(); d = travE.next()) + { + if(!shareVertexEmbeddings) + { + if(getEmbedding(VERTEX, d) == EMBNULL) + embedNewCell(VERTEX, d) ; + if(getEmbedding(VERTEX, phi1(d)) == EMBNULL) + embedNewCell(VERTEX, d) ; + } + + cutEdge(d) ; + travE.mark(d) ; + travE.mark(phi1(d)) ; + +// When importing MR files +// if(embedNewVertices) +// embedNewCell(VERTEX, phi1(d)) ; + } + + //2. split faces - triangular faces + TraversorF travF(*this) ; + for (Dart d = travF.begin(); d != travF.end(); d = travF.next()) + { + Dart old = d ; + if(getDartLevel(old) == getMaxLevel()) + old = phi1(old) ; + + Dart dd = phi1(old) ; + Dart e = phi1(phi1(dd)) ; + splitFace(dd, e) ; + travF.mark(dd) ; + + dd = e ; + e = phi1(phi1(dd)) ; + splitFace(dd, e) ; + travF.mark(dd) ; + + dd = e ; + e = phi1(phi1(dd)) ; + splitFace(dd, e) ; + travF.mark(dd) ; + + travF.mark(e) ; + } + + //3. create inside volumes + setCurrentLevel(getMaxLevel() - 1) ; + TraversorW traW(*this); + for(Dart dit = traW.begin(); dit != traW.end(); dit = traW.next()) + { + Traversor3WV traWV(*this, dit); + + for(Dart ditWV = traWV.begin(); ditWV != traWV.end(); ditWV = traWV.next()) + { + setCurrentLevel(getMaxLevel()) ; + Dart f1 = phi1(ditWV); + Dart e = ditWV; + std::vector v ; + + do + { + v.push_back(phi1(e)); + e = phi2(phi_1(e)); + } + while(e != ditWV); + + splitVolume(v) ; + + //if is not a tetrahedron + unsigned int fdeg = faceDegree(phi2(f1)); + if(fdeg > 3) + { + Dart old = phi2(phi1(ditWV)); + Dart dd = phi1(old) ; + splitFace(old,dd) ; + + Dart ne = phi1(old); + + cutEdge(ne); + + Dart stop = phi2(phi1(ne)); + ne = phi2(ne); + do + { + dd = phi1(phi1(ne)); + + splitFace(ne, dd) ; + + ne = phi2(phi_1(ne)); + dd = phi1(dd); + } + while(dd != stop); + } + + setCurrentLevel(getMaxLevel() - 1) ; + } + } + + setCurrentLevel(getMaxLevel() - 1) ; + for(Dart dit = traW.begin(); dit != traW.end(); dit = traW.next()) + { + Traversor3WV traWV(*this, dit); + + for(Dart ditWV = traWV.begin(); ditWV != traWV.end(); ditWV = traWV.next()) + { + setCurrentLevel(getMaxLevel()) ; + Dart x = phi_1(phi2(phi1(ditWV))); + + if(!isTetrahedron(x)) + { + DartMarkerStore me(*this); + + Dart f = x; + + do + { + Dart f3 = phi3(f); + + if(!me.isMarked(f3)) + { + Dart tmp = phi_1(phi2(phi_1(phi2(phi_1(f3))))); //future voisin par phi2 + + Dart f32 = phi2(f3); + swapEdges(f3, tmp); + + me.markOrbit(EDGE, f3); + me.markOrbit(EDGE, f32); + } + + f = phi2(phi_1(f)); + }while(f != x); + + // When importing MR files + //if(embedNewVertices) + // embedNewCell(VERTEX, x) ; + } + setCurrentLevel(getMaxLevel() - 1) ; + } + } + + popLevel() ; +} + +void Map3MR_PrimalRegular::addNewLevel(bool embedNewVertices) +{ + pushLevel(); + + addLevel(); + setCurrentLevel(getMaxLevel()); + + //create the new level with the old one + for(unsigned int i = m_mrattribs.begin(); i != m_mrattribs.end(); m_mrattribs.next(i)) + { + unsigned int newindex = copyDartLine((*m_mrDarts[m_mrCurrentLevel])[i]) ; // duplicate all darts + (*m_mrDarts[m_mrCurrentLevel])[i] = newindex ; // on the new max level + if(!shareVertexEmbeddings) + (*m_embeddings[VERTEX])[newindex] = EMBNULL ; // set vertex embedding to EMBNULL if no sharing + } + + //subdivision + //1. cut edges + TraversorE travE(*this); + for (Dart d = travE.begin(); d != travE.end(); d = travE.next()) + { + if(!shareVertexEmbeddings) + { + if(getEmbedding(VERTEX, d) == EMBNULL) + embedNewCell(VERTEX, d) ; + if(getEmbedding(VERTEX, phi1(d)) == EMBNULL) + embedNewCell(VERTEX, d) ; + } + + cutEdge(d) ; + travE.mark(d) ; + travE.mark(phi1(d)) ; + +// When importing MR files : activated for DEBUG + if(embedNewVertices) + embedNewCell(VERTEX, phi1(d)) ; + } + + //2. split faces - quadrangule faces + TraversorF travF(*this) ; + for (Dart d = travF.begin(); d != travF.end(); d = travF.next()) + { + Dart old = d; + if(getDartLevel(old) == getMaxLevel()) + old = phi1(old) ; + + Dart dd = phi1(old) ; + Dart next = phi1(phi1(dd)) ; + splitFace(dd, next) ; // insert a first edge + + Dart ne = phi2(phi_1(dd)) ; + cutEdge(ne) ; // cut the new edge to insert the central vertex + travF.mark(dd) ; + +// When importing MR files : activated for DEBUG + if(embedNewVertices) + embedNewCell(VERTEX, phi1(ne)) ; + + dd = phi1(phi1(next)) ; + while(dd != ne) // turn around the face and insert new edges + { // linked to the central vertex + Dart tmp = phi1(ne) ; + splitFace(tmp, dd) ; + travF.mark(tmp) ; + dd = phi1(phi1(dd)) ; + } + travF.mark(ne) ; + } + + //3. create inside volumes + setCurrentLevel(getMaxLevel() - 1) ; + TraversorW traW(*this); + for(Dart dit = traW.begin(); dit != traW.end(); dit = traW.next()) + { + std::vector > subdividedFaces; + subdividedFaces.reserve(64); + Dart centralDart = NIL; + + Traversor3WV traWV(*this, dit); + for(Dart ditWV = traWV.begin(); ditWV != traWV.end(); ditWV = traWV.next()) + { + setCurrentLevel(getMaxLevel()) ; + + Dart e = ditWV; + std::vector v ; + + do + { + v.push_back(phi1(e)); + v.push_back(phi1(phi1(e))); + + if(!Map2::isBoundaryEdge(phi1(e))) + subdividedFaces.push_back(std::pair(phi1(e),phi2(phi1(e)))); + + if(!Map2::isBoundaryEdge(phi1(phi1(e)))) + subdividedFaces.push_back(std::pair(phi1(phi1(e)),phi2(phi1(phi1(e))))); + + e = phi2(phi_1(e)); + } + while(e != ditWV); + + splitSurfaceInVolume(v); + + Dart dd = phi2(phi1(ditWV)); + Dart next = phi1(phi1(dd)) ; + Map2::splitFace(dd, next) ; + + Dart ne = phi2(phi_1(dd)); + Map2::cutEdge(ne) ; + centralDart = phi1(ne); + + dd = phi1(phi1(next)) ; + while(dd != ne) + { + Dart tmp = phi1(ne) ; + Map2::splitFace(tmp, dd) ; + dd = phi1(phi1(dd)) ; + } + + setCurrentLevel(getMaxLevel() - 1) ; + } + + setCurrentLevel(getMaxLevel()) ; + //4 couture des relations precedemment sauvegarde + for (std::vector >::iterator it = subdividedFaces.begin(); it != subdividedFaces.end(); ++it) + { + Dart f1 = phi2((*it).first); + Dart f2 = phi2((*it).second); + + //if(isBoundaryFace(f1) && isBoundaryFace(f2)) + if(phi3(f1) == f1 && phi3(f2) == f2) + sewVolumes(f1, f2, false); + } + embedOrbit(VERTEX, centralDart, getEmbedding(VERTEX, centralDart)); + + setCurrentLevel(getMaxLevel() - 1) ; + } + + //A optimiser + setCurrentLevel(getMaxLevel()-1) ; + TraversorE travE2(*this); + for (Dart d = travE2.begin(); d != travE2.end(); d = travE2.next()) + { + setCurrentLevel(getMaxLevel()) ; + setCurrentLevel(getMaxLevel()-1) ; + embedOrbit(VERTEX, phi1(d), getEmbedding(VERTEX, phi1(d))); + } + setCurrentLevel(getMaxLevel()) ; + + setCurrentLevel(getMaxLevel()-1) ; + TraversorF travF2(*this) ; + for (Dart d = travF2.begin(); d != travF2.end(); d = travF2.next()) + { + setCurrentLevel(getMaxLevel()) ; + embedOrbit(VERTEX, phi2(phi1(d)), getEmbedding(VERTEX, phi2(phi1(d)))); + setCurrentLevel(getMaxLevel()-1) ; + } + setCurrentLevel(getMaxLevel()) ; + + + popLevel() ; +} + +void Map3MR_PrimalRegular::analysis() +{ + assert(getCurrentLevel() > 0 || !"analysis : called on level 0") ; + + decCurrentLevel() ; + + for(unsigned int i = 0; i < analysisFilters.size(); ++i) + (*analysisFilters[i])() ; +} + +void Map3MR_PrimalRegular::synthesis() +{ + assert(getCurrentLevel() < getMaxLevel() || !"synthesis : called on max level") ; + + for(unsigned int i = 0; i < synthesisFilters.size(); ++i) + (*synthesisFilters[i])() ; + + incCurrentLevel() ; +} + +} //end namespace CGoGN + +// setCurrentLevel(getMaxLevel()-1) ; +// //3.1 sauvegarde des relations +// std::vector > subdividedFaces; +// subdividedFaces.reserve(1024); +// TraversorF travF2(*this) ; +// for(Dart dit = travF2.begin(); dit != travF2.end(); dit = travF2.next()) +// { +// setCurrentLevel(getMaxLevel()) ; +// saveRelationsAroundVertex(phi2(phi1(dit)), subdividedFaces); +// setCurrentLevel(getMaxLevel()-1) ; +// } +// +// //3.2 decoupe des faces +// setCurrentLevel(getMaxLevel()) ; +// unsewAroundVertex(subdividedFaces); +// setCurrentLevel(getMaxLevel()-1) ; + +// //3.3 quadrangulation des faces interieurs +// TraversorV travV(*this) ; +// for (Dart dit = travV.begin(); dit != travV.end(); dit = travV.next()) +// { +// setCurrentLevel(getMaxLevel()) ; +// +// Map2::fillHole(phi1(dit)); +// +// Dart old = phi2(phi1(dit)); +// Dart bc = newBoundaryCycle(faceDegree(old)); +// sewVolumes(old, bc, false); +// +// if (isOrbitEmbedded(VERTEX)) +// { +// Dart it = bc; +// do +// { +// copyDartEmbedding(VERTEX, it, phi1(phi3(it))); +// it = phi1(it) ; +// } while(it != bc) ; +// } +// +// Dart dd = phi1(phi1(old)) ; +// splitFace(old,dd) ; +// +// Dart ne = phi1(phi1(old)) ; +// +// cutEdge(ne); +// +// if(embedNewVertices) +// embedNewCell(VERTEX, phi1(ne)) ; +// +// Dart stop = phi2(phi1(ne)); +// ne = phi2(ne); +// do +// { +// dd = phi1(phi1(phi1(ne))); +// +// splitFace(ne, dd) ; +// +// ne = phi2(phi_1(ne)); +// dd = phi1(phi1(dd)); +// } +// while(dd != stop); + + +// quadranguleFace(dit); +// setCurrentLevel(getMaxLevel() - 1) ; +// } + +// //4 couture des relations precedemment sauvegarde +// for (std::vector >::iterator it = subdividedFaces.begin(); it != subdividedFaces.end(); ++it) +// { +// Dart f1 = phi2((*it).first); +// Dart f2 = phi2((*it).second); +// +// if(isBoundaryFace(f1) && isBoundaryFace(f2)) +// { +// std::cout << "plop" << std::endl; +// sewVolumes(f1, f2, false); +// } +// } + + +//void Map3MR_PrimalRegular::splitFaceInSurface(Dart d, Dart e) +//{ +// Map2::splitFace(d, e) ; +// +// if (isOrbitEmbedded(VERTEX)) +// { +// copyDartEmbedding(VERTEX, phi_1(e), d) ; +// copyDartEmbedding(VERTEX, phi_1(d), e) ; +// } +// if (isOrbitEmbedded(FACE)) +// { +// copyDartEmbedding(FACE, phi_1(d), d) ; +// embedNewCell(FACE, e) ; +// copyCell(FACE, e, d) ; +// } +//} +// +//Dart Map3MR_PrimalRegular::cutEdgeInSurface(Dart d) +//{ +// Dart nd = Map2::cutEdge(d) ; +// +// if (isOrbitEmbedded(EDGE)) +// { +// copyDartEmbedding(EDGE, phi2(d), d) ; +// embedNewCell(EDGE, nd) ; +// copyCell(EDGE, nd, d) ; +// } +// +// if(isOrbitEmbedded(FACE)) +// { +// copyDartEmbedding(FACE, nd, d) ; +// Dart e = phi2(nd) ; +// copyDartEmbedding(FACE, phi1(e), e) ; +// } +// +// return nd; +//} + +//void Map3MR_PrimalRegular::saveRelationsAroundVertex(Dart d, std::vector >& vd) +//{ +// //le brin est forcement du niveau cur +// Dart dit = d; +// +// do +// { +// vd.push_back(std::pair(dit,phi2(dit))); +// +// dit = phi2(phi_1(dit)); +// +// }while(dit != d); +//} +// +//void Map3MR_PrimalRegular::unsewAroundVertex(std::vector >& vd) +//{ +// //unsew the edge path +// for(std::vector >::iterator it = vd.begin() ; it != vd.end() ; ++it) +// { +// Dart dit = (*it).first; +// Dart dit2 = (*it).second; +// +// //if(!Map2::isBoundaryFace(dit)) +// //{ +// Map2::unsewFaces(dit); +// +// if(isOrbitEmbedded(VERTEX)) +// { +// copyDartEmbedding(VERTEX, phi2(dit2), dit); +// copyDartEmbedding(VERTEX, phi2(dit), dit2); +// } +// +// if(isOrbitEmbedded(EDGE)) +// { +// +// } +// //} +// } +//} +// +//Dart Map3MR_PrimalRegular::quadranguleFace(Dart d) +//{ +// Dart centralDart = NIL; +// +// Map2::fillHole(phi1(d)); +// +// Dart old = phi2(phi1(d)); +// Dart bc = newBoundaryCycle(faceDegree(old)); +// sewVolumes(old, bc, false); +// +// if (isOrbitEmbedded(VERTEX)) +// { +// Dart it = bc; +// do +// { +// //copyDartEmbedding(VERTEX, it, phi1(phi3(it))); +// embedOrbit(VERTEX, it, getEmbedding(VERTEX, phi1(phi3(it)))); +// it = phi1(it) ; +// } while(it != bc) ; +// } +// +// Dart dd = phi1(phi1(old)) ; +// splitFace(old,dd) ; +// +// Dart ne = phi1(phi1(old)) ; +// +// cutEdge(ne); +// centralDart = phi1(ne); +// +// Dart stop = phi2(phi1(ne)); +// ne = phi2(ne); +// do +// { +// dd = phi1(phi1(phi1(ne))); +// +// splitFace(ne, dd) ; +// +// ne = phi2(phi_1(ne)); +// dd = phi1(phi1(dd)); +// } +// while(dd != stop); +// +// return centralDart; +//} + + + -- GitLab