From 7cf3a062745762ecf9ada1823b5e6d1da13e9725 Mon Sep 17 00:00:00 2001 From: Jonathan Westhues Date: Tue, 20 Jan 2009 21:04:38 -0800 Subject: [PATCH] A very rough set of routines to triangulate by ear clipping. This is O(n^2), not perfectly robust, and the bridge-finding code is particularly bad. But it works, triangulates, and shouldn't ever generate zero-area triangles like gl does. [git-p4: depot-paths = "//depot/solvespace/": change = 1900] --- dsc.h | 5 +- glhelper.cpp | 2 +- polygon.cpp | 69 +++++++++++++ polygon.h | 15 +++ srf/ratpoly.cpp | 13 +-- srf/triangulate.cpp | 231 ++++++++++++++++++++++++++++++++++++++++++++ util.cpp | 21 +++- 7 files changed, 347 insertions(+), 9 deletions(-) diff --git a/dsc.h b/dsc.h index f3f9692c..2fb82e8d 100644 --- a/dsc.h +++ b/dsc.h @@ -47,7 +47,8 @@ public: Vector n2, double d2); static Vector AtIntersectionOfLines(Vector a0, Vector a1, Vector b0, Vector b1, - bool *skew); + bool *skew, + double *pa=NULL, double *pb=NULL); double Element(int i); bool Equals(Vector v, double tol=LENGTH_EPS); @@ -73,6 +74,8 @@ public: Vector ProjectVectorInto(hEntity wrkpl); double DivPivoting(Vector delta); Vector ClosestOrtho(void); + void MakeMaxMin(Vector *maxv, Vector *minv); + bool OutsideAndNotOn(Vector maxv, Vector minv); Point2d Project2d(Vector u, Vector v); }; diff --git a/glhelper.cpp b/glhelper.cpp index 675b1145..aa331fbc 100644 --- a/glhelper.cpp +++ b/glhelper.cpp @@ -181,7 +181,7 @@ void glxFillMesh(int specColor, SMesh *m, DWORD h, DWORD s1, DWORD s2) glBegin(GL_TRIANGLES); } - if(tr->an.EqualsExactly(Vector::From(0, 0, 0))) { + if(1 || tr->an.EqualsExactly(Vector::From(0, 0, 0))) { // Compute the normal from the vertices Vector n = tr->Normal(); glNormal3d(n.x, n.y, n.z); diff --git a/polygon.cpp b/polygon.cpp index ab60e6bb..ae67b189 100644 --- a/polygon.cpp +++ b/polygon.cpp @@ -136,6 +136,66 @@ bool SEdgeList::AssemblePolygon(SPolygon *dest, SEdge *errorAt) { } } +//----------------------------------------------------------------------------- +// Test if the specified edge crosses any of the edges in our list. Two edges +// are not considered to cross if they share an endpoint (within LENGTH_EPS), +// but they are considered to cross if they are coincident and overlapping. +//----------------------------------------------------------------------------- +bool SEdgeList::AnyEdgeCrosses(Vector a, Vector b) { + Vector d = b.Minus(a); + double t_eps = LENGTH_EPS/d.Magnitude(); + + SEdge *se; + for(se = l.First(); se; se = l.NextAfter(se)) { + + Vector dse = (se->b).Minus(se->a); + double tse_eps = LENGTH_EPS/dse.Magnitude(); + + if(a.Equals(se->a) && b.Equals(se->b)) return true; + if(b.Equals(se->a) && a.Equals(se->b)) return true; + + double dist_a = (se->a).DistanceToLine(a, d), + dist_b = (se->b).DistanceToLine(a, d); + + if(fabs(dist_a - dist_b) < LENGTH_EPS) { + // The edges are parallel. + if(fabs(dist_a) > LENGTH_EPS) { + // and not coincident, so can't be interesecting + continue; + } + // The edges are coincident. Make sure that neither endpoint lies + // on the other + double t; + t = ((se->a).Minus(a)).DivPivoting(d); + if(t > t_eps && t < (1 - t_eps)) return true; + t = ((se->b).Minus(a)).DivPivoting(d); + if(t > t_eps && t < (1 - t_eps)) return true; + t = a.Minus(se->a).DivPivoting(dse); + if(t > tse_eps && t < (1 - tse_eps)) return true; + t = b.Minus(se->a).DivPivoting(dse); + if(t > tse_eps && t < (1 - tse_eps)) return true; + // So coincident but disjoint, okay. + continue; + } + + // Lines are not parallel, so look for an intersection. + double t, tse; + bool skew; + Vector pi = Vector::AtIntersectionOfLines(a, b, + se->a, se->b, + &skew, + &t, &tse); + if(skew) continue; + + if(t > t_eps && t < (1 - t_eps) && + tse > tse_eps && tse < (1 - tse_eps)) + { + return true; + } + } + return false; +} + void SContour::AddPoint(Vector p) { SPoint sp; sp.tag = 0; @@ -155,6 +215,13 @@ void SContour::MakeEdgesInto(SEdgeList *el) { } } +void SContour::CopyInto(SContour *dest) { + SPoint *sp; + for(sp = l.First(); sp; sp = l.NextAfter(sp)) { + dest->AddPoint(sp->p); + } +} + Vector SContour::ComputeNormal(void) { Vector n = Vector::From(0, 0, 0); @@ -286,12 +353,14 @@ void SPolygon::FixContourDirections(void) { if(sc->l.n < 1) continue; Vector pt = (sc->l.elem[0]).p; + sc->timesEnclosed = 0; bool outer = true; for(j = 0; j < l.n; j++) { if(i == j) continue; SContour *sct = &(l.elem[j]); if(sct->ContainsPointProjdToNormal(normal, pt)) { outer = !outer; + (sc->timesEnclosed)++; } } diff --git a/polygon.h b/polygon.h index 742ab86d..cd7e4974 100644 --- a/polygon.h +++ b/polygon.h @@ -24,17 +24,25 @@ public: bool AssemblePolygon(SPolygon *dest, SEdge *errorAt); bool AssembleContour(Vector first, Vector last, SContour *dest, SEdge *errorAt); + bool AnyEdgeCrosses(Vector a, Vector b); }; class SPoint { public: int tag; + + static const int UNKNOWN = 0; + static const int NOT_EAR = 1; + static const int EAR = 2; + int ear; + Vector p; }; class SContour { public: int tag; + int timesEnclosed; List l; void AddPoint(Vector p); @@ -45,6 +53,12 @@ public: bool ContainsPointProjdToNormal(Vector n, Vector p); bool AllPointsInPlane(Vector n, double d, Vector *notCoplanarAt); void OffsetInto(SContour *dest, double r); + void CopyInto(SContour *dest); + + bool IsEar(int bp); + bool BridgeToContour(SContour *sc, SEdgeList *el, List *vl); + void ClipEarInto(SMesh *m, int bp); + void UvTriangulateInto(SMesh *m); }; typedef struct { @@ -70,6 +84,7 @@ public: bool IsEmpty(void); Vector AnyPoint(void); void OffsetInto(SPolygon *dest, double r); + void UvTriangulateInto(SMesh *m); }; class STriangle { diff --git a/srf/ratpoly.cpp b/srf/ratpoly.cpp index 1ac80b91..1eff6d74 100644 --- a/srf/ratpoly.cpp +++ b/srf/ratpoly.cpp @@ -482,11 +482,12 @@ void SSurface::TriangulateInto(SShell *shell, SMesh *sm) { } int i, start = sm->l.n; + poly.UvTriangulateInto(sm); + STriMeta meta = { 0, 0x888888 }; - poly.normal = Vector::From(0, 0, 1); - poly.TriangulateInto(sm, meta); for(i = start; i < sm->l.n; i++) { STriangle *st = &(sm->l.elem[i]); + st->meta = meta; st->an = NormalAt(st->a.x, st->a.y); st->bn = NormalAt(st->b.x, st->b.y); st->cn = NormalAt(st->c.x, st->c.y); @@ -562,8 +563,8 @@ SShell SShell::FromExtrusionOf(SBezierLoopSet *sbls, Vector t0, Vector t1) { (ret.surface.FindById(hs1))->trim.Add(&stb1); // The translated curves also trim the surface of extrusion. - (ret.surface.FindById(hsext))->trim.Add(&stb0); - (ret.surface.FindById(hsext))->trim.Add(&stb1); +// (ret.surface.FindById(hsext))->trim.Add(&stb0); +// (ret.surface.FindById(hsext))->trim.Add(&stb1); // And form the trim line Vector pt = sb->Finish(); @@ -586,8 +587,8 @@ SShell SShell::FromExtrusionOf(SBezierLoopSet *sbls, Vector t0, Vector t1) { TrimLine *tlp = &(trimLines.elem[WRAP(i-1, trimLines.n)]); - ss->trim.Add(&(tl->trim)); - ss->trim.Add(&(tlp->trim)); +// ss->trim.Add(&(tl->trim)); +// ss->trim.Add(&(tlp->trim)); } trimLines.Clear(); } diff --git a/srf/triangulate.cpp b/srf/triangulate.cpp index e55a3279..446adf29 100644 --- a/srf/triangulate.cpp +++ b/srf/triangulate.cpp @@ -1,3 +1,234 @@ #include "../solvespace.h" +void SPolygon::UvTriangulateInto(SMesh *m) { + if(l.n <= 0) return; + + normal = Vector::From(0, 0, 1); + + while(l.n > 0) { + FixContourDirections(); + l.ClearTags(); + + // Find a top-level contour, and start with that. Then build bridges + // in order to merge all its islands into a single contour. + SContour *top; + for(top = l.First(); top; top = l.NextAfter(top)) { + if(top->timesEnclosed == 0) { + break; + } + } + if(!top) { + dbp("polygon has no top-level contours?"); + return; + } + + // Start with the outer contour + SContour merged; + ZERO(&merged); + top->tag = 1; + top->CopyInto(&merged); + (merged.l.n)--; + + // List all of the edges, for testing whether bridges work. + SEdgeList el; + ZERO(&el); + top->MakeEdgesInto(&el); + List vl; + ZERO(&vl); + + // And now find all of its holes; + SContour *sc; + for(sc = l.First(); sc; sc = l.NextAfter(sc)) { + if(sc->timesEnclosed != 1) continue; + if(sc->l.n < 1) continue; + + if(top->ContainsPointProjdToNormal(normal, sc->l.elem[0].p)) { + sc->tag = 2; + sc->MakeEdgesInto(&el); + } + } + + bool holesExist, mergedHole; + do { + holesExist = false; + mergedHole = false; + + for(sc = l.First(); sc; sc = l.NextAfter(sc)) { + if(sc->tag != 2) continue; + + holesExist = true; + if(merged.BridgeToContour(sc, &el, &vl)) { + // Merged it in, so we're done with it. + sc->tag = 3; + mergedHole = true; + } else { + // Looks like we can't merge it yet, but that can happen; + // we may have to merge the holes in a specific order. + } + } + if(holesExist && !mergedHole) { + dbp("holes exist, yet couldn't merge one"); + return; + } + } while(holesExist); + + merged.UvTriangulateInto(m); + merged.l.Clear(); + el.Clear(); + vl.Clear(); + + l.RemoveTagged(); + } +} + +bool SContour::BridgeToContour(SContour *sc, + SEdgeList *avoidEdges, List *avoidPts) +{ + int thisp, scp; + + Vector a, b, *f; + for(thisp = 0; thisp < l.n; thisp++) { + a = l.elem[thisp].p; + + for(f = avoidPts->First(); f; f = avoidPts->NextAfter(f)) { + if(f->Equals(a)) break; + } + if(f) continue; + + for(scp = 0; scp < (sc->l.n - 1); scp++) { + b = sc->l.elem[scp].p; + + for(f = avoidPts->First(); f; f = avoidPts->NextAfter(f)) { + if(f->Equals(b)) break; + } + if(f) continue; + + if(avoidEdges->AnyEdgeCrosses(a, b)) { + // doesn't work, bridge crosses an existing edge + } else { + goto haveEdge; + } + } + } + // Tried all the possibilities, didn't find an edge + return false; + +haveEdge: + SContour merged; + ZERO(&merged); + int i, j; + for(i = 0; i < l.n; i++) { + merged.AddPoint(l.elem[i].p); + if(i == thisp) { + // less than or equal; need to duplicate the join point + for(j = 0; j <= (sc->l.n - 1); j++) { + int jp = WRAP(j + scp, (sc->l.n - 1)); + merged.AddPoint((sc->l.elem[jp]).p); + } + // and likewise duplicate join point for the outer curve + merged.AddPoint(l.elem[i].p); + } + } + + // and future bridges mustn't cross our bridge, and it's tricky to get + // things right if two bridges come from the same point + avoidEdges->AddEdge(a, b); + avoidPts->Add(&a); + avoidPts->Add(&b); + + l.Clear(); + l = merged.l; + return true; +} + +bool SContour::IsEar(int bp) { + int ap = WRAP(bp-1, l.n), + cp = WRAP(bp+1, l.n); + + STriangle tr; + ZERO(&tr); + tr.a = l.elem[ap].p; + tr.b = l.elem[bp].p; + tr.c = l.elem[cp].p; + + Vector n = Vector::From(0, 0, -1); + if((tr.Normal()).Dot(n) < LENGTH_EPS) { + // This vertex is reflex, or between two collinear edges; either way, + // it's not an ear. + return false; + } + + // Accelerate with an axis-aligned bounding box test + Vector maxv = tr.a, minv = tr.a; + (tr.b).MakeMaxMin(&maxv, &minv); + (tr.c).MakeMaxMin(&maxv, &minv); + + int i; + for(i = 0; i < l.n; i++) { + if(i == ap || i == bp || i == cp) continue; + + Vector p = l.elem[i].p; + if(p.OutsideAndNotOn(maxv, minv)) continue; + + // A point on the edge of the triangle is considered to be inside, + // and therefore makes it a non-ear; but a point on the vertex is + // "outside", since that's necessary to make bridges work. + if(p.EqualsExactly(tr.a)) continue; + if(p.EqualsExactly(tr.b)) continue; + if(p.EqualsExactly(tr.c)) continue; + + if(tr.ContainsPointProjd(n, p)) { + return false; + } + } + return true; +} + +void SContour::ClipEarInto(SMesh *m, int bp) { + int ap = WRAP(bp-1, l.n), + cp = WRAP(bp+1, l.n); + + STriangle tr; + ZERO(&tr); + tr.a = l.elem[ap].p; + tr.b = l.elem[bp].p; + tr.c = l.elem[cp].p; + m->AddTriangle(&tr); + + // By deleting the point at bp, we may change the ear-ness of the points + // on either side. + l.elem[ap].ear = SPoint::UNKNOWN; + l.elem[cp].ear = SPoint::UNKNOWN; + + l.ClearTags(); + l.elem[bp].tag = 1; + l.RemoveTagged(); +} + +void SContour::UvTriangulateInto(SMesh *m) { + int i; + for(i = 0; i < l.n; i++) { + (l.elem[i]).ear = IsEar(i) ? SPoint::EAR : SPoint::NOT_EAR; + } + + while(l.n > 3) { + for(i = 0; i < l.n; i++) { + if(l.elem[i].ear == SPoint::UNKNOWN) { + (l.elem[i]).ear = IsEar(i) ? SPoint::EAR : SPoint::NOT_EAR; + } + } + for(i = 0; i < l.n; i++) { + if(l.elem[i].ear == SPoint::EAR) { + break; + } + } + if(i >= l.n) { + dbp("couldn't find an ear! fail"); + return; + } + ClipEarInto(m, i); + } + + ClipEarInto(m, 0); // add the last triangle +} diff --git a/util.cpp b/util.cpp index e51275e5..c6b5f343 100644 --- a/util.cpp +++ b/util.cpp @@ -539,6 +539,22 @@ Vector Vector::ClosestOrtho(void) { } } +void Vector::MakeMaxMin(Vector *maxv, Vector *minv) { + maxv->x = max(maxv->x, x); + maxv->y = max(maxv->y, y); + maxv->z = max(maxv->z, z); + + minv->x = min(minv->x, x); + minv->y = min(minv->y, y); + minv->z = min(minv->z, z); +} + +bool Vector::OutsideAndNotOn(Vector maxv, Vector minv) { + return (x > maxv.x + LENGTH_EPS) || (x < minv.x - LENGTH_EPS) || + (y > maxv.y + LENGTH_EPS) || (y < minv.y - LENGTH_EPS) || + (z > maxv.z + LENGTH_EPS) || (z < minv.z - LENGTH_EPS); +} + Vector Vector::AtIntersectionOfPlanes(Vector n1, double d1, Vector n2, double d2) { @@ -552,7 +568,8 @@ Vector Vector::AtIntersectionOfPlanes(Vector n1, double d1, Vector Vector::AtIntersectionOfLines(Vector a0, Vector a1, Vector b0, Vector b1, - bool *skew) + bool *skew, + double *parama, double *paramb) { Vector da = a1.Minus(a0), db = b1.Minus(b0); @@ -567,6 +584,8 @@ Vector Vector::AtIntersectionOfLines(Vector a0, Vector a1, // to solve for da and db double pb = ((a0.Minus(b0)).Dot(dna))/(db.Dot(dna)); double pa = -((a0.Minus(b0)).Dot(dnb))/(da.Dot(dnb)); + if(parama) *parama = pa; + if(paramb) *paramb = pb; // And from either of those, we get the intersection point. Vector pi = a0.Plus(da.ScaledBy(pa));