From 9802b5d1ee0d450f68f58dc5f7ad4199d8333cb1 Mon Sep 17 00:00:00 2001 From: phkahler <14852918+phkahler@users.noreply.github.com> Date: Wed, 1 Jul 2020 16:16:27 -0400 Subject: [PATCH] Better helical triangulation - issue 489. Resolve issue #489 helix has stairsteps. Force helix axis line to 8 segments. Grid triangulation to use a minimum of 4 segments for degree>1. Adds twist dependence for grid triangulation with degree=1. Added a max_dt parameter for PWL creation and use that for helical edges. --- src/srf/ratpoly.cpp | 52 ++++++++++++++++++++++++----------------- src/srf/surface.cpp | 12 +++++++--- src/srf/surface.h | 15 ++++++------ src/srf/triangulate.cpp | 48 ++++++++++++++++++++++++++++++------- 4 files changed, 87 insertions(+), 40 deletions(-) diff --git a/src/srf/ratpoly.cpp b/src/srf/ratpoly.cpp index 7affff02..6498cc4a 100644 --- a/src/srf/ratpoly.cpp +++ b/src/srf/ratpoly.cpp @@ -173,18 +173,18 @@ void SBezier::SplitAt(double t, SBezier *bef, SBezier *aft) const { } } -void SBezier::MakePwlInto(SEdgeList *sel, double chordTol) const { +void SBezier::MakePwlInto(SEdgeList *sel, double chordTol, double max_dt) const { List lv = {}; - MakePwlInto(&lv, chordTol); + MakePwlInto(&lv, chordTol, max_dt); int i; for(i = 1; i < lv.n; i++) { sel->AddEdge(lv[i-1], lv[i]); } lv.Clear(); } -void SBezier::MakePwlInto(List *l, double chordTol) const { +void SBezier::MakePwlInto(List *l, double chordTol, double max_dt) const { List lv = {}; - MakePwlInto(&lv, chordTol); + MakePwlInto(&lv, chordTol, max_dt); int i; for(i = 0; i < lv.n; i++) { SCurvePt scpt; @@ -195,32 +195,42 @@ void SBezier::MakePwlInto(List *l, double chordTol) const { } lv.Clear(); } -void SBezier::MakePwlInto(SContour *sc, double chordTol) const { +void SBezier::MakePwlInto(SContour *sc, double chordTol, double max_dt) const { List lv = {}; - MakePwlInto(&lv, chordTol); + MakePwlInto(&lv, chordTol, max_dt); int i; for(i = 0; i < lv.n; i++) { sc->AddPoint(lv[i]); } lv.Clear(); } -void SBezier::MakePwlInto(List *l, double chordTol) const { +//-------------------------------------------------------------------------------------- +// all variants of MakePwlInto come here. Split a rational Bezier into Piecewise Linear +// segments that don't deviate from the actual curve by more than the chordTol distance. +// max_dt allows to force curves to be split into spans of no more than a certain +// length based on t-parameter. RemoveShortSegments() may delete points when dt <= 0.1 +//-------------------------------------------------------------------------------------- +void SBezier::MakePwlInto(List *l, double chordTol, double max_dt) const { if(EXACT(chordTol == 0)) { // Use the default chord tolerance. chordTol = SS.ChordTolMm(); } + // Never do fewer than three intermediate points for curves; people seem to get + // unhappy when their circles turn into squares, but maybe less + // unhappy with octagons. Now 16-gons. + if (EXACT(max_dt == 0.0)) { + max_dt = (deg == 1) ? 1.0 : 0.25; + } l->Add(&(ctrl[0])); - if(deg == 1) { + // don't split first degee (lines) unless asked to by the caller via max_dt + if((deg == 1) && (max_dt >= 1.0)) { l->Add(&(ctrl[1])); } else { - // Never do fewer than one intermediate point; people seem to get - // unhappy when their circles turn into squares, but maybe less - // unhappy with octagons. - MakePwlInitialWorker(l, 0.0, 0.5, chordTol); - MakePwlInitialWorker(l, 0.5, 1.0, chordTol); + MakePwlInitialWorker(l, 0.0, 0.5, chordTol, max_dt); + MakePwlInitialWorker(l, 0.5, 1.0, chordTol, max_dt); } } -void SBezier::MakePwlWorker(List *l, double ta, double tb, double chordTol) const +void SBezier::MakePwlWorker(List *l, double ta, double tb, double chordTol, double max_dt) const { Vector pa = PointAt(ta); Vector pb = PointAt(tb); @@ -229,16 +239,16 @@ void SBezier::MakePwlWorker(List *l, double ta, double tb, double chordT double d = pm.DistanceToLine(pa, pb.Minus(pa)); double step = 1.0/SS.GetMaxSegments(); - if((tb - ta) < step || d < chordTol) { + if(((tb - ta) < step || d < chordTol) && ((tb-ta) <= max_dt) ) { // A previous call has already added the beginning of our interval. l->Add(&pb); } else { double tm = (ta + tb) / 2; - MakePwlWorker(l, ta, tm, chordTol); - MakePwlWorker(l, tm, tb, chordTol); + MakePwlWorker(l, ta, tm, chordTol, max_dt); + MakePwlWorker(l, tm, tb, chordTol, max_dt); } } -void SBezier::MakePwlInitialWorker(List *l, double ta, double tb, double chordTol) const +void SBezier::MakePwlInitialWorker(List *l, double ta, double tb, double chordTol, double max_dt) const { Vector pa = PointAt(ta); Vector pb = PointAt(tb); @@ -259,13 +269,13 @@ void SBezier::MakePwlInitialWorker(List *l, double ta, double tb, double }); double step = 1.0/SS.GetMaxSegments(); - if( ((tb - ta) < step || d < chordTol) && ((tb-ta) < 0.2) ) { + if( ((tb - ta) < step || d < chordTol) && ((tb-ta) <= max_dt) ) { // A previous call has already added the beginning of our interval. l->Add(&pb); } else { double tm = (ta + tb) / 2; - MakePwlWorker(l, ta, tm, chordTol); - MakePwlWorker(l, tm, tb, chordTol); + MakePwlWorker(l, ta, tm, chordTol, max_dt); + MakePwlWorker(l, tm, tb, chordTol, max_dt); } } diff --git a/src/srf/surface.cpp b/src/srf/surface.cpp index 6afcc791..fbdc5d92 100644 --- a/src/srf/surface.cpp +++ b/src/srf/surface.cpp @@ -746,7 +746,9 @@ void SShell::MakeFromHelicalRevolutionOf(SBezierLoopSet *sbls, Vector pt, Vector sc = {}; sc.isExact = true; sc.exact = sb->TransformedBy(ts, qs, 1.0); - (sc.exact).MakePwlInto(&(sc.pts)); + double max_dt = 0.5; + if (sc.exact.deg > 1) max_dt = 0.25; + (sc.exact).MakePwlInto(&(sc.pts), 0.0, max_dt); // the surfaces already exist so trim with this curve if(j < sections) { @@ -771,7 +773,9 @@ void SShell::MakeFromHelicalRevolutionOf(SBezierLoopSet *sbls, Vector pt, Vector sc = {}; sc.isExact = true; sc.exact = sb->TransformedBy(ts, qs, 1.0); - (sc.exact).MakePwlInto(&(sc.pts)); + double max_dt = 0.5; + if (sc.exact.deg > 1) max_dt = 0.25; + (sc.exact).MakePwlInto(&(sc.pts), 0.0, max_dt); sc.surfA = hs1; // end cap sc.surfB = hs0; // staring cap hSCurve hcb = curve.AddAndAssignId(&sc); @@ -793,7 +797,9 @@ void SShell::MakeFromHelicalRevolutionOf(SBezierLoopSet *sbls, Vector pt, Vector sc.isExact = true; sc.exact = SBezier::From(ss->ctrl[0][0], ss->ctrl[0][1], ss->ctrl[0][2]); sc.exact.weight[1] = ss->weight[0][1]; - (sc.exact).MakePwlInto(&(sc.pts)); + double max_dt = 0.5; + if (sc.exact.deg > 1) max_dt = 0.125; + (sc.exact).MakePwlInto(&(sc.pts), 0.0, max_dt); sc.surfA = revs[j]; sc.surfB = revsp[j]; diff --git a/src/srf/surface.h b/src/srf/surface.h index c82af1ce..d1918b6c 100644 --- a/src/srf/surface.h +++ b/src/srf/surface.h @@ -92,12 +92,12 @@ public: Vector Start() const; Vector Finish() const; bool Equals(SBezier *b) const; - void MakePwlInto(SEdgeList *sel, double chordTol=0) const; - void MakePwlInto(List *l, double chordTol=0) const; - void MakePwlInto(SContour *sc, double chordTol=0) const; - void MakePwlInto(List *l, double chordTol=0) const; - void MakePwlWorker(List *l, double ta, double tb, double chordTol) const; - void MakePwlInitialWorker(List *l, double ta, double tb, double chordTol) const; + void MakePwlInto(SEdgeList *sel, double chordTol=0, double max_dt=0.0) const; + void MakePwlInto(List *l, double chordTol=0, double max_dt=0.0) const; + void MakePwlInto(SContour *sc, double chordTol=0, double max_dt=0.0) const; + void MakePwlInto(List *l, double chordTol=0, double max_dt=0.0) const; + void MakePwlWorker(List *l, double ta, double tb, double chordTol, double max_dt) const; + void MakePwlInitialWorker(List *l, double ta, double tb, double chordTol, double max_dt) const; void MakeNonrationalCubicInto(SBezierList *bl, double tolerance, int depth = 0) const; void AllIntersectionsWith(const SBezier *sbb, SPointList *spl) const; @@ -362,8 +362,9 @@ public: void MakeClassifyingBsp(SShell *shell, SShell *useCurvesFrom); double ChordToleranceForEdge(Vector a, Vector b) const; void MakeTriangulationGridInto(List *l, double vs, double vf, - bool swapped) const; + bool swapped, int depth) const; Vector PointAtMaybeSwapped(double u, double v, bool swapped) const; + Vector NormalAtMaybeSwapped(double u, double v, bool swapped) const; void Reverse(); void Clear(); diff --git a/src/srf/triangulate.cpp b/src/srf/triangulate.cpp index 51a31ddd..d02db791 100644 --- a/src/srf/triangulate.cpp +++ b/src/srf/triangulate.cpp @@ -382,13 +382,24 @@ Vector SSurface::PointAtMaybeSwapped(double u, double v, bool swapped) const { } } +Vector SSurface::NormalAtMaybeSwapped(double u, double v, bool swapped) const { + Vector du, dv; + if(swapped) { + TangentsAt(v, u, &dv, &du); + } else { + TangentsAt(u, v, &du, &dv); + } + return du.Cross(dv).WithMagnitude(1.0); +} + void SSurface::MakeTriangulationGridInto(List *l, double vs, double vf, - bool swapped) const + bool swapped, int depth) const { double worst = 0; // Try piecewise linearizing four curves, at u = 0, 1/3, 2/3, 1; choose // the worst chord tolerance of any of those. + double worst_twist = 1.0; int i; for(i = 0; i <= 3; i++) { double u = i/3.0; @@ -405,16 +416,24 @@ void SSurface::MakeTriangulationGridInto(List *l, double vs, double vf, Vector pm1 = PointAtMaybeSwapped(u, vm1, swapped), pm2 = PointAtMaybeSwapped(u, vm2, swapped); + // 0.999 is about 2.5 degrees of twist over the middle 1/3 V-span. + // we don't check at the ends because the derivative may not be valid there. + double twist = 1.0; + if (degm == 1) twist = NormalAtMaybeSwapped(u, vm1, swapped).Dot( + NormalAtMaybeSwapped(u, vm2, swapped) ); + if (twist < worst_twist) worst_twist = twist; + worst = max(worst, pm1.DistanceToLine(ps, pf.Minus(ps))); worst = max(worst, pm2.DistanceToLine(ps, pf.Minus(ps))); } double step = 1.0/SS.GetMaxSegments(); - if((vf - vs) < step || worst < SS.ChordTolMm()) { + if( ((vf - vs) < step || worst < SS.ChordTolMm()) + && ((worst_twist > 0.999) || (depth > 4)) ) { l->Add(&vf); } else { - MakeTriangulationGridInto(l, vs, (vs+vf)/2, swapped); - MakeTriangulationGridInto(l, (vs+vf)/2, vf, swapped); + MakeTriangulationGridInto(l, vs, (vs+vf)/2, swapped, depth+1); + MakeTriangulationGridInto(l, (vs+vf)/2, vf, swapped, depth+1); } } @@ -432,11 +451,22 @@ void SPolygon::UvGridTriangulateInto(SMesh *mesh, SSurface *srf) { List li, lj; li = {}; lj = {}; - double v = 0; - li.Add(&v); - srf->MakeTriangulationGridInto(&li, 0, 1, /*swapped=*/true); - lj.Add(&v); - srf->MakeTriangulationGridInto(&lj, 0, 1, /*swapped=*/false); + double v[5] = {0.0, 0.25, 0.5, 0.75, 1.0}; + li.Add(&v[0]); + srf->MakeTriangulationGridInto(&li, 0, 1, /*swapped=*/true, 0); + lj.Add(&v[0]); + srf->MakeTriangulationGridInto(&lj, 0, 1, /*swapped=*/false, 0); + + // force 2nd order grid to have at least 4 segments in each direction + if ((li.n < 5) && (srf->degm>1)) { // 4 segments minimun + li.Clear(); + li.Add(&v[0]);li.Add(&v[1]);li.Add(&v[2]);li.Add(&v[3]);li.Add(&v[4]); + } + if ((lj.n < 5) && (srf->degn>1)) { // 4 segments minimun + lj.Clear(); + lj.Add(&v[0]);lj.Add(&v[1]);lj.Add(&v[2]);lj.Add(&v[3]);lj.Add(&v[4]); + } + if ((li.n > 3) && (lj.n > 3)) { // Now iterate over each quad in the grid. If it's outside the polygon, // or if it intersects the polygon, then we discard it. Otherwise we