CN105427362A - 一种快速aiap形状插值算法 - Google Patents

Abstract

本发明公开了一种快速AIAP形状插值算法，包括：1)找出奇异边；2)网格分割；3)设置传播源，传播局部标价；4)初始化序列输出；5)解耦初始形状序列；6)块坐标迭代下降六个步骤。主要解决的问题是如何从给定的两个基于三角形网格的几何形状关键帧中，通过内部插值得到符合物理或视觉经验的过渡序列，并以此为基础，构建一个能够用于变形、运动移植和形状序列编辑的几何处理框架。本发明具有速度更快、矫正错切、适用于二维插值和三维插值、迭代阶段的快速收敛和抗噪特性等优点。

Description

一种快速AIAP形状插值算法

技术领域

本发明涉及计算机图形学和三维动画制作领域，尤其是指一种快速AIAP形状插值算法。

背景技术

几何形状插值是静态几何处理非常重要的一环，是快速获得几何序列的常用手段，也为形状操纵的另一操作形状移植提供依据。在计算机动画设计中，利用关键帧生成动画序列既能减少设计时间，也可降低经济成本。在几何数据压缩中，通过关键帧压缩几何序列也是较为常用的手段。

在目前比较流行的几何形状表示方式中，三角形网格由于其表达和数据结构的简单性而得到广泛的应用。复杂的形状可通过三角形面片的集合来逼近，三角形顶点集可认为是源形状的一个采样。因此，对形状的处理在某种意义下可转化成对三角形面片集合的处理。很多研究机构为了得到更加真实、符合物理或视觉经验的序列，对此做出了大量的研究，从早期的二维折线段插值到目前的三维网格插值，该项技术日趋成熟。尽管如此，现有的关于运动路径的解决方案，依然存在相当的改进空间，例如如何提高插值速度，如何生成合理的大尺度变形路径，如何抗错切等。这些问题在本发明中得到进一步的解决。

具体来说，本发明将从给定的两个基于三角形网格的几何形状关键帧中，通过内部插值得到符合物理或视觉经验的过渡序列并以此为基础，对序列进行变形、运动移植和形状序列编辑的几何处理框架。两个给定的关键帧具有相同的三角形剖分，每个顶点已对应好。在此假设下，基于物体在变形过程中三角形边长在线性变化意义下保持(as-isometric-as-possible,AIAP)的观察，模拟形状过渡的整个过程，未经AIAP处理的顶点坐标线性插值导致人模型的左臂在变形过程中产生明显的长度变化，而经过AIAP处理后，左臂的长度则符合视觉经验。形状插值一般须要解决两个问题，即(网格)形状之间的点对应以及运动路径问题。本发明假设前者已完成，着重于后者。

