林业科学  2019, Vol. 55 Issue (4): 108-121   PDF    
DOI: 10.11707/j.1001-7488.20190411
0

文章信息

夏永杰, 庞勇, 刘鲁霞, 陈博伟, 董斌, 黄庆丰.
Xia Yongjie, Pang Yong, Liu Luxia, Chen Bowei, Dong Bin, Huang Qingfeng.
高精度DEM支持下的多时期航片杉木人工林树高生长监测
Forest Height Growth Monitoring of Cunninghamia lanceolata Plantation Using Multi-Temporal Aerial Photography with the Support of High Accuracy DEM
林业科学, 2019, 55(4): 108-121.
Scientia Silvae Sinicae, 2019, 55(4): 108-121.
DOI: 10.11707/j.1001-7488.20190411

文章历史

收稿日期:2017-03-28
修回日期:2017-11-28

作者相关文章

夏永杰
庞勇
刘鲁霞
陈博伟
董斌
黄庆丰

高精度DEM支持下的多时期航片杉木人工林树高生长监测
夏永杰1, 庞勇1, 刘鲁霞1,2, 陈博伟1, 董斌3, 黄庆丰2     
1. 中国林业科学研究院资源信息研究所 北京 100091;
2. 安徽农业大学林学与园林学院 合肥 230036;
3. 安徽农业大学理学院 合肥 230036
摘要:【目的】集成多时期航片数据和由机载激光雷达数据获取的密集林区数字高程模型,估测多时期杉木人工林冠层高度,并对其生长情况进行定量监测,为多时期航片监测森林生长趋势和评价林地生产力提供可能。【方法】首先基于分类后的激光雷达点云数据获得林下高精度数字高程模型和森林数字表面模型,利用航片数据构建立体像对,通过自动立体匹配算法生成森林冠层的摄影测量数字表面模型,然后借助数字高程模型将2种数字表面模型进行高度归一化,提取研究区多时期森林冠层高度。利用1996、2004年历史航片和2014年数字航片以及激光雷达数据,构建18年内皖南杉木人工林3期森林冠层高度,并对其精度进行分析。【结果】1)由2014年数字航片和激光雷达数据获取的森林冠层高度的R2为0.52,RMSE为1.79 m;2)由2014年数字航片处理得到的森林冠层高度与对应样地实测上层木的平均高验证精度较高,平均绝对误差1.59 m,平均相对误差15%,最大绝对误差3.45 m,最大相对误差30.80%,测量精度85.00%;3)由1996、2004、2014年航片得到3期杉木人工林冠层高度,其增长趋势与树高生长曲线预测趋势一致。【结论】在多山复杂地形条件下,利用航片可准确定量反映山脊向阳面的森林冠层高度变化,但对于山谷阴影处,则会出现冠层高度被低估情况,利用多期航片结合高精度DEM数据可定量反映上层木的冠层高度变化。
关键词:历史航片    密集匹配    激光雷达    数字高程模型    数字表面模型    冠层高度模型    生长监测    
Forest Height Growth Monitoring of Cunninghamia lanceolata Plantation Using Multi-Temporal Aerial Photography with the Support of High Accuracy DEM
Xia Yongjie1, Pang Yong1, Liu Luxia1,2, Chen Bowei1, Dong Bin3, Huang Qingfeng2     
1. Research Institute of Forestry Resource Information Technigues, CAF Beijing 100091;
2. School of Forestry & Landscape Architecture, Anhui Agricultural University Hefei 230036;
3. School of Science, Anhui Agricultural University Hefei 230036
Abstract: 【Objective】This study integrated multi-temporal aerial photographs and DEM derived from airborne LiDAR data to calculate the forest canopy height of Cunninghamia lanceolata and monitor the variation of growth quantitatively.【Method】First of all, high accuracy digital elevation model beneath canopy and forestry digital surface model were constructed based on classified LiDAR point cloud data. Digital surface models were then created by applying an automated stereo-matching algorithm to the scanning copy of aerial photographs. These multi-temporal canopy heights were obtained by subtracting the LiDAR ground elevations from the two kinds of DSM. Using historical aerial photographs of 1996, 2004 and digital aerial photographs, LiDAR data of 2014, multi-temporal CHMs were reconstructed within a period of 18 years, and the accuracy was evaluated and analyzed.【Result】1) The R2 between the canopy height models acquired by LiDAR data and corresponding digital aerial photographs in 2014 is 0.52, and the root mean square error is 1.79 m. 2) Compared with the measurements from field plots, our data showed an accuracy of 85.00% with mean absolute errorand mean relative error of 1.59 m and 15.00%, and the maximum absolute error and maximum relative error of 3.45 m and 30.80% respectively. 3) Combined with the aerial photos of year 1996, 2004 and 2014, these multi-temporal canopy height models of Cunninghamia lanceolata plantation have a similar growth trend to the predicted growth curve.【Conclusion】Based on the results, utilizing aerial photographs can characterize the variation of canopy height in the sunny slope of mountainous terrain. However, for forests located in the valley bottom, the canopy height would be under estimated with aerial photographs. Multi-temporal aerial photographs combining with the high accuracy DEM can reflect the variation of overstory's height, which provide the possibility for monitoring the forest's growth trend and access the forest's productivity.
Key words: historical aerial photographs    dense matching    LiDAR    DEM    DSM    CHM    growth monitoring    

