From 10fc714df918f36dbdeb19179a04be4ca84a399b Mon Sep 17 00:00:00 2001
From: Laurent Rineau <laurent.rineau@cgal.org>
Date: Wed, 13 Dec 2023 14:45:58 +0100
Subject: [PATCH] add options to set the epsilon values (in particular to make
 them 0)

---
 .../CGAL/Conforming_Delaunay_triangulation_3.h   | 10 +++++++---
 .../test/Triangulation_3/cdt_3_from_off.cpp      | 16 +++++++++++++---
 2 files changed, 20 insertions(+), 6 deletions(-)

diff --git a/Triangulation_3/include/CGAL/Conforming_Delaunay_triangulation_3.h b/Triangulation_3/include/CGAL/Conforming_Delaunay_triangulation_3.h
index 69e17508846..122eab37bfd 100644
--- a/Triangulation_3/include/CGAL/Conforming_Delaunay_triangulation_3.h
+++ b/Triangulation_3/include/CGAL/Conforming_Delaunay_triangulation_3.h
@@ -259,6 +259,10 @@ protected:
   }
 
 public:
+  void set_segment_vertex_epsilon(double epsilon) {
+    segment_vertex_epsilon = epsilon;
+  }
+
   Vertex_handle insert(const Point &p, Locate_type lt, Cell_handle c,
                        int li, int lj)
   {
@@ -656,13 +660,12 @@ protected:
         return return_orig_result_point(lambda, orig_pb, orig_pa);
       }
     } else {
-      auto epsilon = 1e-8;
-      if(epsilon > 0) {
+      if(segment_vertex_epsilon > 0) {
         if(!max_bbox_edge_length) {
           update_max_bbox_edge_length();
         }
         auto sq_dist = squared_distance(reference_point, Line{orig_pa, orig_pb});
-        if(sq_dist < CGAL::square(epsilon * *max_bbox_edge_length)) {
+        if(sq_dist < CGAL::square(segment_vertex_epsilon * *max_bbox_edge_length)) {
           std::stringstream ss;
           ss.precision(std::cerr.precision());
           ss << "A constrained segment is too close to a vertex.\n";
@@ -693,6 +696,7 @@ protected:
   Compare_vertex_handle comp = {this};
   Constraint_hierarchy constraint_hierarchy = {comp};
   Bbox_3 bbox{};
+  double segment_vertex_epsilon = 1e-8;
   std::optional<double> max_bbox_edge_length;
   std::map<std::pair<Vertex_handle, Vertex_handle>, Constraint_id> pair_of_vertices_to_cid;
   Insert_in_conflict_visitor insert_in_conflict_visitor = {*this};
diff --git a/Triangulation_3/test/Triangulation_3/cdt_3_from_off.cpp b/Triangulation_3/test/Triangulation_3/cdt_3_from_off.cpp
index 4fddc835d99..2bcdb19c128 100644
--- a/Triangulation_3/test/Triangulation_3/cdt_3_from_off.cpp
+++ b/Triangulation_3/test/Triangulation_3/cdt_3_from_off.cpp
@@ -60,6 +60,8 @@ struct CDT_options
 {
   bool merge_facets = false;
   double ratio = 0.;
+  double vertex_vertex_epsilon = 1e-6;
+  double segment_vertex_epsilon = 1e-8;
   std::string input_filename = CGAL::data_file_path("meshes/mpi.off");
   std::string output_filename{"dump.off"};
   std::string dump_patches_after_merge_filename{};
@@ -81,6 +83,8 @@ Usage: cdt_3_from_off [options] input.off output.off
   --dump-patches-after-merge: dump patches after merging facets
   --dump-patches-borders-prefix: dump patches borders
   --dump-after-conforming: dump mesh after conforming
+  --vertex-vertex-epsilon: epsilon for vertex-vertex min distance (default: 1e-6)
+  --segment-vertex-epsilon: epsilon for segment-vertex min distance (default: 0)
 )";
 }
 
@@ -108,6 +112,12 @@ int main(int argc, char* argv[])
     } else if(arg == "--dump-after-conforming") {
       assert(i + 1 < argc);
       options.dump_after_conforming_filename = argv[++i];
+    } else if(arg == "--vertex-vertex-epsilon") {
+      assert(i + 1 < argc);
+      options.vertex_vertex_epsilon = std::stod(argv[++i]);
+    } else if(arg == "--segment-vertex-epsilon") {
+      assert(i + 1 < argc);
+      options.segment_vertex_epsilon = std::stod(argv[++i]);
     } else if(arg == "--help") {
       help(std::cout);
       return 0;
@@ -253,6 +263,7 @@ auto segment_soup_to_polylines(Range_of_segments&& segment_soup) {
 
 int go(Mesh mesh, CDT_options options) {
   CDT cdt;
+  cdt.set_segment_vertex_epsilon(options.segment_vertex_epsilon);
   auto pmap = get(CGAL::vertex_point, mesh);
 
   auto [patch_id_map, ok] = mesh.add_property_map<face_descriptor, int>("f:patch_id", -1);
@@ -390,11 +401,10 @@ int go(Mesh mesh, CDT_options options) {
     cdt.insert(Point(bbox.xmax() + d_x, bbox.ymax() + d_y, bbox.zmax() + d_z));
   }
   {
-    double espilon = 1e-6;
+    double epsilon = options.vertex_vertex_epsilon;
     auto min_distance = CGAL::approximate_sqrt(std::ranges::min(
         cdt.finite_edges() | std::views::transform([&](auto edge) { return cdt.segment(edge).squared_length(); })));
-    std::cout << "Min distance between vertices: " << min_distance << '\n';
-    if(min_distance < espilon * max_d) {
+    if(min_distance < epsilon * max_d) {
       std::cerr << "ERROR: min distance between vertices is too small\n";
       exit_code = EXIT_FAILURE;
       return exit_code;