Li等[G.Li,L.Yang,S.Wu,W.Tan,X.Chen,andC.Xian,“Planarshapeinterpolationusingrelativevelocityfields,”Computers&Graphics,vol.37,no.5,p.364C375,2013.]提出二维AIAP插值算法，且构建了二维形状序列操纵的框架。该框架中插值算法，可认为是本发明在二维上的一个特殊情况；形状序列编辑的算法，利用了相对速度场来表达运动信息的变化，继而将相对速度场移植至新的目标(编辑对象中，以实现运动移植。虽然Li等人的工作能够成功处理二维的大变形情况，但无法直接推广到三维情况。Kircher等[S.KircherandM.Garland,“Free-formmotionprocessing,”ACMTransactionsonGraphics(TOG),vol.27,no.2,p.12,2008.]利用连接映射(connectionmap)编码三角形对之间的局部位置关系，通过极分解提取旋转分量，并对其进行非线性插值。该文献针对不同参考系，将插值分为绝对插值和相对插值两类，前者对小尺度变形具有良好的插值效果，后者则能处理大尺度变形。虽然形式化描述简洁，但相对插值算法极有可能导致不同程度错切，并且跟插值种子点的选取有较大关系。Kilian等[M.Kilian,N.J.Mitra,andH.Pottmann,“Geometricmodelinginshapespace,”inACMTransactionsonGraphics,vol.26,no.3,2007,p.64.]计算高维形状空间里面点之间的测地线距离来确定插值形状的运动路径。具体地，每个形状都可以将所有的顶点坐标以任意的顺序排列成一个高维列向量，而每个这样的列向量可以看作是高维空间上的一个点。接着，在这个空间上定义一种内积等距的黎曼度量，使之成为黎曼空间。于是，两个三维形状之间的插值便可转化为寻找黎曼空间上两个点之间的测地线。该算法由于形式化比较复杂，且全局地考虑插值序列中所有帧的所有未知量，所以随着网格规模的增大，计算量和内存消耗也非线性地增大，直接导致大模型不能直接插值，需要采用多分辨率表示，这无疑增加额外的计算量。另外，因为采用顶点位置的线性插值来初始化，所以不能处理首尾关键帧之间存在大变形的情况。虽然该方法在一定程度上能解决大变形问题，但须通过手工指定全局旋转的区域才能改善固有的形状全局错切缺陷，而且，随着模型复杂度增加，该全局错切出现的可能性也随之增大。

Lipman等[Y.Lipman,O.Sorkine,D.Levin,andD.Cohen-Or,“Linearrotation-invariantcoordinatesformeshes,”ACMTransactionsonGraphics,vol.24,no.3,pp.479–487,2005.]提出线性旋转不变坐标，即通过建立基于顶点的局部坐标系，以第一和第二离散型的方式编码该顶点一邻域的边向量，最后通过线性方程组重建模型。该方法不能够对付所有的大变形问题(例如相邻三角形夹角，即二面角变化超过180度)。Winkler等[T.Winkler,J.Drieseberg,M.Alexa,andK.Hormann,“Multi-scalegeometryinterpolation,”inComputerGraphicsForum,vol.29,no.2,2010,pp.309–318.]提出基于边的多尺度插值，将模型分为多层结构，自底向上分别对边长和二面角进行线性插值。该方法可解决大变形问题，且其效率是目前除本发明外最高效的算法。但是，该方法的成功需要较好的分层策略，而该策略的执行与模型形状密切相关。Fr等[S.Fr¨ohlichandM.Botsch,“Example-drivendeformationsbasedondiscreteshells,”Comput.Graph.Forum,vol.30,no.8,pp.2246–2257,2011.]的插值量与Winkler等相仿，并考虑了体积，但没有使用分层策略，将插值问题形式化成为一个非线性问题。由于该算法缺乏较好的初始化方法，也未能很理想地处理大变形问题。

发明内容

本发明的目的在于克服现有技术的缺点与不足，提供了一种快速AIAP形状插值算法，从给定的两个基于三角形网格的几何形状关键帧中，通过内部插值得到符合物理或视觉经验的过渡序列。

为实现上述目的，本发明所提供的技术方案为：一种快速AIAP形状插值算法，包括以下步骤：

1)找出奇异边

给定的两个基于三角形网格的几何形状关键帧，M₀和M₁，假定两者具有相同的剖分，基于最优化的标架重建找出所有的奇异边；

2)网格分割

对步骤1)中的三角形网格M₀进行分割，将M₀分割为包含普通块和奇异块两类的子网格；

3)设置传播源，传播局部标价

在步骤2)中的任意一个普通块上放置传播源，以广度优先的方式，利用对旋转角插值后的连接映射，传播边局部标架；

4)初始化序列输出

将每个局部标架的X轴取出，并以特定的网格重构方法插值边长度，最后进行边拼接，输出初始化形状序列；

5)解耦初始形状序列

对步骤4)中的初始形状序列进行解耦，并对称地考虑能量的形式化，计算每条边满足AIAP约束的相对运动场；

6)块坐标迭代下降

对步骤5)中的计算，采用了块坐标下降法，这是一个迭代的过程。迭代结束后，进行边拼接，得最终形状序列。

在步骤1)中，所述的基于最优化的标架重建，其方法为：给定第k帧的连接映射集合，需要求解一个过约束的线性方程组，定义如下的优化问题：

${\mathrm{\ϵ}}_{INI}=\underset{(i,j)\∈E}{\Σ}{w}_{ij}{\mathrm{\ϵ}}_{ijk}---\left(1\right)$

其中，

$\begin{array}{c}{\mathrm{\ϵ}}_{ijk}=\left|\right|F_{ijk}{Q}_{ijk}^0-F_{\widehat{ijk}}|{|}^2+||F_{ijk}{Q}_{ijk}^1-F_{\widehat{ijk}}|{|}^2\\ +\left|\right|F_{ijk}{Q}_{ijk}^2-F_{\widehat{ijk}}|{|}^2+||F_{ijk}{Q}_{ijk}^3-F_{\widehat{ijk}}|{|}^2\end{array}$

