From 1d1e71b52f53ab598ecd785d5fd1c85d070fccc6 Mon Sep 17 00:00:00 2001 From: wmayer Date: Sun, 5 Apr 2015 14:10:03 +0200 Subject: [PATCH] + integrate biarcs algorithm --- src/Mod/Part/App/BSplineCurveBiArcs.cpp | 227 ++++++++++++++++++++++++ src/Mod/Part/App/BSplineCurvePy.xml | 8 + src/Mod/Part/App/BSplineCurvePyImp.cpp | 26 +++ src/Mod/Part/App/CMakeLists.txt | 2 + src/Mod/Part/App/Geometry.cpp | 1 + src/Mod/Part/App/Geometry.h | 9 + src/Mod/Part/App/Tools.cpp | 83 +++++++++ src/Mod/Part/App/Tools.h | 10 ++ 8 files changed, 366 insertions(+) create mode 100644 src/Mod/Part/App/BSplineCurveBiArcs.cpp create mode 100644 src/Mod/Part/App/Tools.cpp diff --git a/src/Mod/Part/App/BSplineCurveBiArcs.cpp b/src/Mod/Part/App/BSplineCurveBiArcs.cpp new file mode 100644 index 000000000..09c943254 --- /dev/null +++ b/src/Mod/Part/App/BSplineCurveBiArcs.cpp @@ -0,0 +1,227 @@ +// This file is released under the BSD license +// +// Copyright (c) 2009, Daniel Heeks +// All rights reserved. +// +// Redistribution and use in source and binary forms, with or without modification, +// are permitted provided that the following conditions are met: +// +// * Redistributions of source code must retain the above copyright notice, this +// list of conditions and the following disclaimer. +// * Redistributions in binary form must reproduce the above copyright notice, this +// list of conditions and the following disclaimer in the documentation and/or +// other materials provided with the distribution. +// * Neither the name of Daniel Heeks nor the names of its contributors may be used +// to endorse or promote products derived from this software without specific prior +// written permission. +// +// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +// AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE +// ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE +// LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +// CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +// SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +// INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +// CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +// ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +// POSSIBILITY OF SUCH DAMAGE. + +#include "PreCompiled.h" +#ifndef _PreComp_ +# include +# include +# include +# include +# include +#endif + +#include "Geometry.h" +#include "Tools.h" + +using Part::GeomBSplineCurve; +using Part::Geometry; + +// Algorithm taken from HeeksCAD +namespace Part { + +bool tangentialArc(const gp_Pnt& p0, const gp_Vec& v0, const gp_Pnt& p1, gp_Pnt& c, gp_Dir& axis) +{ + if (p0.Distance(p1) > Precision::Intersection() && + v0.Magnitude() > Precision::Intersection()){ + gp_Vec v1(p0, p1); + gp_Pnt halfway(p0.XYZ() + v1.XYZ() * 0.5); + gp_Pln pln1(halfway, v1); + gp_Pln pln2(p0, v0); + gp_Lin plane_line; + if (intersect(pln1, pln2, plane_line)) { + gp_Lin l1(halfway, v1); + gp_Pnt p2; + closestPointsOnLines(plane_line, l1, c, p2); + axis = -(plane_line.Direction()); + return true; + } + } + + return false; +} + +class TangentialArc +{ +public: + gp_Pnt m_p0; // start point + gp_Vec m_v0; // start direction + gp_Pnt m_p1; // end point + gp_Pnt m_c; // centre point + gp_Dir m_a; // axis + bool m_is_a_line; + TangentialArc(const gp_Pnt& p0, const gp_Vec& v0, const gp_Pnt& p1) + : m_p0(p0), m_v0(v0), m_p1(p1) + { + // calculate a tangential arc that goes through p0 and p1, with a direction of v0 at p0 + m_is_a_line = !Part::tangentialArc(m_p0, m_v0, m_p1, m_c, m_a); + } + bool isRadiusEqual(const gp_Pnt &p, double tolerance) const + { + if (m_is_a_line) + return true; + + double point_radius = gp_Vec(m_c.XYZ() - p.XYZ()).Magnitude(); + double diff = fabs(point_radius - radius()); + return diff <= tolerance; + } + double radius() const + { + double r0 = gp_Vec(m_p0.XYZ() - m_c.XYZ()).Magnitude(); + double r1 = gp_Vec(m_p1.XYZ() - m_c.XYZ()).Magnitude(); + double r = (r0 + r1)/2; + return r; + } + Geometry* makeArc() const + { + if (m_is_a_line) { + GeomLineSegment* line = new GeomLineSegment(); + line->setPoints(Base::convertTo(m_p0),Base::convertTo(m_p1)); + return line; + } + + gp_Circ c(gp_Ax2(m_c, m_a), radius()); + GC_MakeArcOfCircle arc(c, m_p0, m_p1, true); + GeomArcOfCircle* new_object = new GeomArcOfCircle(); + new_object->setHandle(arc.Value()); + return new_object; + } +}; + +} + +void GeomBSplineCurve::createArcs(double tolerance, std::list& new_spans, + const gp_Pnt& p_start, const gp_Vec& v_start, + double t_start, double t_end, gp_Pnt& p_end, gp_Vec& v_end) const +{ + this->myCurve->D1(t_end, p_end, v_end); + + gp_Pnt p1, p2, p3; + + bool can_do_spline_whole = calculateBiArcPoints(p_start, v_start, p_end, v_end, p1, p2, p3); + + Geometry* arc_object1 = NULL; + Geometry* arc_object2 = NULL; + + if (can_do_spline_whole) { + Part::TangentialArc arc1(p_start, v_start, p2); + Part::TangentialArc arc2(p2, gp_Vec(p3.XYZ() - p2.XYZ()), p_end); + + gp_Pnt p_middle1, p_middle2; + this->myCurve->D0(t_start + ((t_end - t_start) * 0.25), p_middle1); + this->myCurve->D0(t_start + ((t_end - t_start) * 0.75), p_middle2); + + if (!arc1.isRadiusEqual(p_middle1, tolerance) || + !arc2.isRadiusEqual(p_middle2, tolerance)) { + can_do_spline_whole = false; + } + else { + arc_object1 = arc1.makeArc(); + arc_object2 = arc2.makeArc(); + } + } + else { + // calculate_biarc_points failed, just add a line + GeomLineSegment* line = new GeomLineSegment(); + line->setPoints(Base::convertTo(p_start),Base::convertTo(p_end)); + new_spans.push_back(line); + return; + } + + if (can_do_spline_whole) { + new_spans.push_back(arc_object1); + new_spans.push_back(arc_object2); + } + else { + double t_middle = t_start + ((t_end - t_start) * 0.5); + gp_Pnt p_middle; + gp_Vec v_middle; + createArcs(tolerance, new_spans, p_start, v_start, t_start, t_middle, p_middle, v_middle);// recursive + gp_Pnt new_p_end; + gp_Vec new_v_end; + createArcs(tolerance, new_spans, p_middle, v_middle, t_middle, t_end, new_p_end, new_v_end); + } +} + +bool GeomBSplineCurve::calculateBiArcPoints(const gp_Pnt& p0, gp_Vec v_start, + const gp_Pnt& p4, gp_Vec v_end, + gp_Pnt& p1, gp_Pnt& p2, gp_Pnt& p3) const +{ + if (v_start.Magnitude() < Precision::Intersection()) + v_start = gp_Vec(p0, p1); + if (v_end.Magnitude() < Precision::Intersection()) + v_end = gp_Vec(p3, p4); + + v_start.Normalize(); + v_end.Normalize(); + + gp_Vec v = p0.XYZ() - p4.XYZ(); + + double a = 2*(v_start*v_end-1); + double c = v*v; + double b = (v*2)*(v_start+v_end); + if (fabs(a) < Precision::Intersection()) + return false; + + + double d = b*b-4*a*c; + if (d < 0.0) + return false; + + double sd = sqrt(d); + double e1 = (-b - sd) / (2.0 * a); + double e2 = (-b + sd) / (2.0 * a); + if (e1 > 0 && e2 > 0) + return false; + + double e = e1; + if (e2 > e) + e = e2; + if (e < 0) + return false; + + p1 = p0.XYZ() + v_start.XYZ() * e; + p3 = p4.XYZ() - v_end.XYZ() * e; + p2 = p1.XYZ() * 0.5 + p3.XYZ() * 0.5; + + + return true; +} + +std::list GeomBSplineCurve::toBiArcs(double tolerance) const +{ + gp_Pnt p_start; + gp_Vec v_start; + gp_Pnt p_end; + gp_Vec v_end; + this->myCurve->D1(this->myCurve->FirstParameter(), p_start, v_start); + + std::list list; + createArcs(tolerance, list, p_start, v_start, this->myCurve->FirstParameter(), this->myCurve->LastParameter(), p_end, v_end); + return list; +} diff --git a/src/Mod/Part/App/BSplineCurvePy.xml b/src/Mod/Part/App/BSplineCurvePy.xml index 704b25167..6977773ee 100644 --- a/src/Mod/Part/App/BSplineCurvePy.xml +++ b/src/Mod/Part/App/BSplineCurvePy.xml @@ -307,6 +307,14 @@ from the knots table of this B-Spline curve. + + + + Build a list of arcs and lines to approximate the b-spline. + toBiArcs(tolerance) -> list. + + + diff --git a/src/Mod/Part/App/BSplineCurvePyImp.cpp b/src/Mod/Part/App/BSplineCurvePyImp.cpp index 343d853c2..f2e7d2338 100644 --- a/src/Mod/Part/App/BSplineCurvePyImp.cpp +++ b/src/Mod/Part/App/BSplineCurvePyImp.cpp @@ -558,6 +558,7 @@ PyObject* BSplineCurvePy::setPeriodic(PyObject * args) { if (!PyArg_ParseTuple(args, "")) return 0; + std::list new_spans; try { Handle_Geom_BSplineCurve curve = Handle_Geom_BSplineCurve::DownCast (getGeometryPtr()->handle()); @@ -711,6 +712,31 @@ Py::List BSplineCurvePy::getKnotSequence(void) const return list; } +PyObject* BSplineCurvePy::toBiArcs(PyObject * args) +{ + double tolerance = 0.001; + if (!PyArg_ParseTuple(args, "d", &tolerance)) + return 0; + try { + GeomBSplineCurve* curve = getGeomBSplineCurvePtr(); + std::list arcs; + arcs = curve->toBiArcs(tolerance); + + Py::List list; + for (std::list::iterator it = arcs.begin(); it != arcs.end(); ++it) { + list.append(Py::Object((*it)->getPyObject())); + delete (*it); + } + + return Py::new_reference_to(list); + } + catch (Standard_Failure) { + Handle_Standard_Failure e = Standard_Failure::Caught(); + PyErr_SetString(PartExceptionOCCError, e->GetMessageString()); + return 0; + } +} + PyObject* BSplineCurvePy::approximate(PyObject *args) { PyObject* obj; diff --git a/src/Mod/Part/App/CMakeLists.txt b/src/Mod/Part/App/CMakeLists.txt index 068e4eb88..d7c1b8b1d 100644 --- a/src/Mod/Part/App/CMakeLists.txt +++ b/src/Mod/Part/App/CMakeLists.txt @@ -234,6 +234,7 @@ SET(Part_SRCS ${Python_SRCS} AppPart.cpp AppPartPy.cpp + BSplineCurveBiArcs.cpp CrossSection.cpp CrossSection.h Geometry.cpp @@ -252,6 +253,7 @@ SET(Part_SRCS edgecluster.h modelRefine.cpp modelRefine.h + Tools.cpp Tools.h encodeFilename.h OCCError.h diff --git a/src/Mod/Part/App/Geometry.cpp b/src/Mod/Part/App/Geometry.cpp index 9b22f119a..5543682c7 100644 --- a/src/Mod/Part/App/Geometry.cpp +++ b/src/Mod/Part/App/Geometry.cpp @@ -54,6 +54,7 @@ # include # include # include +# include # include # include # include diff --git a/src/Mod/Part/App/Geometry.h b/src/Mod/Part/App/Geometry.h index 1073721bf..5ec2a71c3 100644 --- a/src/Mod/Part/App/Geometry.h +++ b/src/Mod/Part/App/Geometry.h @@ -49,6 +49,7 @@ #include #include #include +#include #include #include @@ -160,6 +161,7 @@ public: std::vector getPoles() const; bool join(const Handle_Geom_BSplineCurve&); void makeC1Continuous(double, double); + std::list toBiArcs(double tolerance) const; // Persistence implementer --------------------- virtual unsigned int getMemSize(void) const; @@ -171,6 +173,13 @@ public: void setHandle(const Handle_Geom_BSplineCurve&); const Handle_Geom_Geometry& handle() const; +private: + void createArcs(double tolerance, std::list& new_spans, + const gp_Pnt &p_start, const gp_Vec &v_start, + double t_start, double t_end, gp_Pnt &p_end, gp_Vec &v_end) const; + bool calculateBiArcPoints(const gp_Pnt& p0, gp_Vec v_start, + const gp_Pnt& p4, gp_Vec v_end, + gp_Pnt& p1, gp_Pnt& p2, gp_Pnt& p3) const; private: Handle_Geom_BSplineCurve myCurve; }; diff --git a/src/Mod/Part/App/Tools.cpp b/src/Mod/Part/App/Tools.cpp new file mode 100644 index 000000000..f16d6f795 --- /dev/null +++ b/src/Mod/Part/App/Tools.cpp @@ -0,0 +1,83 @@ +/*************************************************************************** + * Copyright (c) 2011 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" +#ifndef _PreComp_ +# include +# include +# include +# include +# include +# include +#endif + +#include +#include "Tools.h" + +void Part::closestPointsOnLines(const gp_Lin& lin1, const gp_Lin& lin2, gp_Pnt& p1, gp_Pnt& p2) +{ + // they might be the same point + gp_Vec v1(lin1.Direction()); + gp_Vec v2(lin2.Direction()); + gp_Vec v3(lin2.Location(), lin1.Location()); + + double a = v1*v1; + double b = v1*v2; + double c = v2*v2; + double d = v1*v3; + double e = v2*v3; + double D = a*c - b*b; + double s, t; + + // D = (v1 x v2) * (v1 x v2) + if (D < Precision::Angular()){ + // the lines are considered parallel + s = 0.0; + t = (b>c ? d/b : e/c); + } + else { + s = (b*e - c*d) / D; + t = (a*e - b*d) / D; + } + + p1 = lin1.Location().XYZ() + s * v1.XYZ(); + p2 = lin2.Location().XYZ() + t * v2.XYZ(); +} + +bool Part::intersect(const gp_Pln& pln1, const gp_Pln& pln2, gp_Lin& lin) +{ + bool found = false; + Handle (Geom_Plane) gp1 = new Geom_Plane(pln1); + Handle (Geom_Plane) gp2 = new Geom_Plane(pln2); + + GeomAPI_IntSS intSS(gp1, gp2, Precision::Confusion()); + if (intSS.IsDone()) { + int numSol = intSS.NbLines(); + if (numSol > 0) { + Handle_Geom_Curve curve = intSS.Line(1); + lin = Handle_Geom_Line::DownCast(curve)->Lin(); + found = true; + } + } + + return found; +} diff --git a/src/Mod/Part/App/Tools.h b/src/Mod/Part/App/Tools.h index b48ab6ab6..16bfc4dc6 100644 --- a/src/Mod/Part/App/Tools.h +++ b/src/Mod/Part/App/Tools.h @@ -29,6 +29,9 @@ #include #include +class gp_Lin; +class gp_Pln; + namespace Base { // Specialization for gp_Pnt template <> @@ -83,6 +86,13 @@ private: namespace Part { +PartExport +void closestPointsOnLines(const gp_Lin& lin1, const gp_Lin& lin2, gp_Pnt &p1, gp_Pnt &p2); +PartExport +bool intersect(const gp_Pln& pln1, const gp_Pln& pln2, gp_Lin& lin); +PartExport +bool tangentialArc(const gp_Pnt& p0, const gp_Vec& v0, const gp_Pnt& p1, gp_Pnt& c, gp_Dir& a); + } //namespace Part