联合资源三号与机载LiDAR的林分平均树高估测

邢艳秋 ，张锦绣 ，陈世培 ，高 立 ，关 雷 ，郭慧宇

（1.东北林业大学 森林作业与环境研究中心，黑龙江 哈尔滨 150040；2.黑龙江第三测绘工程院，黑龙江 哈尔滨 150025）

林分平均树高，作为重要的森林测树因子，既能反应林分的生长状况，又能为森林资源的经营与管理提供基础数据。随着航空摄影测量技术与遥感技术的不断革新，新的技术手段和遥感数据在森林资源监测应用中广泛渗透[1-3]。

资源三号（ZY-3）卫星作为我国首颗对地进行立体摄影测量的民用测绘卫星，具有地面测高能力，目前主要用于1∶50 000地形测绘。ZY-3在无控制点情况下，定位精度优于15 m，少量控制点的情况下高程精度优于3 m，平面精度优于4 m[4]。将ZY-3立体像对数据应用于林业研究可以获得大区域连续的森林冠层数字表面模型（digital surface model, DSM），能够有效避免离散点因空间插值形成连续DSM时所引入的误差。与此同时，ZY-3携带的多光谱信息可以用于辅助反映森林植被信息[5]。目前，已有研究人员利用摄影测量手段来提取树高，该类树高提取主要针对单木，且数据源多为航片。如：胥利红利用航片进行解析空中三角测量，完成相对定向和绝对定向，配合影像判读森林区域，所得树高估测模型的决定系数为0.903的树高估测模型，模型精度为87.22%[6]。王佳等以航片为数据源建立可量测的立体模型，进而从模型中提取单木树高，最终所得模型估测精度为93.96%[7]。李明泽等同样以航片为数据源，利用数字摄影测量平台对立木高度进行量测，结果表明，此方法的估测精度为90.92%[8]。虽然上述研究结果表明，基于遥感航片数据能够精确提取森林立木树高，但目前基于ZY-3航空立体影像，且在无控制点的条件下大范围获取林分平均树高的研究相对少见。这可能主要因为在森林冠层相对较密的林区，被动光学信号较难穿透森林冠层到达地面，因此不能获得林下高精度的地形数据[9]。

机载激光雷达（light detection and ranging，LiDAR）作为一种新型主动遥感技术，既能获取地物的水平结构信息，同时又能获得地物的垂直结构信息，因而可以实现区域地形的高精度模拟[10]。如：吴芳等利用机载LiDAR点云数据，通过对植被点数据进行分类移除，最终对植被覆盖相对较高的矿山地面起伏情况进行模拟，生成数字高程模型（digital elevation model, DEM），结果DEM高程误差为0.015 m[11]。尹艳豹等通对激光点云回波数据判别分析进而提取林地DEM，验证得DEM高差总平均值为0.24 m，方差为0.05 m[12]。上述研究结果表明，机载LiDAR点云数据可以生成高精度DEM，弥补ZY-3立体影像在植被覆盖区域提取DEM的不足。同时利用高精度的DEM作为高程基准采集虚拟控制点进行平差能提高ZY-3立体影像数据的定位精度[13]。

因此，本研究通过联合ZY-3数据与机载LiDAR点云数据，以期实现无实测控制点条件下森林林分平均树高的高精度估测，为我国测绘卫星数据在林业结构参数估测研究方面提供基础数据及技术支撑。

1 研究区概况

本研究实验区位于吉林省长春市净月潭国家森林公园，其地处长白山余脉丘陵地带（125°24′～125°28′E，43°45′～43°48′N），高程为175～406 m，森林由人工林与天然次生林组成，总面积达92.52 km2。树种以落叶松Larix gmelinii、樟子松Pinus sylvestris和蒙古栎Quercus mongolica为主。气候为温带半湿润季风气候，区域全年温差约50 ℃，全年平均降水量约654.3 mm。

2 数据及预处理

2.1 野外调查数据

