本发明涉及一种保持局部拓扑不变性的点配准算法。
背景技术：
图像配准技术是计算机视觉、图像处理领域一项非常重要的研究课题，已经成功应用于模式识别、医学图像处理以及数字媒体等多个重要领域；特别是在目前，虚拟现实和现实增强技术已逐渐走出实验室来到人民群众的生活中来，即将成为下一个经济增长点。而配准算法是该技术是否能成功应用的核心关键算法。根据配准对象可分为灰度配准算法和点配准算法；其中，由于点配准算法既可以有效地降低算法复杂度，又可以兼顾配准的精度，具有较高的实际应用价值；因此，点配准技术得到了人们特别地关注；目前，点配准技术按技术线路分为基于分类的方法和基于图匹配的方法。
基于分类的配准算法中主要应用了有限混合模型，特别是高斯混合模型；Chui和Rangarajan通过混合模型将特征点配准问题转化为极大似然估计问题，该方法使用一个确定性退火机制来实现了期望最大化算法，从而可以直接控制点对之间的对应模糊性；一致性点漂移算法将高斯混合模型中的质心，即为两个待配准点集中的模板点集合定义了一个速度场函数，并通过期望最大化算法计算高斯混合模型中的未知参数来实现配准。
Horaud等提出了一个条件期望最大化点配准算法，该方法利用了各向异性协方差模型和条件期望最大化算法来计算刚性及铰接点配准问题，为了增强算法的鲁棒性，这类方法都在混合模型中加入了一个额外的均匀成分；但是这个策略无法同时拟合两个点集中混入的异常点；Sanroma等使用点邻居关系提出了一个刚性点配准算法并且扩展了一个非刚性配准版本。
然而这些方法假设了一个平均混合权重作为了一个不适合的先验概率，实际上使用了简化了的后验概率；没有在混合模型中准确地拟合两个点集中出现的异常点，进而导致降低算法对于异常点的鲁棒性，特别是在两个点集中同时出现异常点时，配准的准确性会显著降低。
基于图匹配的方法是将点的配准问题转化为图的匹配问题，这类方法中，点集中的点作为图中的节点，低于设定阈值的点间的连线作为图的边；通过最大化匹配的边的数目实现配准点集，在计算两个点集合间的相似矩阵时，各增加一行和一列以用来拟合两个点集中的异常点。
Chui等提出了鲁棒性点配准方法，将点的形变和对应关系整合到一个目标函数中，并通过软分配和确定性退火来求解了这个优化问题；Zheng等提出在配准过程中保持局部点之间的邻居关系，并基于形状上下文定义了计算点邻居关系的匹配概率函数。
然而在计算点邻居关系时统计了所有邻域点之间的匹配概率，包含了大量的冗余信息；Lee等使用了标签松弛方法，并且通过离散的角度和半径差量化了点对之间的匹配相关性；但在比较点及其邻域间的拓扑差异时，同样使用了一对多的方式来统计计算相似性概率；这类方法忽略了局部匹配点之间的先验概率，即在计算局部相似性时，忽略了局部点对之间的相似概率，没有建立一个邻域点之间的对应关系，从而降低了配准的准确性。
技术实现要素：
针对现有技术中的上述不足，本发明提供的保持局部拓扑不变性的点配准算法精确地量化了局部拓扑差异，提高了非刚性配准在发生大形变时的准确性，并提高了存在异常点时算法的鲁棒性。
为了达到上述发明目的，本发明采用的技术方案为：提供一种保持局部拓扑不变性的点配准算法，其包括如下步骤：
S1.将点的配准转化为图的匹配；其中，点集中的点作为图中的节点；
S2.定义图中节点的邻域；其中，节点和其邻域整体构成一个局部子图，节点和其邻域中的一个点连接成局部子图中的一条边；
S3.对局部子图进行动态规划，得到局部子图间的匹配概率；
S4.通过局部子图间的匹配概率及松弛标记法对全局匹配概率进行处理，得到点对间的对应关系；
S5.根据点对间的对应关系，通过薄板样条函数对点集合进行弹性形变正则化处理，得到模型参数。
进一步地，S2的具体步骤为：设ST＝{t1,t2...tm}表示含有m个点的点集合T，对于一个给定点ti∈ST，定义ti的邻域为x＝1,2...X；其中，X表示为有X个位于点ti圆形邻域内的点；点ti和其邻域构成了一个关于点ti的局部子图，ti分别和邻域中的点连接成该局部子图的边。
进一步地，邻域的半径为ST中所有点对之间的欧式距离的均值。
进一步地，S3的具体步骤为：将动态规划分为n个阶段，其中，n表示为局部子图中边的数量；对局部子图中的边进行分析处理，得到每一条边的代价；根据代价差计算函数公式，计算出每一阶段任意两条边之间的代价差；根据匹配函数对代价进行规划处理，得到局部子图间的匹配概率。
进一步地，S32的具体步骤为：使用log距离和极坐标角度标定邻域，设中心点ti的距离为0，最远的距离为6，则中心点到邻域点的距离为以30度递进，共划分12个不同的角度值；以表示中心点和邻域点间的角度值，距离线和角度线将邻域范围共划分为72个栅格；以点所在栅格的距离值和角度值构成了一个二元组即为该边的代价。
进一步地，代价差计算公式为其中，参数α和β分别用于调节半径和角度在计算代价差时所占的权重，(ud,σd)和(uθ,σθ)分别为相应半径和角度的均值及方差。
进一步地，局部子图间的匹配概率为：
其中，表示为子图M中的点u对应于子图N中的点v，ω为权重参数。
进一步地，S4的具体步骤为：
利用(M+1)*(N+1)的概率矩阵P来表示各个点对之间的对应概率:
如果点m和点n相匹配，那么概率pmn＝1，否则pmn＝0，此时多出的一行和一列用于拟合点集合中的异常点，同时该矩阵满足：
迭代地更新概率矩阵P，根据当前最大的一个匹配概率，赋予点对间的对应关系。
进一步地，S5的具体步骤为：
通过薄板样条函数进行点集合的弹性形变正则化，在迭代过程中，匹配的点对来估计薄板样条函数的形变模型参数，之后用得到的模型参数来计算模板点集合的下一个形变位置，直到参数收敛到一个全局最优解。
本发明的有益效果为：该保持局部拓扑不变性的点配准算法将点的配准转化为图的匹配，在配准过程中，保持了准确的局部的拓扑关系；并通过动态规划得到局部子图间的拓扑差异，用松弛标记法进行优化处理，得到点对间的对应关系，再通过薄板样条函数来进行点集合的弹性形变正则化；其通过局部子图间的匹配概率来引导两个点集合间匹配的全局概率，在点配准发生形变的过程中，始终保持了节点及其邻域点之间的拓扑稳定性，提高了算法的准确性和鲁棒性。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一种实施例，而不是全部的实施例。基于本发明中的实施例，本邻域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明的保护范围。
为简单起见，以下内容中省略了该技术邻域技术人员所公知的技术常识。
该保持局部拓扑不变性的点配准算法包括以下步骤：
S1.将点的配准转化为图的匹配；其中，点集中的点作为图中的节点。
S2.定义图中节点的邻域；其中，节点和其邻域整体构成一个局部子图，节点和其邻域中的一个点连接成局部子图中的一条边；在具体实施中，将点的配准转化为图的匹配，其点的配准问题转化为两个图中匹配的边的数量的问题。
在实际操作中，设ST＝{t1,t2...tm}表示含有m个点的点集合T，SD＝{d1,d2...dn}表示含有n个点的点集合D，点匹配就是要找到一个匹配函数f:在两个点集合间建立起一种对应关系；当两个点集中点的个数不对等或者点集中存在异常点时，则无法形成一一对应关系；这时在两个点集中分别加入一个伪点点来表示这些无法匹配的点或者异常点，则两个点集合被扩充为和这时对应关系f:则满足一一对应关系，无法匹配的点或者异常点将对应于伪点
当两个点集合发生非刚性形变时，点间的距离变化可以是任意的；然而节点和其邻域间却保持着相对固定的关系；即在非刚性形变下，节点的局部邻域关系或者拓扑结构依然保持着稳定性，基于这一原理，对应的优化函数为：其中，
对于一个给定点ti∈ST，定义ti的邻域为x＝1,2...X；其中，X表示为有X个位于点ti圆形邻域内的点；点ti和其邻域构成了一个关于点ti的局部子图，ti分别和邻域中的点连接成该局部子图的边；邻域的半径为ST中所有点对之间的欧式距离的均值。
S3.对局部子图进行动态规划，得到局部子图间的匹配概率；在具体实施中，将动态规划分为n个阶段，其中，n表示为局部子图中边的数量；对局部子图中的边进行分析处理，得到每一条边的代价；根据代价差计算函数公式，计算出每一阶段任意两条边之间的代价差；根据匹配函数对代价进行规划处理，得到局部子图间的匹配概率。
在实际操作中，使用了log距离和极坐标角度标定邻域，设中心点ti的距离为0，最远的距离为6，则中心点到邻域点的距离为以30度递进，共划分12个不同的角度值；以表示中心点和邻域点间的角度值，距离线和角度线将邻域范围共划分为72个栅格；以点所在栅格的距离值和角度值构成了一个二元组即为该边的代价。
两条边之间的代价差计算公式为其中，参数α和β分别用于调节半径和角度在计算代价差时所占的权重，(ud,σd)和(uθ,σθ)分别为相应半径和角度的均值及方差；此时，局部子图间的匹配概率为：其中，表示为子图M中的点u对应于子图N中的点v，ω为权重参数。
S4.通过局部子图间的匹配概率及松弛标记法对全局匹配概率进行处理，得到点对间的对应关系；在具体实施中，利用(M+1)*(N+1)的概率矩阵P来表示各个点对之间的对应概率:
如果点m和点n相匹配，那么概率pmn＝1，否则pmn＝0，此时多出的一行和一列用于拟合点集合中的异常点，同时该矩阵满足：
迭代地更新概率矩阵P，根据当前最大的一个匹配概率，赋予点对间的对应关系。
S5.根据点对间的对应关系，通过薄板样条函数对点集合进行弹性形变正则化处理，得到模型参数；在具体实施中，通过薄板样条函数进行点集合的弹性形变正则化，在迭代过程中，匹配的点对来估计薄板样条函数的形变模型参数，之后用得到的模型参数来计算模板点集合的下一个形变位置，直到参数收敛到一个全局最优解。
该保持局部拓扑不变性的点配准算法的改进及其具体实施过程如下：
该保持局部拓扑不变性的点配准算法将点配准问题转化为图匹配，点集中的点作为图中的节点，节点和其邻域中的一个点组成图中的一条边；如此，点配准问题就转化为图中匹配的边的数量的问题。
设ST＝{t1,t2...tm}和SD＝{d1,d2...dn}分别表示含有m个点的点集合T和含有n个点的点集合D，点匹配就是要找到一个匹配函数f:在两个点集合间建立起一种对应关系，当两个点集中点的个数不对等或者点集中存在异常点时，则无法形成一一对应关系；这时在两个点集合分别加入一个伪点点来表示这些无法匹配的点或者异常点；则两个点集合被扩充为S'T＝{t1,t2...tM,nil}和这时对应关系f:则满足一一对应关系，无法匹配的点或者异常点将对应于伪点
当两个点集合发生非刚性形变时，点间的距离变化可以是任意的；然而节点和其邻域间却保持着相对固定的关系；即在非刚性形变下，节点的局部邻域关系或者拓扑结构依然保持着稳定性，基于这一原理，对应的优化函数为：其中，目标函数
由于Nm是点m的邻域，同时邻域关系是对称的，即如果i∈Nm，那么m∈Ni。f(m),f(n)表示两个点集合对应点的邻域系统；d(f(m),f(n))度量了这两个邻域系统间的距离，距离为1表示这两个邻域系统匹配，也就是邻域关系或者拓扑结构保持的越好，否则为0。
在非刚性配准过程中，其局部拓扑结构的稳定性保持的越好，配准算法的准确性和鲁棒性也就越高；因为点集合的弹性形变可以是任意的，但是如果节点及其邻域的拓扑关系发生了变化，那么该形变可以认为是无意义的，此时，两个点集合形状之间则失去了对应关系，现有技术的算法中简单地使用了点间的欧式距离来描述邻域关系及其邻域间的差异，此时使用目标函数替代公式其中，同理适用于其余点集合。
当匹配函数f实现了准确的一一对应关系后，那么的最小值为零；无论使用欧式距离或者使用半径角度差，这些方法在量化局部拓扑结构差异时，均忽略了局部邻域点间的对应关系，只是针对所有可能的对应点之间进行了概率求和，此时两个点的邻域间的局部匹配概率则为：qm,n＝∑i∈M∑j∈Nsi,j。
对于一个给定点ti∈ST，定义ti的邻域为x＝1,2...X；其中，X表示为有X个位于点ti圆形邻域内的点；点ti和其邻域构成了一个关于点ti的局部子图，ti分别和邻域中的点连接成该局部子图的边；邻域的半径为ST中所有点对之间的欧式距离的均值。
log距离和极坐标角度标定邻域，设中心点ti的距离为0，最远的距离为6，则中心点到邻域点的距离为以30度递进，共划分12个不同的角度值；以表示中心点和邻域点间的角度值，距离线和角度线将邻域范围共划分为72个栅格；其中，点ti和其邻域构成了一个关于点ti的局部子图，ti分别和邻域中的点连接成该局部子图的边；以点所在栅格的距离值和角度值构成了一个二元组即为该边的代价。
两条边之间的代价差计算公式为
其中，参数α和β分别用于调节半径和角度在计算代价差时所占的权重，(ud,σd)和(uθ,σθ)分别为相应半径和角度的均值及方差；显然，此时两个局部子图中所有相对应边的代价差总和即就是当前点的局部拓扑结构差异；如果局部子图m和局部子图n拓扑结构完全相同的话，那么其对应边的代价差之和为零，这时他们的相似度为最大。
然而，由于缺乏局部对应的先验概率，公式qm,n＝∑i∈M∑j∈Nsi,j不能精确地量化拓扑结构的差异，为了能够得到局部子图中边的对应概率，即局部子图m中的边i对应于局部子图n中的边j的概率，将求局部子图间代价差转化为了一个n阶段最小化目标函数的动态规划；此时qm,n＝Σi∈MΣj∈Nsi,j转化为qm,n＝Min(Σi∈MΣj∈Nsi,j)。
此时，局部子图中边的数量n表示该动态规划的n个阶段，在每一阶段根据代价差计算两条边的代价差，根据目标函数求出最小代价和作为其最优解，而得到的最小代价差之和即就是要求的两个局部子图间的拓扑结构差异。
将局部匹配转化为动态规划时，设某一阶段k为，每匹配一个特征点作为一个阶段；状态变量xK为匹配k个点前的总侯选点xK＝{x1,x2...xn-k+1}；决策变量dk为第k个匹配点的对应候选点；决策允许集合为dk∈xK；状态转移方程为xK+1＝xK-dk；阶段指标为Vk(xk,dk)；递推方程为fk(xk)＝max{Vk(xk,dk)+fk+1(xk+1)}；初始条件为fK(xK)＝VK(xK,dk)；终点条件为f1(x1)＝max{V1(x1,d1)+f2(x2)}，其中f2(x2)可在前第k-1步求得。
通过求解该动态规划问题可以得到一个子图间点的对应关系其表示局部子图M中的点u对应于局部子图N中的点v；通过上述的点对之间的相似性度量和动态规划，可建立局部子图之间的拓扑结构相似先验概率；此时，局部子图间的匹配概率为ω为权重参数。
然后，通过局部子图间的匹配概率及松弛标记法对全局匹配概率进行处理，得到点对间的对应关系；在具体实施中，使用松弛标记法来解这个图匹配的优化问题，利用一个(M+1)*(N+1)的概率矩阵P来表示各个点对之间的对应概率，其率矩阵为：
如果点m和点n相匹配，那么概率pmn＝1，否则pmn＝0，此时多出的一行和一列用于拟合点集合中的异常点，同时该矩阵满足：
基于此概率矩阵，可将目标函数改写为：通过公式迭代地更新概率矩阵P，之后根据当前最大的一个匹配概率，赋予点对间的对应关系。
最后，通过薄板样条函数来进行点集合的弹性形变正则化，在迭代过程中，匹配的点对来估计形变模型参数，然后用得到的模型参数来计算模板点集合的下一个形变位置，直到参数收敛到一个全局最优解。
在实际实验中，先初始化参数，赋初值；之后计算两个集合中点的形状上下文和基于形状上下文的点对之间的匹配概率；再使用概率初始化匹配概率矩阵P并正则化；并根据公式计算局部子图间拓扑结构匹配概率；使用公式更新匹配概率矩阵P，同时对P的行和列进行正则化处理。
接着根据P中点对的匹配概率来计算薄板样条函数的形变参数；再根据薄板样条函数形变参数更新点集的下一个位置；依次迭代次数加1，当到达最大迭代数Imax时或者达到最大匹配点数时，结束，否则重新开始计算两个集合中点的形状上下文和基于形状上下文的点对之间的匹配概率。
在具体实验中，可使用两个点集合-鱼形状点集合和中文字符形状点集合，为测试算法性能，对模型点集合采取不同程度的非刚性形变、形变加部分轮廓点缺失和增加异常点的退化处理；将每种图像由弱到强分为5-6个退化级别，在每个退化级别上根据退化程度各随机生成100个测试点集合；算法的性能误差由配准后两个点集合的平均欧式距离来度量，同时给出误差的方差分布。
该保持局部拓扑不变性的点配准算法考虑了点的局部拓扑关系，将局部点之间的对应关系作为点配准的先验概率；同时，基于局部子图中点对间的位置关系，即局部子图中边的对应关系，并通过动态规划得到局部子图间的拓扑差异，通过局部拓扑结构匹配概率来引导两个点集合间匹配的全局匹配概率，使得在点配准发生形变的过程中，始终保持了节点及其邻域点之间的拓扑稳定性。
对所公开的实施例的上述说明，使本邻域专业技术人员能够实现或使用本发明。对这些实施例的多种修改对本邻域的专业技术人员来说将是显而易见的，本文所定义的一般原理可以在不脱离发明的精神或范围的情况下，在其他实施例中实现。因此，本发明将不会被限制与本文所示的这些实施例，而是要符合与本文所公开的原理和新颖性特点相一致的最宽的范围。