CN107423773A - 三维颅骨的自动配准方法和装置 - Google Patents

Abstract

本申请提供了三维颅骨的自动配准方法和装置。本发明中，首先训练一个边缘类型分类器，然后对待配准三维颅骨进行边缘提取获得三维孔洞边缘，对所有三维孔洞边缘提取形状分布特征，采用边缘类型分类器对三维孔洞边缘自动分类为眼眶边缘、鼻框边缘等，然后对鼻框边缘上的所有边缘点进行PCA以此确定统一的坐标系，从而对待配准三维颅骨的初始位置和姿态进行调整，之后采用CPD算法根据待配准三维颅骨与参考颅骨边缘上的对应关系实现边缘区域的配准，最后采用CSRBF对待配准三维颅骨与参考颅骨配准误差大的区域进行调整，实现待配准三维颅骨与参考颅骨整个颅骨区域的精确配准。

Description

三维颅骨的自动配准方法和装置

技术领域

本申请涉及计算机技术，特别涉及三维颅骨的自动配准方法和装置。

背景技术

通过配准建立三维颅骨模型点对点的精确对应是进行颅骨形态统计分析、颅骨面貌复原的前提和基础。但由于颅骨形状和拓扑结构复杂、个体差异大，加之存在牙齿、部分头骨缺失等情况，使得三维颅骨模型的自动配准一直是一个难题和挑战。

目前，三维颅骨模型的自动配准常用的方法是：先采用迭代最近点(Iterative Closest Point，ICP)算法进行三维颅骨模型的粗配准，在此基础上再进行精细配准。

但是，当两个三维颅骨模型初始位置相差较大时，ICP算法容易陷入局部最优从而无法正确配准。更甚，三维颅骨模型拓扑结构复杂，正侧面存在很多孔洞，而现有ICP算法不能实现这些孔洞的精确配准。

发明内容

本申请提供了三维颅骨的自动配准方法和装置，以实现三维颅骨的精确配准。

本申请提供的技术方案包括：

一种三维颅骨的自动配准方法，包括：

训练一个边缘类型分类器；

利用所述边缘类型分类器自动识别待配准三维颅骨上的各条边缘，所述各条边缘至少包括：眼眶边缘、鼻框边缘、颞骨边缘、上颔骨边缘、下颔骨边缘和其他边缘；

通过对识别出的鼻框边缘上的所有边缘点进行主成分分析PCA实现所述待配准三维颅骨的坐标校正，以调整待配准三维颅骨与参考颅骨之间具有近似的位置和姿态；

对待配准三维颅骨与参考颅骨执行边缘区域的配准，之后对待配准三维颅骨与参考颅骨在边缘区域之外的其他区域进行配准。

一种三维颅骨的自动配准装置，包括：

训练模块，用于训练一个边缘类型分类器；

识别模块，用于利用所述边缘类型分类器自动识别待配准三维颅骨上的各条边缘，所述各条边缘至少包括：眼眶边缘、鼻框边缘、颞骨边缘、上颔骨边缘、下颔骨边缘和其他边缘；

校正模块，用于通过对识别出的鼻框边缘上的所有边缘点进行主成分分析PCA实现所述待配准三维颅骨的坐标校正，以调整待配准三维颅骨与参考颅骨之间具有近似的位置和姿态；

配准模块，用于对待配准三维颅骨与参考颅骨执行边缘区域的配准，之后对待配准三维颅骨与参考颅骨在边缘区域之外的其他区域进行配准。

由以上技术方案可以看出，本发明中，首先训练一个边缘类型分类器，然后对待配准三维颅骨进行边缘提取获得三维孔洞边缘，对所有三维孔洞边缘提取形状分布特征，采用边缘类型分类器对三维孔洞边缘自动分类为眼眶边缘、鼻框边缘等，然后对鼻框边缘上的所有边缘点进行PCA以此确定统一的坐标系从而对待配准三维颅骨的初始位置和姿态进行调整，之后采用CPD算法进行待配准三维颅骨与参考颅骨边缘区域的配准，最后采用CSRBF对待配准三维颅骨与参考颅骨未匹配好的区域进行配准，实现待配准三维颅骨与参考颅骨整个颅骨区域的精确配准。

附图说明

图1为本发明提供的方法流程图；

图2为本发明提供的步骤101实现流程图；

图3为本发明提供的各个类型的边缘示意图；

图4为本发明提供的步骤103实现流程图；

图5a为本发明提供的坐标校正前待配准三维颅骨和参考颅骨位置关系的示意图；

图5b为本发明提供的坐标校正后待配准三维颅骨和参考颅骨位置关系的示意图；

图6为本发明提供的采用CPD算法实现边缘区域配准的流程图；

图7a为本发明提供的用于CPD配准的两个点集的示意图；

图7b为本发明提供的经CPD配准后的两个点集的示意图；

图7c为本发明提供的经CPD配准后的参考颅骨和待配准三维颅骨的配准示意图；

图8a为本发明提供的经过一次CSRBF后参考颅骨和待配准三维颅骨的配准示意图；

