基于无人机载激光雷达点云数据的人工侧柏林单木分割研究

李远航,笪志祥,闫烨琛

(1.黄河勘测规划设计研究院有限公司,河南 郑州 450003;2.水利部黄河流域水治理与水安全重点实验室(筹),河南 郑州 450003;3.天津市水务局,天津 300074;4.天津市地质工程勘测设计院有限公司,天津 300191)

森林是陆地生态建设的重要载体,其内在的植被—水碳耦合和反馈机制,不仅通过与大气进行水分和能量交换,对气候产生重要影响,同时可以促进水源涵养和水土流失控制[1],而且在“碳中和”进程中可以通过森林碳汇等方式有效地消耗CO2[2]。当前,我国人工林经营管理已从追求短期经济效益的粗放管理转变为多种生态服务功能、经济、社会效益相结合以及可持续发展模式,因此全面加强森林资源的经营和管理,有助于提升生态服务功能和效益,抑制全球变暖进程,从而实现生态文明建设的高质量发展。森林资源调查是森林经营管理的重要组成部分和前提,其关键在于对森林结构参数的精确获取和对森林发展情况、趋势的全面考究[3-4]。然而传统的人工野外调查方式具有工作强度大,耗时长,效率低等问题,在面临较大规模的林业调查需求时,工作人员难以保证其精确性;同时传统的无人机拍摄易受天气、光照和地形影响,其拍摄的遥感影像存在着反演和解译结果精确性低等问题[5]。因此提高森林结构参数反演精度和效率,对准确获取森林资源信息和预测森林未来发展情况具有重要意义。

无人机载激光雷达(UAV-LiDAR)作为一种新式的三维点云获取设备,具有成本低、效率高以及能够快速获取高精度点云数据的优点,并且在快速精准地获取林木三维结构、提取和反演森林参数等方面有着巨大优势。基于UAV-LiDAR点云图像的单木分割是单木结构信息提取和参数估测的前提,如估测单木数量、树高、冠幅等,也为后续的森林生物量反演以及森林三维模型构建提供有力支撑[6-7]。目前,基于机载LiDAR点云数据实现单木分割主要有2种思路:其一为基于点云栅格—冠层高度模型(canopy height model,CHM)进行单木分割。于海洋等[8]评估了分水岭分割、邻域增长和特征点决策树以及基于点云的聚类分割算法等4种算法的分割精度,验证了分水岭算法分割单木精确性和适应性最好。陈日强等[9]利用标记控制分水岭分割算法检测和提取果树单木树冠信息,其检测精度和提取准确率分别达到了86.35%和95.03%。其二是利用聚类算法对归一化后的LiDAR点云之间的空间结构关系和属性信息进行聚类分析进而直接进行单木分割[10]。Vega等[11]提出了一种多尺度动态点云分割新方法,通过确定最优树顶位置点作为评价标准进行树冠聚类,正确分割率达82%。刘浩然等[6]设计改进了K-means算法和基于相对间距的阈值分割算法,并评估其应用在人工林的单木分割精度。然而受限于聚类方法所需输入参数较多且复杂,难以大规模推广,相关研究引入了深度机器学习的方法进行模型训练。如Wang等[12]利用快速回归卷积神经网络(faster R-CNN)对橡胶林点云数据进行单木分割,并取得了良好效果。陈鑫鑫[3]采用深度学习框架点云神经网络模型(point net)实现单株树冠分割。但基于深度机器学习的分割方法需要构建以大量实测数据为基础的训练模型,违背了机载激光雷达技术进行单木分割时效率高、成本低的初衷,因此不具备普适性。

通过检索文献发现多数研究的单木分割精度受CHM分辨率限制,难以区分低矮植物;同时处理高分辨率数据时会出现大量过度分割现象,对于高密度林分的树冠边界识别精度不高。目前,前人研究中针对不同空间分辨率下的单木结构检测与信息提取的定量分析研究较少。基于此,本研究采用基于标记控制分水岭分割算法的CHM模型实现单木分割的方法,探讨提取人工侧柏林单木结构信息的可能性和科学性,主要目的为:1)基于点云数据生成的CHM模型,使用标记控制分水岭分割算法实现单木结构信息的检测与提取;2)基于CHM模型对不同空间分辨率的单木结构信息进行精度评价和敏感度分析;3)对不同点云密度的单木分割结果进行影响分析,为基于CHM模型的单木结构信息提取最优空间分辨率选择提供参考。

1 研究区概况及研究方法

1.1 研究区概况

研究区位于天津市蓟州区罗庄子镇,117°21′14″E,40°6′49″N,属暖温带半湿润季风型大陆性气候,四季分明,年平均气温11.5 ℃,年平均降雨量678.6 mm,全年平均相对湿度60%,年平均风速2.2 m·s-1。常年主导风向东南风,由于燕山山脉的屏障作用,加之沟谷幽深,山高坡大,故风力较缓、较弱;无霜期195 d。区域土壤以潮土为主,主要乔木树种为侧柏(图1)。

