【笔记】Libigl从入门到放弃

主要内容Libigl基本知识离散几何参量与算子Libigl基本知识Libigl设计原则没有复杂的数据类型。只使用矩阵和向量；尽可能减少对外部库的依赖只有头文件（也可以做一个静态编译库）函数封装性下载Libiglgit clone https://github.com/libigl/libigl.gitLigigl核心功能只依赖C++标准库和Eigen但是编译起来还是挺费劲的，哈！Mesh表示Lib

界明城

4195人浏览 · 2022-04-30 12:10:11

界明城 · 2022-04-30 12:10:11 发布

主要内容

Libigl基本知识
离散几何参量与算子

Libigl基本知识

Libigl设计原则

没有复杂的数据类型。只使用矩阵和向量；
尽可能减少对外部库的依赖
只有头文件（也可以做一个静态编译库）
函数封装性

下载Libigl

git clone https://github.com/libigl/libigl.git

Ligigl核心功能只依赖C++标准库和Eigen

但是编译起来还是挺费劲的，哈！

Mesh表示

Libigl使用Eigen库编码向量和矩阵

Eigen参考指南: Dense matrix和Spase matrix

一个三角网格编码为一对矩阵：

Eigen::MatrixXd V;
Eigen::MatrixXi F;

V是#N*3的矩阵，每一行存储一个顶点（xyz）；
F是连通关系，每一行存储组成三角形的三个顶点的下标
F中顶点的下标顺序决定了三角形的朝向

这种编码方式的好处是：

节约内存和缓冲
使用下标而不是指针，调试起来很简约
数据很容易复制和序列化

Libigl提供从用网格格式的I/O函数

例如，读取off文件并写出obj文件

igl::readOFF(path, V, F);
igl::writeOBJ(name, V, F);

可视化表面

Libigl使用glfw-based OpenGL 3.2 viewer来可视化surface。同时提供了glfw的相关性质

// Plot a mesh
igl::opengl::glfw::Viewer viewer;
viewer.data().set_mesh(V, F);
viewer.launch();

set_mesh()把网格数据拷贝给viewer，launch()创建一个glfw的窗口
默认相机模式是2-axis，可以改成3-axis trackboall模式：

viewer.core().set_rotation_type(igl::opengl::ViewerCore::ROTATION_TYPE_TRACKBALL);

可以使用标准OpenGL代码来扩展Viewer

使用键鼠进行交互

通过在Viewer中注册事件触发器的回调函数，实现键鼠交互。支持以下几种类型：

bool (*callback_pre_draw)(Viewer& viewer);
bool (*callback_post_draw)(Viewer& viewer);
bool (*callback_mouse_down)(Viewer& viewer, int button, int modifier);
bool (*callback_mouse_up)(Viewer& viewer, int button, int modifier);
bool (*callback_mouse_move)(Viewer& viewer, int mouse_x, int mouse_y);
bool (*callback_mouse_scroll)(Viewer& viewer, float delta_y);
bool (*callback_key_down)(Viewer& viewer, unsigned char key, int modifiers);
bool (*callback_key_up)(Viewer& viewer, unsigned char key, int modifiers);

使用键盘的最简单的例子就是通过按键，在不同网格模型之间切换
在初始化Viewer后，注册回调函数

viewer.callback_key_down = &key_down;

按键事件，在绘制新的mesh前记得要先清楚之前的mesh

bool key_down(igl::opengl::glfw::Viewer& viewer, unsigned char key, int modifier)
{
  if (key == '1')
  {
    viewer.data().clear();
    viewer.data().set_mesh(V1, F1);
    viewer.core.align_camera_center(V1,F1);
  }
  else if (key == '2')
  {
    viewer.data().clear();
    viewer.data().set_mesh(V2, F2);
    viewer.core.align_camera_center(V2,F2);
  }
  return false;
}

标量场的可视化

通过set_colors()可以为顶点或者面附着颜色

viewer.data().set_colors(C);

C是一个#C*3的矩阵，每一行是RGB值；C的行数必须和面片的数目或顶点的数据相等。viewer根据数量自动决定是把颜色加到顶点还是面片上

逐顶点的标量场可以直接使用set_data

viewer.data().set_data(D);

D是#V*1的向量，每个顶点一个值；set_data会在三角形内线性插值出颜色；这个颜色是通过一个colormap查找得到的；默认的colormap是igl::COLOR_MAP_TYPE_VIRIDIS，有21个离散区间；可以自定义colormapset_colormap

Overlays （勉强翻译成图示）

Viewer还支持点、线和文本标签的可视化

viwer.data().add_points(P, Eigen::RowVector3d(r,g,b));

