基于UAV高分影像的陕北人工林结构参数提取

高 飞， 史海静， 税军峰， 张 艳， 郭明航， 温仲明

(1.中国科学院 水利部 水土保持研究所，712100，陕西杨凌；2.西北农林科技大学水土保持研究所，黄土高原土壤侵蚀与旱地农业国家重点实验室，712100，陕西杨凌；3.西北农林科技大学草业与草原学院，712100，陕西杨凌)

以乔木为主的人工林生态系统是我国陕北黄土高原地区防风固沙、气候调节、水源涵养和生态修复的主体，也是该区经济收入稳步发展的自然屏障[1]，黄土高原由于其特殊的地形、气候和土壤条件，人工林生态系统十分脆弱；因此，快速准确地获取该区人工林的结构信息，对于生物量的测定、碳储量的估测以及人工林生态系统健康状况的评估等具有重要意义。而传统的每木检尺的森林结构调查方法，由于其劳动强度大、效率低、精度不高[2]，不能满足大范围、持续性调查的实践需求。

近年来，飞速发展的遥感技术为森林结构参数的获取提供了技术支持。目前，在森林结构参数提取中应用最广的是LiDAR数据。马振宇等[3]利用地基激光雷达实现了天然林区近地面点云数据的精细分类，点云分类精度>89%。Yin等[4]基于无人机的激光雷达数据对红树林进行了树高 (tree height，TH) 和树冠直径 (crown diameter，CD)的测量，其预测精度分别达到34.6%和46.0%。Almeida等[5]利用无人机-激光雷达系统监测了人工林的森林结构，分析了冠层高度、间隙率和叶面积指数之间的相关性。但机载LiDAR 数据获取成本较高，限制了其推广与应用[6]。随着无人机平台的发展、低空空域的开放，价格低廉的无人机平台为林业测绘工作者开辟了一条方便快捷的新途径[7]。Ayana等[8]和Sothe等[9]基于UAV-SfM(structure from motion)技术进行了单株检测、物种分类和碳动力学模拟等森林监测工作。Liang等[10]也使用无人飞行器作为地面测量的替代方法进行森林实地观测。同时，基于计算机视觉中“运动恢复结构”(SfM)算法的无人机影像处理技术不断成熟，无人机遥感测绘弥补了传统卫星遥感影像分辨率低、影像数据质量容易受到云层影响等不足[11]。但限于地区和林分特征差异，仍需对现有方法进行改进或者研发具有特定针对性提取方法，以满足不同地区林分监测需求。

陕北黄土高原地区人工林分稀疏的特点，降低了冠层重叠与遮挡对测量工作造成的困扰[12]，为小型无人机野外作业提供了条件。为此，笔者尝试以搭载RGB高分辨率相机的小型无人机为基础，通过SfM算法获取研究区正射影像 (digital orthophoto model，DOM)及点云数据，探索陕北人工林结构参数(冠幅、树高)、CHM (canopy height model)影像提取的新方法，以期快捷高效的获取陕北人工林的结构参数。

1 试验区概况

试验区位于陕北榆林市靖边县城东南25 km处天赐湾镇内的乔沟湾乡，206省道东侧的人工防护林片区。试验区地理位置E 108°54′54.50″～108°55′4.70″、N 37°24′35.86″～37°24′44.57″，海拔1 587 m。该区域属半干旱内陆性季风气候，年平均气温7.8 ℃，年降水量395.4 mm，主要集中在7、8和9月。光照充足，温差大，气候干燥，区域内主要是人工栽植的杨树(Populussimonii)纯种林以及柠条 (Caraganakorshinskii)灌木丛。实地测量前，在林地样方设置4个地面标识点M1、M2、M3和M4，并按照实测顺序对38棵树木进行了编号为t1～t38，如图1所示。

t1～t38 为38棵研究树木的编号，下同。t1-t38 are the numbering of 38 studied trees，the same below.图1 乔沟湾人工杨树林研究区Fig.1 Study area of artificial poplar forest in Qiaogouwan

2 数据与方法

2.1 地面实测数据获取

通过传统的卷尺测量树冠垂直投影边界的方式观测树冠冠幅，东西向记作东西向冠幅，南北向记作南北向冠幅。采用SNDWAY SW-600A激光测距仪测量试验区的单木树高数据。

2.2 无人机高分影像数据获取

数据获取时间为2019年9月21日，天气晴好、风力较小的时段。数据获取平台为大疆精灵4RTK专业级四旋翼无人机(图2)，主要面向低空摄影测量应用，具备cm级导航定位系统和高性能成像系统。该无人机质量为1.39 kg，最大水平飞行速度为50 km/h，定位模式下最大可倾斜角为25°。无人机内置有效像素2 000万(总像素2 048万)的云台相机，照片最大分辨率为5 472×3 648(3∶2)，具有前后水平60°、垂直±27°，下视前后70°、左右50°的视场角。

