自动化学报  2017, Vol. 43 Issue (4): 622-633   PDF    
基于形状骨架图匹配的文物碎片自动重组方法
张雨禾1, 耿国华1, 魏潇然2, 张靖1, 周明全3     
1. 西北大学信息科学与技术学院 西安 710127;
2. 西北大学新闻传播学院 西安 710127;
3. 北京师范大学信息科学与技术学院 北京 100875
摘要: 为了有效解决文物碎片自动重组中由于断裂部位受损造成几何信息丢失,采用传统几何驱动方法容易失效的问题,本文提出一种基于形状骨架图匹配的文物碎片自动重组方法,将碎片匹配问题转化为碎片表面纹饰中非完整纹元的互补匹配问题.首先,通过提取文物碎片表面特征线得到碎片表面的纹饰信息;然后根据完整纹元的特征确定非完整纹元互补匹配的约束条件,采用视觉骨架剪枝的方法提取完全位于断裂部位的非完整纹元的形状骨架图,基于形状骨架图语法及匹配约束条件判定非完整纹元是否互补匹配;接着,将碎片上非完整纹元的顺序作为上层约束条件,采用基于带剪枝深度优先的搜索方法搜索匹配碎片;最后,以邻接碎片上非完整纹元间公共弦的端点作为邻接约束点,采用最小二乘法计算刚体变换参数得到碎片的初始位置,并采用迭代最近点方法将邻接碎片精确对齐.实验结果表明,该方法能够有效解决断裂部位存在缺损文物碎片的自动重组问题.
关键词: 碎片重组     特征线     表面邻接约束     形状匹配    
Reassembly of Fractured Fragments Based on Skeleton Graphs Matching
ZHANG Yu-He1, GENG Guo-Hua1, WEI Xiao-Ran2, ZHANG Jing1, ZHOU Ming-Quan3     
1. School of Information Science and Technology, Northwest University, Xi'an 710127;
2. School of Journalism and Communication, Northwest University, Xi'an 710127;
3. College of Information Science and Technology, Beijing Normal University, Beijing 100875
Received: 2016-04-15, Accepted: 2016-08-02.
Foundation Item: Supported by National Natural Science Foundation of China (61602380, 61373117, 61673319), Scientiflc Research Project of Shaanxi Provincial Department of Education (16JK2178), Special Research Fund for the Doctoral Program of Higher Education (20136101110019)
Author brief: ZHANG Yu-He   Ph. D. candidate at the School of Information Science and Technology, Northwest University. She received her bachelor degree from Northwest University in 2012. Her research interest covers graphics geometry processing, visualization technology and 3D printing technology;
WEI Xiao-Ran   Ph. D. and lecturer at the School of Journalism and Communication, Northwest University. He received his Ph. D. degree from Northwest University in 2016. His research interest covers graphics geometry processing and 3D printing technology;
ZHANG Jing   Master student at the School of Information Science and Technology, Northwestern University. Her research interest covers computer graphics and visualization;
ZHOU Ming-Quan   Professor at the School of Information Science and Technology, Beijing Normal University. His research interest covers virtual reality and visualization technology, intelligent information processing, database and knowledge base, graphics and image processing
Corresponding author. GENG Guo-Hua   Professor at the School of Information Science and Technology, Northwest University. Her research interest covers graphics image processing and visualization technology. Corresponding author of this paper
Abstract: In order to effectively address the problem that traditional geometry-driven methods fail to reassemble the fragments with incompleteness in fracture surfaces or break-curves, an automatic reassembly method based on skeleton graph of shape matching is proposed, which transfers the matching of fragments into complementary matching of incomplete geometry texels. Firstly, the surface ornamentation information is obtained by extracting the feature lines of the fragments. Secondly, the constraints for matching texels are designed according to the features of complete texels, and a visual skeleton pruning method is utilized to generate the skeletons of the shapes of the incomplete texels and the grammar of skeleton graphs is introduced. Then, the sequence of the incomplete texels is used as the upper constraint and a depth-first searching method with pruning for adjacent fragments is presented. Finally, the endpoints of the common chords of the incomplete texels on the adjacent fragments are used as the initial matching points in order to calculate the initial position of the fragment. After that, the iterated closest point (ICP) method is employed to precisely align the adjacent fragments. Experimental results show that the proposed method can successfully solve the reassembly problem of the fragments with incompleteness in fractured surfaces and break-curves.
Key words: Reassembly of fragments     feature lines     surface adjacency constraints     shape matching    

从20世纪90年代起, 伴随信息技术的快速发展, 计算机辅助文物碎片虚拟拼接技术得到了长足的发展, 对文物保护及复原具有重要意义.相较于传统手工复原方法, 数字化虚拟拼接技术能够不受时间和空间的限制, 在提高文物复原效率的同时, 避免人工复原对文物的二次破坏.

文物数字化虚拟拼接技术, 根据复原过程中是否有专家参与可分为自动复原方法和交互式复原方法.其中, 自动复原方法大多以几何特征为驱动, 通过提取破碎部位的几何特征搜索匹配对, 进而实现碎片拼接.