式中，ε_INT表示在时刻t的插值(嵌入)结果；E表示网格模型的边集，i，j为顶点索引，k表示姿态索引。上式表示，对于姿态k，可以构造关于边(i,j)的正交标架，用三元组表示为F_ijk＝＜μ_ijk,υ_ijk,η_ijk＞其中，和υ_ijk＝η_ijk×μ_ijk·η_ijk是该边对应的两个三角形面积加权所得的法向量。在不引起歧义的情况下，我们使用F_ijk表示包含上述三元组向量的3×3矩阵。

权重W_ij用于减少奇异边对整体能量的影响。奇异边定义方式是：我们首先规定两个相邻三角形的二面角的范围是0度到360度；在插值阶段，想象首尾帧某个特定二面角发生变化，若其变化大于给定阈值，且变化过程中二面角的大小能取180度，则称与此二面角对应的边为奇异边。几何上看，当相邻三角形对从“凹”变“凸”或相反，是奇异边的必要条件。在实验中，我们将阈值设为45度，对奇异边对应的W_ij设为0.001，其他边的为1。

上公式(1)将得到除旋转意义外的唯一解，前提是我们需要制定至少一条边及其局部标架作为锚。通过此方法得出的结果通常会伴有错切，但这无碍最终的结果，AIAP迭代能够将其矫正。

在步骤3)中，所述的边局部标价的传播步骤为：当普通块传播完毕后，找到其相邻的且未被访问的奇异块，利用对旋转角插值后的连接映射，按照公式(1)进行边局部标架优化。在步骤4)中，所述的网格重构方法对边的长度进行了非线性插值：

$\left|\right|e_{ijk}\left|\right|=\sqrt{(1-t_k)\left|\right|e_{ij0}\left|\right|\mathrm{\τ}+t_k\left|\right|e_{ijm}|{|}^{\mathrm{\τ}}}$

式中，e_ijk表示在t_k时刻边e_ij的插值结果，τ表示流形上某个点的切向量，

边拼接操作的定义如下：

$\underset{P_k}{\arg \min }\underset{\left(i,j\right)\∈E}{\mathrm{\Σ}}\left|\right|p_{jk}-p_{ik}-e_{ijk}|{|}^2---\left(2\right)$

式中，P_k表示第k帧，p_ik表示顶点i在t_k时刻的坐标。公式(2)这个线性方程组的稀疏矩阵不会因为不同的帧而改变，因此，采用Cholesky分解稀疏矩阵一次，便能在不同帧采用回代的方式求解。这是非常高效的操作：首先，稀疏矩阵Cholesky分解适用于对称正定矩阵，将目标矩阵分解为上三角和下三角矩阵，在某些特殊情况下(例如对角带宽较小)其操作有接近线性的时间复杂度；其次，回代的时间代价非常低，如果跟Cholesky分解的平摊时间相加，那么重构每一帧的代价也就非常低。由于求解该线性方程组最耗时的操作是Cholesky分解，所以边拼接能够以很低的平均耗时代价，完成大规模网格、多帧插值的计算任务。

迭代算法的最后，需要至少一个的迭代结束条件。本算法采用了迭代次数阈值和能量阈值双结束条件，当任意一个条件满足时，算法结束。

在步骤5)中，所述的对初始序列进行解耦用到的公式为：

${\mathrm{\ϵ}}_{INT}(i,j)=\underset{k=0}{\overset{m-1}{\Σ}}\left|\right|e_{ijk}^{\′}-c_{ijk}{R}_{ijk}{e}_{ijk}|{|}^2$

式中，c_ijk表示平移矩阵，C＝{c_ijk}，R_ijk表示旋转矩阵，R＝{R_ijk}。得益于解耦，原来庞大的线性方程组如今可以解耦成许多独立的方程组，这使得求解大规模网格的差值问题成为可能。另外，由于这些方程组都相互独立，适合使用GPU进行并行计算。

在步骤6)中，所述的块坐标迭代下降中，用到的两个参数分别为C＝{c_ijk}、R＝{R_ijk}，块坐标下降法的具体做法是：先固定边，求解对应的C和R，然后固定C和R，求解边。反复执行直至满足以下任意条件：序列中最大顶点位置的相对(网格包围盒直径的)值小于阈值；迭代次数大于阈值。

本发明与现有技术相比，具有如下优点与有益效果：

