Made initial_point parameter work with any Range

With doc and example (
 example "mesh_3D_image_with_initial_points.cpp" has been renamed to  "mesh_3D_image_with_image_initialization.cpp")

Mesh_3/doc/Mesh_3/CGAL/Mesh_3/parameters.h

@@ -486,9 +486,10 @@ unspecified_type initial_points_generator(const InitialPointsGenerator& generato
  * \ingroup PkgMesh3Parameters
  *
  * The function `parameters::initial_points()` enables the user to
- * specify a `std::vector` of initial points
+ * specify a container model of `Range` of initial points
  * to the mesh generation function `make_mesh_3()`.
- * The initial points vector is of type `std::vector<std::tuple<Weighted_point_3, int, Index>>` where
+ * The initial points `Range` has elements of type
+ * `std::tuple<Weighted_point_3, int, Index>` where
  * `Weighted_point_3` contains the point's position and weight,
  * `int` is the dimension of the minimal dimension subcomplex on which the point lies, and
  * `Index` is the underlying subcomplex index.
@@ -496,7 +497,7 @@ unspecified_type initial_points_generator(const InitialPointsGenerator& generato
  * \tparam MeshDomain model of `MeshDomain_3`
  * \tparam C3t3 model of `MeshComplex_3InTriangulation_3`
  *
- *  @param initial_points a vector containing points of type
+ *  @param initial_points a `Range` containing points of type
  * `std::tuple<C3t3::Triangulation::Geom_traits::Weighted_point_3, int, MeshDomain::Index>`
  *
  * \cgalHeading{Example}
@@ -504,7 +505,7 @@ unspecified_type initial_points_generator(const InitialPointsGenerator& generato
  * \code{.cpp}
  * // Creation of the initial_points vector
  * std::vector<std::tuple<K::Weighted_point_3, int, Mesh_domain::Index>> initial_points;
- * // Mesh generation from labeled image with connexity checks.
+ * // Mesh generation from labeled image with initial points.
  * C3t3 c3t3 = make_mesh_3<c3t3>(domain,
  *                               criteria,
  *                               parameters::initial_points(std::cref(initial_points));//use std::cref to avoid a copy

Mesh_3/doc/Mesh_3/Mesh_3.txt

@@ -737,14 +737,14 @@ without the weights (left, 25563 vertices) and with the weights (right, 19936 ve
 
 \subsubsection Mesh_3DomainsFrom3DImagesWithCustomInitialization Domains From 3D Images, with a Custom Initialization
 
-The example \ref Mesh_3/mesh_3D_image_with_initial_points.cpp is a modification
+The example \ref Mesh_3/mesh_3D_image_with_image_initialization.cpp is a modification
 of \ref Mesh_3/mesh_3D_image.cpp. The goal of that example is to show how
 the default initialization of the triangulation, using random rays, can be
 replaced by a new implementation. In this case, the initialization detects
 all connected components in the 3D segmented image, and inserts points in
 the triangulation for each connected component.
 
-\snippet Mesh_3/mesh_3D_image_with_initial_points.cpp Meshing
+\snippet Mesh_3/mesh_3D_image_with_image_initialization.cpp Meshing
 
 The parameter `CGAL::parameters::initial_points_generator` is used.
 It expects a functor that returns a set of points for the mesh
@@ -785,10 +785,10 @@ Right: the mesh generated after the initialization of all connected components
 \cgalFigureCaptionEnd
 
 Note that the example \ref
-Mesh_3/mesh_3D_image_with_initial_points.cpp also shows how to
+Mesh_3/mesh_3D_image_with_image_initialization.cpp also shows how to
 create a 3D image using the undocumented API of CGAL_ImageIO.
 
-\snippet Mesh_3/mesh_3D_image_with_initial_points.cpp Create the image
+\snippet Mesh_3/mesh_3D_image_with_image_initialization.cpp Create the image
 
 The code of the function `%random_labeled_image()` is in the header file \ref
 Mesh_3/random_labeled_image.h.
@@ -800,6 +800,10 @@ a custom functor that initializes the triangulation.
 A struct `Custom_Initial_points_generator` that places 1D-feature points
 alongside a line is created.
 
+Finally, the exemple \ref Mesh_3/mesh_3D_image_with_initial_points.cpp features
+a point container that initializes the triangulation using the meshing parameter
+`parameters::initial_points()`.
+
 \subsection Mesh_3UsingVariableSizingField Using Variable Sizing Field
 
 \subsubsection Mesh_3SizingFieldasanAnalyticalFunction Sizing Field as an Analytical Function

Mesh_3/doc/Mesh_3/examples.txt

@@ -4,6 +4,7 @@
 \example Mesh_3/mesh_3D_image_with_features.cpp
 \example Mesh_3/mesh_3D_image_with_custom_initialization.cpp
 \example Mesh_3/mesh_3D_image_with_initial_points.cpp
+\example Mesh_3/mesh_3D_image_with_image_initialization.cpp
 \example Mesh_3/mesh_3D_image_with_detection_of_features.cpp
 \example Mesh_3/mesh_3D_image_with_input_features.cpp
 \example Mesh_3/mesh_3D_weighted_image.cpp

Mesh_3/examples/Mesh_3/CMakeLists.txt

@@ -162,6 +162,11 @@ if(TARGET CGAL::CGAL_ImageIO)
     target_link_libraries(mesh_3D_image_with_custom_initialization
                           PUBLIC CGAL::Eigen3_support)
 
+    create_single_source_cgal_program(
+      "mesh_3D_image_with_image_initialization.cpp")
+    target_link_libraries(mesh_3D_image_with_image_initialization
+                          PUBLIC CGAL::Eigen3_support)
+
     create_single_source_cgal_program(
       "mesh_3D_image_with_initial_points.cpp")
     target_link_libraries(mesh_3D_image_with_initial_points
@@ -201,6 +206,7 @@ if(CGAL_ACTIVATE_CONCURRENT_MESH_3 AND TARGET CGAL::TBB_support)
     mesh_3D_image_variable_size
     mesh_3D_image_with_custom_initialization
     mesh_3D_image_with_initial_points
+    mesh_3D_image_with_image_initialization
     mesh_3D_image_with_features
     mesh_3D_image_with_detection_of_features
     mesh_3D_image_with_input_features

Mesh_3/examples/Mesh_3/mesh_3D_image_with_custom_initialization.cpp

@@ -35,6 +35,11 @@ typedef CGAL::Mesh_criteria_3<Tr> Mesh_criteria;
 
 namespace params = CGAL::parameters;
 
+// Custom_Initial_points_generator will put points on the mesh for initialisation.
+// Those points are objects of type std::tuple<Weighted_point_3, int, Index>.
+//   Weighted_point_3 is the point's position and weight,
+//   int is the dimension of the minimal dimension subcomplex on which the point lies,
+//   Index is the underlying subcomplex index.
 struct Custom_Initial_points_generator
 {
   CGAL::Image_3& image_;
@@ -43,20 +48,26 @@ struct Custom_Initial_points_generator
   template <typename OutputIterator, typename MeshDomain, typename C3t3>
   OutputIterator operator()(OutputIterator pts, const MeshDomain& domain, const C3t3& c3t3, int n = 1) const
   {
-    typedef typename C3t3::Triangulation::Geom_traits::Point_3 Point_3;
+    typedef typename C3t3::Triangulation::Geom_traits::Point_3   Point_3;
+    typedef typename C3t3::Triangulation::Geom_traits::Vector_3  Vector_3;
+    typedef typename C3t3::Triangulation::Geom_traits::Segment_3 Segment_3;
 
     typename C3t3::Triangulation::Geom_traits::Construct_weighted_point_3 cwp =
         c3t3.triangulation().geom_traits().construct_weighted_point_3_object();
 
-    // Add points along the segment from
+    // Add points along the segment
-    // (  0.0 50.0 66.66) to
+    Segment_3 segment(Point_3(  0.0, 50.0, 66.66),
-    // (100.0 50.0 66.66)
+                      Point_3(100.0, 50.0, 66.66));
+
+
+    Point_3 source = segment.source();
+    Vector_3 vector = segment.to_vector();
     double edge_size = 5;
-    std::size_t nb = static_cast<int>(100.0 / edge_size);
+    std::size_t nb = static_cast<int>(CGAL::sqrt(segment.squared_length()) / edge_size);
     for (std::size_t i = 1; i < nb; i++)
     {
       *pts++ = std::make_tuple(
-         cwp(Point_3(i*edge_size, 50.0, 66.66), edge_size*edge_size), 1, 0);
+         cwp(source + (i/(double)nb)*vector, edge_size*edge_size), 1, 0);
     }
     return pts;
   }

Mesh_3/examples/Mesh_3/mesh_3D_image_with_image_initialization.cpp

@@ -0,0 +1,60 @@
+#include "random_labeled_image.h"
+#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
+
+#include <CGAL/Mesh_triangulation_3.h>
+#include <CGAL/Mesh_complex_3_in_triangulation_3.h>
+#include <CGAL/Mesh_criteria_3.h>
+
+#include <CGAL/Mesh_3/Construct_initial_points_labeled_image.h>
+
+#include <CGAL/Labeled_mesh_domain_3.h>
+#include <CGAL/make_mesh_3.h>
+#include <CGAL/Image_3.h>
+
+#include <CGAL/SMDS_3/Dump_c3t3.h>
+
+// Domain
+typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
+typedef CGAL::Labeled_mesh_domain_3<K> Mesh_domain;
+
+#ifdef CGAL_CONCURRENT_MESH_3
+typedef CGAL::Parallel_tag Concurrency_tag;
+#else
+typedef CGAL::Sequential_tag Concurrency_tag;
+#endif
+
+// Triangulation
+typedef CGAL::Mesh_triangulation_3<Mesh_domain,CGAL::Default,Concurrency_tag>::type Tr;
+
+typedef CGAL::Mesh_complex_3_in_triangulation_3<Tr> C3t3;
+
+// Criteria
+typedef CGAL::Mesh_criteria_3<Tr> Mesh_criteria;
+
+namespace params = CGAL::parameters;
+
+int main()
+{
+  /// [Create the image]
+  CGAL::Image_3 image = random_labeled_image();
+  /// [Create the image]
+
+  // Domain
+  Mesh_domain domain = Mesh_domain::create_labeled_image_mesh_domain(image);
+
+  // Mesh criteria
+  Mesh_criteria criteria(params::facet_angle(30).facet_size(3).facet_distance(1)
+                         .cell_radius_edge_ratio(3).cell_size(3)
+  );
+
+  /// [Meshing]
+  C3t3 c3t3 = CGAL::make_mesh_3<C3t3>(domain, criteria
+    , params::initial_points_generator(CGAL::Construct_initial_points_labeled_image(image))
+  );
+  /// [Meshing]
+
+  // Output
+  CGAL::dump_c3t3(c3t3, "out");
+
+  return 0;
+}

Mesh_3/examples/Mesh_3/mesh_3D_image_with_initial_points.cpp

@@ -5,17 +5,19 @@
 #include <CGAL/Mesh_complex_3_in_triangulation_3.h>
 #include <CGAL/Mesh_criteria_3.h>
 
-#include <CGAL/Mesh_3/Construct_initial_points_labeled_image.h>
-
 #include <CGAL/Labeled_mesh_domain_3.h>
 #include <CGAL/make_mesh_3.h>
 #include <CGAL/Image_3.h>
 
 #include <CGAL/SMDS_3/Dump_c3t3.h>
 
+
+#include <CGAL/Mesh_3/Detect_features_in_image.h>
+
 // Domain
 typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
-typedef CGAL::Labeled_mesh_domain_3<K> Mesh_domain;
+typedef CGAL::Labeled_mesh_domain_3<K> Image_domain;
+typedef CGAL::Mesh_domain_with_polyline_features_3<Image_domain> Mesh_domain;
 
 #ifdef CGAL_CONCURRENT_MESH_3
 typedef CGAL::Parallel_tag Concurrency_tag;
@@ -35,21 +37,39 @@ namespace params = CGAL::parameters;
 
 int main()
 {
-  /// [Create the image]
+  const std::string fname = CGAL::data_file_path("images/420.inr");
-  CGAL::Image_3 image = random_labeled_image();
+  // Loads image
-  /// [Create the image]
+  CGAL::Image_3 image;
+  if(!image.read(fname)){
+    std::cerr << "Error: Cannot read file " <<  fname << std::endl;
+    return EXIT_FAILURE;
+  }
 
   // Domain
-  Mesh_domain domain = Mesh_domain::create_labeled_image_mesh_domain(image);
+  Mesh_domain domain = Mesh_domain::create_labeled_image_mesh_domain(image
+      , params::features_detector(CGAL::Mesh_3::Detect_features_in_image())
+  );
 
   // Mesh criteria
-  Mesh_criteria criteria(params::facet_angle(30).facet_size(3).facet_distance(1)
+  Mesh_criteria criteria(params::facet_angle(30).facet_size(3).facet_distance(1).edge_size(3)
                          .cell_radius_edge_ratio(3).cell_size(3)
   );
 
+  using Point_3 = K::Point_3;
+  using Weighted_point_3 = K::Weighted_point_3;
+  using Index = Mesh_domain::Index;
+  using Initial_point_t = std::tuple<Weighted_point_3, int, Index>;
+
+  // Creation of the initial_points vector
+  std::vector<Initial_point_t> initial_points = {
+    std::make_tuple(Weighted_point_3(Point_3(30.0, 50.0, 83.33), 30.0), 1, 0),
+    std::make_tuple(Weighted_point_3(Point_3(70.0, 50.0, 83.33), 50.0), 1, 0)
+  };
+
   /// [Meshing]
+  // Mesh generation from labeled image with initial points.
   C3t3 c3t3 = CGAL::make_mesh_3<C3t3>(domain, criteria
-    , params::initial_points_generator(CGAL::Construct_initial_points_labeled_image(image))
+    , params::initial_points(std::cref(initial_points))
   );
   /// [Meshing]
 

Mesh_3/include/CGAL/make_mesh_3.h

@@ -40,9 +40,9 @@ namespace internal {
 
 template < typename C3T3, typename MeshDomain, typename MeshCriteria >
 void
-init_c3t3(C3T3& c3t3, const MeshDomain& domain, const MeshCriteria&,
+add_points_from_generator(C3T3& c3t3, const MeshDomain& domain,
-          const int nb_initial_points,
+                          const int nb_initial_points,
-          const parameters::internal::Initial_points_generator_options<MeshDomain, C3T3>& generator = parameters::internal::Initial_points_generator_generator<MeshDomain, C3T3>()())
+                          const parameters::internal::Initial_points_generator_options<MeshDomain, C3T3>& generator)
 {
   typedef typename C3T3::Triangulation::Geom_traits::Weighted_point_3 Weighted_point_3;
   typedef typename MeshDomain::Index Index;
@@ -70,6 +70,39 @@ init_c3t3(C3T3& c3t3, const MeshDomain& domain, const MeshCriteria&,
   }
 }
 
+template < typename C3T3, typename MeshDomain, typename MeshCriteria >
+void
+init_c3t3(C3T3& c3t3, const MeshDomain& domain, const MeshCriteria&,
+          const int nb_initial_points,
+          const parameters::internal::Initial_points_generator_options<MeshDomain, C3T3>& generator = parameters::internal::Initial_points_generator_generator<MeshDomain, C3T3>()())
+{
+  add_points_from_generator<C3T3, MeshDomain, MeshCriteria>(c3t3, domain, nb_initial_points, generator);
+
+  // If c3t3 initialization is not sufficient (may happen if
+  // the user has not specified enough points ), add some surface points
+  bool need_more_init = c3t3.triangulation().dimension() != 3 || !generator.is_default();
+  if(!need_more_init) {
+    CGAL::Mesh_3::C3T3_helpers<C3T3, MeshDomain> helper(c3t3, domain);
+    helper.update_restricted_facets();
+
+    if (c3t3.number_of_facets() == 0) {
+      need_more_init = true;
+    }
+    else
+    {
+      helper.update_restricted_cells();
+      if(c3t3.number_of_cells() == 0) {
+        need_more_init = true;
+      }
+    }
+  }
+  if(need_more_init) {
+    parameters::internal::Initial_points_generator_options<MeshDomain, C3T3> domain_generator =
+        parameters::internal::Initial_points_generator_generator<MeshDomain, C3T3>()();
+    add_points_from_generator<C3T3, MeshDomain, MeshCriteria>(c3t3, domain, nb_initial_points, domain_generator);
+  }
+}
+
 template < typename EdgeCriteria >
 struct Edge_criteria_sizing_field_wrapper
 {
@@ -471,9 +504,9 @@ C3T3 make_mesh_3(const MeshDomain& domain, const MeshCriteria& criteria, const C
     parameters::internal::Manifold_options manifold_options_param = choose_parameter(get_parameter(np, internal_np::manifold_param), parameters::internal::Manifold_options());
 
     using Initial_points_generator_generator = parameters::internal::Initial_points_generator_generator<MeshDomain, C3T3>;
-    using Initial_points = typename Initial_points_generator_generator::Initial_points;
+    using Value_type = typename Initial_points_generator_generator::Value_type;
-    Initial_points empty_vec;
+    std::vector<Value_type> empty_vec;
-    const Initial_points& initial_points
+    const auto& initial_points
           = choose_parameter(get_parameter_reference(np, internal_np::initial_points_param), empty_vec);
     parameters::internal::Initial_points_generator_options initial_points_generator_options_param =
         Initial_points_generator_generator()

STL_Extension/include/CGAL/STL_Extension/internal/mesh_option_classes.h

@@ -168,6 +168,54 @@ private:
   bool b_;
 };
 
+template <typename Value>
+class Input_const_iterator_interface
+{
+public:
+  virtual ~Input_const_iterator_interface() {}
+  virtual const Value& operator*() = 0;
+  virtual Input_const_iterator_interface<Value>* operator++() = 0;
+  virtual bool operator!=(const Input_const_iterator_interface<Value>* other) const = 0;
+  virtual Input_const_iterator_interface<Value>* clone() = 0;
+};
+
+template <typename Value, typename Iterator>
+struct Input_const_iterator_container
+  : Input_const_iterator_interface<Value>
+{
+  typedef Input_const_iterator_container<Value, Iterator> Self;
+public:
+  Input_const_iterator_container(const Iterator& it) : it_(it) {}
+
+  virtual ~Input_const_iterator_container() {}
+
+  virtual const Value& operator*()
+  {
+    return *it_;
+  }
+
+  virtual Input_const_iterator_interface<Value>* operator++()
+  {
+    ++it_;
+    return this;
+  }
+
+  virtual bool operator!=(const Input_const_iterator_interface<Value>* other) const
+  {
+    const Self* other_casted = dynamic_cast<const Self*>(other);
+    if (other_casted == nullptr)
+      return true;
+    return it_ != other_casted->it_;
+  }
+
+  virtual Input_const_iterator_interface<Value>* clone()
+  {
+    return new Input_const_iterator_container<Value, Iterator>(it_);
+  }
+
+private:
+  Iterator it_;
+};
 
 // options is holding the generator (default or the user's one)
 template <typename MeshDomain, typename C3t3>
@@ -175,10 +223,10 @@ struct Initial_points_generator_options
 {
   typedef typename C3t3::Triangulation::Geom_traits::Weighted_point_3 Weighted_point_3;
   typedef typename MeshDomain::Index Index;
-  typedef typename std::vector<std::tuple<Weighted_point_3, int, Index>> Initial_points;
+  typedef typename std::tuple<Weighted_point_3, int, Index> Value_type;
-  typedef typename std::back_insert_iterator<Initial_points> OutputIterator;
+  typedef typename std::back_insert_iterator<std::vector<Value_type>> OutputIterator;
 
-  template <typename Initial_points_generator>
+  template <typename Initial_points_generator, typename Initial_points>
   Initial_points_generator_options(const Initial_points_generator& generator, const Initial_points& initial_points, bool is_default = false)
     : initial_points_generator_no_number_of_points_(generator)
     , initial_points_generator_(generator)
@@ -186,11 +234,14 @@ struct Initial_points_generator_options
   {
     if (initial_points.size() == 0)
     {
-      initial_points_ = nullptr;
+      begin_it = nullptr;
+      end_it = nullptr;
     }
     else
     {
-      initial_points_ = &initial_points;
+      using Iterator_type = typename Initial_points::const_iterator;
+      begin_it = new Input_const_iterator_container<Value_type, Iterator_type>(initial_points.cbegin());
+      end_it = new Input_const_iterator_container<Value_type, Iterator_type>(initial_points.cend());
     }
   }
 
@@ -208,10 +259,11 @@ struct Initial_points_generator_options
 
   OutputIterator add_initial_points(OutputIterator pts) const
   {
-    if (initial_points_ != nullptr)
+    if (begin_it != nullptr && end_it != nullptr)
     {
-      for (const auto& point_tuple : *initial_points_)
+      Input_const_iterator_interface<Value_type>& it = *(begin_it->clone());
-        *pts++ = point_tuple;
+      for (; it != end_it; ++it)
+        *pts++ = *it;
     }
     return pts;
   }
@@ -221,7 +273,8 @@ struct Initial_points_generator_options
 private:
   const std::function<OutputIterator(OutputIterator&,const MeshDomain&,const C3t3&)> initial_points_generator_no_number_of_points_;
   const std::function<OutputIterator(OutputIterator&,const MeshDomain&,const C3t3&,int)> initial_points_generator_;
-  const Initial_points* initial_points_;
+  Input_const_iterator_interface<Value_type>* begin_it;
+  Input_const_iterator_interface<Value_type>* end_it;
   const bool is_default_;
 };
 
@@ -273,7 +326,7 @@ template <typename MeshDomain, typename C3t3>
 struct Initial_points_generator_generator
 {
   typedef typename CGAL::parameters::internal::Initial_points_generator_options<MeshDomain, C3t3> Initial_points_generator_options;
-  typedef typename Initial_points_generator_options::Initial_points Initial_points;
+  typedef typename Initial_points_generator_options::Value_type     Value_type;
   typedef typename Initial_points_generator_options::OutputIterator OutputIterator;
 
   struct Initial_points_generator_domain_traductor
@@ -322,7 +375,7 @@ struct Initial_points_generator_generator
   template <typename InitialPointsGenerator>
   Initial_points_generator_options operator()(const InitialPointsGenerator& initial_points_generator)
   {
-    Initial_points empty_input_features;
+    std::vector<Value_type> empty_input_features;
     return operator()(initial_points_generator, empty_input_features);
   }