CN105631877A - 基于加权双稀疏约束的非刚性表面配准方法 - Google Patents

Abstract

本发明属于计算机应用领域，为实现对位置和变换同时约束，获得更好的配准结果。为此，本发明采取的技术方案是，基于加权双稀疏约束的非刚性表面配准方法，找到模板点集中的每个点在目标点集中的对应点，应用双稀疏性的L1范数对非刚性配准的能量方程进行约束，从而将具有形变的三维表面数据进行非刚性配准。本发明主要应用于获得三维物体的完整模型。

Description

基于加权双稀疏约束的非刚性表面配准方法

技术领域

本发明属于计算机应用领域，三维物体完整模型的获取。具体讲，涉及基于加权双稀疏约束的非刚性表面配准方法。

背景技术

在计算机图形学和计算机视觉领域，非刚性配准是一个应用于动态三维重建的重要技术。在表面重建、三维物体识别等领域中，此项技术都有着十分重要的应用。近日，如Kinect等深度传感器的应用愈加广泛，但是其所采集到的点云包含太多噪声，得到的扫描数据也只能提供三维表面的一部分信息，必须经过配准融合之后才可以获取物体完整的数字模型。这就使得我们更需要一个对噪声和离群值鲁棒性强、可高效扫描有形变物体的动态场景的非刚性配准算法。

三维形状配准包含刚性配准和非刚性配准。刚性配准致力于找到全局的刚性变换，而非刚性配准则需要找到两个形状之间的一系列局部变换。迄今为止，针对于刚性配准的曲面配准算法得到了长足的发展，算法种类繁多各具特色。其中，由Besl等在1992年提出的ICP(IterativeClosestPoiht，迭代最近点)算法和它的变体依然是刚性配准的主流算法。ICP算法基于三维物体的纯粹几何模型，具有功能强大、精确度高、简单易解等优点。许多研究者也对这一算法进行了详尽的分析和一系列改进。

不同于研究成果较多的刚性配准，非刚性配准指的是从模板点到目标点的运动是非刚性的，它致力于找到模板形状与目标形状之间合适的形变转换。随着扫描设备对三维曲面捕捉能力的逐渐提升，非刚性配准也逐步受到了关注。传统的非刚性配准技术同样是使用迭代最近点与最小二乘法的传统ICP方法。这种方法简单易懂，复杂度相对较低。但是这种利用传统ICP算法来解决非刚性配准问题时，同样对离群值比较敏感。非刚性配准通常被定义为一个优化的问题。大多数方法制订了同时具有位置约束和转换约束的能量方程。其中，位置约束描述了转换后的模板形状和目标形状之间的接近度，转换约束则描述了模型的适合性。大多数工作在位置约束和转换约束中应用的是经典的L2范数(B.Amberg，S.Romdhani，andT.Vetter.Optimalstepnonrigidicpalgorithmsforsurfaceregistration.InIEEEConferenceonComputerVisionandPatternRecognition(CVPR)，pages1-8，2007.)，然而，L2范数形式下的能量方程更容易被噪声和离群值所影响。为了解决这个问题，Yang等(J.Yang，K.Li，K.Li，andY.-K.Lai.Sparsenon-rigidregistrationof3Dshapes.InComputerGraphicsForum，volume34，pages89-99，2015.)提出了稀疏非刚性配准(SNR)的方法，对转换约束应用L1范数模型。然而，他们的位置约束还是基于L2范数的。事实上，对于局部分段的刚性变换，位置误差倾向于集中在局部小区域中。而这是应用L2范数无法良好约束的。

发明内容

为克服现有技术的不足，实现对位置和变换同时约束，获得更好的配准结果。为此，本发明采取的技术方案是，基于加权双稀疏约束的非刚性表面配准方法，找到模板点集中的每个点在目标点集中的对应点，应用双稀疏性的L1范数对非刚性配准的能量方程进行约束，从而将具有形变的三维表面数据进行非刚性配准。

其中，找到模板点集中的每个点在目标点集中的对应点的具体步骤是，

