kindred/create
1
1
Fork 0
Code Issues 44 Pull Requests Actions 1 Packages Projects 9 Releases 2 Wiki Activity

[TD]use Ptolemy's Theorem for bspline to circle conversion

Browse Source
This commit is contained in:
wandererfan
2025-03-21 08:45:48 -04:00
parent dc04f68732
commit 7d511b2bd9
2 changed files with 144 additions and 99 deletions
Download Patch File Download Diff File Expand all files Collapse all files

209
src/Mod/TechDraw/App/Geometry.cpp
View File

@@ -1493,87 +1493,124 @@ TopoDS_Edge GeometryUtils::edgeFromCircleArc(TechDraw::AOCPtr c)
}
//used by DVDim for approximate dims
bool GeometryUtils::isCircle(TopoDS_Edge occEdge)
bool GeometryUtils::isCircle(const TopoDS_Edge& occEdge)
{
double radius;
double radius{0};
Base::Vector3d center;
bool isArc = false;
return GeometryUtils::getCircleParms(occEdge, radius, center, isArc);
}
//! tries to interpret a B-spline edge as a circle. Used by DVDim for approximate dimensions.
//! calculates the curvature of the spline at a number of places and measures the deviation from the average
//! a true circle has constant curvature and would have no deviation from the average.
bool GeometryUtils::getCircleParms(TopoDS_Edge occEdge, double& radius, Base::Vector3d& center, bool& isArc)
//! calculates the radius and center of circles using groups of 4 points on the b-spline. if the
//! groups of 4 points all lie on a circle, we use that circle to get the radius and center.
bool GeometryUtils::getCircleParms(const TopoDS_Edge& occEdge, double& radius, Base::Vector3d& center, bool& isArc)
{
int testCount = 5;
double curveLimit = EWTOLERANCE;
BRepAdaptor_Curve c(occEdge);
Handle(Geom_BSplineCurve) spline = c.BSpline();
double f, l;
f = c.FirstParameter();
l = c.LastParameter();
double parmRange = fabs(l - f);
double parmStep = parmRange/testCount;
std::vector<double> curvatures;
std::vector<gp_Pnt> centers;
gp_Pnt curveCenter;
double sumCurvature = 0;
Base::Vector3d sumCenter, valueAt;
try {
GeomLProp_CLProps prop(spline, f, 3, Precision::Confusion());
constexpr int PointCount{8}; // number of points on the edge to examine (>= 8)
constexpr int TestCount{3}; // number of candidate circles to test
// check only the interior points of the edge
for (int i = 1; i < (testCount - 1); i++) {
prop.SetParameter(parmStep * i);
curvatures.push_back(prop.Curvature());
sumCurvature += prop.Curvature();
prop.CentreOfCurvature(curveCenter);
centers.push_back(curveCenter);
sumCenter += Base::convertTo<Base::Vector3d>(curveCenter);
}
BRepAdaptor_Curve curveAdapt(occEdge);
double firstParam = curveAdapt.FirstParameter();
auto firstPoint = Base::convertTo<Base::Vector3d>(curveAdapt.Value(firstParam));
double lastParam = curveAdapt.LastParameter();
auto lastPoint = Base::convertTo<Base::Vector3d>(curveAdapt.Value(lastParam));
double parmRange = fabs(lastParam - firstParam);
double parmStep = parmRange / PointCount;
std::vector<Base::Vector3d> pointsOnCurve;
for (size_t iPoint = 0; iPoint < PointCount; iPoint++) {
auto iPointMath = static_cast<double>(iPoint);
auto newpoint = curveAdapt.Value(firstParam + iPointMath * parmStep);
pointsOnCurve.push_back(Base::convertTo<Base::Vector3d>(newpoint));
}
catch (Standard_Failure&) {
auto edgeLong = edgeLength(occEdge);
constexpr double LimitFactor{0.001}; // 0.1% not sure about this value
double tolerance = edgeLong * LimitFactor;
isArc = true;
if (firstPoint.IsEqual(lastPoint, tolerance)) {
isArc = false;
}
int passCount{0};
int firstIndex{0};
for (int iTest = 0; iTest < TestCount; iTest++) {
firstIndex++;
auto A = pointsOnCurve.at(firstIndex);
auto B = pointsOnCurve.at(firstIndex + 1);
auto C = pointsOnCurve.at(firstIndex + 2);
auto D = pointsOnCurve.at(firstIndex + 3);
if (pointsAreOnCircle(A, B, C, D, tolerance)) {
passCount++;
}
}
if (passCount != TestCount) {
// at least 1 test failed.
return false;
}
Base::Vector3d avgCenter = sumCenter/ centers.size();
double avgCurve = sumCurvature/ centers.size();
double errorCurve = 0;
// sum the errors in curvature
for (auto& cv: curvatures) {
errorCurve += avgCurve - cv;
// each group of 4 points lies on a circle. Since the groups of 4 overlap, all the points lie
// on the same circle. https://en.wikipedia.org/wiki/Ptolemy%27s_theorem and
// https://math.stackexchange.com/questions/3130053/how-to-check-if-a-set-of-points-in-cartesian-space-could-lie-on-the-circumferenc
// so we can use any three points to make our circle.
auto gPoint0 = Base::convertTo<gp_Pnt>(pointsOnCurve.at(1));
auto gPoint1 = Base::convertTo<gp_Pnt>(pointsOnCurve.at(3));
auto gPoint2 = Base::convertTo<gp_Pnt>(pointsOnCurve.at(5)); //NOLINT readability-magic-numbers
try {
GC_MakeCircle mkCircle(gPoint0, gPoint1, gPoint2);
if (!mkCircle.IsDone()) {
return false;
}
const Handle(Geom_Circle) circleFromParms = mkCircle.Value();
radius = circleFromParms->Circ().Radius();
center = Base::convertTo<Base::Vector3d>(circleFromParms->Circ().Location());
return true;
}
catch (Standard_Failure& err) {
Base::Console().Message("Geo::getCircleParms - failed to make a circle\n");
}
double errorCenter{0};
for (auto& observe : centers) {
auto error = (Base::convertTo<Base::Vector3d>(observe)- avgCenter).Length();
errorCenter += error;
return false;
}
//! returns true if the A, B, C and D all lie on the same circle according to Ptolemy's theorem
//! we can skip the test for same plane, since the points are all on the XY plane(?not true
//! for 3d dims?).
bool GeometryUtils::pointsAreOnCircle(Base::Vector3d A,
Base::Vector3d B,
Base::Vector3d C,
Base::Vector3d D,
double tolerance)
{
auto AB = (B-A).Length();
auto AC = (C-A).Length();
auto AD = (D-A).Length();
auto BC = (C-B).Length();
auto BD = (D-B).Length();
auto CD = (D-C).Length();
auto pieceLength = AB + BC + CD;
auto wholeLength = AD;
if (DU::fpCompare(pieceLength, wholeLength, tolerance)) {
// these points are colinear
return false;
}
// calculate average error in curvature. we are only interested in the magnitude of the error
errorCurve = fabs(errorCurve / curvatures.size());
// calculate the average error in center of curvature
errorCenter = errorCenter / curvatures.size();
auto edgeLong = edgeLength(occEdge);
double centerLimit = edgeLong * 0.01;
isArc = !c.IsClosed();
bool isCircle(false);
if ( errorCurve <= curveLimit &&
errorCenter <= centerLimit) {
isCircle = true;
radius = 1.0/avgCurve;
center = avgCenter;
}
return isCircle;
bool eq1 = DU::fpCompare(AB*CD + AC*BD, AD*BC, tolerance);
bool eq2 = DU::fpCompare(AB*CD + AD*BC, AC*BD, tolerance);
bool eq3 = DU::fpCompare(AC*BD + AD*BC, AB*CD, tolerance);
return eq1 || eq2 || eq3;
}
//! make a circle or arc of circle Edge from BSpline Edge
// Note that the input edge has been inverted by GeometryObject, so +Y points down.
TopoDS_Edge GeometryUtils::asCircle(TopoDS_Edge splineEdge, bool& arc)
TopoDS_Edge GeometryUtils::asCircle(const TopoDS_Edge& splineEdge, bool& arc)
{
double radius{0};
Base::Vector3d center;
@@ -1582,11 +1619,16 @@ TopoDS_Edge GeometryUtils::asCircle(TopoDS_Edge splineEdge, bool& arc)
if (!canMakeCircle) {
throw Base::RuntimeError("GU::asCircle received non-circular edge!");
}
arc = isArc;
gp_Pnt gCenter = Base::convertTo<gp_Pnt>(center);
gp_Dir gNormal{0, 0, 1};
Handle(Geom_Circle) circleFromParms = GC_MakeCircle(gCenter, gNormal, radius);
if (!isArc) {
return BRepBuilderAPI_MakeEdge(circleFromParms);
}
// find the ends of the edge from the underlying curve
BRepAdaptor_Curve curveAdapt(splineEdge);
double firstParam = curveAdapt.FirstParameter();
@@ -1594,40 +1636,42 @@ TopoDS_Edge GeometryUtils::asCircle(TopoDS_Edge splineEdge, bool& arc)
gp_Pnt startPoint = curveAdapt.Value(firstParam);
gp_Pnt endPoint = curveAdapt.Value(lastParam);
if (startPoint.IsEqual(endPoint, EWTOLERANCE)) { //more reliable than IsClosed flag
arc = false;
return BRepBuilderAPI_MakeEdge(circleFromParms);
}
arc = true;
double midRange = (lastParam + firstParam) / 2;
gp_Pnt midPoint = curveAdapt.Value(midRange);
// this should be using circleFromParms as a base instead of points on the original spline.
// could be problems with very small errors?? other versions of GC_MakeArcOfCircle have
// poorly explained parameters.
GC_MakeArcOfCircle mkArc(startPoint, midPoint, endPoint);
auto circleArc = mkArc.Value();
if (!mkArc.IsDone()) {
throw Base::RuntimeError("GU::asCircle failed to create arc");
}
return BRepBuilderAPI_MakeEdge(circleArc);
}
bool GeometryUtils::isLine(TopoDS_Edge occEdge)
bool GeometryUtils::isLine(const TopoDS_Edge& occEdge)
{
BRepAdaptor_Curve c(occEdge);
BRepAdaptor_Curve adapt(occEdge);
Handle(Geom_BSplineCurve) spline = c.BSpline();
double f = c.FirstParameter();
double l = c.LastParameter();
gp_Pnt s = c.Value(f);
gp_Pnt e = c.Value(l);
Handle(Geom_BSplineCurve) spline = adapt.BSpline();
double firstParm = adapt.FirstParameter();
double lastParm = adapt.LastParameter();
auto startPoint = Base::convertTo<Base::Vector3d>(adapt.Value(firstParm));
auto endPoint = Base::convertTo<Base::Vector3d>(adapt.Value(lastParm));
auto edgeLong = edgeLength(occEdge);
bool samePnt = s.IsEqual(e, FLT_EPSILON);
if (samePnt) {
constexpr double LimitFactor{0.001}; // 0.1% not sure about this value
double tolerance = edgeLong * LimitFactor;
if (startPoint.IsEqual(endPoint, tolerance)) {
// either not a line or a zero length line?
return false;
}
Base::Vector3d vs = Base::convertTo<Base::Vector3d>(s);
Base::Vector3d ve = Base::convertTo<Base::Vector3d>(e);
double endLength = (vs - ve).Length();
// in a line the sum of the lengths of the segments should equal the distance
// from start to end
double endPointLength = (endPoint - startPoint).Length();
int low = 0;
int high = spline->NbPoles() - 1;
TColgp_Array1OfPnt poles(low, high);
@@ -1641,15 +1685,12 @@ bool GeometryUtils::isLine(TopoDS_Edge occEdge)
lenTotal += (v2-v1).Length();
}
if (DrawUtil::fpCompare(lenTotal, endLength)) {
return true;
}
return false;
return DrawUtil::fpCompare(lenTotal, endPointLength, tolerance);
}
//! make a line Edge from B-spline Edge
TopoDS_Edge GeometryUtils::asLine(TopoDS_Edge occEdge)
TopoDS_Edge GeometryUtils::asLine(const TopoDS_Edge& occEdge)
{
BRepAdaptor_Curve c(occEdge);

34
src/Mod/TechDraw/App/Geometry.h
View File

@@ -315,7 +315,6 @@ class TechDrawExport BSpline: public BaseGeom
bool isLine();
bool isCircle();
TopoDS_Edge asCircle(bool& isArc);
// void getCircleParms(bool& isCircle, double& radius, Base::Vector3d& center, bool& isArc);
bool intersectsArc(Base::Vector3d p1, Base::Vector3d p2);
std::vector<BezierSegment> segments;
};
@@ -382,23 +381,23 @@ class TechDrawExport Vertex : public TechDraw::Tag
Base::Vector3d point() const { return Base::Vector3d(pnt.x, pnt.y, 0.0); }
void point(Base::Vector3d v){ pnt = Base::Vector3d(v.x, v.y); }
double x() {return pnt.x;}
double y() {return pnt.y;}
double x() const {return pnt.x;}
double y() const {return pnt.y;}
// attribute setters and getters
bool getHlrVisible() { return hlrVisible; }
bool getHlrVisible() const { return hlrVisible; }
void setHlrVisible(bool state) { hlrVisible = state; }
int getRef3d() { return ref3D; }
int getRef3d() const { return ref3D; }
void setRef3d(int ref) { ref3D = ref; }
TopoDS_Vertex getOCCVertex() { return occVertex; }
void setOCCVertex(TopoDS_Vertex newVertex) { occVertex = newVertex; }
bool getCosmetic() { return cosmetic; }
void setOCCVertex(const TopoDS_Vertex& newVertex) { occVertex = newVertex; }
bool getCosmetic() const { return cosmetic; }
void setCosmetic (bool state) { cosmetic = state; }
std::string getCosmeticTag() { return cosmeticTag; }
void setCosmeticTag(std::string t) { cosmeticTag = t; }
bool isCenter() {return m_center;}
void setCosmeticTag(const std::string& t) { cosmeticTag = t; }
bool isCenter() const {return m_center;}
void isCenter(bool state) { m_center = state; }
bool isReference() { return m_reference; }
bool isReference() const { return m_reference; }
void isReference(bool state) { m_reference = state; }
Part::TopoShape asTopoShape(double scale = 1.0);
@@ -447,17 +446,22 @@ class TechDrawExport GeometryUtils
static TopoDS_Edge edgeFromCircle(TechDraw::CirclePtr c);
static TopoDS_Edge edgeFromCircleArc(TechDraw::AOCPtr c);
static bool isCircle(TopoDS_Edge occEdge);
static bool getCircleParms(TopoDS_Edge occEdge, double& radius, Base::Vector3d& center, bool& isArc);
static TopoDS_Edge asCircle(TopoDS_Edge splineEdge, bool& arc);
static bool isLine(TopoDS_Edge occEdge);
static TopoDS_Edge asLine(TopoDS_Edge occEdge);
static bool isCircle(const TopoDS_Edge& occEdge);
static bool getCircleParms(const TopoDS_Edge& occEdge, double& radius, Base::Vector3d& center, bool& isArc);
static TopoDS_Edge asCircle(const TopoDS_Edge& splineEdge, bool& arc);
static bool isLine(const TopoDS_Edge& occEdge);
static TopoDS_Edge asLine(const TopoDS_Edge& occEdge);
static double edgeLength(TopoDS_Edge occEdge);
static TopoDS_Face makePerforatedFace(FacePtr bigCheese, const std::vector<FacePtr>& holesAll);
static std::vector<FacePtr> findHolesInFace(const DrawViewPart* dvp, const std::string& bigCheeseSubRef);
static bool pointsAreOnCircle(Base::Vector3d A,
Base::Vector3d B,
Base::Vector3d C,
Base::Vector3d D,
double tolerance);
};