一种基于点云曲面特征约束的建筑物局部变形分析方法与流程

本发明具体涉及一种基于建筑物点云曲面特征约束的局部变形分析方法。

背景技术：

近年来，为了能够满足城市居民的住房需求，不断涌现出结构复杂的大型建筑物，这也为建筑物的安全埋下了隐患。因此，需要采用快速、高精度的建筑物变形监测手段，对建筑物进行测量，及时获取最新的变形信息。传统三维测量技术主要是高精度地测量目标的某一个或多个离散定位点的三维坐标数据及该点三维特性，该技术仅能测量定位点数据，并且测量不同定位点间的简单几何尺寸。三维激光扫描数据处理技术可以重建目标模型及分析结构特性，并且进行全面的后处理及提取目标结构的复杂几何内容，可以分析目标模型的形变、位移及各部分结构之间的变化关系等。因此，将三维激光扫描技术应用到建筑物变形测量，研究基于三维点云模型的建筑物形体特征的提取原理和方法，具有重要意义。

利用三维激光扫描技术进行建筑物变形监测，首先需要运用合理的方法确定出能够正确反映建筑物变化情况的形体特征。建筑物可以由若干部位或结构组成，而点、线、面等特征又组成了建筑物中的部位与结构。建筑物中的点、线、面特征不仅具有几何约束特征，比如形状、位置、大小等，其内部之间还具有复杂的拓扑关系，而这才是空间实体依存的基础。因此，研究建筑的形体特征就不仅要对建筑物点、线、面对象的位置、形状、长度、大小、面积、体积等约束特征进行研究，还要对建筑物内部以及不同建筑物之间的几何关系特征进行研究，包括它们之间的距离关系、平行关系、垂直关系等。

因此，对于研究建筑物的形体特征及变形分析来说，确定出建筑物的点、线、面特征仅仅只是开始，因为建筑物要符合设计要求，保证建筑物能够安全地运行，其各组成部分之间以及同周围事物之间要满足一定的几何关系，这种几何关系对建筑物各组成部分以及建筑物之间起到约束和校核的作用，而点、线、面等形体特征只是这些几何关系的基础。由于建筑物主要由平面、柱面等二次曲面构成，本发明将三维激光扫描技术应用于建筑物变形监测，提出了基于建筑物点云数据曲面约束关系的局部变形分析方法。

目前，三维激光扫描测量技术进行建筑物的变形监测，已经出现了多种基于点云数据的变形量分析计算方法，主要包括：点对点的直接比较分析法，点云到平面的变形分析法，平面之间的比较分析法，点云数据的特征量比较分析法。对于点对点的变形分析方法，不确定是否能在两个不同的时期采样到完全相同的点，并且从不同角度扫描测量，物体上的点密度可能会有所不同，因此直接的点对点比较是不合适的。基于点到平面或平面之间变形量计算方法，虽然将面特征参数应用到变形监测中，但只是用到了平面特征参数，不能满足大型复杂建筑物的变形监测的需要。基于建筑物整体点云数据特征提取方法，比如主成分分析法，只是利用数学计算方法对整体点云数据的主成分向量进行提取，该方法仍然没有考虑建筑物各曲面之间的约束关系，因此，也不能对建筑物各部位之间的局部空间变形进行判断。除此之外，目前基于建筑物点云数据的曲面特征提取方法，大多只是利用曲率和法失向量，将具有相似特征的点云数据区分开来，来实现点云数据的分割，并没有进行点云数据曲面类型的判断识别，更没有将建筑物曲面之间的约束关系应用到建筑物变形监测中，因此，在建筑物的变形分析中就没有考虑建筑物各部位之间的相对变形，从而会降低建筑物变形监测的精度和可靠性。

技术实现要素：

发明目的：针对现有技术的不足，将建筑物曲面自动识别技术和曲面约束关系判断技术应用到建筑物变形监测中，提出了一种基于点云曲面特征约束的建筑物局部变形分析方法，以提高变形监测精度和可靠性。