激光雷达(light detection and ranging, LiDAR)是近年来国际上发展十分迅速的主动遥感技术,可快速获取真实地表的高精度三维信息,精准估算森林树高、冠层郁闭度和地上生物量等,在森林参数的定量测量和反演上取得了成功应用(庞勇等,2005汪承义,2007),特别是对森林高度和垂直结构的探测能力,具有传统光学遥感数据难以比拟的优势。数字摄影测量利用空三加密处理技术可以得到航空像片中的森林表层绝对高度,但在密集林区却不能获得地面信息(宫鹏等,1999),因此,将LiDAR数据获取的高精度地形与数字摄影测量技术相结合,则可以生成林区冠层高度模型(canopy height model,CHM),反映材积和林地质量。历史航片具有长时间序列、立体成像等特点,且相较同时期卫星像片具备高分辨率的独特优势,可以作为历史森林资源调查的高精度数据源,通过对多时期航片构建立体像对,进行空三加密生成多时期数字表面模型Photo_DSM(the photogrammetric digital surface model),结合LiDAR获取的高精度地形,则可以得到多时期CHM数据,对森林资源调查以及森林生长预测将产生积极作用。

由立体航空影像获得的森林冠层表面绝对高程与由LiDAR数据获得的林下地形相结合提取林分冠层高度和单木冠层高度,试验证明LiDAR数据获得的LiDAR_CHM与2种数据结合生成的Photo_CHM之间具有较高的相关性(St-Onge et al., 2011; Véga et al., 2013)。Véga等(2009)基于1946—2003年的历史航片与近期LiDAR覆盖数据,构建了1946—2003年的多时期CHM,并与森林生长方程结合求算立地指数,取得了良好效果。Bohlin等(2016)利用多时相(生长季和落叶季)航空像片构建CHM,结合航空影像光谱信息对不同季节落叶林的生长情况进行变化监测,并以此信息进行分类,分类精度得到提高。Ginzler等(2015)基于2007—2012年夏季航片,按照研究区地形差异采用不同匹配策略对瑞士全国航片数据进行立体匹配,进而更新全国森林冠层高度,同时也进行了单木和林分尺度上的落叶林、针叶林树高精度比较。Nurminen等(2015)对研究区森林树高进行估测,得出通过现代摄影测量技术处理68年间的航片,可以恢复时间序列的森林高度信息。

宫鹏等(1999)将历史航片像对生成的DSM进行编辑加工生成DEM,从已获得的DSM中减去相应的DEM得到了消除地面起伏的树冠表面三维模型,并比较了25年间橡树(Quercus palustris)草原变化, 使融合三维表面数据和多光谱数据以提高遥感信息提取的精度成为可能。王佳等(2011)利用VirtuoZo软件对数字航片构建立体模型,测量研究区50株样木高度,与全站仪实测树高进行比较,结果显示2种方式获取的树高很高相关性。李明泽等(2009)利用航片及相应的外方位元素建立立体像对模型,测量研究区32株样木高度,并与超声波测高器实地测量的样木高度进行对比,结果表明2组数据差异不显著,基于数字摄影测量的立木高度实际精度为90.92%。

目前,国内研究主要是基于平缓地形的研究区,编辑由航片生成的DSM,进而得到DEM或构建立体像对,求得单木顶部坐标和底部坐标,二者相减得到树高,这种方法虽然可以得到较高精确度的树高,但是需要获取树木底部高程,在密集林区或复杂地形的山区需要大量手工编辑工作。而由LiDAR数据获取的高精度地形结合多时期航片提取历史树高的研究工作还很少开展。本研究对复杂地形区域的多时期历史航片与最新数字航片分别匹配立体像对,生成多期DSM,与LiDAR数据相结合生成多期CHM,并估测林分尺度上的树高变化,以期为多时期航片监测森林生长趋势和评价林地生产力提供可能。

1 研究区概况与数据来源 1.1 研究区概况

研究区位于安徽省黄山市休宁县,中心地理坐标为118°2′E、29°7′N,海拔150 ~360 m。气候为亚热带季风湿润气候,四季分明,春秋短、夏冬长,热量丰富,雨水充沛,日照时数和日照百分率偏低,云雾多、湿度大,年平均气温15.5~16.4 ℃,降水量1 395~1 702 mm。地形以山地为主,中低部山地土壤主要为黄壤,土层较厚,肥力较高,有利于树木、茶叶等生长。主要森林类型为亚热带常绿落叶阔叶混交林,主要树种有杉木(Cunninghamia lanceolata)、毛竹(Phyllostachys edulis)等。

1.2 机载LiDAR数据和数字航片获取

采用中国林业科学研究院的CAF-LiCHy系统(Pang et al., 2016)获取LiDAR数据和数字航片数据,获取时间为2014年10月5日,飞行平台为塞斯纳208B。LiDAR传感器型号为LMS-Q680i,其主要技术参数如下:1)波长为1 550 nm;2)扫描速率每秒138线;3)扫描角±30°;4)发射频率400 kHz;5)点密度4~5 m-2;6)采样间隔1 ns。本研究使用的LiDAR点云数据,每个激光返回点包含三维坐标(x, y, z)、强度(intensity)、回波次数(return number)等信息。数字相机型号为Hasselblad H4D-60,航向和旁向重叠率分别为60%和30%。

1.3 历史航片获取

历史航片数据由瑞士威特公司生产的传统框幅式量测胶片相机RC-10采集,像幅为23 cm×23 cm,拍摄时间分别为1996年5月26日和2004年12月5日。原始航片为黑白胶片,影像采集清晰度高、反差适中、色调均匀。对原始胶片影像进行底片扫描(900 dpi),得到空间分辨率1 m的数字化影像, 文件大小为60~70 M,扫描影像行列大小为8 218× 8 218。东西方向飞行,影像航向重叠率为64%~70%,旁向重叠率为31%~42%,每时期2张航片覆盖研究区域。具体参数如表 1所示。

表 1 航空相机和像片参数 Tab.1 Technical parameters of the aerial camera and photographs
1.4 样地数据获取

研究区总面积为48 000 m2,1992年进行采伐迹地更新造林,统一种植杉木幼苗,至2014年林分平均年龄为22年。在研究区选择16块样地(20 m ×30 m),采用胸径尺测定1.3 m处胸径,起测径级为5 cm,记录胸径、树高、枝下高、冠幅等信息。利用差分GPS(differential GPS,DGPS)定位样地4个角点坐标,样地分布如图 1所示。

