the most stable version up to the moment, max 3 hanging pfaces

Kinetic_shape_reconstruction/examples/Kinetic_shape_reconstruction/CMakeLists.txt

@@ -28,7 +28,7 @@ if(CGAL_FOUND)
       set(targets
       # kinetic_2d_example
       kinetic_precomputed_shapes_example
-      # kinetic_reconstruction_example
+      kinetic_reconstruction_example
       # kinetic_random_shapes_example
       )
 

Kinetic_shape_reconstruction/include/CGAL/KSR_3/Data_structure.h

@@ -401,6 +401,11 @@ public:
         CGAL_assertion(sp_idx_2 != KSR::no_element());
 
         pairs.push_back(std::make_pair(std::make_pair(sp_idx_1, sp_idx_2), false));
+
+        // Makes results much better! ??
+        // Probably because it gives more available intersections between planes
+        // that is the same as increasing k. Is it good? No! Is it correct?
+        pairs.push_back(std::make_pair(std::make_pair(sp_idx_2, sp_idx_1), false));
         ++num_1_intersected;
         if (m_verbose) {
           std::cout << "pair 1: " << std::to_string(sp_idx_1) << "/" << std::to_string(sp_idx_2) << std::endl;
@@ -2581,6 +2586,14 @@ public:
         pv_prev, pv_next, pvertex, ivertex, is_open, reverse,
         iedges, new_crossed, new_pvertices);
     }
+
+    // TODO: IDEA:
+    // What about using potental but go pair by pair and check each potential pface
+    // and add pfaces even if they are between two other blocked pfaces in order to
+    // avoid holes? However, in this case, I need to change implementation of the traverse_iedges()!
+    // It should directly accept potential and I also need to convert crossed to potential
+    // in order to be able to use it as well. In this case, I also do not need to traverse
+    // open case from both sides, one side is enough.
   }
 
   void traverse_iedges(
@@ -2623,8 +2636,8 @@ public:
       const auto& iedge_ip = iedges[ip].first;
 
       bool is_occupied_iedge, is_bbox_reached;
-      // std::tie(is_occupied_iedge, is_bbox_reached) = this->is_occupied(pvertex, ivertex, iedge_i);
-      std::tie(is_occupied_iedge, is_bbox_reached) = this->is_occupied(pvertex, iedge_i);
+      std::tie(is_occupied_iedge, is_bbox_reached) = this->is_occupied(pvertex, ivertex, iedge_i);
+      // std::tie(is_occupied_iedge, is_bbox_reached) = this->is_occupied(pvertex, iedge_i);
 
       const bool is_limit_line = this->is_limit_line(pvertex, iedge_i, is_occupied_iedge);
 
@@ -3941,7 +3954,7 @@ public:
       }
     } while (!quit);
 
-    std::size_t stop_value = 1;
+    std::size_t stop_value = 3;
     if (num_removed_pfaces >= stop_value) {
       std::cout << "* number of removed hanging pfaces: " << num_removed_pfaces << std::endl;
     }

Kinetic_shape_reconstruction/include/CGAL/Kinetic_shape_reconstruction_3.h

@@ -1110,12 +1110,12 @@ private:
         CGAL_assertion(nfaces[0] == pface);
 
         bool is_occupied_iedge_1, is_bbox_reached_1;
-        // std::tie(is_occupied_iedge_1, is_bbox_reached_1) = m_data.collision_occured(pvertex, iedge);
-        std::tie(is_occupied_iedge_1, is_bbox_reached_1) = m_data.is_occupied(pvertex, iedge);
+        std::tie(is_occupied_iedge_1, is_bbox_reached_1) = m_data.collision_occured(pvertex, iedge);
+        // std::tie(is_occupied_iedge_1, is_bbox_reached_1) = m_data.is_occupied(pvertex, iedge);
 
         bool is_occupied_iedge_2, is_bbox_reached_2;
-        // std::tie(is_occupied_iedge_2, is_bbox_reached_2) = m_data.collision_occured(pother, iedge);
-        std::tie(is_occupied_iedge_2, is_bbox_reached_2) = m_data.is_occupied(pother, iedge);
+        std::tie(is_occupied_iedge_2, is_bbox_reached_2) = m_data.collision_occured(pother, iedge);
+        // std::tie(is_occupied_iedge_2, is_bbox_reached_2) = m_data.is_occupied(pother, iedge);
 
         const bool is_limit_line_1 = m_data.is_limit_line(pvertex, iedge, is_occupied_iedge_1);
         const bool is_limit_line_2 = m_data.is_limit_line(pother , iedge, is_occupied_iedge_2);