2012年10月对研究区进行了野外数据采集。根据林分的密度不同，在进行野外数据采集时布设了半径分别为10 m和15 m的圆形样地。利用Trimble GeoXT6000 GPS对样地中心以及样地内所有活立木进行定位，定位精度约0.5 m。利用超声波测高仪对样地内所有单木树高使用进行3次测量，取3次量测均值作为单木树高的实测值。同时利用胸径尺对样地内的所有单木胸径进行3次测量，同样取测量的均值作为胸径的实测值。最终在研究区共采集37个样方（如图1所示），其中包括34个针叶林（23个樟子松和11个落叶松）和3个阔叶林（2个蒙古栎和1个杨树）。

因遥感林分平均树高估测研究中多采用胸高断面积加权平均高[14]，因此本研究同样采用胸高断面积加权平均高作为林分平均树高，其具体计算公式如式（1）所示：

式中，HL为第L（L=1，2，…，n）个样地的实测树高，hi为第i株立木的树高，gi为第i株立木的胸高断面积，n为样地内立木株数。

2.2 ZY-3影像数据

ZY-3为太阳同步卫星，高度约505km，回归周期59 d，重访周期5 d。本研究使用的ZY-3影像数据于2013年5月22日采集得到。

三线阵立体影像由3台具有立体视角的光学相机沿卫星飞行方向推扫成像，分别是前视22°和后视22°全色延时积分成像（TDI CCD）相机（分辨率3.5 m，幅宽52 km），下视全色TDI CCD相机（分辨率2.1 m，幅宽51 km）。谱段范围0.5～0.8 μm，侧摆能力±32°。

多光谱影像由一台分辨率5.8 m，幅宽51 km， 摆 能 力 ±32°， 含有蓝（0.45～0.52 μm）、绿（0.52～ 0.59 μm）、红（0.63～0.69 μm）和近红外（0.77～0.89 μm）四个通道的多光谱相机采集。

图 1 研究区域Fig.1 Study area

2.3 机载LiDAR点云数据

本研究使用的机载LiDAR离散点云数据是由Leica ALS70型传感器于2012年5月31日采集得到。传感器的飞行高度为560 m，扫描频率40.3 Hz，激光脉冲波长1 064 nm，光束发散角0.22 mrad，光斑直径0.28 m，点云密度约10 个/m2。

设置点云数据的强度范围为0～255，建立对应的灰度图像，提取与林区DEM“同源”的正射影像（分辨率1 m），作为ZY-3数据定向处理的平面控制基准。

3 研究方法

在林分平均树高估测过程中，需要解决的几个关键问题：①利用机载LiDAR数据提取DEM；②利用ZY-3数据提取DSM；③建立林分平均树高估测模型并评价模型精度。具体的技术流程如图2所示。

3.1 DEM提取方法

本研究获取的原始LiDAR数据通常包含部分噪声点，因此在使用之前需对其进行去噪处理。利用不规则三角格网（Triangulated irregular network，TIN）算法将去噪后的LiDAR点云数据分为地面点及非地面点，并利用TIN算法将离散地面点生成能够有利描述地形细节的不规则三角网，最后内插生成数据结构简单，便于使用和管理的规则格网林区DEM（分辨率1 m）。

图 2 技术流程图Fig.2 Flow chart of the method

3.2 DSM提取方法及高程精度评价

3.2.1 DSM提取方法

对立体影像进行绝对定向，相对定向，平差处理，最终利用多类像对组合提取包含地面林分树高的林区DSM。

为完成绝对定向，以林区DEM为高程基准，“同源”正射影像为平面基准，采用快速傅里叶相位变换匹配法进行“虚拟控制点”采集，根据ZY-3的RPC参数文件建立有理函数模型，使用“虚拟控制点”将立体影像进行移动、转动和放缩，建立像面点坐标与地面点坐标对应关系。在控制点采集过程中，控制点搜索半径不小于参考点如图3a所示相对于三视影像上同名点各自偏移距离的最大值，本研究中参考点相对于前视，下视和后视影像上同名点偏移距离分别如图3b，图3c和图3d所示，最终确定搜索半径为200 m，同时设置搜索拒绝阈值为0.75。

图 3 同名点偏移情况Fig.3 Homonymy point offsets

在绝对定向基础上，以林区DEM为高程控制基准，采用快速傅里叶相位变换匹配法进行连接点采集，以光栅格网法进行采样，还原立体像对中左右影像在成像时的相对位置，获得所有成对的同名光线相交，完成相对定向。在连接点采集过程中，设置搜索半径为100 m，拒绝阈值0.76。同名点采集完成后利用最小二乘原理进行平差，减小像点与实际地面点对应的差异。