发明内容：本发明所述的一种基于建筑物点云曲面特征约束的局部变形分析方法，包括以下步骤：

(1)对预先获取的建筑物密集点云数据进行预处理和数据分割，得到建筑物各表面的点云数据；

(2)对建筑物各表面点云数据进行曲面特征识别，判断建筑物表面的曲面类型；

(3)根据建筑物的曲面类型，对点云数据进行曲面拟合计算，得到建筑物各表面的曲面特征参数；

(4)建筑物曲面特征约束关系识别，即根据建筑物曲面特征参数，对建筑物曲面的尺寸约束和结构约束进行计算；

(5)根据建筑物曲面特征约束关系，判断建筑物曲面之间的空间位置关系，并分析曲面之间的局部变形。

进一步地，所述步骤(1)包括以下步骤：

(11)对获取的点云数据进行配准、噪声剔除等初步处理；

(12)利用建筑物点云层次聚类分割方法对点云数据进行分割，将相距较远又比较密集的点云数据块进行初始分割，同时，结合点云数据的光谱特征和几何特征，对经初始分割后的点云数据块进行细分，直到每个点云数据块都具有单一的几何特征，得到分割后的点云数据。

进一步地，所述步骤(2)包括以下步骤：

(21)用一般二次曲面的表达式对点云数据进行拟合，确定出各数据点到其拟合曲面的距离，统计出最小的几何距离误差；将给定的阈值与最小几何距离误差值进行比较，若阈值小于最小几何距离误差，则该曲面为nubrs曲面；否则，则转入(22)；

(22)根据点云数据微分几何量估算方法，计算点云数据中各点的法矢和曲率，若各点的法矢n为一固定数值c，并且曲率k为零或近似为零，则该点云数据区域判断为平面；否则，转入(23)；

(23)计算点云数据中各点的高斯映射和法曲率映射，若高斯映射为一维分布，并且其对应的法曲率映射ka等于常数，kb＝0，则识别为圆柱面；否则，转入(24)；

(24)若高斯映射为二维均匀分布，并且其法曲率映射|k|min＝|k|max＝c＝常数，则识别为球面。

进一步地，步骤(3)所述曲面拟合计算包括平面拟合、球面拟合、圆柱面拟合和圆锥面拟合。

进一步地，所述步骤(4)包括以下步骤：

(41)两个平面之间的约束识别：