图 1 野外调查样地分布 Fig. 1 The spatial distribution of field plots
2 研究方法

基于多时期航片构建立体像对,通过密集匹配得到摄影测量点云,并进行光束法区域网平差,生成反映树冠表层绝对高度的Photo_DSM。由于激光作为机载LiDAR传感器的信息传播脉冲,能部分穿透稀疏冠层或冠层间隙到达地面返回回波信息,因此可准确获取林下地形信息。对采集的LiDAR点云数据分类滤波,内插生成高精度LiDAR_DSM(digital surface model generate from LiDAR data)和LiDAR_DEM(digital elevation model generate from LiDAR data)。以LiDAR_DEM为坐标基准,对Photo_DSM转换坐标系,对应像元配准。LiDAR_DEM表示林下地形的起伏情况,由于研究区地质条件稳定,地形随时间变化极其微小,所以本研究将其作为不变量,分别用研究区的多期DSM(Photo_DSM, LiDAR_DSM)与LiDAR_DEM做差。将多期Photo_DSM归一化到同一高程平面,从而生成时间序列的树冠高度分布影像。本研究的技术路线如图 2所示。

图 2 技术路线 Fig. 2 Technical route
2.1 航片处理 2.1.1 历史航片处理

历史航片主要采用VirtuoZo软件进行处理(Zhang et al., 1996)。由于框幅式量测相机拍摄过程中机载设备限制,无法同时采集拍摄瞬间的航片外方位元素(xyz 3个瞬时姿态角度)数据,因此不能直接通过外方位元素进行空中三角测量,需要地面控制点进行单片空间后方交会,解算出单张航片的外方位元素。本研究中,LiDAR_DSM保留了充足的地物纹理、形状和大小等信息,同时能保证精确的地面坐标信息,满足提取控制点要求,所以地面控制点从LiDAR_DSM上进行选取。控制点选取道路拐角、房屋边角等人工建筑。

选取框幅式历史航片控制点主要遵循以下原则:1)研究时段内稳定,无变化,在LiDAR_DSM上能同时获取;2)在影像中易于寻找和辨识(Micheletti et al., 2015);3)分布均匀,控制点数量满足精度要求。本研究在1996和2004年历史航片中共选取56个控制点,主要分布于村落及林间道路、田埂,在郁闭度较高的林区,由于森林遮盖无法寻找精确控制点,所以对于空中三角测量、解算航片外方位元素会存在一定误差。

运用VirtuoZo软件,首先通过控制点对历史航片进行单片空间后方交会,解算共线条件方程,求出拍摄瞬间2张航片的外方位元素;然后构建航片立体像对,重采样生成左右核线影像,后进行密集匹配。密集匹配是将核线影像进行低通滤波,增大采样间隔,得到一个像元总数逐渐变小、分辨率逐渐降低的影像金字塔,然后基于特征点匹配算法和局部多点松弛影像匹配算法(Zhang et al., 1996),由低分辨率到高分辨率对核线影像逐级匹配,同时利用高频信息进行精确相关(张剑清, 2011)。得到核线影像之间的密集匹配点后,再利用航片外方位元素进行空间前方交会,解算共线条件方程获得同名像点的三维坐标(xyz)。本研究中,由于密集森林纹理性较差,山体阴影对影像清晰度也有影响,在密集匹配时出现了“孔洞”现象,而对于此现象,软件会自行进行“孔洞”插值,因此不可避免地造成了高程精度降低。同时,由于飞行高度较高,在纹理缺乏区、地形遮蔽区通常会出现比研究区实际地形或高或低的匹配结果,因此需要人工调整等高线高度使匹配结果与视差影像显示的实际地形起伏相一致(杨战辉等, 1998)。本研究对冠层高度进行研究,需要调整等高线与冠层表面起伏一致,重采样可以得到空间分辨率为1 m的Photo_DSM。

2.1.2 数字航片处理

采用Pix4DMapper软件对2014年数字航片进行处理(St-Onge et al., 2015)。由于在拍摄数字航片的瞬间,IMU同时记录了航片的外方位元素,所以本研究结合外方位元素对数字航片直接进行密集匹配与光束法区域网平差等处理,得到研究区的Photo_DSM。鉴于后续处理过程中尺度需要统一,本研究将Photo_DSM的空间分辨率重采样到1 m。

2.2 LiDAR数据处理

采用TerraSolid软件对LiDAR点云数据进行处理。LiDAR点云数据包括地面点、植被点等,一般根据激光点之间的空间关系进行分类。分类的基本过程为:首先去除明显高于或低于地表物体的噪声点;然后提取地面点;最后分离植被点和非植被点。

处理过程主要包含以下3步:1)去除噪声点,噪声点一般为孤立点,包括一些明显高于地表物体和明显低于地表的激光点,根据绝对高程或设定阈值去除明显的异常点;2)点云分类,将点云分为地面点和非地面点,本研究区中主要地物为森林,故提取非地面点为森林反射的激光雷达脉冲点;3) LiDAR_DEM和LiDAR_DSM生成,将分类后的地面点和非地面点分别采用TIN内插算法和最大高程插值法生成空间分辨率为1 m的LiDAR_DEM和LiDAR_DSM。

2.3 影像配准 2.3.1 坐标系统一

影像间精确配准的关键是历史航片的绝对定向,控制点的选取尤为重要。LiDAR因其高密度点云、高精度地理位置,可以提供纹理信息丰富、位置精确的LiDAR_DSM。控制点以LiDAR_DSM为基准影像进行选取,结合历史航片处理得到Photo_DSM,后将Photo_DSM转换到与LiDAR_DSM相同的坐标系中(WGS84 UTM 50N)。因LiDAR_DSM和LiDAR_DEM是从相同的点云数据源经分类后插值得到的,故二者在坐标位置上完全一致。通过坐标系转换,本研究将Photo_DSM和LiDAR_DEM统一到相同坐标系中。

