来源：第三维度
    作者：孙守迁+, 徐爱国, 黄琦, 王鑫
    单位：浙江大学计算机科学与技术学院

    摘要: 提出了一种快速处理三维服装仿真中角色与服装碰撞的方法.该方法能够满足交互式实时仿真环境的需求.在预处理阶段,根据蒙皮动画的特点,将角色的几何形状以球和简化凸包等简单碰撞体进行估计.在实时模拟阶段,这些碰撞体跟随骨架运动并代替角色模型网格完成与虚拟服装之间的碰撞处理.此外,为了能够快速计算碰撞响应信息,该方法还利用外围映射机制进一步开发了相交测试的空间局部性.实验结果表明,应用该方法可以有效避免衣片与角色模型之间的相互穿透,同时大幅度地减少碰撞处理计算量.实时仿真系统对于复杂服装网格仍然保持了较高的模拟帧速率.

    近年来,基于物理的服装仿真动画技术受到了国内外学者的广泛关注.一方面,它可以生成以往无法比拟的真实感模拟效果;另一方面,计算机硬件性能的不断提高也促使这项技术逐步趋于实用化.与其相关的研究工作主要涉及3个方面的内容,即布料物理建模、动力学微分方程数值求解以及可变形体碰撞检测和响应.其中,碰撞检测和响应是实时仿真的一个重要计算瓶颈.布料作为高度可变形体与虚拟角色模型发生大量接触碰撞,设计一种高效、精确的检测算法难度极大,特别是在兼顾系统稳定性和实时性的情况下.

    已有的一些方法更侧重于碰撞检测和响应的精确性,其中比较重要的一类是以层次包围体加速进行几何基元相交测试[1-5],可选择的包围体类型涉及轴向包围盒(axis-aligned bounding boxes,简称AABB)、有向包围盒(oriented bounding box,简称OBB)以及离散有向多面体(k-discrete orientation polytopes,简称k-DOP)等.显然,对于可变形体,需要在模拟过程中频繁更新层次数据结构,这使得形状简单的包围体反而表现出更高的处理效率.但从总体上看,该类方法对于一个实时模拟系统来说计算量依然过大.另一类加速方法基于空间剖分.Zhang[6]采用哈希表索引来降低剖分体素化信息的存储需求,并得到了比层次包围体方法更快的碰撞处理速度.然而,由于索引基元可能会在相邻离散检测时间点内跨越多个体素,对于大步长模拟,该方法并不稳定,文献[7]对它的应用仅限制在1/90s,所以很难满足实时系统的需求.Teschner[8]对包括哈希函数、哈希表的大小以及剖分单元尺寸在内的一些参数进行了分析优化,并且采用四面体网格作为碰撞处理对象,虽然作者表示算法可以扩展到三角形网格,但没有给出关于模拟时间步长以及可能造成的效率损失等问题的详细讨论.

    基于图像空间的方法充分利用图形绘制硬件来减少计算负担.Vassilev[9]首先将该思想用于人物角色与服装的碰撞检测,这类方法不需要进行任何预处理,并且可以利用图形硬件的可编程能力进一步加速.然而在文献[10]中Heidelberger指出,由于需要进行多次绘制来提取不同检测对象的深度信息,包括频繁的显存读取操作,硬件加速策略更适合于大规模场景,对于几千个面片的低复杂度仿真环境,其整体性能并非总是优于软件计算.

    文献[11,12]的研究工作与本文相关度较大.Rudomin[11]用椭球体对角色模型进行估计,将变形网格之间的碰撞接触问题转换为衣片质点相对于椭球基元的穿透测试,极大地提高了系统整体性能.仿真过程中,衣片质点先随着距离最近的椭球体移动,然后再进行物理模拟和碰撞处理,保证了系统的稳定性.然而由此带来的问题也很明显:没有发生碰撞的质点受椭球体带动并非基于物理模拟,这在某些情况下会产生较大的几何效应,特别是对于宽松的服装,将显著影响仿真效果.为此,Oh[12]提出了相对模拟的思想,服装作为一个整体仅随角色模型平移,除此之外的质点运动由后续的物理模拟和碰撞共同决定,这使得几何效应减弱并提高了模拟真实感,但在保证系统的稳定性方面要涉及更多的计算.从角色模型估计的方式来看,Oh的方法与文献[11]类似,主要是采用两端截断椭球体和球体作为拟合基元.一方面需要用户交互调整这些基元的描述参数来获得较好的拟合质量,另一方面,由于球体自身的几何特性,原始角色模拟的细节部分将被过度平滑,特别是在拟合基元数目比较少的情况下,通常会造成角色形体特征的严重失真,例如过于臃肿,进而降低最终的仿真质量.

    本文提出的方法同样基于碰撞体估计,但我们充分考虑了蒙皮动画的特性以及碰撞体对角色模型的拟合能力,进一步的算法优化还包括针对系统实时性的凸包约束简化和预计算外围映射体.另外,在保证系统稳定性方面,我们通过引入质点几何约束机制来避免文献[12]中的复杂求交计算,而与文献[11]相比,该策略产生的几何效应更小.

    1 物理模拟系统

    布料物理建模方法典型分为两类,即连续力学模型和粒子系统.其中,粒子系统将力相互作用关系局部化,计算量相对较小,是实时仿真系统中普遍采用的模型.对于动力学微分方程数值求解问题,显式方法由于受系统稳定性的限制,通常只能取较小的时间步长,隐式方法虽然可以进行大步长积分,但要求解规模较大的稀疏线性方程组,最后,完整的模型系统还要高效处理衣片与环境对象(例如人物角色)的接触碰撞问题.

    1.1 布料物理建模

    为了保证系统的实时性,质点弹簧模型作为最简单的粒子系统是合适的选择.另外,考虑到任意衣片形状的生成以及多个衣片之间的有效缝合拼接,以三角网格作为布料的几何载体更为方便.图1是建立在三角网格上的布料质点弹簧结构,网格顶点对应物理模型质点,边定义为拉伸剪切弹簧,两个相邻三角形非公共顶点的连接线段定义为弯曲弹簧.系统包含内力和外力,其中,外力计算如下:

    分别对应重力、阻尼力(理论上为内力,为计算方便,作为外力引入)和流体粘滞力(如风力).

