CN112229385A - 一种用3个三维平面基元直接求解扫描点云地理化参数的方法 - Google Patents

Abstract

一种用3个三维平面基元直接求解扫描点云地理化参数的方法，包括如下步骤：步骤1：以三维平面转换方程为基础方程，以扫描仪坐标系中三维平面基元的姿态参数和位置参数为观测值，求解未知的地理化参数。由于平面转换方程中地理化的姿态参数和位置参数没有同时在一个方程出现，因此解算面基元的姿态参数的方程为第一组，以面基元的姿态参数为观测值，解算旋转矩阵中三个独立参数；其余方程为第二组，解算地理化的位置参数。步骤2：获取地理化参数的方差阵，建立地理化参数误差模型。本发明可得到完整点云地理化数学模型、地理化参数及其方差阵，为研究和利用地理化后点云质量评估提供支持。

Description

一种用3个三维平面基元直接求解扫描点云地理化参数的 方法

技术领域

本发明涉及三维平面不同坐标系之间的转换，用3个三维平面基元直接求解扫描点云地理化参数及其精度评定的方法，属于工程测量领域。

背景技术

人们总是习惯性的将客观对象放在以铅垂线为基准线、以北方向为基本方向的地理空间中去观察、表达和研究，这种地理空间坐标系就是测绘领域的地理参考坐标系，也称为工程测量坐标系。不同传感器平台获取的点云，需要将其坐标系转换到统一的地理坐标系中，这个过程在点云处理理论中称为点云地理化，其关键是地理化参数求解。学术界针对点云地理化提出了大量点云地理化方法和地理化精度的评定方法。直接求解法在工程测量领域很常见，它是用最少的已知点直接计算测站仪器中心在地理坐标系中的3维位置坐标和姿态的方法。在三维激光扫描测量中，用于解算扫描仪在地理坐标系中6个地理化参数的基本元素称为基元，按其几何特征可分为点基元、线基元和面基元。点基元就是传统意义的控制点，三维激光扫描技术出现后，线基元和面基元作为点云拼接/地理化的基本单元，已然成为地理化参数/坐标转换参数的基元。然而用线基元和面基元解算地理化参数的方法还没有完善，例如没有线基元和面基元的地理化参数直接求解方法，本专利就是发明一种用3个三维平面基元求解扫描点云地理化参数的方法，本发明的意义在于：

(1)在受地形限制的某些情况下进行设站扫描，扫描视场内没有已知点基元、或只已知一个或两个点基元，此时无法解算6个地理化参数，可寻求扫描物体上的已经过地理化的点云与本站点云公共的三维平面，来解算地理化参数。

(2)全站仪是一种高精度三维坐标测量仪器，用于无合作目标法测定已知平面特征参数，本专利方法更容易使扫描仪与全站仪进行联合测量。

(3)直接解算地理化参数的结果可作为平差法计算6个地理化参数时的可靠初始值。

(4)能获取地理化参数的方差阵，使建立地理化参数的误差模型更容易。

在建筑物等目标的扫描点云中有大量的平面对象，所以本专利在三维激光扫描的数据处理和应用中有重要意义。

发明内容

本发明的目的是提供用3个三维平面基元直接求解扫描点云地理化参数的方法，包括如下步骤：

步骤1：以三维平面转换方程为基础方程，以扫描仪坐标系中三维平面基元的姿态参数和位置参数为观测值，求解未知的地理化参数。由于平面转换方程中地理化的姿态参数和位置参数没有同时在一个方程出现，因此解算面基元的姿态参数的方程为第一组，以面基元的姿态参数为观测值，解算旋转矩阵中三个独立参数；其余方程为第二组，解算地理化的位置参数。在三维笛卡尔o-xyz坐标系内(如扫描仪坐标系)，三维平面E分3类，用于表达物体的正面、侧面和水平面，其方程为

式中E_x代表第一类不平行于x轴的三维平面，E_y代表第二类不平行于y轴的三维平面，E_z代表第三类不平行于z轴三维平面，a_i,b_i,c_i(i＝x,y,z)为三维平面的特征参数。已知地理化参数

称T_S＝[X_S Y_S Z_S]^T为地理化的位置参数；

ω、κ为地理化的姿态参数，由其构成旋转矩阵

为了方便，

设反对称矩阵

其元素a、b、c相互独立，R表达为

式中，I为三阶单位阵，Δ＝1+a²+b²+c²。我们已经推导出3类三维平面反转换模型为

式中Q＝[A_i,B_i,C_i]^T,i＝x,y,z是E在工程测量坐标系O-XYZ中的3类平面特征参数。

由上面3式可见，姿态相关参数求解方程与位置参数求解方程是独立的，可分组解算。以三维平面姿态观测值a_i,b_i,i＝x,y,z为第一组解算参数，以反对称矩阵的元素a、b、c为未知数，用最小二乘先解未知数a、b、c，然后计算R、

ω、κ。根据(3)-(6)式得到方程

M₁T_P＝N₁ (7)

式中

按最小二乘平差可解得

以三维平面位置观测值c_i,i＝x,y,z为第二组解算参数，用于解算地理化位置参数X_S、Y_S、Z_S，按下式计算

式中，

步骤2：按最小二乘平差方法计算单位权中误差

可得T_P的方差阵Δ_P为