1)找到模板点集中的每个点在目标点集中的对应点，为了保证准确性，利用高维特征之间的相似性进行排序，只取前面一部分更为可靠的对应点，找到特殊点的对应关系之后，再根据特殊点的对应关系，用k-d树最近邻搜索的方法来找目标点集所有剩余点的对应点；

2)将非刚性配准问题建模：

E(X；f)＝E_data(X；f)+αE_smooth(X)(1)

设为模板template点集，其中包含N个点，为目标target点集，其中包含M个点；其中，v_i是模板点集中的第i个点，u_i是目标点集中与之相对的对应点，f：{1...N}→{1...M}表示从模板点集到目标点集的对应索引映射，定义X_i是对于点v_i的一个3×4转换矩阵，那么则是一组非刚性变换矩阵，定义包含N个待解的转换矩阵，给定映射f，目的是找到将模板点集转换到目标点集的转换矩阵X；

E_data(X；f)和E_smooth(X)分别是数据项和平滑项，α是调整这两项的权重系数，其中，数据项E_data(X；f)衡量了各个点位置的准确性，平滑项E_smooth(X)保证临近的点拥有尽量相似的变换，在获得了较好的初值情况下，迭代最近点算法可以得到很好的收敛性，将上述方程重新描述，

$\begin{array}{cc}min& \left|\right|W(VX-{\tilde{U}}_f)|{|}_1+\mathrm{\α}|\left|BX\right|{|}_1\end{array}---\left(2\right)$

其中，W是权重项，对于每一个点，都有其对应的权重w_i； 是笛卡尔坐标系下的u_f(i)，u_f(i)是目标点集中与模板点集中第i个点相对应的点；引入矩阵其中，B的第i行ε为模板点集中边的个数，v为模板点集中顶点的个数，在接下来的步骤中，定义稀疏矩阵K，K的行数即为ε，K的列数即为v；

3)找到模板点集中的每个点在目标点集中的对应点：

用SHOT方法，找到从模板点集到目标点集的可靠对应点，此时所找到点的w_i设为1，其余设为0，第一次迭代结束后，除了之前的可靠对应点，再次寻找每个点的距离最近点，如果两个对应点之间的距离小于某个阈值，则认为是可靠的，此时所找到点w_i设为1，其余设为0。

应用双稀疏性的L1范数对非刚性配准的能量方程进行约束具体步骤是，

1)定义稀疏矩阵K，得到要解的最优化方程：

$\begin{array}{cc}min& \left|\right|W(VX-{\tilde{U}}_f)|{|}_1+\mathrm{\α}|\left|BX\right|{|}_1\end{array}$

1.1)对于模板点集中的每一个点，设其变换参数是一个3×4的变换矩阵X_i，模板中一共有n个对应点，则待解的矩阵是一个4n×3的矩阵X：[X₁…X_n]^T，每个模板点集中顶点的坐标为v_i，设v_i＝[x，y，z，1]^T，目标点集中与之相对应的点为u_i，w_i表示匹配的可靠性，如果在目标点集中没有找到与模板点集中的点v_i相对应的点，则将w_i设为0，否则设为1，

$W=diag\left(\begin{array}{ccc}\sqrt{w_1},& ...,& \sqrt{w_N}\end{array}\right)---\left(4\right)$

其中 是笛卡尔坐标系下的u_f(i)，u_f(i)是目标点集中与模板点集中第i个点相对应的点；

1.2)定义稀疏矩阵K，K的行数为模板点集中边的个数ε，K的列数为模板点集中顶点的个数v。如果第r条边连接第i个和第j个定点，则k_r，i＝1，k_r，j＝-1，

$E_{smooth}\left(X\right)={\mathrm{\Σ}}_{v_i\∈V}{\mathrm{\Σ}}_{v_j\∈N_i}\left|\right|X_i{v}_i-X_j{v}_i|{|}_1---\left(5\right)$

2)利用增广拉格朗日方法进行最终求解

利用增广拉格朗日方法进行最终求解具体步骤是，用SuiteSparse库来表示稀疏矩阵、进行LDL分解、并求解稀疏矩阵的线性方程：

总式：

