polyhedron helper functions

include/igl/copyleft/cgal/convex_hull.cpp

@@ -40,6 +40,8 @@ IGL_INLINE void igl::copyleft::cgal::convex_hull(
       }
       Polyhedron_3 poly;
       CGAL::convex_hull_3(points.begin(),points.end(),poly);
+      // https://stackoverflow.com/a/27535366/148668
+      // It always should be.
       assert(poly.is_pure_triangle() && "Assuming CGAL outputs a triangle mesh");
       polyhedron_to_mesh(poly,W,G);
       break;

include/igl/copyleft/cgal/join_coplanar_neighboring_facets.cpp

@@ -0,0 +1,35 @@
+#include "join_coplanar_neighboring_facets.h"
+#include <CGAL/Polyhedron_3.h>
+
+template <typename Polyhedron>
+IGL_INLINE void igl::copyleft::cgal::join_coplanar_neighboring_facets(
+  Polyhedron & poly)
+{
+  using Plane_3 = typename Polyhedron::Traits::Plane_3;
+  for(auto iter = poly.edges_begin(); iter != poly.edges_end();)
+  {
+    auto e = iter++;
+    auto f1 = e->facet();
+    auto f2 = e->opposite()->facet();
+    auto p1 = Plane_3(
+      f1->halfedge()->vertex()->point(), 
+      f1->halfedge()->next()->vertex()->point(), 
+      f1->halfedge()->next()->next()->vertex()->point());
+    auto p2 = Plane_3(
+      f2->halfedge()->vertex()->point(), 
+      f2->halfedge()->next()->vertex()->point(), 
+      f2->halfedge()->next()->next()->vertex()->point());
+    if(p1 == p2)
+    {
+      poly.join_facet(e);
+    }
+  }
+}
+
+
+#ifdef IGL_STATIC_LIBRARY
+// Explicit template instantiation
+#include <CGAL/Polyhedron_items_with_id_3.h>
+#include <CGAL/Simple_cartesian.h>
+template void igl::copyleft::cgal::join_coplanar_neighboring_facets<CGAL::Polyhedron_3<CGAL::Simple_cartesian<double>, CGAL::Polyhedron_items_with_id_3, CGAL::HalfedgeDS_default, std::allocator<int>>>(CGAL::Polyhedron_3<CGAL::Simple_cartesian<double>, CGAL::Polyhedron_items_with_id_3, CGAL::HalfedgeDS_default, std::allocator<int>>&);
+#endif

include/igl/copyleft/cgal/join_coplanar_neighboring_facets.h

@@ -0,0 +1,20 @@
+#ifndef IGL_COPYLEFT_CGAL_JOIN_COPLANAR_NEIGHBORING_FACETS_H
+#define IGL_COPYLEFT_CGAL_JOIN_COPLANAR_NEIGHBORING_FACETS_H
+#include "../../igl_inline.h"
+namespace igl
+{
+  namespace copyleft
+  {
+    namespace cgal
+    {
+      template <typename Polyhedron>
+      IGL_INLINE void join_coplanar_neighboring_facets(
+        Polyhedron & poly);
+    }
+  }
+}
+
+#ifndef IGL_STATIC_LIBRARY
+#include "join_coplanar_neighboring_facets.cpp"
+#endif
+#endif

include/igl/copyleft/cgal/polyhedron_to_mesh.cpp

@@ -8,21 +8,21 @@
 #include "polyhedron_to_mesh.h"
 #include <CGAL/Polyhedron_3.h>
 #include "assign_scalar.h"
+#include <unordered_map>
 