1、速度更快，分别体现在初始化阶段和非线性优化的迭代阶段。在初始化阶段使用局部标架传播为主，求解局部局域线性方程组为辅的方法。简单地讲，本算法先找出首尾两帧中二面角变化特殊的三角形对，并将三角形对的公共边记录下来，作为局部线性方程组的未知量，其他的边则使用局部标架传播的方法来求解。于是，与Kircher等人的方法和Kilian等人的方法比较，本算法无须求解庞大的全局线性方程组，而同时保证二流形表面不产生面片翻转，继而能够高效地得到符合后续迭代要求的初始化解。在非线性优化的迭代阶段，由于采用快坐标下降法，原本的非线性问题可转化为相继迭代两个二次优化问题；另外，两个二次优化问题解耦后，原本庞大的、由所有模型边序列所构成的重建方程组可转化为多组独立的、针对单条边序列的方程组，显然，这是非常适合采用GPU并行加速。无论是初始化还是非线性优化阶段，模型顶点位置的重建运算，几乎都只与其拓扑结构相关，而在本算法中拓扑结构始终不变，因此，算法中的所有顶点重建方程组系数矩阵只需一次分解，便可多次重用。

2、矫正错切。使用连接映射的方法，几乎都会产生错切效果，而错切本身也可以是因为非保长的运动路径所导致。从实验中发现，模型序列发生错切与局部标架的重建方程组的未知量个数有密切关系，未知量越多，越可能发生错切。因此，本发明的初始化算法能够有效地减少错切的风险。其次，迭代阶段能够找出保持边长的相对速度场，也就能在一定程度上矫正初始化阶段所带来的错切。

3、适用于二维插值和三维插值。Yang等人的方法主要针对二维的情况，其算法无法直接地推广到三维，包括初始化算法和迭代算法(前者须重新提出，后者须要修改)。而本发明则适用于二维和三维的情况，具有良好的维度兼容性。

4、迭代阶段的快速收敛和抗噪特性。对于大部分小尺度变形的模型，迭代阶段即使不使用效率最高的顶点位置的线性插值，只采用第一帧的姿态作为中间帧的初始化，也能够在5次迭代以内收敛。对于大部分大尺度变形的模型，初始化后迭代次数也能控制在8次以内。另外，作为迭代算法的一个特殊效果，当模型序列的中间帧(不包括首尾两帧)被少量噪声污染，通过迭代，可将噪声去除。

附图说明

图1为本发明的实验狮子模型序列。首尾是给定的两个关键帧，中间的是内插结果。

图2为本发明实验的左右两臂分别经历不同尺度的变形示意图，左臂是小尺度，右臂是大尺度。左边的模型是使用简单的顶点坐标线性插值得到的结果，右边的模型则是在左边结果的基础上对边长进行了尽量保持等距(AIAP)的操作而得到的结果。

图3为本发明展示的变形示意图。

图4为本发明的实验的运动移植示意图。第一行的首尾和第二行的第一帧是给定的帧。第一行是插值结果，第二行是参考第一行的运动所生成的序列。

图5为本发明的实验形状序列编辑示意图。第一行是给定的运动序列，第二行的第四帧是修改第一行对应帧的姿态后的结果，第二行的其他帧是编辑后的序列。

图6为本发明的实验中间帧出现了三角形翻转示意图。

图7为本发明的实验错切示意图。第一行的结果错切明显，第二行则经过AIAP迭代后的效果。

图8为本发明算法流程图。

图9为本发明算法展示的相邻边局部标架之间存在着连接映射示意图。

图10为本发明展示的直接传播标架所造成的错误。最右边的图的虚线所示的边是经检测的到的奇异边。

图11为本发明的初始化算法流程图。

图12为本发明使用二维网格序列鹅和海豚插值的结果。

图13为本发明使用带有材质的二维网格序列插值蛇和海豚的插值结果。

图14为本发明插值得到的棒的弯曲轨迹。

图15为本发明实验对称变形的插值结果。

图16为本发明实验粗糙狮子模型的大变形插值结果。

图17为本发明实验抗噪结果图。第一行为噪受声污染的序列，从第二行开始，由上到下依次是迭代2次，4次和6次的结果。

具体实施方式

下面结合具体实施例对本发明作进一步说明。