Fig.1 The structure of mass-spring model
图1 质点弹簧模型结构

    在该模型中,内力并没有显式计算,而是作为一个相对约束条件通过弹簧迭代校正过程体现.给定一段弹簧的当前状态,计算总的校正向量:

    xi和xj为弹簧校正前两端质点的位置,为弹簧自由状态下的长度.然后计算校正后的质点目标位置:

    对应质点i和j的逆质量(质量的倒数).该校正过程很直观,而且通过将其中一个质点的逆质量置零就可以方便地对它进行约束与文献[13]单纯校正超弹簧效应不同,为了模拟弹簧的完整性质,该模型对弹簧的过度压缩也进行了校正.

    1.2 数值积分方法

    由第1.1节可知,本文的布料物理模型只计算外力,弹簧内力作为相对约束引入.由于外力在数值上通常比较小,因此大步长显式积分方法也适用,但考虑到碰撞引入的非连续性振荡在极端情况下也会造成求解系统的崩溃,这里选用更为稳定的Verlet积分器:

    其中,h为时间步长.该积分器计算简单,只涉及质点位置,没有速度项的显式求解,因而也不存在速度和位置的同步问题,系统能够在较大程度上保持稳定.对于质点位置积分的误差阶数为4,该精度对于一个实时模拟系统来说已经足够了.另外,部分外力的计算需要用到质点速度项,可以由位置信息进行估计.

    1.3 碰撞检测与响应

    在本文的系统中,碰撞检测只针对于布料质点与环境几何体的穿透,碰撞响应作为绝对约束引入,采用几何响应策略,并且仅考虑质点位置的调整.以球碰撞对象为例,质点发生穿透当且仅当下式成立:

    其中,x为质点位置,ccenter和r分别对应球心和半径.碰撞响应点计算如下:

    这里,我们假设碰撞为完全非弹性,即将质点平移到球体表面.采用这种碰撞处理方法效率很高,但存在比较明显的缺陷.如图2(a)所示,质点平移到球体表面却无法避免边的穿透,特别是当布料网格比较稀疏时.另外,质点过于靠近碰撞体表面也容易引发深度缓存精度不够的问题,进而影响模拟的绘制质量.一种有效的解决办法是对几何实体以厚度阈值δ进行膨胀,在其表面形成一层薄的包围体,并将检测边界外推到该包围体.如图2(b)所示,增加球体碰撞检测半径避免了边的穿透和相关的一些绘制问题.

