pull out cdt code into seperate functions

Former-commit-id: d568bd30a28dcfe25f4de8c19f03b68ee0d9c416

include/igl/copyleft/cgal/SelfIntersectMesh.h

@@ -72,9 +72,9 @@ namespace igl
           typedef CGAL::Segment_2<Kernel>  Segment_2; 
           typedef CGAL::Triangle_2<Kernel> Triangle_2; 
           // 2D Constrained Delaunay Triangulation types
-          typedef CGAL::Triangulation_vertex_base_2<Kernel>  TVB_2;
-          typedef CGAL::Constrained_triangulation_face_base_2<Kernel> CTFB_2;
-          typedef CGAL::Triangulation_data_structure_2<TVB_2,CTFB_2> TDS_2;
+          //typedef CGAL::Triangulation_vertex_base_2<Kernel>  TVB_2;
+          //typedef CGAL::Constrained_triangulation_face_base_2<Kernel> CTFB_2;
+          //typedef CGAL::Triangulation_data_structure_2<TVB_2,CTFB_2> TDS_2;
           typedef CGAL::Exact_intersections_tag Itag;
           // Axis-align boxes for all-pairs self-intersection detection
           typedef std::vector<Triangle_3> Triangles;

include/igl/copyleft/cgal/insert_into_cdt.cpp

@@ -0,0 +1,69 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+// 
+// Copyright (C) 2016 Alec Jacobson
+// 
+// This Source Code Form is subject to the terms of the Mozilla Public License 
+// v. 2.0. If a copy of the MPL was not distributed with this file, You can 
+// obtain one at http://mozilla.org/MPL/2.0/.
+#include "insert_into_cdt.h"
+#include <CGAL/Point_3.h>
+#include <CGAL/Segment_3.h>
+#include <CGAL/Triangle_3.h>
+
+template <typename Kernel>
+IGL_INLINE void igl::copyleft::cgal::insert_into_cdt(
+  const CGAL::Object & obj,
+  const CGAL::Plane_3<Kernel> & P,
+  CGAL::Constrained_triangulation_plus_2<
+    CGAL::Constrained_Delaunay_triangulation_2<
+      Kernel,
+      CGAL::Triangulation_data_structure_2<
+        CGAL::Triangulation_vertex_base_2<Kernel>,
+        CGAL::Constrained_triangulation_face_base_2< Kernel>
+      >,
+      CGAL::Exact_intersections_tag
+    > 
+  > 
+  & cdt)
+{
+  typedef CGAL::Point_3<Kernel>    Point_3;
+  typedef CGAL::Segment_3<Kernel>  Segment_3; 
+  typedef CGAL::Triangle_3<Kernel> Triangle_3; 
+  typedef CGAL::Plane_3<Kernel>    Plane_3;
+
+  if(const Segment_3 *iseg = CGAL::object_cast<Segment_3 >(&obj)) 
+  {
+    // Add segment constraint
+    cdt.insert_constraint( P.to_2d(iseg->vertex(0)),P.to_2d(iseg->vertex(1)));
+  }else if(const Point_3 *ipoint = CGAL::object_cast<Point_3 >(&obj)) 
+  {
+    // Add point
+    cdt.insert(P.to_2d(*ipoint));
+  } else if(const Triangle_3 *itri = CGAL::object_cast<Triangle_3 >(&obj)) 
+  {
+    // Add 3 segment constraints
+    cdt.insert_constraint( P.to_2d(itri->vertex(0)),P.to_2d(itri->vertex(1)));
+    cdt.insert_constraint( P.to_2d(itri->vertex(1)),P.to_2d(itri->vertex(2)));
+    cdt.insert_constraint( P.to_2d(itri->vertex(2)),P.to_2d(itri->vertex(0)));
+  } else if(const std::vector<Point_3 > *polyp = 
+      CGAL::object_cast< std::vector<Point_3 > >(&obj)) 
+  {
+    const std::vector<Point_3 > & poly = *polyp;
+    const size_t m = poly.size();
+    assert(m>=2);
+    for(size_t p = 0;p<m;p++)
+    {
+      const size_t np = (p+1)%m;
+      cdt.insert_constraint(P.to_2d(poly[p]),P.to_2d(poly[np]));
+    }
+  }else {
+    throw std::runtime_error("Unknown intersection object!");
+  }
+}
+
+#ifdef IGL_STATIC_LIBRARY
+// Explicit template specialization
+#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
+template void igl::copyleft::cgal::insert_into_cdt<CGAL::Epick>(CGAL::Object const&, CGAL::Plane_3<CGAL::Epick> const&, CGAL::Constrained_triangulation_plus_2<CGAL::Constrained_Delaunay_triangulation_2<CGAL::Epick, CGAL::Triangulation_data_structure_2<CGAL::Triangulation_vertex_base_2<CGAL::Epick, CGAL::Triangulation_ds_vertex_base_2<void> >, CGAL::Constrained_triangulation_face_base_2<CGAL::Epick, CGAL::Triangulation_face_base_2<CGAL::Epick, CGAL::Triangulation_ds_face_base_2<void> > > >, CGAL::Exact_intersections_tag> >&);
+template void igl::copyleft::cgal::insert_into_cdt<CGAL::Epeck>(CGAL::Object const&, CGAL::Plane_3<CGAL::Epeck> const&, CGAL::Constrained_triangulation_plus_2<CGAL::Constrained_Delaunay_triangulation_2<CGAL::Epeck, CGAL::Triangulation_data_structure_2<CGAL::Triangulation_vertex_base_2<CGAL::Epeck, CGAL::Triangulation_ds_vertex_base_2<void> >, CGAL::Constrained_triangulation_face_base_2<CGAL::Epeck, CGAL::Triangulation_face_base_2<CGAL::Epeck, CGAL::Triangulation_ds_face_base_2<void> > > >, CGAL::Exact_intersections_tag> >&);
+#endif