@@ -1146,6 +1146,13 @@ private:
         } else {
 
           if (m_debug) std::cout << "- free, any k, continue" << std::endl;
+          if (
+            m_data.is_occupied(pvertex, iedge).first ||
+            m_data.is_occupied(pother , iedge).first) {
+
+            CGAL_assertion_msg(false,
+            "ERROR: TWO PVERTICES SNEAK TO THE OTHER SIDE EVEN WHEN WE HAVE A POLYGON!");
+          }
           // CGAL_assertion_msg(false, "TODO: FREE, ANY K, CONTINUE!");
         }
 
@@ -1222,8 +1229,8 @@ private:
     CGAL_assertion(nfaces[0] == pface);
 
     bool is_occupied_iedge, is_bbox_reached;
-    // std::tie(is_occupied_iedge, is_bbox_reached) = m_data.collision_occured(pvertex, iedge);
-    std::tie(is_occupied_iedge, is_bbox_reached) = m_data.is_occupied(pvertex, iedge);
+    std::tie(is_occupied_iedge, is_bbox_reached) = m_data.collision_occured(pvertex, iedge);
+    // std::tie(is_occupied_iedge, is_bbox_reached) = m_data.is_occupied(pvertex, iedge);
 
     const bool is_limit_line = m_data.is_limit_line(pvertex, iedge, is_occupied_iedge);
 

Kinetic_shape_reconstruction/test/Kinetic_shape_reconstruction/kinetic_3d_test_all.cpp