图2 大疆精灵4RTK四旋翼无人机Fig.2 DJI Phantom 4 RTK UAV

此次实验的飞行高度为80 m，飞行路线的航向重叠度和旁向重叠度分别为80%和70%，设置的重叠度满足航空摄影测量的要求。本次飞行共拍摄454张航空影像，覆盖面积为0.177 km2。无人机航线和相机拍摄点位置如图3所示。

图3 无人机航线和相机拍摄点位置Fig.3 UVA route and camera shooting position

2.3 研究方法

2.3.1 总体技术流程 基于影像拼接及三维建模软件Agisoft Photoscan Pro、点云数据处理软件Leica Cyclone 9、以及图像处理软件eCognition Developer 6.4和ArcGIS10.3，获取本次试验所需要的三维点云与DOM、数字表面模型 (digital surface model，DSM)和数字高程模型 (digital elevation model，DEM)等数据。整个技术过程(图4)采取全自动处理结合目视解译的方式进行，在快速生成研究区三维点云和正射影像的同时，也准确获取人工林单木冠幅及树高等森林结构信息。

图4 基于UAV高分影像的陕北人工林结构参数提取技术流程Fig.4 Technical flow of structural parameters extraction of plantation in northern Shanxi based on UAV high-resolution image

2.3.2 单木冠幅的提取 结合试验区面积较小，无人机RGB高分影像分辨率高，只包含红、绿、蓝3色可见波段的特点，首先利用eCogintion多尺度分割功能对图像进行分割，分割尺度参数设置为48，形状指数为0.1，紧密度为0.5。在不同尺度上对影像进行分割，可避免采用同一尺度分割影像容易造成的“分割不足”或“分割过度”现象[13]。由于研究区为人工杨树造林区，没有人工建筑的干扰，因此将分类系统设为树冠和地面2类，然后以可见光植被指数VI′ (vegetation index)、亮度 (brightness)、形状指数 (shape index)、长宽比 (length/width)为分类依据进行分类，接着手动创建训练样本用灰色代表土壤、绿色代表冠幅、黑色代表沟壑和树木阴影，执行分类后，自动分类结果如图5a和b所示。最后将绿色冠幅层的矢量文件导入ArcGIS中，结合试验区DOM目视修改矢量冠幅细节，使分类的矢量冠幅最大程度的契合DOM图像，计算并提取冠幅径长和冠幅面积数值(以17号树t17为例，如图5c和d所示)。由于NDVI(normalized difference vegetation index)的计算需要近红外波段的参与，对于只有红、绿和蓝3个波段的RGB影像，本研究通过计算得到试验区可见光植被指数VI′[14]，作为分类参照依据。

第一个项目是外科楼自然冷源利用，即在过渡季和冬季，利用自然冷源（冷却塔将冷却水降温）给空调冷冻水降温，供应ICU病房和手术室各房间的制冷需求。在这个项目上，医院也是首次采用合同能源管理的方式进行，节能效果十分显著，年节电49.5%，年节天然气2.1万立方米。该项目入选住房和城乡建设部优秀案例，纳入《建筑合同能源管理项目案例精选》。

图5 单木冠幅的提取流程Fig.5 Extraction process of crown width of an individual tree

VI′=(2g-r-b)-(1.4r-b)。

(1)

(2)

(3)

(4)

式中：VI′为可见光植被指数；r、g、b分别为红、绿、蓝3色波段占3色波段中的比例；R、G和B分别为每个RGB通道的图像的实际像素值，pixels/inch。

2.3.3 冠层高度模型的提取 CHM体现的是森林的冠层高度分布信息，是在点云数据中提取树高的基础影像数据，本研究通过2种方式获取了试验区的CHM。一种是在Photoscan中通过点云分类直接分离地面点和植被点，从而生成试验区的CHM。另外一种是通过照片拼接得出试验区的DOM，点云分类后导出生成试验区的数字表面模型和数字高程模型(图6)，在ArcGIS中通过DSM和DEM最邻近内插法重采样后相减计算间接获得的CHM。

图6 间接获取冠层高度模型所需图像Fig.6 Indirectly obtain the image required for canopy height model

2.3.4 单木树高的提取 传统的利用局部最大值算法结合基于无人机影像生成的CHM能够较有效地提取树木顶点[15-17]，这种方法提取树高是可行的，但树高的提取精度受CHM的精度和分辨率影响较大。因此，研究利用Leica Cyclone三维点云数据图模块对样地内的树木高度进行测量。新建点云数据库后，在点云数据库中生成研究区相应的三维立体模型，截取样方区域(图7)；利用Reference plane功能对样方内单木进行切割处理，切换视角，选择地面特征点作为水平基准面建立空三坐标系(图8)；调整侧视角度，切准树木顶点，树木顶点与每株树木周围地面点的相对高差即为提取树高。

