CN105486249A - 三维扫描数据的自适应底面消除方法 - Google Patents

Abstract

本发明公开了一种三维扫描数据的自适应底面消除方法，包括以下步骤：S1、将三维扫描数据通过坐标转换投影到近似投影平面上；S2、根据近似投影平面上的数据形成按矩阵形式排列的数据；将数据进行区域划分，并为每一个分区数据建立局部坐标系，将每个分区的数据在对应的局部坐标系进行转换；S3、针对每一个分区，利用差值参量，根据其局部坐标系中每个数据的坐标值计算在拓印基准面上对应的坐标值；S4、对步骤S3处理得到的拓印基准面上的数据投影到平面上，之后进行去噪。本发明的方法在不降低原始三维数据精度的前提下，对全局的三维数据在不转化为面绘制模型数据的基础上直接进行处理。在处理过程中避免了多余的人工参与对数据的干扰。

Description

三维扫描数据的自适应底面消除方法

技术领域

本发明属于书法拓印领域，更具体涉及一种三维扫描数据的自适应底面消除方法。

背景技术

三维扫描可以应用到许多领域，例如制造业、逆向工程、医学、文化产业等。三维扫描仪器按原理不同分为不同类别，但采集的都是点云数据。现有技术中，点云数据的处理均在三维空间中进行，处理方法基本分为两类：

1、先利用专业软件将点云重构为曲面模型，再导入到建模正向处理软件中进行模型编辑。

2、根据点云数据行为特征，利用逆向工程类软件重构为曲面，切换点云数据和曲面进行编辑。

第一种处理方法被广泛应用于文化产业，第二种处理方法被广泛应用于工业领域。由于三维扫描在工业制造逆向工程领域起步较早，因此第二种处理方法中逆向工程类软件对点云处理通常具备较强的曲面重构处理能力，常用的逆向软件有Imageware、GeomagicStudio、RapidForm等。现有技术中的两种处理方法在编辑之前，需要从专业软件中先进行数据精简，再转为面绘制模型，这就有可能丢失了对处理结果至关重要的细节信息，以致后期的模型结果失真。

以对钟碑文的扫描为例，在进行钟碑刻扫描数据的去底面处理时，需要根据以下几个步骤进行处理：

【步骤一】在逆向软件中打开扫描数据，将点云数据调整至底面与操作平面呈垂直状态。

【步骤二】利用点云选择工具，选择底面点云数据进行删除操作。

【步骤三】将剩余的数据转为obj或其他通用面绘制通用格式导出。

【步骤四】将面绘制模型导入建模软件中，对剩余没有切割彻底的模型面片进行二次筛选删除。

上述方法在处理钟碑刻扫描数据时，运用现有的工具对钟碑刻扫描数据去除底面的工作中，需要降低点云的数据量进行格式转化，以满足不同的处理软件对原始数据的要求，在此过程中不但增加了大量处理时间以及冗余的数据信息，而且丢失了对处理结果至关重要的细节信息，以致后期的模型结果失真。

另外，现有技术中在去底面的工作中，需要操作者对精简数据后的模型数据在主观意识下进行编辑，这要求操作者具备一定的经验和三维软件的操作技能，导致每个人的处理结果是不一样的，不能完全复现钟碑刻的表面文字雕饰信息。

发明内容

本发明要解决的技术问题是如何自动分离三维扫描数据中底面背景部分的数据，从而获得拓印基准面上的数据，同时不会由于精简数据造成信息丢失。

为了解决上述技术问题，本发明提供了一种三维扫描数据的自适应底面消除方法，所述方法包括以下步骤：

S1、将三维扫描数据通过坐标转换投影到近似投影平面上；

S2、根据所述近似投影平面上的数据形成按矩阵形式排列的数据；

将所述数据进行区域划分，并为每一个分区数据建立局部坐标系，将每个分区的数据在对应的局部坐标系进行转换；

S3、针对每一个分区，利用差值参量，根据其局部坐标系中每个数据的坐标值计算在拓印基准面上对应的坐标值，实现将所述步骤S2处理得到的数据投影到拓印基准面，完成数据分割；

S4、对所述步骤S3处理得到的所述拓印基准面上的数据投影到平面上，之后进行去噪。

优选地，所述步骤S4还包括以下步骤：遍历投影深度参量进行反低通滤波，其中所述投影深度为对应数据点的深度信息。

优选地，所述步骤S1中利用如下公式将三维扫描数据通过坐标转换投影到所述近似投影平面上：

$\left(\begin{array}{c}{x}^{\′}\\ y^{\′}\\ z^{\′}\end{array}\right)=\left(\begin{array}{ccc}\frac{a_1}{A}& \frac{a_2}{A}& \frac{a_3}{A}\\ \frac{b_1}{B}& \frac{b_2}{B}& \frac{b_3}{B}\\ \frac{c_1}{C}& \frac{c_2}{C}& \frac{c_3}{C}\end{array}\right)\left(\left(\begin{array}{c}x\\ y\\ z\end{array}\right)-\left(\begin{array}{c}{x}_0\\ y_0\\ z_0\end{array}\right)\right)$

式中，(x,y,z)^T为所述三维扫描数据中一个元素的坐标，(x′,y′,z′)^T为所述近似投影平面上对应元素的坐标，

$\left(\begin{array}{c}A\\ B\\ C\end{array}\right)=\left(\begin{array}{c}\sqrt{a_1^2+a_2^2+a_3^2}\\ \sqrt{b_1^2+b_2^2+b_3^2}\\ \sqrt{c_1^2+c_2^2+c_3^2}\end{array}\right)$

其中，c＝(c₁,c₂,c₃)^T

b＝b-b×c×c

$a=b\×c=\left|\begin{array}{ccc}i& j& k\\ b_1& b_2& b_3\\ c_1& c_2& c_3\end{array}\right|=(b_2{c}_3-b_3{c}_2,b_1{c}_3-b_3{c}_1,b_1{c}_2-b_2{c}_1)$

a、b、c为所述坐标转换的单位坐标矢量，(x₀，y₀，z₀)为原点；

所述近似投影平面为：

c₁x+c₂y+c₃z＝0。

