From 6647c9e0b37a72db757b08428329f5fa8cb56cc5 Mon Sep 17 00:00:00 2001
From: RonaldoCMP <rcmpersiano@gmail.com>
Date: Fri, 4 Jun 2021 00:57:50 +0100
Subject: [PATCH 1/3] Speeding up of hull3d

---
 hull.scad            | 49 ++++++++++++++++++++++----------------------
 tests/test_hull.scad |  8 ++++----
 2 files changed, 29 insertions(+), 28 deletions(-)

diff --git a/hull.scad b/hull.scad
index 1bffeb7..4c32f26 100644
--- a/hull.scad
+++ b/hull.scad
@@ -81,8 +81,8 @@ function _backtracking(i,points,h,t,m,all) =
 
 // clockwise check (2d)
 function _is_cw(a,b,c,all) = 
-		all ? cross(a-c,b-c)<=EPSILON*norm(a-c)*norm(b-c) :
-		cross(a-c,b-c)<-EPSILON*norm(a-c)*norm(b-c);
+    all ? cross(a-c,b-c)<=EPSILON*norm(a-c)*norm(b-c) :
+    cross(a-c,b-c)<-EPSILON*norm(a-c)*norm(b-c);
 
 
 // Function: hull2d_path()
@@ -111,7 +111,7 @@ function hull2d_path(points, all=false) =
          ip = sortidx(points) )
     // lower hull points
     let( lh = 
-            [ for( 	i = 2,
+            [ for(   i = 2,
                     k = 2, 
                     h = [ip[0],ip[1]]; // current list of hull point indices 
                   i <= n;
@@ -120,7 +120,7 @@ function hull2d_path(points, all=false) =
                     i = i+1
                  ) if( i==n ) h ][0] )
     // concat lower hull points with upper hull ones
-    [ for( 	i = n-2,
+    [ for(   i = n-2,
             k = len(lh), 
             t = k+1,
             h = lh; // current list of hull point indices 
@@ -129,7 +129,7 @@ function hull2d_path(points, all=false) =
             h = [for(j=[0:1:k-2]) h[j], if(i>0) ip[i]],
             i = i-1
          ) if( i==-1 ) h ][0] ;
-			 
+       
 
 function _hull_collinear(points) =
     let(
@@ -207,25 +207,33 @@ function _hull3d_iterative(points, triangles, planes, remaining, _i=0) =
     let (
         // pick a point
         i = remaining[_i],
+        // evaluate the triangle plane equations at point i
+        planeq_val = planes*[each points[i], -1],
         // find the triangles that are in conflict with the point (point not inside)
-        conflicts = _find_conflicts(points[i], planes),
-        // for all triangles that are in conflict, collect their halfedges
+        conflicts = [ for (i = [0:1:len(planes)-1]) if (planeq_val[i]>0) i ],
+        // collect the halfedges of all triangles that are in conflict 
         halfedges = [ 
-            for(c = conflicts, i = [0:2]) let(
-                j = (i+1)%3
-            ) [triangles[c][i], triangles[c][j]]
+            for(c = conflicts, i = [0:2])
+                [triangles[c][i], triangles[c][(i+1)%3]]
         ],
         // find the outer perimeter of the set of conflicting triangles
         horizon = _remove_internal_edges(halfedges),
-        // generate a new triangle for each horizon halfedge together with the picked point i
-        new_triangles = [ for (h = horizon) concat(h,i) ],
-        // calculate the corresponding plane equations
-        new_planes = [ for (t = new_triangles) plane3pt_indexed(points, t[0], t[1], t[2]) ]
+        // generate new triangles connecting point i to each horizon halfedge vertices
+        tri2add = [ for (h = horizon) concat(h,i) ],
+        // add tria2add and remove conflict triangles
+        new_triangles = 
+            concat( tri2add,
+                    [ for (i = [0:1:len(planes)-1]) if (planeq_val[i]<=0) triangles[i] ] 
+                  ),
+        // add the plane equations of new added triangles and remove the plane equations of the conflict ones
+        new_planes = 
+            concat( [ for (t = tri2add) plane3pt_indexed(points, t[0], t[1], t[2]) ],
+                    [ for (i = [0:1:len(planes)-1]) if (planeq_val[i]<=0) planes[i] ] 
+                  )
     ) _hull3d_iterative(
         points,
-        //  remove the conflicting triangles and add the new ones
-        concat(list_remove(triangles, conflicts), new_triangles),
-        concat(list_remove(planes, conflicts), new_planes),
+        new_triangles,
+        new_planes,
         remaining,
         _i+1
     );
@@ -238,13 +246,6 @@ function _remove_internal_edges(halfedges) = [
 ];
 
 
-function _find_conflicts(point, planes) = [
-    for (i = [0:1:len(planes)-1])
-        if (in_front_of_plane(planes[i], point))
-            i
-];
-
-
 function _find_first_noncoplanar(plane, points, i=0) = 
     (i >= len(points) || !points_on_plane([points[i]],plane))? i :
     _find_first_noncoplanar(plane, points, i+1);
diff --git a/tests/test_hull.scad b/tests/test_hull.scad
index b294d5a..396a2b9 100644
--- a/tests/test_hull.scad
+++ b/tests/test_hull.scad
@@ -87,12 +87,12 @@ module test_hull() {
                  [-4.81079, -0.967785, -10.4726], [-0.949023, 23.1441, -2.08208], [16.1256, -8.2295, -24.0113], [6.45274,
                  -7.21416, 23.1409], [22.8274, 1.07038, 19.1756], [-10.6256, -10.0112, -6.12274], [6.29254, -7.81875,
                  -24.4037], [22.8538, 8.78163, -6.82567], [-1.96142, 19.1728, -1.726]];
-    assert_equal(hull(rand25_3d),[[21, 29, 11], [29, 21, 20], [21, 14, 20], [20, 14, 26], [15, 0, 28], [13, 29, 31], [0, 15,
+    assert_equal(sort(hull(rand25_3d)),sort([[21, 29, 11], [29, 21, 20], [21, 14, 20], [20, 14, 26], [15, 0, 28], [13, 29, 31], [0, 15,
                                  31], [15, 9, 31], [9, 24, 31], [24, 13, 31], [28, 0, 31], [11, 29, 35], [29, 13, 35], [15,
                                  21, 36], [9, 15, 36], [24, 9, 36], [13, 24, 36], [15, 28, 37], [28, 26, 37], [28, 31, 39],
                                  [31, 29, 39], [14, 21, 42], [21, 15, 42], [26, 14, 42], [15, 37, 42], [37, 26, 42], [29, 20,
                                  43], [39, 29, 43], [20, 26, 43], [26, 28, 43], [21, 13, 44], [13, 36, 44], [36, 21, 44],
-                                 [21, 11, 45], [11, 35, 45], [13, 21, 45], [35, 13, 45], [28, 39, 47], [39, 43, 47], [43, 28, 47]]);
+                                 [21, 11, 45], [11, 35, 45], [13, 21, 45], [35, 13, 45], [28, 39, 47], [39, 43, 47], [43, 28, 47]]));
 
     /*  // Inconsistently treats coplanar faces: sometimes face center vertex is included in output, sometimes not
     test_cube_3d = [for(x=[1:3], y=[1:3], z=[1:3]) [x,y,z]];
@@ -178,12 +178,12 @@ module test_hull3d_faces() {
                  [-4.81079, -0.967785, -10.4726], [-0.949023, 23.1441, -2.08208], [16.1256, -8.2295, -24.0113], [6.45274,
                  -7.21416, 23.1409], [22.8274, 1.07038, 19.1756], [-10.6256, -10.0112, -6.12274], [6.29254, -7.81875,
                  -24.4037], [22.8538, 8.78163, -6.82567], [-1.96142, 19.1728, -1.726]];
-    assert_equal(hull(rand25_3d),[[21, 29, 11], [29, 21, 20], [21, 14, 20], [20, 14, 26], [15, 0, 28], [13, 29, 31], [0, 15,
+    assert_equal(sort(hull(rand25_3d)), sort([[21, 29, 11], [29, 21, 20], [21, 14, 20], [20, 14, 26], [15, 0, 28], [13, 29, 31], [0, 15,
                                  31], [15, 9, 31], [9, 24, 31], [24, 13, 31], [28, 0, 31], [11, 29, 35], [29, 13, 35], [15,
                                  21, 36], [9, 15, 36], [24, 9, 36], [13, 24, 36], [15, 28, 37], [28, 26, 37], [28, 31, 39],
                                  [31, 29, 39], [14, 21, 42], [21, 15, 42], [26, 14, 42], [15, 37, 42], [37, 26, 42], [29, 20,
                                  43], [39, 29, 43], [20, 26, 43], [26, 28, 43], [21, 13, 44], [13, 36, 44], [36, 21, 44],
-                                 [21, 11, 45], [11, 35, 45], [13, 21, 45], [35, 13, 45], [28, 39, 47], [39, 43, 47], [43, 28, 47]]);
+                                 [21, 11, 45], [11, 35, 45], [13, 21, 45], [35, 13, 45], [28, 39, 47], [39, 43, 47], [43, 28, 47]]));
 }
 test_hull3d_faces();
 

From 6793563aec8bb74b316cf1d32ab341ecfda1ded3 Mon Sep 17 00:00:00 2001
From: RonaldoCMP <rcmpersiano@gmail.com>
Date: Fri, 4 Jun 2021 13:07:01 +0100
Subject: [PATCH 2/3] Allows a tolerance in the computation of plane equations

---
 hull.scad | 10 +++++++---
 1 file changed, 7 insertions(+), 3 deletions(-)

diff --git a/hull.scad b/hull.scad
index 4c32f26..814ce6c 100644
--- a/hull.scad
+++ b/hull.scad
@@ -208,9 +208,11 @@ function _hull3d_iterative(points, triangles, planes, remaining, _i=0) =
         // pick a point
         i = remaining[_i],
         // evaluate the triangle plane equations at point i
+//			  xx=[for(i=[0:len(planes)-1],p=[planes[i]]) if(len(p)!=4) echo(i=i,len(p))0 ],//echo([each points[i], -1]),
+//        planeq_val = [for(p=planes) p*[each points[i], -1]],
         planeq_val = planes*[each points[i], -1],
         // find the triangles that are in conflict with the point (point not inside)
-        conflicts = [ for (i = [0:1:len(planes)-1]) if (planeq_val[i]>0) i ],
+        conflicts = [ for (i = [0:1:len(planes)-1]) if (planeq_val[i]>EPSILON) i ],
         // collect the halfedges of all triangles that are in conflict 
         halfedges = [ 
             for(c = conflicts, i = [0:2])
@@ -220,15 +222,17 @@ function _hull3d_iterative(points, triangles, planes, remaining, _i=0) =
         horizon = _remove_internal_edges(halfedges),
         // generate new triangles connecting point i to each horizon halfedge vertices
         tri2add = [ for (h = horizon) concat(h,i) ],
+//				w=[for(t=tri2add) if(collinear(points[t[0]],points[t[1]],points[t[2]])) echo(t)],
         // add tria2add and remove conflict triangles
         new_triangles = 
             concat( tri2add,
-                    [ for (i = [0:1:len(planes)-1]) if (planeq_val[i]<=0) triangles[i] ] 
+                    [ for (i = [0:1:len(planes)-1]) if (planeq_val[i]<=EPSILON) triangles[i] ] 
                   ),
+//				y=[for(t=tri2add) if([]==plane3pt_indexed(points,t[0],t[1],t[2])) echo(tri2add=t,pts=[for(ti=t) points[ti]])],
         // add the plane equations of new added triangles and remove the plane equations of the conflict ones
         new_planes = 
             concat( [ for (t = tri2add) plane3pt_indexed(points, t[0], t[1], t[2]) ],
-                    [ for (i = [0:1:len(planes)-1]) if (planeq_val[i]<=0) planes[i] ] 
+                    [ for (i = [0:1:len(planes)-1]) if (planeq_val[i]<=EPSILON) planes[i] ] 
                   )
     ) _hull3d_iterative(
         points,

From 60ad2839d0e7ddea969b598a7f12ff28f23c0725 Mon Sep 17 00:00:00 2001
From: RonaldoCMP <rcmpersiano@gmail.com>
Date: Fri, 4 Jun 2021 13:09:47 +0100
Subject: [PATCH 3/3] Cosmetic

---
 hull.scad | 6 +++---
 1 file changed, 3 insertions(+), 3 deletions(-)

diff --git a/hull.scad b/hull.scad
index 814ce6c..1dc147c 100644
--- a/hull.scad
+++ b/hull.scad
@@ -208,7 +208,7 @@ function _hull3d_iterative(points, triangles, planes, remaining, _i=0) =
         // pick a point
         i = remaining[_i],
         // evaluate the triangle plane equations at point i
-//			  xx=[for(i=[0:len(planes)-1],p=[planes[i]]) if(len(p)!=4) echo(i=i,len(p))0 ],//echo([each points[i], -1]),
+//        xx=[for(i=[0:len(planes)-1],p=[planes[i]]) if(len(p)!=4) echo(i=i,len(p))0 ],//echo([each points[i], -1]),
 //        planeq_val = [for(p=planes) p*[each points[i], -1]],
         planeq_val = planes*[each points[i], -1],
         // find the triangles that are in conflict with the point (point not inside)
@@ -222,13 +222,13 @@ function _hull3d_iterative(points, triangles, planes, remaining, _i=0) =
         horizon = _remove_internal_edges(halfedges),
         // generate new triangles connecting point i to each horizon halfedge vertices
         tri2add = [ for (h = horizon) concat(h,i) ],
-//				w=[for(t=tri2add) if(collinear(points[t[0]],points[t[1]],points[t[2]])) echo(t)],
+//        w=[for(t=tri2add) if(collinear(points[t[0]],points[t[1]],points[t[2]])) echo(t)],
         // add tria2add and remove conflict triangles
         new_triangles = 
             concat( tri2add,
                     [ for (i = [0:1:len(planes)-1]) if (planeq_val[i]<=EPSILON) triangles[i] ] 
                   ),
-//				y=[for(t=tri2add) if([]==plane3pt_indexed(points,t[0],t[1],t[2])) echo(tri2add=t,pts=[for(ti=t) points[ti]])],
+//        y=[for(t=tri2add) if([]==plane3pt_indexed(points,t[0],t[1],t[2])) echo(tri2add=t,pts=[for(ti=t) points[ti]])],
         // add the plane equations of new added triangles and remove the plane equations of the conflict ones
         new_planes = 
             concat( [ for (t = tri2add) plane3pt_indexed(points, t[0], t[1], t[2]) ],