2.3.2 对应像元配准

在进行影像对应像元配准时,如何在Photo_DSM和LiDAR_DEM上寻找到充足的控制点是关键。由于本研究区主要地物为森林,Photo_DSM是森林冠层上表面的绝对高度,表现为光滑曲面,而LiDAR_DEM描述的是地面高程,没有丰富的地表纹理信息,因此在二者之间寻找控制点几乎不可能,需要一种新途径解决上述问题。在生成Photo_DSM过程中,也生成了对应区域的正射影像,正射影像是航片经Photo_DSM逐像元微分纠正处理得到的,所以正射影像像元位置与Photo_DSM像元位置对应一致且含有丰富的光谱纹理信息(吕丽萍, 2006)。而LiDAR_DSM和LiDAR_DEM像元坐标位置一致,同时含有丰富的地面信息。所以,本研究选取LiDAR_DSM为基准影像,以正射影像为配准媒介,利用控制点对Photo_DSM进行配准,将2014年LiDAR_DSM分别与1996、2004和2014年正射影像进行配准,进而得到控制点对3期Photo_DSM进行纠正,使影像间逐像元一一对应。

2.4 CHM生成

利用LiDAR数据和航片数据估测研究区CHM,具体公式如下:

$ \text{Photo }\!\!\_\!\!\text{ CHM=Photo }\!\!\_\!\!\text{ DSM}-\text{LiDAR }\!\!\_\!\!\text{ DEM;} $ (1)
$ \text{LiDA}{{\text{R}}_{-}}\text{CHM}=\text{LiDA}{{\text{R}}_{-}}\text{DSM}-\text{LiDA}{{\text{R}}_{-}}\text{DEM}\text{。} $ (2)

式中:LiDAR_DSM、Photo_DSM分别为利用LiDAR数据、航片数据经处理后得到的森林冠层顶部绝对高度;LiDAR_DEM为由LiDAR数据经处理得到的林下地表绝对高度。

2.5 检验方法 2.5.1 绝对定向与相对配准精度分析

控制点对历史航片的绝对定向影响较大,直接决定Photo_DSM的水平和高程精度。设dx、dy、dz为控制点平差解算后在物方坐标系下xyz3个方向的误差,以上述参量对绝对定向精度进行评价。

由于本研究重点为估测树高,因此LiDAR_DEM需与Photo_DSM精确相对配准。设Error x、Error y分别为配准点在xy方向相对于最大似然位置的误差,RMS为均方根,以上述参量对相对配准精度进行评价。

2.5.2 树高精度验证

本研究缺少1996、2004年样地实测树高和对应年份的LiDAR数据,故从以下2个角度对Photo_CHM精度进行验证。

1) 对2014年Photo_CHM和2014年LiDAR_CHM进行回归分析,计算2014年Photo_CHM对2014年实测数据的误差和测量精度,根据2014年Photo_CHM的精度,间接说明采用相似密集匹配算法得到的1996、2004年Photo_CHM的精度和可靠性。

2) 结合16块样地的Photo_CHM在1996、2004和2014年的变化趋势与杉木树高生长曲线的一致性,进一步检验1996、2004年Photo_CHM的可靠性。

2.5.3 Photo_CHM与样地树高精度分析

数字航片只能获取上层木的分布情况,中低层木因被上层木遮挡无法在航片上显示,因此验证Photo_CHM与样地树高精度时需要提取上层木树高。结合正射影象,对样地中可被判读的单木树冠进行目视解译并统计数量,其总数约占立木总数1/3,同时还需要满足冠幅较大易识别的特点。结合样地实测数据设置条件如下:1)树高选取,以样地为单元对树高由大到小进行排列,选取前1/3树高;2)冠幅选取,按照东西冠幅和南北冠幅同时满足大于等于2.5 m的标准样地树高选取实测树高。

为方便进行Photo_CHM与样地树高之间的精度比较,选择绝对误差(absolute error, AE)和相对误差(relative error, RE)进行统计分析,其计算公式(樊仲谋, 2015)如下:

$ \text{AE}=X-{{X}_{0}};$ (3)
$ \text{RE}=\frac{\text{AE}}{{{X}_{0}}}。$ (4)

式中:X为立体像对提取的树高;X0为样地调查获取的树高。

3 结果与分析 3.1 多时期CHM生成及变化分析

按照上述方法对数据进行处理,得到1996、2004和2014年的正射影像和Photo_CHM影像如图 345所示。

图 3 1996年正射影像(a)和Photo_CHM影响(b) Fig. 3 The orthoimage(a) and Photo_CHM image(b) in 1996
图 4 2004年正射影像(a)和Photo_CHM影像(b) Fig. 4 The orthoimage(a) and Photo_CHM image(b) in 2004
图 5 2014年正射影像(a)和Photo_CHM影像(b) Fig. 5 The orthoimage(a) and Photo_CHM image(b) in 2014

Photo_CHM影像展示了研究区内每个像元处的树冠表面高度。空间分辨率为1 m,一个像元的面积为1 m2,通常一株成熟杉木的冠幅为9~16 m2,所以Photo_CHM可显示杉木单木的树冠高度情况。单木冠层高度可从树冠侧视图(图 6a)读取,单木树冠所占像元数量可从树冠俯视图(图 6b)读取。两图结合,可以反映单木树冠在不同像元位置的高程。

图 6 树冠侧视图(a)与树冠俯视图(b) Fig. 6 Side view of canopy(a) and aerial view of canopy(b)

