【libigl】Libigl_Tutorials -第四章(共9节）

第四章 形状变形

现代基于网格的形状变形方法满足手柄（网格上选定的顶点或区域）处的用户变形约束，并将这些手柄变形平滑地传播到形状的其余部分，而不会删除或扭曲细节。Libigl 提供了各种最先进的变形技术的实现，范围从基于二次网格的能量最小化器到蒙皮方法，再到非线性弹性启发的技术。

一、 双调和变形 Biharmonic Deformation

2000 年至 2010 年间的研究产生了一系列技术，将基于手柄的形状变形问题转化为二次能量最小化问题或等效于线性偏微分方程的解。

这些技术有很多风格，但典型的子集是那些考虑双拉普拉斯方程解的技术，即双调和函数 ¹¹ 。此四阶偏微分方程在边界条件中提供了足够的灵活性，以确保句柄约束（在细化的限制中） C1�1 连续性 ¹⁵ 。

 双调和曲面 ¶

让我们首先通过考虑双调和曲面来开始讨论双调和变形。我们将随意将双调和曲面定义为其位置函数相对于某些初始参数化是双调和的曲面：

Δ2x′=0Δ2�′=0

并受到一些句柄约束，概念化为“边界条件”：

x′b=xbc.��′=���.

其中 x′�′ 是表面上点的未知 3D 位置。所以我们要求每个空间坐标函数的双拉普拉斯为零。

在 libigl 中，可以通过 igl::harmonic 并设置 k=2�=2 （双谐波）来解决双调和问题：

// U_bc contains deformation of boundary vertices b
igl::harmonic(V,F,b,U_bc,2,U);

这会产生一个插入手柄约束的平滑表面，但表面上的所有原始细节都将被平滑掉。最明显的是，如果原始表面还不是双调和的，那么给予所有手柄恒等变形（将它们保持在静止位置）将不会再现原始表面。相反，结果将是对这些手柄位置进行插值的双调和曲面。

因此，我们可以得出结论，这不是一种直观的形状变形技术。

双调和变形场 ¶

现在我们知道变形技术的一个有用属性是“静止姿势再现”：不对手柄应用变形也不应对形状应用变形。

为了通过构造来保证这一点，我们可以使用变形场（即位移） d� 而不是直接使用位置 x� 。然后变形位置可以恢复为

x′=x+d.�′=�+�.

对手柄约束的变形场进行插值的平滑变形场 d� 将施加平滑变形形状 x′�′ 。当然，我们考虑双调和变形场：

Δ2d=0Δ2�=0

受到相同的手柄约束，但根据边界（手柄）处的隐含变形场重写：

db=xbc−xb.��=���−��.

我们再次可以将 igl::harmonic 与 k=2�=2 一起使用，但这一次求解变形场，然后恢复变形位置：

// U_bc contains deformation of boundary vertices b
D_bc = U_bc - igl::slice(V,b,1);
igl::harmonic(V,F,b,D_bc,2,D);
U = V+D;

双调和变形示例将雕像的头部变形为双调和表面（顶部）并使用双调和位移（底部）。

与“微分坐标”和拉普拉斯曲面编辑的关系 ¶

双调和函数（无论是位置还是位移）是双拉普拉斯方程的解，也是“拉普拉斯能量”的最小化器。例如，对于位移 d� ，能量读取

∫S∥Δd∥2dA,∫�‖Δ�‖2��,

我们定义 ΔdΔ� 来简单地应用拉普拉斯坐标。

通过 Laplace(-Beltrami) 算子的线性，我们可以根据原始位置 x� 和未知位置 x′=x−d�′=�−� 重新表达该能量：

∫S∥Δ(x′−x)∥2dA=∫S∥Δx′−Δx)∥2dA.∫�‖Δ(�′−�)‖2��=∫�‖Δ�′−Δ�)‖2��.

在 Sorkine 等人的早期工作中，数量 Δx′Δ�′ 和 ΔxΔ� 被称为“微分坐标” ²¹ 。因此，它们的变形（没有线性化旋转）相当于双调和变形场。

#include <igl/colon.h> // 包含libigl库中的colon函数，用于生成等差数列
#include <igl/harmonic.h> // 包含libigl库中的harmonic函数，用于生成调和映射
#include <igl/readOBJ.h> // 包含libigl库中的readOBJ函数，用于读取.obj文件
#include <igl/readDMAT.h> // 包含libigl库中的readDMAT函数，用于读取.dmat文件
#include <igl/opengl/glfw/Viewer.h> // 包含libigl库配合glfw库使用的Viewer组件
#include <algorithm> // 包含STL算法库
#include <iostream> // 包含输入输出流库