文物碎片主要分为薄壁类碎片和非薄壁类碎片.针对薄壁类文物碎片, 自动复原方法主要基于空间轮廓线的匹配. Üçoluk等[1]将轮廓曲线上顶点的曲率和挠率作为匹配的特征向量, 然后通过特征串匹配的方法搜索匹配对. Cooper等[2]利用贝叶斯方法估计碎片的分布, 并利用断裂处的轮廓线进行碎片匹配和对齐.樊少荣等[3]同样采用基于轮廓线匹配的方法进行碎片拼接, 但其不同在于该方法在区分内轮廓线和外轮廓线的基础上实现碎片拼接, 因此针对非薄壁类文物碎片同样适用. Oxholm等[4]定义的特征串由轮廓曲线上顶点的曲率、挠率和颜色组成, 然后采取最长公共子串实现轮廓线的匹配, 求解匹配对.

非薄壁类碎片的匹配问题是计算机视觉、模式识别和机器智能中一个颇具挑战的问题[5], 针对非薄壁类碎片的匹配, 往往采用基于断裂面几何特征的方法. Winkelbach等[6]采用随机采样的方法, 利用顶点的法向和曲率等特征, 实现断裂面的拼合. Huang等[7]通过计算点云模型上顶点的积分不变量值刻画碎块断裂面的尖锐程度, 进而对断裂面进行分割, 然后基于向前搜索技术和表面一致性的约束方法将断裂面进行两两匹配及融合.随后, Chen等[8]利用曲面中心点、局部曲面的类型及二维直方图刻画特征点处的局部曲面描述子, 通过几何散列算法对碎片进行匹配.王坚等[9]将曲面匹配问题转化为图论中的最大权团搜索问题并采用自旋图确定匹配关系.李姬俊男等[5]提出一种基于空间曲面特征优化的匹配算法, 通过计算模型表面点体积积分不变量形成匹配约束簇, 然后定义3类空间几何一致性约束, 并采用最大独立集方法对非正确匹配对进行消除.

交互式的虚拟文物复原方法则大多需要考古学家和文物修复工作者根据自己的经验来指导复原过程. Mellado等[10]提出了一种基于反馈式的交互复原方法, 通过用户确定碎片间的初始位置, 为后续的拼接过程提供了可信的初始输入值. Palmas等[11]提出了一种融合几何特征及用户经验的交互式复原方法.

现有几何驱动的方法, 无论是基于空间轮廓曲线或基于断裂面的方法, 都要求碎片的断裂部位较完整, 能够提供足够的几何信息用于匹配对搜索.兵马俑碎片由于其陶类质地, 在倒塌时遭到地面撞击, 断裂部位普遍存在缺损, 且由于其体积较大, 破碎时会产生大量碎片.因此, 现有几何驱动的方法往往容易失效.然而, 采用交互式的方法虽然可以从一定程度上解决该问题, 但是依靠人工的方式往往效率较低, 当碎片数量较多时, 即使是具有一定经验的文物修复人员, 也需要耗费大量的时间和精力确定碎片间的邻接关系.

针对上述问题, 本文提出一种基于形状骨架图匹配的文物碎片自动重组方法.该方法的核心思想是将碎片匹配问题转化为碎片表面纹饰中非完整纹元的互补匹配问题:首先, 通过提取碎片的特征线获得碎片表面的纹饰信息; 然后通过分析纹饰中完整纹元的特征定义非完整纹元互补匹配的约束条件, 提取断裂部位非完整纹元的形状骨架图并定义相关语法, 从而判定非完整纹元是否互补匹配; 最后, 将纹饰中互补匹配纹元的顺序作为上层约束, 给出邻接碎片的搜索策略, 最终得到邻接碎片并将其拼合.

1 表面纹饰提取及相关术语

本文将碎片表面纹饰的完整性作为邻接碎片的约束条件, 从而搜索匹配碎片对.因此, 本节首先通过提取碎片的特征线获得碎片表面的纹饰信息.

1.1 表面纹饰提取方法

虽然兵马俑碎片断裂处存在缺损, 但其表面的纹饰信息则普遍较完整, 能够为碎片的重组提供足够的信息, 且大部分文物复原人员也利用该信息进行人工修复.因此, 本文首先提取碎片表面的纹饰信息.

本文算法的处理对象为利用三维激光扫描仪扫描得到的兵马俑碎片点云模型$P=\left \{ p_{i} \right \}\in \boldsymbol{R}^{3}, i=1, 2, \cdots, N$.采用文献[12]中的特征线提取方法提取兵马俑碎片特征线, 从而获得碎片的纹饰信息.由于文献[12]中的特征线提取方法在生成特征线时, 采用了基于$L_{1}$中值的方法生成特征线, 将$L_{1}$中值局部化, 利用式 (1) 重建特征点的形状:

$ \begin{eqnarray}\label{equ1} & Q=\arg \min\sum_{i\in I}{}\sum_{j\in J}{}\left \| x_{i}-c_{j} \right \|\times\\&\theta (\left \| q_{i}-c_{j} \right \|)+ R(X, Q)& \end{eqnarray}\\[-6mm] $ (1)

其中, 第一项为$L_{1}$中值局部化项, 第二项为惩罚项, 防止构成特征线的点距离较近, 且使其呈线性; 最后将形状重建后的特征点进行连接得到特征线.因此, 采用该方法生成的特征线则更准确, 基本步骤如下:

步骤1. 给定点云模型$p$, 利用文献[13]中的方法估计点云上顶点的法向.