图8b为本发明提供的经过五次CSRBF后参考颅骨和待配准三维颅骨的配准示意图；

图9为本发明提供的装置结构图。

具体实施方式

为了使本发明的目的、技术方案和优点更加清楚，下面结合附图和具体实施例对本发明进行详细描述。

参见图1，图1为本发明提供的方法流程图。如图1所示，该流程可包括以下步骤：

步骤101，训练一个边缘类型分类器。

步骤102，利用所述边缘类型分类器自动识别待配准三维颅骨上的各条边缘，各条边缘至少包括：眼眶边缘、鼻框边缘、颞骨边缘、上颔骨边缘、下颔骨边缘和其他边缘。

步骤103，通过对识别出的鼻框边缘上的所有边缘点进行主成分分析(Principle Components Analysis，PCA)实现待配准三维颅骨的坐标校正，以调整待配准三维颅骨与参考颅骨之间具有近似的位置和姿态。

步骤104，对待配准三维颅骨与参考颅骨执行边缘区域的配准，之后对待配准三维颅骨与参考颅骨在边缘区域之外的其他区域进行配准。

通过步骤101至步骤104，最终实现了待配准三维颅骨与参考颅骨精确配准。

下面对图1所示的各个步骤进行具体描述：

参见图2，图2为本发明提供的步骤101实现流程图。在执行图2所示流程之前，预先会设置一个数据库，该数据库中包含多个待训练的三维颅骨模型。

如图2所示，该流程可包括以下步骤：

步骤201，对数据库中的每一个三维颅骨模型进行边缘提取，以提取出边界边、三维孔洞边缘、边缘点。

这里，待训练的三维颅骨模型都是采用三角形或四边形网格模型表示。当待训练的三维颅骨模型采用三角形网格模型表示时，边界边就是只被一个三角形网格使用的边，当待训练的三维颅骨模型采用四边形网格模型表示时，边界边就是只被一个四边形网格使用的边。

这里，三维孔洞边缘则是由多个不同的边界边首尾相连组成的，在实际应用中，从一个三维颅骨模型中提取出的三维孔洞边缘不止一个。边缘点是指边界边两端的端点。

步骤202，标注每一个待训练三维颅骨模型中的每一条三维孔洞边缘的类型。

作为本发明的一个实施例，步骤202中，可手工标注。

步骤203，从每一三维孔洞边缘提取形状分布特征。

鉴于形状分布算法简单，并具有很强的不变性和鲁棒性，作为一个优选实施例，本步骤203可具体采用D2形状分布函数从每一三维孔洞边缘提取形状分布特征，具体为：

步骤a1，采用统一的距离范围，并把距离范围平均分为N段；

步骤a2，计算三维孔洞边缘中每两个边缘点之间的距离，确定该计算出的距离落在哪一段距离范围内，则该距离范围的数量增加1；

步骤a3，对各段距离范围内的数量进行归一化，获得统计直方图。即完成了从三维孔洞边缘提取形状分布特征的操作。通过步骤a1至步骤a3，能够保证最终每一三维孔洞边缘可采用一个N维的向量表示。

步骤204，对所有三维孔洞边缘，基于提取出的形状分布特征以及标注的边缘类型，采用支持向量机(Support Vector Machines，SVM)对所有三维孔洞边缘进行分类训练得到边缘类型分类器，三维孔洞边缘至少分为以下6类：眼眶边缘、鼻框边缘、颞骨边缘、上颔骨边缘、下颔骨边缘以及其他边缘。

在本发明中，初始默认所有三维孔洞边缘类型为其他边缘，然后手工标注其中的眼眶边缘、鼻框边缘、颞骨边缘、上颌骨边缘和下颌骨边缘，并计算所有边缘的形状分布特征，从而获得了一个边缘类型分类训练样本集合。基于该训练样本集合，采用SVM得到一个边缘类型分类器。这里，以数据库中包含331个待训练的三维颅骨模型为例，则在实验中，本申请发明人通过这331个待训练的三维颅骨模型共5673条三维孔洞边缘的形状分布特征和对应的边缘类型采用SVM进行训练，从而获得了一个比较完整、健全的边缘类型分类器。

至此，完成图2所示的流程。

在步骤102中，当输入待配准三维颅骨时，首先按照如上述步骤201的方法从待配准三维颅骨中自动提取所有边界边、三维孔洞边缘、边缘点，之后按照步骤203的方式对每一条三维孔洞边缘计算形状分布特征，最后根据各条三维孔洞边缘提取出的形状分布特征并利用训练得到的边缘类型分类器就会自动识别各条三维孔洞边缘的边缘类型。这里，边缘类型至少为：眼眶边缘、鼻框边缘、颞骨边缘、上颔骨边缘、下颔骨边缘、以及其他边缘。图3示出了各个类型的边缘。在下文中，一般不涉及其他边缘，所以其他边缘不做考虑。

下面对步骤103进行描述：

