#ifndef HALFEDGEQUADCLEAN_H #define HALFEDGEQUADCLEAN_H #include #include #include #include #include #include namespace vcg { namespace tri { /*! * \brief The class provides methods for detecting doublets and singlets and removing them * */ template class HalfedgeQuadClean { protected: typedef typename MeshType::VertexPointer VertexPointer; typedef typename MeshType::EdgePointer EdgePointer; typedef typename MeshType::HEdgePointer HEdgePointer; typedef typename MeshType::FacePointer FacePointer; typedef typename MeshType::VertexIterator VertexIterator; typedef typename MeshType::EdgeIterator EdgeIterator; typedef typename MeshType::HEdgeIterator HEdgeIterator; typedef typename MeshType::FaceIterator FaceIterator; /*! Adds to a queue all faces on the 1-ring of a given vertex * * \param q Queue of faces * \param vp Pointer to the vertex * */ static void add_faces(queue &q, VertexPointer vp) { vector faces = HalfEdgeTopology::get_incident_faces(vp); for(typename vector::iterator fi = faces.begin(); fi != faces.end(); ++fi) { if(*fi) q.push(*fi); } } /*! Removes doublets from all the faces into a queue until the queue empties * * \param m Mesh * \param faces Set of all faces modified by the removals * \param q Queue of faces to clean * */ static void remove_doublets(MeshType &m, set &faces, queue &q) { while(!q.empty()) { FacePointer fp = q.front(); q.pop(); if( !fp->IsD() ) { faces.insert(remove_doublet(m,fp, &q)); } } } /*! Removes a doublet and all the other ones that the removal may have generated * * \param m Mesh * \param fp Pointer to the face to clean from doublets * \param Queue of faces to check for new doublets * * \return The new face generated after doublet removal */ static FacePointer remove_doublet(MeshType &m, FacePointer fp, queue *q = NULL) { vector hedges = HalfEdgeTopology::find_doublet_hedges_quad(fp); assert(hedges.size() <= 4); switch(hedges.size()) { // No doublets case 0: return NULL; // A single doublet case 1: if(q) { add_faces(*q, hedges[0]->HNp()->HVp()); add_faces(*q, hedges[0]->HPp()->HVp()); } return HalfEdgeTopology::doublet_remove_quad(m, hedges[0]->HVp()); // Two doublets on the same face case 2: { if(q) { add_faces(*q, hedges[0]->HNp()->HVp()); add_faces(*q, hedges[0]->HPp()->HVp()); } FacePointer fp1 = HalfEdgeTopology::doublet_remove_quad(m, hedges[0]->HVp()); // Removal of the doublet may generate a singlet if(HalfEdgeTopology::is_singlet_quad(fp1)) { HEdgePointer hp = HalfEdgeTopology::singlet_remove_quad(m, fp1); if(hp) { if(q) { if(hp->HFp()) q->push(hp->HFp()); if(hp->HOp()->HFp()) q->push(hp->HOp()->HFp()); } int valence1, valence2; valence1 = HalfEdgeTopology::vertex_valence(hp->HVp()); valence2 = HalfEdgeTopology::vertex_valence(hp->HOp()->HVp()); // Check if the removal of the singlet generated other singlets, then iteratively remove them while(valence1 == 1 || valence2 == 1) { assert(! (valence1 == 1 && valence2 == 1)); FacePointer singlet_pointer; if(valence1 == 1 ) singlet_pointer = hp->HFp(); else singlet_pointer = hp->HOp()->HFp(); hp = HalfEdgeTopology::singlet_remove_quad(m, singlet_pointer); if(!hp) break; if(q) { if(hp->HFp()) q->push(hp->HFp()); if(hp->HOp()->HFp()) q->push(hp->HOp()->HFp()); } valence1 = HalfEdgeTopology::vertex_valence(hp->HVp()); valence2 = HalfEdgeTopology::vertex_valence(hp->HOp()->HVp()); } } } return fp1; } // Four doublets: simply remove one of the two faces case 4: HalfEdgeTopology::remove_face(m,fp->FHp()->HOp()->HFp()); return fp; default: assert(0); } } public: /*! Removes doublets from all the faces on the 1-ring of the given vertices * * \param m Mesh * \param Set of all faces modified by the removals * \param vertices Vector of vertices that will be * */ static void remove_doublets(MeshType &m, set &faces, vector vertices) { queue q; for(typename vector::iterator vi = vertices.begin(); vi != vertices.end(); ++vi) { vector inc_faces = HalfEdgeTopology::get_incident_faces(*vi); for(typename vector::iterator fi = inc_faces.begin(); fi != inc_faces.end(); ++fi) { if(*fi) if( !((*fi)->IsD()) ) q.push(*fi); } } remove_doublets(m, faces, q); } /*! Removes doublets from all the faces on the 1-ring of the given vertices * * \param m Mesh * \param vertices Vector of vertices that will be * */ static void remove_doublets(MeshType &m, vector vertices) { set faces; remove_doublets(m,faces,vertices); } /*! Removes doublets from a set of faces * * \param m Mesh * \param set of faces to clean * */ static void remove_doublets(MeshType &m, set &faces) { queue q; for(typename set::iterator fi = faces.begin(); fi != faces.end(); ++fi) { if(*fi) if( !((*fi)->IsD()) ) q.push(*fi); } remove_doublets(m, faces, q); } /*! Removes all doublets from a mesh * * \param m Mesh * * \return Number of doublets removed */ static int remove_doublets(MeshType &m) { int count; int removed = 0; do { count = 0; for(FaceIterator fi = m.face.begin(); fi != m.face.end(); ++fi) { if( !((*fi).IsD()) ) { if(remove_doublet(m, &(*fi))) count++; } } removed += count; }while(count != 0); return removed; } /*! Removes all singlets from a mesh * * \param m Mesh * * \return Number of singlets removed */ static int remove_singlets(MeshType &m) { int removed = 0; int count; do { count = 0; for(FaceIterator fi = m.face.begin(); fi != m.face.end(); ++fi) { if( !((*fi).IsD()) ) { if( HalfEdgeTopology::is_singlet_quad(&(*fi)) ) { HalfEdgeTopology::singlet_remove_quad(m, &(*fi)); count++; } } } removed += count; }while(count != 0); return removed; } /*! Checks if a mesh has singlets * * \param m Mesh * * \return Value indicating whether mesh has singlets */ static bool has_singlets(MeshType &m) { for(FaceIterator fi = m.face.begin(); fi != m.face.end(); ++fi) { if( !((*fi).IsD()) ) { if( HalfEdgeTopology::is_singlet_quad(&(*fi)) ) return true; } } return false; } /*! Checks if a mesh has doublets * * \param m Mesh * * \return Value indicating whether mesh has doublets */ static bool has_doublets(MeshType &m) { for(FaceIterator fi = m.face.begin(); fi != m.face.end(); ++fi) { if( !((*fi).IsD()) ) { if( HalfEdgeTopology::has_doublet_quad(&(*fi)) ) return true; } } return false; } /*! Performs rotation of selected edges computing the best rotation * * * \param m Mesh * \param hedges Vector of halfedges representing the edges to rotate * \param faces Set of modified faces (every modified face will be inserted into this set) * */ template static void flip_edges(MeshType &m, vector &hedges, set &faces) { for(typename vector::iterator hi = hedges.begin(); hi != hedges.end(); ++hi) { // edge must be shared by two faces if((*hi)->HFp() && (*hi)->HOp()->HFp()) { if(!(*hi)->IsD() && !(*hi)->HFp()->IsD() && !(*hi)->HOp()->HFp()->IsD()) { typename PriorityType::FlipType fliptype = PriorityType::best_flip( *hi ); if(fliptype != PriorityType::NO_FLIP) { vector vertices; // Record old vertices for future removal of doublets vertices.push_back((*hi)->HVp()); vertices.push_back((*hi)->HOp()->HVp()); // Add modified faces into the set of faces faces.insert((*hi)->HFp()); faces.insert((*hi)->HOp()->HFp()); bool cw = (fliptype == PriorityType::CW_FLIP); HalfEdgeTopology::edge_rotate_quad((*hi),cw); // After a rotation doublets may have generated remove_doublets(m, faces, vertices); } } } } } /*! Performs edge rotations on the entire mesh computing the best rotation for each edge * * * \param m Mesh * */ template static int flip_edges(MeshType &m) { int count; int ret=0; do { count = 0; for(typename MeshType::EdgeIterator ei = m.edge.begin(); ei != m.edge.end(); ++ei) { if( !(ei->IsD()) ) { HEdgePointer hp = ei->EHp(); if(hp->HFp() && hp->HOp()->HFp()) { typename PriorityType::FlipType fliptype = PriorityType::best_flip( hp ); if(fliptype != PriorityType::NO_FLIP) { vector vertices; // Record old vertices for future removal of doublets vertices.push_back(hp->HVp()); vertices.push_back(hp->HOp()->HVp()); bool cw = (fliptype == PriorityType::CW_FLIP); HalfEdgeTopology::edge_rotate_quad(hp,cw); // After a rotation doublets may have generated remove_doublets(m, vertices); count++; } } } } ret+=count; }while(count != 0); return ret; } }; /*! * \brief Generic priority for edge rotations * */ template class EdgeFlipPriority { public: typedef typename MeshType::HEdgePointer HEdgePointer; /// Possible types of rotation enum FlipType { NO_FLIP, CW_FLIP, CCW_FLIP}; /*! * Computes the best rotation to perform * * \param hp Pointer to an halfedge representing the edge to rotate * * \return The best type of rotation */ static FlipType best_flip( HEdgePointer hp); }; /*! * \brief Priority based on maximizing vertices regularity * */ template class VertReg: public EdgeFlipPriority { public: typedef typename MeshType::VertexPointer VertexPointer; typedef typename MeshType::EdgePointer EdgePointer; typedef typename MeshType::HEdgePointer HEdgePointer; typedef typename MeshType::FacePointer FacePointer; typedef EdgeFlipPriority Base; typedef typename Base::FlipType FlipType; /// Default Constructor VertReg(){} ~VertReg(){} /*! * Computes the best rotation to perform for maximizing vertices regularity * * \param hp Pointer to an halfedge representing the edge to rotate * * \return The best type of rotation */ static FlipType best_flip( HEdgePointer hp) { assert(hp); assert(!hp->IsD()); vector vps1 = HalfEdgeTopology::getVertices(hp->HFp(), hp); vector vps2 = HalfEdgeTopology::getVertices(hp->HOp()->HFp(), hp->HOp()); valarray valences(6); /* Indices of vertices into vector valences: 3-------2 | | | f1 | | | 0-------1 | | | f2 | | | 4-------5 */ // Compute valencies of vertices for(int i=0; i<4; i++) valences[i] = HalfEdgeTopology::vertex_valence(vps1[i]) - 4 ; valences[4] = HalfEdgeTopology::vertex_valence(vps2[2]) - 4; valences[5] = HalfEdgeTopology::vertex_valence(vps2[3]) - 4; // Vector containing sums of the valencies in all the possible cases: no rotation, ccw rotation, cw rotation vector sums; // positions: // sums[0]: now (No rotation) // sums[1]: flipping ccw // sums[2]: flipping cw // No rotation sums.push_back( pow(valences, 2.0).sum() ); // CCW valences[0]--; valences[1]--; valences[2]++; valences[4]++; sums.push_back( pow(valences, 2.0).sum() ); // CW valences[2]--; valences[4]--; valences[3]++; valences[5]++; sums.push_back( pow(valences, 2.0).sum() ); if( sums[2]<= sums[1] && sums[2]< sums[0] ) return Base::CW_FLIP; else if( sums[1]< sums[2] && sums[1]< sums[0] ) return Base::CCW_FLIP; return Base::NO_FLIP; } }; /*! * \brief Priority based on minimizing homeometry * */ template class Homeometry: public EdgeFlipPriority { public: typedef typename MeshType::VertexPointer VertexPointer; typedef typename MeshType::EdgePointer EdgePointer; typedef typename MeshType::HEdgePointer HEdgePointer; typedef typename MeshType::FacePointer FacePointer; typedef EdgeFlipPriority Base; typedef typename Base::FlipType FlipType; /// Default Constructor Homeometry(){} ~Homeometry(){} /*! * Computes the best rotation to perform for minimizing the distance from homeometry * * \param hp Pointer to an halfedge representing the edge to rotate * * \return The best type of rotation */ static FlipType best_flip( HEdgePointer hp) { assert(hp); assert(!hp->IsD()); vector face1 = HalfEdgeTopology::getVertices(hp->HFp(), hp); vector face2 = HalfEdgeTopology::getVertices(hp->HOp()->HFp(), hp->HOp()); // Vector containing sums of the valencies in all the possible cases: no rotation, ccw rotation, cw rotation vector sums; // positions: // sums[0]: now (No rotation) // sums[1]: flipping ccw // sums[2]: flipping cw // No rotation sums.push_back( distance_from_homeometry(face1, face2, 0) ); // CCW face1[1] = face2[2]; face2[1] = face1[2]; sums.push_back( distance_from_homeometry(face1, face2, 1) ); // CW face1[2] = face2[3]; face2[2] = face1[3]; sums.push_back( distance_from_homeometry(face1, face2, 2) ); if( sums[2]<= sums[1] && sums[2]< sums[0] ) return Base::CW_FLIP; else if( sums[1]< sums[2] && sums[1]< sums[0] ) return Base::CCW_FLIP; return Base::NO_FLIP; } protected: /*! * Computes the area of a quad * * \param vertices Vector of the four vertices of the quad * * \return Area of the quad */ static float area(vector &vertices) { assert(vertices.size() == 4); float tri1 = Norm( (vertices[1]->cP() - vertices[0]->cP()) ^ (vertices[2]->cP() - vertices[0]->cP()) ); float tri2 = Norm( (vertices[2]->cP() - vertices[0]->cP()) ^ (vertices[3]->cP() - vertices[0]->cP()) ); return (tri1+tri2) / 2; } /*! * Computes the distance of two faces from being homeometirc * * \param face1 Vector of vertices belonging to the first face * \param face2 Vector of vertices belonging to the second face * \param i Index of the edge to compute * * \return The computed homeometry */ static float distance_from_homeometry(vector &face1, vector &face2, int i) { // Ideal edge length float mu = sqrt( (area(face1) + area(face2)) / 2 ); // Length of the i-th edge (the edge changed with a rotation) float edge_length = Distance( face1[i]->cP(), face1[i+1]->cP() ); // Length of the diagonals valarray diagonals(4); diagonals[0] = Distance( face1[0]->cP(), face1[2]->cP() ); diagonals[1] = Distance( face1[1]->cP(), face1[3]->cP() ); diagonals[2] = Distance( face2[0]->cP(), face2[2]->cP() ); diagonals[3] = Distance( face2[1]->cP(), face2[3]->cP() ); // Ideal diagonal length float ideal_diag_length = SQRT_TWO*mu; float sum_diagonals = pow(diagonals - ideal_diag_length, 2.0).sum(); return (pow (edge_length - mu , static_cast(2.0)) + sum_diagonals); } }; } } #endif // HALFEDGEQUADCLEAN_H