步骤2. 采用基于泊松分布的区域生长法, 将点云划分为不同的曲面片, 同时位于多个曲面片相交区域的点判定为潜在特征点, 同时储存所有潜在特征点所属区域信息.

步骤3. 对潜在特征点所属区域信息进行分析, 将位于两个曲面相交区域的点判定为边点, 将位于多个曲面相交区域的点判定为角点, 其余点为噪声点.并将具有相同区域信息的点聚为一类, 然后通过角点的区域信息, 得到特征线的连接性信息.

步骤4. 采用基于$L_1$中值的方法对每个聚类中的点进行形状重建, 得到特征线, 并根据特征线的连接性信息, 连接特征线端点和相应角点, 特征线提取结束.

根据该方法提取得到的兵马俑碎片表面纹饰如图 1所示.

图 1 G10-19号俑部分碎片表面纹饰特征线提取结果 Figure 1 The results of feature lines extraction from fragments of Warriors G10-19
1.2 表面纹饰描述相关术语

为了构造能有效表示碎片断裂处纹饰的形状骨架图及相关语法, 本节对相关术语进行定义.碎片表面完整的矩形块为纹元, 即构成纹饰的最小元素, 如图 1 (b) 所示.每个纹元含有4个顶点, 每个顶点称为纹元的度, 那么, 一个完整纹元的度为4.组成特征线的点在文献[12]中被称为边点, 选取与边点相连接的$k$个点为其局部邻域点, 然后对局部邻域点进行主成分分析 (Principle component analysis, PCA), 若较大的两个特征值$\lambda _{1}, \lambda _{2}$对应的特征向量${\pmb v_{1}}, {\pmb v_{2}}$近似垂直, 则该点为纹元的顶点, 如图 1所示.根据文献[12]中的区域信息分析, 本文将角点定义为断裂点, 即:纹元在断裂处的端点.特征线上其他的非断裂点及非纹元顶点则被称为边缘点.

当两个断裂点之间或一个断裂点与边缘点之间由一条特征线连接且中间不含有其他断裂点时, 本文称这两个点是连通的, 如图 1 (e) 中点$p, q$.连接两个连通点的直线段称为断裂弦, 连接两个连通点的特征线则为断裂线, 如图 1 (e) 所示.

一个非完整纹元上至少含有1个断裂点, 如图 1所示.当一个非完整纹元上含有至少一对连通的断裂点时, 本文称该纹元完全位于碎片的断裂部位, 例如图 1 (e) 所示, 图 1 (d)(f) 中所示纹元则不完全位于断裂部位.

2 匹配纹元对搜索

碎片表面纹饰中纹元的完整性可以为碎片提供邻接约束信息.文物碎片拼接问题是典型的形状匹配问题, 属于局部匹配中的互补匹配[14].因此, 本章首先分析纹元的破碎形状, 然后根据完整纹元的特征定义非完整纹元互补匹配的相关约束; 受文献[15-16]的启发, 本文提取完全处于断裂部位的非完整纹元的形状骨架图, 并定义其语法; 最后, 借助于非完整纹元的形状骨架图, 设计基于形状骨架图的互补纹元匹配判定方法.

2.1 纹元破碎形状分析

完整纹元的度为4, 因此, 纹元破碎形状大致可分为4种: 1) 度为1与度为3的非完整纹元, 如图 2 (a) 所示; 2) 2个度为2的非完整纹元, 如图 2 (b) 所示; 3) 2个度为1的非完整纹元及度为2的非完整纹元, 如图 2 (c)(d) 所示; 4) 4个度为1的非完整纹元, 如图 2 (e)(f) 所示.这4种情况中, 均存在邻接碎片共弦的情况, 即:两个邻接碎片的断裂弦相同且该断裂弦由两个断裂点构成; 但在第3) 种和第4) 种情况中还会出现T型分裂, 即:两个邻接碎片的公共断裂弦由连通的断裂点与边缘点构成, 如图 2 (d)(f) 所示.

图 2 4种纹元破碎情况 Figure 2 Four styles of fractured texels
2.2 邻接约束条件

Palmas等[11]提出了约束的概念, 通过引入约束来描述碎片断裂部位的特征以确定碎片间的空间关系.本节对碎片纹饰中纹元的完整性判定因素进行分析, 并根据这些因素确定非完整纹元互补匹配的约束条件.

1) 完整纹元的度

根据第1.2节中的定义, 一个完整纹元的度应为4, 可定义约束1如下:

约束1.  对任一完整纹元$T$:

$ Degree(T)=4 $ (2)

其中, $Degree(\cdot )$为纹元$T$的度.

2) 匹配度

当非完整纹元拼接构成一个完整纹元时, 根据图 2中所示纹元破碎的4种情况, 本文在搜索匹配块时, 两个匹配块的度之和应满足约束2:

约束2.  对非完整纹元上的断裂弦$T_{1}, T_{2}$:

$ Degree(T_{1})+Degree(T_{2})\leq4 $ (3)

3) 共弦性

邻接的两个非完整纹元必定会有相同的一条断裂弦, 如图 2中黑色虚线所示, 因此, 邻接的两个非完整纹元需满足约束3:

约束3.  对邻接的两个非完整纹元$T_{1}, T_{2}$上的断裂弦$l_{T_{1}}, l_{T_{2}}$:

