Add Surface_mesh_skeletonization (with polyhedron 3 implementation)

CMakeLists.txt

@@ -209,6 +209,7 @@ if ( CGAL_FOUND AND ${CGAL_MAJOR_VERSION} GREATER 4 )
     if (EIGEN3_FOUND)
       include(CGAL_Eigen3_support)
       add_subdirectory(SWIG_CGAL/Polygon_mesh_processing)
+      add_subdirectory(SWIG_CGAL/Surface_mesh_skeletonization)
     else()
       message(STATUS "Eigen 3.2 or later was not found, Polygon_mesh_processing bindings will not be available")
     endif(EIGEN3_FOUND)

SWIG_CGAL/Surface_mesh_skeletonization/CGAL_Surface_mesh_skeletonization.i

@@ -0,0 +1,103 @@
+// ------------------------------------------------------------------------------
+// ------------------------------------------------------------------------------
+
+%define AW3_DOCSTRING
+"SWIG wrapper for the CGAL 3D Surface_mesh_skeletonization provided under the GPL-3.0+ license"
+%enddef
+%module (package="CGAL", docstring=AW3_DOCSTRING) CGAL_Surface_mesh_skeletonization
+
+%include "SWIG_CGAL/common.i"
+Decl_void_type()
+
+SWIG_CGAL_add_java_loadLibrary(CGAL_Surface_mesh_skeletonization)
+SWIG_CGAL_package_common()
+
+%include <std_string.i>
+%include "SWIG_CGAL/Common/Iterator.h"
+%include "SWIG_CGAL/Kernel/typedefs.h"
+%import  "SWIG_CGAL/Common/Macros.h"
+%import  "SWIG_CGAL/Kernel/CGAL_Kernel.i"
+%include  "SWIG_CGAL/Common/Output_iterator_wrapper.h"
+%include  "SWIG_CGAL/Surface_mesh_skeletonization/structs.h" // include so structs are recognized in python
+
+//include files
+%{
+  #include <SWIG_CGAL/Kernel/typedefs.h>
+  #include <SWIG_CGAL/Kernel/Point_3.h>
+  #include <SWIG_CGAL/Polyhedron_3/all_includes.h>
+  #include <SWIG_CGAL/Point_set_3/all_includes.h>
+  #include <SWIG_CGAL/Surface_mesh_skeletonization/impl_polyhedron_3.h>
+  #include <SWIG_CGAL/Common/Output_iterator_wrapper.h>
+%}
+
+%types(Point_3*,Point_3); //needed so that the identifier SWIGTYPE_p_Point_3 is generated
+
+%pragma(java) jniclassimports=%{
+  import CGAL.Kernel.Point_3;
+  import java.util.Iterator;
+  import java.util.Collection;
+  import CGAL.Polyhedron_3.Polyhedron_3;
+%}
+
+%pragma(java) moduleimports=%{
+  import CGAL.Kernel.Point_3;
+  import java.util.Iterator;
+  import java.util.Collection;
+  import CGAL.Polyhedron_3.Polyhedron_3;
+%}
+
+%define Integer_output_iterator boost::function_output_iterator<Container_writer<int,int> > %enddef
+SWIG_CGAL_output_iterator_typemap_in(Integer_output_iterator,int,Integer,int,swig_types[0],"Ljava/lang/Integer;")
+
+//import definitions of Polyhedron objects
+%import "SWIG_CGAL/Polyhedron_3/CGAL_Polyhedron_3.i"
+
+//import Polyhedron_3 wrapper types
+SWIG_CGAL_import_Polyhedron_3_SWIG_wrapper
+
+//import STL vector
+%include "std_vector.i"
+
+// output types
+typedef Point_3 Point_3_wrapper;
+typedef std::vector<Point_3_wrapper> Polyline_wrapper;
+
+// main function
+void surface_mesh_skeletonization(const Polyhedron_3_SWIG_wrapper &polyhedron, std::vector<Polyline_wrapper> &output_skeleton, 
+                                  std::vector<Polyline_wrapper> &output_correspondence, bool debug=false);
+
+void surface_mesh_skeletonization(const Polyhedron_3_SWIG_wrapper &polyhedron, std::vector<Polyline_wrapper> &output_skeleton, 
+                                  std::vector<Polyline_wrapper> &output_correspondence, const local_remeshing_parameters& params1, 
+                                  const algorithm_termination_parameters& params2, const vertex_motion_parameters& params3, bool debug=false);
+
+// get default parameters - allow setting them later
+local_remeshing_parameters get_local_remeshing_parameters(const Polyhedron_3_SWIG_wrapper &polyhedron);
+algorithm_termination_parameters get_algorithm_termination_parameters(const Polyhedron_3_SWIG_wrapper &polyhedron);
+vertex_motion_parameters get_vertex_motion_parameters(const Polyhedron_3_SWIG_wrapper &polyhedron);
+
+// inline implementations - wrap get_data & call implementation from impl_polyhedron_3.h 
+%{
+  void surface_mesh_skeletonization(const Polyhedron_3_SWIG_wrapper &polyhedron,
+    std::vector<Polyline_wrapper> &output_skeleton, std::vector<Polyline_wrapper> &output_correspondence, bool debug = false) {
+    skeletonization(polyhedron.get_data(), output_skeleton, output_correspondence, debug);
+  }
+
+  void surface_mesh_skeletonization(const Polyhedron_3_SWIG_wrapper &polyhedron, std::vector<Polyline_wrapper> &output_skeleton, 
+                                  std::vector<Polyline_wrapper> &output_correspondence, const local_remeshing_parameters& params1, 
+                                  const algorithm_termination_parameters& params2, const vertex_motion_parameters& params3, bool debug = false){
+    skeletonization_advanced(polyhedron.get_data(), output_skeleton, output_correspondence, 
+                             params1, params2, params3, debug);
+  }
+
+  local_remeshing_parameters get_local_remeshing_parameters(const Polyhedron_3_SWIG_wrapper &polyhedron) {
+      return get_local_remeshing_parameters(polyhedron.get_data());
+  }
+  algorithm_termination_parameters get_algorithm_termination_parameters(const Polyhedron_3_SWIG_wrapper &polyhedron) {
+      return get_algorithm_termination_parameters(polyhedron.get_data());
+  }
+  vertex_motion_parameters get_vertex_motion_parameters(const Polyhedron_3_SWIG_wrapper &polyhedron) {
+      return get_vertex_motion_parameters(polyhedron.get_data());
+  }
+
+%}
+

