diff --git a/src/Mod/Mesh/App/CMakeLists.txt b/src/Mod/Mesh/App/CMakeLists.txt index a0d547c3b5..fa355f18c6 100644 --- a/src/Mod/Mesh/App/CMakeLists.txt +++ b/src/Mod/Mesh/App/CMakeLists.txt @@ -88,6 +88,8 @@ SET(Core_SRCS Core/Triangulation.h Core/Trim.cpp Core/Trim.h + Core/TrimByPlane.cpp + Core/TrimByPlane.h Core/tritritest.h Core/Visitor.cpp Core/Visitor.h diff --git a/src/Mod/Mesh/App/Core/Elements.cpp b/src/Mod/Mesh/App/Core/Elements.cpp index 2422d45cf8..e56fe310c9 100644 --- a/src/Mod/Mesh/App/Core/Elements.cpp +++ b/src/Mod/Mesh/App/Core/Elements.cpp @@ -268,6 +268,25 @@ bool MeshGeomEdge::IntersectWithLine (const Base::Vector3f &rclPt, return dist2 + dist3 <= dist1 + eps; } +bool MeshGeomEdge::IntersectWithPlane (const Base::Vector3f &rclPt, + const Base::Vector3f &rclDir, + Base::Vector3f &rclRes) const +{ + float dist1 = _aclPoints[0].DistanceToPlane(rclPt, rclDir); + float dist2 = _aclPoints[1].DistanceToPlane(rclPt, rclDir); + + // either both points are below or above the plane + if (dist1 * dist2 >= 0.0f) + return false; + + Base::Vector3f u = _aclPoints[1] - _aclPoints[0]; + Base::Vector3f b = rclPt - _aclPoints[0]; + float t = b.Dot(rclDir) / u.Dot(rclDir); + rclRes = _aclPoints[0] + t * u; + + return true; +} + void MeshGeomEdge::ProjectPointToLine (const Base::Vector3f &rclPoint, Base::Vector3f &rclProj) const { diff --git a/src/Mod/Mesh/App/Core/Elements.h b/src/Mod/Mesh/App/Core/Elements.h index 716fb007a9..fe3f43827a 100644 --- a/src/Mod/Mesh/App/Core/Elements.h +++ b/src/Mod/Mesh/App/Core/Elements.h @@ -168,6 +168,10 @@ public: * with the edge. The intersection must be inside the edge. If there is no intersection false is returned. */ bool IntersectWithLine (const Base::Vector3f &rclPt, const Base::Vector3f &rclDir, Base::Vector3f &rclRes) const; + /** Calculates the intersection point of the plane defined by the base \a rclPt and the direction \a rclDir + * with the edge. The intersection must be inside the edge. If there is no intersection false is returned. + */ + bool IntersectWithPlane (const Base::Vector3f &rclPt, const Base::Vector3f &rclDir, Base::Vector3f &rclRes) const; /** * Calculates the projection of a point onto the line defined by the edge. The caller must check if * the projection point is inside the edge. diff --git a/src/Mod/Mesh/App/Core/TrimByPlane.cpp b/src/Mod/Mesh/App/Core/TrimByPlane.cpp new file mode 100644 index 0000000000..49b81d11b0 --- /dev/null +++ b/src/Mod/Mesh/App/Core/TrimByPlane.cpp @@ -0,0 +1,168 @@ +/*************************************************************************** + * Copyright (c) 2019 Werner Mayer * + * * + * This file is part of the FreeCAD CAx development system. * + * * + * This library is free software; you can redistribute it and/or * + * modify it under the terms of the GNU Library General Public * + * License as published by the Free Software Foundation; either * + * version 2 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 Library General Public License for more details. * + * * + * You should have received a copy of the GNU Library General Public * + * License along with this library; see the file COPYING.LIB. If not, * + * write to the Free Software Foundation, Inc., 59 Temple Place, * + * Suite 330, Boston, MA 02111-1307, USA * + * * + ***************************************************************************/ + +#include "PreCompiled.h" +#include + +#include "TrimByPlane.h" +#include "Grid.h" +#include "Iterator.h" + +using namespace MeshCore; + +MeshTrimByPlane::MeshTrimByPlane(MeshKernel &rclM) + : myMesh(rclM) +{ +} + +MeshTrimByPlane::~MeshTrimByPlane() +{ +} + +void MeshTrimByPlane::CheckFacets(const MeshFacetGrid& rclGrid, const Base::Vector3f& base, const Base::Vector3f& normal, + std::vector &trimFacets, std::vector& removeFacets) const +{ + // Go through the grid and check for each cell if its bounding box intersects the plane. + // If the box is completely below the plane all facets will be kept, if it's above the + // plane all triangles will be removed. + std::vector checkElements; + MeshGridIterator clGridIter(rclGrid); + for (clGridIter.Init(); clGridIter.More(); clGridIter.Next()) { + Base::BoundBox3f clBBox3d = clGridIter.GetBoundBox(); + if (clBBox3d.IsCutPlane(base, normal)) { + // save all elements in checkElements + clGridIter.GetElements(checkElements); + } + else if (clBBox3d.CalcPoint(0).DistanceToPlane(base, normal) > 0.0f) { + // save all elements in removeFacets + clGridIter.GetElements(removeFacets); + } + } + + // remove double elements + std::sort(checkElements.begin(), checkElements.end()); + checkElements.erase(std::unique(checkElements.begin(), checkElements.end()), checkElements.end()); + + trimFacets.reserve(checkElements.size()/2); // reserve some memory + for (auto it = checkElements.begin(); it != checkElements.end(); ++it) { + MeshGeomFacet clFacet = myMesh.GetFacet(*it); + if (clFacet.IntersectWithPlane(base, normal)) { + trimFacets.push_back(*it); + removeFacets.push_back(*it); + } + else if (clFacet._aclPoints[0].DistanceToPlane(base, normal) > 0.0f) { + removeFacets.push_back(*it); + } + } + + // remove double elements + std::sort(removeFacets.begin(), removeFacets.end()); + removeFacets.erase(std::unique(removeFacets.begin(), removeFacets.end()), removeFacets.end()); +} + +void MeshTrimByPlane::CreateOneFacet(const Base::Vector3f& base, const Base::Vector3f& normal, unsigned short shift, + const MeshGeomFacet& facet, std::vector& trimmedFacets) const +{ + unsigned short nul = shift % 3; + unsigned short one = (shift + 1) % 3; + unsigned short two = (shift + 2) % 3; + + Base::Vector3f p1, p2; + MeshGeomEdge edge; + + edge._aclPoints[0] = facet._aclPoints[nul]; + edge._aclPoints[1] = facet._aclPoints[one]; + edge.IntersectWithPlane(base, normal, p1); + + edge._aclPoints[0] = facet._aclPoints[nul]; + edge._aclPoints[1] = facet._aclPoints[two]; + edge.IntersectWithPlane(base, normal, p2); + + MeshGeomFacet create; + create._aclPoints[0] = facet._aclPoints[nul]; + create._aclPoints[1] = p1; + create._aclPoints[2] = p2; + trimmedFacets.push_back(create); +} + +void MeshTrimByPlane::CreateTwoFacet(const Base::Vector3f& base, const Base::Vector3f& normal, unsigned short shift, + const MeshGeomFacet& facet, std::vector& trimmedFacets) const +{ + unsigned short nul = shift % 3; + unsigned short one = (shift + 1) % 3; + unsigned short two = (shift + 2) % 3; + + Base::Vector3f p1, p2; + MeshGeomEdge edge; + + edge._aclPoints[0] = facet._aclPoints[nul]; + edge._aclPoints[1] = facet._aclPoints[two]; + edge.IntersectWithPlane(base, normal, p1); + + edge._aclPoints[0] = facet._aclPoints[one]; + edge._aclPoints[1] = facet._aclPoints[two]; + edge.IntersectWithPlane(base, normal, p2); + + MeshGeomFacet create; + create._aclPoints[0] = facet._aclPoints[nul]; + create._aclPoints[1] = facet._aclPoints[one]; + create._aclPoints[2] = p1; + trimmedFacets.push_back(create); + + create._aclPoints[0] = facet._aclPoints[one]; + create._aclPoints[1] = p2; + create._aclPoints[2] = p1; + trimmedFacets.push_back(create); +} + +void MeshTrimByPlane::TrimFacets(const std::vector& trimFacets, const Base::Vector3f& base, + const Base::Vector3f& normal, std::vector& trimmedFacets) +{ + trimmedFacets.reserve(2 * trimFacets.size()); + for (auto it = trimFacets.begin(); it != trimFacets.end(); ++it) { + MeshGeomFacet facet = myMesh.GetFacet(*it); + float dist1 = facet._aclPoints[0].DistanceToPlane(base, normal); + float dist2 = facet._aclPoints[1].DistanceToPlane(base, normal); + float dist3 = facet._aclPoints[2].DistanceToPlane(base, normal); + + // only one point below + if (dist1 < 0.0f && dist2 > 0.0f && dist3 > 0.0f) { + CreateOneFacet(base, normal, 0, facet, trimmedFacets); + } + else if (dist1 > 0.0f && dist2 < 0.0f && dist3 > 0.0f) { + CreateOneFacet(base, normal, 1, facet, trimmedFacets); + } + else if (dist1 > 0.0f && dist2 > 0.0f && dist3 < 0.0f) { + CreateOneFacet(base, normal, 2, facet, trimmedFacets); + } + // two points below + else if (dist1 < 0.0f && dist2 < 0.0f && dist3 > 0.0f) { + CreateTwoFacet(base, normal, 0, facet, trimmedFacets); + } + else if (dist1 > 0.0f && dist2 < 0.0f && dist3 < 0.0f) { + CreateTwoFacet(base, normal, 1, facet, trimmedFacets); + } + else if (dist1 < 0.0f && dist2 > 0.0f && dist3 < 0.0f) { + CreateTwoFacet(base, normal, 2, facet, trimmedFacets); + } + } +} diff --git a/src/Mod/Mesh/App/Core/TrimByPlane.h b/src/Mod/Mesh/App/Core/TrimByPlane.h new file mode 100644 index 0000000000..0da728355f --- /dev/null +++ b/src/Mod/Mesh/App/Core/TrimByPlane.h @@ -0,0 +1,66 @@ +/*************************************************************************** + * Copyright (c) 2019 Werner Mayer * + * * + * This file is part of the FreeCAD CAx development system. * + * * + * This library is free software; you can redistribute it and/or * + * modify it under the terms of the GNU Library General Public * + * License as published by the Free Software Foundation; either * + * version 2 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 Library General Public License for more details. * + * * + * You should have received a copy of the GNU Library General Public * + * License along with this library; see the file COPYING.LIB. If not, * + * write to the Free Software Foundation, Inc., 59 Temple Place, * + * Suite 330, Boston, MA 02111-1307, USA * + * * + ***************************************************************************/ + +#ifndef MESHTRIM_BY_PLANE_H +#define MESHTRIM_BY_PLANE_H + +#include +#include + +namespace MeshCore { + +/** + * Trim the the facets in 3D with a plane + * \author Werner Mayer + */ +class MeshExport MeshTrimByPlane +{ +public: + MeshTrimByPlane(MeshKernel& mesh); + ~MeshTrimByPlane(); + +public: + /** + * Checks all facets for intersection with the plane and writes all touched facets into the vector + */ + void CheckFacets(const MeshFacetGrid& rclGrid, const Base::Vector3f& base, const Base::Vector3f& normal, + std::vector& trimFacets, std::vector& removeFacets) const; + + /** + * The facets from \a trimFacets will be trimmed or deleted and \a trimmedFacets holds the newly generated facets + */ + void TrimFacets(const std::vector& trimFacets, const Base::Vector3f& base, + const Base::Vector3f& normal, std::vector& trimmedFacets); + +private: + void CreateOneFacet(const Base::Vector3f& base, const Base::Vector3f& normal, unsigned short shift, + const MeshGeomFacet& facet, std::vector& trimmedFacets) const; + void CreateTwoFacet(const Base::Vector3f& base, const Base::Vector3f& normal, unsigned short shift, + const MeshGeomFacet& facet, std::vector& trimmedFacets) const; + +private: + MeshKernel& myMesh; +}; + +} //namespace MeshCore + +#endif //MESHTRIM_BY_PLANE_H diff --git a/src/Mod/Mesh/App/Mesh.cpp b/src/Mod/Mesh/App/Mesh.cpp index 5b70871f23..842e259825 100644 --- a/src/Mod/Mesh/App/Mesh.cpp +++ b/src/Mod/Mesh/App/Mesh.cpp @@ -50,6 +50,7 @@ #include "Core/SetOperations.h" #include "Core/Triangulation.h" #include "Core/Trim.h" +#include "Core/TrimByPlane.h" #include "Core/Visitor.h" #include "Core/Decimation.h" @@ -60,7 +61,7 @@ using namespace Mesh; float MeshObject::Epsilon = 1.0e-5f; -TYPESYSTEM_SOURCE(Mesh::MeshObject, Data::ComplexGeoData); +TYPESYSTEM_SOURCE(Mesh::MeshObject, Data::ComplexGeoData) MeshObject::MeshObject() { @@ -1053,6 +1054,21 @@ void MeshObject::trim(const Base::Polygon2d& polygon2d, this->_kernel.AddFacets(triangle); } +void MeshObject::trim(const Base::Vector3f& base, const Base::Vector3f& normal) +{ + MeshCore::MeshTrimByPlane trim(this->_kernel); + std::vector trimFacets, removeFacets; + std::vector triangle; + + MeshCore::MeshFacetGrid meshGrid(this->_kernel); + trim.CheckFacets(meshGrid, base, normal, trimFacets, removeFacets); + trim.TrimFacets(trimFacets, base, normal, triangle); + if (!removeFacets.empty()) + this->deleteFacets(removeFacets); + if (!triangle.empty()) + this->_kernel.AddFacets(triangle); +} + MeshObject* MeshObject::unite(const MeshObject& mesh) const { MeshCore::MeshKernel result; diff --git a/src/Mod/Mesh/App/Mesh.h b/src/Mod/Mesh/App/Mesh.h index cd151d9ed0..0f0e4edc75 100644 --- a/src/Mod/Mesh/App/Mesh.h +++ b/src/Mod/Mesh/App/Mesh.h @@ -223,6 +223,7 @@ public: float fMinEps = 1.0e-2f, bool bConnectPolygons = false) const; void cut(const Base::Polygon2d& polygon, const Base::ViewProjMethod& proj, CutType); void trim(const Base::Polygon2d& polygon, const Base::ViewProjMethod& proj, CutType); + void trim(const Base::Vector3f& base, const Base::Vector3f& normal); //@} /** @name Selection */ diff --git a/src/Mod/MeshPart/Gui/Command.cpp b/src/Mod/MeshPart/Gui/Command.cpp index 5d0af4d1b0..3435b67599 100644 --- a/src/Mod/MeshPart/Gui/Command.cpp +++ b/src/Mod/MeshPart/Gui/Command.cpp @@ -101,8 +101,8 @@ void CmdMeshPartTrimByPlane::activated(int) msgBox.setIcon(QMessageBox::Question); msgBox.setWindowTitle(qApp->translate("MeshPart_TrimByPlane","Trim by plane")); msgBox.setText(qApp->translate("MeshPart_TrimByPlane","Select the side you want to keep.")); - QPushButton* inner = msgBox.addButton(qApp->translate("MeshPart_TrimByPlane","Inner"), QMessageBox::ActionRole); - QPushButton* outer = msgBox.addButton(qApp->translate("MeshPart_TrimByPlane","Outer"), QMessageBox::ActionRole); + QPushButton* inner = msgBox.addButton(qApp->translate("MeshPart_TrimByPlane","Below"), QMessageBox::ActionRole); + QPushButton* outer = msgBox.addButton(qApp->translate("MeshPart_TrimByPlane","Above"), QMessageBox::ActionRole); QPushButton* split = msgBox.addButton(qApp->translate("MeshPart_TrimByPlane","Split"), QMessageBox::ActionRole); msgBox.setDefaultButton(inner); msgBox.exec(); @@ -123,65 +123,41 @@ void CmdMeshPartTrimByPlane::activated(int) return; } - Base::Placement plm = static_cast(plane.front())->Placement.getValue(); - Base::Vector3d normal(0,0,1); - plm.getRotation().multVec(normal, normal); - Base::Vector3d up(-1,0,0); - plm.getRotation().multVec(up, up); - Base::Vector3d view(0,1,0); - plm.getRotation().multVec(view, view); - - Base::Vector3d base = plm.getPosition(); - Base::Rotation rot(Base::Vector3d(0,0,1), view); - Base::Matrix4D mat; - rot.getValue(mat); - Base::ViewProjMatrix proj(mat); + Base::Placement plnPlacement = static_cast(plane.front())->Placement.getValue(); openCommand("Trim with plane"); std::vector docObj = Gui::Selection().getObjectsOfType(Mesh::Feature::getClassTypeId()); for (std::vector::iterator it = docObj.begin(); it != docObj.end(); ++it) { + Base::Vector3d normal(0,0,1); + plnPlacement.getRotation().multVec(normal, normal); + Base::Vector3d base = plnPlacement.getPosition(); + Mesh::MeshObject* mesh = static_cast(*it)->Mesh.startEditing(); - Base::BoundBox3d bbox = mesh->getBoundBox(); - double len = bbox.CalcDiagonalLength(); - // project center of bbox onto plane and use this as base point - Base::Vector3d cnt = bbox.GetCenter(); - double dist = (cnt-base)*normal; - base = cnt - normal * dist; - proj.setTransform(Base::Matrix4D()); + // Apply the inverted mesh placement to the plane because the trimming is done + // on the untransformed mesh data + Base::Placement meshPlacement = mesh->getPlacement(); + meshPlacement.invert(); + meshPlacement.multVec(base, base); + meshPlacement.getRotation().multVec(normal, normal); - Base::Vector3d p1 = base + up * len; - Base::Vector3d p2 = base - up * len; - Base::Vector3d p3 = p2 + normal * len; - Base::Vector3d p4 = p1 + normal * len; - p1 = proj(p1); - p2 = proj(p2); - p3 = proj(p3); - p4 = proj(p4); - - // must be set after getting the transformed polygon points - proj.setTransform(mesh->getTransform()); - - Base::Polygon2d polygon2d; - polygon2d.Add(Base::Vector2d(p1.x, p1.y)); - polygon2d.Add(Base::Vector2d(p2.x, p2.y)); - polygon2d.Add(Base::Vector2d(p3.x, p3.y)); - polygon2d.Add(Base::Vector2d(p4.x, p4.y)); + Base::Vector3f plnBase = Base::convertTo(base); + Base::Vector3f plnNormal = Base::convertTo(normal); if (role == Gui::SelectionRole::Inner) { - mesh->trim(polygon2d, proj, Mesh::MeshObject::INNER); + mesh->trim(plnBase, plnNormal); static_cast(*it)->Mesh.finishEditing(); } else if (role == Gui::SelectionRole::Outer) { - mesh->trim(polygon2d, proj, Mesh::MeshObject::OUTER); + mesh->trim(plnBase, -plnNormal); static_cast(*it)->Mesh.finishEditing(); } else if (role == Gui::SelectionRole::Split) { Mesh::MeshObject copy(*mesh); - mesh->trim(polygon2d, proj, Mesh::MeshObject::INNER); + mesh->trim(plnBase, plnNormal); static_cast(*it)->Mesh.finishEditing(); - copy.trim(polygon2d, proj, Mesh::MeshObject::OUTER); + copy.trim(plnBase, -plnNormal); App::Document* doc = (*it)->getDocument(); Mesh::Feature* fea = static_cast(doc->addObject("Mesh::Feature")); fea->Label.setValue((*it)->Label.getValue());