由于三维颅骨可通过不同的设备或软件获得，因此不同三维颅骨的初始位置和姿态可能存在较大差异，为了消除这样的差异，在本发明中，所有的三维颅骨均采用统一的坐标系，基于此，本申请就需要对待配准三维颅骨进行坐标校正，以满足该统一的坐标系，具体可参见图4所示的流程：

参见图4，图4为本发明提供的步骤103实现流程图。如图4所示，该流程可包括以下步骤：

步骤401，对待配准三维颅骨的鼻框边缘上的所有边缘点坐标组成的矩阵进行PCA，计算前三个最大特征值对应的特征向量(分别记作u1，u2，u3)。

步骤402，确定待配准三维颅骨的坐标系：X轴的正方向为u1的反方向，Z轴的正方向为u2的方向，Y轴的正方向为u3的方向，坐标原点为鼻框边缘中心点。

这里，鼻框边缘中心点的坐标轴值(x、y、z)通过以下方式确定：

鼻框边缘中心点的X坐标轴值x为鼻框边缘上所有边缘点的X坐标轴值的平均值；

鼻框边缘中心点的Y坐标轴值y为鼻框边缘上所有边缘点的Y坐标轴值的平均值；

鼻框边缘中心点的Z坐标轴值z为鼻框边缘上所有边缘点的Z坐标轴值的平均值。

也就是说，鼻框边缘中心点实质是鼻框边缘上所有边缘点的中心点。

上述步骤402之所以如此定义坐标系，原因是尽可能和人类学中常用的法兰克福坐标系相近，当然，本领域其他技术人员也可采用其他定义坐标系的方式，本发明并不再一一具体限定。

需要说明的是，由于u₁、u₂和u₃的正负方向和三维颅骨所在的坐标系相关，因此由它们确定的X轴、Y轴、Z轴可能与参考颅骨的X轴、Y轴、Z轴正负方向不一致，在本发明中，需要对正负方向不一致的坐标轴进行调整，从而保证待配准三维颅骨经坐标校正后和参考颅骨的姿态和位置相似。比如，当发现已确定的X轴的正负方向与参考颅骨的X轴的正负方向不一致，则调整已确定的X轴的正负方向，当发现已确定的Y轴的正负方向与参考颅骨的Y轴的正负方向不一致，则调整已确定的Y轴的正负方向，当发现已确定的Z轴的正负方向与参考颅骨的Z轴的正负方向不一致，则调整已确定的Z轴的正负方向。

具体的调整过程为：

步骤b1，计算矢量v＝p_nosecenter-p_skullcenter，其中p_nosecenter为鼻框边缘中心点，p_skullcenter为三维颅骨中心点；

步骤b2，判断v是否与已确定的Y轴的夹角大于90°即判断v与已确定的Y轴方向的点积是否小于0，如果是，则Y轴反向，即将已确定的Y轴正方向反向旋转180°，这样保证Y轴正向是朝外的；

步骤b3，判断以下条件：下颌骨边缘中心点的Z轴坐标值z1值<上颌骨边缘中心点的Z轴坐标值z2值<鼻框边缘中心点的Z轴坐标值z3值<眼眶边缘中心点的Z轴坐标值z4值是否满足，如果不满足，则Z轴和X轴反向，即将已确定的Z轴、X轴的正方向反向旋转180°。

通过上述定义的坐标系即实现了待配准三维颅骨的坐标校正，以调整待配准三维颅骨与参考颅骨之间具有近似的位置和姿态。图5a和5b分别显示了坐标校正前后待配准三维颅骨和参考颅骨之间的位置关系，其中深色显示的是参考颅骨，浅色显示的是待配准三维颅骨。从图5a和5b可知，即使待配准三维颅骨和参考颅骨两个颅骨的初始位置和姿态相差很大，通过坐标校正后，可得到相近的姿态和位置,从而实现了待配准三维颅骨和参考颅骨的预配准。

在本发明中，待配准三维颅骨经过坐标校正后，可进一步根据鼻框边缘与眼眶边缘位置之间的相互关系，确定出左右方位，具体为：比较眼眶边缘中心点的X轴坐标值与鼻框边缘中心点的X轴坐标值，如果前者大于后者，则眼眶边缘为右眼眶边缘，反之，眼眶边缘为左眼眶边缘；

采用同样的方法也可确定颞骨边缘和上颌骨边缘的左右方位，具体为：

比较颞骨边缘中心点的X轴坐标值与鼻框边缘中心点的X轴坐标值，如果前者大于后者，则颞骨边缘为右颞骨边缘，反之，颞骨边缘为左颞骨边缘；

比较上颔骨边缘中心点的X轴坐标值与鼻框边缘中心点的X轴坐标值，如果前者大于后者，则上颔骨边缘为右边的上颔骨边缘，反之，上颔骨边缘为左边的上颔骨边缘。

下面描述步骤104中对待配准三维颅骨与参考颅骨执行边缘区域的配准(简称边缘区域配准)的具体实现：

