基于地物形态特征的行道树点云渐进提取方法

1.本发明涉及车载lidar数据处理技术领域，尤其涉及一种基于地物形态特征的行道树点云渐进提取方法。


背景技术：

2.车载lidar技术可以高效、精确地获取地表以及建筑物等人工地物的高密度三维离散点云数据，并广泛应用于三维重建、城市规划与管理、车辆导航、无人驾驶、灾害应急与评估等方面，其中点云分类提取起着关键作用。基于三维激光扫描数据的树木及属性信息提取研究已经日趋成熟，但在对大面积场景下进行数据扫描时，车载lidar系统具有相对的弱势。因此，有很多关于行道树的研究已经逐渐转移到车载激光扫描数据处理上。
3.martin针对点云中树木点云之间互不遮挡的这种特殊情况，对单株行道树的属性进行参数化，最后利用这些参数对行道树进行了建模；yu等人利用单点标记过程获得树木大致位置及几何形状,已标记点云的后验分布采用贝叶斯模型进行获得，同时采用可逆跳转马尔可夫蒙特卡洛算法(r.jmcmc)模拟进行实验，最终得到树木的最佳位置，具有一定实用性；王东等人提出结合随机森林和区域分割的方法实现行道树的提取，引入随机森林算法提取行道树的点云并通过dbscan聚类和属性信息提取方法提取出单株行道和行道树的属性信息，该方法具有适应性，但是其精度还有待进一步提高；张焱等人引入随机森林算法，通过多维特征向量构建、特征选择、随机森林模型建立和随机森林分类等步骤从地物点云中准确分离出行道树点云，但在面对较复杂的环境时，就会有相对的局限性。上述行道树点云的提取算法没有较好的普适性，没有考虑到部分复杂环境下存在有相距树干很近的低矮物难以剔除的情况，影响对行道树点云高精度的获得。


技术实现要素：