Fig.2 Avoiding edge penetration with inflation of object
图2 对象膨胀避免边穿透

    2 基于简化凸包的碰撞处理

    以几何形状简单的包围体代替原始模型网格进行相交测试可以极大地提高系统性能,为此引入简化凸包.首先,在三维空间所有类型的包围体中,凸包相对于原始模型的紧致性最好;其次,对于非凸模型,可以通过有限的体分解来获取更好的拟合性能;第三,增加凸包约束简化,算法能够在模拟质量和速度之间进行有效平衡;最后,凸包自身的特殊性质使一些耗时的几何包含和相交测试过程得到简化,并且还可以借助预计算进一步加速.

    2.1 三维凸包

    凸包作为一种基本结构,在很多领域都有重要应用,其更广泛意义上的定义是在n维空间.本文主要利用三维凸包来处理非精确碰撞问题,因此有n=3.为了下面叙述清楚起见,首先引入相关的基本概念和性质定义.

    定义1. 三维空间中由若干个平面多边形围成的几何图形称为多面体.围成多面体的各个多边形叫做多面体的面,两个面的公共边叫做多面体的棱边,邻接棱边的公共点为多面体的顶点,连接不在同一个面上的两个顶点的线段为多面体的对角线.

    定义2. 把多面体的任意一个面展开成平面,如果其他面都位于这个平面的同一侧,则把这种特殊的多面体称为凸多面体.

    性质1. 凸多面体的所有对角线均包含在该凸多面体内,所有棱边都为凸棱边.

    性质2. 凸多面体中任意两个相邻面所形成的二面角是凸的(<π),与某一顶点关联的棱边之间形成的平面角度之和小于2π.

    性质3. 包含在凸多面体内的所有空间点都位于其任意一个面的同一侧,即负侧.

    定义3. 给定三维空间的N个点p1,p2,…,pn,凸包是包含该点集的最小凸多面体,可表示为

    本文主要处理三角网格凸包,即凸包的每个面均为三角形.另外,我们还假设给定的空间点集中至少有4个点不共面,这样可以保证其上的凸包总是具备体的性质,不会出现到三维点、线以及面的退化.

    有多种方法可以计算三维凸包,例如QuickHull算法[14],凸包的生成过程仅包含在模型预处理阶段.

    2.2 凸包约束简化

    给定原始模型网格,我们以定义在其顶点集合上的三维凸包作为碰撞替代体.考虑到碰撞处理过程中的边穿透问题,同样需要对凸包进行膨胀.简单的方法是,对每个凸包顶点以厚度阈值δ沿其各自法向量向外平移.然而,膨胀后得到的多面体并不能保持凸性,要由新顶点位置重新计算凸包.

    凸包面基元数目直接影响算法效率,为了进一步提高系统的性能,本节引入凸包网格简化来有效平衡碰撞处理的精度与速度.这里,我们采用二次误差面删除简化算法,单步操作如图3所示.

