extend mesh repair functions

This commit is contained in:
wmayer
2019-03-13 01:39:27 +01:00
parent 51aeb5867a
commit 4e018c506b
14 changed files with 289 additions and 56 deletions

View File

@@ -1090,43 +1090,37 @@ void MeshAlgorithm::CheckFacets(const MeshFacetGrid& rclGrid, const Base::ViewPr
// Precompute the polygon's bounding box
Base::BoundBox2d clPolyBBox = rclPoly.CalcBoundBox();
// Falls true, verwende Grid auf Mesh, um Suche zu beschleunigen
if (bInner)
{
// if true use grid on mesh to speed up search
if (bInner) {
BoundBox3f clBBox3d;
BoundBox2d clViewBBox;
std::vector<unsigned long> aulAllElements;
// Iterator fuer die zu durchsuchenden B-Boxen des Grids
// iterator for the bounding box grids
MeshGridIterator clGridIter(rclGrid);
// alle B-Boxen durchlaufen und die Facets speichern
for (clGridIter.Init(); clGridIter.More(); clGridIter.Next())
{
for (clGridIter.Init(); clGridIter.More(); clGridIter.Next()) {
clBBox3d = clGridIter.GetBoundBox();
clViewBBox = clBBox3d.ProjectBox(&fixedProj);
if (clViewBBox.Intersect(clPolyBBox))
{
// alle Elemente in AllElements sammeln
if (clViewBBox.Intersect(clPolyBBox)) {
// collect all elements in aulAllElements
clGridIter.GetElements(aulAllElements);
}
}
// doppelte Elemente wieder entfernen
// remove duplicates
std::sort(aulAllElements.begin(), aulAllElements.end());
aulAllElements.erase(std::unique(aulAllElements.begin(), aulAllElements.end()), aulAllElements.end());
Base::SequencerLauncher seq( "Check facets", aulAllElements.size() );
Base::SequencerLauncher seq("Check facets", aulAllElements.size());
for (it = aulAllElements.begin(); it != aulAllElements.end(); ++it)
{
for (it = aulAllElements.begin(); it != aulAllElements.end(); ++it) {
bNoPointInside = true;
clGravityOfFacet.Set(0.0f, 0.0f, 0.0f);
MeshGeomFacet rclFacet = _rclMesh.GetFacet(*it);
for (int j=0; j<3; j++)
{
for (int j=0; j<3; j++) {
clPt2d = fixedProj(rclFacet._aclPoints[j]);
clGravityOfFacet += clPt2d;
if ((clPolyBBox.Contains(Base::Vector2d(clPt2d.x, clPt2d.y)) &&
rclPoly.Contains(Base::Vector2d(clPt2d.x, clPt2d.y))) ^ !bInner)
{
if (clPolyBBox.Contains(Base::Vector2d(clPt2d.x, clPt2d.y)) &&
rclPoly.Contains(Base::Vector2d(clPt2d.x, clPt2d.y))) {
raulFacets.push_back(*it);
bNoPointInside = false;
break;
@@ -1134,30 +1128,25 @@ void MeshAlgorithm::CheckFacets(const MeshFacetGrid& rclGrid, const Base::ViewPr
}
// if no facet point is inside the polygon then check also the gravity
if (bNoPointInside == true)
{
clGravityOfFacet *= 1.0f/3.0f;
if (bNoPointInside) {
clGravityOfFacet *= 1.0f/3.0f;
if ((clPolyBBox.Contains(Base::Vector2d(clGravityOfFacet.x, clGravityOfFacet.y)) &&
rclPoly.Contains(Base::Vector2d(clGravityOfFacet.x, clGravityOfFacet.y))) ^ !bInner)
raulFacets.push_back(*it);
if (clPolyBBox.Contains(Base::Vector2d(clGravityOfFacet.x, clGravityOfFacet.y)) &&
rclPoly.Contains(Base::Vector2d(clGravityOfFacet.x, clGravityOfFacet.y)))
raulFacets.push_back(*it);
}
seq.next();
}
}
// Dreiecke ausserhalb schneiden, dann alles durchsuchen
else
{
// When cutting triangles outside then go through all elements
else {
Base::SequencerLauncher seq("Check facets", _rclMesh.CountFacets());
for (clIter.Init(); clIter.More(); clIter.Next())
{
for (int j=0; j<3; j++)
{
for (clIter.Init(); clIter.More(); clIter.Next()) {
for (int j=0; j<3; j++) {
clPt2d = fixedProj(clIter->_aclPoints[j]);
if ((clPolyBBox.Contains(Base::Vector2d(clPt2d.x, clPt2d.y)) &&
rclPoly.Contains(Base::Vector2d(clPt2d.x, clPt2d.y))) ^ !bInner)
{
rclPoly.Contains(Base::Vector2d(clPt2d.x, clPt2d.y))) == false) {
raulFacets.push_back(clIter.Position());
break;
}
@@ -1809,6 +1798,16 @@ void MeshRefPointToFacets::RemoveNeighbour(unsigned long pos, unsigned long face
_map[pos].erase(facet);
}
void MeshRefPointToFacets::RemoveFacet(unsigned long facetIndex)
{
unsigned long p0, p1, p2;
_rclMesh.GetFacetPoints(facetIndex, p0, p1, p2);
_map[p0].erase(facetIndex);
_map[p1].erase(facetIndex);
_map[p2].erase(facetIndex);
}
//----------------------------------------------------------------------------
void MeshRefFacetToFacets::Rebuild (void)

View File

@@ -379,6 +379,7 @@ public:
Base::Vector3f GetNormal(unsigned long) const;
void AddNeighbour(unsigned long, unsigned long);
void RemoveNeighbour(unsigned long, unsigned long);
void RemoveFacet(unsigned long);
protected:
void SearchNeighbours(const MeshFacetArray& rFacets, unsigned long index, const Base::Vector3f &rclCenter,

View File

@@ -26,6 +26,7 @@
#ifndef _PreComp_
# include <algorithm>
# include <map>
# include <queue>
#endif
#include "Degeneration.h"
@@ -505,15 +506,97 @@ bool MeshFixDegeneratedFacets::Fixup()
return true;
}
unsigned long MeshFixDegeneratedFacets::RemoveEdgeTooSmall (float fMinEdgeLength, float fMinEdgeAngle)
bool MeshRemoveSmallEdges::Fixup()
{
typedef std::pair<unsigned long, int> FaceEdge; // (face, edge) pair
typedef std::pair<float, FaceEdge> FaceEdgePriority;
MeshTopoAlgorithm topAlg(_rclMesh);
MeshRefPointToFacets vf_it(_rclMesh);
const MeshFacetArray &rclFAry = _rclMesh.GetFacets();
const MeshPointArray &rclPAry = _rclMesh.GetPoints();
rclFAry.ResetInvalid();
rclPAry.ResetInvalid();
rclPAry.ResetFlag(MeshPoint::VISIT);
std::size_t facetCount = rclFAry.size();
std::priority_queue<FaceEdgePriority,
std::vector<FaceEdgePriority>,
std::greater<FaceEdgePriority> > todo;
for (std::size_t index = 0; index < facetCount; index++) {
for (int i=0; i<3; i++) {
const MeshFacet& facet = rclFAry[index];
const Base::Vector3f& p1 = rclPAry[facet._aulPoints[i]];
const Base::Vector3f& p2 = rclPAry[facet._aulPoints[(i+1)%3]];
float distance = Base::Distance(p1, p2);
if (distance < fMinEdgeLength) {
unsigned long facetIndex = static_cast<unsigned long>(index);
todo.push(std::make_pair(distance, std::make_pair(facetIndex, i)));
}
}
}
bool removedEdge = false;
while (!todo.empty()) {
FaceEdge faceedge = todo.top().second;
todo.pop();
// check if one of the face pairs was already processed
if (!rclFAry[faceedge.first].IsValid())
continue;
// the facet points may have changed, so check the current distance again
const MeshFacet& facet = rclFAry[faceedge.first];
const Base::Vector3f& p1 = rclPAry[facet._aulPoints[faceedge.second]];
const Base::Vector3f& p2 = rclPAry[facet._aulPoints[(faceedge.second+1)%3]];
float distance = Base::Distance(p1, p2);
if (distance >= fMinEdgeLength)
continue;
// collect the collapse-edge information
EdgeCollapse ce;
ce._fromPoint = rclFAry[faceedge.first]._aulPoints[faceedge.second];
ce._toPoint = rclFAry[faceedge.first]._aulPoints[(faceedge.second+1)%3];
ce._removeFacets.push_back(faceedge.first);
unsigned long neighbour = rclFAry[faceedge.first]._aulNeighbours[faceedge.second];
if (neighbour != ULONG_MAX)
ce._removeFacets.push_back(neighbour);
std::set<unsigned long> vf = vf_it[ce._fromPoint];
vf.erase(faceedge.first);
if (neighbour != ULONG_MAX)
vf.erase(neighbour);
ce._changeFacets.insert(ce._changeFacets.begin(), vf.begin(), vf.end());
if (topAlg.IsCollapseEdgeLegal(ce)) {
topAlg.CollapseEdge(ce);
for (auto it : ce._removeFacets) {
vf_it.RemoveFacet(it);
}
for (auto it : ce._changeFacets) {
vf_it.RemoveNeighbour(ce._fromPoint, it);
vf_it.AddNeighbour(ce._toPoint, it);
}
removedEdge = true;
}
}
if (removedEdge) {
topAlg.Cleanup();
_rclMesh.RebuildNeighbours();
}
return true;
#if 0
unsigned long ulCtLastLoop, ulCtFacets = _rclMesh.CountFacets();
const MeshFacetArray &rclFAry = _rclMesh.GetFacets();
const MeshPointArray &rclPAry = _rclMesh.GetPoints();
MeshFacetArray::_TConstIterator f_beg = rclFAry.begin();
// repeat until no facet van be removed
// repeat until no facet can be removed
do {
MeshRefPointToFacets clPt2Facets(_rclMesh);
@@ -523,7 +606,7 @@ unsigned long MeshFixDegeneratedFacets::RemoveEdgeTooSmall (float fMinEdgeLength
std::set<std::pair<unsigned long, unsigned long> > aclPtDelList;
MeshFacetIterator clFIter(_rclMesh), clFN0(_rclMesh), clFN1(_rclMesh), clFN2(_rclMesh);
MeshFacetIterator clFIter(_rclMesh);
for (clFIter.Init(); clFIter.More(); clFIter.Next()) {
MeshGeomFacet clSFacet = *clFIter;
Base::Vector3f clP0 = clSFacet._aclPoints[0];
@@ -537,26 +620,26 @@ unsigned long MeshFixDegeneratedFacets::RemoveEdgeTooSmall (float fMinEdgeLength
unsigned long ulP1 = clFacet._aulPoints[1];
unsigned long ulP2 = clFacet._aulPoints[2];
if ((Base::Distance(clP0, clP1) < fMinEdgeLength) ||
(clE20.GetAngle(-clE12) < fMinEdgeAngle)) {
if (Base::Distance(clP0, clP1) < fMinEdgeLength) {
// delete point P1 on P0
aclPtDelList.insert(std::make_pair
(std::min<unsigned long>(ulP1, ulP0), std::max<unsigned long>(ulP1, ulP0)));
clFIter.SetFlag(MeshFacet::INVALID);
}
else if ((Base::Distance(clP1, clP2) < fMinEdgeLength) ||
(clE01.GetAngle(-clE20) < fMinEdgeAngle)) {
else if (Base::Distance(clP1, clP2) < fMinEdgeLength) {
// delete point P2 on P1
aclPtDelList.insert(std::make_pair
(std::min<unsigned long>(ulP2, ulP1), std::max<unsigned long>(ulP2, ulP1)));
clFIter.SetFlag(MeshFacet::INVALID);
}
else if ((Base::Distance(clP2, clP0) < fMinEdgeLength) ||
(clE12.GetAngle(-clE01) < fMinEdgeAngle)) {
else if (Base::Distance(clP2, clP0) < fMinEdgeLength) {
// delete point P0 on P2
aclPtDelList.insert(std::make_pair
(std::min<unsigned long>(ulP0, ulP2), std::max<unsigned long>(ulP0, ulP2)));
clFIter.SetFlag(MeshFacet::INVALID);
}
}
#if 0
// remove points, fix indices
for (std::set<std::pair<unsigned long, unsigned long> >::iterator pI = aclPtDelList.begin();
pI != aclPtDelList.end(); ++pI) {
@@ -588,7 +671,7 @@ unsigned long MeshFixDegeneratedFacets::RemoveEdgeTooSmall (float fMinEdgeLength
}
}
}
#endif
ulCtLastLoop = _rclMesh.CountFacets();
_rclMesh.RemoveInvalids();
}
@@ -596,7 +679,8 @@ unsigned long MeshFixDegeneratedFacets::RemoveEdgeTooSmall (float fMinEdgeLength
_rclMesh.RebuildNeighbours();
return ulCtFacets - _rclMesh.CountFacets();
return ulCtFacets > _rclMesh.CountFacets();
#endif
}
// ----------------------------------------------------------------------

View File

@@ -324,16 +324,38 @@ public:
* Removes degenerated facets.
*/
bool Fixup ();
/**
* Removes all facets with an edge smaller than \a fMinEdgeLength without leaving holes or gaps
* in the mesh. Returns the number of removed facets.
*/
unsigned long RemoveEdgeTooSmall (float fMinEdgeLength = MeshDefinitions::_fMinPointDistance,
float fMinEdgeAngle = MeshDefinitions::_fMinEdgeAngle);
private:
float fEpsilon;
};
/**
* The MeshRemoveSmallEdges class tries to fix degenerations by removing small edges.
* @see MeshFixDegeneratedFacets
* @author Werner Mayer
*/
class MeshExport MeshRemoveSmallEdges : public MeshValidation
{
public:
/**
* Construction.
*/
MeshRemoveSmallEdges (MeshKernel &rclM, float fMinEdgeLen = MeshDefinitions::_fMinPointDistance)
: MeshValidation(rclM), fMinEdgeLength(fMinEdgeLen) { }
/**
* Destruction.
*/
~MeshRemoveSmallEdges () { }
/**
* Removes all facets with an edge smaller than \a fMinEdgeLength without leaving holes or gaps
* in the mesh.
*/
bool Fixup ();
private:
float fMinEdgeLength;
};
/**
* The MeshEvalDeformedFacets class searches for deformed facets. A facet is regarded as deformed
* if an angle is < 30 deg or > 120 deg.

View File

@@ -1214,3 +1214,22 @@ bool MeshGeomFacet::IsPointOfSphere(const MeshGeomFacet& rFacet) const
return false;
}
float MeshGeomFacet::AspectRatio() const
{
Base::Vector3f center;
float radius = CenterOfCircumCircle(center);
float minLength = FLOAT_MAX;
for (int i=0; i<3; i++) {
const Base::Vector3f& p0 = _aclPoints[i];
const Base::Vector3f& p1 = _aclPoints[(i+1)%3];
float distance = Base::Distance(p0, p1);
if (distance < minLength) {
minLength = distance;
}
}
if (minLength == 0)
return FLOAT_MAX;
return radius / minLength;
}

View File

@@ -484,6 +484,9 @@ public:
* If one of the facet's points is inside the sphere true is returned, otherwise false.
*/
bool IsPointOfSphere(const MeshGeomFacet& rFacet) const;
/** The aspect ratio is the ratio of the radius of the circum-circle to the shortest edge length.
*/
float AspectRatio() const;
protected:
Base::Vector3f _clNormal; /**< Normal of the facet. */

View File

@@ -953,6 +953,35 @@ bool MeshTopoAlgorithm::CollapseEdge(unsigned long ulFacetPos, unsigned long ulN
return true;
}
bool MeshTopoAlgorithm::IsCollapseEdgeLegal(const EdgeCollapse& ec) const
{
std::vector<unsigned long>::const_iterator it;
for (it = ec._changeFacets.begin(); it != ec._changeFacets.end(); ++it) {
MeshFacet f = _rclMesh._aclFacetArray[*it];
if (!f.IsValid())
return false;
// ignore the facet(s) at this edge
if (f.HasPoint(ec._fromPoint) && f.HasPoint(ec._toPoint))
continue;
MeshGeomFacet tria1 = _rclMesh.GetFacet(f);
f.Transpose(ec._fromPoint, ec._toPoint);
MeshGeomFacet tria2 = _rclMesh.GetFacet(f);
if (tria1.GetNormal() * tria2.GetNormal() < 0.0f)
return false;
}
if (!_rclMesh._aclPointArray[ec._fromPoint].IsValid())
return false;
if (!_rclMesh._aclPointArray[ec._toPoint].IsValid())
return false;
return true;
}
bool MeshTopoAlgorithm::CollapseEdge(const EdgeCollapse& ec)
{
std::vector<unsigned long>::const_iterator it;

View File

@@ -79,7 +79,7 @@ public:
* Swaps the common edge of two adjacent facets even if the operation might
* be illegal. To be sure that this operation is legal, check either with
* IsSwapEdgeLegal() or ShouldSwapEdge() before.
* An illegel swap edge operation can produce non-manifolds, degenerated
* An illegal swap edge operation can produce non-manifolds, degenerated
* facets or it might create a fold on the surface, i.e. geometric overlaps
* of several triangles.
*/
@@ -116,6 +116,12 @@ public:
* by three facets is supposrted.
*/
bool CollapseVertex(const VertexCollapse& vc);
/**
* Checks whether a collapse edge operation is legal, that is fulfilled if none of the
* adjacent facets flips its normal. If this operation is legal
* true is returned, false is returned if this operation is illegal.
*/
bool IsCollapseEdgeLegal(const EdgeCollapse& ec) const;
/**
* Collapses the common edge of two adjacent facets. This operation removes
* one common point of the collapsed edge and the facets \a ulFacetPos and

View File

@@ -98,5 +98,17 @@ The two angles are given in radian.
</Documentation>
<Parameter Name="Area" Type="Float"/>
</Attribute>
<Attribute Name="AspectRatio" ReadOnly="true">
<Documentation>
<UserDocu>The aspect ratio of the facet</UserDocu>
</Documentation>
<Parameter Name="AspectRatio" Type="Float"/>
</Attribute>
<Attribute Name="CircumCircle" ReadOnly="true">
<Documentation>
<UserDocu>The center and radius of the circum-circle</UserDocu>
</Documentation>
<Parameter Name="CircumCircle" Type="Tuple"/>
</Attribute>
</PythonExport>
</GenerateModel>

View File

@@ -29,6 +29,7 @@
#include <Mod/Mesh/App/FacetPy.cpp>
#include <Base/VectorPy.h>
#include <Base/GeometryPyCXX.h>
using namespace Mesh;
@@ -220,6 +221,35 @@ Py::Float FacetPy::getArea(void) const
return Py::Float(tria.Area());
}
Py::Float FacetPy::getAspectRatio(void) const
{
FacetPy::PointerType face = this->getFacetPtr();
if (!face->isBound()) {
return Py::Float(-1.0);
}
const MeshCore::MeshKernel& kernel = face->Mesh->getKernel();
MeshCore::MeshGeomFacet tria = kernel.GetFacet(face->Index);
return Py::Float(tria.AspectRatio());
}
Py::Tuple FacetPy::getCircumCircle(void) const
{
FacetPy::PointerType face = this->getFacetPtr();
if (!face->isBound()) {
return Py::None();
}
const MeshCore::MeshKernel& kernel = face->Mesh->getKernel();
MeshCore::MeshGeomFacet tria = kernel.GetFacet(face->Index);
Base::Vector3f center;
float radius = tria.CenterOfCircumCircle(center);
Py::Tuple tuple(2);
tuple.setItem(0, Py::Vector(center));
tuple.setItem(1, Py::Float(radius));
return tuple;
}
PyObject *FacetPy::getCustomAttributes(const char* /*attr*/) const
{
return 0;

View File

@@ -1128,6 +1128,15 @@ void MeshObject::refine()
this->_segments.clear();
}
void MeshObject::removeSmallEdges(float length)
{
unsigned long count = _kernel.CountFacets();
MeshCore::MeshRemoveSmallEdges eval(_kernel, length);
eval.Fixup();
if (_kernel.CountFacets() < count)
this->_segments.clear();
}
void MeshObject::optimizeTopology(float fMaxAngle)
{
MeshCore::MeshTopoAlgorithm topalg(_kernel);
@@ -1411,8 +1420,8 @@ void MeshObject::validateDeformations(float fMaxAngle, float fEps)
{
unsigned long count = _kernel.CountFacets();
MeshCore::MeshFixDeformedFacets eval(_kernel,
Base::toRadians(30.0f),
Base::toRadians(120.0f),
Base::toRadians(15.0f),
Base::toRadians(150.0f),
fMaxAngle, fEps);
eval.Fixup();
if (_kernel.CountFacets() < count)

View File

@@ -255,6 +255,7 @@ public:
/** @name Topological operations */
//@{
void refine();
void removeSmallEdges(float);
void optimizeTopology(float);
void optimizeEdges();
void splitEdges();

View File

@@ -136,7 +136,12 @@ mesh.write(Stream=file,Format='STL',[Name='Object name',Material=colors])</UserD
<UserDocu>Remove a list of facet indices from the mesh</UserDocu>
</Documentation>
</Methode>
<Methode Name="getInternalFacets" Const="true">
<Methode Name="removeSmallEdges">
<Documentation>
<UserDocu>Remove all edges that are smaller than a given length</UserDocu>
</Documentation>
</Methode>
<Methode Name="getInternalFacets" Const="true">
<Documentation>
<UserDocu>Builds a list of facet indices with triangles that are inside a volume mesh</UserDocu>
</Documentation>

View File

@@ -1254,6 +1254,19 @@ PyObject* MeshPy::refine(PyObject *args)
Py_Return;
}
PyObject* MeshPy::removeSmallEdges(PyObject *args)
{
float length;
if (!PyArg_ParseTuple(args, "f", &length))
return NULL;
PY_TRY {
getMeshObjectPtr()->removeSmallEdges(length);
} PY_CATCH;
Py_Return;
}
PyObject* MeshPy::optimizeTopology(PyObject *args)
{
float fMaxAngle=-1.0f;