$\begin{array}{c}L(X,C,A,Y_1,Y_2,{\mathrm{\&mu;}}_1,{\mathrm{\&mu;}}_2)\\ =\mathrm{\&alpha;}\left|\right|A|{|}_1+|\left|C\right|{|}_1+<Y_{1}C-W_D\left(VX-{\tilde{U}}_f\right)>+\frac{{\mathrm{\&mu;}}_1}{2}\left|\right|C-W_D\left(VX-{\tilde{U}}_f\right)|{|}_F^2\\ +<Y_2,A-W_{S}BX>+\frac{{\mathrm{\&mu;}}_2}{2}\left|\right|A-W_{S}BX|{|}_F^2\end{array}$

其中， $C-problem:C^{k+1}=\arg \min L_{{\mathrm{\&mu;}}^k}(C,A^K,X^k,Y_1^K,Y_2^K),A-problem:$

$A^{k+1}=\arg \min L_{{\mathrm{\&mu;}}^k}(A,C^K,X^k,Y_1^K,Y_2^K),X-problem:X^{k+1}=\arg \min L_{{\mathrm{\&mu;}}^k}(A^{K+1},C^{K+1},Y_1^K,Y_2^K),$

其中，(μ₁，μ₂)是正的常数，(Y₁，Y₂)是拉格朗日乘子，<·，·>表示将两个矩阵看成长向量的内积。在增广拉格朗日解法的框架下，(μ₁，μ₂)和(Y₁，Y₂)可以有效更新，对变量C、A与X进行迭代最小化，更新拉格朗日乘子Y，最终得到变换矩阵X。

本发明的特点及有益效果是：

本发明用迭代最近点的算法计算模板形状和目标形状之间的变形关系。针对有形变的非刚性数据进行配准。相比于L1范数，L2范数更容易受到噪声和离群值的影响，基于此，在迭代最近点的基础上，我们提出了双稀疏的非刚性配准算法，也就是说对于数据项和平滑项都用L1范数来进行稀疏约束。这样可以更好的识别正常数据及离群值。它具有以下特点：

1、简单易懂，复杂度相对较低，易于实现。

2、利用SHOT方法采用高维特征进行特征点的寻找，这样找到的特征点可靠性更高。

3、L2范数形式下的能量方程更容易被噪声和离群值所影响。约束性最好的L0范数具有非凸性，这使得求解变得非常困难。所以我们采用L0范数的最优凸近似L1范数进行约束，L1范数最小化是凸优化问题，可以进行线性方程的求解。

4、由于三维数据点集数据庞大，用SuiteSparse库来表示稀疏矩阵、进行LDL分解、并求解稀疏矩阵的线性方程。

附图说明：

本发明上述的和/或附加的方面和优点从下面结合附图对实施例的描述中将变得明显和容易理解：

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

图2为应用加权双稀疏约束方法在跳跃数据库上的配准结果。蓝色的点云代表重叠变形后的目标形状，灰色的点云代表模板形状，在重建的点云上，配准误差以颜色图表征；

图2(a)为模板三维模型(上)和目标三维模型(下)；

图2(b)为应用加权双稀疏约束方法的配准结果。蓝色的点云代表重叠变形后的目标形状，灰色的点云代表模板形状，在重建的点云上，配准误差以颜色图表征。

具体实施方式

本发明根据扫描系统提供的三维表面数据，采用加权双稀疏约束，对有形变的两个三维物体表面进行配准，从而可以获得三维物体的完整模型。，请审核看可否按这样的文本递交申请。

本发明克服了已有方法的不足，采用双L1范数对位置和变换同时约束，获得了更好的配准结果。具体来说，找到模板点集中的每个点在目标点集中的对应点，应用双稀疏性的L1范数对非刚性配准的能量方程进行约束，从而将具有形变的三维表面数据进行非刚性配准。具体方法包括以下步骤：

1)找到模板点集中的每个点在目标点集中的对应点。为了保证准确性，我们利用高维特征之间的相似性进行排序，只取前面一部分更为可靠的对应点。找到特殊点的对应关系之后，再根据特殊点的对应关系，用k-d树最近邻搜索的方法来找目标点集所有剩余点的对应点。

2)将非刚性配准问题建模：

E(X；f)＝E_data(X；f)+αE_smooth(X)(1)