图7 点云三维立体模型样方区域Fig.7 Quadrat region by the point cloud 3D model

图8 点云模型切割后单木树高提取Fig.8 Point cloud model to extract tree height after cutting

3 结果与分析

3.1 单木冠幅的提取结果

由于实地无法直接获取单木冠幅的面积数值，为便于提取精度评价，试验获取单木冠幅的南北径长及东西径长数值。表1为冠幅东西、南北径长的具体提取值和实测值以及提取面积值。

表1 冠幅提取值与实测值

在图像多尺度分割后的基础上，将遥感混淆矩阵分类精度评价方法[18]用于单木冠幅信息提取精度评价中，主要评价指标为总体精度(overall accuracy，OA)和Kappa系数。通过对遥感影像分类的精度(表2)计算得出影像总体分类精度为96%，Kappa系数为0.74，分类质量和冠幅提取均属于很好的水平。

表2 eCognition分类精度评价表

总体精度和Kappa系数的具体计算公式：

(5)

(6)

(7)

通过对实测及提取的冠幅东西径长、南北径长以及东西和南北方向的均值之间作简单的线性回归分析发现，实测值与提取值之间均呈较强的线性关系。如图9a、b和c所示，东西向冠幅径实测值与提取值之间的R2达到0.90；南北向冠幅径实测值与提取值之间的R2达到0.91；东西、南北两个方向冠幅径长的均值拟合后的R2也达到0.95，这表明，基于面向对象多尺度分割后自定义分类指标的无人机高分影像树冠信息提取方法可行，且精度满足陕北人工林生长状况快速评价与动态遥感监测的需求。

图9 冠幅及其均值、树高散点图Fig.9 Scatter plots of crown width, mean crown width and tree height

3.2 CHM和单木树高的提取结果

在Photoscan中通过点云分类，分离地面点和地上植被点，直接提取生成的CHM(图10a)，其最高高度为1 534.76 m，最低高度为1 451.71 m，最大高度差为83.05 m。通过点云分类生成DSM和DEM后，通过ArcGIS间接处理后得到的CHM(图10b)，其最高高度为15.83 m，最低高度为-8.33 m，正向最大高度差为15.83 m。由于点云分类不够准确，存在着错分类、漏分类等情况，另外原始数据中仍存在着大量的噪声点，因此导致Photoscan导出的CHM的高度达到82.64 m，远远超过正常的植被层高度。此外，虽然间接计算处理后得到的CHM在精度上有所提高，但其精度依旧无法用于进一步获取试验区准确的树高数据。这是因为在DSM格网尺寸重采样到DEM栅格尺寸的过程中，通过重采样DEM数据中误差公式RMSE计算，始终存在着约1/7的误差值，这会导致计算生成CHM存在一定程度上的误差[19]。实验中5 cm的DEM影像网尺寸缩小至0.5 cm，无论采用双线性内插法、最近邻内插法和三次立方卷积内插法进行重采样，通过重采样DEM数据中误差公式RMSE，式(8)计算，其生成的DEM的误差分别约在0.61 cm、1.11 cm和0.56 cm左右[20]。另外在影像中的山谷阴影处会出现高度被低估情况，仍需要结合高精度的DEM数据才能准确反映出冠层高度变化，但通过RGB影像去获取纯净的DEM就目前来说是不现实的。

图10 2种方式获取的研究区CHMFig.10 Study area CHM obtained in two ways

(8)

式中：R为均方根误差RMSE；Zk为重采样DEM数据的格网高程，m；ak为对应的基准DEM数据高程，m(k=1，2，…，n)。

因此，基于无人机影像产生的点云数据在Leica Cyclone中，将样地中目视识别的单木根据每个单木树冠中心最高点的位置，提取相应的单木树高，计算了实测树高与点云提取树高的误差与相对误差，其中最高树高为10.52 m，最低树高3.78 m，符合间接计算获得CHM的高差范围，结果如表3所示。图9d为样地提取的树高与实测树高的散点图和对应线性拟合结果，从图中我们可以看到提取的单木树高与实测值之间具有明显的线性关系，决定系数R2约为0.8。

表3 树高提取值与实测值

4 讨论