include/igl/copyleft/cgal/insert_into_cdt.h

@@ -0,0 +1,59 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+// 
+// Copyright (C) 2016 Alec Jacobson
+// 
+// This Source Code Form is subject to the terms of the Mozilla Public License 
+// v. 2.0. If a copy of the MPL was not distributed with this file, You can 
+// obtain one at http://mozilla.org/MPL/2.0/.
+#ifndef IGL_COPYLEFT_CGAL_INSERT_INTO_CDT_H
+#define IGL_COPYLEFT_CGAL_INSERT_INTO_CDT_H
+#include "../../igl_inline.h"
+
+#include <CGAL/Plane_3.h>
+#include <CGAL/Constrained_Delaunay_triangulation_2.h>
+#include <CGAL/Constrained_triangulation_plus_2.h>
+
+namespace igl
+{
+  namespace copyleft
+  {
+    namespace cgal
+    {
+      // Given a current 2D constrained Delaunay triangulation (cdt), insert a
+      // 3D "object" (e.g., resulting from intersecting two triangles) into the
+      // cdt, by projecting it via the given plane (P) and adding appropriate
+      // constraints.
+      //
+      // Inputs:
+      //   obj  CGAL::Object representing a vertex, segment, or (convex)
+      //     polygon. All vertices should lie on the plane P. If not, then this
+      //     adds the _projection_ of this object to the cdt and that might not
+      //     be what you wanted to do.
+      //   P  plane obj lies on and upon which the cdt is being performed
+      //   cdt  current CDT, see output
+      // Outputs:
+      //   cdt  CDT updated to contain constraints for the given object
+      // 
+      template <typename Kernel>
+      IGL_INLINE void insert_into_cdt(
+        const CGAL::Object & obj,
+        const CGAL::Plane_3<Kernel> & P,
+        CGAL::Constrained_triangulation_plus_2<
+          CGAL::Constrained_Delaunay_triangulation_2<
+            Kernel,
+            CGAL::Triangulation_data_structure_2<
+              CGAL::Triangulation_vertex_base_2<Kernel>,
+              CGAL::Constrained_triangulation_face_base_2< Kernel>
+            >,
+            CGAL::Exact_intersections_tag
+          > 
+        > 
+        & cdt);
+    }
+  }
+}
+
+#ifndef IGL_STATIC_LIBRARY
+#  include "insert_into_cdt.cpp"
+#endif
+#endif