设为模板(template)点集，其中包含N个点。为目标(target)点集，其中包含M个点。其中，v_i是模板点集中的第i个点，u_i是目标点集中与之相对的对应点。f：{1...N}→{1...M}表示从模板点集到目标点集的对应索引映射。定义X_i是对于点v_i的一个3×4转换矩阵，那么则是一组非刚性变换矩阵。定义包含N个待解的转换矩阵。给定了映射f，算法的目的是找到将模板点集转换到目标点集的转换矩阵X。

E_data(X；f)和E_smooth(X)分别是数据项和平滑项，α是调整这两项的权重系数。其中，数据项E_data衡量了各个点位置的准确性，平滑项E_smooth保证临近的点拥有尽量相似的变换。在获得了较好的初值情况下，迭代最近点算法可以得到很好的收敛性，将上述方程重新描述，

$\begin{array}{cc}min& \left|\right|W(VX-{\tilde{U}}_f)|{|}_1+\mathrm{\&alpha;}|\left|BX\right|{|}_1\end{array}---\left(2\right)$

其中矩阵的构建，方法见4)。

其中，W是权重项，对于每一个点，都有其对应的权重w_i， 是笛卡尔坐标系下的u_f(i)，u_f(i)是目标点集中与模板点集中第i个点相对应的点；引入矩阵其中，B的第i行|ε|为模板点集中边的个数，|v|为模板点集中顶点的个数。在接下来的步骤中，我们定义了稀疏矩阵K，K的行数即为|ε|，K的列数即为|v|。具体设定见3)。

3)找到模板点集中的每个点在目标点集中的对应点。

用SHOT方法，找到从模板点集到目标点集的可靠对应点，此时所找到点的w_i设为1，其余设为0。第一次迭代结束后，除了之前的可靠对应点，再次寻找每个点的距离最近点，如果两个对应点之间的距离小于某个阈值，则认为是可靠的，此时所找到点w_i设为1，其余设为0。

4)定义稀疏矩阵K，得到要解的最优化方程：

$\begin{array}{cc}min& \left|\right|W(VX-{\tilde{U}}_f)|{|}_1+\mathrm{\&alpha;}|\left|BX\right|{|}_1\end{array}$

41)对于模板点集中的每一个点，设其变换参数是一个3×4的变换矩阵X_i，模板中一共有n个对应点，则待解的矩阵是一个4n×3的矩阵X：[X₁…X_n]^T。每个模板点集中顶点的坐标为v_i，设v_i＝[x，y，z，1]^T，目标点集中与之相对应的点为u_i，w_i表示匹配的可靠性，如果在目标点集中没有找到与模板点集中的点v_i相对应的点，则将w_i设为0，否则设为1。

$W=diag\left(\begin{array}{ccc}\sqrt{w_1},& ...,& \sqrt{w_N}\end{array}\right)---\left(4\right)$

其中 是笛卡尔坐标系下的u_f(i)，u_f(i)是目标点集中与模板点集中第i个点相对应的点

42)定义稀疏矩阵K，K的行数为模板点集中边的个数|ε|，K的列数为模板点集中定点的个数|v|。如果第r条边连接第i个和第j个定点，则k_r，i＝1，k_r，j＝-1。

$E_{smooth}\left(X\right)={\mathrm{\&Sigma;}}_{v_i\&Element;V}{\mathrm{\&Sigma;}}_{v_j\&Element;N_i}\left|\right|X_i{v}_i-X_j{v}_i|{|}_1---\left(5\right)$

5)利用增广拉格朗日方法进行最终求解。用SuiteSparse库来表示稀疏矩阵、进行LDL分解、并求解稀疏矩阵的线性方程。

52)总式：

$\begin{array}{c}L(X,C,A,Y_1,Y_2,{\mathrm{\&mu;}}_1,{\mathrm{\&mu;}}_2)\\ =\mathrm{\&alpha;}\left|\right|A|{|}_1+|\left|C\right|{|}_1+<Y_{1}C-W_D\left(VX-{\tilde{U}}_f\right)>+\frac{{\mathrm{\&mu;}}_1}{2}\left|\right|C-W_D\left(VX-{\tilde{U}}_f\right)|{|}_F^2\\ +<Y_2,A-W_{S}BX>+\frac{{\mathrm{\&mu;}}_2}{2}\left|\right|A-W_{S}BX|{|}_F^2\end{array}$