本实施例所述的快速AIAP形状插值算法，其具体情况如下：从给定的两个基于三角形网格的几何形状关键帧中，通过内部插值得到符合物理或视觉经验的过渡序列，图1的狮子序列便是内插结果。以此为基础，本发明提出快速AIAP形状插值算法能够用于变形(如图3所示)、运动移植(如图4所示)和形状序列编辑(如图5所示)的几何处理框架。两个给定的关键帧具有相同的三角形剖分，每个顶点已对应好。在此假设下，基于物体在变形过程中三角形边长在线性变化意义下保持(as-isometric-as-possible,AIAP)的观察，模拟形状过渡的整个过程，如图2所示，未经AIAP处理的顶点坐标线性插值导致人模型的左臂在变形过程中产生明显的长度变化，而经过AIAP处理后，左臂的长度则符合视觉经验。形状插值一般须要解决两个问题，即(网格)形状之间的点对应以及运动路径问题。本发明的算法可以避免出现如图6所示的大变形情况，并且能够避免出现图7所示的错切情况的出现。

首先，进行初始化阶段。在初始化阶段使用局部标架传播为主，求解局部局域线性方程组为辅的方法。简单地讲，本算法先找出首尾两帧中二面角变化特殊的三角形对，并将三角形对的公共边记录下来，作为局部线性方程组的未知量，其他的边则使用局部标架传播的方法来求解。本算法无须求解庞大的全局线性方程组，而同时保证二流形表面不产生面片翻转，继而能够高效地得到符合后续迭代要求的初始化解。初始化结果将对后面的优化有着非常重要的影响。如果初始化结果出现三角形反转或退化的情形，那么后续的优化将很有可能无法纠错，导致算法失败；即使能够通过某些方法避免上述问题，初始解也位于一个非常庞大的解空间中，如果初始化出现严重扭曲，必然要求后续的AIAP迭代次数增多；一个比较好的初始化结果，同时初始化的计算代价小，对提高插值效率有着重要的影响。

然后，进行非线性优化的迭代阶段。在非线性优化的迭代阶段，由于采用快坐标下降法，原本的非线性问题可转化为相继迭代两个二次优化问题；另外，两个二次优化问题解耦后，原本庞大的、由所有模型边序列所构成的重建方程组可转化为多组独立的、针对单条边序列的方程组。无论是初始化还是非线性优化阶段，模型顶点位置的重建运算，几乎都只与其拓扑结构相关，而在本算法中拓扑结构始终不变，因此，算法中的所有顶点重建方程组系数矩阵只需一次分解，便可多次重用。

为了结合优化和传播的优势，提出一种基于两者的混合重建算法。首先，在预处理阶段，按照定义检测出所有的奇异边。接着，将源三维模型按边分割为连通块的集合G，包含两种不同类型的块：普通块和奇异块。普通块是边界全是奇异边的最大连通子网格，而奇异块则为内部边全为奇异边、边界全为普通边的最大连通子网格。然后，给定一个传播种子(即某条边，假设放在某个普通快内)的标架，则进行标架传播，若当前块已遍历完毕，则将当前块的边界作为相邻奇异块的边界条件，求解奇异块对应的线性方程组；当奇异块求解完毕后，在转向普通块的传播。反复这个过程，就能完成整个网格局部标架集合的求解。为了能进一步减少奇异块中的未知量，我们通过简单的混合，把邻边标架皆为已知的奇异边标架解析地求解出来。

本实施例上述的快速AIAP形状插值算法，算法流程如图8所示，包括以下步骤：

1)找出奇异边

给定的两个基于三角形网格的几何形状关键帧，M₀和M₁，假定两者具有相同的剖分，基于最优化的标架重建找出所有的奇异边。基于最优化的标架重建，其方法为：给定第k帧的连接映射集合，需要求解一个过约束的线性方程组，定义如下的优化问题：

${\mathrm{\ϵ}}_{INI}=\underset{(i,j)\∈E}{\Σ}{w}_{ij}{\mathrm{\ϵ}}_{ijk}---\left(1\right)$

其中，

$\begin{array}{c}{\mathrm{\ϵ}}_{ijk}=\left|\right|F_{ijk}{Q}_{ijk}^0-F_{\widehat{ijk}}|{|}^2+||F_{ijk}{Q}_{ijk}^1-F_{\widehat{ijk}}|{|}^2\\ +\left|\right|F_{ijk}{Q}_{ijk}^2-F_{\widehat{ijk}}|{|}^2+||F_{ijk}{Q}_{ijk}^3-F_{\widehat{ijk}}|{|}^2\end{array}$