优选地，所述坐标转换的单位坐标矢量c的计算方法如下：

计算所述三维扫描数据到所述近似投影平面的距离平方和；

建立对该距离平方和进行最小二乘法拟合的公式，并利用高斯消元法求解公式得到坐标矢量c。

优选地，所述步骤S1之前还包括以下步骤：

利用三维扫描仪按照预定行数和预定列数进行扫描得到所述三维扫描数据。

优选地，所述步骤S2中形成按矩阵形式排列的数据包括以下步骤：

S21、在所述近似投影平面上的数据中找到最下面的扫描点、最上面的扫描点、最左面的扫描点以及最右面的扫描点，并以四个扫描点作为顶点形成四边形；

S22、以所述四边形的底边作为一条边，以最上面的扫描点和最右面的扫描点作为顶点形成一个矩形；

S23、利用所述矩形对对应的所述近似投影平面上的数据进行剪裁，去掉所述矩形之外的数据；

S24、计算所述步骤S23得到的数据的宽度和高度；其中所述宽度为所述矩形内数据的列数，所述高度为所述矩形内数据的行数；

S25、对于所述矩形内每一个数据，根据其三维坐标确定其在所述按矩阵形式排列的数据中的行和列，并建立其深度信息与其在所述按矩阵形式排列的数据中位置的英映射关系。

优选地，所述步骤S25中，对于所述矩形内每一个数据，其在所述按矩阵形式排列的数据中的行、列以及所述映射关系利用下面公式计算：

$i=height*\frac{y_m-y_{\min }}{y_{\max }-y_{\min }+\mathrm{\ϵ}}$

$j=width*\frac{x_m-x_{\min }}{x_{\max }-x_{\min }+\mathrm{\ϵ}}$

D(i,j)＝z'_m

m＝1,...,numpoint

式中，i表示行，j表示列，D(i，j)＝z_m为所述映射关系，m为所述矩形内数据总个数，(x_m，y_m，z_m)为数据在所述近似投影平面上的位置坐标。

优选地，所述步骤S23还包括确定列数的子步骤，具体为：

对于所述矩形内数据中两个相邻的数据点，计算其在X方向上的距离是否大于列距离阈值，若是则所述两个相邻的数据点为两列。

优选地，所述步骤S25之还包括去除所述图像矩阵中空缺点步骤：

以所述空缺点为中心，计算其周边多个数据点的深度值，并以计算得到的多个深度值的平均值作为所述空缺点的深度值。

优选地，所述步骤S2中所述局部坐标系的坐标轴为：

$e^x={(\frac{(x_2-x_1)}{w},\frac{(y_2-y_1)}{w})}^T$

$e^y={(\frac{(x_3-x_2)}{h},\frac{(y_3-y_2)}{h})}^T$

式中，w为所述矩形的宽，h为所述矩形的高，(x₁,y₁)、(x₂,y₂)、(x₃,y₃)、(x₄,y₄)为所述矩形的四个顶点

所述局部坐标系的变量为：

$\mathrm{\ξ}=\frac{2(x_2-x_1)(x-x_0)}{w}+\frac{2(y_2-y_1)(y-y_0)}{w}$

$\mathrm{\η}=\frac{2(x_3-x_2)(x-x_0)}{h}+\frac{2(y_3-y_2)(y-y_0)}{h}$

其中，

$w=\sqrt{{(x_2-x_1)}^2+{(y_2-y_1)}^2}=\sqrt{{(x_3-x_4)}^2+{(y_3-y_4)}^2}$

$h=\sqrt{{(x_3-x_2)}^2+{(y_3-y_2)}^2}=\sqrt{{(x_4-x_1)}^2+{(y_4-y_1)}^2}.$

优选地，所述步骤S3中，利用下面公式计算局部坐标系下数据在所述拓印基准面上的坐标：

式中，(X_k,Y_k,Z_k){k＝0:8}为所述插值参量，为三维扫描数据在对应局部坐标系中的坐标值。

优选地，所述插值参量的计算包括以下步骤：

计算局部坐标系下数据到所述拓印基准面的距离平方和；

利用最小二乘法，按照所述差值参量使所述平方和最小计算得到所述插值参量。

本发明提供了一种三维扫描数据的自适应底面消除方法，本发明的方法在不降低原始三维数据精度的前提下，对全局的三维数据在不转化为面绘制模型数据的基础上直接进行处理。并且在处理过程方面，避免了多余的人工参与对数据的干扰，降低操作者在三维模型编辑处理方面的准备知识。同时本发明的方法可以在短时间内去除钟碑刻底面部分数据，得到较好的分离结果，能对任意扫描仪的通用.asc格式钟碑刻点云直接进行处理。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为本发明的一个较佳实施例的三维扫描数据的自适应底面消除方法流程图；

图2为本发明原始的三维扫描数据示意图；

图3为本发明中步骤S2处理后的效果示意图；

图4为本发明中步骤S3处理后的效果示意图；

图5为本发明中步骤S4处理后的效果示意图；

图6为本发明中步骤S3中差值分割示意图。

具体实施方式

下面结合附图和实施例对本发明作进一步详细描述。以下实施例用于说明本发明，但不能用来限制本发明的范围。

一种三维扫描数据的自适应底面消除方法，如图1所示，所述方法包括以下步骤：

S1、将三维扫描数据通过坐标转换投影到近似投影平面上；

S2、根据所述近似投影平面上的数据形成按矩阵形式排列的数据；

将所述数据进行区域划分，并为每一个分区数据建立局部坐标系，将每个分区的数据在对应的局部坐标系进行转换；

S3、针对每一个分区，利用差值参量，根据其局部坐标系中每个数据的坐标值计算在拓印基准面上对应的坐标值，实现将所述步骤S2处理得到的数据投影到拓印基准面，完成数据分割；

S4、对所述步骤S3处理得到的所述拓印基准面上的数据投影到平面上，之后进行去噪。

上述方法在不降低原始三维数据精度的前提下，对全局的三维数据在不转化为面绘制模型数据的基础上直接进行处理。并且在处理过程方面，避免了多余的人工参与对数据的干扰，降低操作者在三维模型编辑处理方面的准备知识。同时本发明的方法可以在短时间内去除钟碑刻底面部分数据，得到较好的分离结果，能对任意扫描仪的通用.asc格式钟碑刻点云文件直接进行处理。