-template <
-  typename Polyhedron,
-  typename DerivedV,
-  typename DerivedF>
-IGL_INLINE void igl::copyleft::cgal::polyhedron_to_mesh(
-  const Polyhedron & poly,
-  Eigen::PlainObjectBase<DerivedV> & V,
-  Eigen::PlainObjectBase<DerivedF> & F)
+namespace
 {
-  V.resize(poly.size_of_vertices(),3);
-  F.resize(poly.size_of_facets(),3);
-  typedef typename Polyhedron::Vertex_const_iterator Vertex_iterator;
-  std::map<Vertex_iterator,size_t> vertex_to_index;
+
+  template <
+    typename Polyhedron,
+    typename DerivedV>
+  inline void polyhedron_to_mesh_helper_vertices(
+    const Polyhedron & poly,
+    Eigen::PlainObjectBase<DerivedV> & V,
+    std::unordered_map<typename Polyhedron::Vertex_const_iterator,size_t> & vertex_to_index)
   {
+    using namespace igl::copyleft::cgal;
+    V.resize(poly.size_of_vertices(),3);
     size_t v = 0;
     for(
       typename Polyhedron::Vertex_const_iterator p = poly.vertices_begin();
@@ -36,8 +36,23 @@ IGL_INLINE void igl::copyleft::cgal::polyhedron_to_mesh(
       v++;
     }
   }
+}
+
+template <
+  typename Polyhedron,
+  typename DerivedV,
+  typename DerivedF>
+IGL_INLINE void igl::copyleft::cgal::polyhedron_to_mesh(
+  const Polyhedron & poly,
+  Eigen::PlainObjectBase<DerivedV> & V,
+  Eigen::PlainObjectBase<DerivedF> & F)
+{
+  std::unordered_map<typename Polyhedron::Vertex_const_iterator,size_t> vertex_to_index;
+  polyhedron_to_mesh_helper_vertices(poly,V,vertex_to_index);
+
+  F.resize(poly.size_of_facets(),3);
   {
-    size_t f = 0;
+    int f = 0;
     for(
       typename Polyhedron::Facet_const_iterator facet = poly.facets_begin();
       facet != poly.facets_end();
@@ -47,24 +62,74 @@ IGL_INLINE void igl::copyleft::cgal::polyhedron_to_mesh(
         facet->facet_begin();
       // Facets in polyhedral surfaces are at least triangles.
       assert(CGAL::circulator_size(he) == 3 && "Facets should be triangles");
-      size_t c = 0;
+      int c = 0;
       do {
         //// This is stooopidly slow
         // F(f,c) = std::distance(poly.vertices_begin(), he->vertex());
         F(f,c) = vertex_to_index[he->vertex()];
         c++;
-      } while ( ++he != facet->facet_begin());
+      } while (++he != facet->facet_begin());
       f++;
     }
   }
+
+}
+
+template <
+  typename Polyhedron,
+  typename DerivedV,
+  typename DerivedI,
+  typename DerivedC
+  >
+IGL_INLINE void igl::copyleft::cgal::polyhedron_to_mesh(
+  const Polyhedron & poly,
+  Eigen::PlainObjectBase<DerivedV> & V,
+  Eigen::PlainObjectBase<DerivedI> & I,
+  Eigen::PlainObjectBase<DerivedC> & C)
+{
+  // I think even the unordered_map can be avoided if we const expr detect that 
+  // the Polyhedron uses items with id and then use that id as the index.
+  //
+  // It would... But the caller needs to not only use Polyhedron_items_with_id_3
+  // but needs to loop over vertices and set v->id() = i++ before calling this
+  // function (or this function needs to do that but then poly can't be a const
+  // input anymore).
+  //
+  // Probably there should be a special version of this function where the user
+  // can promise that ids have been set.
+  std::unordered_map<typename Polyhedron::Vertex_const_iterator,size_t> vertex_to_index;
+  polyhedron_to_mesh_helper_vertices(poly,V,vertex_to_index);
+  C.resize(poly.size_of_facets() + 1);
+  I.resize(poly.size_of_halfedges());
+  {
+    int f = 0;
+    C(0) = 0;
+    int i = 0;
+    for (auto fit = poly.facets_begin(); fit != poly.facets_end(); ++fit, ++f)
+    {
+      auto he = fit->halfedge();
+      const auto he0 = he;
+      int deg = 0;
+      do
+      {
+        //I(i++) = std::distance(poly.vertices_begin(), he->vertex());
+        I(i++) = vertex_to_index[he->vertex()];
+        deg++;
+        he = he->next();
+      } while (he != he0);
+      C(f+1) = C(f) + deg;
+    }
+  }
 }
 
 #ifdef IGL_STATIC_LIBRARY
-// Explicit template instantiation
-// generated by autoexplicit.sh
 #include <CGAL/Simple_cartesian.h>
 #include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
 #include <CGAL/Polyhedron_items_with_id_3.h>
+// Explicit template instantiation
+// generated by autoexplicit.sh
+template void igl::copyleft::cgal::polyhedron_to_mesh<CGAL::Polyhedron_3<CGAL::Simple_cartesian<double>, CGAL::Polyhedron_items_with_id_3, CGAL::HalfedgeDS_default, std::allocator<int>>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>>(CGAL::Polyhedron_3<CGAL::Simple_cartesian<double>, CGAL::Polyhedron_items_with_id_3, CGAL::HalfedgeDS_default, std::allocator<int>> const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1>>&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1>>&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1>>&);
+// generated by autoexplicit.sh
 template void igl::copyleft::cgal::polyhedron_to_mesh<CGAL::Polyhedron_3<CGAL::Epick, CGAL::Polyhedron_items_3, CGAL::HalfedgeDS_default, std::allocator<int> >, Eigen::Matrix<double, -1, 2, 0, -1, 2>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(CGAL::Polyhedron_3<CGAL::Epick, CGAL::Polyhedron_items_3, CGAL::HalfedgeDS_default, std::allocator<int> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 2, 0, -1, 2> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
 template void igl::copyleft::cgal::polyhedron_to_mesh<CGAL::Polyhedron_3<CGAL::Epick, CGAL::Polyhedron_items_3, CGAL::HalfedgeDS_default, std::allocator<int> >, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(CGAL::Polyhedron_3<CGAL::Epick, CGAL::Polyhedron_items_3, CGAL::HalfedgeDS_default, std::allocator<int> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
 template void igl::copyleft::cgal::polyhedron_to_mesh<CGAL::Polyhedron_3<CGAL::Simple_cartesian<double>,CGAL::Polyhedron_items_with_id_3, CGAL::HalfedgeDS_default, std::allocator<int> >, Eigen::Matrix<double, -1, -1, 0,-1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(CGAL::Polyhedron_3<CGAL::Simple_cartesian<double>,CGAL::Polyhedron_items_with_id_3, CGAL::HalfedgeDS_default, std::allocator<int> > const&,Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0,-1, -1> >&);

include/igl/copyleft/cgal/polyhedron_to_mesh.h

@@ -16,12 +16,13 @@ namespace igl
   {
     namespace cgal
     {
-      /// Convert a CGAL Polyhedron to a mesh (V,F)
+      /// Convert a CGAL Polyhedron (assumed to be all triangles) to a mesh (V,F)
       ///
       /// @tparam Polyhedron  CGAL Polyhedron type (e.g. Polyhedron_3)
       /// @param[in] poly  cgal polyhedron
       /// @param[out] V  #V by 3 list of vertex positions
       /// @param[out] F  #F by 3 list of triangle indices
+      ///
       template <
         typename Polyhedron,
         typename DerivedV,
@@ -30,6 +31,26 @@ namespace igl
         const Polyhedron & poly,
         Eigen::PlainObjectBase<DerivedV> & V,
         Eigen::PlainObjectBase<DerivedF> & F);
+      /// Convert a CGAL Polyhedron to a polygon mesh (V,I,C)
+      ///
+      /// @tparam Polyhedron  CGAL Polyhedron type (e.g. Polyhedron_3)
+      /// @param[in] poly  cgal polyhedron
+      /// @param[out] V  #V by 3 list of vertex positions
+      /// @param[out] I  #I vectorized list of polygon corner indices into rows of some matrix V
+      /// @param[out] C  #polygons+1 list of cumulative polygon sizes so that C(i+1)-C(i) =
+      ///     size of the ith polygon, and so I(C(i)) through I(C(i+1)-1) are the
+      ///     indices of the ith polygon
+      template <
+        typename Polyhedron,
+        typename DerivedV,
+        typename DerivedI,
+        typename DerivedC
+        >
+      IGL_INLINE void polyhedron_to_mesh(
+        const Polyhedron & poly,
+        Eigen::PlainObjectBase<DerivedV> & V,
+        Eigen::PlainObjectBase<DerivedI> & I,
+        Eigen::PlainObjectBase<DerivedC> & C);
     }
   }
 }

tests/include/igl/copyleft/cgal/convex_hull.cpp

@@ -0,0 +1,20 @@
+#include <test_common.h>
+#include <igl/copyleft/cgal/convex_hull.h>
+
+TEST_CASE(
+  "igl_copyleft_cgal_convex_hull: cube",
+  "[igl/copyleft/cgal/]")
+{
+  Eigen::MatrixXd V(8,3);
+  V << 0,0,0,
+       1,0,0,
+       1,1,0,
+       0,1,0,
+       0,0,1,
+       1,0,1,
+       1,1,1,
+       0,1,1;
+  Eigen::MatrixXi F;
+  igl::copyleft::cgal::convex_hull(V,F);
+  REQUIRE(F.rows() == 12); // cube has 12 triangles
+}

tests/include/igl/copyleft/cgal/polyhedron_to_mesh.cpp

@@ -0,0 +1,58 @@
+#include <test_common.h>
+#include <igl/copyleft/cgal/polyhedron_to_mesh.h>
+#include <igl/copyleft/cgal/mesh_to_polyhedron.h>
+#include <igl/copyleft/cgal/join_coplanar_neighboring_facets.h>
+#include <CGAL/Simple_cartesian.h>
+#include <CGAL/Polyhedron_3.h>
+#include <CGAL/Polyhedron_items_with_id_3.h>
+
+TEST_CASE(
+  "igl_copyleft_cgal_polyhedron_to_mesh: cube",
+  "[igl/copyleft/cgal/]")
+{
+  Eigen::MatrixXd V(8,3);
+  V << 0,0,0,
+       1,0,0,
+       1,1,0,
+       0,1,0,
+       0,0,1,
+       1,0,1,
+       1,1,1,
+       0,1,1;
+  Eigen::MatrixXi F(12,3);
+  F <<
+    0,1,2,
+    0,2,3,
+    4,7,5,
+    5,7,6,
+    0,4,1,
+    1,4,5,
+    1,5,2,
+    2,5,6,
+    2,6,3,
+    3,6,7,
+    3,7,0,
+    0,7,4;
+  CGAL::Polyhedron_3<
+    CGAL::Simple_cartesian<double>, 
+    CGAL::Polyhedron_items_with_id_3, 
+    CGAL::HalfedgeDS_default, 
+    std::allocator<int> > 
+    poly;
+  REQUIRE( igl::copyleft::cgal::mesh_to_polyhedron(V,F,poly));
+  Eigen::MatrixXd V_out;
+  Eigen::MatrixXi F_out;
+  igl::copyleft::cgal::polyhedron_to_mesh(poly,V_out,F_out);
+  REQUIRE( V_out.rows() == 8 );
+  REQUIRE( F_out.rows() == 12 );
+  test_common::assert_eq(V,V_out);
+  test_common::assert_eq(F,F_out);
+  igl::copyleft::cgal::join_coplanar_neighboring_facets(poly);
+  REQUIRE( poly.size_of_facets() == 6 );
+  Eigen::VectorXi I,C;
+  igl::copyleft::cgal::polyhedron_to_mesh(poly,V_out,I,C);
+  REQUIRE( V_out.rows() == 8 );
+  REQUIRE( C.size() == 6+1 );
+  auto sizes = (C.tail(C.size()-1) - C.head(C.size()-1)).eval();
+  REQUIRE((sizes.array() == 4).all());
+}