Mesh Deformation Transfer

Deformation Transfer For Triangle Mesh

概述

问题描述

对于两个具有一定视觉相似度的三角网格：原网格 S0 和目标网格 T0，根据原网格已知的变形序列 S={S0,…,S|S|}，生成目标网格的对应的变形序列T={T0,…,T|S|}。

banner

问题假设

1. S 和 T 中的网格分别具有相同的拓扑结构，但两个集合对应网格之间不要求。
2. S0 与 T0 应当具有一定的视觉相似度，并且相关点对通过人为标记的方式体现。

基本思想

1. 根据人为标定的标记点，计算由 S0 到 T0 的三角面片对应关系 M={(s0,t0),…,(s|M|,t|M|)}。
2. 对 Si,i∈{1,…,|S|} ，由于 Si 与 S0 具有相同的拓扑，可以根据对应关系M，将Si每个三角面的变换作用到T0，加上一些约束条件，得到变形后的 Ti 。

三角面变换

在三维空间中，对于三角面Fi=[v0,v1,v2] 到另一个三角面~Fi=[˜v0,˜v1,˜v2] 的仿射变换可以分解为线型部分Qi 和非线性部分 di。求解该仿射变换需要用四对点到点的关系，对每个三角面引入第四个点： v3=v0+(v1−v0)×(v2−v0)/√|(v1−v0)×(v2−v0)|

线性部分Qi=˜ViVi−1，其中Vi=[v1−v0v2−v0v3−v0]

对应关系计算

根据标定点，将 S0 在保持拓扑不变的前提下变形为 T0。

根据三角面变换中的定义，将S0中每个三角面的第四个点加到顶点序列之后： x=(v0,v1,…,vn−1⏟原始顶点,vn,…,vn+m−1⏟新增顶点)

其中n为S0中原始顶点的个数，m为三角面个数，x∈R3×(n+m)。

通过改动x，来实现S0到T0的变形，具体表现为定义损失函数，最小二乘法搜索最优解： min˜x(wSES,wIEI,wCEC)约束标记点对中 原网格顶点与目标网格顶点相同

损失项	表达式	备注
平滑性 (smoothness)	ES(x)=∑Qi∑Qj∈adj(Qi)‖Qi−Qj‖2F	对每一个三角面，周围三角面的形变应该尽量与之相似
不变性 (identity)	EI(x)=∑Qi‖Qi−I‖2F	每个三角面的形变尽可能小
最近点损失	EC(x;c0,…,cn)=∑i‖vi−ci‖2F	原网格的每个顶点都应该尽可能贴近与目标网格的最近点

当S0在保持拓扑不变的情况下变形为接近 T0后，计算两个网格三角面之间的对应关系M。

对S0中的每个三角面，在T0中至少存在一个三角面与之存在最近关系，反之亦然。所以两个网格中三角面的对应关系为多对多，并非映射。

形变迁移

在有原网格与目标网格三角面的对应关系M, 以及原网格 初始网格S0 与形变网格Si 的每个三角面形变关系后，我们可以直接将S0 到Si的形变迁移到T0 上，得到Ti：

minQi+di|M|−1∑j=0‖Ssj−Ttj‖2F约束 Qjvi+dj=Qkvi+dk,∀i,∀j,k∈adj(vi)

这里的Ssj表述为网格Si中标号为sj的三角形的变换Qsj，Ttj同理。

对于实际求解，可以将上述对三角面变换的逼近转换到对顶点变换的逼近。

由于三角面变换逼近时，可能会出现边缘撕裂的情况，所以需要添加三角形的邻域约束。而转用顶点表达时，通过同一网格类型拓扑不变的假设，蕴含了变换后的网格不会出现边缘撕裂。

具体的顶点公式推导在后面。

细节推导

三角面变换

以上描述的优化函数大多都是用三角面的形变表示，而优化目标是顶点序列x，需要进行推导将三角面的形变转换为关于顶点序列 x 的表达式，即。 Qi=˜ViVi−1展开顶点序列 x→Qi=xˆV−1i

计算展开x Qi3×3=˜ViVi−1=[˜vi1−˜vi0˜vi2−˜vi0˜vi3−˜vi0]Vi−1=[˜vi1˜vi2˜vi3]Vi−1−[˜vi0˜vi0˜vi0]Vi−1=[˜v0˜v1…˜vn+m−1][…000…100…010…001…]Vi−1−[˜v0˜v1…˜vn+m−1][…111…000…000…000…]Vi−1=[˜v0˜v1…˜vn+m−1][…−1−1−1…100…010…001…]Vi−1=xˆV−1i