进一步地，所述步骤S4中，利用反低通滤波进行去噪。将经过步骤S3处理后得到的拓印基准面投影到平面上，形成平面字帖。对于平面碑刻的情况，这种投影可以直接将拓印基准面的底面直接投影到平面上。对于钟鼎铭文等凸版刻字，在投影时要考虑近似曲面的曲率变化和近似曲面上浮尺寸的影响。同时，因底面不平整造成的斑点、底面或字模损坏造成的缺陷等干扰信息，需要降噪处理。在步骤S4中通过反低通滤波算法来消除近似曲面的曲率变化对投影的影响以及降噪处理。图5显示了经过反低通滤波滤除底面噪点和畸变影响后的图像。

进一步地，所述步骤S1中利用如下公式将将三维扫描数据通过坐标转换投影到所述近似投影平面上：

$\left(\begin{array}{c}{x}^{\′}\\ y^{\′}\\ z^{\′}\end{array}\right)=\left(\begin{array}{ccc}\frac{a_1}{A}& \frac{a_2}{A}& \frac{a_3}{A}\\ \frac{b_1}{B}& \frac{b_2}{B}& \frac{b_3}{B}\\ \frac{c_1}{C}& \frac{c_2}{C}& \frac{c_3}{C}\end{array}\right)\left(\left(\begin{array}{c}x\\ y\\ z\end{array}\right)-\left(\begin{array}{c}{x}_0\\ y_0\\ z_0\end{array}\right)\right)$

式中，(x,y,z)^T为所述三维扫描数据中一个元素的坐标，(x′,y′,z′)^T为所述近似投影平面上对应元素的坐标，

$\left(\begin{array}{c}A\\ B\\ L\end{array}\right)=\left(\begin{array}{c}\sqrt{a_1^2+a_2^2+a_3^2}\\ \sqrt{b_1^2+b_2^2+b_3^2}\\ \sqrt{c_1^2+c_2^2+c_3^2}\end{array}\right)$

其中，c＝(c₁,c₂,c₃)^T

b＝b-b×c×c

$a=b\×c=\left|\begin{array}{ccc}i& j& k\\ b_1& b_2& b_3\\ c_1& c_2& c_3\end{array}\right|=(b_2{c}_3-b_3{c}_2,b_1{c}_3-b_3{c}_1,b_1{c}_2-b_2{c}_1)$

a、b、c为所述坐标转换的坐标矢量，(x₀，y₀，z₀)为原点；

所述近似投影平面为：

c₁x+c₂y+c₃z＝0。

所述坐标转换的坐标矢量c的计算方法如下：

(1)计算所述三维扫描数据到所述近似投影平面的距离平方和；

(2)建立对该距离平方和进行最小二乘法拟合的公式，并利用高斯消元法求解公式得到坐标矢量c。

进一步地，所述步骤S1之前还包括以下步骤：

利用三维扫描仪按照预定行数和预定列数进行扫描得到所述三维扫描数据。

进一步地，所述步骤S2中形成按矩阵形式排列的数据包括以下步骤：

S21、在所述近似投影平面上的数据中找到最下面的扫描点、最上面的扫描点、最左面的扫描点以及最右面的扫描点，并以四个扫描点作为顶点形成四边形；

S22、以所述四边形的底边作为一条边，以最上面的扫描点和最右面的扫描点作为顶点形成一个矩形；

S23、利用所述矩形对对应的所述近似投影平面上的数据进行剪裁，去掉所述矩形之外的数据；

S24、计算所述步骤S23得到的数据的宽度和高度；其中所述宽度为所述矩形内数据的列数，所述高度为所述矩形内数据的行数；

S25、对于所述矩形内每一个数据，根据其三维坐标确定其在所述按矩阵形式排列的数据中的行和列，并建立其深度信息与其在所述按矩阵形式排列的数据中位置的英映射关系。

所述步骤S25中，对于所述矩形内每一个数据，其在所述按矩阵形式排列的数据中的行、列以及所述映射关系利用下面公式计算：

$i=height*\frac{y_m-y_{\min }}{y_{\max }-y_{\min }+\mathrm{\ϵ}}$

$j=width*\frac{x_m-x_{\min }}{x_{\max }-x_{\min }+\mathrm{\ϵ}}$

D(i,j)＝z'_m

m＝1,...,numpoint

式中，i表示行，j表示列，D(i，j)＝z_m为所述映射关系，m为所述矩形内数据总个数，(x_m，y_m，z_m)为数据在所述近似投影平面上的位置坐标。

所述步骤S23还包括确定列数的子步骤，具体为：

对于所述矩形内数据中两个相邻的数据点，计算其在X方向上的距离是否大于列距离阈值，若是则所述两个相邻的数据点为两列。

所述步骤S25之后还包括去除所述图像矩阵中空缺点步骤：

以所述空缺点为中心，计算其周边多个数据点的深度值，并以计算得到的多个深度值的平均值作为所述空缺点的深度值。

所述步骤S2中所述局部坐标系的坐标轴为：

$e^x={(\frac{(x_2-x_1)}{w},\frac{(y_2-y_1)}{w})}^T$

$e^y={(\frac{(x_3-x_2)}{h},\frac{(y_3-y_2)}{h})}^T$

式中，w为所述矩形的宽，h为所述矩形的高，(x₁,y₁)、(x₂,y₂)、(x₃,y₃)、(x₄,y₄)为所述矩形的额四个顶点

所述局部坐标系的变量为：

$\mathrm{\ξ}=\frac{2(x_2-x_1)(x-x_0)}{w}+\frac{2(y_2-y_1)(y-y_0)}{w}$

$\mathrm{\η}=\frac{2(x_3-x_2)(x-x_0)}{h}+\frac{2(y_3-y_2)(y-y_0)}{h}$

其中，

$w=\sqrt{{(x_2-x_1)}^2+{(y_2-y_1)}^2}=\sqrt{{(x_3-x_4)}^2+{(y_3-y_4)}^2}$