根据空中三角测量原理，利用“虚拟控制点”求出地面加密点的物方三维坐标，将三类立体像对进行组合，基于离散加密点并以“丘陵”地貌类型为约束，提取前下视DSM，下后视DSM和前后视DSM，通过自适应方法将三类DSM最优化合成三视DSM。

3.2.2 DSM高程精度评价方法

理想情况下，在无地物遮挡的地面点，DEM与DSM的高程应保持一致。由于ZY-3数据与机载LiDAR数据的类型、成像时间、分辨率及处理方法等差异，加之地形起伏变化，实验结果不可能达到理想情况。同时，林区DSM与高程基准林区DEM的差值为图上量测树高值，该值的准确性直接影像树高估测模型精度，故需要对DSM数据进行高程精度评价，选择高程精度最佳的DSM进行平均树高估测。在无覆盖的城区道路面提取同名点，同时，为保证研究区高程精度分析结果的准确性，在林区宽阔无植被覆盖的道路上也选取少量的同名点进行高程比较。求同名点DSM与DEM的差值，计算DSM数据高程均方根误差（Root Mean Square Error, RMSE），高程平均正误差值和高程平均负误差值。RMSE和高程平均正、负误差绝对值越小，说明DSM高程精度越好。

3.3 植被指数提取方法

ZY-3多光谱影像以正射校正后的下视影像为参考影像进行几何校正，将校正后的多光谱影像进行波段组合计算，提取10个植被指数如下所示：修正土壤调整植被指数，比差植被指数，修正单比植被指数，修正三角植被指数，非线性指数，优化土壤调整植被指数，归一化差值植被指数，绿色归一化植被指数，三角植被指数，红外百分比植被指数。

3.4 平均树高估测模型建立与评价

3.4.1 基于回归分析方法建立平均树高估测模型

以图上树高量测值为自变量，实测树高值为因变量进行一元回归分析，建立初始树高估测模型。图上树高量测值如式（2）所示：

式中，Hmea为图上树高量测值，Hdsm为DSM算数平均值，Hdem为DEM算数平均值。

在初始树高估测模型的基础上，对各项植被指数进行T检验，认为T值对应概率小于0.05的植被指数引入回归模型是有意义的，并将该指数作为回归分析自变量，最后得到改进树高估测模型如式（3）所示：

式中，Hpre为树高改进预测值，Hmea为树高图上量测值，Vk为植被指数值，k为植被指数个数（1≤k≤10），a，bk和c为待求解的改进树高预测模型的回归系数。

3.4.2 利用DSM高程误差修正树高估测模型

ZY-3提取的DSM相对于高程基准存在高程误差。在DSM中地面高程高于基准DEM的区域ZY-3量测高度值会高于实测值，而在DSM中地面高程低于基准DEM的区域ZY-3量测高度值会低于实测值。为了减小DSM高程误差对预测结果的影响，根据高程精度评价结果，利用DSM高程误差对回归树高估测模型的预测值进行修正，在树高预测值低于实测树高的样地加上DSM高程负误差绝对值，在树高预测值高于实测树高的样地减去DSM高程正误差绝对值。进而得到修正初始树高估测模型与修正改进树高估测模型。

3.4.3 模型评价

26个样方数据用于建立树高估测模型，剩余的11个样方数据用于模型精度评价，其中模型拟合优度用AdjR2来评价，预测精度用P值评价。P值公式如式（4）所示：

式中，H为野外实测树高，hˆ为模型预测树高。

4 结果与分析

4.1 DEM提取结果

从机载LiDAR点云中提取了分辨率为1 m的DEM图像，高程范围为180.9～385.8 m。经目视分析，DEM细节清晰，在森林覆盖区域图像上直接反应林下地面起伏。各验证样地的DEM平均值如表1所示。

4.2 DSM提取结果及高程精度评价

绝对定向过程中共采集控制点59个，控制点X，Y方向上的均方根误差分别为1.27 m和1.29 m。相对定向过程中采集二度和三度连接点共275对，连接点X，Y方向上的均方根误差分别为0.53 m和0.88 m。