根据(9)式，用方差传播律，T_S的方差阵Δ_S为

附图说明

图1是某大型城市建筑物图；

图2是传统方法和本发明方法得到的控制点地理化后三维点位误差图。

具体实施方式

下面结合实例、附图对本发明作进一步说明。

为某大型城市建筑物高精度建模而进行的三维激光扫描点云数据地理化。使用Riegl VZ-1000三维激光扫描仪，扫描时单站扫描精度为20千分度，点云抽稀后单站2500万点左右，抽稀点云分辨率2cm，共扫描18站；在建筑物附近4个标准控制点上布设了球形标靶(标准标靶)，另外在建筑物周围布设无控制点的球形标靶(检测标靶)，用于检测地理化精度。地理参考平面坐标系选用CGCS2000坐标系，实验区中央子午线为117°，高程坐标系选用黄海85高程坐标系，投影方式为高斯投影，以下统一称为工程测量坐标系。建筑物周围有4个标准控制点，只有第一站S1站能扫描到3个标准标靶和3个检测标靶，可用三维坐标转换方程进行地理化；第二站S2只能扫描到1个标准标靶和3个检测标靶，可用三维坐标转换方程进行地理化(与第一站相同)。以S2第二站的地理化说明本发明的实施。

步骤1：S1和S2的重叠平面P1、P2和P3，作为地理化三维平面基元，如图1所示，在S1地理化后的点云中提取上述3个平面，拟合平面在工程测量坐标系中的特征参数Q＝[A_i,B_i,C_i]^T,i＝x,y,z，在S2扫描点云中提取上述3个平面，拟合其在扫描仪坐标系的平面特征参数q＝[a_i,b_i,c_i]^T,i＝x,y,z，平面拟合误差均小于2mm。以三维平面姿态观测值a_i,b_i,i＝x,y,z为第一组解算参数，以反对称矩阵的元素a、b、c为未知数，用最小二乘先解未知数a、b、c，然后计算R、

ω、κ。根据(3)-(6)式得到方程

M₁T_P＝N₁ (7)

式中

按最小二乘平差可解得

以三维平面位置观测值c_i,i＝x,y,z为第二组解算参数，用于解算地理化位置参数X_S、Y_S、Z_S，按下式计算

式中，

步骤2：按最小二乘平差方法计算单位权中误差

可得T_P的方差阵Δ_P为

根据(9)式，用方差传播律，T_S的方差阵Δ_S为

用本专利方法和传统标靶ICP点云拼接方法相比较，本专利有如下优势：

1)解决只有少量标准控制点基元时点云的地理化问题。传统地理化方法每站都需要不少于3个标准控制点，而本专利方法借助物体三组平面进行地理化，可以有效减少标准控制点的要求，提高效率。

2)本专利方法地理化后三维点位误差小3mm，如图2所示，高于传统标靶ICP点云拼接方法。

3)本专利地理化参数解算方法的精度除了与扫描仪本身的精度有关，还与平面的平整度有关，平整度越高，平面拟合误差越小，平面特征参数的精度就越高。

Claims (1)

1.一种用3个三维平面基元直接求解扫描点云地理化参数的方法，包括如下步骤：

步骤1：以三维平面转换方程为基础方程，以扫描仪坐标系中三维平面基元的姿态参数和位置参数为观测值，求解未知的地理化参数。由于平面转换方程中地理化的姿态参数和位置参数没有同时在一个方程出现，因此解算面基元的姿态参数的方程为第一组，以面基元的姿态参数为观测值，解算旋转矩阵中三个独立参数；其余方程为第二组，解算地理化的位置参数。在三维笛卡尔o-xyz坐标系内(如扫描仪坐标系)，三维平面E分3类，用于表达物体的正面、侧面和水平面，其方程为

式中E_x代表第一类不平行于x轴的三维平面，E_y代表第二类不平行于y轴的三维平面，E_z代表第三类不平行于z轴三维平面，a_i,b_i,c_i(i＝x,y,z)为三维平面的特征参数。已知地理化参数

称T_S＝[X_S Y_S Z_S]^T为地理化的位置参数；

ω、κ为地理化的姿态参数，由其构成旋转矩阵

为了方便，

设反对称矩阵

其元素a、b、c相互独立，R表达为

式中，I为三阶单位阵，Δ＝1+a²+b²+c²。我们已经推导出3类三维平面反转换模型为

式中Q＝[A_i,B_i,C_i]^T,i＝x,y,z是E在工程测量坐标系O-XYZ中的3类平面特征参数。

由上面3式可见，姿态相关参数求解方程与位置参数求解方程是独立的，可分组解算。以三维平面姿态观测值a_i,b_i,i＝x,y,z为第一组解算参数，以反对称矩阵的元素a、b、c为未知数，用最小二乘先解未知数a、b、c，然后计算R、

ω、κ。根据(3)-(6)式得到方程

M₁T_P＝N₁ (7)

式中

按最小二乘平差可解得

以三维平面位置观测值c_i,i＝x,y,z为第二组解算参数，用于解算地理化位置参数X_S、Y_S、Z_S，按下式计算

式中，

步骤2：按最小二乘平差方法计算单位权中误差

可得T_P的方差阵Δ_P为

根据(9)式，用方差传播律，T_S的方差阵Δ_S为

Images