$h=\sqrt{{(x_3-x_2)}^2+{(y_3-y_2)}^2}=\sqrt{{(x_4-x_1)}^2+{(y_4-y_1)}^2}.$

进一步地，所述步骤S3中，利用下面公式计算局部坐标系下数据在所述拓印基准面上的坐标：

式中，(X_k,Y_k,Z_k){k＝0:8}为插值参量，为三维扫描数据在对应局部坐标系中的坐标值。

所述插值参量的计算包括以下步骤：

(1)计算局部坐标系下数据到所述拓印基准面的距离平方和；

(2)利用最小二乘法，按照所述差值参量使所述平方和最小计算得到所述插值参量。

进一步地，所述近似曲面为所述三维数据的二次函数，或者为三维数据的二次差值函数。

下面结合具体实施例进行说明。

本实施例的方法包括以下步骤：

【步骤一】将三维扫描数据投影到近似投影平面上：三维扫描得到的原始数据(即三维扫描数据)是一组未被表达的3D点云，如图2所示。需要将三维数据通过相应的坐标转换投影到一个近似投影平面上。坐标转换包括以下步骤：

设原始数据具有m个数据点，其中心点为：

$x_0=\left(\underset{i=0}{\overset{m-1}{\Σ}}{x}_i\right)/m$

$y_0=\left(\underset{i=0}{\overset{m-1}{\Σ}}{y}_i\right)/m---\left(1\right)$

$z_0=\left(\underset{i=0}{\overset{m-1}{\Σ}}{z}_i\right)/m$

式中，(x_i、y_i、z_i)为原始点云数据的坐标。

设原始坐标系中的坐标为(x,y,z)^T，变换坐标系中的坐标为(x′,y′,z′)^T，转换坐标系的单位坐标矢量和原点分别为：

$e^x={(e_1^x,e_2^x,e_3^x)}^T$

$\begin{array}{c}{e}^y={(e_1^y,e_2^y,e_3^y)}^T\\ e^z={(e_1^z,e_2^z,e_3^z)}^T\end{array}---\left(2\right)$

u＝(x₀,y₀,z₀)^T

则转换后的数据坐标与原始数据空间坐标的转换关系可以表示为：

$\left(\begin{array}{c}{x}^{\′}\\ y^{\′}\\ z^{\′}\end{array}\right)=\left(\begin{array}{ccc}{e}_1^x& e_2^x& e_3^x\\ e_1^y& e_2^y& e_2^y\\ e_1^z& e_2^z& e_3^z\end{array}\right)\left(\left(\begin{array}{c}x\\ y\\ z\end{array}\right)-\left(\begin{array}{c}{x}_0\\ y_0\\ z_0\end{array}\right)\right)---\left(3\right)$

做如下定义：

$e^x={(\frac{a_1}{\sqrt{a_1^2+a_2^2+a_3^2}},\frac{a_2}{\sqrt{a_1^2+a_2^2+a_3^2}},\frac{a_3}{\sqrt{a_1^2+a_2^2+a_3^2}})}^T$

$e^y={(\frac{b_1}{\sqrt{b_1^2+b_2^2+b_3^2}},\frac{b_2}{\sqrt{b_1^2+b_2^2+b_3^2}},\frac{b_3}{\sqrt{b_1^2+b_2^2+b_3^2}})}^T---\left(4\right)$

$e^z={(\frac{c_1}{\sqrt{c_1^2+c_2^2+c_3^2}},\frac{c_2}{\sqrt{c_1^2+c_2^2+c_3^2}},\frac{c_3}{\sqrt{c_1^2+c_2^2+c_3^2}})}^T$

u＝(x₀,y₀,z₀)^T

其中，e^z为近似投影平面的法向量。

若设

$\left(\begin{array}{c}A\\ B\\ L\end{array}\right)=\left(\begin{array}{c}\sqrt{a_1^2+a_2^2+a_3^2}\\ \sqrt{b_1^2+b_2^2+b_3^2}\\ \sqrt{c_1^2+c_2^2+c_3^2}\end{array}\right)---\left(5\right)$

则转换方程为：

$\left(\begin{array}{c}{x}^{\′}\\ y^{\′}\\ z^{\′}\end{array}\right)=\left(\begin{array}{ccc}\frac{a_1}{A}& \frac{a_2}{A}& \frac{a_3}{A}\\ \frac{b_1}{B}& \frac{b_2}{B}& \frac{b_3}{B}\\ \frac{c_1}{C}& \frac{c_2}{C}& \frac{c_3}{C}\end{array}\right)\left(\left(\begin{array}{c}x\\ y\\ z\end{array}\right)-\left(\begin{array}{c}{x}_0\\ y_0\\ z_0\end{array}\right)\right)---\left(6\right)$

设近似投影平面方程为：

c₁x+c₂y+c₃z＝0(7)

则原始数据中的点到近似投影平面的距离平方和为：

$\Π=\underset{i=1}{\overset{m-1}{\Σ}}{z}^{\′2}=\underset{i=0}{\overset{m-1}{\Σ}}{\left(c_1\right(x_i-x_0)+c_2(y_i-y_0)+c_3(z_i-z_0\left)\right)}^2---\left(8\right)$

对该距离平方和进行最小二乘法拟合，公式为：

$\frac{\∂\Π}{\∂c_1}=\underset{i=0}{\overset{m-1}{\Σ}}2\left(c_1\right(x_i-x_0)+c_2(y_i-y_0)+c_3(z_i-z_0\left)\right)(x_i-x_0)=0$

$\frac{\∂\Π}{\∂c_2}=\underset{i=0}{\overset{m-1}{\Σ}}2\left(c_1\right(x_i-x_0)+c_2(y_i-y_0)+c_3(z_i-z_0\left)\right)(y_i-y_0)=0---\left(9\right)$

$\frac{\∂\Π}{\∂c_3}=\underset{i=1}{\overset{m-1}{\Σ}}2\left(c_1\right(x_i-x_0)+c_2(y_i-y_0)+c_3(z_i-z_0\left)\right)(z_i-z_0)=0$