SWIG_CGAL/Surface_mesh_skeletonization/CMakeLists.txt

@@ -0,0 +1,8 @@
+SET (LIBSTOLINKWITH CGAL_Kernel_cpp CGAL::Eigen3_support)
+if (TBB_FOUND)
+  set(LIBSTOLINKWITH ${LIBSTOLINKWITH} TBB::tbb TBB::tbbmalloc Threads::Threads)
+endif()
+# Modules
+ADD_SWIG_CGAL_JAVA_MODULE   ( Surface_mesh_skeletonization ${LIBSTOLINKWITH} )
+ADD_SWIG_CGAL_PYTHON_MODULE ( Surface_mesh_skeletonization ${LIBSTOLINKWITH} )
+ADD_SWIG_CGAL_RUBY_MODULE   ( Surface_mesh_skeletonization ${LIBSTOLINKWITH} )

SWIG_CGAL/Surface_mesh_skeletonization/README.txt

@@ -0,0 +1,8 @@
+See https://github.com/CGAL/cgal-swig-bindings/wiki/User_package
+
+Files list:
+CGAL_Surface_mesh_skeletonization.i - Swig interface file.
+structs.h - Added classes/structs. For now, contains the algorithm parameters' structs.
+common.h - template functions, used by impl*.h files. Contains the algorithm running & output printing functions.
+impl_polyhedron_3.h - algorithm implementation
+CMakeLists.txt - how to build package. Note: requires CGAL::Eigen3_support in order to compile.

SWIG_CGAL/Surface_mesh_skeletonization/common.h