式中，ε_INT表示在时刻t的插值(嵌入)结果；E表示网格模型的边集，i，j为顶点索引，k表示姿态索引。上式表示，对于姿态k，可以构造关于边(i,j)的正交标架，用三元组表示为F_ijk＝＜μ_ijk,υ_ijk,η_ijk＞其中，和υ_ijk＝η_ijk×μ_ijk·η_ijk是该边对应的两个三角形面积加权所得的法向量。在不引起歧义的情况下，我们使用F_ijk表示包含上述三元组向量的3×3矩阵

权重W_ij用于减少奇异边对整体能量的影响。奇异边定义方式是：我们首先规定两个相邻三角形的二面角的范围是0度到360度；在插值阶段，想象首尾帧某个特定二面角发生变化，若其变化大于给定阈值，且变化过程中二面角的大小能取180度，则称与此二面角对应的边为奇异边。几何上看，当相邻三角形对从“凹”变“凸”或相反，是奇异边的必要条件。在实验中，我们将阈值设为45度，对奇异边对应的W_ij设为0.001，其他边的为1。

上公式(1)将得到除旋转意义外的唯一解，前提是我们需要制定至少一条边及其局部标架作为锚。通过此方法得出的结果通常会伴有错切，但这无碍最终的结果，AIAP迭代能够将其矫正。

2)网格分割

对找出奇异边时使用的三角形网格M₀进行分割，将M₀分割为包含普通块和奇异块两类的子网格。

3)设置传播源，传播局部标价

在完成分割后的网格上的任意一个普通块上放置传播源，以广度优先的方式，利用对旋转角插值后的连接映射，传播边局部标架。连接映射的示意图如图9所示。边局部标价的传播步骤为：当普通块传播完毕后，找到其相邻的且未被访问的奇异块，利用对旋转角插值后的连接映射，按照公式(1)进行边局部标架优化。如果采用不正确的标价传播策略，对某些模型做初始化会产生如图10所示的错误。合理的传播方式是：当奇异块优化完毕后，找到期相邻的且未被访问的普通块，进行类似的边局部脚架传播。反复执行直至所有块都已访问完毕。图11描述了标价传播的算法流程。

4)初始化序列输出

将每个局部标架的X轴取出，并以特定的网格重构方法插值边长度，最后进行边拼接，输出初始化形状序列。网格重构方法对边的长度进行了非线性插值：

$\left|\right|e_{ijk}\left|\right|=\sqrt{(1-t_k)\left|\right|e_{ij0}\left|\right|\mathrm{\τ}+t_k\left|\right|e_{ijm}|{|}^{\mathrm{\τ}}}$

边拼接操作的定义如下：

式中，e_ijk表示在t_k时刻边e_ij的插值结果，τ表示流形上某个点的切向量，边拼接操作的定义如下：

$\underset{P_k}{\arg \min }\underset{\left(i,j\right)\∈E}{\mathrm{\Σ}}\left|\right|p_{jk}-p_{ik}-e_{ijk}|{|}^2---\left(2\right)$

式中，P_k表示第k帧，p_ik表示顶点i在t_k时刻的坐标。公式(2)这个线性方程组的稀疏矩阵不会因为不同的帧而改变，因此，采用Cholesky分解稀疏矩阵一次，便能在不同帧采用回代的方式求解。这是非常高效的操作：首先，稀疏矩阵Cholesky分解适用于对称正定矩阵，将目标矩阵分解为上三角和下三角矩阵，在某些特殊情况下(例如对角带宽较小)其操作有接近线性的时间复杂度；其次，回代的时间代价非常低，如果跟Cholesky分解的平摊时间相加，那么重构每一帧的代价也就非常低。由于求解该线性方程组最耗时的操作是Cholesky分解，所以边拼接能够以很低的平均耗时代价，完成大规模网格、多帧插值的计算任务。

迭代算法的最后，需要至少一个的迭代结束条件。本算法采用了迭代次数阈值和能量阈值双结束条件，当任意一个条件满足时，算法结束。

5)解耦初始形状序列

对初始形状序列进行解耦，并对称地考虑能量的形式化，计算每条边满足AIAP约束的相对运动场。对初始序列进行解耦用到的公式为：

${\mathrm{\ϵ}}_{INT}(i,j)=\underset{k=0}{\overset{m-1}{\Σ}}\left|\right|e_{ijk}^{\′}-c_{ijk}{R}_{ijk}{e}_{ijk}|{|}^2$

式中，c_ijk表示平移矩阵，C＝{c_ijk}，R_ijk表示旋转矩阵，R＝{R_ijk}。得益于解耦，原来庞大的线性方程组如今可以解耦成许多独立的方程组，这使得求解大规模网格的差值问题成为可能。另外，由于这些方程组都相互独立，适合使用GPU进行并行计算。