double bc_frac = 1.0; // 控制边界条件变化的因子
double bc_dir = -0.03; // 控制bc_frac增减的方向和步长
bool deformation_field = false; // 控制是否显示变形场
Eigen::MatrixXd V,U,V_bc,U_bc; // 定义顶点矩阵V和U，边界条件下的V_bc和U_bc
Eigen::VectorXd Z; // 定义向量Z
Eigen::MatrixXi F; // 定义面矩阵F
Eigen::VectorXi b; // 定义控制边界的索引向量b


bool pre_draw(igl::opengl::glfw::Viewer & viewer) // 定义预绘制函数
{
  using namespace Eigen; // 使用Eigen命名空间
  // Determine boundary conditions
  if(viewer.core().is_animating) // 如果观察者处于动画状态
  {
    bc_frac += bc_dir; // 更新bc_frac
    bc_dir *= (bc_frac>=1.0 || bc_frac<=0.0?-1.0:1.0); // 如果bc_frac达到边界，则改变增减方向
  }


  const MatrixXd U_bc_anim = V_bc+bc_frac*(U_bc-V_bc); // 计算边界条件动画的变形位置
  if(deformation_field) // 如果显示变形场
  {
    MatrixXd D; // 定义变形场矩阵D
    MatrixXd D_bc = U_bc_anim - V_bc; // 计算边界条件下的变形场
    igl::harmonic(V,F,b,D_bc,2,D); // 通过调和映射计算D
    U = V+D; // 更新顶点位置
  }else // 如果不显示变形场
  {
    igl::harmonic(V,F,b,U_bc_anim,2.,U); // 通过调和映射直接计算顶点位置
  }
  viewer.data().set_vertices(U); // 更新观察者中的顶点位置
  viewer.data().compute_normals(); // 重新计算法线
  return false; // 返回false，表示不中断当前的绘制流程
}


bool key_down(igl::opengl::glfw::Viewer &viewer, unsigned char key, int mods) // 定义按键响应函数
{
  switch(key) // 根据按键进行不同的操作
  {
    case ' ': // 如果按下空格键
      viewer.core().is_animating = !viewer.core().is_animating; // 切换动画状态
      return true;
    case 'D': // 如果按下D键
    case 'd': // 如果按下d键
      deformation_field = !deformation_field; // 切换是否显示变形场
      return true;
  }
  return false; // 返回false，表示不拦截剩余的按键事件
}


int main(int argc, char *argv[]) // 程序主函数
{
  using namespace Eigen; // 使用Eigen命名空间
  using namespace std; // 使用std命名空间
  igl::readOBJ(TUTORIAL_SHARED_PATH "/decimated-max.obj",V,F); // 读取obj文件，获取顶点矩阵V和面矩阵F
  U=V; // 将U初始化为V，即初始变形为0
  // S(i) = j: j<0 (vertex i not in handle), j >= 0 (vertex i in handle j)
  VectorXi S; // 定义选择向量S，用于标识某顶点是否是控制柄的一部分
  igl::readDMAT(TUTORIAL_SHARED_PATH "/decimated-max-selection.dmat",S); // 读取.dmat文件获取选择向量S
  igl::colon<int>(0,V.rows()-1,b); // 使用colon函数填充b，从0到V的行数-1
  b.conservativeResize(stable_partition( b.data(), b.data()+b.size(),
   [&S](int i)->bool{return S(i)>=0;})-b.data()); // 根据选择向量S过滤控制边界索引b


  // Boundary conditions directly on deformed positions
  U_bc.resize(b.size(),V.cols()); // 调整U_bc的尺寸，用于存储边界条件的变形位置
  V_bc.resize(b.size(),V.cols()); // 调整V_bc的尺寸，用于存储边界条件下的初始位置
  for(int bi = 0;bi<b.size();bi++) // 遍历边界索引b
  {
    V_bc.row(bi) = V.row(b(bi)); // 设置V_bc为V在边界索引b处的行
    switch(S(b(bi))) // 根据选择向量S的值进行不同的操作
    {
      case 0:
        // Don't move handle 0
        U_bc.row(bi) = V.row(b(bi)); // 如果是控制柄0，U_bc不变化
        break;
      case 1:
        // move handle 1 down
        U_bc.row(bi) = V.row(b(bi)) + RowVector3d(0,-50,0); // 如果是控制柄1，U_bc向下移动
        break;
      case 2:
      default:
        // move other handles forward
        U_bc.row(bi) = V.row(b(bi)) + RowVector3d(0,0,-25); // 如果是其它控制柄，U_bc向前移动
        break;
    }
  }


  // Pseudo-color based on selection
  MatrixXd C(F.rows(),3); // 定义颜色矩阵C
  RowVector3d purple(80.0/255.0,64.0/255.0,255.0/255.0); // 紫色
  RowVector3d gold(255.0/255.0,228.0/255.0,58.0/255.0); // 金色
  for(int f = 0;f<F.rows();f++) // 遍历所有的面
  {
    if( S(F(f,0))>=0 && S(F(f,1))>=0 && S(F(f,2))>=0) // 如果面的所有顶点都是控制柄的一部分
    {
      C.row(f) = purple; // 设置该面的颜色为紫色
    }else
    {
      C.row(f) = gold; // 否则设置为金色
    }
  }


  // Plot the mesh with pseudocolors
  igl::opengl::glfw::Viewer viewer; // 创建一个Viewer
  viewer.data().set_mesh(U, F); // 设置Viewer的网格数据
  viewer.data().show_lines = false; // 不显示网格线
  viewer.data().set_colors(C); // 设置面的颜色
  viewer.core().trackball_angle = Eigen::Quaternionf(sqrt(2.0),0,sqrt(2.0),0); // 设置观察者的视角
  viewer.core().trackball_angle.normalize(); // 正规化四元数来表示旋转
  viewer.callback_pre_draw = &pre_draw; // 设置预绘制回调函数
  viewer.callback_key_down = &key_down; // 设置按键响应回调函数
  //viewer.core().is_animating = true; // 初始化时设置为动画状态，该行代码被注释
  viewer.core().animation_max_fps = 30.; // 设置动画的最大帧率为30fps
  cout<<
    "Press [space] to toggle deformation."<<endl<<
    "Press 'd' to toggle between biharmonic surface or displacements."<<endl; // 输出控制指令
  viewer.launch(); // 启动Viewer
}