@@ -128,18 +128,6 @@ int main (const int argc, const char** argv) {
   assert(run_test("data/stress-test-4/test-4-rnd-polygons-4-6.off" , ks, num_iters, num_tests));
   assert(run_test("data/stress-test-4/test-5-rnd-polygons-6-4.off" , ks, num_iters, num_tests));
   assert(run_test("data/stress-test-4/test-6-rnd-polygons-5-6.off" , ks, num_iters, num_tests));
-
-  // Edge case tests.
-  assert(run_test("data/edge-case-test/test-flat-bbox-xy.off", ks, num_iters, num_tests)); // flat bbox / 2 coplanar in XY
-  assert(run_test("data/edge-case-test/test-flat-bbox-xz.off", ks, num_iters, num_tests)); // flat bbox / 2 coplanar in XZ
-  assert(run_test("data/edge-case-test/test-flat-bbox-yz.off", ks, num_iters, num_tests)); // flat bbox / 2 coplanar in YZ
-  assert(run_test("data/edge-case-test/test-2-polygons.off"  , ks, num_iters, num_tests)); // edge touch
-  assert(run_test("data/edge-case-test/test-4-polygons.off"  , ks, num_iters, num_tests)); // edge touch / 2 coplanar
-  assert(run_test("data/edge-case-test/test-5-polygons.off"  , ks, num_iters, num_tests)); // edge touch / vertex touch / 2 coplanar
-
-  // All below: unlikely working data sets!
-
-  // Stress tests 4.
   assert(run_test("data/stress-test-4/test-7-rnd-polygons-7-6.off" , ks, num_iters, num_tests));
   assert(run_test("data/stress-test-4/test-8-rnd-polygons-7-8.off" , ks, num_iters, num_tests));
   assert(run_test("data/stress-test-4/test-9-rnd-polygons-12-4.off", ks, num_iters, num_tests));
@@ -149,7 +137,13 @@ int main (const int argc, const char** argv) {
   assert(run_test("data/real-data-test/building-b-15squares-15planes.off" , ks, num_iters, num_tests));
 
   // Edge case tests.
-  assert(run_test("data/edge-case-test/test-20-polygons.off", ks, num_iters, num_tests)); // 2 overlap and coplanar
+  assert(run_test("data/edge-case-test/test-flat-bbox-xy.off", ks, num_iters, num_tests)); // flat bbox / 2 coplanar in XY
+  assert(run_test("data/edge-case-test/test-flat-bbox-xz.off", ks, num_iters, num_tests)); // flat bbox / 2 coplanar in XZ
+  assert(run_test("data/edge-case-test/test-flat-bbox-yz.off", ks, num_iters, num_tests)); // flat bbox / 2 coplanar in YZ
+  assert(run_test("data/edge-case-test/test-2-polygons.off"  , ks, num_iters, num_tests)); // edge touch
+  assert(run_test("data/edge-case-test/test-4-polygons.off"  , ks, num_iters, num_tests)); // edge touch / 2 coplanar
+  assert(run_test("data/edge-case-test/test-5-polygons.off"  , ks, num_iters, num_tests)); // edge touch / vertex touch / 2 coplanar
+  assert(run_test("data/edge-case-test/test-20-polygons.off" , ks, num_iters, num_tests)); // 2 overlap and coplanar
 
   std::cout << std::endl << "--OUTPUT STATS:" << std::endl;
   std::cout << "* number of iterations per test: " << num_iters << std::endl;

Kinetic_shape_reconstruction/todo.md

@@ -1,10 +1,8 @@
 * Use the exact tag in the Delaunay triangulation.
 * After converting exact to inexact, we may have equal points in the bbox faces. We should fix that. Can we?
 * Can we accelerate initializer when using exact kernel?
-* In the merge pvertices event, check for parallel cases not only for the cropping case but for the other cases, too!
-* I should not use pface_of_pvertex() but rather find the correct face that will intersect or not the corresponding iedge and if I should decrement k for it then I do, otherwise I stop/continue. I am pretty sure that pface_of_pvertex() is wrong especially for the open case, back || front subcase.
-
 
+QUESTIONS:
 1. The first step: we compute N full line arrangements. That gives an intersection graph. But to get it, we need to intersect 3D planes. Right?
 2. How do they handle cases with multiple coplanar polygons?
 3. Can we avoid kinetic completely? What if we use labeling instead?
@@ -12,13 +10,8 @@
 5. It looks like for open case, we never insert new polygons. Or better, do we need open case at all? Maybe it should be reassigned into closing case! Because, in the original paper, they do not have such a case.
 6. Graph-cut based surface extraction does not guarantee a 2-manifold. What if we simply stop when we meet a roof polygon or another wall polygon? In this case, we do not guarantee convex polyhedra but who cares?
 7. Do we use the trick from the paper when inserting only the first three events (corresponding to the three closest lines) of an active vertex in the queue instead of all possible events?
-8. What about multiple simultaneous collisions between k primitives at the same time? Do we handle that? Do we need to add a random perturbation as in the original code?
-9. How to do the subdivision and merging efficiently?
-10. How should we choose the pface whose k is updated? Can it differ for different event types?
-11. Do we need to add extra triangles as in my experimental code to handle k?
 
-
-QUESTIONS:
+ANSWERS:
 1. Do we need to intersect all 3D planes in order to find N full line arrangements?
 - Yes.
 
@@ -32,16 +25,4 @@ QUESTIONS:
 - No, it is super difficult.
 
 5. When two polygons intersect at the very beginning, does it count as an intersection? Can we squeeze through the whole or not?
-- Both are ok.
-
-Ideas:
-Do the functions collision_occured/is_occupied work correctly? Not sure.
-Should I use uniform k with respect to the support plane rather than pface? I think yes.
-Try to use for pvertex->iedge events the function version: is_occupied(pvertex, ivertex, iedge).
-Study again the case with 7 and 12 input polygons. They have 1 and 3 hanging pfaces respectively for k = 1. 27 k = 1 for building b and 19 k = 2 for the same building.
-Try to simplify the case as much as possible such that choice of pface does not change anything, we do not have many events until the bad case happens, we are sure that we correctly passed all intermediate events, use only the case k = 1 so that we do not propagate.
-Try to do again what I already tried in the experimental code but for the new refactored case.
-Try to use strict version of the function is_occupied() everywhere.
-What if I extend the initially intersected polygons until they fill the whole cell, then hanging pfaces for the case of the test with 7 and 12 polygons (test-8 and test-9) should not happen!
-Try to see what JP is doing in more detail and compare his steps with ours!
-Where I should use which version of is_occupied()?
+- Both are ok.