在本发明中，采用一致点漂移(Coherent Point Drift，CPD)算法实现上述边缘区域配准。CPD算法是一种基于概率的点集非刚性配准算法，被视为基于速度场运动一致性约束的最大似然估计问题，其通过变分的方法表达运动一致性约束并得到规则化最大似然估计。在CPD算法中，一个点集作为高斯混合模型的内核，另一个点集作为高斯混合模型的数据，采用确定性退火的EM算法优化求解最大似然估计过程的参数，从而找到两个点集的对应关系与非刚性变换，即内核点集向数据点集配准。这种方法可以估计复杂的非线性非刚性变换，而且在含有噪声和溢出点的情况下有较强的鲁棒性。因此可用于颅骨有缺损的情况。

图6为本发明提供的采用CPD算法实现边缘区域配准的流程图。如图6所示，该流程可包括：

步骤601，将待配准三维颅骨上除其他边缘之外的边缘上的所有边缘点组成数据点集，将参考颅骨上满足条件的边缘上的顶点组成内核点集；所述条件为：类型与待配准三维颅骨上边缘的类型相同。

在本发明采用的CPD算法进行边缘区域配准中，两个点集由待配准三维颅骨和参考颅骨各自的边缘点构成，但并不是所有的边缘点。由于待配准三维颅骨可能存在缺损，因此根据待配准三维颅骨中包含的边缘对应地在参考颅骨上选择出相同类型的边缘，将待配准三维颅骨上各类边缘上的顶点组成数据点集，而将参考颅骨上选择出的边缘上的顶点组成内核点集。

为避免距离边缘区域较远的颅顶在配准中产生错误的变形，优选地，作为本发明的一个实施例，本发明中，在参考颅骨顶部手工标定了一些顶点，对于其中每个顶点，求其在待配准颅骨上的最近点，然后把这些最近点添加到数据点集中。而对数据点集中每一顶点，计算其在参考颅骨上的最近点，将计算出的最近点添加至所述内核点集。

由于待配准颅骨可能存在缺损，因此，优选地，作为本发明的一个实施例，可对数据点集中新增加的每一个顶点，依此求其在参考颅骨上的最近点，之后把这些最近点增加到内核点集中，以尽可能避免上述的缺陷。

步骤602，基于数据点集和内核点集并采用CPD算法计算非刚性变换，将该变换作用至参考颅骨的各个顶点以实现经过CPD变换后的参考颅骨与待配准三维颅骨在相同类型的边缘区域配准。

这里，CPD算法已是比较成熟的算法，采用CPD算法进行非刚性变换计算的方式是CPD算法的一种实现，比较成熟，不再赘述。

至此，通过步骤601至步骤603实现了采用CPD算法进行边缘区域的配准。图7a示出了用于CPD配准的两个点集，图7b示出了的是经CPD配准后的两个点集；图7a、7b中黑色显示的为内核点集，灰色显示的是数据点集；图7c示出了CPD配准后的待配准三维颅骨与参考颅骨，在图7c中，深色显示的是参考颅骨，浅色显示的是待配准三维颅骨。从图7c可以看出经过CPD变形后的参考颅骨在共有边缘处与待配准颅骨能够很好地匹配上，但在距离边缘较远的区域，则未能完全匹配上。

下面对步骤104中对待配准三维颅骨与参考颅骨在边缘区域之外的其他区域进行配准(简称全局配准)进行描述：

为了进一步对待配准三维颅骨与参考颅骨中所有颅骨区域进行精确配准，本发明采用具有紧支撑的径向基函数(Compact Support Radial BasisFunctions,CSRBF)实现上述的全局配准。

优选地，这里的径向基函数可定义为Wendland CSRBF，简称Wendland函数(通过ψ_d,k(r)表示)。Wendland函数，其要求每一特征点在三维空间中的左右域为一个半径可调的球体，并只对待配准三维颅骨与参考颅骨之间未匹配的区域进行匹配，维持已经配准的区域。

当给定空间维数d，平滑度C^2k(R)以及欧式距离r时，ψ_d,k(r)可表示为：

其中，

为截断多项式，

为积分运算，在公式1中执行k次。

从公式1可以看出，ψ_d,k(r)只有在r≤1时有效，有效范围缩放至a，a在下文进行了定义，这里不再赘述。在本发明中，缩放后的ψ_d,k(r)的数学属性保持不变，即：

ψ_a(r)＝ψ(r/a)。

对于3维空间，k＝0,1,2的情况，ψ_d,k(r)分别表示如下：

在本发明中，采用ψ_3,1(r)作为CSRBF中的径向基函数，则采用CSRBF对经过边缘区域配准后的参考颅骨上的顶点进行变形具体为：

利用公式2对经过边缘区域配准后的参考颅骨上的顶点进行变形：

其中为一个顶点，为顶点到特征点p_i＝(xi,yi,zi)之间的欧式距离，m为特征点的个数，α_i(i＝1,...m)为未知权重，这些权重可通过以下公式2a把参考颅面模型上的每一个特征点p_i一一映射到待配准颅面模型的对应特征点q_i求解得到：

u(p_i)＝q_i,i＝1,2...m 公式2a

公式2a表示成矩阵形式，则如公式3所示：

Kα＝ΔQ 公式3