$ \left | length(l_{T_{1}}) -length(l_{T_{2}})\right |\leq\tau $ (4)

其中, $length(\cdot )$为断裂弦长.断裂弦由连通的两个断裂点或断裂点与边缘点形成. $\tau $度量了两个弦的长度差.

4) 凹凸互补性

文物碎片拼接问题是形状匹配中的互补匹配, 那么邻接的两个非完整纹元的断裂线也应为凹凸互补的.因此, 本文采用断裂线上点的积分不变量进行约束.

设断裂线$C=\partial D$, 计算$C$上的点$p$的面积积分不变量, 其积分不变量[7, 17]定义如式 (5) 所示:

$ A_{r}(p)=\int_{p+r}l_{D}(x)d_{x} $ (5)

其中, 圆$D$的半径为$r$, $A_{r}(p)$为该非完整纹元在曲线段$p\pm r$区域与圆$D$相交部分的面积.那么, 邻接的两个非完整纹元上断裂线积分不变量需满足:

约束4.  对邻接的两个非完整纹元断裂线$C_{i}, C_{j} $上连续的对应点$p_{i}, p_{j}$满足:

$ 0\leq\pi length(r)^{2}-(A_{r}(p_{i})+A_{r}(p_{j}))\leq\gamma $ (6)

其中, $\left | p_{i}-p_{i+1} \right |=r, \left | p_{j}-p_{j+1} \right |=r$.

2.3 形状骨架图及其语法定义

相互匹配的非完整纹元应满足第2.2节中定义的相关约束.然而, 对于T型分裂纹元则很难确定断裂弦及断裂线; 另外, 当纹元的度小于3时, 仅仅根据该纹元的度, 很难确定其属于哪一种破碎情况.

因此, 本节提取完全位于碎片断裂处的非完整纹元的形状骨架图, 根据其潜在形状信息搜索匹配非完整纹元对.

首先, 需要将非完整纹元的特征线进行降维.由于兵马俑碎片局部纹元的弯曲较小, 因此, 本文采用PCA方法对非完整纹元的特征线进行降维, 得到二维纹元形状.

然后, 提取非完整纹元的形状骨架.由于兵马俑碎片是由人工制作的陶制文物碎片, 完全处于断裂部位的非完整纹元的特征线为曲线但并不平滑, 那么, 本文在对匹配对进行判定时, 需要考虑该非完整纹元的全局形状信息, 忽略局部形状凹凸点.因此, 本文采用文献[18]中提出的视觉骨架剪枝的方法, 提取非完整纹元的形状骨架.

最后, 本节对形状骨架图的语法进行定义.根据文献[16]中按照奇点及其邻接奇点最大圆的半径变化趋势将奇点分为4类: 1) 第1类奇点处最大圆的半径单调变化; 2) 第2类奇点处最大圆的半径为局部最小值; 3) 第3类奇点处最大圆的半径不变化; 4) 第4类奇点处最大圆的半径为局部最大值.为了判定相匹配非完整纹元, 本文对形状骨架图语法进行定义, 具体如下:

定义1.  形状骨架图语法 (The skeleton graph grammar, SGG) $G=(V, \Sigma, R, S)$:

1) $V=\left \{ 1, 2, 3, 4\right\}$, 骨架点类型;

2) $\Sigma =\left \{ \Phi, \Psi, \Theta \right\}$, 终结符;

3) $S =\#_{i}$, 起始符;

4) $R=\left \{ R_{1}, R_{2}, R_{3}, R_{4}, R_{5}\right\}$, 语法规则.

其中, $V=\left \{ 1, 2, 3, 4\right\}$分别对应上述4种不同类型的骨架点; $\Phi$表示该骨架的终结点为顶点, $\Psi$表示骨架的终结点为断裂点, $\Theta$表示骨架的终结点为边缘点. $\#_{i}$表示同一碎片上第$i$个非完整纹元的奇点图.语法规则$R=\left \{ R_{1}, R_{2}, R_{3}, R_{4}, R_{5}\right\}$定义如图 3所示.图 1中部分非完整纹元的形状骨架图及语法如图 4所示.

图 3 形状骨架图语法规则 Figure 3 Grammatical rules of skeleton graphs
图 4 部分非完整纹元奇点图 Figure 4 Skeleton graphs of some incomplete texels

得到非完整纹元的形状骨架图后, 可根据终结符$\Phi$的个数得到该纹元的度, 同时根据终结点间的连通性确定断裂弦及断裂线.

2.4 基于形状骨架图的纹元互补匹配判定

本文算法将邻接碎片的匹配问题转化为碎片表面纹饰的匹配问题, 即:当且仅当两个碎片断裂部位的非完整纹元序列一一匹配时, 该碎片为邻接碎片, 可进行重组, 流程如图 5所示.

图 5 邻接碎片搜索流程 Figure 5 Searching pipeline of adjacent fragments

在提取到能够反映碎片断裂处非完整纹元的形状骨架后, 需要根据其形状骨架语法及约束条件设计相应的非完整纹元互补匹配算法, 从而搜索邻接碎片, 完成重组.本节设计非完整纹元互补匹配判定方法.