4.本发明要解决的技术问题是针对上述现有技术的不足，提供一种基于地物形态特征的行道树点云渐进提取方法。
5.对提取行道树点云环境的考虑可以有效提高算法的普适性，直接影响行道树点云的提取精度。因此本文通过提高提取算法的普适性和提取精度：对各类地物的形态特征进行深入分析，根据面积、高差、投影、形态四个方面的形态特征渐进提取行道树；对初始提取的行道树进行优化，采用行道树点云二值化方法对树冠以下点云进行优化，提取小面积对象，得到最后的行道树点云数据，在很大程度上提高了行道树点云的提取精度。
6.为解决上述技术问题，本发明提供一种基于地物形态特征的行道树点云渐进提取方法，包括以下步骤：
7.步骤1：采用kd树方法(k-dimensional tree)对原始车载lidar点云数据进行粗差剔除，对去噪后的点云采用布料模拟滤波(cloth simulation filtering，csf)得到地面点云和非地面点云。
8.步骤2：通过对栅格图像八邻域聚类进而实现离散的非地面点云的聚类分割，得到
n个点云块，作为数据处理单元。
9.步骤2.1：非地面点云栅格化；
10.步骤2.1.1：依据地物点云的覆盖范围建立xoy二维平面格网，公式如下所示，投影点云到格网平面；
[0011][0012][0013]
式中，x
max
、x
min
、y
max
、y
min
为点云覆盖范围内x方向和y方向上的最大值和最小值；gridsize为格网单元大小，gridsize＝0.45米；m、n分别为x方向和y方向上的格网单元的数目；int(
·
)为求整函数。
[0014]
步骤2.1.2：将每个格网内点云数据的最低高程值赋给每个格网，对空格网赋值为0，实现地物点云数据的栅格化；
[0015]
步骤2.1.3：对每个格网内的点云进行记录存储；
[0016]
步骤2.2：统计每个网格中点的数目n
p
，判断其是否大于给定的数目阈值tn，如果n
p
》tn，则将该格网保留，记为非空格网，反之，则将包含点的数目小于或等于数目阈值tn的格网视为空格网，赋值为0；
[0017]
步骤2.3：将对应栅格影像上非空格网处的矩阵值赋值为-1，设定索引变量i，并设置初始值i＝1；
[0018]
步骤2.4：随机选取a中任一非空格网为初始中心格网，赋值i给中心格网；
[0019]
步骤2.5：探测中心格网八邻域属性值为-1的格网，将其属性值更新为i，步骤2.4和该步骤中处理过的非空格网后续不再进行探测；
[0020]
步骤2.6：对于新探测到的非0格网，继续八邻域探测，重复该过程直至找不到非空格网，至此一个聚类单元处理完毕；
[0021]
步骤2.7：令i＝i+1，重复步骤2.4到步骤2.6，，直到所有格网点云完成八邻域探测，得到n个独立的聚类单元。
[0022]
步骤3：提取每一个聚类单元的特征，具体包括高程相关特征、投影面积、投影跨度、投影面积比及其它特征等。
[0023]
对每一个聚类单元分别进行如下特征提取：
[0024]
(1)聚类单元高程相关特征包括格网化点云的归一化高程、高程最大值、高程最小值和高程差；
[0025]
(2)投影面积sj为单元中的点云投影到二维xoy平面上的面积；
[0026]
(3)投影跨度k
x
和ky分别为每一个聚类单元在格网x方向和y方向上投影的最大跨度；
[0027]
其中，x方向跨度为投影单元对应的x方向上格网编号的最大值与最小值的差值；y方向跨度为投影单元对应的y方向上单元格网编号的最大值和最小值的差值；
[0028]
(4)投影面积比sr为每一个聚类单元中整体点云在xoy平面上的投影面积s1和部分点云在xoy平面投影面积s2的比值；
[0029]
(5)每一个聚类单元的面积sm以及包含的点云个数，其中面积sm为每个聚类单元所
占平面格网的数目。
[0030]
步骤4：基于对各类地物形态特征分析的基础上，依据面积、高差、投影、形态四个方面的形态特征渐进提取初始行道树点云。
[0031]
步骤4.1：对各类地物的形态特征进行分析；
[0032]
步骤4.1.1：地面信息存在于整个环境中，局部范围内基本表现为平面，在整个环境中则是有起伏的，且地面点的高程是最小的；
[0033]
步骤4.1.2：建筑物分布于道路的两侧，垂直于地面，建筑物的高程在场景信息中最高，其在投影平面的形状为长条形状并且其截面的投影形状与整个建筑物的投影形状相差不大，投影面积大于其它地物；
[0034]
步骤4.1.3：行道分布在道路与建筑物之间，从高程来看的话，其高低于建筑物、路灯和电线杆等地物，高于其它的地物，行道树的投影形状呈现为簇状，投影面积比建筑物投影面积小，大于其它地物的投影面积，且行道树的树干部分的投影面积与整体的投影面积相差较大；
[0035]
步骤4.1.4：车辆分布在路面或者道路两侧，高程低于行道树、建筑物、路灯等地物，高于低矮植被，投影面积较小；
[0036]
步骤4.1.5：路分布在行道树之间或者道路两侧，高程低于建筑物，高于其它地物，为5～6米，投影形状呈不规则形状，具有很好的辨识性，投影面积较小，通常不大于1～2平方米；
[0037]
步骤4.1.6：低矮植被分布在行道树的周围，高程高于地面，低于其它地物，其最大高程同树干点的最高点高程一样，投影形状为不规则形状，投影面积较小。
[0038]
步骤4.2：二次噪声筛选：统计每个单元的面积sm(m＝1
…
n)，此处面积为聚类单元所占格网的数目，同时设定面积阈值sy，判断聚类单元的面积是否大于面积阈值，如果存在聚类单元面积小于面积阈值sy，还需计算聚类单元中包含的点云数目，当点数小于数目阈值tm时，则判定为二次噪声，进行滤除；
[0039]
步骤4.3：低矮地物筛选：计算的最大高差值，设置高差阈值th，如果最大高差值小于高差阈值th，则将该聚类单元判定为低矮地物，并删除掉该聚类单元；如果最大高差值大于设置的高差阈值th，则保存该聚类单元；
[0040]
步骤4.4：建筑物筛选：设定x和y方向上的跨度阈值t
kx
和t
ky
，如果聚类单元在x和y方向上的跨度小于给定的跨度阈值t
kx
和t
ky
时，则将聚类单元视为杆状地物，反之则为建筑物；
[0041]
步骤4.5：杆状地物筛选：设置阈值根据聚类单元的投影面积设置面积阈值s
p
，如果某个聚类单元的投影面积大于所设置的面积阈值s
p
，则判定该聚类单元为树木；若投影面积小于阈值s
p
，则判定为其他杆状地物。
[0042]
步骤5：采用二值化算法对初始行道树点云进行优化，得到最终的行道树结果；
[0043]
步骤5.1：选择初始行道树点云数据中1.5米以下的点云数据，对其格网化，格网大小为0.1米，记录每个格网内的点云数据；
[0044]
步骤5.2：统计并判断所有格网是否为空，若格网为空值则赋值为0；若格网内有点云数据，则给该格网赋值为1，生成二值化栅格图像；
[0045]
步骤5.3：对二值图像进行聚类，设定索引变量j并初始j＝1；
[0046]
步骤5.4：随机选二值栅格图像中任一非0格网为初始中心格网，赋值j给中心格网；
[0047]
步骤5.5：探测中心格网八邻域属性值为非0的格网，将非0格网的属性值更新为j，步骤5.4和该步骤中处理过的非空格网后续不再进行探测；
[0048]
步骤5.6：对于新探测到的非0格网，继续八邻域探测，重复该过程直至找不到非空格网，至此一个聚类单元处理完毕；
[0049]
步骤5.7：令j＝j+1，重复步骤5.4到步骤5.6，直到所有格网完成八邻域探测，得到m个独立的聚类单元。
[0050]
步骤5.8：判断上述结果中每个聚类单元面积是否小于给定面积阈值t0，设置t0＝0.01平方米，对应格网单元的数目为10，如果聚类单元面积小于t0，则认为其对应的点云为树干点云；如果聚类单元面积大于t0，则认为其对应的点云为非行道树点云，对其进行剔除。
[0051]
采用上述技术方法所产生的有益效果在于：
[0052]
(1)本发明提供了基于地物形态特征的行道树点云渐进提取方法，为行道树点云高精度的提取提供了技术支持，在深入分析常见各类地物形态特征的基础上，充分利用地物形态特征对行道树点云进行快速高效的提取；
[0053]
(2)本发明在基于地物形态特征渐进提取行道树点云后，依据点云二值化算法，提取出小面积区域对象，实现行道树点云优化，使行道树点云提取精度极大程度上地提高。
附图说明
[0054]
图1为本发明具体实施方式中基于地物形态特征的行道树点云渐进提取方法流程图；
[0055]
图2为本发明具体实施方式中点云块效果图；
[0056]
图3为本发明具体实施方式中常见地物投影形态图；
[0057]
图4为本发明具体实施方式中实验数据的实验区域图；
[0058]
图5为本发明具体实施方式中实验数据的实验结果图。
具体实施方式
[0059]
下面结合附图和实施例，对本发明的具体实施方式作进一步详细描述。以下实施例用于说明本发明，但不用来限制本发明的范围。
[0060]
本发明提供一种基于地物形态特征的行道树点云渐进提取方法，如图1所示，包括如下步骤：
[0061]
步骤1：采用kd树方法对原始车载lidar点云数据进行粗差剔除，对去噪后的点云采用布料模拟滤波得到地面点云和非地面点云。
[0062]
步骤2：通过对栅格图像八邻域聚类进而实现离散的非地面点云的聚类分割，得到n个点云块，作为数据处理单元。
[0063]
步骤2.1：非地面点云栅格化；
[0064]
步骤2.1.1：依据地物点云的覆盖范围建立xoy二维平面格网，公式如下所示，投影点云到格网平面；
[0065][0066][0067]
式中，x
max
、x
min
、y
max
、y
min
为点云覆盖范围内x方向和y方向上的最大值和最小值；gridsize为格网单元大小，gridsize＝0.45米；m、n分别为x方向和y方向上的格网单元的数目；int(
·
)为求整函数。
[0068]
步骤2.1.2：将每个格网内点云数据的最低高程值赋给每个格网，对空格网赋值为0，实现地物点云数据的栅格化；
[0069]
步骤2.1.3：对每个格网内的点云进行记录存储；
[0070]
步骤2.2：统计每个网格中点的数目n
p
，判断其是否大于给定的数目阈值tn，本发明设置为3，如果n
p
》tn，则将该格网保留，记为非空格网，反之，则将包含点的数目小于或等于数目阈值tn的格网视为空格网，赋值为0；
[0071]
步骤2.3：将对应栅格影像上非空格网处的矩阵值赋值为-1，设定索引变量i，并设置初始值i＝1；
[0072]
步骤2.4：随机选取a中任一非空格网为初始中心格网，赋值i给中心格网；
[0073]
步骤2.5：探测中心格网八邻域属性值为-1的格网，将其属性值更新为i，步骤2.4和该步骤中处理过的非空格网后续不再进行探测；
[0074]
步骤2.6：对于新探测到的非0格网，继续八邻域探测，重复该过程直至找不到非空格网，至此一个聚类单元处理完毕；
[0075]
步骤2.7：令i＝i+1，重复步骤2.4到步骤2.6，，直到所有格网点云完成八邻域探测，得到n个独立的聚类单元，效果如图2所示。
[0076]
步骤3：提取每一个聚类单元数据处理单元的特征，具体包括高程相关特征、投影面积、投影跨度、投影面积比及其它特征等。
[0077]
对每一个聚类单元分别进行如下特征提取：
[0078]
(1)聚类单元数据处理单元高程相关特征包括格网化点云的归一化高程、高程最大值、高程最小值和高程差；
[0079]
(2)投影面积sj为单元中的点云投影到二维xoy平面上的面积；
[0080]
(3)投影跨度k
x
和ky分别为每一个聚类单元在格网x方向和y方向上投影的最大跨度；
[0081]
其中，x方向跨度为投影单元对应的x方向上格网编号的最大值与最小值的差值；y方向跨度为投影单元对应的y方向上单元格网编号的最大值和最小值的差值；
[0082]
(4)投影面积比sr为每一个聚类单元中整体点云在xoy平面上的投影面积s1和部分点云在xoy平面投影面积s2的比值；
[0083]
(5)每一个聚类单元的面积sm以及包含的点云个数，其中面积sm为每个聚类单元所占平面格网的数目。
[0084]
步骤4：基于对各类地物形态特征分析的基础上，依据面积、高差、投影、形态四个方面的形态特征渐进提取初始行道树点云。
[0085]
步骤4.1：对各类地物的形态特征进行分析，常见各类地物投影形态如图3所示；
[0086]
步骤4.1.1：地面信息存在于整个环境中，局部范围内表现为平面，在整个环境中则是有起伏的，且地面点的高程是最小的；
[0087]
步骤4.1.2：建筑物分布于道路的两侧，垂直于地面，建筑物的高程在场景信息中最高，其在投影平面的形状为长条形状并且其截面的投影形状与整个建筑物的投影形状相差不大，投影面积大于其它地物；
[0088]
步骤4.1.3：行道树分布在道路与建筑物之间，从高程来看的话，其高程低于建筑物、路灯和电线杆等地物，高于其它的地物，行道树的投影形状呈现为簇状，投影面积比建筑物投影面积小，大于其它地物的投影面积，且行道树的树干部分的投影面积与整体的投影面积相差较大；
[0089]
步骤4.1.4：车辆分布在路面或者道路两侧，高程低于行道树、建筑物、路灯等地物，高于低矮植被，投影面积较小；
[0090]
步骤4.1.5：路灯分布在行道树之间或者道路两侧，高程低于建筑物，高于其它地物，为5～6米，投影形状呈不规则形状，具有很好的辨识性，投影面积较小，通常不大于1～2平方米；
[0091]
步骤4.1.6：低矮植被分布在行道树的周围，高程高于地面，低于其它地物，其最大高程同树干点的最高点高程一样，投影形状为不规则形状，投影面积较小。
[0092]
步骤4.2：二次噪声筛选：统计每个聚类单元的面积sm(m＝1
…
n)，此处面积为聚类单元所占格网的数目，同时设定面积阈值sy，判断聚类单元的面积是否大于面积阈值，如果存在聚类单元面积小于面积阈值sy，还需计算聚类单元中包含的点云数目，当点数小于数目阈值tm时，则判定为二次噪声，进行滤除，本实施例中tm＝10；
[0093]
步骤4.3：低矮地物筛选：计算的最大高差值，设置高差阈值th，本实施例设置为1.55米，如果最大高差值小于高差阈值th，则将该聚类单元判定为低矮地物，并删除掉该聚类单元；如果最大高差值大于设置的高差阈值th，则保存该聚类单元；
[0094]
步骤4.4：建筑物筛选：设定x和y方向上的跨度阈值t
kx
和t
ky
，本实施例设置跨度阈值t
kx
和t
ky
分别为20米、40米，如果聚类单元在x和y方向上的跨度小于给定的跨度阈值t
kx
和t
ky
时，则将聚类单元视为杆状地物，反之则为建筑物；
[0095]
步骤4.5：杆状地物筛选：根据聚类单元的投影面积设置面积阈值s
p
，本实施例设置为1.5平方米，如果某个聚类单元的投影面积大于所设置的面积阈值s
p
，则判定该聚类单元为树木；若投影面积小于阈值s
p
，则判定为其他杆状地物。
[0096]
步骤5：采用二值化算法对初始行道树点云进行优化，得到最终的行道树结果；
[0097]
步骤5.1：选择初始行道树点云数据中1.5米以下的点云数据，对其格网化，格网大小为0.1米，记录每个格网内的点云数据；
[0098]
步骤5.2：统计并判断所有格网是否为空，若格网为空值则赋值为0；若格网内有点云数据，则给该格网赋值为1，生成二值化栅格图像；
[0099]
步骤5.3：对二值图像进行聚类，设定索引变量j并初始j＝1；
[0100]
步骤5.4：随机选二值栅格图像中任一非0格网为当前格网，赋值i给当前格网；
[0101]
步骤5.5：探测中心格网八邻域属性值为非0的格网，将非0格网的属性值更新为j，步骤5.4和该步骤中处理过的非空格网后续不再进行探测；
[0102]
步骤5.6：对于新探测到的非0格网，继续八邻域探测，重复该过程直至找不到非空
格网，至此一个聚类单元处理完毕；
[0103]
步骤5.7：令j＝j+1，重复步骤5.4到步骤5.6，直到所有格网完成八邻域探测，得到m个独立的聚类单元。
[0104]
步骤5.8：判断上述结果中每个聚类单元面积是否小于给定面积阈值t0，本发明设置t0＝0.01平方米，对应格网单元的数目为10，如果聚类单元面积小于t0，则认为其对应的点云为树干点云；如果聚类单元面积大于t0，则认为其对应的点云为非行道树点云，对其进行剔除。
[0105]
本实施例所用实验数据集为大学校园实测的河南焦作地区车载lidar点云数据，数据位于校园区域，且采集时间为盛夏季节，因此植被繁茂。实验区域中包含直角形建筑物、道路、行道树、低矮地物，空中杂物等地物，地势较平坦，包含点云数据1101958个，数据显示如图4所示。
[0106]
实验在处理器intel(r)core(tm)i5-6200u cpu@2.30ghz 2.40ghz,
×
64位,16g内存和win 10操作系统的计算机上进行，借助matlab 2020a软件平台编辑和运行实验算法，并基于标准数据，评价本发明的点云提取精度。同时为验证本发明算法的有效性，通过分析得到本发明算法的提取精度评价为：第一类误差为3.6％，第二类误差为8.7％，总误差为5.0％，kappa系数为0.8646，树木提取精度为95％；证明在本实施例中本发明的基于地物形态特征的行道树点云渐进提取方法具有较好的效果。图5为本实例的最终提取效果图。
[0107]
最后应说明的是：以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分或者全部技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明权利要求所限定的范围。