其中K为m×m的矩阵，矩阵中的元素k_ij＝ψ_a,3,1||p_i-q_i||，α＝(α₁,α₂...α_m)^T，ΔQ＝(q₁-p₁,q₂-p₂...q_m-p_m)^T。

权重α_i(i＝1,...m)的值一旦确定后，就可以根据公式2对经过边缘区域配准后的参考颅骨上的顶点进行变形，使它们更贴切地匹配到待配准三维颅骨上。在本发明中，CSRBF可执行多次，直到执行次数超过用户设定的阈值，或最大配准误差小于用户设定的另一阈值。

在本发明中，特征点{p_i|i＝1,2,...,m}和{q_i|i＝1,2,...,m}是自动生成的，它们分别分布在参考颅骨和待配准三维颅骨配准误差大的区域上，具体通过以下方式生成：

因为待配准三维颅骨可能存在缺失，因此在求特征点时，首先对待配准三维颅骨上的每一个顶点，求解它在参考颅骨上的最近点作为对应点，三维颅骨上的每一个顶点与其在参考颅骨上的对应点组成一个点对，将各点对中两个点之间的距离看作配准误差，然后对各点对按照配准误差从大到小的顺序进行排序，之后从排序后的第一个点对开始采样，就能得到特征点{p_i|i＝1,2,...,m}和{q_i|i＝1,2,...,m}。并且，为了使特征点分布更合理，在对排序后的点对进行采样时，需满足如下公式4示出的距离限制条件：

||p_i-p_j||＞0.5a,p_i≠p_j,p_j∈P 公式4

其中a为ψ_a,3,1(r)的支持范围，a＝3.66Δ，其中Δ为所有对应点对在X,Y,Z方向轴上的最大位移；P＝{p_i|i＝1,2,...,m}，为参考颅骨上的CSRBF特征点集合。

采用公式4这一限制，能够保证特征点P＝{p_i|i＝1,2,...,m}不会聚集在某一处，对应地，{q_i|i＝1,2,...,m}也不会聚集在某一处。

图8a和8b分别示出了对图7c中两个颅骨进一步进行1次CSRBF和5次CSRBF配准后的结果。在图8a和8b中，深色显示的是参考颅骨，浅色显示的是待配准三维颅骨。可以看出，CSRBF能保持边缘区域的配准结果，同时对配准误差大的区域的配准效果进行调整，使得最终两个颅骨很好地匹配到了一起。

以上对本发明提供的方法进行了描述。

下面对本发明提供的装置进行描述：

参见图9，图9为本发明提供的装置结构图。如图9所示，该装置可包括：

训练模块，用于训练一个边缘类型分类器；

识别模块，用于利用所述边缘类型分类器自动识别待配准三维颅骨上的各条边缘，所述各条边缘至少包括：眼眶边缘、鼻框边缘、颞骨边缘、上颔骨边缘、下颔骨边缘和其他边缘；

校正模块，用于通过对识别出的鼻框边缘上的所有边缘点进行主成分分析PCA实现所述待配准三维颅骨的坐标校正，以调整待配准三维颅骨与参考颅骨之间具有近似的位置和姿态；

配准模块，用于对待配准三维颅骨与参考颅骨执行边缘区域的配准，之后对待配准三维颅骨与参考颅骨在边缘区域之外的其他区域进行配准。

优选地，在本发明中，训练模块训练一个边缘类型分类器包括：

对数据库中的每一个三维颅骨模型进行边缘提取，以提取出边界边、三维孔洞边缘、边缘点；其中，数据库中的每一个三维颅骨模型是采用三角形或四边形网格模型表示，边界边是指三维颅骨模型上只被一个三角形或四边形网格使用的边，三维孔洞边缘是由不同边界边首尾相连组成的，边缘点是指边界边两端的端点；

标注每一个三维颅骨模型中每一条三维孔洞边缘的类型；

从每一三维孔洞边缘提取形状分布特征；

对所有三维孔洞边缘，基于提取出的形状分布特征和标注的边缘类型并采用支持向量机SVM进行分类训练得到边缘类型分类器，三维孔洞边缘至少分为以下6类：眼眶边缘、鼻框边缘、颞骨边缘、上颔骨边缘、下颔骨边缘以及其他边缘。

优选地，在本发明中，识别模块通过以下步骤利用边缘类型分类器自动识别待配准三维颅骨上的各条边缘：

对待配准三维颅骨进行边缘提取，以提取出边界边、三维孔洞边缘、边缘点；其中，待配准三维颅骨是采用三角形或四边形网格模型表示，边界边是指待配准三维颅骨上只被一个三角形或四边形网格使用的边，三维孔洞边缘是由不同边界边首尾相连组成的，边缘点是指边界边上的两个端点；

从提取出的每一三维孔洞边缘提取形状分布特征；

针对每一三维孔洞边缘，根据从该三维孔洞边缘提取出的形状分布特征并利用所述边缘类型分类器自动识别该三维孔洞边缘的边缘类型，边缘类型至少为：眼眶边缘、鼻框边缘、颞骨边缘、上颔骨边缘、下颔骨边缘和其他边缘。