根据第2.1节中对破碎纹元的分析可知, 纹元度越大, 可匹配纹元的情况越少, 反之, 可匹配纹元的情况随着纹元度的减小而增加.例如:当非完整纹元的度为3时, 仅可能与度为1的非完整纹元匹配; 当非完整纹元的度为2时, 由于存在T型分裂, 因此其可与度为1或2的非完整纹元匹配; 同理, 当非完整纹元的度为1时, 与度为1、2或3的非完整纹元均可匹配.因此, 在搜索匹配纹元对时, 本文以度较大的纹元作为基准, 分别根据第2.2节中的约束条件进行判定, 如无特殊说明, 本节所指的纹元均为完全处于断裂部位的非完整纹元, 判定方法如算法1所示, 具体判定步骤如下:

步骤1. 奇点图的终结符$\Phi$的个数可以有效表示非完整纹元的度, 此时利用约束2进行判定.例如:图 4 (a) 中的纹元可能与图 4 (b), (d), (e) 中的纹元匹配; 图 4 (c) 中的纹元也可能与图 4 (b), (d), (e) 中的纹元匹配; 但图 4 (a) 中的纹元与图 4 (c) 中的纹元则不匹配.

步骤2. 由于纹元匹配为形状匹配中的互补型, 那么纹元在破碎时, 断裂点或边缘点应同时出现在两个邻接的非完整纹元中, 即:构成断裂弦的点类型应相同.因此, 同类型终结符才可进行匹配.那么, 图 4 (a) 中的纹元可能与图 4 (b)(e) 中的纹元匹配; 图 4 (c) 中的纹元仅可能与图 4 (d) 中的纹元匹配.

步骤3. 计算形状骨架图中所有断裂弦的长度及断裂线的积分不变量, 然后利用约束3及约束4对匹配对进行判定, 此时, 图 4 (a) 中的纹元仅可能与图 4 (b) 中的纹元匹配, 图 4 (c) 中的纹元与图 4 (d) 中的纹元匹配.

步骤4. 计算匹配后纹元的度, 若度为4, 则该纹元匹配结束; 若度不为4, 则将拼接后碎片当作一个整体, 生成其形状骨架图, 然后转向步骤1, 直到搜索到匹配纹元.

需要注意的是, 形状骨架图及面积积分不变量虽具有旋转、平移及尺度不变性, 但弦长约束仅具有旋转、平移不变性.因此, 输入的碎片应固定原尺度, 不能进行放大及缩小操作.

算法1.非完整纹元互补匹配判定

输入. 两个非完整纹元$T_{1}, T_{2}$

输出. 匹配纹元对

While(!约束1)   //不满足约束1

{

 if ($Degree(T_{1})+Degree(T_{2})\leq 4$)  //约束2

  {

    if ($\Sigma _{T_{1}}=\Sigma _{T_{2}}$)  //语法图中的终结符相同

     {

      if (约束3  & &  约束4)

      {

        $T_{1}\bigodot T_{2}$   //非完整纹元匹配

      }

     }

  }

}

3 邻接碎片搜索及拼合

部分文物碎片表面纹饰断裂处含有多个纹元, 因此, 基于断裂处非完整纹元的互补匹配, 本节首先给出邻接碎片搜索策略, 然后对邻接碎片进行拼合.

3.1 邻接碎片搜索

上述章节介绍了邻接碎片上非完整纹元匹配约束条件及其判定方法, 但不能保证非完整纹元互补匹配对的唯一性, 非完整纹元间可能出现一对多的匹配.然而, 邻接碎片的非完整纹元间需要同时满足连续性, 形成自我约束, 如图 6所示.因此, 本节定义两个邻接碎片间的约束:

图 6 碎片纹元顺序示意图 Figure 6 Illustration for the sequences of texels on fragments

约束5.  两个碎片$F_{1}, F_{2}$邻接, 当且仅当其上的非完整纹元$T_{i}, T_{j}$满足

$ \left \{ T_{i-1}, T_{i}, \cdots, T_{n}\right\}\bigodot \left \{ T_{j-1}, T_{j}, \cdots, T_{m}\right\} $ (7)

其中, $\bigodot $表示一一匹配.

纹元顺序可通过计算断裂点的邻域得到, 即:

$ T_{i+1}=\left \{p_{T_{i+1}}\subset Neighbor(p_{T_{i}})\right \} $ (8)

其中, $p_{T_i}$$T_i$的断裂点, $Neighbor(p_{T_{i}})$$p_{T_i}$的球邻域, 如图 7所示.当碎片处于俑体中间部位时, 纹元可能会形成环路, 如图 1 (b) 所示.

图 7 纹元顺序计算示意图 Figure 7 Illustration for the computation of the texels sequences

当给定一组文物碎片时, 首先选择一个碎片作为基准搜索匹配碎片, 当该碎片上非完整纹元较少时, 即使满足判定准则, 也可能会出现错误匹配的情况, 反之, 若该碎片上非完整纹元较多时, 会产生大量冗余计算.因此, 本文将碎片$F_i$选作基准碎片, 若其满足:

$ \begin{split} {T_{i-1}, T_i, \cdots, T_n}\subset F_{i}\\ T_{i-1} \neq T_{n}~~~~~~~~~~\\ n\geq 3~~~~~~~~~~~~~~~\\ Degree(T_j)=3~~~~~~~\\ \end{split} $ (9)