include/igl/copyleft/cgal/projected_cdt.cpp

@@ -0,0 +1,81 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+// 
+// Copyright (C) 2016 Alec Jacobson
+// 
+// This Source Code Form is subject to the terms of the Mozilla Public License 
+// v. 2.0. If a copy of the MPL was not distributed with this file, You can 
+// obtain one at http://mozilla.org/MPL/2.0/.
+#include "projected_cdt.h"
+#include "insert_into_cdt.h"
+#include "assign_scalar.h"
+#include "../../list_to_matrix.h"
+template <typename Kernel, typename Index>
+IGL_INLINE void igl::copyleft::cgal::projected_cdt(
+  const std::vector<CGAL::Object> & objects,
+  const CGAL::Plane_3<Kernel> & P,
+  std::vector<CGAL::Point_3<Kernel> >& vertices,
+  std::vector<std::vector<Index> >& faces)
+{
+  typedef CGAL::Point_3<Kernel>    Point_3;
+  typedef CGAL::Segment_3<Kernel>  Segment_3; 
+  typedef CGAL::Triangle_3<Kernel> Triangle_3; 
+  typedef CGAL::Plane_3<Kernel>    Plane_3;
+  typedef CGAL::Triangulation_vertex_base_2<Kernel>  TVB_2;
+  typedef CGAL::Constrained_triangulation_face_base_2<Kernel> CTFB_2;
+  typedef CGAL::Triangulation_data_structure_2<TVB_2,CTFB_2> TDS_2;
+  typedef CGAL::Exact_intersections_tag Itag;
+  typedef CGAL::Constrained_Delaunay_triangulation_2<Kernel,TDS_2,Itag> CDT_2;
+  typedef CGAL::Constrained_triangulation_plus_2<CDT_2> CDT_plus_2;
+  CDT_plus_2 cdt;
+  for(const auto & obj : objects) insert_into_cdt(obj,P,cdt);
+  // Read off vertices of the cdt, remembering index
+  std::map<typename CDT_plus_2::Vertex_handle,Index> v2i;
+  size_t count=0;
+  for (
+    auto itr = cdt.finite_vertices_begin();
+    itr != cdt.finite_vertices_end(); 
+    itr++) 
+  {
+    vertices.push_back(P.to_3d(itr->point()));
+    v2i[itr] = count;
+    count++;
+  }
+  // Read off faces and store index triples
+  for (
+    auto itr = cdt.finite_faces_begin();
+    itr != cdt.finite_faces_end(); 
+    itr++)
+  {
+    faces.push_back( 
+      { v2i[itr->vertex(0)], v2i[itr->vertex(1)], v2i[itr->vertex(2)] });
+  }
+}
+
+template < typename Kernel, typename DerivedV, typename DerivedF>
+IGL_INLINE void igl::copyleft::cgal::projected_cdt(
+  const std::vector<CGAL::Object> & objects,
+  const CGAL::Plane_3<Kernel> & P,
+  Eigen::PlainObjectBase<DerivedV> & V,
+  Eigen::PlainObjectBase<DerivedF> & F)
+{
+  std::vector<CGAL::Point_3<Kernel> > vertices;
+  std::vector<std::vector<typename DerivedF::Scalar> > faces;
+  projected_cdt(objects,P,vertices,faces);
+  V.resize(vertices.size(),3);
+  for(int v = 0;v<vertices.size();v++)
+  {
+    for(int d = 0;d<3;d++)
+    {
+      assign_scalar(vertices[v][d], V(v,d));
+    }
+  }
+  list_to_matrix(faces,F);
+}
+
+#ifdef IGL_STATIC_LIBRARY
+// Explicit template specialization
+// generated by autoexplicit.sh
+template void igl::copyleft::cgal::projected_cdt<CGAL::Epick, long>(std::::vector<CGAL::Object, std::::allocator<CGAL::Object> > const&, CGAL::Plane_3<CGAL::Epick> const&, std::::vector<CGAL::Point_3<CGAL::Epick>, std::::allocator<CGAL::Point_3<CGAL::Epick> > >&, std::::vector<std::::vector<long, std::::allocator<long> >, std::::allocator<std::::vector<long, std::::allocator<long> > > >&);
+// generated by autoexplicit.sh
+template void igl::copyleft::cgal::projected_cdt<CGAL::Epeck, long>(std::::vector<CGAL::Object, std::::allocator<CGAL::Object> > const&, CGAL::Plane_3<CGAL::Epeck> const&, std::::vector<CGAL::Point_3<CGAL::Epeck>, std::::allocator<CGAL::Point_3<CGAL::Epeck> > >&, std::::vector<std::::vector<long, std::::allocator<long> >, std::::allocator<std::::vector<long, std::::allocator<long> > > >&);
+#endif