图 345中绿色表示高冠层森林的分布区域,红色表示低冠层森林的分布区域,蓝色表示CHM为负值的区域。但在实际研究区中,不存在CHM为负值的情况,分析其原因主要是地形阴影、林区冠层光滑、树木密集、地形遮挡等因素造成匹配算法在密集匹配像点时出现“孔洞”现象,导致Photo_DSM被低估(St-Onge et al., 2011; Lisein et al., 2013; Harwin et al., 2012),且地面控制点误差也会使Photo_DSM低于对应区域的LiDAR_DEM(智长贵等, 2006)。从图中可看出,不同颜色区域色彩连贯均一,树冠细节处体现极少,分析其原因主要是太阳光不能穿透密集树冠,无林下地表返回光线,航空相机接收不到含有地表信息的光谱信号;同时,航空像片在原始分辨率下不能体现树冠内部细节和树冠间孔隙(White et al., 2013)。

将LiDAR_DSM与LiDAR_DEM做差,得到空间分辨率为1 m的LiDAR_CHM影像如图 7所示。因为本研究中LiDAR点云密度为4~5 m-2,所以LiDAR_CHM既包含充足的冠层表面信息,又包含丰富的冠层间隙和林下信息。

图 7 研究区阴影分布和等高线显示(a)及2014年LiDAR_CHM分布(b) Fig. 7 The shadow distribution with contour display in research area(a) and the distribution of LiDAR_CHM (b) in 2014
3.2 绝对定向精度评价

利用地面控制点对历史航片进行绝对定向后的精度显示(部分控制点和总体精度),结果如表 23所示。从表中可看出,部分控制点经平差解算后,dx、dy、dz分别存在0~3 m的误差,主要原因是在对历史航片的立体像对取控制点时,纹理性较差及分辨率较低导致立体像对间寻找控制点存在误差,但总体满足误差在3 m内。

表 2 1996年控制点精度显示 Tab.2 The accuracy display of control point in 1996
表 3 2004年控制点精度显示 Tab.3 The accuracy display of control point in 2004
3.3 相对配准精度评价

LiDAR_DEM与Photo_DSM相对配准精度如表 45所示(部分配准点)。从表中可看出,LiDAR_DEM与1996、2004年Photo_DSM相对配准总体精度在1个像素左右,基本满足要求。

表 4 相对配准精度(LiDAR_DEM与1996年Photo_DSM) Tab.4 The relative accuracy of registration (LiDAR_DEM and Photo_DSM in 1996)
表 5 相对配准精度(LiDAR_DEM与2004年Photo_DSM) Tab.5 The relative accuracy of registration (LiDAR_DEM and Photo_DSM in 2004)
3.4 Photo_CHM与LiDAR_CHM的相关性分析

本研究使用LiDAR_CHM作为验证数据对Photo_CHM进行验证。激光雷达具有主动发射脉冲的属性,可以穿透稀疏冠层到达地面,归一化后的LiDAR点云数据可以反映整个冠层的垂直分布情况;而航片是以地表物体反射的可见光作为信息媒介,框幅式相机只能反映密集森林冠层上表面的分布情况(图 8)。所以,Photo_CHM与LiDAR_CHM在稀疏冠层部分的误差较大。

图 8 林区LiDAR_DSM(a)、Photo_DSM(b)成像示意 Fig. 8 Diagrammatic drawing of LiDAR_DSM(a) and Photo_DSM(b)

为了减小Photo_CHM与LiDAR_CHM在稀疏冠层部分的误差,本研究首先以几组不同搜索窗口分别计算Photo_CHM的均值(Photo_CHMmean)、最大值(Photo_CHMmax)、LiDAR点云平均高(LiDAR_CHMmean)、75%(LiDAR_CHMP75)和99%(LiDAR_CHMP99)分位数处树高,然后采用ENVI软件对Photo_CHM提取的不同CHM与LiDAR数据提取的不同CHM进行相关性分析,结果如表 6所示。从表中可看出,3×3搜索窗口Photo_CHMmean与LiDAR_CHMP99相关性最大,相关系数为0.862。空间分辨率为3 m的LiDAR_CHMP99对于体现杉木的单木树高效果良好,符合研究区实际情况。

表 6 不同CHM提取及相关系数 Tab.6 Relative analyze of extracted different CHMs by search windows in different sizes
3.5 树高精度分析

以2014年LiDAR_CHM和样地实测数据对2014年Photo_CHM进行精度分析,并结合树高生长曲线对1996、2004年Photo_CHM进行精度验证。

3.5.1 2014年Photo_CHM与2014年LiDAR_CHM回归性分析

由上文可知,Photo_CHMmean与LiDAR_CHMP99的相关系数最大。在此基础上,对Photo_CHMmean与LiDAR_CHMP99进行回归性分析如图 9a所示,并以2 m为分类间隔,对CHM(LiDAR_CHM、Photo_CHM)按大小分类后统计相对树高分布如图 9b所示。从图中可知,Photo_CHM与LiDAR_CHM线性回归的R2=0.52,RMSE=1.79 m。

图 9 样地内LiDAR_CHM与Photo_CHM相关性分析(a)及其高度相对分布(b) Fig. 9 Relationship between LiDAR_CHM and Photo_CHM(a) and the relative distribution of height class(b)
3.5.2 2014年Photo_CHM对2014年实测数据的误差和测量精度

将Photo_CHM得到的16块样地的平均冠层高度(CHMmean)与实测树高的上层木平均高进行对比,计算误差和测量精度,结果如表 7所示。由表 7可知,Photo_CHM得到的16块样地的平均冠层高度(CHMmean)与实测上层木平均高之间的最大绝对误差为3.45 m,平均绝对误差为1.59 m,最大相对误差为30.80%,平均相对误差为15.00%,测量精度为85.00%。

表 7 树高提取值与实测数据对比 Tab.7 The aerial photography measurement of tree height comparison with the field measured height
3.5.3 基于杉木树高生长曲线的1996、2004年Photo_CHM精度验证

研究区杉木树高生长曲线(汪少俊, 2015)如下:

$ h=23.0463/[1+5.508\text{6EXP}(-0.0798X)]。$ (5)

式中:h为杉木林的算术平均高;X为杉木树龄。