即:非完整纹元个数大于等于3、未形成环路, 且度为3的非完整纹元个数最多.

在搜索该基准碎片的邻接匹配碎片时, 采用带剪枝的深度优先搜索策略搜索匹配对.为加快搜索, 在该基准碎片$F_i$上的所有非完整纹元中, 首先选择度最大的纹元$T_i$作为待匹配纹元, 将与其相匹配的最优纹元$T_j$所在碎片$F_j$, 作为邻接匹配对.然后再对邻域纹元$T_{i-1}$$T_{i+1}$$T_{j+1}$$T_{j-1}$进行判定, 若匹配, 则继续判定邻域纹元是否匹配; 否则, 重新搜索与纹元$T_i$匹配的次优匹配纹元$T_t$所在碎片$F_t$.最后非完整纹元匹配数最多的碎片即为最优邻接碎片对.

算法步骤如算法2所示.

算法2.多碎片匹配对搜索

输入. 一组兵马俑铠甲碎片点云模型

输出. 匹配对集合

While (exist (未拼接碎片))

{

  1) 选取满足式 (9) 的碎片$F_i$作为基准碎片

   选取未匹配纹元$T_i\subset F_i$$Degree(T_i)$最大

  2) 搜索非完整纹元$T_i\subset F_i$的匹配纹元$T_j\subset F_j$

  While ($!\{T_i, \cdots, T_t\}\bigodot \{T_j, \cdots, T_n\}$)

   {

    3) 分别判定$T_{i-1} \subset F_i$, $T_{i+1}\subset F_i$

    与$T_{j-1}\subset F_j, T_{j+1}\subset F_j$是否匹配

    if ($\{T_i, \cdots, T_t\}\bigodot \{T_j, \cdots, T_n\}$)

     $i=i+1$$i=i-1$

     $j=j+1$$j=j-1$

    转向3)

    else

    转向2)

   }

}

3.2 碎片拼合

当确定匹配碎片对$F_i, F_j$后, 由于$F_i, F_j$的位置和方向各不相同, 因此需要根据断裂点或边缘点间的对应关系, 估计两碎片的初始位置.

本文在判定匹配对时, 采用了共弦约束 (约束3), 那么, 邻接的两个碎片在对齐时应保证弦对齐.不共线的3个点可确定一个平面, 本文选取至少2个不共线的弦, 由于碎片边缘可能会存在破损, 因此计算刚体变换的弦均应由两个连通的断裂点构成.

假设构成这些弦的断裂点分别构成点集$A, B$, $A$是由$B$中的点经过旋转和平移得到的, 如式 (10) 所示:

$ {\pmb p}_{j}=R {\pmb p}_{i}+{\pmb t} $ (10)

其中, ${\pmb p}_{j}$为属于集合$A$中的点构成的列向量, ${\pmb p}_{i}$为属于集合$B$中的点构成的列向量, $R$为旋转矩阵, ${\pmb t}$为平移向量.

已知弦端点的坐标及对应关系, 采用最小二乘法[19]估计旋转平移矩阵, 使得式 (11) 最小, 然后采用文献[20]中给出的迭代算法进行求解.

$ E=\sum\limits_{i = 1}^N{\left \|{\pmb p}_{j}-(R {\pmb p}_{i}+{\pmb t})\right \|}^2 $ (11)

得到邻接碎片的初始位置后, 采用迭代最近点方法 (Iterated closest point, ICP) [21]对邻接碎片精确对齐.

4 实验结果与分析

本文算法在Intel Core i7 3.33 GHz的CPU、16 GB内存的PC机、Visual studio 2010及PCL环境上进行了实现, 实验数据均为兵马俑1号坑第三次挖掘出土的部分兵马俑碎片点云模型.点云模型由非专业人员利用Creaform VIU手持式扫描仪扫描获得, 扫描分辨率为3.91 mm, 扫描得到的点云未经过任何去噪处理, 点云中均含有实际噪声.碎片在断裂处均存在破损, 断裂处几何特征缺失较严重.

4.1 本文算法结果

图 8 $\sim$ 10所示为本文算法在部分邻接碎片上的实验结果.图 8所示为G10-19号俑中部分碎片的拼接结果, 图 9图 10所示分别为碎片、碎片特征线及其拼接结果.实验结果表明, 本文算法能够有效判定匹配纹元, 从而搜索邻接碎片对, 并将其拼合.

图 8 G10-19号俑部分邻接碎片拼接结果 Figure 8 Reassembly results of some adjacent fragments of Warriors G10-19
图 9 G10-13号俑部分邻接碎片拼接结果 Figure 9 Reassembly results of some adjacent fragments of Warriors G10-13
图 10 G10-52号俑部分邻接碎片拼接结果 Figure 10 Reassembly results of some adjacent fragments of Warriors G10-52

另外, 本文利用文献[12]的方法提取文物碎片表面特征线, 该特征线提取算法经过了噪声点判定及$L_1$几何重建两步, 具有较强的鲁棒性, 如图 9图 10所示的特征线提取结果.由于本文后续步骤均在所提取出的特征线上实现, 因此, 算法对给定的含有噪声的点云模型具有较好的抗噪性.

