基于单站三维激光扫描数据的林木胸径提取方法

2019-08-08柯敏，解静，刘昔

浙江林业科技 2019年2期

柯敏，解静，刘昔

（天津水运工程勘察设计院天津市水运工程测绘技术重点实验室，天津 300456）

胸径又称干径，指乔木主干离地表面1.3 m处的直径，断面畸形时，测取最大值和最小值的平均值[1]。在绿化工程投资审计过程中，林木胸径大小作为审核重点之一，其精度直接影响工程造价[2]。传统的胸径测量方法主要是用标准测量工具，如测树胸径尺、游标卡尺等，在地面向上1.3 m处测量树木的直径。该方法效率低下，数据获取周期长，时效性差，且劳动强度较大，不能满足实际调查需求[3-4]。

为快速获取林木胸径数据，近年来国内技术人员尝试利用包括三维激光扫描技术在内的许多技术进行试验研究。在技术应用方面，关晓平等[5]利用电子全站仪，应用三角高程法及截面圆法分别计算树高及胸径，成果精度可靠。曹先革等[6]通过分析影响地面三维激光扫描精度的因素及误差，提出外业时扫描仪距离被扫描物体之间的距离尽量保持在50 m以内，点云数据拼接应优先选择标靶拼接。在算法研究方面，樊丽等[7]利用二维激光扫描仪获取单木胸径处点云数据，分别采用最小二乘圆拟合算法、Nlinfit算法、Taubin算法、Pratt算法对点云数据进行拟合处理，得出基于最小方差迭代圆拟合算法结果精度最高，算法稳定。孙浩等[8]利用二维激光扫描5种不同管径的PVC管，应用几何法与拟合法，结合4种角度补偿算法计算管径并分析精度，表明拟合法在所有计算方法中精度与稳定性最高。刘伟乐等[9]利用三维激光扫描技术对少量杨属Populus植物进行多站扫描、配准后分析，提出林木胸径点云最理想的厚度为1～2 cm，采用快速凸包算法提取胸径效果最好。王祺等[10]对采用三维激光扫描技术获得的树木胸径点云数据，利用蚁群算法和B样条曲线拟合技术计算树木胸径，相对传统算法精度取得较大提高。王宁宁等[11]利用三维扫描技术获取高密度人工林中杨树的冠形特征数据，通过提取胸径处截面面积，利用圆面积公式计算胸径值并进行回归分析，结果表明测量精度可靠。李丹等[12]利用多站TSL数据通过配准、过滤、抽稀等处理，生成灰度影像，结合Hough变换和椭圆拟合方法提取单木胸径等参数，并与外业实际测量数据对比，结果一致性较高。

上述研究中基于点云数据进行林木胸径提取主要包括以下特点：（1）多为利用标靶装置，通过对多站激光扫描数据进行拼接后，提取胸径处较为完整的截面数据进行胸径提取。（2）研究林木样本数量较少，数据量较小，因而可利用相关算法快速进行最优解研究，保证林木胸径提取精度。对于城市绿化林木，主要位于道路两侧，与森林、公园等林木呈片状分布不同，主要呈条带状密集分布。因而，城市绿化林木胸径调查工作具有林木数量多、密度大、通视条件差、调查精度要求高等特点，林木胸径提取难度较大。如采用上述方法进行胸径提取，容易造成标靶数量设置过多、利用率低，从而影响调查工作进度，同时对于获取的海量数据进行提取及处理，也比较复杂及困难。

针对绿化林木的特点，通过现场试验，采用单站获取的三维激光扫描数据，利用最小二乘圆拟合算法进行胸径快速提取，并对提取精度进行分析统计，拟解决传统作业方法效率低、时效差、周期长的问题，优化多站数据拼接作业方法数据量庞大、处理难度高的问题，实现快速提取绿化林木胸径的目的，以期为城市绿化林木调查统计分析奠定基础。

1 试验方法与数据采集

1.1 采集方式选取

目前，搭载三维激光扫描设备的平台有多种，常见的包括机载测量平台、车载测量平台、地面设站测量平台及背包测量平台。根据三维激光扫描技术定位方式及扫描距离、精度等技术指标的不同，目前可应用于绿化林木测量的平台主要包括机载测量平台、背包测量平台及设站测量平台，其定位方式及技术指标见下表1。其中表1所列精度为激光扫描仪相对于某一特定距离条件下精度，对于脉冲式扫描仪一般标定距离为100 m。

表1 激光扫描测量平台参数Table1 Parameter comparison table of laser scanning measurement platform

三维激光扫描测量误差主要包括：定位误差、测距误差、测角误差、集成误差及其他误差[13]。机载测量平台及背包测量平台属于移动测量平台，其位置及姿态处于不断运动变化过程中，平台内部集成了惯导系统，受惯导精度、定位精度及系统集成精度等因素综合影响，测量精度均低于设站测量平台。为满足绿化林木胸径测量精度要求，本试验选用点位精度更高的地面设站式三维激光测量平台[13-15]。

1.2 现场数据采集

