diff --git a/src/Mod/Mesh/App/Core/Elements.h b/src/Mod/Mesh/App/Core/Elements.h index fbea5359e..64f53c479 100644 --- a/src/Mod/Mesh/App/Core/Elements.h +++ b/src/Mod/Mesh/App/Core/Elements.h @@ -70,6 +70,14 @@ public: /** MeshEdge just a pair of two point indices */ typedef std::pair MeshEdge; +struct MeshExport EdgeCollapse +{ + unsigned long _fromPoint; + unsigned long _toPoint; + std::vector _removeFacets; + std::vector _changeFacets; +}; + /** * The MeshPoint class represents a point in the mesh data structure. The class inherits from * Vector3f and provides some additional information such as flag state and property value. @@ -250,6 +258,10 @@ public: * Decrement the index for each corner point that is higher than \a ulIndex. */ inline void Decrement (unsigned long ulIndex); + /** + * Checks if the facets references the given point index. + */ + inline bool HasPoint(unsigned long) const; /** * Replaces the index of the neighbour facet that is equal to \a ulOrig * by \a ulNew. If the facet does not have a neighbourt with this index @@ -826,6 +838,17 @@ inline void MeshFacet::Decrement (unsigned long ulIndex) if (_aulPoints[2] > ulIndex) _aulPoints[2]--; } +inline bool MeshFacet::HasPoint(unsigned long ulIndex) const +{ + if (_aulPoints[0] == ulIndex) + return true; + if (_aulPoints[1] == ulIndex) + return true; + if (_aulPoints[2] == ulIndex) + return true; + return false; +} + inline void MeshFacet::ReplaceNeighbour (unsigned long ulOrig, unsigned long ulNew) { if (_aulNeighbours[0] == ulOrig) diff --git a/src/Mod/Mesh/App/Core/TopoAlgorithm.cpp b/src/Mod/Mesh/App/Core/TopoAlgorithm.cpp index 11ff5c764..c0aa18656 100644 --- a/src/Mod/Mesh/App/Core/TopoAlgorithm.cpp +++ b/src/Mod/Mesh/App/Core/TopoAlgorithm.cpp @@ -25,8 +25,8 @@ #ifndef _PreComp_ # include -# include -# include +# include +# include #endif #include @@ -235,102 +235,102 @@ void MeshTopoAlgorithm::OptimizeTopology(float fMaxAngle) } } } - -// Cosine of the maximum angle in triangle (v1,v2,v3) -static float cos_maxangle(const Base::Vector3f &v1, - const Base::Vector3f &v2, - const Base::Vector3f &v3) -{ - float a = Base::Distance(v2,v3); - float b = Base::Distance(v3,v1); - float c = Base::Distance(v1,v2); - float A = a * (b*b + c*c - a*a); - float B = b * (c*c + a*a - b*b); - float C = c * (a*a + b*b - c*c); - return 0.5f * std::min(std::min(A,B),C) / (a*b*c); // min cosine == max angle -} - -static float swap_benefit(const Base::Vector3f &v1, const Base::Vector3f &v2, - const Base::Vector3f &v3, const Base::Vector3f &v4) -{ - Base::Vector3f n124 = (v4 - v2) % (v1 - v2); - Base::Vector3f n234 = (v3 - v2) % (v4 - v2); - if ((n124 * n234) <= 0.0f) - return 0.0f; // avoid normal flip - - return std::max(-cos_maxangle(v1,v2,v3), -cos_maxangle(v1,v3,v4)) - - std::max(-cos_maxangle(v1,v2,v4), -cos_maxangle(v2,v3,v4)); -} - -float MeshTopoAlgorithm::SwapEdgeBenefit(unsigned long f, int e) const -{ - const MeshFacetArray& faces = _rclMesh.GetFacets(); - const MeshPointArray& vertices = _rclMesh.GetPoints(); - - unsigned long n = faces[f]._aulNeighbours[e]; - if (n == ULONG_MAX) - return 0.0f; // border edge - - unsigned long v1 = faces[f]._aulPoints[e]; - unsigned long v2 = faces[f]._aulPoints[(e+1)%3]; - unsigned long v3 = faces[f]._aulPoints[(e+2)%3]; - unsigned short s = faces[n].Side(faces[f]); - if (s == USHRT_MAX) { - std::cerr << "MeshTopoAlgorithm::SwapEdgeBenefit: error in neighbourhood " - << "of faces " << f << " and " << n << std::endl; - return 0.0f; // topological error - } - unsigned long v4 = faces[n]._aulPoints[(s+2)%3]; - if (v3 == v4) { - std::cerr << "MeshTopoAlgorithm::SwapEdgeBenefit: duplicate faces " - << f << " and " << n << std::endl; - return 0.0f; // duplicate faces - } - return swap_benefit(vertices[v2], vertices[v3], - vertices[v1], vertices[v4]); -} -typedef std::pair FaceEdge; // (face, edge) pair -typedef std::pair FaceEdgePriority; +// Cosine of the maximum angle in triangle (v1,v2,v3) +static float cos_maxangle(const Base::Vector3f &v1, + const Base::Vector3f &v2, + const Base::Vector3f &v3) +{ + float a = Base::Distance(v2,v3); + float b = Base::Distance(v3,v1); + float c = Base::Distance(v1,v2); + float A = a * (b*b + c*c - a*a); + float B = b * (c*c + a*a - b*b); + float C = c * (a*a + b*b - c*c); + return 0.5f * std::min(std::min(A,B),C) / (a*b*c); // min cosine == max angle +} + +static float swap_benefit(const Base::Vector3f &v1, const Base::Vector3f &v2, + const Base::Vector3f &v3, const Base::Vector3f &v4) +{ + Base::Vector3f n124 = (v4 - v2) % (v1 - v2); + Base::Vector3f n234 = (v3 - v2) % (v4 - v2); + if ((n124 * n234) <= 0.0f) + return 0.0f; // avoid normal flip + + return std::max(-cos_maxangle(v1,v2,v3), -cos_maxangle(v1,v3,v4)) - + std::max(-cos_maxangle(v1,v2,v4), -cos_maxangle(v2,v3,v4)); +} + +float MeshTopoAlgorithm::SwapEdgeBenefit(unsigned long f, int e) const +{ + const MeshFacetArray& faces = _rclMesh.GetFacets(); + const MeshPointArray& vertices = _rclMesh.GetPoints(); + + unsigned long n = faces[f]._aulNeighbours[e]; + if (n == ULONG_MAX) + return 0.0f; // border edge + + unsigned long v1 = faces[f]._aulPoints[e]; + unsigned long v2 = faces[f]._aulPoints[(e+1)%3]; + unsigned long v3 = faces[f]._aulPoints[(e+2)%3]; + unsigned short s = faces[n].Side(faces[f]); + if (s == USHRT_MAX) { + std::cerr << "MeshTopoAlgorithm::SwapEdgeBenefit: error in neighbourhood " + << "of faces " << f << " and " << n << std::endl; + return 0.0f; // topological error + } + unsigned long v4 = faces[n]._aulPoints[(s+2)%3]; + if (v3 == v4) { + std::cerr << "MeshTopoAlgorithm::SwapEdgeBenefit: duplicate faces " + << f << " and " << n << std::endl; + return 0.0f; // duplicate faces + } + return swap_benefit(vertices[v2], vertices[v3], + vertices[v1], vertices[v4]); +} + +typedef std::pair FaceEdge; // (face, edge) pair +typedef std::pair FaceEdgePriority; void MeshTopoAlgorithm::OptimizeTopology() { - // Find all edges that can be swapped and insert them into a - // priority queue - const MeshFacetArray& faces = _rclMesh.GetFacets(); - unsigned long nf = _rclMesh.CountFacets(); - std::priority_queue todo; - for (unsigned long i = 0; i < nf; i++) { - for (int j = 0; j < 3; j++) { - float b = SwapEdgeBenefit(i, j); - if (b > 0.0f) - todo.push(std::make_pair(b, std::make_pair(i, j))); - } - } - - // Edges are sorted in decreasing order with respect to their benefit - while (!todo.empty()) { - unsigned long f = todo.top().second.first; - int e = todo.top().second.second; - todo.pop(); - // Check again if the swap should still be done - if (SwapEdgeBenefit(f, e) <= 0.0f) - continue; - // OK, swap the edge - unsigned long f2 = faces[f]._aulNeighbours[e]; - SwapEdge(f, f2); - // Insert new edges into queue, if necessary - for (int j = 0; j < 3; j++) { - float b = SwapEdgeBenefit(f, j); - if (b > 0.0f) - todo.push(std::make_pair(b, std::make_pair(f, j))); - } - for (int j = 0; j < 3; j++) { - float b = SwapEdgeBenefit(f2, j); - if (b > 0.0f) - todo.push(std::make_pair(b, std::make_pair(f2, j))); - } - } + // Find all edges that can be swapped and insert them into a + // priority queue + const MeshFacetArray& faces = _rclMesh.GetFacets(); + unsigned long nf = _rclMesh.CountFacets(); + std::priority_queue todo; + for (unsigned long i = 0; i < nf; i++) { + for (int j = 0; j < 3; j++) { + float b = SwapEdgeBenefit(i, j); + if (b > 0.0f) + todo.push(std::make_pair(b, std::make_pair(i, j))); + } + } + + // Edges are sorted in decreasing order with respect to their benefit + while (!todo.empty()) { + unsigned long f = todo.top().second.first; + int e = todo.top().second.second; + todo.pop(); + // Check again if the swap should still be done + if (SwapEdgeBenefit(f, e) <= 0.0f) + continue; + // OK, swap the edge + unsigned long f2 = faces[f]._aulNeighbours[e]; + SwapEdge(f, f2); + // Insert new edges into queue, if necessary + for (int j = 0; j < 3; j++) { + float b = SwapEdgeBenefit(f, j); + if (b > 0.0f) + todo.push(std::make_pair(b, std::make_pair(f, j))); + } + for (int j = 0; j < 3; j++) { + float b = SwapEdgeBenefit(f2, j); + if (b > 0.0f) + todo.push(std::make_pair(b, std::make_pair(f2, j))); + } + } } void MeshTopoAlgorithm::DelaunayFlip(float fMaxAngle) @@ -891,6 +891,27 @@ bool MeshTopoAlgorithm::CollapseEdge(unsigned long ulFacetPos, unsigned long ulN return true; } +bool MeshTopoAlgorithm::CollapseEdge(const EdgeCollapse& ec) +{ + std::vector::const_iterator it; + for (it = ec._removeFacets.begin(); it != ec._removeFacets.end(); ++it) { + MeshFacet& f = _rclMesh._aclFacetArray[*it]; + f.SetInvalid(); + } + + for (it = ec._changeFacets.begin(); it != ec._changeFacets.end(); ++it) { + MeshFacet& f = _rclMesh._aclFacetArray[*it]; + + // The neighbourhood might be broken from now on!!! + f.Transpose(ec._fromPoint, ec._toPoint); + } + + _rclMesh._aclPointArray[ec._fromPoint].SetInvalid(); + + _needsCleanup = true; + return true; +} + bool MeshTopoAlgorithm::CollapseFacet(unsigned long ulFacetPos) { MeshFacet& rclF = _rclMesh._aclFacetArray[ulFacetPos]; diff --git a/src/Mod/Mesh/App/Core/TopoAlgorithm.h b/src/Mod/Mesh/App/Core/TopoAlgorithm.h index b0262d5f7..598804de9 100644 --- a/src/Mod/Mesh/App/Core/TopoAlgorithm.h +++ b/src/Mod/Mesh/App/Core/TopoAlgorithm.h @@ -40,6 +40,8 @@ namespace MeshCore { class AbstractPolygonTriangulator; +struct EdgeCollapse; + /** * The MeshTopoAlgorithm class provides several algorithms to manipulate a mesh. * It supports various mesh operations like inserting a new vertex, swapping the @@ -113,11 +115,11 @@ public: * Collapses the common edge of two adjacent facets. This operation removes * one common point of the collapsed edge and the facets \a ulFacetPos and * \a ulNeighbour from the data structure. - * @note If \a ulNeighbour is the neighbour facet on the i-th side then the - * i-th point is removed whereat i is 0, 1 or 2. If the other common point - * should be removed then CollapseEdge() should be invoked with transposed - * arguments of \a ulFacetPos and \a ulNeighbour, i.e. CollapseEdge - * ( \a ulNeighbour, \a ulFacetPos ). + * @note If \a ulNeighbour is the neighbour facet on the i-th side of + * \a ulFacetPos then the i-th point is removed whereas i is 0, 1 or 2. + * If the other common point should be removed then CollapseEdge() + * should be invoked with swapped arguments of \a ulFacetPos and + * \a ulNeighbour, i.e. CollapseEdge( \a ulNeighbour, \a ulFacetPos ). * * @note The client programmer must make sure that this is a legal operation. * @@ -132,6 +134,10 @@ public: * must take care not to use such elements. */ bool CollapseEdge(unsigned long ulFacetPos, unsigned long ulNeighbour); + /** + * Convenience function that passes already all needed information. + */ + bool CollapseEdge(const EdgeCollapse& ec); /** * Removes the facet with index \a ulFacetPos and all its neighbour facets. * The three vertices that are referenced by this facet are replaced by its @@ -144,7 +150,7 @@ public: * inconsistent stage. To make the structure consistent again Cleanup() should * be called. * The reason why this cannot be done automatically is that it would become - * quite slow if a lot of edges should be collapsed. + * quite slow if a lot of facets should be collapsed. * * @note While the mesh structure has invalid elements the client programmer * must take care not to use such elements.