图 11 $\sim$ 13所示为本文算法在一组整俑碎片上的实验结果, 由于兵马俑铠甲正面与背面的纹饰存在极大的相似性, 因此在判定邻接碎片对时复杂度更高.实验数据如表 1所示.实验结果表明, 本文算法在成组碎片上, 仍然能够获得较好的结果.如图 13中方框区域所示, 当碎片上仅有小部分缺失时, 本文算法仍能获得较好的结果.

图 11 G10-54号俑部分邻接碎片拼接结果 Figure 11 Reassembly results of some adjacent fragments of Warriors G10-54
图 12 G10-39号俑部分邻接碎片拼接结果 Figure 12 Reassembly results of some adjacent fragments of Warriors G10-39
图 13 G10-36号俑部分邻接碎片拼接结果 Figure 13 Reassembly results of some adjacent fragments of Warriors G10-36
表 1 兵马俑碎片模型数据 Table 1 Data of fragment models of warriors

本文将图 1中3号碎片裁剪为薄壁点云后, 特征提取结果如图 14所示, 断裂处的轮廓线被漏提取.这是由于本文算法基于文献[12]中的特征线提取方法, 针对由曲面相交形成的特征具有较好的效果, 而对于薄壁类模型的轮廓点提取效果较差, 造成断裂点漏提取的问题, 从而影响后续非完整纹元互补匹配对的识别等.因此, 本文算法针对薄壁类文物碎片模型, 效果较差.

图 14 薄壁碎片模型特征线提取结果 Figure 14 Feature lines extraction of a shell fragment model
4.2 运行时间

