+ integrate biarcs algorithm

This commit is contained in:
wmayer 2015-04-05 14:10:03 +02:00
parent d3af1ea85a
commit 1d1e71b52f
8 changed files with 366 additions and 0 deletions

View File

@ -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 <gp_Ax2.hxx>
# include <gp_Pln.hxx>
# include <gp_Circ.hxx>
# include <GC_MakeArcOfCircle.hxx>
# include <Geom_BSplineCurve.hxx>
#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<Base::Vector3d>(m_p0),Base::convertTo<Base::Vector3d>(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<Geometry*>& 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<Base::Vector3d>(p_start),Base::convertTo<Base::Vector3d>(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<Geometry*> 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<Geometry*> list;
createArcs(tolerance, list, p_start, v_start, this->myCurve->FirstParameter(), this->myCurve->LastParameter(), p_end, v_end);
return list;
}

View File

@ -307,6 +307,14 @@ from the knots table of this B-Spline curve.</UserDocu>
</UserDocu>
</Documentation>
</Methode>
<Methode Name="toBiArcs">
<Documentation>
<UserDocu>
Build a list of arcs and lines to approximate the b-spline.
toBiArcs(tolerance) -> list.
</UserDocu>
</Documentation>
</Methode>
<Methode Name="join">
<Documentation>
<UserDocu>

View File

@ -558,6 +558,7 @@ PyObject* BSplineCurvePy::setPeriodic(PyObject * args)
{
if (!PyArg_ParseTuple(args, ""))
return 0;
std::list<Geometry*> 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<Geometry*> arcs;
arcs = curve->toBiArcs(tolerance);
Py::List list;
for (std::list<Geometry*>::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;

View File

@ -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

View File

@ -54,6 +54,7 @@
# include <GeomLProp_CLProps.hxx>
# include <GeomLProp_SLProps.hxx>
# include <gp.hxx>
# include <gp_Ax2.hxx>
# include <gp_Circ.hxx>
# include <gp_Elips.hxx>
# include <gp_Hypr.hxx>

View File

@ -49,6 +49,7 @@
#include <TopoDS_Shape.hxx>
#include <gp_Ax1.hxx>
#include <gp_Dir.hxx>
#include <list>
#include <Base/Persistence.h>
#include <Base/Vector3D.h>
@ -160,6 +161,7 @@ public:
std::vector<Base::Vector3d> getPoles() const;
bool join(const Handle_Geom_BSplineCurve&);
void makeC1Continuous(double, double);
std::list<Geometry*> 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<Geometry*>& 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;
};

View File

@ -0,0 +1,83 @@
/***************************************************************************
* Copyright (c) 2011 Werner Mayer <wmayer[at]users.sourceforge.net> *
* *
* 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 <gp_Pln.hxx>
# include <gp_Lin.hxx>
# include <Geom_Plane.hxx>
# include <GeomAPI_IntSS.hxx>
# include <Geom_Line.hxx>
# include <Precision.hxx>
#endif
#include <Base/Vector3D.h>
#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;
}

View File

@ -29,6 +29,9 @@
#include <gp_Dir.hxx>
#include <gp_XYZ.hxx>
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