其中， $C-problem:C^{k+1}=\arg \min L_{{\mathrm{\&mu;}}^k}(C,A^K,X^k,Y_1^K,Y_2^K),A-problem:$

$A^{k+1}=\arg \min L_{{\mathrm{\&mu;}}^k}(A,C^K,X^k,Y_1^K,Y_2^K),X-problem:X^{k+1}=\arg \min L_{{\mathrm{\&mu;}}^k}(A^{K+1},C^{K+1},Y_1^K,Y_2^K),$

其中，(μ₁，μ₂)是正的常数，(Y₁，Y₂)是拉格朗日乘子，<·，·>表示将两个矩阵看成长向量的内积。在增广拉格朗日解法的框架下，(μ₁，μ₂)和(Y₁，Y₂)可以有效更新。对变量C、A与X进行迭代最小化，更新拉格朗日乘子Y，最终可以得到我们需要的变换矩阵X。

下面结合附图和具体实施方式进一步详细说明本发明。

本发明用迭代最近点的算法计算模板形状和目标形状之间的变形关系。在迭代最近点的基础上，找到模板点集中的每个点在目标点集中的对应点。应用双稀疏性的L1范数对非刚性配准的能量方程进行约束，从而将具有形变的三维表面数据进行非刚性配准。

本发明的主要框架迭代最近点有两个步骤，首先，通过上一次迭代结果预测模板点集与目标点集之间的对应关系，应用迭代最近点的方法找到二者之间的映射。然后，建立双稀疏性约束的能量方程，解优化方程得到最优解X。将X代入第一步骤中重新找对应点，这样迭代多次之后便可以得到一个最优解。

1)为模板点集中的每个点找到在目标点集中的对应点。用SHOT方法为模板点集中的每个点找到在目标点集中的对应点，为了确保准确性，我们将高维特征之间的相似性进行排序，只取前面一部分更为可靠的对应点。找到特殊点的对应关系之后，再根据特殊点的对应关系，用k-d树最近邻搜索的方法来找目标点集所有剩余点的对应点。

2)将非刚性配准问题建模：

E(X；f)＝E_data(X；f)+αE_smooth(X)(1)

设为模板(template)点集，其中包含N个点。为目标(target)点集，其中包含M个点。其中，v_i是模板点集中的第i个点，u_i是目标点集中与之相对的对应点。f：{1...N}→{1...M}表示从模板点集到目标点集的对应索引映射。定义X_i是对于点v_i的一个3×4转换矩阵，那么则是一组非刚性变换矩阵。定义包含N个待解的转换矩阵。给定了映射f，算法的目的是找到将模板点集转换到目标点集的转换矩阵X。

E_data(X；f)和E_smooth(X)分别是数据项和平滑项，α是调整这两项的权重系数。其中，数据项E_data衡量了各个点位置的准确性，平滑项E_smooth保证临近的点拥有尽量相似的变换。在获得了较好的初值情况下，迭代最近点算法可以得到很好的收敛性，将上述方程重新描述，

$\begin{array}{cc}min& \left|\right|W(VX-{\tilde{U}}_f)|{|}_1+\mathrm{\&alpha;}|\left|BX\right|{|}_1\end{array}---\left(2\right)$

其中矩阵的构建，方法见4)。

其中，W是权重项，对于每一个点，都有其相对应的权重w_i； 是笛卡尔坐标系下的u_f(i)，u_f(i)是目标点集中与模板点集中第i个点相对应的点；引入矩阵其中，B的第i行|ε|为模板点集中边的个数，|v|为模板点集中顶点的个数。在接下来的步骤中，我们定义了稀疏矩阵K，K的行数即为|ε|，K的列数即为|v|。具体设定见3)。

3)找到模板点集中的每个点在目标点集中的对应点。

31)用SHOT方法，找到从模板点集到目标点集的前20％的可靠对应点，此时所找到点的w_i设为1，其余设为0。第一次迭代结束后，除了前20％的可靠对应点，再次寻找每个点的距离最近点，如果两个对应点之间的距离小于数据总长度的1/20，则认为是可靠的，此时所找到点w_i设为1，其余设为0。