此代码是一个使用libigl库进行几何处理和视图渲染的C++程序。程序的主要功能是读取OBJ格式的网格模型，计算网格的调和映射并实现交互式动画效果。用户可以通过按下空格键来开始或停止动画，通过按下'd'键来切换是否显示变形场。程序还能够根据用户的选择对网格的某些部分进行伪彩色渲染。总的来说，这个程序演示了如何使用libigl库处理网格和实现基于变形的视觉效果。

二、多谐波变形 Polyharmonic Deformation

我们可以通过考虑拉普拉斯算子的不同幂来概括双调和变形，从而产生一系列以下形式的偏微分方程：

Δkd=0.Δ��=0.

与 k∈1,2,3,…�∈1,2,3,… 。 k� 的选择决定了手柄处的连续性级别。特别是， k=1�=1 暗示边界处的 C0�0 ， k=2�=2 暗示 C1�1 ， k=3�=3 暗示 C2�2 意味着 Ck−1��−1 。

int k = 2;// or 1,3,4,...
igl::harmonic(V,F,b,bc,k,Z);

PolyharmonicDeformation 示例将平坦域（左）变形为凸块，作为各种 k� 谐波偏微分方程的解。

#include <igl/colon.h> // 包含libigl库中的colon函数，用于生成向量中等自增的序列
#include <igl/harmonic.h> // 包含libigl库中的harmonic函数，用于调和映射计算
#include <igl/readOBJ.h> // 包含libigl库中的readOBJ函数，用于读取OBJ格式的文件
#include <igl/opengl/glfw/Viewer.h> // 包含libigl库和glfw库相关的Viewer类
#include <algorithm> // 包含算法相关的标准库函数，比如stable_partition
#include <iostream> // 包含cin、cout等标准输入输出流的函数和对象


double z_max = 1.0; // 控制网格在z轴方向上的变形比例
double z_dir = -0.03; // 控制z_max增减的方向和步长
int k = 2; // 调和映射的次数，可以通过按键修改
bool resolve = true; // 标志是否需要重新解调和映射问题
Eigen::MatrixXd V,U; // 定义原始顶点矩阵V，变形后的顶点矩阵U
Eigen::VectorXd Z; // 用于存放调和映射计算结果的Z轴方向的值
Eigen::MatrixXi F; // 面矩阵F
Eigen::VectorXi b; // 边界顶点索引向量b
Eigen::VectorXd bc; // 边界条件向量bc，决定每个边界顶点的变形方式


bool pre_draw(igl::opengl::glfw::Viewer & viewer) // 定义预绘制时的回调函数
{
  using namespace Eigen; // 使用Eigen命名空间
  if(resolve) // 如果标志位resolve为true，则说明需要重新计算调和映射
  {
    // 解调和映射问题，输入为网格的顶点和面，边界顶点索引，边界条件，映射次数
    igl::harmonic(V,F,b,bc,k,Z);
    resolve = false; // 重置resolve标志位，避免下一次绘制时重复计算
  }
  U.col(2) = z_max*Z; // 将Z轴方向的变形应用到网格上
  viewer.data().set_vertices(U); // 更新viewer的顶点位置
  viewer.data().compute_normals(); // 重新计算法线
  if(viewer.core().is_animating) // 如果viewer正在进行动画模式
  {
    z_max += z_dir; // 更新z_max
    // 当z_max达到[0,1]的边界时，改变z_dir的方向
    z_dir *= (z_max>=1.0 || z_max<=0.0?-1.0:1.0);
  }
  return false; // 返回false，表示预绘制回调结束
}


