一种地面GPS辅助的LiDAR数据高程差的误差修正方法与流程

本发明属于测绘科学与技术领域，涉及一种地面GPS辅助的机载激光点云误差修正方法，主要应用于高精度DEM/DSM生产、城市三维建模、高精度地貌特征的拾取等领域。

背景技术：

机载LiDAR(Light Detection And Ranging)系统是一种快速采集地表信息的新兴测量技术，在硬质开阔地面的高程精度可高达15cm，平面精度优于30cm，常用来大面积地获取密集且精确的地形和地物数据。目前已广泛地应用于电力巡线、公路选线、森林参数估计、3D城市建模等领域，机载LiDAR系统是由激光扫描仪、GPS/IMU(Inertial MeasurementUnit)惯性导航装置、计算机控制导航系统、存储设备等硬件集成在一起的一个复杂的测距系统，因此，许多系统内在的因素制约着LiDAR点云的平面和高程精度。

在LiDAR点云数据获取过程中，常存在多种误差源，如激光扫描仪测距误差、扫描镜控制误差、激光扫描仪与POS设备时间同步误差、激光扫描仪与POS设备的偏心以及安置角误差、POS观测值系统误差等。在LiDAR数据处理的实际生产中发现，点云的高程误差一般在0.3m至1.0m之间，目前还没有系统的修正点云数据高程误差的方法。

技术实现要素：

本发明的目的在于提供一种地面GPS辅助的LiDAR数据高程差的误差修正方法，可以有效的对机载LiDAR数据误差进行修正。

为达到上述目的，本发明采用了以下技术方案：

步骤1获取包含针对同一目标区的K个航带的激光雷达点云数据，作为原始点云数据，其中K≥3，对原始点云数据进行滤波，剔除植被和建筑物，只留下有效的地表点云数据；

步骤2格网化处理：用四边形网格对经过步骤1处理后的点云数据进行划分；

步骤3经过格网化处理后，将K组点云数据转换到一个绝对坐标下，分别对K组点云数据，在形状为长方形的航带上，在长度和宽度上每间隔拾取一个点，将这个点的高程值eij作为矩阵元素值，其中，i、j分别是点云取点的行数和列数，这样，可以由K组点云数据分别组装成一个点云高程值矩阵，即，可以组装出第一组，第二组，第三组，……第K组点云数据的点云高程值矩阵E1,E2,E3,…,EK；

步骤4计算K组点云高程值矩阵的均值高程值矩阵Eμ，然后分别计算各组点云高程值矩阵El(l＝1,2,……，K)与均值高程值矩阵Eμ的绝对差参数ηl，对ηl进行排序，即则第p1组为最佳点云数据；

步骤5利用高精度GPS在目标区内静态测量一定数量的观测点，观测点的数量根据被测的地形进行选择，地势越平坦则需要布置的观测点越多，且相邻观测点的间距通常不应超过500米。从各个GPS观测点测得的数据值中分别提取高程值，组装成GPS高程值矩阵GGPS，同时从最佳点云数据中，提取与GPS高程值矩阵GGPS平面坐标对应的观测点的点云数据高程值，组装成最优点云数据高程值参数矩阵Qp；

步骤6计算GPS高程值矩阵GGPS与最优点云数据高程值参数矩阵Qp元素差值的均值μp，将Qp的各个元素加上μp，推算出新的参数矩阵Q′p；计算GPS高程值矩阵GGPS与参数矩阵Q′p元素差值的方差σp；

步骤7若σp≤σMAX，则该组点云数据修正后符合误差要求，可以使用；若σp>σMAX，则说明此组点云数据不符合使用要求，应进行舍弃，以下一组点云数据(依步骤4给出的排序)为最优点云数据，转步骤8；

步骤8，重复步骤6、7，直至选出可以使用的首组点云数据(即最合适的点云数据)。

假设第pi组为最合适的点云数据，则将计算出的GPS高程值矩阵GGPS与第pi组点云数据高程值参数矩阵元素差值的均值以代数计算的形式累加到第pi组点云数据每一个激光脚点的Z向三维坐标数据上。

本发明具有如下有益效果：

本发明通过多次反复扫描被测对象，在地面使用GPS高精度定位等方式，使机载LiDAR点云数据中的高程值误差得到有效估计和预算，从而对DEM数据以及以后三维模型的绘制打下了良好的基础。本发明采用GPS数据与点云数据的误差绝对值以及方差大小作为判断点云数据是否合格的判别条件，通过计算参数对点云数据进行修正，使得修正后的点云数据更加的符合实际。

附图说明

图1是本发明实施例的流程图；

图2是GPS地面观测点在航带上的布置形式，图2中：1为GPS布置位点。

具体实施方式

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

本发明的实施步骤如下，参见图1：

步骤1，采用机载激光雷达由先至后对测区进行K个架次飞行后，获得包含K个航带的激光雷达点云数据，为防止无效点云数据组干扰，其中K≥3；每一个航带的激光雷达点云数据中均包括多个激光脚点的测量数据，且每一个激光脚点的测量数据均包括该激光脚点的(X，Y，Z)三维坐标数据和测量数据获取时间；