6)块坐标迭代下降

在解耦初始形状序列是用到的计算，采用了块坐标下降法，这是一个迭代的过程。迭代结束后，进行边拼接，得最终形状序列。块坐标迭代下降中，用到的两个参数分别为C＝{c_ijk}、R＝{R_ijk}，块坐标下降法的具体做法是：先固定边，求解对应的C和R，然后固定C和R，求解边。反复执行直至满足以下任意条件：序列中最大顶点位置的相对(网格包围盒直径的)值小于阈值；迭代次数大于阈值。

综上所述，在采用以上方案后，本发明提出了一种新的、基于三角形网格的快速形状插值算法，可使生成的形状序列在线性变化意义下尽可能地保持边长(AIAP)，即通过寻找保持边长的速度场来解决两个关键帧之间的运动路径问题。本形状插值算法的技术特点在于：

a)基于三角形网格的边，利用块坐标下降算法将非线性的AIAP最优化问题转化为简单的二次优化问题，然后使用分治策略对未知量进行解耦，从而高效地求解。

b)针对三维网格的大变形问题，采用基于边的连接映射(connectionmap)技术来初始化非线性AIAP问题的解。为了提高初始解的求解速度，提出一种新的、基于连接映射插值的优化与传播混合算法。该算法基于以下观察：边的正交标架在奇异边集中的区域传播时会不稳定，传播错误累计明显，需要最优化策略来解决。该算法不仅高效，而且能够非常方便地移植到其他已有的非线性形状插值算法中，为其提供初始化解。

本发明经过实验证明其可行性，能广泛应用于各种模型。图12、图13分别展示了为本发明在二维情况下的插值结果以及为本发明对带材质的二维序列插值的结果。图14、图15、图16分别展示了本发明对弯曲棒的插值结果、对对称变形的插值结果以及对大变形序列的插值结果。图17表明了本发明具有很好的抗噪特性。从结果可以看出，本发明插值速度快，鲁棒性强，能够很好地处理大变形的情况。

以上所述之实施例子只为本发明之较佳实施例，并非以此限制本发明的实施范围，故凡依本发明之形状、原理所作的变化，均应涵盖在本发明的保护范围内。

Claims (5)

1.一种快速AIAP形状插值算法，其特征在于：包括以下步骤：

1)找出奇异边

给定的两个基于三角形网格的几何形状关键帧，M₀和M₁，假定两者具有相同的剖分，基于最优化的标架重建找出所有的奇异边；

2)网格分割

对步骤1)中的三角形网格M₀进行分割，将M₀分割为包含普通块和奇异块两类的子网格；

3)设置传播源，传播局部标价

在步骤2)中的任意一个普通块上放置传播源，以广度优先的方式，利用对旋转角插值后的连接映射，传播边局部标架；

4)初始化序列输出

将每个局部标架的X轴取出，并以特定的网格重构方法插值边长度，最后进行边拼接，输出初始化形状序列；

5)解耦初始形状序列

对步骤4)中的初始形状序列进行解耦，并对称地考虑能量的形式化，计算每条边满足AIAP约束的相对运动场；

6)块坐标迭代下降

对步骤5)中的计算，采用了块坐标下降法，这是一个迭代的过程，迭代结束后，进行边拼接，以得最终形状序列。

2.根据权利要求1所述的一种快速AIAP形状插值算法，其特征在于：在步骤1)中，所述的基于最优化的标架重建，其方法为：给定第k帧的连接映射集合，需要求解一个过约束的线性方程组，定义如下的优化问题：

${\mathrm{\ϵ}}_{INI}=\underset{(i,j)\∈E}{\Σ}{w}_{ij}{\mathrm{\ϵ}}_{ijk}---\left(1\right)$

其中，

$\begin{array}{c}{\mathrm{\ϵ}}_{ijk}=\left|\right|F_{ijk}{Q}_{ijk}^0-F_{\widehat{ijk}}|{|}^2+||F_{ijk}{Q}_{ijk}^1-F_{\widehat{ijk}}|{|}^2\\ +\left|\right|F_{ijk}{Q}_{ijk}^2-F_{\widehat{ijk}}|{|}^2+||F_{ijk}{Q}_{ijk}^3-F_{\widehat{ijk}}|{|}^2\end{array}$