4)定义稀疏矩阵K，得到要解的最优化方程：

$\begin{array}{cc}min& \left|\right|W(VX-{\tilde{U}}_f)|{|}_1+\mathrm{\&alpha;}|\left|BX\right|{|}_1\end{array}$

41)对于模板点集中的每一个点，设其变换参数是一个3×4的变换矩阵X_i，模板中一共有n个对应点，则待解的矩阵是一个4n×3的矩阵X：[X₁…X_n]^T。每个模板点集中顶点的坐标为v_i，设v_i＝[x，y，z，1]^T，目标点集中与之相对应的点为u_i，w_i表示匹配的可靠性，如果在目标点集中没有找到与模板点集中的点v_i相对应的点，则将w_i设为0，否则设为1。

$W=diag\left(\begin{array}{ccc}\sqrt{w_1},& ...,& \sqrt{w_N}\end{array}\right)---\left(4\right)$

其中 是笛卡尔坐标系下的u_f(i)，u_f(i)是目标点集中与模板点集中第i个点相对应的点

42)定义稀疏矩阵K，K的行数为模板点集中边的个数|ε|，K的列数为模板点集中顶点的个数|v|。如果第r条边连接第i个和第j个定点，则k_r，i＝1，k_r，i＝-1。引入矩阵其中，B的第i行 $b_{i:}:=k_{i:}\&CircleTimes;v_i^T$

$E_{smooth}\left(X\right)={\mathrm{\&Sigma;}}_{v_i\&Element;V}{\mathrm{\&Sigma;}}_{v_j\&Element;N_i}\left|\right|X_i{v}_i-X_j{v}_i|{|}_1---\left(5\right)$

综合以上两式，要解的最优化方程为：

$\begin{array}{cc}min& \left|\right|W(VX-{\tilde{U}}_f)|{|}_1+\mathrm{\&alpha;}|\left|BX\right|{|}_1\end{array}$

5)利用增广拉格朗日方法进行最终求解。

51)用SuiteSparse库来表示稀疏矩阵、进行LDL分解、并求解稀疏矩阵的线性方程。

52)总式：

$\begin{array}{c}L(X,C,A,Y_1,Y_2,{\mathrm{\&mu;}}_1,{\mathrm{\&mu;}}_2)\\ =\mathrm{\&alpha;}\left|\right|A|{|}_1+|\left|C\right|{|}_1+<Y_{1}C-W_D\left(VX-{\tilde{U}}_f\right)>+\frac{{\mathrm{\&mu;}}_1}{2}\left|\right|C-W_D\left(VX-{\tilde{U}}_f\right)|{|}_F^2\\ +<Y_2,A-W_{S}BX>+\frac{{\mathrm{\&mu;}}_2}{2}\left|\right|A-W_{S}BX|{|}_F^2\end{array}$

其中， $C-problem:C^{k+1}=\arg \min L_{{\mathrm{\&mu;}}^k}(C,A^K,X^k,Y_1^K,Y_2^K),A-problem:$

$A^{k+1}=\arg \min L_{{\mathrm{\&mu;}}^k}(A,C^K,X^k,Y_1^K,Y_2^K),X-problem:X^{k+1}=\arg \min L_{{\mathrm{\&mu;}}^k}(A^{K+1},C^{K+1},Y_1^K,Y_2^K),$

其中，(μ₁，μ₂)是正的常数，(Y₁，Y₂)是拉格朗日乘子，<·，·>表示将两个矩阵看成长向量的内积。在增广拉格朗日解法的框架下，(μ₁，μ₂)和(Y₁，Y₂)可以有效更新。对辅助变量C、A与变量X进行迭代最小化，更新拉格朗日乘子Y，最终可以得到我们需要的变换矩阵X。α根据每次迭代中E_data项与E_smooth项的值来选取，调整二者之间的比重。

Claims (4)

1.一种基于加权双稀疏约束的非刚性表面配准方法，其特征是，找到模板点集中的每个点在目标点集中的对应点，应用双稀疏性的L1范数对非刚性配准的能量方程进行约束，从而将具有形变的三维表面数据进行非刚性配准。

