diff --git a/.gitattributes b/.gitattributes index dbadd85852f..57f9f214914 100644 --- a/.gitattributes +++ b/.gitattributes @@ -3776,6 +3776,8 @@ Surface_mesher/examples/Surface_mesher/skip_vcproj_auto_generation -text Surface_mesher/find_debug_macros -text Surface_mesher/find_non_debug_macros -text Surface_mesher/test/Surface_mesher/combined_spheres.cin -text +Surface_modeling/examples/Surface_modeling/Deform_mesh.cpp -text +Surface_modeling/examples/Surface_modeling/Deform_mesh.h -text Surface_modeling/examples/Surface_modeling/example.cpp -text Surface_modeling/examples/Surface_modeling/k-ring_BGL.cpp -text Surface_reconstruction_points_3/demo/Surface_reconstruction_points_3/GlSplat/cmake/FindGLEW.cmake -text diff --git a/Surface_modeling/examples/Surface_modeling/Deform_mesh.cpp b/Surface_modeling/examples/Surface_modeling/Deform_mesh.cpp new file mode 100644 index 00000000000..d85a86848c0 --- /dev/null +++ b/Surface_modeling/examples/Surface_modeling/Deform_mesh.cpp @@ -0,0 +1,134 @@ +#include "Deform_mesh.h" + + +Deform_mesh::Deform_mesh(Polyhedron &P) +{ + polyhedron = P; +} + +Deform_mesh::~Deform_mesh(void) +{ +} + +void Deform_mesh::region_of_interest(Vertex_iterator begin, Vertex_iterator end, size_t k) +{ + roi.clear(); + for (Vertex_iterator vit = begin; vit != end; vit ++) + { + Vertex_handle handle = vit; + roi.push_back(handle); + } + roi.push_back(end); + + int idx_lv = 0; // pointing the neighboring vertices on current level + int idx_lv_end; + + for (size_t lv = 0; lv < k; lv++) + { + idx_lv_end = roi.size(); + for ( ;idx_lv < idx_lv_end; idx_lv++ ) + { + Vertex_handle vh = roi[idx_lv]; + HV_circulator wc = vh->vertex_begin(), done(wc); + do { + Vertex_handle wh = wc->opposite()->vertex(); + vector ::iterator result = find(roi.begin(), roi.end(), wh); + if (result == roi.end()) + { + roi.push_back(wh); + } + ++wc; + }while(wc != done); + } + } +} + +void Deform_mesh::handles(Vertex_iterator begin, Vertex_iterator end) +{ + hdl.clear(); + for (Vertex_iterator vit = begin; vit != end; vit ++) + { + Vertex_handle handle = vit; + hdl.push_back(handle); + } + hdl.push_back(end); +} + +void Deform_mesh::preprocess() +{ + CGAL_TRACE("Calls preprocess()\n"); + + double time_init = clock(); + + double duration_assembly = 0.0; + double duration_prefactor = 0.0; + + // get #variables + unsigned int nb_variables = polyhedron.size_of_vertices(); + + CGAL_TRACE(" Creates matrix...\n"); + // Assemble linear system A*X=B + Taucs_solver_traits::Matrix A(nb_variables); // matrix is symmetric definite positive + Taucs_solver_traits::Vector X(nb_variables), B(nb_variables); + + assemble_laplacian(A, "uni"); + + duration_assembly = (clock() - time_init)/CLOCKS_PER_SEC; + CGAL_TRACE(" Creates matrix: done (%.2lf s)\n", duration_assembly); + CGAL_TRACE(" Pre-factorizing linear system...\n"); + + // Pre-factorizing the linear system A*X=B + time_init = clock(); + double D; + if(!solver.pre_factor(A, D)) + return; + duration_prefactor = (clock() - time_init)/CLOCKS_PER_SEC; + + CGAL_TRACE(" Pre-factorizing linear system: done (%.2lf s)\n", duration_prefactor); + +} + +void Deform_mesh::assemble_laplacian(Taucs_solver_traits::Matrix& A, string type) +{ + map idx; + int index = 0; + for (Vertex_iterator vit = polyhedron.vertices_begin(); vit != polyhedron.vertices_end(); vit++) + { + idx[vit] = index; + index++; + } + + for (Vertex_iterator vit = polyhedron.vertices_begin(); vit != polyhedron.vertices_end(); vit++) + { + Vertex_handle vi = vit; + double diagonal = 0; + int idx_i = idx[vi]; + HV_circulator wc = vi->vertex_begin(), done(wc); + do { + Vertex_handle vj = wc->opposite()->vertex(); + double wij = 1; + if (type == "cot") // cotangent Laplacian weights + { + ; + } + int idx_j = idx[vj]; + A.set_coef(idx_i, idx_j, -wij, true); // off-diagonal coefficient + diagonal += wij; + ++wc; + }while(wc != done); + // diagonal coefficient + A.set_coef(idx_i, idx_i, diagonal, true); + } + + // handle constraints + for (int i = 0; i < hdl.size(); i++) + { + int idx_i = idx[hdl[i]]; + A.set_coef(idx_i, idx_i, 1.0); + } +} + +void Deform_mesh::operator ()(Point p, Vector v) +{ + p = p - v; +} \ No newline at end of file diff --git a/Surface_modeling/examples/Surface_modeling/Deform_mesh.h b/Surface_modeling/examples/Surface_modeling/Deform_mesh.h new file mode 100644 index 00000000000..18adfb3a703 --- /dev/null +++ b/Surface_modeling/examples/Surface_modeling/Deform_mesh.h @@ -0,0 +1,57 @@ +#pragma once + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +typedef CGAL::Simple_cartesian Kernel; +typedef Kernel::Vector_3 Vector; +typedef Kernel::Point_3 Point; +typedef CGAL::Polyhedron_3 Polyhedron; + +typedef Polyhedron::Vertex_const_handle Vertex_handle; +typedef Polyhedron::Vertex_iterator Vertex_iterator; +typedef Polyhedron::Halfedge_around_vertex_const_circulator HV_circulator; + +using namespace std; +using namespace CGAL; + +class Deform_mesh +{ +private: + Polyhedron polyhedron; // target mesh + vector roi; + vector hdl; + vector dsplc; // displacement of handles + + Taucs_solver_traits solver; + +public: + // The constructor gets the Polyhedron that we will model + Deform_mesh(Polyhedron &P); + + // Release ressources + ~Deform_mesh(void); + + // The region of interest and the handles are a set of vertices + void region_of_interest(Vertex_iterator begin, Vertex_iterator end, size_t k); + + void handles(Vertex_iterator begin, Vertex_iterator end); + + // Before we can model we have to do some precomputation + void preprocess(); + + // Assemble Laplacian matrix A of linear system A*X=B + void assemble_laplacian(Taucs_solver_traits::Matrix& A, string type); + + // The operator will be called in a real time loop from the GUI. + void operator()(Point p, Vector v); + +};