最终从ZY-3数据提取了4个DSM数据（分辨率3.5 m），在全区共选取了168个同名点对DSM数据进行高程精度评价。研究表明，四类DSM数据的高程正负误差分布均表现为无规律性，其中三视DSM具有最小高程均方根误差1.433 m，最小高程差方差1.871，最小高程正误差0.851 m和最小高程负误差-1.175 m。经目视分析，三视DSM中建筑物轮廓清晰，道路线条流畅，DSM过渡光滑，很少有异常剧烈的突起和凹陷，能够从DSM中再现平面地物起伏的形态。认为三视DSM为林分平均树高估测的最佳森林冠顶数字表面模型。各验证样地的DSM平均值如表1所示。

4.3 平均树高模型结果及精度评价

经一元回归建模后得到初始树高估测模型，如式（5）所示：

经多元回归建模后得到改进树高估测模型，如式（6）所示。10个植被指数经t检验后，仅有优化土壤调整植被指数（Optimized Soil Adjusted Vegetation Index, OSAVI）与修正土壤调整植被指 数（Modified Soil Adjusted Vegetation Index,MSAVI）对应的T检验值概率显著小于0.05，适合被引入改进树高估测模型。这两个植被指数均由对植被敏感的红波段与近红外波段组合计算得到，能增加植被信号的动态范围，同时进一步最小化土壤背景的影响。

表1 实验结果†Table 1 Experiment results m

利用11个验证样地来验证模型精度。树高实测值，图上量测值，4个模型的预测值以及精度如表1所示。4个模型对应的树高预测趋势及验证样地树高误差按模型精度由小到大排列，如图4所示。

初始树高估测模型精度为85.48%，其中，针叶林精度为86.46%，阔叶林精度为75.70%。改进树高估测模型精度为87.47%，其中，针叶林精度为88.58%，阔叶林精度为76.41%。修正初始树高估测模型精度提高到92.10%，其中，针叶林精度提高到93.32%，阔叶林精度提高到79.93%。修正改进树高估测模型精度提高到93.29%，其中，针叶林精度提高到94.55%，阔叶林精度提高到80.63%。实验结果表明，经高程误差修正的改进树高估测模型为最佳树高估测模型。

图4 树高预测趋势及树高误差Fig.4 Predicting trends of mean canopy height and errors

最佳模型对针叶林树高的预测精度较高，误差约为1 m，而对阔叶林树高的预测精度相对较小，误差约为5 m。考虑到阔叶林样地数量较少，树高估测模型包含较少的阔叶林高度和光谱信息，所以导致阔叶林树高的预测精度较低。

野外样地实测数据是在10月采集的，而用于提取树高的ZY-3数据是在5月收集的，林分在不同时间其生长状态也不同，数据采集的时间差会导致预测树高都要略低于实测树高，这个现象是合理的。

5 结论与讨论

本研究以吉林省长春市净月潭国家森林公园林区为研究区，基于ZY-3影像数据、机载LiDAR数据和野外实测数据进行林分平均树高估测研究，所得结论主要如下：

（1）从机载LiDAR点云中提取的“同源”高精度DEM与正射影像数据，作为ZY-3定向基准可在无控制点条件下提高其定位测高能力。

（2）在ZY-3立体影像提取的多个DSM中，经目视对比分析及高程精度评价分析，确定三视DSM的高程精度优于二视DSM。

（3）经高程误差修正后的改进树高估测模型精度最高，模型的调整决定系数AdjR2=0.913，预测精度为93.29%。

相比王佳，李明泽，胥利红等学者利用航片数据进行对点的森林立木树高提取研究，本研究利用ZY-3航空立体影像与机载LiDAR点云数据，在无控制点条件下进行对面的林分平均树高提取，既保证了较高的树高模型估测精度，又实现了较大连续区域快速的林分平均树高估测，亦成功实验了将测绘数据应用到林分测高。由于本研究中阔叶林样本数量较少，树高估测模型对该类林分建立的高度预测能力较针叶林林分弱，因此在下一步的研究中将加大阔叶林样的比例，甚至增加针阔混交类型的林分，提高估测模型对各类林分平均树高的估测精度以及模型自身的普适能力。