传统的基于像元提取冠幅的方法，充分利用影像中的光谱信息，但仍存在分类精度低，结果错误多，空间数据大量冗余的缺点[21]。Huang等[22]在利用标记控制分水岭分割算法对桂花 (Osmanthusfragrans)和罗汉松(Podocarpusmacrophyllus)无人机影像的冠幅提取研究中发现，分水岭算法对噪声和图像中较大的光谱变化较为敏感，会导致影像过分分割的现象产生。而基于面向对象提取冠幅的方法利用了影像的光谱和空间信息，通过图像分割得到若干同质像元组成的对象，可以消除大量噪声得到精度更高的信息量，能很好契合高分影像冠幅的提取处理。本研究通过试验区DOM影像，并利用eCogintion多尺度分割的方法，引入可见光植被指数VI′作为分类依据，提取出冠幅的径长和面积数据，总体精度达到96%，冠幅分类的Kappa系数为0.74，分类结果较好。

对于单木树高参数的提取，目前大多基于LiDAR点云数据在计算机中构建出森林CHM后，通过点云滤波找出每株单木的最高点进行单木识别，实现单木树高的自动提取。该方法获取树高是可行的，但其容易受到CHM分辨率、点云滤波精度、单木识别精度的限制[23]。Wallace等[24]研究发现CHM点密度的增加(从5个点/m2增加到50个点/m2)会导致算法的遗漏率显著提高(高达8%)，这说明单木识别算法、点密度对CHM提取单木树高的准确性有着极为重要的影响。因此，笔者以Cyclone中点云数据图模块对研究区单木树高进行直接提取，无需在CHM中进行单木识别操作，避免CHM精度及识别精度的限制，并且在与实测树高数据进行比对后，提取的单木树高与实测值之间具有明显的线性关系，决定系数R2为0.8。由于本研究区域较小，对于大范围的人工造林地，还需进一步的点云数据处理算法的改进，以实现自动快速的获取大范围区域内的单木树高提取。

对于不同林地情况，人工林与天然林中林木的生长形态与冠层结构存在着明显的差异，无人机森林结构参数提取方法适应性也不尽相同。Brieger等[25]基于高分辨率无人机的摄影测量点云对西伯利亚天然林的树冠直径和树高进行了提取，计算得出实测与提取树高与冠幅的简单线性回归曲线平均R2分别为0.77和0.46，并阐述了该方法的有效性取决于林分结构的复杂性。Karpina等[26]基于无人机对人工树木的生长进行了自动测量以及生物量估算，其中树高测量的精度约为5 cm。本文利用无人机对陕北黄土高原地区人工林结构参数进行了提取，计算得出实测与提取树高与冠幅的简单线性回归曲线平均R2分别为0.80和0.95，所提出的森林结构参数获取方法在陕北地区人工林具有很好的适用性。刘鲁霞等[27]以云南省普洱市天然林与杉木人工林为研究对象，利用地基激光雷达(terrestrial laser scanning，TLS)扫描的点云数据提取了样地的单木胸径与树高，发现天然林单木树高估测结果为R2=0.77、RMSE=1.46 m；人工林单木树高估测结果为R2=0.94、RMSE=0.96 m，人工林在精度方面明显要优于天然林。尽管Tian等[28]将TLS和基于无人机图像的点云集成技术相结合，提出了一种新的树高提取方法，该方法明显提高了CHM的精度，但是Wallace等[29]采用机载激光扫描(airborne laser scanning，ALS)和SfM两种遥感技术，分别从小型多旋翼无人机平台上获取了林地三维结构信息的研究发现，SfM摄影测量技术在获取越来越密集的冠层覆盖下的地形表面方面表现的确不如ALS。这2种方式都只能描述冠层结构简单、覆盖相对较低的区域的地形表面和冠层特性，无法在结构复杂的天然林中获取森林结构参数，这也有待更多学者进一步的研究。

5 结论

1)研究通过面向对象多尺度分割结合可见光植被参数VI′的方法，分割尺度参数设置为48，形状指数为0.1，紧密度为0.5进行冠幅提取，然后结合Cyclone中点云数据图模块，在树底部与地面接触处为原点建立与树干垂直的空三坐标系，对样地内单木树高进行直接提取，提取了陕北地区人工杨树林的单木冠幅与树高，与实测数据线性拟合后整体精度都较高。

2)研究发现无人机影像直接产生CHM分辨率较低，Cyclone点云数据图模块树高提取一定程度上规避了这一缺点，可实现快速准确单木信息获取的目的。但还需对整个处理算法展开进一步的改进，以便实现自动快速的获取大范围区域内的林分结构参数提取。

3)本研究提出的单木冠幅、树高以及CHM获取方法十分契合陕北地区人工林林分多为纯林，冠层交叉遮挡较少，且株行距较大等的特点，该方法可为当下陕北地区日趋自动化的林业调查工作提供一定的参考。