一般的目标函数

对于一般情况，期望在保持原网格拓扑不变的同时满足每个三角面的变形目标Ci，写作： |M|−1∑i=0‖Qi−Ci‖2F=|M|−1∑i=0‖xˆV−1i−Ci‖2F=|M|−1∑i=0‖(xˆV−1i−Ci)T‖2F=|M|−1∑i=0‖(ˆV−1i)TxT−CTi‖2F=‖((ˆV−10)T(ˆV−11)T…(ˆV−1|M|−1)T)xT−(CT0CT1…CT|M|−1)‖2F=‖AxT−b‖2F

即目标函数变为 E(x)=‖AxT−b‖2F 的形式，使用最小二乘法^[1]，得到 ∂E(x)∂x=∂(xTATAx−2bTAx+bTb)∂x=2ATAx−2ATb=0ATAx=ATb

具体的目标函数

平滑性损失 ES(x)=∑Qi∑Qj∈adj(Qi)‖Qi−Qj‖2F=‖((ˆV−10−ˆV−1j0)T(ˆV−10−ˆV−1j1)T…(ˆV−10−ˆV−1j|adj(Qi)|−1)T(ˆV−11−ˆV−1j0)T…)xT−0‖2F=‖ASxT−bS‖2F

其中AS∈R3q×(n+m)， bS∈R3q×3，q=∑i|adj(Qi)|。

不变性损失 EI(x)=∑Qi‖Qi−I‖2F=‖((ˆV−10)T(ˆV−11)T…(ˆV−1m−1)T)xT−(II…I)‖=‖AIxT−bI‖2F

其中AI∈R3m×(n+m)，bI∈R3m×3。

最近点损失 EC(x;c0,…,cn)=∑i‖vi−ci‖2=‖IxT−(CT0CT1…CTm+n−1)‖2F=‖ACxT−bC‖2F

其中AC∈R3m×(n+m)，bC∈R3m×3。

为了统一使用变量xT，需要有n+m个点的形式，但实际上只需要计算前n个顶点的最近点。

形变迁移损失 EQ(x)=|M|−1∑j=0‖Ssj−Ttj‖2F=|M|∑j=1‖STsj−ˆV−1tjTxT‖2F=‖−((ˆV−1t0)T(ˆV−1t1)T…(ˆV−1t|M|−1)T)xT+(STs0STs1…STs|M|−1)‖2F=‖AQxT−bQ‖2F

其中AQ∈R3|M|×(n+m)，bQ∈R3|M|×3，|M|是对应关系中元素的个数，满足|M|≥m。

实验

注意事项

对应关系求解: 标记点约束

在求解minxES、minxEI和minxEQ时，需要约束S0和T0对应的标记点相同： ATAxT=ATb(Au)TAu(xu)T+(Am)TAm(˜xm)T=ATb(Au)TAu(xu)T=ATb−(Am)TAm(˜xm)T

对于S0，Au、xu分别为未标记点的系数矩阵以及顶点序列，Am、xm为已标记点的系数矩阵以及顶点序列，˜xm为T0中已标记点的顶点序列。

再求解出xu后再根据约束条件xm=˜xm 计算出结果 ˜x。

对应关系求解: 目标函数组合与迭代

对应关系计算时，可以将相关目标函数写在一起： ˜x=argminx‖(wSASwIAIwCAC)xT−(wSbSwIbIwCbC)‖2F

其中wS=1.0，wI=0.001，wC=[0,1…,5000]。

稀疏线性方程组求解^[2]

1. 直接求解
   • LU分解：Ax=b⟹LUx=b
   • Cholesky分解^[3]：Ax=b⟹LTLx=b
2. 间接求解
   • Jacobi method / Gauss-Seidel method

优化

1. 稀疏矩阵的计算：乘法、转置、切片、拼接和索引。
2. 稀疏方程组求解。
3. 计算EC时，最近点的计算。
4. 计算对应关系M时，最近三角面重心以及法线夹角的计算。

计算流程

计算形变迁移，本质上在优化四个目标函数：ES、EI、EC和EQ。其中系数矩阵几乎都包含ˆV−1的计算。

参考

1. https://eeweb.engineering.nyu.edu/iselesni/lecture_notes/least_squares/least_squares_SP.pdf ↩︎
2. https://zhuanlan.zhihu.com/p/479913328 ↩︎
3. https://zhuanlan.zhihu.com/p/112091443 ↩︎