图1 研究区概况图及采集点云数据Fig.1 Overview of the study area and collected point cloud data

1.2 数据获取

本研究采用的是经纬M300RTK无人机搭载禅思LI雷达镜头开展研究工作,主要包括研究区勘察、像控点布设及测量、航线规划及机载激光雷达系统测绘、室内数据处理分析等环节,具体设备参数如表1。通过前期收集资料及实地踏勘后,对样地中的树木数量、树高及冠幅等信息进行测定,其中样地中共有树木97株。之后进行航线规划,对测区进行了机载雷达测量,数据采集时间为2022年8月10日,天气状况良好,飞行高度80 m,激光旁向重叠率50%,作业速度8 m·s-1,平均点密度为981.31点·m-2。

表1 LiDAR参数Table 1 LiDAR parameters

1.3 数据处理

1.3.1 预处理 本次研究主要利用大疆智图(DJI Terra)软件进行航测数据处理生成数字正射影像(digital orthophoto map,DOM)及完成点云数据预处理,对机载激光雷达航摄获取的点云数据进行解算,并进行各航带数据的拼接工作,最后导出三维可视化激光点云数据[13-15];将预处理后的点云数据加载到LiDAR360 UAV(Green Valley Co.,中国)软件中进行处理,其中包括点云裁剪、重采样、点云去除噪点(图2A)、滤波(图2B)、点云归一化(图2C)及图像增强(图2D)等。为了验证不同点云密度对点云聚类分割的影响,首先对点云数据进行密度重分类,即从原始密度100%分别抽稀到50%、20%和10%。点云去噪采用基于高斯去噪的方式将点云中的偏离点进行剔除,减小数字高程模型(DEM)制作过程中因噪点产生的影响。基于简单的物理过程模拟的地面滤波方法(CSF)采用地面滤波处理进行地面点的提取,将点云中地面点与非地面点分离。点云归一化是为了消除地形起伏对DEM值的影响,从而提高后续点云分割和回归分析的准确性。

图2 经过预处理的图像Fig.2 The pictures after pre-process

1.3.2 CHM模型生成 采用的是改进的渐进加密三角网滤波算法[16](improved progressive TIN densification,IPTD)对点云中的地面点和非地面点进行自动分类,通过种子点生成一个稀疏的三角网,然后通过迭代处理逐层加密,直至将所有地面点分类完毕。由于LiDAR扫描的点云密度较大,难以对地面点进行识别,因此采用基于中位地面点和分离低点的地形分类方式进行平滑处理。数字高程模型(DEM)是通过有限的地形高程数据实现对地面地形的数字化模拟,是通过不规则三角网插值的方法进行模型生成。数字地表模型(digital surface model,DSM)反映了地表建筑物、桥梁和树木等高度的地面高程模型,是采用反距离权重插值算法生成的。冠层高度模型(CHM)是表示地表包括建筑物、树木和杂草等物体高度的高度模型,通过数字表面模型与数字高程模型的差值得到。

1.3.3 基于标记控制分水岭分割算法的单木分割 本研究主要采用基于标记控制分水岭分割算法[17](marked-controlled watered segmentation,MCWS)进行CHM单木分割,这是一种基于拓扑理论的数字形态分割算法,其基本理论是将CHM的高点处视为山峰,低点处视为山谷,如果用水填充,不同山谷的水将开始汇合。为了避免这种情况,在水汇合的地方建立屏障,即不同的树冠边界,这些屏障将决定分割的结果。为了能够精确识别低矮植物冠层特征,提高高密度下林木分割的精确度,并且为了防止过分割,采用不同空间分辨率(0.1、0.3、0.5和1.0 m)来进行林木冠层分割,从而获取单木位置、树高、冠幅直径、冠幅面积和树木边界。同时为了减少CHM图像中产生的“孔洞”或“空值”,采用反距离权重插值(inverse distance weight,IDW)对“孔洞”进行平滑化填充处理,并对CHM图像进行增强效果处理,便于分离不同单木树冠[18]。

1.4 精度分析

1.4.1 林木位置和数量分析 采用3个指标对实测林木与基于高分辨率正射影像与激光雷达点云三维形态目视解译数据之间的空间和数量关系进行衡量,主要通过以下公式开展精度分析。

(1)

(2)

(3)

式中:P代表单木分割精确率,表示正确检测单木数量占整个检测结果的比例;R代表单木分割查全率,表示正确检测单木数量占真实树顶数量的比例;F代表错分和漏分的总精度,是精确率和查全率的综合结果;TP是被正确分割的单木数量;FP是被错误分割的单木数量;FN是被遗漏的单木数量。