本文算法针对成组兵马俑碎片的重组时间如表 2所示, 其中, 纹饰特征线提取为Stage ${\rm{\#}}$1, 本文在提取特征线时, 每次同时生成4个碎片的特征线; 形状骨架图生成与互补匹配对判定为Stage ${\rm{\#}}$2, 每次同时生成4个非完整纹元的形状骨架图; 碎片拼合为Stage ${\rm{\#}}$3.实验表明, 随着碎片数量的增加, 运行时间也会随之增加, 这是因为纹饰特征线提取为本文算法最耗时的部分, 需要对点云上的所有点进行遍历; 另外, 非完整纹元数量会直接影响匹配对搜索时间.

表 2 本文算法运行时间 (s) Table 2 Execute time of the proposed algorithm (s)
4.3 对比实验

图 15所示分别为采用基于轮廓线方法[3]及基于断裂面方法[5]的碎片重组结果, 实验数据分别为G10-19号俑, G10-54号俑及G10-36号俑部分碎片.其中, G10-19号俑中断裂面存在较严重的缺损, G10-54号俑中碎片轮廓线缺损较严重, G10-36号俑中碎片存在缺失.

图 15 传统算法实验结果 Figure 15 Experimental results of traditional methods

如果两个断裂面能够正确匹配, 它们对齐拼接后不应该发生明显的碰撞[22], 那么当两个碎片精确配准之后, 它们断裂面的拼合基本上达到严丝合缝, 不会发生明显的重叠交错, 即相互渗透.如图 15 (a) 所示为采用文献[3]中基于轮廓线方法对碎片进行重组的结果, 左下图中的线分别为3个碎片断裂面的轮廓线, 可以看出, 当断裂处存在缺损时, 断裂面的轮廓线无法较好地匹配, 因此会出现渗透现象, 如图 15 (a) 中方框标注所示.而采用文献[5]中基于断裂面的拼合方法, 由于断裂面存在缺损, 其几何信息并不完整, 不但会造成渗透现象, 并且断裂处表面的纹饰结构并不连续, 如图 15 (b) 中方框标注所示; 当碎片存在小部分缺失时, 则会导致错误的重组结果, 如图 15 (c) 所示.实验结果表明, 相较于传统基于轮廓线或基于断裂面的方法, 本文方法针对断裂处存在缺损碎片具有更好的效果.

5 总结

本文将文物碎片重组问题转化为碎片表面纹饰中非完整纹元的匹配问题, 采用基于形状骨架图匹配的方法, 在断裂部位几何特征缺失的情况下引入表面纹饰邻接约束, 有效地解决了断裂部位受损碎片的自动重组问题.该方法有效避免了直接利用模型上的几何特征而造成的噪声敏感、采样敏感及计算量过大等问题.

为了能够更准确地描述模型间的匹配关系, 综合考虑多种约束条件来扩展特征描述符将是一个有意义的方向.

参考文献
1 Üçoluk G, Toroslu I H. Automatic reconstruction of broken 3-D surface objects. Computers & Graphics, 1999, 23 (4): 573–582.
2 Cooper D B, Willis A, Andrews S, Baker J, Cao Y, Han D J, Kang K B, Kong W X, Leymarie F F, Orriols X, Velipasalar S, Vote E L, Joukowsky M S, Kimia B B, Laidlaw D H, Mumford D. Bayesian pot-assembly from fragments as problems in perceptual-grouping and geometric-learning. In:Proceedings of the 16th International Conference on Pattern Recognition. Quebec City, Canada:IEEE, 2002. 297-302
3 Fan Shao-Rong, Ru Shao-Feng, Zhou Ming-Quan, Geng Guo-Hua. A method of extracting feature contour from triangular mesh surface model of fractured solid. Journal of Computer-Aided Design & Computer Graphics, 2005, 17 (9): 2003–2009.
( 樊少荣, 茹少峰, 周明全, 耿国华. 破碎刚体三角网格曲面模型的特征轮廓线提取方法. 计算机辅助设计与图形学学报, 2005, 17 (9): 2003–2009. )
4 Oxholm G, Nishino K. Reassembling thin artifacts of unknown geometry. In:Proceedings of the 12th International Conference on Virtual Reality, Archaeology and Cultural Heritage. Aire-la-Ville, Switzerland:Eurographics Association Press, 2011. 49-56
5 Li Ji-Jun-Nan, Geng Guo-Hua, Zhou Ming-Quan, Kang Xin-Yue. Surface feature optimization for virtual matching of relic fragments. Journal of Computer-Aided Design & Computer Graphics, 2014, 26 (12): 2149–2154.
( 李姬俊男, 耿国华, 周明全, 康馨月. 文物碎块虚拟拼接中的表面特征优化. 计算机辅助设计与图形学学报, 2014, 26 (12): 2149–2154. )
6 Winkelbach S, Rilk M, Schönfelder C, Wahl F M. Fast random sample matching of 3D fragments. In:Proceedings of the 26th DAGM Symposium on Pattern Recognition, Lecture Notes in Computer Science. Tübingen, Germany:Springer, 2004. 129-136
7 Huang Q X, Flöry S, Gelfand N, Hofer M, Pottmann H. Reassembling fractured objects by geometric matching. ACM Transactions on Graphics, 2006, 25 (3): 569–578. DOI:10.1145/1141911
8 Chen H, Bhanu B. 3D free-form object recognition in range images using local surface patches. Pattern Recognition Letters, 2007, 28 (10): 1252–1262. DOI:10.1016/j.patrec.2007.02.009
9 Wang Jian, Zhou Lai-Shui. Surface rough matching algorithm based on maximum weight clique. Journal of Computer-Aided Design & Computer Graphics, 2008, 20 (2): 167–173.
( 王坚, 周来水. 基于最大权团的曲面粗匹配算法. 计算机辅助设计与图形学学报, 2008, 20 (2): 167–173. )
10 Mellado N, Reuter P, Schlick C. Semi-automatic geometry-driven reassembly of fractured archeological objects. In:Proceedings of the 11th International Conference on Virtual Reality, Archaeology and Cultural Heritage. Aire-la-Ville, Switzerland:Eurographics Association Press, 2010. 33-38
11 Palmas G, Pietroni N, Cignoni P, Scopigno R. A computer-assisted constraint-based system for assembling fragmented objects. In:Proceedings of the 2013 Digital Heritage International Congress. Marseille, France:IEEE, 2013. 529-536
12 Zhang Y H, Geng G H, Wei X R, Zhang S L, Li S S. A statistical approach for extraction of feature lines from point clouds. Computers & Graphics, 2016, 56 (5): 31–45.
13 Hoppe H, DeRose T, Duchamp T, McDonald J, Stuetzle W. Surface reconstruction from unorganized points. ACM SIGGRAPH Computer Graphics, 1992, 26 (2): 71–78. DOI:10.1145/142920
14 Fan Shao-Rong. Research on Fractured Objects Complementary Shape Matching and Aligning[Ph.D. dissertation], Northwest University, China, 2005
( 樊少荣. 破碎刚体互补形状匹配与拼接方法研究[博士学位论文], 西北大学, 中国, 2005 )
15 Sebastian T B, Klein P N, Kimia B B. Recognition of shapes by editing their shock graphs. IEEE Transactions on Pattern Analysis & Machine Intelligence, 2004, 26 (5): 550–571.
16 Siddiqi K, Shokoufandeh A, Dickenson S J, Zucker S W. Shock graphs and shape matching. In:Proceedings of the 6th International Conference on Computer Vision. Bombay, India:IEEE, 1998. 222-229
17 Wei Xiao-Ran. Extraction of Boundary from Relic Three-Dimensional Model[Master dissertation], Northwest University, China, 2011
( 魏潇然. 三维文物模型边界特征提取[硕士学位论文], 西北大学, 中国, 2011 )
18 Bai X, Latecki L J, Liu W Y. Skeleton pruning by contour partitioning with discrete curve evolution. IEEE Transactions on Pattern Analysis & Machine Intelligence, 2007, 29 (3): 449–462.
19 Arun K S, Huang T S, Blostein S D. Least-squares fitting of two 3-D point sets. IEEE Transactions on Pattern Analysis & Machine Intelligence, 1987, PAMT-9 (5): 698–700.
20 Huang T S, Blostein S D, Margerum E A. Least-squares estimation of motion parameters from 3-D point correspondences. In:Proceedings of the 1986 IEEE Conference Computer Vision and Pattern Recognition. Miami Beach, USA:IEEE, 1986. 112-115
21 Besl P J, McKay H D. A method for registration of 3-D shapes. IEEE Transactions on Pattern Analysis & Machine Intelligence, 1992, 14 (2): 239–256.
22 Li Qun-Hui. Research on Fractured Slid Recovery Based on Bracture Surfaces Matching[Ph.D. dissertation], Northwest University, China, 2013
( 李群辉. 基于断裂面匹配的破碎刚体复原研究[博士学位论文], 西北大学, 中国, 2013 )