假定已初步提取出两个平面的相关参数p1(q0,q1,q2,q3)，p2(q'0,q'1,q'2,q'3)，两平面的法矢为：n1＝(q0,q1,q2)，n2＝(q'0,q'1,q'2)，两个平面的夹角∠p1,p2＝a，且

①两个平面共面：则相应参数应满足：

其中，ε和d为阈值，ε的取值范围为0.95-1.05，d的取值小于平面拟合误差；

②两个面平行：则平面参数需要满足：

|cosa|∈ε

其中，ε取值为0.95-1.05；

③两个平面垂直：则平面参数必须满足：

|cosa|∈ε

其中，ε取值范围为0-0.05；

(42)平面与圆柱面之间的约束识别：

假设平面p的相关参数为：(q0,q1,q2,q3)，柱面cy的相关参数为：(c0,c1,c2,c3,c4,c5,c6)，则p与cy之间的约束关系包括：

①cy轴线满足在平面p上约束时有：

②p与cy之间满足相切约束时有：

③p与cy之间满足垂直约束时有：

其中，δ为阈值，δ的取值范围为0-0.05；

(43)圆柱面与圆柱面之间的约束识别：

假定两个圆柱面的参数分别为：cy1:(c0,c1,c2,c3,c4,c5,c6)，cy2:(c'0,c'1,c'2,c'3,c'4,c'5,c'6)，则两个圆柱面同轴时应该满足：

其中，δ为阈值，δ的取值范围为0-0.05；

(44)平面与圆锥面之间的约束识别：

假定锥面方程的相关参数为co：(o0,o1,o2,o3,o4,o5,o6)，则平面与锥面满足相切约束时有：

其中，δ为阈值，δ的取值范围为0-0.05；

(45)球面与圆柱面之间的约束识别：

假定球面的相关参数为s：(s0,s1,s2,s3)，并令：

则球心过柱面的轴线时应该满足：

其中，δ为阈值，δ的取值范围为0-0.05。

有益效果：与现有技术相比，本发明的有益效果：本发明可以在大量散乱的建筑物点云数据中，对曲面特征和参数进行快速有效地识别和提取，结合建筑物的曲面约关系，将建筑物曲面类型识别技术和曲面约束关系判断技术应用到建筑物变形分析中，对建筑物不同部位的局部变形情况进行判断，从多个局部角度去分析和评价建筑物的安全性，提高了建筑物变形分析的精确程度和可靠性。

附图说明

图1为本发明的流程图；

图2为点云数据的曲面特征识别流程图；

图3为球面参数化示意图；

图4为圆柱面参数化示意图；

图5为圆锥面参数化示意图；

图6为桥梁扫描的柱面点云图；

图7为柱面a-e的点云数据图；

图8为圆柱面扫描平面图；

图9为圆柱面的点云图；

图10为圆柱面拟合模型图。

具体实施方式

下面结合附图对本发明的技术方案作进一步的说明。如图1所示，本发明具体包括以下步骤：

1、对预先获取建筑物的密集点云数据进行预处理和数据分割，得到建筑物各表面的点云数据。

首先对获取的点云数据进行配准、噪声剔除等初步处理，然后，利用建筑物点云层次聚类分割方法对点云数据进行分割，将相距较远又比较密集的点云数据块进行初始分割，同时，结合点云数据的光谱特征和几何特征，对经初始分割后的点云数据块进行细分，直到每个点云数据块都具有单一的几何特征，得到分割后的点云数据。

2、对建筑物各表面点云数据进行曲面特征识别，即在对点云数据进行分割后，计算离散点云曲面的法矢和曲率特性，结合高斯映射和法曲率映射来识别建筑物表面的曲面类型。具体包括以下步骤，如图2所示。

(1)用一般二次曲面的表达式对点云数据进行拟合，确定出各数据点到其拟合曲面的距离，统计出最小的几何距离误差。将给定的阈值与最小几何距离误差值进行比较，若阈值小于最小几何距离误差，则该曲面为nubrs曲面。否则，则转入第(2)步。

(2)根据点云数据微分几何量估算方法，计算点云数据中各点的法矢和曲率，若各点的法矢(n)为一固定数值(c)，并且曲率k为零或近似为零，则该点云数据区域判断为平面。否则，转入第(3)步。

(3)计算点云数据中各点的高斯映射和法曲率映射，若高斯映射为一维分布，并且其对应的法曲率映射ka等于常数，kb＝0。，则识别为圆柱面。否则，转入第(4)步。

(4)若高斯映射为二维均匀分布，并且其法曲率映射|k|min＝|k|max＝c＝常数，则识别为球面。

3、建筑物曲面特征参数提取：根据曲面类型识别的结果，对点云数据采用最小二乘原理进行曲面拟合，得到其具体的几何参数，通过参数提取得出特征曲面。

基于点云数据的曲面特征参数提取包括两步，即二次曲面拟合和特征参数计算：

(1)二次曲面拟合

①平面拟合

假设平面方程为：

ax+by+cz+d＝0

式中，a、b、c、d为平面方程系数。

为方便构造误差方程，令上式转换为：

z＝ax+by+c

则点云数据中第i个点的误差方程为：

vi＝axi+byi+c-zi

根据最小二乘原理，可得法方程如下式所示：

式中，n为观测点个数，按上式即可计算出平面参数a、b、c。

②球面拟合

假设球面方程为：

(x-x0)²+(y-y0)²+(z-z0)²＝r²

其中，x0、y0、z0、r为球面参数，o＝(x0,y0,z0)表示球面中心，r表示球面半径。对空间点pi到球面(球心为o，半径为r)的距离为：

为了方便对球心和半径进行计算，先对球面进行参数化，如图3所示，pn为球表面上离坐标原点o0距离最近的点，pn点是连接坐标原点和球心之后与球表面的一个交点；其中n为单位向量，从原点指向球心；ρ表示n方向上从原点离球表面上最近点的代数距离；k表示点pn处的最大曲率值，球面半径可以用1/k表示。n用球坐标可以表示如下：

球面坐标可以表示为：

可以得到球面拟合的残差为：

利用最小二乘法求解残差平方和的最小值来确定球面参数。

③圆柱面拟合

假设圆柱面方程为：

[(x-x0)m-(y-y0)l]²+[(y-y0)n-(z-z0)m]²+[(z-z0)l-(x-x0)n]²＝r²

其中，x0,y0,z0,m,n,l,r为圆柱面参数，o(x0,y0,z0)表示圆柱面轴线上任意一点，a＝(m,n,l)表示圆柱面轴线方向的单位向量，r表示圆柱面半径。

如图4所示，对柱面进行参数化，pi为圆柱表面上离坐标原点o0距离最近的点，pi是原点到旋转轴的投影线与圆柱面的一个交点；a表示圆柱面的轴线方向；a和n均为单位向量，显然n为圆柱面的一个法向量，则a·n＝0；k为点pi处的最大曲率值，这里1/k可以表示圆柱面半径。

实际局部圆柱面的两端面圆心c1、c2和圆柱面半径r，这些参数的转化如下：

r＝1/k

为了得到两端面的圆心c1、c2，将所有的数据点投影到圆柱面的旋转轴上，然后由投影参数t可以确定两端点，也就是圆心c1、c2。各数据点到旋转轴上的投影参数为：

然后将得到的参数ti(i＝1,…,n)排序，最大的与最小的分别为两端点，即圆柱面两端面圆圆心c1、c2对应的投影参数，还原可得c1、c2两点：

旋转轴的方向为a，为了确定旋转轴的位置，将原点投影到旋转轴上，(实际上是沿着方向n投影)即点o'＝(x'0,y'0,z'0)：

这样就得到了有着常见的几何参数的圆柱面。

可以得到圆柱面拟合的残差为：

利用最小二乘法求解残差平方和的最小值来确定圆柱面参数。

④圆锥面拟合

假设圆锥面方程为：

[(x-x0)m-(y-y0)l]²+[(y-y0)n-(z-z0)m]²+[(z-z0)l-(x-x0)n]²＝[(x-x0)²+(y-y0)²+(z-z0)²]cos²θ

其中，x0,y0,z0,m,n,l,θ为圆锥面参数，p＝(x0,y0,z0)表示圆锥面顶点，a＝(m,n,l)表示圆锥面轴线方向的单位向量，θ表示圆锥面的半锥顶角。

如图5所示，圆锥面参数化，pi为圆锥面上离原点距离最近的点(对于圆锥面来说，是原点垂直于最近一条母线的投影线与圆锥面的一个交点)；k表示在点pi处圆锥面的最大曲率值，a表示圆柱面的旋转轴线方向，a和n均为单位向量，与圆柱面不同，在圆锥面中a和n并不垂直，因此将a也用球坐标表示为：

a＝(cosσsinτ,sinσsinτ,cosτ)

可以得到圆锥面拟合的残差为：

利用最小二乘法求解残差平方和的最小值来确定圆锥面参数。

(2)曲面特征参数提取

①平面方程：假设平面方程为

ax+by+cz+d＝0

其中a,b,c,d为平面参数，则n＝(a,b,c)为平面法向量特征。

②球面方程：假设球面方程为

(x-x0)²+(y-y0)²+(z-z0)²-r²＝0

其中：(x0,y0,z0)是球面中心；r是球面半径。该方程一般形式可表示为：

f(x,y,z)＝x²+y²+z²+ax+by+cz+d＝0

两种表达方法的系数对应关系为：

上式中(x0,y0,z0)、r等为球面对应的特征参数。

③圆柱面方程：假设柱面方程为

(x-x0)²+(y-y0)²+(z-z0)²-(nx(x-x0)+ny(y-y0)+nz(z-z0))²-r²＝0

其中：(x0,y0,z0)是轴线通过的固定点、(nx,ny,nz)是轴线方向单位矢量、r是圆柱面半径等为圆柱面特征参数。

④圆锥方程：假设圆锥面的方程为

[(x-x0)²+(y-y0)²+(z-z0)²]cos²α-(nx(x-x0)+ny(y-y0)+nz(z-z0))²＝0

其中：(x0,y0,z0)是圆锥顶点、(nx,ny,nz)是轴线方向单位矢量、α是圆锥半顶角等为圆锥面的特征参数。

根据以上方法求出建筑物曲面特征参数，还可以进一步计算对应的面积或者体积等特征信息，确定出建筑物面的几何尺寸特征。

4、建筑物曲面特征约束关系识别：根据建筑物曲面特征参数，对建筑物曲面的尺寸约束和结构约束进行计算。具体包括以下步骤：

(1)两个平面之间的约束识别

假定已初步提取出两个平面的相关参数p1(q0,q1,q2,q3)，p2(q'0,q'1,q'2,q'3)，两平面的法矢为：n1＝(q0,q1,q2)，n2＝(q'0,q'1,q'2)，两个平面的夹角∠p1,p2＝a，且

①两个平面共面：则相应参数应满足

其中，ε和d为阈值，ε的取值范围为0.95-1.05，d的取值小于平面拟合误差。

②两个面平行：则平面参数需要满足

|cosa|∈ε

其中，ε取值为0.95-1.05。

③两个平面垂直：则平面参数必须满足

|cosa|∈ε

其中，ε取值范围为0-0.05。

(2)平面与圆柱面之间的约束识别

假设平面p的相关参数为：(q0,q1,q2,q3)，柱面cy的相关参数为：(c0,c1,c2,c3,c4,c5,c6)，则p与cy之间的约束关系包括：

①cy轴线满足在平面p上约束时有：

②p与cy之间满足相切约束时有：

③p与cy之间满足垂直约束时有：

其中，δ为阈值，δ的取值范围为0-0.05。

(3)圆柱面与圆柱面之间的约束识别

假定两个圆柱面的参数分别为：cy1:(c0,c1,c2,c3,c4,c5,c6)，cy2:(c'0,c'1,c'2,c'3,c'4,c'5,c'6)，则两个圆柱面同轴时应该满足：

其中，δ为阈值，δ的取值范围为0-0.05。

(4)平面与圆锥面之间的约束识别

假定锥面方程的相关参数为co：(o0,o1,o2,o3,o4,o5,o6)，则平面与锥面满足相切约束时有：

其中，δ为阈值，δ的取值范围为0-0.05。

(5)球面与圆柱面之间的约束识别

假定球面的相关参数为s：(s0,s1,s2,s3)，并令：

则球心过柱面的轴线时应该满足：

其中，δ为阈值，δ的取值范围为0-0.05。

最后，根据建筑物曲面特征约束关系，判断建筑物曲面之间的空间位置关系，并分析各曲面之间的局部变形。

实例一：桥梁柱面变形分析计算

如图6所示，为了对桥梁柱面之间的约束关系进行检测，检查桥梁构件的施工质量，利用美国faro公司focus^3d三维激光扫描仪对桥梁进行了三维扫描测量，在桥梁的侧面选取5个面a、b、c、d、e，其对应的点云数据如图7所示。为了对桥梁施工质量进行检测，查阅桥梁设计资料知：a、b、c、d、e五个面设计为平面，并且个面之间的设计约束关系如表1所示。

表1各平面约束关系表(设计)

对获取的点云数据进行了配准、噪声剔除和数据分割等处理，得到分割后的点云数据如图7所示。

所选取的桥梁柱面点云数据数量分别为a面38410个，b面42180个，c面29124个，d面39821个，e面39090个。经计算a-e面的点云数据曲率|k|min都小于1.9362×10^-4，|k|max都小于4.3252×10^-4，即各曲面点云数据的曲率几乎为零，判断为平面。

计算出各平面的法向矢量为：

pa＝[0.509,0.861,-0.015]pb＝[0.867,-0.498,0.000]

pc＝[0.503,0.864,-0.004]pd＝[0.534,-0.290,-0.794]

pe＝[0.524,0.852,-0.018]

以a面和b面之间的约束关系计算为例：p1＝[0.509,0.861,-0.015]，p2＝[0.867,-0.498,0.000]。

则|cosa|＝0.0125

其余个平面之间的约束识别要素如表2所示。

表2各平面约束识别要素|cosa|计算表

①平面a和平面b之间的约束关系识别：

|cosaab|＝0.0125＜0.05

判断平面a和平面b相互垂直。

②平面a和平面c之间的约束关系识别：

|cosaac|＝1

判断平面a和平面b相互平行。

③平面b和平面d之间的约束关系识别：

|cosabd|＝0.6074

得：

abd＝52.59°

判断平面b和平面d的夹角为52.59度。

同理，可以识别出其余各平面之间约束关系。

由以上各平面之间约束关系的计算值知，平面a、c、e相互平行，平面b、d与平面a、c、e垂直，平面b和平面d成52.59度的夹角，这与已知各平面之间的约束关系基本一致。说明该桥梁施工完成各平面之间的约束关系同设计值相符，满足设计要求。

实例二：建筑物承重柱子同轴检测

为了对房屋柱子施工质量进行检测，选取如图8所示的房屋承重柱子，该柱子分为两部分，下半部分c2为柱子的基础部分，支撑着上半部分c1。根据房屋设计资料知，下半部分c2的半径比上半部分c1大，两个柱子的设计中心轴线在同一条直线上，即两个柱子设计同轴。如图9所示，为了对房屋柱子进行同轴检测，选取了圆柱面点云数据中的399057个数据点。

对获取的点云数据进行了配准、噪声剔除和数据分割等处理，得到分割后的点云数据，获取圆柱面拟合模型，如图10所示。

经计算其法矢的高斯映射为圆，并且其对应的曲率|k|min均值为0，|k|max的均值为1.67，其方差为0.0021，判断该区域为圆柱面。

建筑物曲面特征参数提取如表3、表4所示：

表3圆柱面几何参数计算表

表4圆柱面约束识别参数计算表

根据圆柱面曲面约束关系特征计算公式：

和

将以上数据代入公式，得：

tz0＝|c'3c4-c3c'4|＝|1.642×5.083-1.538×5.126|＝0.0077

tz1＝|c'3c5-c3c'5|＝|1.642×(-0.083)-1.538×(-0.076)|＝0.019

tz2＝|c'3c6-c3c'6|＝|1.642×3.431-1.538×3.659|＝0.0062

cz0＝||c0|-|c'0||＝0.0013

cz1＝||c1|-|c'1||＝0.0018

cz2＝||c2|-|c'2||＝0.000

分析以上计算数据，可知tz0、tz1、tz2cz0、cz1、cz2的值小于0.05(拟合误差)，说明圆柱c1和c2的轴线方向差异很小，即两轴线基本同向，并且两圆柱具有共同的距离远点最近的点，因此，两圆柱同轴。

根据实例一和实例二的数据处理结果可以看出，利用本发明中建筑物曲面约束关系识别算法可以很好地判断建筑物中各种曲面之间的约束关系，为后续建筑物形体结构的判断和变形分析提供理论依据。