对原始点云数据进行滤波，剔除植被和建筑物，只留下有效的地表点云数据，数据包括该激光脚点的(X，Y，Z)三维坐标数据和测量数据获取时间；

步骤2，规则格网化处理：采用数据处理器对所处理点云数据进行格网划分，即用多个尺寸相同的方形网格对所述点云数据进行划分，且划分后所形成格网的行数和列数分别为：

其中，YMAX为所处理点云数据中Y轴坐标值最大的激光脚点的Y轴坐标值，YMIN为所处理点云数据中Y轴坐标值最小的激光脚点的Y轴坐标值，XMAX为所处理点云数据中X轴坐标值最大的激光脚点的X轴坐标值，XMIN为所处理点云数据中X轴坐标值最小的激光脚点的X轴坐标值；L为所述方形网格的边长，a为所处理点云数据的密度；LMAX为划分网格的行数，RMAX划分网格的列数；

步骤3，分别对K组点云数据，在形状为长方形的航带上，在长度和宽度上每隔拾取一个点，将这个点的高程值eij作为矩阵元素值，其中，i、j分别是点云取点的行数和列数。这样，可以组装成一个高程值矩阵：

其中，

区域内的地形采样点可以任意标记在网格的顶点、边上或者网格内，如果地形采样点在网格的边上或者网格内，可以通过线性插值的方法获得该采样点的高程值。这样航带上的每一次测量都可以生成一个高程值矩阵；

步骤4，将K组点云数据转换到一个绝对坐标下，按照步骤3所述的方法，可以组装出第二组，第三组，……第K组点云数据的点云高程值矩阵E2,E3,…,EK，将K组点云数据转换到一个绝对坐标下，可以保证不同点云高程值矩阵的元素El(i,j)对应点的平面一样，不会出现平面偏差；

步骤5，计算K组点云高程值矩阵的均值高程值矩阵Eμ，作为航带点云数据是否合格的判别参数：

分别计算各组点云高程值矩阵El与均值高程值矩阵Eμ的绝对差参数ηl：

式中，ηl为第l组点云高程值矩阵El与均值高程值矩阵Eμ之差的绝对值，Eμ(i,j)为均值高程值矩阵Eμ第i行第j列的元素，El(i,j)为第l组点云高程值矩阵El第i行第j列的元素；

比较ηl值的大小，对ηl进行排序，即则第p1组为最佳点云数据；

步骤6，利用测量精度可达mm级的高精度GPS在测区内静态测量20个观测点。每个观测点的持续测量时间＞100min，有效观测卫星总数＞5个。这些观测点原则上均匀分布在被测带宽内呈四边形的20个地点上，参见图2；

为了避免地形起伏引起的2次高程误差的影响，选取的GPS观测位置可以根据实际情况进行微调，尽量为视野开阔，地势平坦，周边地表无茂密植被。每一个GPS点的测量数据均包括观测点的(X，Y，Z)三维坐标数据和测量数据获取时间；为保证测量精度，可多次测量取平均值；

为了使GPS测量数据与点云高程值矩阵元素值进行直接对比，应使得GPS测量点的平面坐标在点云绝对坐标内呈整数值，这样可以使得GPS测量点与点云高程值矩阵某一元素的平面坐标重合，有利于直接对比高程值；

从20个GPS观测点测得的数据值中提取高程值，组装成地面GPS测量的固定点高程值矩阵GGPS：

从第p1组点云数据中，提取与地面GPS测量的固定点高程值矩阵GGPS平面坐标对应的观测点的点云数据高程值，组装成最优点云数据高程值参数矩阵

步骤7，计算地面GPS测量的固定点高程值矩阵GGPS与最优点云数据高程值参数矩阵元素差值的均值

式中，为地面GPS测量的固定点高程值矩阵GGPS与最优点云数据高程值参数矩阵元素差值的均值，GGPS(i,j)为地面GPS测量的固定点高程值矩阵GGPS第i行第j列的元素，为最优点云数据高程值参数矩阵第i行第j列的元素；

将的各个元素加上地面GPS测量的固定点高程值矩阵GGPS与最优点云数据高程值参数矩阵元素差值的均值得到新的元素值即：

这样可以得到新的参数矩阵

计算地面GPS测量的固定点高程值矩阵GGPS与最优点云数据高程值参数矩阵元素差值的方差

式中，为地面GPS测量的固定点高程值矩阵GGPS与最优点云数据高程值参数矩阵元素差值的方差，为最优点云数据高程值参数矩阵第i行第j列的元素；

步骤8，若则修正后的点云数据符合误差要求，可以使用；若则说明此组点云数据不符合使用要求，应进行舍弃，重新依据步骤5给出的点云数据排序选取下一组点云数据(例如，P1的下一组为P2，依次类推)为最佳点云数据；重新选取最佳点云数据后，按照步骤6、7进行处理，直至选出最合适的一组点云数据，假设最合适的点云数据为第pi组，则将计算出的以代数计算的形式累加到第pi组点云数据每一个激光脚点的Z向三维坐标数据上，即可得到修正后的点云数据。

σMAX的取值取决于对测区精度的要求，精度要求越高，则σMAX越小，σMAX可根据实际情况进行选取。