CatmullClark.cpp

#include <CGAL/Simple_cartesian.h>

#include <CGAL/Surface_mesh.h>
#include <CGAL/boost/graph/graph_traits_Surface_mesh.h>
#include <CGAL/boost/graph/copy_face_graph.h>

#include <CGAL/subdivision_method_3.h>

#include <iostream>
#include <fstream>
#include <vector>

typedef CGAL::Simple_cartesian<double>         Kernel;
typedef CGAL::Surface_mesh<Kernel::Point_3>    PolygonMesh;
typedef boost::graph_traits<PolygonMesh>::halfedge_descriptor halfedge_descriptor;
typedef boost::graph_traits<PolygonMesh>::face_descriptor face_descriptor;


int main(int argc, char** argv) {

  PolygonMesh pmesh, orig;
  std::ifstream in(argv[1]);
  if(in.fail()) {
    std::cerr << "Could not open input file " << argv[1] << std::endl;
    return 1;
  }
  in >> orig;

  // Make a copy of the mesh.
  // The indices of vertices/halfedges/faces are identical for the two
  CGAL::copy_face_graph(orig, pmesh);

  // Add a face property in the mesh that will be subdivided
  // Each face in the subdivided mesh will store the face in the original mesh 
  PolygonMesh::Property_map<face_descriptor,face_descriptor> f2f;
  bool found;
  boost::tie(f2f, found) = pmesh.property_map<face_descriptor,face_descriptor>("f:2f");
 
  // We store the non-border halfedges before we subdivide
  std::vector<halfedge_descriptor> he;
  for(halfedge_descriptor hd : halfedges(pmesh)){
    if(! is_border(hd, pmesh)){
      he.push_back(hd);
    }
  }
  // Now we subdivide 
  CGAL::Subdivision_method_3::CatmullClark_subdivision(pmesh);


  // and we can exploit that
  // (1) we know that the old halfedges are those with the target vertex on the same halfedge before subdivison
  // (2) that the indices of halfedges in orig have a correspondence in pmesh
  for(halfedge_descriptor hd : he){
    put(f2f, face(hd, pmesh), face(hd,orig));
  }

  return 0;
}