将杉木生长曲线以树龄为自变量、树高为因变量进行展示。利用1996、2004和2014年Photo_CHM,计算3期Photo_CHM中每块样地的CHMmean,最终得到16块样地的Photo_CHM变化趋势线与杉木生长曲线叠加展示如图 10所示。

图 10 理论树高与样地平均树高变化趋势 Fig. 10 The change trend of the mean tree height in the field spot and theoretical tree height

图 10中蓝色曲线为杉木树高生长曲线,可见树高随着树龄增加呈二次线性增大。其他不同颜色曲线是由3期Photo_CHM得到的16组Photo_CHM变化趋势线,当树龄为22年(2014年)时,Photo_CHM与LiDAR_CHM回归性较好,且16块样地的CHMmean对样地实测树高误差较小,测量精度较高,可知2014年Photo_CHM是准确的。同时,16组Photo_CHM变化趋势线与杉木树高生长曲线变化趋势类似,这也从树龄与树高的关系角度侧面证明了1996和2004年Photo_CHM具有较高的准确度。

3.6 CHM变化

图 1112中黄色为森林CHM逐渐增加区域,其变化趋势符合杉木生长规律;红色和绿色分别为森林CHM降低和CHM增长过快区域,其变化趋势不符合样地实际情况和杉木生长规律。分析其原因主要是因为有1期或2期历史航片数据在该区域生成的CHM异常,结合正射影像和阴影变化图(图 8a),可知CHM变化异常区域主要位于研究区东部山脊处以及山谷,该区域地形起伏大,传感器视角易被遮挡,部分地面信息因此缺失;同时,起伏的地形会形成阴影,使纹理信息缺失严重,进而影响密集匹配效果。杉木为喜阳树种,喜欢潮湿和风力较弱的环境,阳光充足可加快杉木成长(袁晓红等, 2012)。山谷区域因常年阴影,导致Photo_CHM被低估;山脊处因是裸地且面积较小,所以无明显的高程变化。分别统计16块样地的CHM平均值,得到不同样地及对应CHM变化图,如图 13所示。

图 11 1996—2004年CHM变化 Fig. 11 The variation of CHM from 1996 to 2004
图 12 2004—2014年CHM变化 Fig. 12 The variation of CHM from 1996 to 2004
图 13 不同样地及对应CHM Fig. 13 Different field plots and corresponding CHMs
4 讨论

本研究集成多时期航片数据和由机载LiDAR数据获取的密集林区数字高程模型,估测多时期杉木人工林冠层高度,并对其生长情况进行定量监测。结果发现,Photo_CHM提取精度较高,历史航片可以作为一种数据源估测森林冠层高度;但是Photo_CHM对于森林冠层分布的刻画不如LiDAR_CHM精细(St-Onge et al., 2008)。尽管如此,历史航片作为LiDAR数据无法替代的历史数据源在估测森林冠层高度上仍具有较高的精度和应用潜力。

4.1 密集匹配的影响因素

本研究的关键技术环节为密集森林区域的密集匹配。在匹配过程中,由于太阳、传感器、地物之间的几何关系引起航片间的差异,导致生成的Photo_CHM与LiDAR_CHM之间差异很大。如果在一张影像中出现某地物点,但是由于阴影或遮挡等因素,在另一张影像上没有该地物点,则密集匹配会失败(St-Onge et al., 2008),因此1个像点位置(xyz)的确定至少需要2张航片中同时出现该点(St-Onge et al., 2015)。

匹配算法是密集匹配的另外一个重要影响因素,本研究对研究区冠层高度进行平滑处理,抹平了突出的单木冠顶,因此抬高了冠层间隙高度(Véga et al., 2013)。

本研究历史航片是以测绘为目的获取的,航向重叠度较低,因此造成立体像对间纹理变化较大,在以密集匹配为基础获得森林冠层高度时,对Photo_DSM影响较大。1期历史树高通过对立体像对的处理获得,但是由于没有多余观测值及较低的航向重叠度,导致高程精度较低(张永军等, 2005)。

4.2 阴影分布对密集匹配的影响

St-Onge等(2008)分别对地形平缓和地形起伏较大区域生成了Photo_CHM,在像元尺度上对Photo_CHM与LiDAR_CHM的相关性进行比较,地形平缓和地形起伏较大区域的R2在3 m空间分辨率下分别为0.79和0.46。本研究在3 m空间分辨率下的R2为0.52,介于二者之间,并且地形更加复杂。St-Onge等(2008)虽然得出像元阴影会降低R2,但是并没有将阴影对Photo_CHM的影响做进一步讨论。

当太阳以一定的高度角和方位角照射时,航片中向阳面与背阴面的光谱差异很大,山谷和山脊在航片上的差异也很明显,这会对立体像对的密集匹配造成挑战,进而影响Photo_CHM,具体如图 14所示。

图 14a为相同位置处(图 6b图 8b)Photo_CHM和LiDAR_CHM的横截面显示,蓝色曲线表示LiDAR_CHM,红色曲线表示Photo_CHM,A1、A2、A3为山脊背阴面阴影区域,B1、B2、B3、B4为向阳面。可以看出,LiDAR_CHM表现冠层细节的能力较Photo_CHM强。由图 14b可知,LiDAR_CHM相对于Photo_CHM在冠层高度上整体偏高。

图 14 LiDAR_CHM和Photo_CHM在林区的横截面显示(a)及其高度相对分布(b) Fig. 14 The transect display of LiDAR_CHM and Photo_CHM (a) and the relative distribution of height class(b)

从整体上看,除B2部分区域外,向阳面的Photo_CHM较背阴面的Photo_CHM和LiDAR_CHM更加吻合。对于阴影区域,A1、A2和A3部分区域出现明显的冠层高度被低估情况,其他区域可准确刻画冠层高度,这说明阴影会导致Photo_CHM被明显低估。结合密集匹配算法的特点可以推断出,利用航空像片估测平缓地区单一树种的上层木树高具有更好的效果。