@@ -0,0 +1,143 @@
+#ifndef COMMON_H
+#define COMMON_H
+
+#include <CGAL/boost/graph/split_graph_into_polylines.h>
+#include "structs.h"
+
+// All functions here are template to allow in the future additional skeleton & mesh types (see examples in CGAL cpp files)
+
+////////////////////////////////////////////
+// Advanced algorithm running
+////////////////////////////////////////////
+
+template <typename TSkeletonizaion>
+void set_params(TSkeletonizaion& mcs,  const local_remeshing_parameters& params1) {
+  if (mcs.max_triangle_angle() != params1.max_triangle_angle) {
+      mcs.set_max_triangle_angle(params1.max_triangle_angle);
+  }
+  if (mcs.min_edge_length() != params1.min_edge_length) {
+      mcs.set_min_edge_length(params1.min_edge_length);
+  }
+}
+
+template <typename TSkeletonizaion>
+void set_params(TSkeletonizaion& mcs,  const algorithm_termination_parameters& params2) {
+  if (mcs.max_iterations() != params2.max_iterations) {
+      mcs.set_max_iterations(params2.max_iterations);
+  }
+  if (mcs.area_variation_factor() != params2.area_variation_factor) {
+      mcs.set_area_variation_factor(params2.area_variation_factor);
+  }
+}
+
+template <typename TSkeletonizaion>
+void set_params(TSkeletonizaion& mcs,  const vertex_motion_parameters& params3) {
+  if (mcs.quality_speed_tradeoff() != params3.quality_speed_tradeoff) {
+      mcs.set_quality_speed_tradeoff(params3.quality_speed_tradeoff);
+  }
+  if (mcs.medially_centered_speed_tradeoff() != params3.medially_centered_speed_tradeoff) {
+      mcs.set_medially_centered_speed_tradeoff(params3.medially_centered_speed_tradeoff);
+  }
+  if (mcs.is_medially_centered() != params3.is_medially_centered) {
+      mcs.set_is_medially_centered(params3.is_medially_centered);
+  }
+}
+
+template <typename TSkeleton, typename TMesh> 
+void run(const TMesh& tmesh, TSkeleton& skeleton, const local_remeshing_parameters& params1,
+         const algorithm_termination_parameters& params2, const vertex_motion_parameters& params3, bool debug = false) {
+  CGAL::Mean_curvature_flow_skeletonization<TMesh> mcs(tmesh);
+
+  set_params(mcs, params1);
+  set_params(mcs, params2);
+  set_params(mcs, params3);
+
+  /*// 1. Contract the mesh by mean curvature flow.
+  mcs.contract_geometry();
+
+  // 2. Collapse short edges and split bad triangles.
+  mcs.collapse_edges();
+  mcs.split_faces();
+
+  // 3. Fix degenerate vertices.
+  mcs.detect_degeneracies();
+
+  // Perform the above three steps in one iteration.
+  mcs.contract();
+  */
+  // Iteratively apply step 1 to 3 until convergence.
+  mcs.contract_until_convergence();
+
+  // Convert the contracted mesh into a curve skeleton and
+  // get the correspondent surface points
+  mcs.convert_to_skeleton(skeleton);
+}
+
+
+////////////////////////////////////////////
+// Output writing
+////////////////////////////////////////////
+
+template <typename TSkeleton, typename Tvertex_descriptor, typename Tout>
+struct Display_polylines{ 
+  const TSkeleton& skeleton;
+  std::vector<Tout>& out;
+
+  Display_polylines(const TSkeleton& skeleton, std::vector<Tout>& out)
+    : skeleton(skeleton), out(out)
+  {}
+
+  void start_new_polyline(){
+      out.push_back(Tout());
+  }
+  void add_node(const Tvertex_descriptor& v){
+    out.back().push_back(skeleton[v].point);
+  }
+  void end_polyline() {
+  }
+};
+
+template <typename TSkeleton, typename TMesh, typename Tout> 
+void write_results(const TSkeleton& skeleton, const TMesh& tmesh, std::vector<Tout>& output, std::vector<Tout>& output_correspondence, bool debug) {
+
+  if (debug) {
+    std::cout << "Number of vertices of the skeleton: " << boost::num_vertices(skeleton) << "\n";
+    std::cout << "Number of edges of the skeleton: " << boost::num_edges(skeleton) << "\n";
+  }
+
+  Display_polylines<TSkeleton, typename TSkeleton::vertex_descriptor, Tout> display(skeleton,output);
+  CGAL::split_graph_into_polylines(skeleton, display);
+
+  // Output skeleton points and the corresponding surface points
+  for (const auto & v : CGAL::make_range(vertices(skeleton)))
+    for (const auto & vd : skeleton[v].vertices)
+    {
+      output_correspondence.push_back(Tout());
+      output_correspondence.back().push_back(skeleton[v].point);
+      output_correspondence.back().push_back(get(CGAL::vertex_point, tmesh, vd));
+    }
+}
+
+////////////////////////////////////////////
+// Get algorithm parameters (allow setting them)
+////////////////////////////////////////////
+
+template <typename TMesh> 
+local_remeshing_parameters get_local_remeshing_parameters(const TMesh& tmesh) {
+    CGAL::Mean_curvature_flow_skeletonization<TMesh> mcs(tmesh);
+    return local_remeshing_parameters{mcs.max_triangle_angle(), mcs.min_edge_length()};
+}
+
+template <typename TMesh> 
+algorithm_termination_parameters get_algorithm_termination_parameters(const TMesh& tmesh) {
+    CGAL::Mean_curvature_flow_skeletonization<TMesh> mcs(tmesh);
+    return algorithm_termination_parameters{mcs.max_iterations(), mcs.area_variation_factor()};
+}
+
+template <typename TMesh> 
+vertex_motion_parameters get_vertex_motion_parameters(const TMesh& tmesh) {
+    CGAL::Mean_curvature_flow_skeletonization<TMesh> mcs(tmesh);
+    return vertex_motion_parameters{mcs.quality_speed_tradeoff(), mcs.is_medially_centered(), mcs.medially_centered_speed_tradeoff()};
+}
+
+#endif

SWIG_CGAL/Surface_mesh_skeletonization/impl_polyhedron_3.h