include/igl/copyleft/cgal/projected_cdt.h

@@ -0,0 +1,60 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+// 
+// Copyright (C) 2016 Alec Jacobson
+// 
+// This Source Code Form is subject to the terms of the Mozilla Public License 
+// v. 2.0. If a copy of the MPL was not distributed with this file, You can 
+// obtain one at http://mozilla.org/MPL/2.0/.
+#ifndef IGL_COPYLEFT_CGAL_PROJECTED_CDT_H
+#define IGL_COPYLEFT_CGAL_PROJECTED_CDT_H
+#include "../../igl_inline.h"
+#include <Eigen/Core>
+#include <CGAL/Plane_3.h>
+#include <CGAL/Point_3.h>
+#include <vector>
+namespace igl
+{
+  namespace copyleft
+  {
+    namespace cgal
+    {
+      // Given a list of objects (e.g., resulting from intersecting a triangle
+      // with many other triangles), construct a constrained Delaunay
+      // triangulation on a given plane (P), by inersting constraints for each
+      // object projected onto that plane.
+      //
+      // Inputs:
+      //   objects  list of objects. This should lie on the given plane (P),
+      //     otherwise they are added to the cdt _after_ their non-trivial
+      //     projection
+      //   P  plane upon which all objects lie and upon which the CDT is
+      //     conducted
+      // Outputs:
+      //   vertices  list of vertices of the CDT mesh _back on the 3D plane_
+      //   faces  list of list of triangle indices into vertices
+      //   
+      template <typename Kernel, typename Index>
+      IGL_INLINE void projected_cdt(
+        const std::vector<CGAL::Object> & objects,
+        const CGAL::Plane_3<Kernel> & P,
+        std::vector<CGAL::Point_3<Kernel> >& vertices,
+        std::vector<std::vector<Index> >& faces);
+      // Outputs:
+      //   V  #V by 3 list of vertices of the CDT mesh _back on the 3D plane_,
+      //     **cast** from the number type of Kernel to the number type of
+      //     DerivedV
+      //   F  #F by 3 list of triangle indices into V
+      template < typename Kernel, typename DerivedV, typename DerivedF>
+      IGL_INLINE void projected_cdt(
+        const std::vector<CGAL::Object> & objects,
+        const CGAL::Plane_3<Kernel> & P,
+        Eigen::PlainObjectBase<DerivedV> & V,
+        Eigen::PlainObjectBase<DerivedF> & F);
+    }
+  }
+}
+
+#ifndef IGL_STATIC_LIBRARY
+#  include "projected_cdt.cpp"
+#endif
+#endif

include/igl/copyleft/cgal/remesh_intersections.cpp

