- added tutorial example LSCM + helper functions

Former-commit-id: d8dd6bcf8b25c2b5628433c5e78f8e9ea50c4278

include/igl/areamatrix.cpp

@@ -0,0 +1,53 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+//
+// Copyright (C) 2014 Daniele Panozzo <daniele.panozzo@gmail.com>
+//
+// 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 "areamatrix.h"
+#include <vector>
+
+// Bug in unsupported/Eigen/SparseExtra needs iostream first
+#include <iostream>
+#include <unsupported/Eigen/SparseExtra>
+
+#include <igl/boundary_vertices_sorted.h>
+
+template <typename DerivedV, typename DerivedF, typename Scalar>
+IGL_INLINE void igl::areamatrix(
+  const Eigen::PlainObjectBase<DerivedV> & V,
+  const Eigen::PlainObjectBase<DerivedF> & F,
+  Eigen::SparseMatrix<Scalar>& A)
+{
+  using namespace Eigen;
+  using namespace std;
+
+	SparseMatrix<Scalar> aux (V.rows() * 2, V.rows() * 2);
+	SparseMatrix<Scalar> auxT(V.rows() * 2, V.rows() * 2);
+
+	vector<Triplet<Scalar> > auxTripletList;
+	vector<Triplet<Scalar> > auxTTripletList;
+
+  Eigen::VectorXi bnd;
+  igl::boundary_vertices_sorted(V,F,bnd);
+
+	for(int k = 0; k < bnd.size(); k++)
+  {
+		int i = bnd[k];
+		int j = k + 1 == bnd.size() ? bnd[0] : bnd[k+1];
+		auxTripletList.push_back(Triplet<Scalar>(i+V.rows(), j, -0.5));
+		auxTripletList.push_back(Triplet<Scalar>(i, j+V.rows(), 0.5));
+		auxTTripletList.push_back(Triplet<Scalar>(j, i+V.rows(), -0.5));
+		auxTTripletList.push_back(Triplet<Scalar>(j+V.rows(), i, 0.5));
+	}
+
+	aux.setFromTriplets(auxTripletList.begin(), auxTripletList.end());
+	auxT.setFromTriplets(auxTTripletList.begin(), auxTTripletList.end());
+	A = (aux + auxT)/2;
+}
+
+#ifndef IGL_HEADER_ONLY
+// Explicit template specialization
+// generated by autoexplicit.sh
+#endif

include/igl/areamatrix.h

@@ -0,0 +1,44 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+//
+// Copyright (C) 2014 Daniele Panozzo <daniele.panozzo@gmail.com>
+//
+// 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_AREAMATRIX_H
+#define IGL_AREAMATRIX_H
+#include "igl_inline.h"
+
+#include <Eigen/Dense>
+#include <Eigen/Sparse>
+
+namespace igl
+{
+  // Constructs the symmetric area matrix A, s.t.
+  // [V.col(0)' V.col(1)'] * A * [V.col(0); V.col(1)] is the area of the planar mesh (V,F).
+  // Note: (V,F) must be a genus-0 mesh, with a single boundary
+  //
+  // Templates:
+  //   DerivedV  derived type of eigen matrix for V (e.g. derived from
+  //     MatrixXd)
+  //   DerivedF  derived type of eigen matrix for F (e.g. derived from
+  //     MatrixXi)
+  //   Scalar  scalar type for eigen sparse matrix (e.g. double)
+  // Inputs:
+  //   V  #V by 2 list of mesh vertex positions
+  //   F  #F by 3 list of mesh faces (must be triangles)
+  // Outputs:
+  //   A  #Vx2 by #Vx2 area matrix
+  //
+  template <typename DerivedV, typename DerivedF, typename Scalar>
+  IGL_INLINE void areamatrix(
+    const Eigen::PlainObjectBase<DerivedV> & V,
+    const Eigen::PlainObjectBase<DerivedF> & F,
+    Eigen::SparseMatrix<Scalar>& A);
+}
+
+#ifdef IGL_HEADER_ONLY
+#  include "areamatrix.cpp"
+#endif
+
+#endif

include/igl/kronecker_product.cpp