$\underset{i=0}{\overset{m-1}{\Σ}}(x_i-x_0)(x_i-x_0)c_1+\underset{i=0}{\overset{m-1}{\Σ}}(x_i-x_0)(y_i-y_0)c_2+\underset{i=0}{\overset{m-1}{\Σ}}(x_i-x_0)(z_i-z_0)e_3=0$

$\underset{i=0}{\overset{m-1}{\Σ}}(y_i-y_0)(x_i-x_0)c_1+\underset{i=0}{\overset{m-1}{\Σ}}(y_i-y_0)(y_i-y_0)c_2+\underset{i=0}{\overset{m-1}{\Σ}}(y_i-y_0)(z_i-z_0)c_3=0---\left(10\right)$

$\underset{i=0}{\overset{m-1}{\Σ}}(z_i-z_0)(x_i-x_0)c_1+\underset{i=0}{\overset{m-1}{\Σ}}(z_i-z_0)(y_i-y_0)c_2+\underset{i=0}{\overset{m-1}{\Σ}}(z_i-z_0)(z_i-z_0)c_3=0$

$\underset{i=0}{\overset{m-1}{\Σ}}(x_i-x_0)(x_i-x_0)c_1+\underset{i=0}{\overset{m-1}{\Σ}}(x_i-x_0)(y_i-y_0)c_2+\underset{i=0}{\overset{m-1}{\Σ}}(x_i-x_0)(z_i-z_0)c_3=0$

$\underset{i=0}{\overset{m-1}{\Σ}}(y_i-y_0)(x_i-x_0)c_1+\underset{i=0}{\overset{m-1}{\Σ}}(y_i-y_0)(y_i-y_0)c_2+\underset{i=0}{\overset{m-1}{\Σ}}(y_i-y_0)(z_i-z_0)c_3=0$

$\underset{i=0}{\overset{m-1}{\Σ}}(z_i-z_0)(x_i-x_0)c_1+\underset{i=0}{\overset{m-1}{\Σ}}(z_i-z_0)(y_i-y_0)c_2+\underset{i=0}{\overset{m-1}{\Σ}}(z_i-z_0)(z_i-z_0)c_3=0---\left(11\right)$

整理化简得到求解的目标方程为：

$\left(\begin{array}{ccc}\underset{0}{\overset{m-1}{\Σ}}{(x_i-x_0)}^2& \underset{0}{\overset{m-1}{\Σ}}(x_i-x_0)(y_i-y_0)& \underset{0}{\overset{m-1}{\Σ}}(x_i-x_0)(z_i-z_0)\\ \underset{0}{\overset{m-1}{\Σ}}(x_i-x_0)(y_i-y_0)& \underset{0}{\overset{m-1}{\Σ}}{(y_i-y_0)}^2& \underset{0}{\overset{m-1}{\Σ}}(y_i-y_0)(z_i-z_0)\\ \underset{0}{\overset{m-1}{\Σ}}(x_i-x_0)(z_i-z_0)& \underset{0}{\overset{m-1}{\Σ}}(y_i-y_0)(z_i-z_0)& \underset{0}{\overset{m-1}{\Σ}}{(z_i-z_0)}^2\end{array}\right)\left(\begin{array}{c}{c}_1\\ c_2\\ c_3\end{array}\right)=\left(\begin{array}{c}0\\ 0\\ 0\end{array}\right)---\left(12\right)$

使用高斯消元法求解该方程，并做归一化处理得到：c＝(c₁,c₂,c₃)^T。此即为近似投影平面的法向量，也是转换坐标系的z'轴方向矢量。

考虑三维扫描仪的特性是按列扫描的，y′轴方向矢量的求解只要简单地用前N个数据计算优势方向即为y′的方向。具体做法是：

DoitoN

$b=b+\frac{x_{i+1}-x_i}{\left|\right|x_{i+1}-x_i\left|\right|}$

EndDo(13)

temp＝b·c

b＝b-temp×c

$b=\frac{b}{\left|\right|b\left|\right|}$

x_i+1、x_i为第i+1、i个原始数据的横坐标，利用公式i取值从1到N依次循环累计算得到b，作为b的初始值，再利用公式temp＝b·c、b＝b-temp×c、计算得到作为坐标矢量的b。本实施例中将b进行了单一化处理。

x′轴方向矢量可以通过已经建立的两个坐标轴方向矢量叉乘获得：

$a=b\×c=\left|\begin{array}{ccc}i& j& k\\ b_1& b_2& b_3\\ c_1& c_2& c_3\end{array}\right|=(b_2{c}_3-b_3{c}_2,b_1{c}_3-b_3{c}_1,b_1{c}_2-b_2{c}_1)---\left(14\right)$

进一步归一化：

本实施例中，对坐标矢量均进行了归一化处理。至此，转换后的坐标系的单位坐标矢量和原点(x₀、y₀、z₀)都已确定，原始数据的近似投影平面也唯一确定。根据公式(6)的坐标转换关系，可以计算出转换后的数据空间坐标。图3展示了经过坐标转换后的点云数据在近似投影平面上的二维投影图像。

【步骤二】在近似投影平面上建立局部坐标系。扫描得到的三维扫描数据既可以拼装成一个完整的三维立体模型，也可以进行分区域处理，重构其局部信息。一般地说，任何钟鼎碑碣都会以不同的区域铭刻不同的文字，在每一个区域的数据，可以假设处在同一个曲面或平面上。这样的曲面或平面数据量很大，大部分情况要通过多次采集才能获得区域内的全部三维信息。每次采集的数据可以称为一帧数据。在一帧数据内，采集设备是按固定的行、列数量采集的，比如M×N个采集点。这种规则排列的数据，要划分为有限个单元(K×L)，只需均匀分割一帧数据，使得每个单元内的数据量为(m＝M/K、n＝N/L)个即可。但由【步骤一】得到的近似投影平面其内部或边界可能有空隙，需要对点云数据进行均匀分布处理，使帧数据规则排列形成严格的二维矩阵形式，并且矩阵元素的下标编号与其空间坐标一一对应。下面对具体步骤进行说明：

首先需要根据原始数据构建广义矩形。如果全部数据在空间某个角度观察时表现为一个平均意义上的近似矩形形状，就可以在这个平面上构造一个包含全部数据点的(面积)最小矩形，这个矩形可以作为平面坐标外框架。构建矩形的方法是：