式中，ε_INT表示在时刻t的插值结果；E表示网格模型的边集，i，j为顶点索引，k表示姿态索引；上式表示，对于姿态k，可以构造关于边(i,j)的正交标架，用三元组表示为F_ijk＝＜μ_ijk,υ_ijk,η_ijk＞其中，和υ_ijk＝η_ijk×μ_ijk·η_ijk是该边对应的两个三角形面积加权所得的法向量；在不引起歧义的情况下，使用F_ijk表示包含上述三元组向量的3×3矩阵；

权重W_ij用于减少奇异边对整体能量的影响；奇异边定义方式是：首先规定两个相邻三角形的二面角的范围是0度到360度；在插值阶段，想象首尾帧某个特定二面角发生变化，若其变化大于给定阈值，且变化过程中二面角的大小能取180度，则称与此二面角对应的边为奇异边；几何上看，当相邻三角形对从“凹”变“凸”或相反，是奇异边的必要条件；在实验中，将阈值设为45度，对奇异边对应的W_ij设为0.001，其他边的为1；

上公式(1)将得到除旋转意义外的唯一解，前提是需要制定至少一条边及其局部标架作为锚；通过此方法得出的结果通常会伴有错切，但这无碍最终的结果，AIAP迭代能够将其矫正；

在步骤3)中，所述的边局部标价的传播步骤为：当普通块传播完毕后，找到其相邻的且未被访问的奇异块，利用对旋转角插值后的连接映射，按照公式(1)进行边局部标架优化；当奇异块优化完毕后，找到期相邻的且未被访问的普通块，进行类似的边局部脚架传播；反复执行直至所有块都已访问完毕。

3.根据权利要求1所述的一种快速AIAP形状插值算法，其特征在于：在步骤4)中，所述的网格重构方法对边的长度进行了非线性插值：

$\left|\right|e_{ijk}\left|\right|=\sqrt{(1-t_k)\left|\right|e_{ij0}\left|\right|\mathrm{\τ}+t_k\left|\right|e_{ijm}|{|}^{\mathrm{\τ}}}$

式中，e_ijk表示在t_k时刻边e_ij的插值结果，τ表示流形上某个点的切向量，

边拼接操作的定义如下：

$\underset{P_k}{\arg \min }\underset{\left(i,j\right)\∈E}{\mathrm{\Σ}}\left|\right|p_{jk}-p_{ik}-e_{ijk}|{|}^2---\left(2\right)$

式中，P_k表示第k帧，p_ik表示顶点i在t_k时刻的坐标；公式(2)这个线性方程组的稀疏矩阵不会因为不同的帧而改变，因此，采用Cholesky分解稀疏矩阵一次，便能在不同帧采用回代的方式求解；首先，稀疏矩阵Cholesky分解适用于对称正定矩阵，将目标矩阵分解为上三角和下三角矩阵，在对角带宽较小的情况下，其操作有接近线性的时间复杂度；其次，回代的时间代价非常低，如果跟Cholesky分解的平摊时间相加，那么重构每一帧的代价也就非常低；由于求解该线性方程组最耗时的操作是Cholesky分解，所以边拼接能够以很低的平均耗时代价，完成大规模网格、多帧插值的计算任务；

迭代算法的最后，需要至少一个的迭代结束条件；本算法采用了迭代次数阈值和能量阈值双结束条件，当任意一个条件满足时，算法结束。

4.根据权利要求1所述的一种快速AIAP形状插值算法，其特征在于：在步骤5)中，所述的对初始序列进行解耦用到的公式为：

${\mathrm{\ϵ}}_{INT}(i,j)=\underset{k=0}{\overset{m-1}{\Σ}}\left|\right|e_{ijk}^{\′}-c_{ijk}{R}_{ijk}{e}_{ijk}|{|}^2$

式中，c_ijk表示平移矩阵，C＝{c_ijk}，R_ijk表示旋转矩阵，R＝{R_ijk}；得益于解耦，原来庞大的线性方程组如今可以解耦成许多独立的方程组，这使得求解大规模网格的差值问题成为可能；另外，由于这些方程组都相互独立，适合使用GPU进行并行计算。

5.根据权利要求1所述的一种快速AIAP形状插值算法，其特征在于：在步骤6)中，所述的块坐标迭代下降中，用到的两个参数分别为C＝{c_ijk}、R＝{R_ijk}，块坐标下降法的具体做法是：先固定边，求解对应的C和R，然后固定C和R，求解边；反复执行直至满足以下任意条件：序列中最大顶点位置的相对值小于阈值；迭代次数大于阈值。