bool key_down(igl::opengl::glfw::Viewer &viewer, unsigned char key, int mods) // 定义按键被按下时的回调函数
{
  switch(key) // 根据不同的按键进行不同的操作
  {
    case ' ': // 空格键
      // 切换viewer的动画模式状态
      viewer.core().is_animating = !viewer.core().is_animating;
      break;
    case '.': // '.'键，增加调和映射的次数k
      k++;
      k = (k>4?4:k); // 限制k的最大值为4
      resolve = true; // 需要重新计算调和映射
      break;
    case ',': // ','键，减少调和映射的次数k
      k--;
      k = (k<1?1:k); // 限制k的最小值为1
      resolve = true; // 需要重新计算调和映射
      break;
  }
  return true; // 返回true，表示按键事件处理完毕
}


int main(int argc, char *argv[]) // 程序的主函数
{
  using namespace Eigen; // 使用Eigen命名空间
  using namespace std; // 使用std命名空间
  igl::readOBJ(TUTORIAL_SHARED_PATH "/bump-domain.obj",V,F); // 读取OBJ格式文件
  U=V; // 将顶点矩阵U初始化为V
  // Find boundary vertices outside annulus
  typedef Matrix<bool,Dynamic,1> VectorXb; // 定义布尔型向量
  VectorXb is_outer = (V.rowwise().norm().array()-1.0)>-1e-15; // 计算外环顶点
  VectorXb is_inner = (V.rowwise().norm().array()-0.15)<1e-15; // 计算内环顶点
  VectorXb in_b = is_outer.array() || is_inner.array(); // 合并内外环顶点作为边界顶点
  igl::colon<int>(0,V.rows()-1,b); // 使用colon函数初始化边界顶点索引
  b.conservativeResize(stable_partition( b.data(), b.data()+b.size(),
   [&in_b](int i)->bool{return in_b(i);})-b.data()); // 根据in_b排列
  bc.resize(b.size(),1); // 根据b的大小调整边界条件的大小
  for(int bi = 0;bi<b.size();bi++) // 遍历所有边界顶点
  {
    // 根据顶点是外环还是内环给出边界条件
    bc(bi) = (is_outer(b(bi))?0.0:1.0);
  }




  // Pseudo-color based on selection
  MatrixXd C(F.rows(),3); // 定义面的颜色矩阵C
  RowVector3d purple(80.0/255.0,64.0/255.0,255.0/255.0); // 紫色向量
  RowVector3d gold(255.0/255.0,228.0/255.0,58.0/255.0); // 金色向量
  for(int f = 0;f<F.rows();f++) // 遍历所有的面
  {
    // 根据面的顶点是否在边界中给出颜色
    if( in_b(F(f,0)) && in_b(F(f,1)) && in_b(F(f,2)))
    {
      C.row(f) = purple; // 边界顶点的面设为紫色
    }else
    {
      C.row(f) = gold; // 非边界顶点的面设为金色
    }
  }


  // Plot the mesh with pseudocolors
  igl::opengl::glfw::Viewer viewer; // 创建Viewer对象
  viewer.data().set_mesh(U, F); // 设置网格数据
  viewer.data().show_lines = false; // 不显示网格线
  viewer.data().set_colors(C); // 设置面的颜色
  viewer.core().trackball_angle = Eigen::Quaternionf(0.81,-0.58,-0.03,-0.03); // 设置初始的视角
  viewer.core().trackball_angle.normalize(); // 归一化以确保它是个有效的旋转
  viewer.callback_pre_draw = &pre_draw; // 设置预绘制时的回调函数
  viewer.callback_key_down = &key_down; // 设置按键响应回调函数
  viewer.core().is_animating = true; // 启动时自动开始动画
  viewer.core().animation_max_fps = 30.; // 最大帧率限制为30FPS
  cout<<
    "Press [space] to toggle animation."<<endl<<
    "Press '.' to increase k."<<endl<<
    "Press ',' to decrease k."<<endl; // 输出提示信息
  viewer.launch(); // 启动viewer
}

这段代码是一个调和映射演示程序，利用了libigl库来进行表面曲线编辑。这个程序读取一个OBJ格式的网格，寻找边界顶点，并基于这些边界顶点执行调和映射。用户可以通过按键修改调和映射的参数，以及启动和停止动画。这个程序具体展示了如何实现基于调和映射的曲面变形。

三、有界双调和权重Bounded Biharmonic Weights