@@ -0,0 +1,64 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+//
+// Copyright (C) 2014 Daniele Panozzo <daniele.panozzo@gmail.com>
+//
+// 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 "kronecker_product.h"
+
+// Bug in unsupported/Eigen/SparseExtra needs iostream first
+#include <iostream>
+#include <unsupported/Eigen/SparseExtra>
+
+template <typename Scalar>
+IGL_INLINE Eigen::SparseMatrix<Scalar> igl::kronecker_product(
+  const Eigen::SparseMatrix<Scalar> & A,
+  const Eigen::SparseMatrix<Scalar> & B)
+{
+  using namespace Eigen;
+  using namespace std;
+
+  // Convert B in triplets format
+  MatrixXd B_triplets(B.nonZeros(),3);
+  int count = 0;
+  for (int k=0; k<B.outerSize(); ++k)
+    for (SparseMatrix<double>::InnerIterator it(B,k); it; ++it)
+      B_triplets.row(count++) << it.row(), it.col(), it.value();
+
+  MatrixXd C_triplets(B_triplets.rows()*A.nonZeros(),3);
+  count = 0;
+  for (int k=0; k<A.outerSize(); ++k)
+    for (SparseMatrix<double>::InnerIterator it(A,k); it; ++it)
+    {
+      int i = it.row();
+      int j = it.col();
+      double v = it.value();
+
+      MatrixXd B_triplets_copy = B_triplets;
+      B_triplets_copy.col(0) = B_triplets_copy.col(0).array() + double(B.rows()*i);
+      B_triplets_copy.col(1) = B_triplets_copy.col(1).array() + double(B.cols()*j);
+      B_triplets_copy.col(2) = B_triplets_copy.col(2).array() * v;
+
+      C_triplets.block(count*B_triplets.rows(),0,
+                       B_triplets.rows(), B_triplets.cols()) = B_triplets_copy;
+
+      count++;
+    }
+
+  typedef Eigen::Triplet<double> T;
+  std::vector<T> triplets;
+  triplets.reserve(C_triplets.rows());
+  
+  for(unsigned i=0; i<C_triplets.rows(); ++i)
+    triplets.push_back(T(C_triplets(i,0),C_triplets(i,1),C_triplets(i,2)));
+  SparseMatrix<Scalar> C(A.rows()*B.rows(),A.cols()*B.cols());
+  C.setFromTriplets(triplets.begin(),triplets.end());
+  
+  return C;
+}
+
+#ifndef IGL_HEADER_ONLY
+// Explicit template specialization
+// generated by autoexplicit.sh
+#endif

include/igl/kronecker_product.h

@@ -0,0 +1,35 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+//
+// Copyright (C) 2014 Daniele Panozzo <daniele.panozzo@gmail.com>
+//
+// 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_KRONECKERPRODUCT_H
+#define IGL_KRONECKERPRODUCT_H
+#include "igl_inline.h"
+
+#include <Eigen/Dense>
+#include <Eigen/Sparse>
+
+namespace igl
+{
+  // Computes the Kronecker product between sparse matrices A and B.
+  //
+  // Inputs:
+  //   A  #M by #N sparse matrix
+  //   B  #P by #Q sparse matrix
+  // Outputs:
+  //      #M*#P by #N*#Q sparse matrix
+  //
+  template <typename Scalar>
+  IGL_INLINE Eigen::SparseMatrix<Scalar> kronecker_product(
+    const Eigen::SparseMatrix<Scalar> & A,
+    const Eigen::SparseMatrix<Scalar> & B);
+}
+
+#ifdef IGL_HEADER_ONLY
+#  include "kronecker_product.cpp"
+#endif
+
+#endif

include/igl/lscm.cpp

@@ -0,0 +1,83 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+//
+// Copyright (C) 2014 Daniele Panozzo <daniele.panozzo@gmail.com>
+//
+// 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 "lscm.h"
+
+// Bug in unsupported/Eigen/SparseExtra needs iostream first
+#include <iostream>
+#include <unsupported/Eigen/SparseExtra>
+
+#include <igl/areamatrix.h>
+#include <igl/cotmatrix.h>
+#include <igl/kronecker_product.h>
+#include <igl/boundary_vertices_sorted.h>
+#include <igl/min_quad_with_fixed.h>
+
+IGL_INLINE Eigen::MatrixXd igl::lscm(
+  const Eigen::MatrixXd& V,
+  const Eigen::MatrixXi& F,
+  const Eigen::VectorXi& b,
+  const Eigen::MatrixXd& bc)
+{
+  using namespace Eigen;
+  
+  // Assemble the area matrix (note that A is #Vx2 by #Vx2)
+  SparseMatrix<double> A;
+  igl::areamatrix(V,F,A);
+
+  // Assemble the cotan laplacian matrix
+  SparseMatrix<double> L;
+  igl::cotmatrix(V,F,L);
+
+  SparseMatrix<double> I2(2,2);
+  I2.insert(0,0) = 1;
+  I2.insert(1,1) = 1;
+
+  SparseMatrix<double> L_flat = igl::kronecker_product(I2,L);
+
+  VectorXi b_flat = b;
+  MatrixXd bc_flat = bc;
+  
+  if (b.size() < 2)
+  {
+    // if no boundary conditions are provided, fix two boundary points
+    Eigen::VectorXi bnd;
+    igl::boundary_vertices_sorted(V,F,bnd);
+
+    int v1 = bnd(0);
+    int v2 = bnd(round(bnd.size()/2));
+
+    b_flat.resize(4);
+    b_flat << v1, v1+1, v2, v2+1;
+    bc_flat.resize(4,1);
+    bc_flat << 0, 1, 0, 0;
+  }
+
+  // Minimize the LSCM energy
+  SparseMatrix<double> Q = -0.5*L_flat + A;
+  const VectorXd B_flat = VectorXd::Zero(V.rows()*2);
+  igl::min_quad_with_fixed_data<double> data;
+  igl::min_quad_with_fixed_precompute(Q,b_flat,SparseMatrix<double>(),true,data);
+
+  MatrixXd W_flat;
+  if(!min_quad_with_fixed_solve(data,B_flat,bc_flat,VectorXd(),W_flat))
+    assert(false);
+
+  MatrixXd V_uv(V.rows(),2);
+  for (unsigned i=0;i<V.rows();++i)
+  {
+    V_uv(i,0) = W_flat(i);
+    V_uv(i,1) = W_flat(i+V.rows());
+  }
+
+  return V_uv;
+}
+
+#ifndef IGL_HEADER_ONLY
+// Explicit template specialization
+// generated by autoexplicit.sh
+#endif