@@ -0,0 +1,79 @@
+#include <CGAL/Polyhedron_3.h>
+#include <CGAL/extract_mean_curvature_flow_skeleton.h>
+#include <SWIG_CGAL/Kernel/Point_3.h>
+#include <fstream>
+#include <vector>
+
+#include "structs.h"
+#include "common.h"
+
+// output types
+typedef Point_3 Point_3_wrapper;
+typedef std::vector<Point_3_wrapper> Polyline_wrapper;
+
+// input types
+typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel;
+typedef CGAL::Polyhedron_3<Kernel, CGAL::Polyhedron_items_with_id_3> Polyhedron;
+
+// algorithm types
+typedef CGAL::Mean_curvature_flow_skeletonization<Polyhedron> Skeletonization;
+typedef Skeletonization::Skeleton                             Skeleton;
+
+// define a double-precision representation of seconds
+using fsecs = std::chrono::duration<double, std::chrono::seconds::period>;
+
+int skeletonization(const Polyhedron & tmesh, std::vector<Polyline_wrapper> & output_skeleton, 
+                    std::vector<Polyline_wrapper> & output_correspondence, bool debug = false)
+{
+  if (!CGAL::is_triangle_mesh(tmesh))
+  {
+    std::cout << "Input geometry is not triangulated." << std::endl;
+    return EXIT_FAILURE;
+  }
+
+  auto begin = std::chrono::steady_clock::now();
+
+  Skeleton skeleton;
+  CGAL::extract_mean_curvature_flow_skeleton(tmesh, skeleton);
+
+  auto end = std::chrono::steady_clock::now();
+
+  write_results(skeleton, tmesh, output_skeleton, output_correspondence, debug=debug);
+
+  if (debug) {
+    std::cout << "Skeletonization took: " << std::chrono::duration_cast<fsecs>(end - begin).count() << " seconds" << std::endl;
+  }
+
+  return EXIT_SUCCESS;
+}
+
+int skeletonization_advanced(const Polyhedron & tmesh, std::vector<Polyline_wrapper> & output_skeleton, 
+                             std::vector<Polyline_wrapper> & output_correspondence, const local_remeshing_parameters& params1, 
+			     const algorithm_termination_parameters& params2, const vertex_motion_parameters& params3, bool debug = false)
+{
+  if (!CGAL::is_triangle_mesh(tmesh))
+  {
+    std::cout << "Input geometry is not triangulated." << std::endl;
+    return EXIT_FAILURE;
+  }
+
+  if (debug) {
+    std::cout << "Parameters: " << params1.__str__() << params2.__str__() << params3.__str__() << std::endl;
+  }
+
+  auto begin = std::chrono::steady_clock::now();
+
+  Skeleton skeleton;
+  run(tmesh, skeleton, params1, params2, params3, debug);
+
+  auto end = std::chrono::steady_clock::now();
+
+  write_results(skeleton, tmesh, output_skeleton, output_correspondence, debug);
+
+  if (debug) {
+    std::cout << "Skeletonization took: " << std::chrono::duration_cast<fsecs>(end - begin).count() << " seconds" << std::endl;
+  }
+
+  return EXIT_SUCCESS;
+}
+

SWIG_CGAL/Surface_mesh_skeletonization/structs.h

@@ -0,0 +1,41 @@
+#ifndef STRUCTS_H
+#define STRUCTS_H
+
+#include <sstream>
+struct local_remeshing_parameters
+{
+  double max_triangle_angle;
+  double min_edge_length;
+
+  std::string __str__() const { // match the python __str__ saved word
+    std::ostringstream oss(std::ostringstream::out);
+    oss << "{'max_triangle_angle':" << (*this).max_triangle_angle << ", 'min_edge_length':" <<  (*this).min_edge_length << "}";
+    return oss.str();
+  }
+};
+
+struct algorithm_termination_parameters
+{
+  std::size_t max_iterations;
+  double area_variation_factor;
+  std::string __str__() const { // match the python __str__ saved word
+    std::ostringstream oss(std::ostringstream::out);
+    oss << "{'max_iterations':" << (*this).max_iterations << ", 'area_variation_factor':" <<  (*this).area_variation_factor << "}";
+    return oss.str();
+  }
+};
+
+struct vertex_motion_parameters
+{
+  double quality_speed_tradeoff;
+  bool is_medially_centered;
+  double medially_centered_speed_tradeoff;
+  std::string __str__() const { // match the python __str__ saved word
+    std::ostringstream oss(std::ostringstream::out);
+    oss << "{'quality_speed_tradeoff':" << (*this).quality_speed_tradeoff << ", 'is_medially_centered':" <<  (*this).is_medially_centered 
+        << ", 'medially_centered_speed_tradeoff':" <<  (*this).medially_centered_speed_tradeoff  << "}";
+    return oss.str();
+  }
+};
+
+#endif