1.4.2 树冠提取分析 通过对比实测林木与目视解译数据之间的空间和数量关系,来评估两者空间匹配精度。按照前人研究[9]所得,一般将两者识别精度分成6类[9],分别为优识别(Pr),实测林木树冠和提取树冠的重叠面积超过两者自身面积的50%;良识别(Gr),即实测林木树冠和提取树冠的重叠面积仅超过某一者面积的50%;遗漏(Mi),即实测林木树冠和提取树冠的重叠面积低于50%;合并(Me),即实测林木树冠覆盖多个提取树冠,且至少2个提取树冠与实测树冠的重叠面积均超过提取树冠面积的50%;分解(Sp),实测林木树冠覆盖多个提取树冠,且至少2个实测树冠与提取树冠的重叠面积均超过实测树冠面积的50%;错误识别(Wr),树冠被错误提取,即FP。

基于以上定义,通过准确率(AR,%)、欠分割误差(OE,%)和过分割误差(CE,%)总体评价树冠提取精度,其计算公式如下

(4)

(5)

(6)

式中:NPm、NGm、NMi、NMe、NSp、NWr分别为优识别、良识别、遗漏、合并、分解和错误识别下树冠的数量,N为以上所有识别树冠的总和。

1.4.3 林木树高和冠幅精度分析 为了评估单木分割后各林木结构特征的准确率,采用线性回归的方式对实测林木和参考林木进行验证,选用决定系数(coefficient of determination,R2)和均方根误差(root mean square error,RMSE)分别评价实测林木和参考林木之间的相关性和误差。

2 结果与分析

2.1 数字高度模型(CHM)的生成

由图3可知,经过实地林分调查,本研究设置最低树高为1 m作为区分人工侧柏林树冠和低矮杂草、土壤的高度阈值,筛选树冠高度信息,剔除非树冠像素信息,从而生成冠层高度模型。

图3 数字地形模型(A)、数字表面模型(B)和数字高度模型(C)Fig.3 Digital Elevation Model (A),Digital Surface Model (B) and Digital Height Model (C)

2.2 不同空间分辨率单木分割及精度检验

本研究通过选取空间分辨率为0.1 m×0.1 m、0.3 m×0.3 m、0.5 m×0.5 m和1.0 m×1.0 m分析不同空间分辨率的CHM信息,分别绘制了单木分割结果(图4a、图4b、图4c、图4d)、效果(图4e、图4f、图4g、图4h)及精度分析图(图5)。P均超过了0.75,其中0.3 m×0.3 m空间分辨率下的精确率最高,其次为0.5 m×0.5 m分辨率。0.1 m×0.1 m分辨率的F仅低于0.3 m×0.3 m分辨率。1.0 m×1.0 m分辨率下的单木分割P、R、F均最低,出现了部分单木未被检测的现象,遗漏分割现象较明显。

图4 不同空间分辨率下的单木分割效果Fig.4 Single tree segmentation effects under different spatial resolutions

图5 不同空间分辨率下的单木分割精度分析Fig.5 Single tree segmentation precision analysis under different spatial resolutions

2.3 树冠提取及精度分析

通过将不同空间分辨率下的提取树冠和实际树冠进行轮廓叠加分析发现,两者在轮廓边缘上不能完全重叠(图6、表2),其树冠提取准确率随空间分辨率的增大呈现出先增高后降低的趋势,当空间分辨率为0.3 m×0.3 m时,准确率达到最高。欠分割误差随空间分辨率的增大而逐渐增高,分辨率为1.0 m×1.0 m时欠分割误差最大为23%,过分割误差则呈现相反的趋势,分辨率为0.1 m×0.1 m时过分割误差最大为24%。

表2 不同空间分辨率下提取树冠与实际树冠的精度分析Table 2 Precision analysis of extracted crown and referenced crown under different spatial resolutions

图6 不同空间分辨率下提取树冠与实际树冠的轮廓叠加Fig.6 Contours of the extracted crown and the actual crown be superposed under different spatial resolutions

为了验证基于CHM方法提取树冠和树高信息的精确性,本研究对实测数据集和提取数据集进行线性拟合分析,结果如图7所示。其拟合结果R2均>0.6,随空间分辨率增大呈现先增高后降低的趋势,表现为0.3 m×0.3 m空间分辨率下的R2最高为0.818 9和0.837 1,1.0 m×1.0 m分辨率下的R2最低为0.617 9和0.688 6。同时RMSE呈相反趋势,表现为0.3 m×0.3 m空间分辨率下的RMSE最低为3.847 5和3.691 5,1.0 m×1.0 m分辨率下的RMSE最高为4.824 6和5.157 6。

图7 提取数据集与实测数据集的线性拟合结果Fig.7 Linear fitting results of the extracted dataset and the actual dataset