优选地，在本发明中，校正模块通过对识别出的鼻框边缘上的所有边缘点进行主成分分析PCA实现所述待配准三维颅骨的坐标校正包括：

对鼻框边缘上各个边缘点坐标组成的矩阵进行主成分分析PCA，计算前三个最大特征值对应的特征向量(分别记作u1，u2，u3)；

确定待配准三维颅骨的坐标系：X轴的正方向为u1的反方向，Z轴的正方向为u2的方向，Y轴的正方向为u3的方向，坐标原点为鼻框边缘中心点。

优选地，在本发明中，当待配准三维颅骨的坐标系中X轴、Y轴、Z轴的正负方向与参考颅骨的坐标系中X轴、Y轴、Z轴的正负方向不一致时，校正模块进一步执行以下步骤：

计算矢量v＝p_nosecenter-p_skull_center，其中，p_nosecenter为鼻框边缘中心点，p_skullcenter为颅骨边缘中心点；

判断v与已确定的Y轴的夹角是否大于90°，如果是，将已确定的Y轴正方向反向旋转180°；

判断以下条件：下颌骨边缘中心点的Z轴坐标值z1值<上颌骨边缘中心点的Z轴坐标值z2值<鼻框边缘中心点的Z轴坐标值z3值<眼眶边缘中心点的Z轴坐标值z4值是否满足，如果否，则将已确定的Z轴、X轴的正方向反向旋转180°。

优选地，在本发明中，校正模块进一步执行以下步骤：

比较眼眶边缘中心点的X轴坐标值与鼻框边缘中心点的X轴坐标值，如果前者大于后者，则眼眶边缘为右眼眶边缘，反之，眼眶边缘为左眼眶边缘；

比较颞骨边缘中心点的X轴坐标值与鼻框边缘中心点的X轴坐标值，如果前者大于后者，则颞骨边缘为右颞骨边缘，反之，颞骨边缘为左颞骨边缘；

比较上颔骨边缘中心点的X轴坐标值与鼻框边缘中心点的X轴坐标值，如果前者大于后者，则上颔骨边缘为右边的上颔骨边缘，反之，上颔骨边缘为左边的上颔骨边缘。

优选地，在本发明中，配准模块对待配准三维颅骨与参考颅骨执行边缘区域的配准包括：

将待配准三维颅骨上除其他边缘之外的边缘上的所有边缘点组成数据点集；

将参考颅骨上满足条件的边缘上的所有边缘点组成内核点集；所述条件为：类型与待配准三维颅骨上组成数据点集的边缘的类型相同；

基于所述数据点集和所述内核点集采用一致点漂移CPD算法计算非刚性变换，将该变换作用至参考颅骨的各个顶点以实现经过CPD变换后的参考颅骨与待配准三维颅骨在相同类型的边缘区域配准。

优选地，本发明中，所述配准模块对待配准三维颅骨与参考颅骨在边缘区域之外的其他区域进行配准包括：

采用具有紧支撑的径向基函数CSRBF对经过边缘区域配准后的参考颅骨上的顶点进行变形，以使变形后的参考颅骨匹配到待配准三维颅骨；

所述径向基函数CSRBF定义为Wendland函数，由ψ_d,k(r)表示，以下公式1示出了ψ_d,k(r)；

为截断多项式，

为积分运算，在公式1中执行k次；

ψ_d,k(r)在r≤1时有效，有效范围缩放至a，ψ_a(r)＝ψ(r/a)；对于3维空间，和k＝0,1,2的情况，ψ_d,k(r)分别表示如下：

其中，当径向基函数为ψ_3,1(r)时，采用CSRBF对经过边缘区域配准后的参考颅骨上的顶点进行变形包括：

利用公式2对经过边缘区域配准后的参考颅骨上的顶点进行变形：

其中，为参考颅骨上的一个顶点，为顶点到特征点p_i＝(xi,yi,zi)之间的欧式距离，m为特征点的个数，α_i(i＝1,...m)为未知权重，通过把参考颅骨上的每一个特征点p_i一一映射到待配准三维颅骨上的对应特征点q_i求解得到，表示成矩阵形式如公式3所示：

Kα＝ΔQ 公式3

其中K为m×m的矩阵，矩阵中的元素k_ij＝ψ_a,3,1||p_i-q_i||，α＝(α₁,α₂...α_m)^T，ΔQ＝(q₁-p₁,q₂-p₂...q_m-p_m)^T；

特征点{p_i|i＝1,2,...,m}和{q_i|i＝1,2,...,m}是自动生成的，分别分布在参考颅骨和待配准三维颅骨配准误差大的区域上，具体通过以下方式生成：首先对待配准三维颅骨上的每一个顶点，求解它在参考颅骨上的最近点作为对应点，三维颅骨上的每一个顶点与其在参考颅骨上的对应点组成一个点对，将各点对中两个点之间的距离看作配准误差，然后对各点对按照配准误差从大到小的顺序进行排序，之后从排序后的第一个点对开始采样，就能得到特征点{p_i|i＝1,2,...,m}和{q_i|i＝1,2,...,m}，并且，在对排序后的点对进行采样时，需满足如下公式4示出的距离限制条件：