2.如权利要求1所述的基于加权双稀疏约束的非刚性表面配准方法，其特征是，找到模板点集中的每个点在目标点集中的对应点的具体步骤是，

1)找到模板点集中的每个点在目标点集中的对应点，为了保证准确性，利用高维特征之间的相似性进行排序，只取前面一部分更为可靠的对应点，找到特殊点的对应关系之后，再根据特殊点的对应关系，用k-d树最近邻搜索的方法来找目标点集所有剩余点的对应点；

2)将非刚性配准问题建模：

E(X；f)＝E_data(X；f)+αE_smooth(X)(1)

设为模板template点集，其中包含N个点，为目标target点集，其中包含M个点；其中，v_i是模板点集中的第i个点，u_i是目标点集中与之相对的对应点，f：{1...N}→{1...M}表示从模板点集到目标点集的对应索引映射，定义X_i是对于点v_i的一个3×4转换矩阵，那么则是一组非刚性变换矩阵，定义包含N个待解的转换矩阵，给定映射f，目的是找到将模板点集转换到目标点集的转换矩阵X；

E_data(X；f)和E_smooth(X)分别是数据项和平滑项，α是调整这两项的权重系数，其中，数据项E_data(X；f)衡量了各个点位置的准确性，平滑项E_smooth(X)保证临近的点拥有尽量相似的变换，在获得了较好的初值情况下，迭代最近点算法可以得到很好的收敛性，将上述方程重新描述，

其中，W是权重项，对于每一个点，都有其对应的权重w_i； 是笛卡尔坐标系下的u_f(i)，u_f(i)是目标点集中与模板点集中第i个点相对应的点；引入矩阵其中，B的第i行ε为模板点集中边的个数，v为模板点集中顶点的个数，在接下来的步骤中，定义稀疏矩阵K，K的行数即为ε，K的列数即为v；

3)找到模板点集中的每个点在目标点集中的对应点：

用SHOT方法，找到从模板点集到目标点集的可靠对应点，此时所找到点的w_i设为1，其余设为0，第一次迭代结束后，除了之前的可靠对应点，再次寻找每个点的距离最近点，如果两个对应点之间的距离小于某个阈值，则认为是可靠的，此时所找到点w_i设为1，其余设为0。

3.如权利要求1所述的基于加权双稀疏约束的非刚性表面配准方法，其特征是，应用双稀疏性的L1范数对非刚性配准的能量方程进行约束具体步骤是，

1)定义稀疏矩阵K，得到要解的最优化方程：

1.1)对于模板点集中的每一个点，设其变换参数是一个3×4的变换矩阵X_i，模板中一共有n个对应点，则待解的矩阵是一个4n×3的矩阵X:[X₁…X_n]^T，每个模板点集中顶点的坐标为v_i，设v_i＝[x,y,z,1]^T，目标点集中与之相对应的点为u_i，w_i表示匹配的可靠性，如果在目标点集中没有找到与模板点集中的点v_i相对应的点，则将w_i设为0，否则设为1，

其中 是笛卡尔坐标系下的u_f(i)，u_f(i)是目标点集中与模板点集中第i个点相对应的点；

1.2)定义稀疏矩阵K，K的行数为模板点集中边的个数ε，K的列数为模板点集中顶点的个数v。如果第r条边连接第i个和第j个定点，则k_r,i＝1,k_r,j＝-1,

2)利用增广拉格朗日方法进行最终求解。

4.如权利要求3所述的基于加权双稀疏约束的非刚性表面配准方法，其特征是，利用增广拉格朗日方法进行最终求解具体步骤是，用SuiteSparse库来表示稀疏矩阵、进行LDL分解、并求解稀疏矩阵的线性方程：

总式：

其中，

其中，(μ₁,μ₂)是正的常数，(Y₁,Y₂)是拉格朗日乘子，<·，·>表示将两个矩阵看成长向量的内积。在增广拉格朗日解法的框架下，(μ₁，μ₂)和(Y₁，Y₂)可以有效更新，对变量C、A与X进行迭代最小化，更新拉格朗日乘子Y，最终得到变换矩阵X。