2.4 点云密度对单木分割结果影响分析

为了验证点云密度对单木分割效果的影响,将点云数据从原始密度100%(图8a、图8e)分别抽稀到50%(图8b、图8f)、20%(图8c、图8g)和10%(图8d、图8h),本研究将处理目标的空间分辨率设置为0.3 m×0.3 m。由表3可以发现,随着点云密度的降低,单木数量的识别精度不断降低,精度值随之下降,其最高值为100%的点云密度,其正确识别的树为85株,F为89%,精度最低的识别结果为10%点云密度,正确识别的单木下降至62株,其F为73%。同时漏分和错分的单木的数量则呈增加趋势,当点云密度为100%时,其漏分和错分的单木数量分别为12株和8株,当点云密度为10%时,其漏分和错分的单木数量最高分别为35株和10株。

表3 不同点云密度下单木分割结果的精度分析Table 3 The precision analysis of single tree segmentation results under different point cloud density

3 结论与讨论

3.1 结论

本研究探索了使用无人机载激光雷达技术获取单木分割信息的科学性与可行性,分析了空间分辨率和点云密度对单木分割和树冠提取结果的影响并评估了精确性。

1)基于无人机载激光雷达技术生成冠层高度模型可以准确分割人工侧柏林和提取树冠信息,总体上单木分割精确率均>75%,重叠严重或矮小的林木难以被检测到;树冠轮廓提取精确率均>65%,树冠之间的重叠、交叉会增加单木分割的难度;通过对实测树冠直径、树高与提取数据的回归分析,其决定系数均>0.6,均方根误差均<5.2。

2)选择合适的空间分辨率有助于提高单木分割精度,当分辨率为0.3 m×0.3 m时单木分割和树冠提取结果均为最优。其中单木分割精确率为91%,树冠提取精确率为75%,回归分析的决定系数为0.818 9和0.837 1,均方根误差为3.847 5和3.691 5。随着空间分辨率降低,单木分割精度和树冠提取精度均出现下降,出现欠分割现象;随着空间分辨率增大,单木分割数量甚至高于真实的林木数量,导致单木和树冠出现被多次分割的过分割现象。

3)单木检测与树冠提取的结果还受到点云密度的影响。本研究发现随着点云密度的降低,单木数量的识别精度不断降低,精度值随之下降,同时漏分和错分的单木的数量则呈增加趋势。当点云密度为100%时,其正确识别的树为85株,F值为89%,当点云密度为10%时,正确识别的单木下降至62株,其F值为73%。

3.2 讨论

3.2.1 空间分辨率对单木分割结果影响分析 随着林业调查对测量精度要求的不断提高,推动了无人机激光雷达技术的更新发展。在应用于不同场景、不同尺度下的林业数据采集时,可调节的空间分辨率有助于更高的图像精度和细节识别。本研究通过分析空间分辨率0.1 m×0.1 m、0.3 m×0.3 m、0.5 m×0.5 m和1.0 m×1.0 m下的CHM信息和树冠提取信息发现,空间分辨率为0.3 m×0.3 m时,单木分割和树冠提取精度最高,该现象与陈日强等[9]、Yin等[19]提取果树单木树冠信息的结果一致。当空间分辨率过小时,单木分割精度和树冠提取精度均出现下降,原因可能是CHM处理过于平滑,图像细节难以有效识别,局部最高点模糊,导致遗漏分割的树冠信息增多,出现欠分割现象。当空间分辨率过大时,CHM信息存在过度处理的可能,导致冠层高度表面出现很多与真实冠层高度不相符的凹凸坑[20-22],单木分割数量甚至高于真实的林木数量,在图像细节被扩大的同时,图像噪声进一步被放大,导致单木和树冠出现被多次分割的过分割现象。

3.2.2 点云密度对单木分割结果影响分析 树冠边界的确定不仅受到空间分辨率的影响,还受到点云密度的影响,树冠之间的重叠、交叉使得树冠边界变得模糊,同样增加单木分割的难度[23-24]。通过对不同点云密度下的分割结果进行分析发现,单木分割精度随着点云密度的减少总体呈下降趋势,由89%下降到73%,点云密度越小,错分漏分现象越明显,于海洋等[6]通过设置不同的点云密度验证点云分割结果的差异性表明了随着点云密度的降低,样地中被正确识别的单木数量逐步降低。由于人工林分种植有较高的郁闭度和林分密度,且株行距较小,同时缺乏相应的抚育管理措施,各个林木的树冠往往重叠严重,导致难以明确地区分树冠边界。点云密度过小时,图像细节识别能力较低,难以准确区分密集林木的树冠边界,容易将重叠部分划分为同一个林木,导致错分漏分现象严重,单木分割准确率也会降低。