1、在近似投影平面上的数据中，依次找到最下面的点A、最右面的点B、最上面的点C和最左面的点D，以ABCD为顶点建立一个四边形。

2、以AB作为矩形的一条边，以C或D作为矩形的顶点，建立一个矩形。其顶点分别是(x₁,y₁)、(x₂,y₂)、(x₃,y₃)、(x₄,y₄)。

3、利用得到的矩形对数据进行剪裁，丢弃矩形之外的数据。

其次，确定数据帧的宽度和高度，这可以由数据点矩阵排列的行数和列数确定。由于三维扫描仪是按列扫描得到的数据，因此只需要判断前后两个数据点在y’方向上的距离是否大于某个阈值(即距离阈值)。若超过该阈值，则增加统计的列数值，否则说明这两个数据点在同一列。当所有矩形内的数据点扫描完毕，则列数就是该帧数据的宽度，同时统计过程中数据帧具有最多数据点数那一列的点数就是数据帧的高度。

根据计算出的数据帧的宽度和高度，建立近似投影平面上数据坐标到二维矩阵(图像矩阵)下标的映射，设近似投影平面上数据点坐标为(x_m,y_m,z'_m)，二维矩阵的每个元素表示为D(i,j)，共有numpoint个数据点，则坐标映射公式为：

$i=height*\frac{y_m-y_{\min }}{y_{\max }-y_{\min }+\mathrm{\ϵ}}$

$j=width*\frac{x_m-x_{\min }}{x_{\max }-x_{\min }+\mathrm{\ϵ}}---\left(15\right)$

D(i,j)＝z'm

m＝1,...,numpoint

通过上面的额计算得到的图像矩阵中存在空隙，即有些矩阵元素的值是没有意义的深度值，本实施例使用均值滤波算法进行修正：以空缺点为中心，计算周边八个数据点有意义的深度值的平均值作为该点的深度值，如果周围的数据点都是空缺点，则该点也是空缺点。经过上述处理，得到了按严格矩阵排列的规则排列图像矩阵，并且图像矩阵的每个矩阵元素的下标编号都与数据点在近似投影平面上的坐标具有对应关系。

在此基础上，为减少运算量并提高局部精度，同一个区域还可以进一步划分为若干个(准矩形)单元。如果单元(即分区)的尺寸较小，在单元内部的底面可以近似的看作一个平面，则处在单元上的字模曲面可以用平面函数表达。将图像矩阵划分为若干个单元表达后，之后为每一个单元建立局部坐标：

单元内的3D数据可以按横纵两个方向进行局部编号。

{D_ij(x,y,z):i＝0:m-1；j＝0:n-1}(16)

局部坐标影射为：

{D(ξ_i,η_j)(x,y,z)|-1≤ξ＜1,-1≤η＜1}

其中， ${\mathrm{\ξ}}_i=\frac{i-m/2}{m/2};{\mathrm{\η}}_j=\frac{i-n/2}{n/2}---\left(17\right)$

或： $\begin{array}{c}{x}_{ij}=x({\mathrm{\ξ}}_i,{\mathrm{\η}}_j)\\ \{y_{ij}=y({\mathrm{\ξ}}_i,{\mathrm{\η}}_j)|i=0:m-1;j=0:n-1\}\\ z_{ij}=z({\mathrm{\ξ}}_i,{\mathrm{\η}}_j)\end{array}---\left(18\right)$

考虑帧数据存储格式特点，比较普遍的存储方式是一位数组格式将图像矩阵中数据进行存储。

在M×N个数据组成的一维数组D’中，在第k×l(0<k<M；0<l<N)个点上，截取其中m×n(0<m<M-k；0<n<N-l)阵列数据D的方式是：

{D(x_i,y_i,z_i)＝D′(x_j,y_j,z_j),j＝((k+i/m)×M+l+mod(i/m))|i＝0:m×n-1}

(19)

局部坐标下的可以表示为：

{D(ξ_i,η_i)(x,y,z)|-1≤ξ＜1,-1≤η＜1}

其中， ${\mathrm{\&xi;}}_i=\frac{\mathrm{mod}(i/m)-m/2}{m/2};{\mathrm{\&eta;}}_i=\frac{i/m-n/2}{n/2}---\left(20\right)$

或： $\begin{array}{c}{x}_i=x({\mathrm{\&xi;}}_i,{\mathrm{\&eta;}}_i)\\ \{y_i=y({\mathrm{\&xi;}}_i,{\mathrm{\&eta;}}_i)|i=0:m\&times;n-1\}\\ z_i=z({\mathrm{\&xi;}}_i,{\mathrm{\&eta;}}_i)\end{array}---\left(21\right)$

针对图像矩阵的矩阵形式，以矩形(x₁,y₁)(x₂,y₂)(x₃,y₃)(x₄,y₄)为单元，构建局部坐标系。为计算方便，将四个节点按照左下角、右下角、右上角、左上角的顺序排列。

首先确定局部坐标系的原点：

$u={(x_0,y_0,z_0)}^T={(\frac{x_1+x_2+x_3+x_4}{4},\frac{y_1+y_2+y_3+y_4}{4},0)}^T---\left(22\right)$

矩形的宽和高为：

$\begin{array}{c}w=\sqrt{{(x_2-x_1)}^2+{(y_2-y_1)}^2}=\sqrt{{(x_3-x_4)}^2+{(y_3-y_4)}^2}\\ h=\sqrt{{(x_3-x_2)}^2+{(y_3-y_2)}^2}=\sqrt{{(x_4-x_1)}^2+{(y_4-y_1)}^2}\end{array}---\left(23\right)$

则局部坐标系坐标轴为：

$\begin{array}{c}{e}^x={(e_1^x,e_2^x)}^T={(\frac{(x_2-x_1)}{w},\frac{(y_2-y_1)}{w})}^T\\ e^y={(e_1^y,e_2^y)}^T={(\frac{(x_3-x_2)}{h},\frac{(y_3-y_2)}{h})}^T\end{array}---\left(24\right)$

局部变量定义为：

$\begin{array}{c}\mathrm{\&xi;}=(x-x_0,y-y_0)e^x=\frac{2(x_2-x_1)(x-x_0)}{w}+\frac{2(y_2-y_1)(y-y_0)}{w}\\ \mathrm{\&eta;}=(x-x_0,y-y_0)e^y=\frac{2(x_3-x_2)(x-x_0)}{h}+\frac{2(y_3-y_2)(y-y_0)}{h}\end{array}---\left(25\right)$

这样，原始的三维局部数据就变换成了在局部坐标表达下的二维数据。利用同样的方法可将近似投影平面上的数据投影到对应的局部坐标系中。

【步骤三】拓印基准面的构造与数据切割。从空间角度考虑，从底面到字模之间，往往存在一个过渡层。在提取字模信息时，应该把过渡层信息去除，所采取的办法是不以底面为切割曲面，而是采用曲面拟合的方式拟合出钟鼎碑刻表面的近似曲面，并将拟合曲面从底面高度处上浮一定距离(通过点到近似曲面的距离实现的)形成拓印基准面，在拓印基准面的高度上完成数据切割(切割其实就是以距离为分割项，只提取上浮平面的一侧)，同时也完成了底面消除。

如果不考虑尺寸影响，钟鼎碑刻表面可以看作是比较平滑的曲面，可以用三维空间的二次函数表达。为了保证算法的普适性，也可以统一使用二次插值函数描述。

完成了数据局部化后(即步骤二)，三维扫描数据在单元内部的插值可以使用9节点等参数单元插值(以下为曲面拟合的步骤)。9个插值节点的分布，单元节点编号及其局部坐标如6图所示。相应节点的形函数可以表示为：

设在单元节点上的三维数据为(X_k,Y_k,Z_k){k＝0:8}，则单元内任意一点的三维数据为：

其中(X_k,Y_k,Z_k){k＝0:8}为插值参量，差值参量通过特定的算法求解。在几何意义上，为单元节点上的三维数据。

在选定的单元内部，三维数据集合可以表示为：

(x_ij,y_ij,z_ij){i＝0:m；j＝0:n}(28)

使用最小二乘法求解插值参量：

三维数据集合上的数据点到插值曲面的距离平方和为：

$\&Pi;=\underset{i=0,j=0}{\overset{m,n}{\&Sigma;}}\{{(x_{ij}-x({\mathrm{\&xi;}}_i,{\mathrm{\&eta;}}_j\left)\right)}^2+{(y_{ij}-y({\mathrm{\&xi;}}_i,{\mathrm{\&eta;}}_j\left)\right)}^2+{(z_{ij}-z({\mathrm{\&xi;}}_i,{\mathrm{\&eta;}}_j\left)\right)}^2\}---\left(29\right)$

根据最小二乘法，最优的插值参量使上述平方和最小，即：

$\frac{\&part;\&Pi;}{X_k}=0$

$\begin{array}{cc}\frac{\&part;\&Pi;}{Y_k}=0& k=0:8\end{array}---\left(30\right)$

$\frac{\&part;\&Pi;}{Y_k}=0$

经过适当演算，建立关于插值参量的方程：

AX＝B(31)

其中：

X＝(X_k,Y_k,Z_k)^T{k＝0:8}

求解方程(31)，可以得到参量X的数值。利用方程(27)，就得到了字模曲面在本单元上的拓印基准面的方程表达式。

一个单元求解完成后，其相邻单元相邻节点上的数值可以采用已经求解的数值，也可以另外求解。如果是每个单元之间单独求解，则可能出现单元之间细微的差异，在用近似解拼接图形时，就可能出现一定的“裂缝”。因此，在求解相邻单元时，应该采用递次求解的办法，每一个新单元只求解未知的节点数据。就会避免裂缝现象的出现。图4显示了经过曲面拟合和底面消除处理后的字模曲面投影图像。

【步骤四】噪点滤除与字模投影。本步骤将经过步骤三处理后得到的拓印基准面投影为平面字帖。在本步骤中，通过反低通滤波来消除拓印基准面的曲率变化对投影的影响以及降噪处理。图5显示了经过反低通滤除底面噪点和畸变影响后的平面字帖图像。

本文的方法，在不降低原始数据精度的前提下，对全局的三维扫描数据在不转化为面绘制模型数据的基础上直接进行处理，提高了处理精度。早处理过程方面，只通过参数控制计算的精度，避免了多余的人工参与对数据的干扰，降低操作者在三维模型编辑处理方面的准备知识。同时本发明的方法可以在短时间内去除钟碑刻底面部分数据，得到较好的分离结果，能对任意扫描仪的通用.asc格式钟碑刻点云文件直接进行处理。

以上实施方式仅用于说明本发明，而非对本发明的限制。尽管参照实施例对本发明进行了详细说明，本领域的普通技术人员应当理解，对本发明的技术方案进行各种组合、修改或者等同替换，都不脱离本发明技术方案的精神和范围，均应涵盖在本发明的权利要求范围当中。

Claims (12)

1.一种三维扫描数据的自适应底面消除方法，其特征在于，所述方法包括以下步骤：

S1、将三维扫描数据通过坐标转换投影到近似投影平面上；

S2、根据所述近似投影平面上的数据形成按矩阵形式排列的数据；

将所述数据进行区域划分，并为每一个分区数据建立局部坐标系，将每个分区的数据在对应的局部坐标系进行转换；

S3、针对每一个分区，利用差值参量，根据其局部坐标系中每个数据的坐标值计算在拓印基准面上对应的坐标值，实现将所述步骤S2处理得到的数据投影到拓印基准面，完成数据分割；

S4、对所述步骤S3处理得到的所述拓印基准面上的数据投影到平面上，之后进行去噪。

2.根据权利要求1所述的方法，其特征在于，所述步骤S4还包括以下步骤：遍历投影深度进行反低通滤波，其中所述投影深度为对应数据点的深度信息。

3.根据权利要求1所述的方法，其特征在于，所述步骤S1中利用如下公式将三维扫描数据通过坐标转换投影到所述近似投影平面上：

$\left(\begin{array}{c}{x}^{\&prime;}\\ y^{\&prime;}\\ z^{\&prime;}\end{array}\right)=\left(\begin{array}{ccc}\frac{a_1}{A}& \frac{a_2}{A}& \frac{a_3}{A}\\ \frac{b_1}{B}& \frac{b_2}{B}& \frac{b_3}{B}\\ \frac{c_1}{C}& \frac{c_2}{C}& \frac{c_3}{C}\end{array}\right)(\left(\begin{array}{c}x\\ y\\ z\end{array}\right)-\left(\begin{array}{c}{x}_0\\ y_0\\ z_0\end{array}\right))$

式中，(x,y,z)^T为所述三维扫描数据中一个元素的坐标，(x′,y′,z′)^T为所述近似投影平面上对应元素的坐标，

$\left(\begin{array}{c}A\\ B\\ C\end{array}\right)=\left(\begin{array}{c}\sqrt{a_1^2+a_2^2+a_3^2}\\ \sqrt{b_1^2+b_2^2+b_3^2}\\ \sqrt{c_1^2+c_2^2+c_3^2}\end{array}\right)$

其中，c＝(c₁,c₂,c₃)^T

b＝b-b×c×c

$a=b\&times;c=\left|\begin{array}{ccc}i& j& k\\ b_1& b_2& b_3\\ c_1& c_2& c_3\end{array}\right|=(b_2{c}_3-b_3{c}_2,b_1{c}_3-b_3{c}_1,b_1{c}_2-b_2{c}_1)$

a、b、c为所述坐标转换的单位坐标矢量，(x₀，y₀，z₀)为原点；

所述近似投影平面为：

c₁x+c₂y+c₃z＝0。

4.根据权利要求1所述的方法，其特征在于，所述坐标转换的单位坐标矢量c的计算方法如下：

计算所述三维扫描数据到所述近似投影平面的距离平方和；

建立对该距离平方和进行最小二乘法拟合的公式，并利用高斯消元法求解公式得到坐标矢量c。

5.根据权利要求1所述的方法，其特征在于，所述步骤S1之前还包括以下步骤：

利用三维扫描仪按照预定行数和预定列数进行扫描得到所述三维扫描数据。

6.根据权利要求5所述的方法，其特征在于，所述步骤S2中形成按矩阵形式排列的数据包括以下步骤：

S21、在所述近似投影平面上的数据中找到最下面的扫描点、最上面的扫描点、最左面的扫描点以及最右面的扫描点，并以四个扫描点作为顶点形成四边形；

S22、以所述四边形的底边作为一条边，以最上面的扫描点和最右面的扫描点作为顶点形成一个矩形；

S23、利用所述矩形对对应的所述近似投影平面上的数据进行剪裁，去掉所述矩形之外的数据；

S24、计算所述步骤S23得到的数据的宽度和高度；其中所述宽度为所述矩形内数据的列数，所述高度为所述矩形内数据的行数；

S25、对于所述矩形内每一个数据，根据其三维坐标确定其在所述按矩阵形式排列的数据中的行和列，并建立其深度信息与其在所述按矩阵形式排列的数据中位置的英映射关系。

7.根据权利要求6所述的方法，其特征在于，所述步骤S25中，对于所述矩形内每一个数据，其在所述按矩阵形式排列的数据中的行、列以及所述映射关系利用下面公式计算：

$i=height*\frac{y_m-y_{\min }}{y_{\max }-y_{\min }+\mathrm{\&epsiv;}}$

$j=width*\frac{x_m-x_{\min }}{x_{\max }-x_{\min }+\mathrm{\&epsiv;}}$

D(i,j)＝z'_m

m＝1,...,numpoint

式中，i表示行，j表示列，D(i，j)＝z_m为所述映射关系，m为所述矩形内数据总个数，(x_m，y_m，z_m)为数据在所述近似投影平面上的位置坐标。

8.根据权利要求6所述的方法，其特征在于，所述步骤S23还包括确定列数的子步骤，具体为：

对于所述矩形内数据中两个相邻的数据点，计算其在X方向上的距离是否大于列距离阈值，若是则所述两个相邻的数据点为两列。

9.根据权利要求6所述的方法，其特征在于，所述步骤S25之后还包括去除所述图像矩阵中空缺点步骤：

以所述空缺点为中心，计算其周边多个数据点的深度值，并以计算得到的多个深度值的平均值作为所述空缺点的深度值。

10.根据权利要求6所述的方法，其特征在于，所述步骤S2中所述局部坐标系的坐标轴为：

$e^x={(\frac{(x_2-x_1)}{w},\frac{(y_2-y_1)}{w})}^T$

$e^y={(\frac{(x_3-x_2)}{h},\frac{(y_3-y_2)}{h})}^T$

式中，w为所述矩形的宽，h为所述矩形的高，(x₁,y₁)、(x₂,y₂)、(x₃,y₃)、(x₄,y₄)为所述矩形的四个顶点

所述局部坐标系的变量为：

$\mathrm{\&xi;}=\frac{2(x_2-x_1)(x-x_0)}{w}+\frac{2(y_2-y_1)(y-y_0)}{w}$

$\mathrm{\&eta;}=\frac{2(x_3-x_2)(x-x_0)}{h}+\frac{2(y_3-y_2)(y-y_0)}{h}$

其中，

$w=\sqrt{{(x_2-x_1)}^2+{(y_2-y_1)}^2}=\sqrt{{(x_3-x_4)}^2+{(y_3-y_4)}^2}$

$h=\sqrt{{(x_3-x_2)}^2+{(y_3-y_2)}^2}=\sqrt{{(x_4-x_1)}^2+{(y_4-y_1)}^2}.$

11.根据权利要求1所述的方法，其特征在于，所述步骤S3中，利用下面公式计算局部坐标系下数据在所述拓印基准面上的坐标：

式中，(X_k,Y_k,Z_k){k＝0:8}为所述插值参量，为三维扫描数据在对应局部坐标系中的坐标值。

12.根据权利要求11所述的方法，其特征在于，所述插值参量的计算包括以下步骤：

计算局部坐标系下数据到所述拓印基准面的距离平方和；

利用最小二乘法，按照所述差值参量使所述平方和最小计算得到所述插值参量。