发明内容
为了解决现有技术存在的上述问题,本发明目的在于提供一种基于单木分割的地面和机载LiDAR点云配准方法。
本发明所采用的技术方案为:
一种基于单木分割的地面和机载LiDAR点云配准方法包括以下步骤:
步骤一、将获取的地面LiDAR点云数据和机载LiDAR点云数据进行LiDAR点云预处理;
步骤二、对LiDAR点云预处理所得的三维激光点云数据进行高程投影转换,获取二维高程投影影像;
步骤三、利用Hough变换圆检测法从二维高程投影图像中检测出单木数据;
步骤四、利用种子区域增长法提取出二维高程投影图像中的单木数据,得单木点图像;
步骤五、将步骤四所得的单木点图像转换成对应的三维点云数据,转换后得地面单木点云数据和机载单木点云数据;
步骤六、基于单木中心的ICP配准参数计算;
步骤七、将步骤五所得的地面单木点云数据和机载单木点云数据进行配准,并对配准结果进行精度评价。
进一步地,步骤一中,LiDAR点云预处理包括点云去噪处理和点云格式转换处理。
进一步地,步骤二中,高程投影转换中,将三维激光点云数据在二维XOY平面的覆盖范围内以一定的采样间隔d划分为规则的网格;
每个网格的高程值为落在对应网格范围内的所有点云数据的最大高程值,然后将得到的最大高程值按照公式(1)变换成灰度值Gi;
其中,[Xmin,Xmax]和[Ymin,Ymax]分别为二维XOY平面的斜对角的端值,Gi为格网的灰度值,Zi为格网的高程值,Zmax和Zmin分别为点云的最大高程值和最小高程值,Width为影像宽度,影像灰度范围为[0,255]。
进一步地,步骤三中,检测出的单木数据包括单木中心的二维坐标和包含单木的圆形区域数据集。
进一步地,步骤四中,以步骤三中检测得到的单木中心的二维坐标为种子区域增长的起始点,设置区域增长灰度阈值和区域增长方式,提取出二维高程投影图像中的单木数据,得单木点图像。
进一步地,步骤五中,将步骤四中提取到的单木点图像根据每个像素的坐标,从像素高程值对应的点云的三维坐标中找到与之对应的三维坐标。
进一步地,步骤六中,将步骤五所得的地面单木点云数据和机载单木点云数据进行ICP配准,计算ICP配准参数。
进一步地,步骤七中,以配准后的机载单木点云数据、迭代过程中的误差和精细配准变换的矩阵作为精度评价的指标。
本发明的有益效果为:
1)在点云配准过程中,采用了基于单木分割的ICP算法进行精细配准,减少了用于计算配准参数的点云数量,在保证精度的前提下,大大提高了算法的效率。
2)在单木分割过程中首先采用了点云高程投影生成灰度图像的方法,这一方法不仅提高了点云的可读性,而且为后续处理提供方便。
3)使用霍夫变换圆检测法和种子区域增长法检测和提取出高程投影图像中的单木,快速并且准确。
具体实施方式
下面结合附图及具体实施例对本发明做进一步阐释。
本实施例提供一种基于单木分割的地面和机载LiDAR点云配准方法,包括以下步骤:
步骤一、将获取的地面LiDAR点云数据和机载LiDAR点云数据进行LiDAR点云预处理,LiDAR点云预处理包括点云去噪处理和点云格式转换处理。具体包括以下步骤:
步骤1.1,点云去噪处理:将地面LiDAR点云数据和机载LiDAR点云数据中的噪声点去除。
步骤1.2,点云格式转换处理:本申请采用的原始数据是同一地区(含有多棵单木)的地面LiDAR点云数据和机载LiDAR点云数据,数据格式为.las文件格式。
因获取的地面LiDAR点云数据和机载LiDAR点云数据均为.las文件格式,将.las文件格式转换为只包含点云三维坐标(X,Y,Z)信息的文本文件格式(删除多余的字段,滤除其中冗余信息,只保留点云的三维坐标),最后保存成txt格式便于程序读取。
步骤二、对LiDAR点云预处理所得的三维激光点云数据进行高程投影转换,获取二维高程投影影像。
高程投影转换中,将三维激光点云数据在二维XOY平面的覆盖范围内以一定的采样间隔d划分为规则的网格;
每个网格的高程值为落在对应网格范围内的所有点云数据的最大高程值,然后将得到的最大高程值按照公式(1)变换成灰度值Gi;从而得到一幅点云高程投影图像。点云高程投影不仅提高了点云数据内容的易读性,而且为后续处理提供了便利。
其中,[Xmin,Xmax]和[Ymin,Ymax]分别为二维XOY平面的斜对角的端值,Gi为格网的灰度值,Zi为格网的高程值,Zmax和Zmin分别为点云的最大高程值和最小高程值,Width为影像宽度,影像灰度范围为[0,255]。
高程投影转换中,将地面LiDAR点云数据和机载LiDAR点云数据分别转换到高程投影空间中,得到一组二维高程投影影像,该过程使用到的特征指标有:
①采样间隔
二维高程投影影像为二维影像,采样间隔d由影像宽度Width和点云范围[Xmin,Xmax]确定。
在本实施例中,采样间隔d取0.25米。
②投影影像灰度范围
在本实施例中,投影影像灰度范围[a,b]取[0,255]。
步骤三、利用Hough变换圆检测法从二维高程投影图像中检测出单木数据;检测出的单木数据包括单木中心的二维坐标和包含单木的圆形区域数据集。具体包括以下步骤:
步骤3.1,单木中心检测,利用改进后的Hough变换圆检测法来实现圆形检测,通过边缘检测、计算圆周法向量的交点等实现单木中心检测,得单木中心的二维坐标。
步骤3.2,单木圆检测,由单木中心的二维坐标推导出单木圆半径,得包含单木的圆形区域数据集。
步骤3.2.1,计算某一个圆心到所有圆周线的距离,该圆心所对应的圆的半径的值包含在这些距离中。
步骤3.2.2,设定两个阈值,分别定义为最大半径和最小半径,然后使某一个圆心到所有圆周线的距离在最大半径和最小半径之间,设置阈值的目的是为了避免半径过大或者过小,因为投影图像中单木的树冠覆盖面半径往往都是在一定范围内的。
步骤3.2.3,对筛选出来的距离进行排序。
步骤3.2.4,找到距离相同的那些值,并计算出相同值的数量。
步骤3.2.5,设定一个阈值,只有在距离相同值的数量大于该阈值时,才认为该距离是该圆心对应的圆半径。
步骤3.2.6,重复步骤3.2.1-步骤3.2.5,直至找到所有圆心对应的半径。
利用Hough变换圆检测法从二维高程投影图像中检测出单木数据,该过程使用到的特征指标有:
①圆心之间的最小距离
为避免在同一棵单木中检测到多个单木中心,需设置单木圆心之间的最小距离,本实施例中,设置的圆心之间的最小距离为10个像素。
②Canny算子边缘检测高阈值
本实施例中,此阈值为50个灰度值。
③圆心检测和半径推导的共有阈值
本实施例中,此阈值为11个像素。
④最大半径和最小半径
为了保证检测到的圆形都为单木,最大半径和最小半径分别为4个像素和10个像素。
⑤圆内像素灰度均值阈值
为了避免一些错误检测,圆内像素灰度均值阈值设置为100。
步骤四、利用种子区域增长法提取出二维高程投影图像中的单木数据,得单木点图像;,提取出二维高程投影图像中的单木数据,得单木点图像。具体包括以下步骤:
步骤4.1,以步骤三中检测得到的单木中心的二维坐标为种子区域增长的起始点;
步骤4.2,设置合适的区域增长灰度阈值和区域增长方式;
步骤4.3,比较种子点与周围点的灰度差值,把所有满足约束的点包含进来,并选择新的种子点。
步骤4.4,以新的种子点进行区域增长,直至再也没有满足条件的像素可以包括进来。
步骤4.5,选择新的起始种子点,重复步骤4.1~4.4,直到提取出所有的单木。
步骤4.2中使用到的特征指标有:
①区域增长方式
为保证能够提取出完整的单木,采用八邻域区域增长的方式。
②区域增长灰度阈值
为避免提取到非单木点和保证能够提取出完整的单木,区域增长灰度阈值设置为50。
步骤五、将步骤四所得的单木点图像转换成对应的三维点云数据,转换后得地面单木点云数据和机载单木点云数据;具体包括以下步骤:
步骤5.1,根据高程投影记录下每个网格的高程值对应的点云的三维坐标。
步骤5.2,将提取到的单木点图像根据每个像素的坐标,从像素高程值对应的点云的三维坐标中找到与之对应的三维坐标,从而实现二维灰度图像到三维点云数据的还原。
步骤六、基于单木中心的ICP配准参数计算;将步骤五所得的地面单木点云数据和机载单木点云数据进行ICP配准,计算ICP配准参数。具体包括以下步骤:
步骤6.1,对于目标点云中的每个点,对应匹配参考点云中的最近点;
步骤6.2,求解使目标点云和对应的参考点云的均方根(RMS)最小的刚体变换,计算旋转参数和平移参数;
步骤6.3,使用上一步获得的变换矩阵来变换目标点云;
步骤6.4,迭代(重新关联点),直到满足迭代终止条件(可以是迭代次数或者是误差小于某一给定的阈值)。这里的误差最小:可以是相邻两次均方根差的绝对值小于某一限差。
步骤6.5,得到最终的ICP配准参数。
步骤6.4中,设置的迭代终止条件为:迭代次数不大于20次和相邻两次均方根差的绝对值小于0.001。
步骤七、将步骤五所得的地面单木点云数据和机载单木点云数据进行配准,并对配准结果进行精度评价。具体包括以下步骤:
步骤7.1,以初始配准后的机载单木点云数据为目标点云数据,以地面单木点云数据为参考点云数据,以用初始配准参数进行坐标变换后得到的机载单木点云数据的最高点为控制点,以步骤六得到的ICP配准参数进行点云精细配准。
步骤7.2,配准精度评价,以配准后的机载单木点云数据、迭代过程中的误差和精细配准变换的矩阵作为精度评价的指标。
本发明不局限于上述可选的实施方式,任何人在本发明的启示下都可得出其他各种形式的产品。上述具体实施方式不应理解成对本发明的保护范围的限制,本发明的保护范围应当以权利要求书中界定的为准,并且说明书可以用于解释权利要求书。