include/igl/lscm.h

@@ -0,0 +1,44 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+//
+// Copyright (C) 2014 Daniele Panozzo <daniele.panozzo@gmail.com>
+//
+// 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_LSCM_H
+#define IGL_LSCM_H
+#include "igl_inline.h"
+
+#include <Eigen/Dense>
+#include <Eigen/Sparse>
+
+namespace igl
+{
+  // Compute a Least-squares conformal map parametrization following the algorithm
+  // presented in: Spectral Conformal Parameterization,
+  //               Patrick Mullen, Yiying Tong, Pierre Alliez and Mathieu Desbrun
+  // Note: (V,F) must be a genus-0 mesh, with a single boundary
+  //
+  // Inputs:
+  //   V  #V by 3 list of mesh vertex positions
+  //   F  #F by 3 list of mesh faces (must be triangles)
+  //   b  #b boundary indices into V
+  //   bc #b by #W list of boundary values
+// Outputs:
+  //   UV #V by 2 list of 2D mesh vertex positions in UV space
+  //
+  // Note: if b and bc are empty, lscm automatically removes the null space
+  //       by fixing two farthest points on the boundary
+
+  IGL_INLINE Eigen::MatrixXd lscm(
+    const Eigen::MatrixXd& V,
+    const Eigen::MatrixXi& F,
+    const Eigen::VectorXi& b,
+    const Eigen::MatrixXd& bc);
+}
+
+#ifdef IGL_HEADER_ONLY
+#  include "lscm.cpp"
+#endif
+
+#endif

tutorial/502_LSCMParam/CMakeLists.txt

@@ -0,0 +1,11 @@
+cmake_minimum_required(VERSION 2.6)
+project(502_LSCMParam)
+
+include("../CMakeLists.shared")
+
+set(SOURCES
+${PROJECT_SOURCE_DIR}/main.cpp
+)
+
+add_executable(${CMAKE_PROJECT_NAME} ${SOURCES} ${SHARED_SOURCES})
+target_link_libraries(${CMAKE_PROJECT_NAME} ${SHARED_LIBRARIES})

tutorial/502_LSCMParam/main.cpp

@@ -0,0 +1,53 @@
+#include <igl/readOFF.h>
+#include <igl/viewer/Viewer.h>
+
+#include <igl/lscm.h>
+
+Eigen::MatrixXd V;
+Eigen::MatrixXi F;
+Eigen::MatrixXd V_uv;
+
+bool key_down(igl::Viewer& viewer, unsigned char key, int modifier)
+{
+
+  if (key == '1')
+    // Plot the 3D mesh
+    viewer.set_mesh(V,F);
+  else if (key == '2')
+    // Plot the mesh in 2D using the UV coordinates as vertex coordinates
+    viewer.set_mesh(V_uv,F);
+
+  viewer.compute_normals();
+
+  return false;
+}
+
+int main(int argc, char *argv[])
+{
+  using namespace Eigen;
+  using namespace std;
+
+  // Load a mesh in OFF format
+  igl::readOFF("../shared/camelhead.off", V, F);
+
+  // LSCM parametrization
+  V_uv = igl::lscm(V,F,Eigen::VectorXi(),Eigen::MatrixXd());
+
+  // Scale the uv
+  V_uv *= 5;
+
+  // Plot the mesh
+  igl::Viewer viewer;
+  viewer.set_mesh(V, F);
+  viewer.set_uv(V_uv);
+  viewer.callback_key_down = &key_down;
+
+  // Disable wireframe
+  viewer.options.show_lines = false;
+
+  // Draw checkerboard texture
+  viewer.options.show_texture = true;
+
+  // Launch the viewer
+  viewer.launch();
+}