绘制点，坐标是P的每一行，颜色统一是rgb；点的大小可以统一设置viewer.data().point_size

viewer.data().add_edges(P1,P2,Eigen::RowVector3d(r,g,b));

绘制线，连接P1和P2的对应行的坐标，P1和P2必须是同大小

viewer.data().add_label(p, str);

绘制标签，在位置p，字符串是str

例子105：使用Eigen快速找到模型的bounding box

Viewer菜单

Libigl使用Dear ImGui实现了默认菜单；也可以自制新的界面

例子106：自制菜单

多个网格

Libigl的Viewer提供绘制多个网格的功能。通过viewer.selected_data_index选择网格；默认index是0

例子107：绘制多个网格，使用键盘切换

多个视图

Libigl的Viewer提供基本的多视点绘制网格的功能

通过View::append_core()可以添加view core。最多可以添加31个cores，每个core指派了一个unsigned int的id，通过Viewer::core(id)可以访问到
如果有多于一个的view core，用户需要指定每个视窗的大小、位置。还需要明确当窗口大小变化时，这个视窗怎么跟着变

举个栗子：左右俩窗口

viewer.callback_post_resize = [&](igl::opengl::glfw::Viewer &v, int w, int h) {
  v.core( left_view).viewport = Eigen::Vector4f(0, 0, w / 2, h);
  v.core(right_view).viewport = Eigen::Vector4f(w / 2, 0, w - (w / 2), h);
  return true;
};

可以通过Viewer::selected_core_index()获取当前鼠标指向的视窗；

通过viewer.data(mesh_id).set_visible()可以控制每个mesh的可见性

Viewer Guizmos（视图小组件）

Viewer和ImGuizmo集成，为网格提供widgets

例子109：提供一个变换widget

Msh Viewer

Libigl可以读取混合网格（.msh文件格式）

MatCaps

Matcaps（material captures），也叫环境贴图，是一种简单的基于图像的绘制技术

例子111：环境贴图

离散几何参量与算子

Normals

表面法向有很多种计算方法。

例子201：计算，并可视化法向

#include <igl/read_triangle_mesh.h>
#include <igl/opengl/glfw/Viewer.h>
#include <igl/per_vertex_normals.h>
#include <igl/per_face_normals.h>
#include <igl/per_corner_normals.h>
#include <iostream>

Eigen::MatrixXd V;  // 顶点
Eigen::MatrixXi F;  // 面片

Eigen::MatrixXd N_vertices;
Eigen::MatrixXd N_faces;
Eigen::MatrixXd N_corners;

// 自定义的数据路径
const std::string DATA_PATH = "/Users/wangxc/CLionProjects/libigl-tutorial-data";

// 按键时间
bool key_down(igl::opengl::glfw::Viewer& viewer, unsigned char key, int modifier)
{
    switch(key)
    {
        case '1':
            viewer.data().set_normals(N_faces);
            return true;
        case '2':
            viewer.data().set_normals(N_vertices);
            return true;
        case '3':
            viewer.data().set_normals(N_corners);
            return true;
        default:
            break;
    }
    return false;
}

int main(int argc, char *argv[])
{
    igl::read_triangle_mesh(argc > 1 ? argv[1] : DATA_PATH + "/fandisk.off", V, F);

    // Compute per-face normals
    igl::per_face_normals(V, F, N_faces);
    // Compute per-vertex normals
    igl::per_vertex_normals(V, F, N_vertices);
    // Compute per-corner normals
    igl::per_corner_normals(V, F, 20,N_corners);

    // Plot the mesh
    igl::opengl::glfw::Viewer viewer;
    viewer.callback_key_down = &key_down;
    viewer.data().show_lines = false;
    viewer.data().set_mesh(V, F);
    viewer.data().set_normals(N_faces);
    std::cout <<
        "Press '1' for per-face normals." << std::endl <<
        "Press '2' for per-vertex normals." << std::endl <<
        "Press '3' for per-corner normals." << std::endl;
    viewer.launch();
}

逐面片per-face

逐面片法向容易出现分段constant，得到的表面法向不平滑，揭示了其内在的离散化

void per_face_normals(vertices, faces, z, normals);
// z: 一个3维向量，当三角面片退化的时候，指定的法向，例如(1,1,1)
void per_face_normals(vertices, faces, normals);
// 默认z=(0,0,0)
void per_face_normals_stable(vertics, faces, normals);
// 面片顺序不会影响输出
void per_face_normals(V, I, C, N, VV, FF, J);
// Inputs:
// 	V: vertex positions
// 	I: #I vectorized list of polygon corner indices into rows of some matrix V
// 	C: #polygons+1 list of cumulative polygon sizes. C(i+1)-C(i) = size of the ith polygon; I(C(i))到I(C(i+1)-1)是第i个多边形的下标
// Outputs:
//	N: per-face normals
//	VV: #I+#polygon * 3, list of auxiliary triangle mesh vertex positions
// 	FF: #I * 3, list of triangle indices into rows of VV
//	J: #I, list of indices into original polygons