||p_i-p_j||＞0.5a,p_i≠p_j,p_j∈P 公式4

其中a为ψ_a,3,1(r)的支持范围，a＝3.66Δ，其中Δ为所有对应点对在X,Y,Z方向轴上的最大位移；P＝{p_i|i＝1,2,...,m}，为参考颅骨上的特征点集合。

以上对本发明提供的装置进行了描述。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明保护的范围之内。

Claims (10)

1.一种三维颅骨的自动配准方法，其特征在于，该方法包括：

训练一个边缘类型分类器；

利用所述边缘类型分类器自动识别待配准三维颅骨上的各条边缘，所述各条边缘至少包括：眼眶边缘、鼻框边缘、颞骨边缘、上颔骨边缘、下颔骨边缘和其他边缘；

通过对识别出的鼻框边缘上的所有边缘点进行主成分分析PCA实现所述待配准三维颅骨的坐标校正，以调整待配准三维颅骨与参考颅骨之间具有近似的位置和姿态；

对待配准三维颅骨与参考颅骨执行边缘区域的配准，之后对待配准三维颅骨与参考颅骨在边缘区域之外的其他区域进行配准。

2.根据权利要求1所述的方法，其特征在于，所述训练一个边缘类型分类器包括：

对数据库中的每一个三维颅骨模型进行边缘提取，以提取出边界边、三维孔洞边缘、边缘点；其中，数据库中的每一个三维颅骨模型是采用三角形或四边形网格模型表示，边界边是指三维颅骨模型上只被一个三角形或四边形网格使用的边，三维孔洞边缘是由不同边界边首尾相连组成的，边缘点是指边界边两端的端点；

标注每一个三维颅骨模型中每一条三维孔洞边缘的类型；

从每一三维孔洞边缘提取形状分布特征；

对所有三维孔洞边缘，基于提取出的形状分布特征和标注的边缘类型并采用支持向量机SVM进行分类训练得到边缘类型分类器，三维孔洞边缘至少分为以下6类：眼眶边缘、鼻框边缘、颞骨边缘、上颔骨边缘、下颔骨边缘以及其他边缘。

3.根据权利要求1所述的方法，其特征在于，所述利用边缘类型分类器自动识别待配准三维颅骨上的各条边缘包括：

对待配准三维颅骨进行边缘提取，以提取出边界边、三维孔洞边缘、边缘点；其中，待配准三维颅骨是采用三角形或四边形网格模型表示，边界边是指待配准三维颅骨上只被一个三角形或四边形网格使用的边，三维孔洞边缘是由不同边界边首尾相连组成的，边缘点是指边界边上的两个端点；

从提取出的每一三维孔洞边缘提取形状分布特征；

针对每一三维孔洞边缘，根据从该三维孔洞边缘提取出的形状分布特征并利用所述边缘类型分类器自动识别该三维孔洞边缘的边缘类型，边缘类型至少为：眼眶边缘、鼻框边缘、颞骨边缘、上颔骨边缘、下颔骨边缘和其他边缘。

4.根据权利要求1所述的方法，其特征在于，所述通过对识别出的鼻框边缘上的所有边缘点进行主成分分析PCA实现所述待配准三维颅骨的坐标校正包括：

对鼻框边缘上各个边缘点坐标组成的矩阵进行主成分分析PCA，计算前三个最大特征值对应的特征向量(分别记作u1，u2，u3)；

确定待配准三维颅骨的坐标系：X轴的正方向为u1的反方向，Z轴的正方向为u2的方向，Y轴的正方向为u3的方向，坐标原点为鼻框边缘中心点。

5.根据权利要求4所述的方法，其特征在于，当待配准三维颅骨的坐标系中X轴、Y轴、Z轴的正负方向与参考颅骨的坐标系中X轴、Y轴、Z轴的正负方向不一致时，该方法进一步包括：

计算矢量v＝p_nosecenter-p_skullcenter，其中，p_nosecenter为鼻框边缘中心点，p_skullcenter为三维颅骨中心点；

判断v与已确定的Y轴的夹角是否大于90°，如果是，将已确定的Y轴正方向反向旋转180°；

判断以下条件：下颌骨边缘中心点的Z轴坐标值z1值<上颌骨边缘中心点的Z轴坐标值z2值<鼻框边缘中心点的Z轴坐标值z3值<眼眶边缘中心点的Z轴坐标值z4值是否满足，如果否，则将已确定的Z轴、X轴的正方向反向旋转180°。

6.根据权利要求4所述的方法，其特征在于，该方法进一步包括：

比较眼眶边缘中心点的X轴坐标值与鼻框边缘中心点的X轴坐标值，如果前者大于后者，则眼眶边缘为右眼眶边缘，反之，眼眶边缘为左眼眶边缘；

比较颞骨边缘中心点的X轴坐标值与鼻框边缘中心点的X轴坐标值，如果前者大于后者，则颞骨边缘为右颞骨边缘，反之，颞骨边缘为左颞骨边缘；