Fig.3 Simplification of convex hull mesh
图3 凸包网格的简化

    除了一般的简化操作有效性检验,如防止面翻转以外,针对凸包网格还需定义额外的两个约束测试条件:

    (1) 一致性保持.要保证简化结果仍为凸多面体,如图3所示,单步三角形删除操作仅影响多面体网格的一个局部,如果相关棱边都为凸棱边,则根据性质1,可以断定该操作是保持凸性的.为此,引入单个棱边的凸凹性检测规则:对于多面体的一条棱边,由性质2,判定为凸当且仅当其相邻两个三角形所形成的二面角小于π.

    检测仅在简化影响区域的局部进行,涉及与新顶点vnew相关的所有三角形包含的棱边.

    (2) 体包含约束.简化凸包必须确保对原始凸包以及网格模型的包含性.比较简单的办法是根据性质3,测试被删除面T的3个顶点与vnew相关的所有三角形面的空间位置关系,满足包含性当且仅当它们位于所有这些平面之上或负侧.

    显然,较大的δ值会得到更高的简化倍率.另外,用户还可以通过设定剩余三角形数目下限来控制简化结果.

    2.3 点在凸包内的快速判别

    判断布料质点穿透要做点在凸包内的测试,可以先在预处理阶段基于统计学方法[15]计算凸包网格的一个OBB进行快速排除,如果测试点不在OBB中,则它也同样位于凸包之外.OBB本身作为表面法向量两两垂直的长方体,能够将其表示为大小一致的AABB绕特定轴旋转的结果,以此,我们定义一个凸包局部坐标空间,原点为OBB中心ocenter,3个坐标轴则由全局基向量经过上述AABB的同等旋转变换得到.

    对于通过OBB测试的质点,需要以凸包自身作进一步的判别.选择OBB的最长轴方向为法向量nproj,过点ocenter定义一个投影平面Pproj,将凸包网格面按照其各自法向量ni与nproj的点积结果分类:ni⋅nproj>0所有面的对应三角形顶点沿-nproj方向投影到Pproj,保持原有的拓扑关系并在Pproj上得到一个平面坐标系下的二维网格,定义为上平面网格.ni⋅nproj<0所有面的对应三角形顶点沿nproj方向投影到Pproj,相应二维网格定义为下平面网格,上、下平面网格具有相同的封闭凸多边形轮廓.ni⋅nproj=0所有面不作任何处理.将测试点位置沿nproj投影,在Pproj上得到的平面投影点设为xproj,如果xproj不包含在上、下平面网格的投影凸多边形轮廓内,那么,测试点位于凸包之外;否则,分别确定包含点xproj的上、下平面网格三角形,这是典型的平面点定位问题.本文使用Kirkpatrick提出的层次剖分方法[16],它的时间复杂度为O(logn)级.由定位的上下平面网格三角形可以索引三维凸包中对应的两个面,最后,测试点在凸包内当且仅当该点同时位于这两个凸包展开平面的负侧.

    2.4 外围体映射加速碰撞响应点计算

    需要对穿透质点作正确的碰撞响应,这里,我们引入响应射线的处理方法.预先确定某个空间源点,由此引一条通过穿透质点当前位置的射线,将其与凸包面距离最近的交点作为碰撞响应点,即质点平移到凸包表面的目标位置.以响应射线遍历凸包所有面进行相交测试无疑会消耗大量的处理时间,为此,本文通过预计算外围映射体对该过程进行了优化.另外,随着映射体类型选取的不同,响应射线的源点定位也有差异.

    映射球.对于形状与球体接近的凸包,计算以ocenter为中心、包含所有顶点的包围球,并构造相应的球面坐标系(r,θ,φ),有θ∈[0,2π)和φ∈[0,π].如图4所示,θ和φ分别在其各自定义域内等间距角度取值对该球体进行经纬剖分,每一个剖分体定义为一个映射单元,设在θ和φ上形成的区间数分别为m和n,则得到一个m×n大小的映射单元表,每个表项包含一个相关三角形链表ltris.对凸包的所有面进行测试,如果某个三角形与映射单元体相交或包含在内,则插入到对应表项的ltris中.

Fig.4 Map sphere and the index table of triangles
图4 映射球和三角形索引表

    给定一个已经变换到局部坐标空间的质点位置x(xlocal,ylocal,zlocal),按下式计算球面坐标参数:

    以ocenter为源点,并以方向定义一条响应射线Rdˆresponse,其中,

    由θx和φx确定映射单元表项,计算Rresponse与表项ltris中某个三角形的交点即为所求碰撞响应点.

    映射圆柱.对于形状较为狭长的凸包,以其OBB长轴为主轴,计算包含所有顶点的有向包围圆柱并构造相应的圆柱坐标系(r,θ,z),有θ∈[0,2π)和z∈[zmin,zmax].如图5所示,对θ和z分别在其各自定义域内等间距角度和距离取值对圆柱体进行剖分,同样得到一个映射单元表,计算每一个表项的相关三角形链表ltris.

Fig.5 Map cylinder and the index table of triangles
图5 映射圆柱及三角形索引表

    计算质点局部位置x(xlocal,ylocal,zlocal)的圆柱坐标参数:

    响应射线源点为x在圆柱坐标系z轴上的投影点,方向如下式:

    由θx和zx确定映射单元表项,碰撞响应点的计算与映射球情况类似.当相交三角形不唯一时,取距离最zprojx近的交点.