图 9b图 14b的主要区别是Photo_CHM在冠层高度上整体偏高,出现这种情况的原因在于本研究区16块样地中4、6和12号样地各有一部分位于阴影处,其余都位于向阳面,存在Photo_CHM被低估的情况较少。结合图 13可知,7号样地Photo_CHM明显高于LiDAR_CHM,这是样地内Photo_CHM整体偏高的主要原因,因此结合坡度等因素对其分析是下一步的研究方向。

4.3 限制性因素及误差分析

因为历史航片的空间分辨率为1 m,所以如果匹配时产生几个像元间的匹配误差,则会对森林CHM产生很大影响。且历史航片是以胶片为存储介质的,影像质量随时间下降,纹理性低、对比度差,也会影响密集匹配效果。本研究历史航片是由黑白胶片扫描所得,由于胶片感光性限制等原因,造成航片中地物色彩单一,不利于寻找地物控制点,进而影响到了航片的绝对定向。

历史航片拍摄过程中,由于飞行平台高度较高,导致影像分辨率低,在获取小区域的精确冠层高度过程中会产生尺度误差。当摄影比例尺为1:20 000时,会有1~1.5 m高程误差(Micheletti et al., 2015),本研究历史航片的摄影比例尺为1:32 000和1:35 000,高程误差在2 m左右。

由于较高的飞行高度及起伏较大的地形,对历史航片进行密集匹配后,在研究区通常会出现比实际地形或高或低的匹配结果,需要人工进行编辑。人工编辑作为一种误差源也会对匹配结果产生影响,进而增加大范围森林冠层高度的监测难度。因此,需要发展一种适应于密集森林、低纹理航片的精确匹配算法,以便对大区域森林进行高效率密集匹配。

另外,受较长的时间跨度和拍摄任务多样性等因素影响,历史航片获取的季节时相往往不一致。本研究2期历史航片是不同季节拍摄获取的,1996年航片拍摄时间为5月,2004年航片拍摄时间为12月,太阳方位角和高度角显著不同,使得叶片反射光面积和颜色存在差异(蔡文峰等, 2010)。

拍摄视角也会对Photo_DSM估测造成影响。航空摄影是中心投影方式,在1996年历史航片研究区位于上边缘,在2004年历史航片中则接近中心投影点,反映在航片中则是研究区位置偏移及高度错估,对研究区Photo_DSM高程产生影响(陈德良, 1998)。

5 结论

1) 数字航片结合外方位数据,获得冠层表面Photo_DSM,利用LiDAR可穿透冠层到达地面的特殊性质获得地形数据,既可以得到准确冠层高度,又可以作为长期观测的数据分析基础。

2) 对于研究区多山的限制性因素,利用历史航片能准确定量反映山脊向阳面的森林冠层高度变化,对于山脊阴影和山谷处,则会出现冠层高度被低估的情况。太阳高度角和方位角作为影响阴影分布的重要因素,需要进行考虑和分析。

3) 本研究的关键技术环节为密集森林区域的密集匹配和影像配准,降低其可弥补误差,估测森林冠层高度会出现更好的结果。