@@ -8,6 +8,7 @@
 //
 #include "remesh_intersections.h"
 #include "assign_scalar.h"
+#include "projected_cdt.h"
 #include "../../get_seconds.h"
 #include "../../unique.h"
 
@@ -162,75 +163,27 @@ IGL_INLINE void igl::copyleft::cgal::remesh_intersections(
       std::vector<Point_3>& vertices,
       std::vector<std::vector<Index> >& faces) -> void 
     {
+      std::vector<CGAL::Object> objects;
+
       CDT_plus_2 cdt;
       // insert each face into a common cdt
       for (const auto& fid : involved_faces)
       {
         const auto itr = offending.find(fid);
         const auto& triangle = T[fid];
-        cdt.insert_constraint(P.to_2d(triangle[0]), P.to_2d(triangle[1]));
-        cdt.insert_constraint(P.to_2d(triangle[1]), P.to_2d(triangle[2]));
-        cdt.insert_constraint(P.to_2d(triangle[2]), P.to_2d(triangle[0]));
+        objects.push_back(CGAL::make_object(triangle));
         if (itr == offending.end())
         {
           continue;
         }
         for (const auto& index_obj : itr->second) 
         {
-            //const auto& ofid = index_obj.first;
-            const auto& obj = index_obj.second;
-            if(const Segment_3 *iseg = CGAL::object_cast<Segment_3 >(&obj)) {
-                // Add segment constraint
-                cdt.insert_constraint(
-                        P.to_2d(iseg->vertex(0)),P.to_2d(iseg->vertex(1)));
-            }else if(const Point_3 *ipoint = CGAL::object_cast<Point_3 >(&obj)) {
-                // Add point
-                cdt.insert(P.to_2d(*ipoint));
-            } else if(const Triangle_3 *itri = CGAL::object_cast<Triangle_3 >(&obj)) {
-                // Add 3 segment constraints
-                cdt.insert_constraint(
-                        P.to_2d(itri->vertex(0)),P.to_2d(itri->vertex(1)));
-                cdt.insert_constraint(
-                        P.to_2d(itri->vertex(1)),P.to_2d(itri->vertex(2)));
-                cdt.insert_constraint(
-                        P.to_2d(itri->vertex(2)),P.to_2d(itri->vertex(0)));
-            } else if(const std::vector<Point_3 > *polyp = 
-                    CGAL::object_cast< std::vector<Point_3 > >(&obj)) {
-                //cerr<<REDRUM("Poly...")<<endl;
-                const std::vector<Point_3 > & poly = *polyp;
-                const Index m = poly.size();
-                assert(m>=2);
-                for(Index p = 0;p<m;p++)
-                {
-                    const Index np = (p+1)%m;
-                    cdt.insert_constraint(P.to_2d(poly[p]),P.to_2d(poly[np]));
-                }
-            }else {
-                throw std::runtime_error("Unknown intersection object!");
-            }
+          //const auto& ofid = index_obj.first;
+          const auto& obj = index_obj.second;
+          objects.push_back(obj);
         }
       }
-      // Read off vertices of the cdt, remembering index
-      std::map<typename CDT_plus_2::Vertex_handle,Index> v2i;
-      size_t count=0;
-      for (
-        auto itr = cdt.finite_vertices_begin();
-        itr != cdt.finite_vertices_end(); 
-        itr++) 
-      {
-        vertices.push_back(P.to_3d(itr->point()));
-        v2i[itr] = count;
-        count++;
-      }
-      // Read off faces and store index triples
-      for (
-        auto itr = cdt.finite_faces_begin();
-        itr != cdt.finite_faces_end(); 
-        itr++)
-      {
-        faces.push_back( 
-          { v2i[itr->vertex(0)], v2i[itr->vertex(1)], v2i[itr->vertex(2)] });
-      }
+      projected_cdt(objects,P,vertices,faces);
     };
 
     // Given p on triangle indexed by ori_f, add point to list of vertices return index of p.