试验区为天津市市政绿化区域，位于滨海新区海滨大道临港工业区段西侧，整个区域呈带状分布，南北长约1 km，东西宽约20 m。区域内林木众多、分布密集，其中西侧为乔木，东侧为灌木，本次仅对乔木进行数据采集。乔木间距约为2.5 m×2.5 m，树种为‘107杨’Populus×euramericana‘74/76’及刺槐Robinia pseudoacacia，种植时间均为2013年4月，所用苗木树龄约为3 a。乔木自西向东呈排状分布，其中西侧3排为‘107杨’，其余5排为刺槐[16]。

2018年4月，用澳大利亚Maptek公司的Maptek I-Site 8820型三维激光扫描仪、美国Trimble公司的Trimble R8 GPS双频接收机采集数据。根据测区林木密集程度，从测区北侧开始，设置每测站间距为50 m，通过逐次设站扫描获取数据，外业三维激光扫描作业步骤如下：

（1）将三维激光扫描仪架设在设计站位，对仪器进行整平，将GPS RTK天线安装在扫描仪指定位置。

（2）利用BHCORS提供的网络RTK服务，分别获取扫描设站点及后视点的三维坐标信息。

（3）将三维激光扫描仪照准后视点，并将扫描设站点及后视点坐标信息输入三维激光扫描仪控制手簿中，开始扫描测量。扫描区域角度设置为360°，扫描质量设置为中档。

（4）本站扫描完成后，检查数据完整性，如有缺失则补充扫描，如没有则搬站至下一站位扫描，直至完成整个测区扫描作业。

1.3 数据预处理

激光点云数据预处理包括两部分，一是对各测站扫描数据进行裁切、拼接，二是对点云数据进行过滤处理。

由于采用单站扫描点云数据进行胸径测量，为保证扫描点云密度及精度，以外业扫描相邻两站点连线的中垂线为界线，与绿化林木东西两侧边界构成的封闭区域为单站点云数据裁切范围，利用处理软件剔除范围以外点云数据，得到尺寸为50 m×20 m的矩形单站点云数据。然后将各测站经裁切后的点云数据，利用设站点坐标进行拼接，得到整个实验区域点云数据。

对试验区域点云数据进行地面点分类，剔除地表以上点云数据，获得地面点云数据。以地面点云数据构建数字地面模型（DTM，Digital Terrain Model），以DTM依据，截取其垂直方向上胸径处指定厚度的点云数据，过滤其他点云数据。胸径处点云数据的截取厚度分为0 cm，2 cm，4 cm及6 cm共4个等级，提取数量不同的点云数据，分析胸径测量的精度，结果见图1。由图1可知，截取厚度为0 cm且距离超过30 m时，因点云数量不够无法计算胸径值，其余截取厚度情况下均可计算胸径值。根据扫描距离、点云数量及胸径计算精度，胸径点云最适宜的厚度为4 cm。

图1 胸径点云数据不同截取厚度Figure1 Different thickness values for intercepted DBH point cloud data of DBH point cloud data

2 胸径计算与精度统计与分析

2.1 最小二乘圆拟合算法

假设（x0，y0）为圆心坐标，R为圆半径，则圆方程展开后可表示为：

令a=-2x0，b=-2y0，c=x02+y02-R2，则式（1）可以变换为：

求解出参数a，b，c，便可得到圆心坐标及半径值。假设圆弧上有限的采集点坐标为（xi，yi），采集点至圆心的距离为di，差值为ΔI=di2-R2（i=1，…，n）则有：

令F（a，b，c）为△i2之和，则有：

最小二乘圆拟合可表示为求F（a，b，c）最小。对a，b，c分别求偏导，令偏导等于0，分析比较极值点的函数值可求得最小值，公式如下：

求解式（5）可计算出参数a，b，c值，进而得到圆心坐标（x0，y0）及圆半径。

2.2 林木胸径计算

将提取的原胸径点云数据导入Auto CAD软件中，剔除其他干扰点后导出保存为DWG格式文件。利用LISP语言编制基于最小二乘圆拟合算法的程序，对各株林木胸径点云数据进行计算、统计后，得到各株林木胸径处圆心坐标及胸径值。胸径提取见图2。

图2 胸径提取示意图Figure2 The schematic diagram of DBH calculation

2.3 精度统计与分析

为验证本方法计算的林木胸径的准确性，对实验区域林木胸径现场抽样进行人工测量，将三维激光扫描点云提取的与现场测量的林木胸径值进行比对，并记录相关检测结果。依据规范要求[17]，共随机抽测104株，其中‘107杨’39株，刺槐65株，基本覆盖整个实验区域，胸径计算值与现场量测值差值统计见表2。

在胸径计算值与量测值对比的104组数据中，差值大于2 cm的有5株，＞1～2 cm的17株，≤1 cm的82株，粗差率不足5%，误差分布基本合理。现场对林木进行查看，发现该5株粗差林木均在西侧树皮严重缺失，因点云数据不完整造成差值较大。其余林木差值均小于2 cm，主要由仪器测量误差及胸径算法精度造成。不考虑上述5株林木因树皮缺失造成数据差值较大，对剩余99株林木胸径计算值按树种分别进行线性回归分析，结果见图3。