参考文献(References)
蔡文峰, 李凤日. 2010. 基于VrituoZo系统对林木冠幅信息的提取. 森林工程, 26(2): 4-7.
(Cai W F, Li F R. 2010. Extraction of tree' crown diameter information based on VirtuoZo system. Forest Engineering, 26(2): 4-7. DOI:10.3969/j.issn.1001-005X.2010.02.002 [in Chinese])
陈德良. 1998. 航空摄影测量像片的几何关系分析. 福建农业大学学报, 27(2): 206-210.
(Chen D L. 1998. Analysis of geometric relationships of aerial photographs. Journal of FuJian Agricultural University, 27(2): 206-210. [in Chinese])
樊仲谋. 2015.摄影测树原理与技术方法研究.北京: 北京林业大学博士学位论文.
(Fan Z M. 2015. Photogrammetric principles and techniques for tree measurement. Beijing: PhD thesis of Beijing Forestry University.) http://d.wanfangdata.com.cn/Thesis/Y2850928
宫鹏, 梅雪良, Biging G S, 等. 1999. 利用数字摄影测量探测橡树草原变化. 遥感学报, 3(4): 285-289.
(Gong P, Mei X L, Biging G S, et al. 1999. Monitoring oak woodland change using digital photogrammetry. Journal of Remote Sensing, 3(4): 285-289. [in Chinese])
李明泽, 范文义, 张元元. 2009. 基于全数字摄影测量的林分立木高度量测. 北京林业大学学报, 31(2): 74-79.
(Li M Z, Fan W Y, Zhang Y Y. 2009. Measuring heights of standing trees based on digital photogrammetry. Journal of Beijing Forestry University, 31(2): 74-79. DOI:10.3321/j.issn:1000-1522.2009.02.011 [in Chinese])
吕丽萍. 2006. 双像法与单片法生成DOM的比较. 工程地球物理学报, 3(1): 45-48.
(Lü L P. 2006. Comparison of DOM production between stereo pair plotting and single photodifferential correction. Chinese Journal of Engineering Geophysics, 3(1): 45-48. DOI:10.3969/j.issn.1672-7940.2006.01.008 [in Chinese])
庞勇, 李增元, 陈尔学, 等. 2005. 激光雷达技术及其在林业上的应用. 林业科学, 41(3): 129-136.
(Pang Y, Li Z Y, Chen E X, et al. 2005. LiDAR remote sensing technology and its application in forestry. Scientia Silvae Sinicae, 41(3): 129-136. DOI:10.3321/j.issn:1001-7488.2005.03.022 [in Chinese])
汪承义. 2007.基于航空激光雷达数据的建筑物重建技术研究.北京: 中国科学院研究生院博士学位论文.
(Wang C Y. 2007. Research on building reconstruction based on airborne LiDAR Data. Beijing: PhD thesis of Beijing Graduate University of Chinese Academy of Sciences.) http://d.wanfangdata.com.cn/Thesis/Y1143486
汪少俊. 2015. 绩溪县杉木生长规律研究. 安徽林业科技, 41(3): 12-16.
(Wang S J. 2015. Analysis on the growth rhythm of Cunninghamia lanceolata in Jixi County. Anhui Forestry Science and Technology, 41(3): 12-16. DOI:10.3969/j.issn.2095-0152.2015.03.004 [in Chinese])
王佳, 冯仲科. 2011. 航空数字摄影测量对林分立木测高及精度分析. 测绘科学, 36(6): 77-79.
(Wang J, Feng Z K. 2011. Stand tree height measured by digital photogrammetry. Science of Surveying and Mapping, 36(6): 77-79. [in Chinese])
杨战辉, 张力. 1998. 用VirtuoZo数字摄影测量工作站生产DEM、DOM的试验. 测绘通报, (11): 10-11, 30.
(Yang Z H, Zhang L. 1998. Experiment of producing DEM, DOM with VirtuoZo digital photogrammetry workstation. Bulletin of Surveying and Mapping, (11): 10-11, 30. DOI:10.3969/j.issn.0494-0911.1998.11.003 [in Chinese])
袁晓红, 李际平. 2012. 杉木人工林南北坡向树高胸径生长曲线研究. 西北林学院学报, 27(2): 180-183.
(Yuan X H, Li J P. 2012. Height, DBH growth model of fir artificial forest on northern and southern slopes. Journal of Northwest Forestry University, 27(2): 180-183. DOI:10.3969/j.issn.1001-7461.2012.02.36 [in Chinese])
张剑清. 2011. 摄影测量学. 武汉: 武汉大学出版社, 150-171.
(Zhang J Q. 2011. Photogrameetry. Wuhan: Wuhan University Press, 150-171. [in Chinese])
张永军, 张勇. 2005. 大重叠度影像的相对定向与前方交会精度分析. 武汉大学学报:信息科学版, 30(2): 126-130.
(Zhang Y J, Zhang Y. 2005. Analysis of precision of relative orientation and forward intersection with high-overlap images. Editorial Board of Geomatics and Information Science of Wuhan University, 30(2): 126-130. [in Chinese])
智长贵, 秦学军, 韩爱惠. 2006. 基于光束法空中三角测量理论的树高测量方法研究. 林业资源管理, 8(4): 92-95.
(Zhi C G, Qin X J, Han A H. 2006. Study on measurement of tree height based on theory of bundle aerotriangulation. Forest Resources Management, 8(4): 92-95. DOI:10.3969/j.issn.1002-6622.2006.04.021 [in Chinese])
Bohlin J, Wallerman J, Fransson J E S. 2016. Deciduous forest mapping using change detection of multi-temporal canopy height models from aerial images acquired at leaf-on and leaf-off conditions. Scandinavian Journal of Forest Research, 31(5): 1-27.
Ginzler C, Hobi M. 2015. Countrywide stereo-image matching for updating digital surface models in the framework of the Swiss national forest inventory. Remote Sensing, 7(4): 4343-4370. DOI:10.3390/rs70404343
Harwin S, Lucieer A. 2012. Assessing the accuracy of georeferenced point clouds produced via multi-view stereopsis from unmanned aerial vehicle (UAV) imagery. Remote Sensing, 4(6): 1573-1599. DOI:10.3390/rs4061573
Lisein J, Pierrot-Deseilligny M, Bonnet S, et al. 2013. A photogrammetric workflow for the creation of a forest canopy height model from small unmanned aerial system imagery. Forests, 4(4): 922-944. DOI:10.3390/f4040922
Micheletti N, Lane S N, Chandler J H. 2015. Application of archival aerial photogrammetry to quantify climate forcing of alpine landscapes. The Photogrammetric Record, 30(15): 143-165.
Nurminen K, Litkey P, Honkavaara E, et al. 2015. Automation aspects for the georeferencing of photogrammetric aerial image archives in forested scenes. Remote Sensing, 7(2): 1565-1593. DOI:10.3390/rs70201565
Pang Y, Li Z Y, Ju H B, et al. 2016. LiCHy:the CAF's LiDAR, CCD and hyperspectral integrated airborne observation system. Remote Sensing, 8(5): 398. DOI:10.3390/rs8050398
St-Onge B, Audet F A, Bégin J. 2015. Characterizing the height structure and composition of a boreal forest using an individual tree crown approach applied to photogrammetric point clouds. Forests, 6(11): 3899-3922.
St-Onge B, Jumelet J, Cobello M, et al. 2011. Measuring individual tree height using a combination of stereo photogrammetry and LiDAR. Canadian Journal of Forest Research, 34(10): 2122-2130.
St-Onge B, Véga C, Fournier R A, et al. 2008. Mapping canopy height using a combination of digital stereo photogrammetry and LiDAR. International Journal of Remote Sensing, 29(11): 3343-3364. DOI:10.1080/01431160701469040
Véga C, St-Onge B. 2009. Mapping site index and age by linking a time of canopy height models with growth curves. Forest Ecology and Management, 257(3): 951-959. DOI:10.1016/j.foreco.2008.10.029
Véga C, St-Onge B. 2013. Height growth reconstruction of a boreal forest canopy over a period of 58 years using a combination of photogrammetric and LiDAR models. Remote Sensing of Environment, 112(4): 1784-1794.
White J C, Wulder M A, Vastaranta M, et al. 2013. The utility of image-based points clouds for forest inventory:a comparison with airborne laser scanning. Forests, 4(3): 518-536. DOI:10.3390/f4030518
Zhang J Q, Zhang Z X, Shen W M, et al. 1996. VirtuoZo digital photogrammetry system and its theoretical foundation and key algorithms. International Archives of Photogrammetry and Remote Sensing,, XXXI(part B2): 424-429.