per_face_normals
- 预设N的大小
- 遍历F的每个row（face），计算每个row里的第1个顶点减去第0个顶点的标量；第2个顶点减去第0个顶点的标量
- 叉乘
- 求模
- 判断如果模为0（退化），赋予z（可预设为z=(0,0,0)）；否则，归一化，输出结果

逐顶点per-vertex

可以在顶点计算法向，然后内插三角形内部法向（Phong Shading），得到的结果比较平滑。大多数逐顶点法向是计算相邻面法向的平均值，区别在于加权方式：uniform还是area-based，还是angle-based，后两者更符合离散化一些

半边结构的计算方式更有效一些

void per_vertex_normals(V, F, weigtingType, N)
{
  // Input: 
  // 	V: #Vx3
  //	F: #Fx3
  //	weightingType: uniform-0, area-1, angle-2, default-3,
  // Output:
  // 	N: #Vx3
}
void per_vertex_normals(V, F, N); // 省去weighting

void per_vertex_normals(V, F, weighting, FN, N)
{
  // Inputs:
  //	FN: #Fx3 matrix of face normals
}
void per_vertex_normals(V, F, FN, N)

单看函数有两类，一类是输入顶点、面片、权重方式（或忽略），输出法向；一类是输入顶点、面片、权重方式、面片法向，输出法向
per_vertex_normals(V, F, weighting, N)
- 先计算per_face_normal，返回PFN
- 然后调用per_vertex_normals(V, F, weighting, PFN, N)
- 具体计算是，遍历每个面，对面上每个顶点，加权法向（每个面参与3个顶点的累加，两层for循环，效率还行）

per-corner

大致理解per-corner法向，是计算每个面的法向量，然后再计算法向量之间的夹角。若大于某个阈值（比如20°），就记录下来

per-face normals 数量是per-vertex normals的两倍；per-corner normals的数量是per-face normals的3倍

高斯曲率

高斯曲率是主曲率的乘积；但是高斯曲率是intrinsic度量。

在连续表面，高斯曲率等于主曲率的乘积；在离散情况下，“离散高斯曲率”在一个三角网格上的一种定义方式是基于顶点的angular deficit：
$k_{G}(v_i)=2\pi -\sum_{j \in N(i)}\theta_{ij}$

template <typename DerivedV, typename DerivedF, typename DeriveddK>
IGL_INLINE void gaussian_curvature(
	const Eigen::MatrixBase<DerivedV>& V,
	const Eigen::MatrixBase<DerivedF>& F,
	Eigen::PlainObjectBase<DerivedK>& K);

调用IGL的internal_angles计算内角和

template <typename DerivedV, typename DerivedF, typename DerivedK>
IGL_INLINE void internal_angles(
  const Eigen::MatrixBase<DerivedV>& V,
  const Eigen::MatrixBase<DerivedF>& F,
  Eigen::PlainObjectBase<DerivedK> & K);

内角和的计算（每个顶点的，累加工作丢在了gaussian_curvature.cpp里面
- 如果是纯三角形，调用igl::internal_angles_using_squared_edge_lengths
- 如果还包含了非三角形面片，则通过计算夹角atant2的方式

IGL的这种按照面来遍历，然后实际上累积更改的是顶点的代码方式，很有意思。刚开始看还没转过弯来

曲率方向

平均曲率：主曲率的均值
$H=\frac{1}{2}(k_1+k_2)$
一种提取平均曲率额方法是检查应用到表面位置的Laplace-Beltrami operator。这个结果也叫做平均曲率法向：
$-\Delta \bold{x}=H\bold{n}$
另一种比较鲁邦的方法是通过quadric fitting决定主曲率和方向。

计算主曲率的方法

template <
  typename DerivedV,
  typename DerivedF,
  typename DerivedPD1,
  typename DerivedPD2,
  typename DerivedPV1,
  typename DerivedPV2>
IGL_INLINE void principal_curvature(
  const Eigen::MatrixBase<DerivedV>& V,
  const Eigen::MatrixBase<DerivedF>& F,
  Eigen::PlainObjectBase<DerivedPD1>& PD1,
  Eigen::PlainObjectBase<DerivedPD2>& PD2,
  Eigen::PlainObjectBase<DerivedPV1>& PV1,
  Eigen::PlainObjectBase<DerivedPV2>& PV2,
  unsigned radius = 5,
  bool useKring = true);