More poking at Booleans; generate the unsplit intersection curves

for planes against planes.

[git-p4: depot-paths = "//depot/solvespace/": change = 1907]

Makefile

@@ -42,6 +42,7 @@ SSOBJS   = $(OBJDIR)\solvespace.obj \
 SRFOBJS =  $(OBJDIR)\ratpoly.obj \
            $(OBJDIR)\triangulate.obj \
            $(OBJDIR)\boolean.obj \
+           $(OBJDIR)\surfinter.obj \
 
 
 RES = $(OBJDIR)\resource.res

dsc.h

@@ -75,6 +75,8 @@ public:
     double DivPivoting(Vector delta);
     Vector ClosestOrtho(void);
     void MakeMaxMin(Vector *maxv, Vector *minv);
+    static bool BoundingBoxesDisjoint(Vector amax, Vector amin,
+                                      Vector bmax, Vector bmin);
     bool OutsideAndNotOn(Vector maxv, Vector minv);
     Point2d Project2d(Vector u, Vector v);
 };

srf/boolean.cpp

@@ -10,10 +10,8 @@ void SShell::MakeFromDifferenceOf(SShell *a, SShell *b) {
 
 SCurve SCurve::MakeCopySplitAgainst(SShell *against) {
     SCurve ret;
-    ZERO(&ret);
-
-    ret.isExact = isExact;
-    ret.exact = exact;
+    ret = *this;
+    ZERO(&(ret.pts));
 
     Vector *p;
     for(p = pts.First(); p; p = pts.NextAfter(p)) {
@@ -127,6 +125,29 @@ void SShell::CopySurfacesTrimAgainst(SShell *against, SShell *into,
     }
 }
 
+void SShell::MakeIntersectionCurvesAgainst(SShell *agnst, SShell *into) {
+    SSurface *sa;
+    for(sa = agnst->surface.First(); sa; sa = agnst->surface.NextAfter(sa)) {
+        SSurface *sb;
+        for(sb = surface.First(); sb; sb = surface.NextAfter(sb)) {
+            // Intersect every surface from our shell against every surface
+            // from agnst; this will add zero or more curves to the curve
+            // list for into.
+            sa->IntersectAgainst(sb, into);
+        }
+    }
+}
+
+void SShell::CleanupAfterBoolean(void) {
+    SSurface *ss;
+    for(ss = surface.First(); ss; ss = surface.NextAfter(ss)) {
+        (ss->orig).Clear();
+        (ss->inside).Clear();
+        (ss->onSameNormal).Clear();
+        (ss->onFlipNormal).Clear();
+    }
+}
+
 void SShell::MakeFromBoolean(SShell *a, SShell *b, int type) {
     // Copy over all the original curves, splitting them so that a
     // piecwise linear segment never crosses a surface from the other
@@ -136,9 +157,14 @@ void SShell::MakeFromBoolean(SShell *a, SShell *b, int type) {
 
     // Generate the intersection curves for each surface in A against all
     // the surfaces in B
+    a->MakeIntersectionCurvesAgainst(b, this);
 
     // Then trim and copy the surfaces
     a->CopySurfacesTrimAgainst(b, this, type, true);
     b->CopySurfacesTrimAgainst(a, this, type, false);
+
+    // And clean up the piecewise linear things we made as a calculation aid
+    a->CleanupAfterBoolean();
+    b->CleanupAfterBoolean();
 }
 

srf/ratpoly.cpp

@@ -545,6 +545,18 @@ void SSurface::ClosestPointTo(Vector p, double *u, double *v) {
     }
 }
 
+void SSurface::GetAxisAlignedBounding(Vector *ptMax, Vector *ptMin) {
+    *ptMax = Vector::From(VERY_NEGATIVE, VERY_NEGATIVE, VERY_NEGATIVE);
+    *ptMin = Vector::From(VERY_POSITIVE, VERY_POSITIVE, VERY_POSITIVE);
+
+    int i, j;
+    for(i = 0; i <= degm; i++) {
+        for(j = 0; j <= degn; j++) {
+            (ctrl[i][j]).MakeMaxMin(ptMax, ptMin);
+        }
+    }
+}
+
 void SSurface::MakeEdgesInto(SShell *shell, SEdgeList *sel, bool asUv) {
     STrimBy *stb;
     for(stb = trim.First(); stb; stb = trim.NextAfter(stb)) {
@@ -700,13 +712,19 @@ void SShell::MakeFromExtrusionOf(SBezierLoopSet *sbls, Vector t0, Vector t1,
 
             // The translated curves trim the flat top and bottom surfaces.
             (surface.FindById(hs0))->trim.Add(&stb0);
+            (curve.FindById(hc0))->surfA = hs0;
+
             (surface.FindById(hs1))->trim.Add(&stb1);
+            (curve.FindById(hc1))->surfA = hs1;
 
             // The translated curves also trim the surface of extrusion.
             stb0 = STrimBy::EntireCurve(this, hc0, true);
             (surface.FindById(hsext))->trim.Add(&stb0);
+            (curve.FindById(hc0))->surfB = hsext;
+
             stb1 = STrimBy::EntireCurve(this, hc1, false);
             (surface.FindById(hsext))->trim.Add(&stb1);
+            (curve.FindById(hc1))->surfB = hsext;
 
             // And form the trim line
             Vector pt = sb->Finish();
@@ -734,6 +752,9 @@ void SShell::MakeFromExtrusionOf(SBezierLoopSet *sbls, Vector t0, Vector t1,
             ss->trim.Add(&stb);
             stb = STrimBy::EntireCurve(this, tlp->hc, false);
             ss->trim.Add(&stb);
+
+            (curve.FindById(tl->hc))->surfA = ss->h;
+            (curve.FindById(tlp->hc))->surfB = ss->h;
         }
         trimLines.Clear();
     }

srf/surface.h

@@ -82,13 +82,17 @@ class SCurve {
 public:
     hSCurve         h;
 
-    hSCurve         newH; // when merging with booleans
+    hSCurve         newH;       // when merging with booleans
+    bool            interCurve; // it's a newly-calculated intersection
 
     bool            isExact;
     SBezier         exact;
 
     List<Vector>    pts;
 
+    hSSurface       surfA;
+    hSSurface       surfB;
+
     static SCurve FromTransformationOf(SCurve *a, Vector t, Quaternion q);
     SCurve MakeCopySplitAgainst(SShell *against);
 
@@ -126,6 +130,12 @@ public:
 
     List<STrimBy>   trim;
 
+    // The trims broken down into piecewise linear segments.
+    SEdgeList       orig;
+    SEdgeList       inside;
+    SEdgeList       onSameNormal;
+    SEdgeList       onFlipNormal;
+
     static SSurface FromExtrusionOf(SBezier *spc, Vector t0, Vector t1);
     static SSurface FromPlane(Vector pt, Vector u, Vector v);
     static SSurface FromTransformationOf(SSurface *a, Vector t, Quaternion q, 
@@ -134,11 +144,13 @@ public:
     SSurface MakeCopyTrimAgainst(SShell *against, SShell *shell,
                                  int type, bool opA);
     void TrimFromEdgeList(SEdgeList *el);
+    void IntersectAgainst(SSurface *b, SShell *into);
 
     void ClosestPointTo(Vector p, double *u, double *v);
     Vector PointAt(double u, double v);
     void TangentsAt(double u, double v, Vector *tu, Vector *tv);
     Vector NormalAt(double u, double v);
+    void GetAxisAlignedBounding(Vector *ptMax, Vector *ptMin);
 
     void TriangulateInto(SShell *shell, SMesh *sm);
     void MakeEdgesInto(SShell *shell, SEdgeList *sel, bool asUv);
@@ -162,7 +174,9 @@ public:
     void MakeFromBoolean(SShell *a, SShell *b, int type);
     void CopyCurvesSplitAgainst(SShell *against, SShell *into);
     void CopySurfacesTrimAgainst(SShell *against, SShell *into, int t, bool a);
+    void MakeIntersectionCurvesAgainst(SShell *against, SShell *into);
     void MakeEdgeListUseNewCurveIds(SEdgeList *el);
+    void CleanupAfterBoolean(void);
 
     void MakeFromCopyOf(SShell *a);
     void MakeFromTransformationOf(SShell *a, Vector trans, Quaternion q);

srf/surfinter.cpp

@@ -0,0 +1,56 @@
+#include "solvespace.h"
+
+void SSurface::IntersectAgainst(SSurface *b, SShell *into) {
+    Vector amax, amin, bmax, bmin;
+    GetAxisAlignedBounding(&amax, &amin);
+    b->GetAxisAlignedBounding(&bmax, &bmin);
+
+    if(Vector::BoundingBoxesDisjoint(amax, amin, bmax, bmin)) {
+        // They cannot possibly intersect, no curves to generate
+        return;
+    }
+
+    if(degm == 1 && degn == 1 && b->degm == 1 && b->degn == 1) {
+        // Plane-plane intersection, easy; result is a line
+        Vector pta = ctrl[0][0], ptb = b->ctrl[0][0];
+        Vector na = NormalAt(0, 0), nb = b->NormalAt(0, 0);
+        na = na.WithMagnitude(1);
+        nb = nb.WithMagnitude(1);
+
+        Vector d = (na.Cross(nb));
+
+        if(d.Magnitude() < LENGTH_EPS) {
+            // parallel planes, no intersection
+            return;
+        }
+
+        Vector inter = Vector::AtIntersectionOfPlanes(na, na.Dot(pta),
+                                                      nb, nb.Dot(ptb));
+
+        // The intersection curve can't be longer than the longest curve
+        // that lies in both planes, which is the diagonal of the shorter;
+        // so just pick one, and then give some slop, not critical.
+        double maxl = ((ctrl[0][0]).Minus(ctrl[1][1])).Magnitude();
+
+        Vector v;
+        SCurve sc;
+        ZERO(&sc);
+        sc.surfA = h;
+        sc.surfB = b->h;
+        v = inter.Minus(d.WithMagnitude(2*maxl));
+        sc.pts.Add(&v);
+        v = inter.Plus(d.WithMagnitude(2*maxl));
+        sc.pts.Add(&v);
+        sc.interCurve = true;
+       
+        // Now split the line where it intersects our existing surfaces
+        SCurve split = sc.MakeCopySplitAgainst(into);
+        sc.Clear();
+
+        into->curve.AddAndAssignId(&split);
+    }
+
+    // need to implement general numerical surface intersection for tough
+    // cases, just giving up for now
+}
+

util.cpp

@@ -555,6 +555,17 @@ bool Vector::OutsideAndNotOn(Vector maxv, Vector minv) {
            (z > maxv.z + LENGTH_EPS) || (z < minv.z - LENGTH_EPS);
 }
 
+bool Vector::BoundingBoxesDisjoint(Vector amax, Vector amin,
+                                   Vector bmax, Vector bmin)
+{
+    int i;
+    for(i = 0; i < 3; i++) {
+        if(amax.Element(i) < bmin.Element(i) - LENGTH_EPS) return true;
+        if(amin.Element(i) > bmax.Element(i) + LENGTH_EPS) return true;
+    }
+    return false;
+}
+
 Vector Vector::AtIntersectionOfPlanes(Vector n1, double d1,
                                       Vector n2, double d2)
 {