张翔,刘洋,玉山,苏日娜,阿茹汗
(内蒙古师范大学 地理科学学院(内蒙古自治区遥感与地理信息系统重点实验室),呼和浩特 010022)
森林是陆地生态系统的主体,作为全球气候系统的重要自然资源,对实现“绿色碳库”有着重要贡献。森林资源在陆地生态系统碳循环中具有调节气候、涵养水源、固碳释氧、养分循环和维持生物多样性的重要作用[1-2]。森林结构相关参数是描述森林环境的重要指标[3],可以为相关地类划分、林分类型划分、小班区划提供数据基础,也是实现“数字林业”和“精准林业”的必然要求[4]。在森林资源调查中, 利用遥感影像提取树顶点和监测树高,可以准确获取森林单木信息,为开展大面积森林材积量评价提供有利的数据支持。
监测森林参数和有效管理森林资源是当前林业工作者主要的工作内容。传统地面实测是森林调查工作中最常用的方法,在条件良好且地形相对平坦的森林地区,通常需要使用目测法[5]和经纬仪[6]并结合测树学相关理论对活立木单木进行树高测定。但是,在地形、坡度和林分密度等因素复杂的林区,观察点的水平距离或位置的复杂性给确定树高带来了很大的困难。然而传统地面实测难以在短时期内快速、精准且大面积地进行森林调查,这对实现森林资源的快速评估较为不利。近年来,随着遥感技术的不断发展,激光雷达为森林进行大面积估测提供了新的途径。由于无人机激光雷达冠层高度模型(CHM)数据具有三维属性特征、测量精度高、能够快速获取森林垂直结构信息的能力,被广泛应用于森林结构参数 (地上生物量、郁闭度、植被覆盖度等) 的遥感估测研究中。然而,很多研究表明,激光雷达点云数据密度低,会对相关树高[7-8]、冠幅[9]和生物量[10-11]造成低估。
机载激光雷达在地面一定高度处从空中向下扫描,当激光雷达点云数据足够密时,树冠特征明显,就能够准确提取森林单木的位置。国内外研究学者对此进行了相关研究。穆喜云等[12]通过无人机机载激光雷达数据提取内蒙古大兴安岭森林结构参数。皋厦等[13]利用无人机机载激光雷达和多光谱数据判断单木分割方法进行单木位置和冠幅提取。Tejada等[14]和Liang等[15]研究认为激光雷达数据能够快速、准确提取森林结构参数,这为森林空间分布研究提供参考。
使用激光雷达数据搜索窗口的树冠识别算法,可以有效检测树冠顶点和树冠边缘特征,这将有利于确定森林单木位置。李响等[16]、耿林等[17]采用CHM数据设置局部最大值窗口搜索树顶点,能够提高单木树顶点探测提取精度,为单木树高提取确定位置。Holmgren等[18]、Zhen等[19]使用激光雷达数据,利用窗口局部最大值算法,提取树顶点位置的树高。Wulder等[20]和Wang等[21]通过不同分辨率CHM图像进行平滑滤波处理,以消除CHM图像噪点,结合目视解译方法,在一定树冠区域内,提取CHM局部最高值点作为树顶点,提取树顶点位置的树高。目前,大多数研究使用局部最大值的方法[22]从CHM确定树顶,但是往往带来过识别。由于CHM图像表面较为粗糙,分水岭变化结果通常与实际树冠有所差异,因此,结合形态学运算的性质,针对每个分水岭水槽进行图像开运算优化,优化后的结果将最大程度上还原真实的树冠形状,优化完成后的水槽被认为是一个独立的单木。
本研究利用无人机机载激光雷达数据的CHM 图像,进行中值滤波和均值滤波处理,解释CHM不同分辨率与滤波窗口大小对森林树顶点提取的共同作用,结合形态学运算的性质,并对单木尺度树顶点提取进行精度评价,最后与实测树高进行对比分析,验证无人机激光雷达数据提取单木树高的可行性。为此,本研究选取内蒙古兴安盟阿尔山市杜拉尔林场为实验区,实现单木尺度森林资源结构参数的精准评估,为后续森林蓄积量和生物量估算研究提供数据参考和理论支撑。
实验区地理位置在内蒙古自治区兴安盟阿尔山市天池镇西北部,海拔805~1 495 m,海拔起伏变化较为明显,位于寒温带大陆性季风气候区,受山地局部气候影响,年平均气温1.48 ℃,平均年降水量437 mm[23]。主要植被类型为落叶松 (Larixgmelinii)、白桦 (Betulaplatyphylla)和山杨 (Populusdavidiana),总面积4.99万hm2,森林覆盖率41%。 在整个杜拉尔林场实验区中,东北部山势险峻,中南部沟谷宽阔,地势平坦,林分郁闭度较高。选取该实验区56 m×50 m样地正射影像图作为参考,对单木位置关系进行目视解译,将单木位置标注于无人机多光谱正射影像中。试验样地如图1所示。
图1 试验样地Fig.1 Experimental area
1.2.1 地面实测数据
2021年7月11日,实地调查时以研究区正射影像为基础,现场选取不同地形坡位的天然林和人工林典型地块,设置了 4 个采样点,样方面积为25 m×25 m,使用激光测距仪对样点周围解译到的树木进行树高测量,共计获得46株样本树高数据,并利用实时动态测量技术RTK(中海达iRTK 2)对每个采样点进行了定位。
1.2.2 无人机多光谱影像数据
2021年7月11日,无人机多光谱影像数据来源于飞马公司D2000产品,相机型号为SONY a6000,飞行高度383 m,航向重叠率80%,旁向重叠率60%。有效像素2 430万,传感器尺寸23.5 mm×15.6 mm,25 mm定焦,并配备IMU惯性导航系统。在天气晴朗无风的条件下,采用D-CAM2 000航测模块对实验区进行正射影像,影像空间分辨率为0.2 m。
1.2.3 无人机激光雷达数据
研究区无人机机载激光雷达数据获取时间于2021年7月15日,共计1个飞行航带,飞行高为300 m。该飞行航带激光雷达点云密度为平均70个/m2。由中海达长续航六旋翼无人机 long120无人机机构,该旋翼无人机搭载 AP15(X)惯性导航系统的Riegl VUX-1激光雷达传感器获取。该Riegl VUX-1传感器通过近红外激光束脉冲测距原理,可提供高达500 000脉冲/s的测量速率和330°视场角。
本研究基于无人机飞马公司D2 000,D-CAM2 000航测模块对研究区进行正射影像,通过目视解译无人机多光谱数据,从中寻找参考树顶点。利用无人机Riegl VUX-1激光雷达数据,在LiDAR360软件中进行激光点云数据去噪、归一化和滤波分类,构建不同空间分辨率冠层高度模型(CHM),对相关CHM栅格数据进行中值、均值滤波设置模板与相关权重系数进行相乘,提取相关激光雷达树高数据,并进行图像处理,站点位置利用 RTK(中海达iRTK 2)进行打点定位。整理每个地面实测树高数据,以对相关数据进行精度验证,实现激光雷达遥感数据的森林树顶点确定与树高精度验证,相关操作步骤如图2所示。
图2 技术路线图Fig.2 Technology roadmap
利用LiDAR360软件对无人机激光雷达数据进行滤波和归一化处理,将激光雷达数据进行点云自动分类,从中分离出数字表面模型(DSM)和数字高程模型(DEM)进行相减,提取每一个像元坐标所对应的不同空间分辨率CHMs数据(0.1 、0.2 、0.3 、0.4 m),根据现场实际情况,样地中最低的树高大致为2 m,将CHMs图像高度阈值设置为2 m[24],目的是将森林植被类型与其他植被类型分开,如图3(a)—图3(d)所示。
图3 CHM不同空间分辨率Fig.3 CHM different spatial resolutions
以人工目视解译得到的树顶点位置为参考,将参考树顶点的0.8 m缓冲区的位置内有且只有一个提取到的树顶点,则为提取正确。如果不存在,则为提取遗漏。若在缓冲区内存在大于等于2个参考树顶点时,最接近参考树顶点的为提取正确点,其余参考树顶点都是提取错误点。另外,提取树顶点不在参考树顶点的缓冲区内时,也认为是提取错误点。将提取正确、提取遗漏和提取错误树顶点结果进行单木尺度的相关精度验证,计算公式分别为
(1)
(2)
(3)
式中:TP为正确提取的树顶点总数;FN为遗漏的树顶点总数;FP为错误提取的树顶点总数;R为召回率,表示树顶提取正确点与目视解译中得到树顶点的比例;P为准确率,表示树顶提取正确点与所有提取树顶点的比例;F为测度,是对召回率和准确性的完整描述,F值越高,精度越高。
为验证本研究所用方法提取树高的可行性,用决定系数R2对拟合度进行评价,R2值越大表示提取树高与实测树高拟合效果越好。使用均方根误差(RMSE)对树高提取精度进行评价,RMSE值越小,表示提取树高与实测树高越接近。
结合无人机多光谱正射影像图作为参考地图进行目视解译,确定样地树木参考树顶点,共计确定参考树顶点为184个。根据表1—表4,分别将无人机激光雷达数据产生不同分辨率的CHM图像进行中值滤波和均值滤波设置,利用局部最大值算法,提取森林树顶点。从这些表可以看出,树顶点提取会受到CHM影像分辨率、窗口大小和滤波类型的共同作用。当CHM空间分辨率一致时,随着均值滤波和中值滤波窗口大小的不断增加,提取森林树顶点数量将会不断减小。当均值滤波和中值滤波窗口相同时,CHM空间分辨率越低,从中提取到的提取森林树顶点数将会越少。
表1 基于CHM 0.1 m得到的树顶点数(窗口大小)
表2 基于CHM 0.2 m得到的树顶点数(窗口大小)
表3 基于CHM 0.3 m得到的树顶点数(窗口大小)
表4 基于CHM 0.4 m得到的树顶点数(窗口大小)Tab.4 Number of tree vertices extracted based on CHM 0.4 m (window size)
根据不同空间分辨率CHM图像(图4),结合目视解译的局部最大值算法,在单木尺度中提取树顶点,并对单木尺度下树顶点提取进行精度评价。在利用不同空间分辨率CHM图像提取的树顶点与参考树顶点相同或基本接近,见表5。其中,树木遗漏点为19~29株,树木错误点为13~18株,召回率为82.84%~87.97%,准确率88.61%~91.45%,F测度为85.63%~89.68%。其中, CHM0.4的F测度最高为89.68%,树木提取正确点为139株,遗漏点为19株,错误点为13株。
从CHM单木尺度提取精度中建立CHM不同空间分辨率图像,选择46个树顶提取正确点提取树高,结合实际测量的树高对CHM进行精度验证。根据图5可知,CHM不同空间分辨率图像会对单木尺度树高提取产生一定的影响。从中发现F测度最高的CHM0.4图像提取树高与实测树高之间具有良好的线性回归关系,拟合程度较高(R2=0.95,P<0.01)。从精度验证来看,RMSE=0.91 m,说明利用无人机激光雷达数据计算得到的CHM0.4图像可以有效对单木树高进行计算,数据精度高。
表5 单木尺度树顶点提取精度评价Tab.5 Evaluation of vertex extraction accuracy of single tree scale
图5 CHM不同空间分辨率不同树高提取结果与精度验证Fig.5 Extraction results and accuracy verification of CHM with different spatial resolutions and different tree heights
影响单木尺度下树提取精度的主要因素是CHM。CHM的空间分辨率直接决定了单木尺度下树顶点提取的位置,并且较高分辨率的CHM可以直接影响单木分离,这与广大研究结果较为一致。刘江俊等[24]利用不同空间分辨率CHM提取森林树顶点,研究发现CHM空间分辨率越低,提取到的树顶点数就越少。白明雄等[25]研究认为,使用无人机高分辨率可见光影像,利用分水岭分割算法对CHM进行单木分割及树高提取是可行的,与实地测量得到的树高值进行精度验证,R2=0.893, RMSE=1.23 m, 估测精度为87.58%。在本研究中,由于研究区选取的样地影像面积较小,地势较为平坦,森林郁闭度较高,树种较为单一,这对森林树顶点和树高的提取创造出一定的客观条件。
在无人机正射影像中,正射影像数据要比卫星遥感图像几何形变小,大部分树冠和冠层间隙匹配良好,适合评估森林冠层范围。另外,在利用激光雷达数据CHM评估森林树顶点和树高时,无法避免冠层间隙中小部分无效值的CHM像元没有被激光雷达点云数据填充[26]。利用局部最大值方法提取单木冠层高分辨率CHM图像时,CHM图像包含多个树顶点像元,在一定程度上会出现树顶点提取遗漏和错误,这会对森林树顶点识别和单木尺度分割造成误差,但并不影响树顶点识别和单木树高估测的结果。
在利用无人机多光谱影像提取树顶点时,使用不同固定窗口的局部最大值搜索算法,可以有效滤去小于最小树高的树顶点。Walsworth等[27]研究发现在具有代表典型森林条件的显著噪声和阴影的合成数据中,滤波图像处理技术能够更好地识别和描绘树冠。Li等[28]将加利福尼亚州内华达山脉的针叶混交林单木分割算法应用于树木分割,根据召回率、精确度和F测度对激光雷达数据结果进行评估,结果表明,森林树顶点召回率86%,森林树顶点正确率为94%,总体F测度为90%。这与激光雷达CHM数据广大研究结果较为一致,说明使用激光雷达数据提取森林树顶点具有较高的精度。
在本研究中,基于CHM不同空间分辨率,均值和中值平滑算法均会改变栅格数据边界内的所有像元值,将相邻领域内之间的像元差异变小。不同的平滑算法,也会对树顶点位置、数量和对应位置的数值产生影响。对于CHM栅格数据,未来的研究,需要采用一种邻域平滑算法来判断当前像元与领域像元的关系,以提高单木识别精度。
利用无人机激光雷达数据进行预处理,通过单木识别算法,从不同空间分辨率的CHM影像中得到树冠顶点位置,采用均值滤波和中值滤波遥感图像处理技术,对同空间分辨率的CHM影像进行处理,采用局部最大值目视解译法,从该位置中提取森林树木顶点和单木树高,研究结果表明:
(1)利用CHM图像提取森林树顶点时,会受到CHM空间分辨率和图像处理窗口大小的综合作用。当CHM空间分辨率较高(0.1 m),滤波窗口为3×3时,树木提取点接近184,较为接近无人机激光雷达影像森林树顶点的实际情况。
(2)利用无人机激光雷达不同分辨率CHM图像数据,使用局部最大值方法可以有效在单木尺度条件下提取树顶点。随着CHM图像空间分辨率不断增加,单木尺度中提取树木点将会不断减少。
(3)从无人机激光雷达数据产生的CHM图像提取46株树木正确点与实地测量的树高结果进行精度验证。结果表明,激光雷达产生的CHM对于森林树高的提取有较高的精度(RMSE=0.91 m)和较强的拟合性(R2=0.95)。