比较上颔骨边缘中心点的X轴坐标值与鼻框边缘中心点的X轴坐标值，如果前者大于后者，则上颔骨边缘为右边的上颔骨边缘，反之，上颔骨边缘为左边的上颔骨边缘。

7.根据权利要求1所述的方法，其特征在于，对待配准三维颅骨与参考颅骨执行边缘区域的配准包括：

将待配准三维颅骨上除其他边缘之外的边缘上的所有边缘点组成数据点集；

将参考颅骨上满足条件的边缘上的所有边缘点组成内核点集；所述条件为：类型与待配准三维颅骨上组成数据点集的边缘的类型相同；

基于所述数据点集和所述内核点集采用一致点漂移CPD算法计算非刚性变换，将该变换作用至参考颅骨的各个顶点以实现经过CPD变换后的参考颅骨与待配准三维颅骨在相同类型的边缘区域配准。

8.根据权利要求1所述的方法，其特征在于，对待配准三维颅骨与参考颅骨在边缘区域之外的其他区域进行配准包括：

步骤1，对待配准三维颅骨上的每一个顶点，求解它在参考颅骨上的最近点作为对应点，这两个顶点组成一个点对，并将它们之间的欧式距离作为配准误差；

步骤2，对所有点对按照配准误差从大到小的顺序进行排序，之后从排序后的第一个点对开始采样，得到用于具有紧支撑的径向基函数CSRBF算法的特征点{p_i|i＝1,2,...,m}和{q_i|i＝1,2,...,m}，其中，m为特征点个数，{p_i|i＝1,2,...,m}位于参考颅骨上，{q_i|i＝1,2,...,m}位于待配准三维颅骨上；其中，在采样时，所述特征点{p_i|i＝1,2,...,m}满足以下距离条件：两个特征点之间的距离大于0.5a，其中a＝3.66Δ，为CSRBF算法的支持范围，Δ为所有点对在XYZ方向轴上的最大位移；

步骤3，基于所述特征点{p_i|i＝1,2,...,m}和特征点{q_i|i＝1,2,...,m}，采用CSRBF算法计算非刚性变换，将该变换作用至参考颅骨的各个顶点以实现对配准误差大的区域进行调整，使两个颅骨在所有区域匹配上；

其中，步骤1至步骤3允许执行多次，直到执行次数超过用户设定的阈值或最大配准误差小于用户设定的另一阈值。

9.一种三维颅骨的自动配准装置，其特征在于，该装置包括：

训练模块，用于训练一个边缘类型分类器；

识别模块，用于利用所述边缘类型分类器自动识别待配准三维颅骨上的各条边缘，所述各条边缘至少包括：眼眶边缘、鼻框边缘、颞骨边缘、上颔骨边缘、下颔骨边缘和其他边缘；

校正模块，用于通过对识别出的鼻框边缘上的所有边缘点进行主成分分析PCA实现所述待配准三维颅骨的坐标校正，以调整待配准三维颅骨与参考颅骨之间具有近似的位置和姿态；

配准模块，用于对待配准三维颅骨与参考颅骨执行边缘区域的配准，之后对待配准三维颅骨与参考颅骨在边缘区域之外的其他区域进行配准。

10.根据权利要求9所述的装置，其特征在于，所述训练模块通过以下步骤训练一个边缘类型分类器：

对数据库中的每一个三维颅骨模型进行边缘提取，以提取出边界边、三维孔洞边缘、边缘点；其中，数据库中的每一个三维颅骨模型是采用三角形或四边形网格模型表示，边界边是指三维颅骨模型上只被一个三角形或四边形网格使用的边，三维孔洞边缘是由不同边界边首尾相连组成的，边缘点是指边界边两端的端点；

标注每一个三维颅骨模型中每一条三维孔洞边缘的类型；

从每一三维孔洞边缘提取形状分布特征；

对所有三维孔洞边缘，基于提取出的形状分布特征和标注的边缘类型并采用支持向量机SVM进行分类训练得到边缘类型分类器，三维孔洞边缘至少分为以下6类：眼眶边缘、鼻框边缘、颞骨边缘、上颔骨边缘、下颔骨边缘以及其他边缘；

所述识别模块通过以下步骤利用边缘类型分类器自动识别待配准三维颅骨上的各条边缘：

对待配准三维颅骨进行边缘提取，以提取出边界边、三维孔洞边缘、边缘点；其中，待配准三维颅骨是采用三角形或四边形网格模型表示，边界边是指待配准三维颅骨上只被一个三角形或四边形网格使用的边，三维孔洞边缘是由不同边界边首尾相连组成的，边缘点是指边界边上的两个端点；

提取每一三维孔洞边缘的形状分布特征；

针对每一三维孔洞边缘，根据从该三维孔洞边缘提取出的形状分布特征并利用所述边缘类型分类器自动识别该三维孔洞边缘的边缘类型，边缘类型至少为：眼眶边缘、鼻框边缘、颞骨边缘、上颔骨边缘、下颔骨边缘和其他边缘。