CN104849722A - 一种LiDAR波形综合特征的单木识别方法 - Google Patents

Abstract

本发明公开了一种LiDAR波形综合特征的单木识别方法，借助机载小光斑全波形LiDAR传感器进行数据采集；LiDAR波形数据预处理：单木定位和冠幅提取；基于发射能量及传感器与地物的距离信息对LiDAR波形数据进行校正；构建体元框架并进行LiDAR波形的结构化分解；提取单木的波形特征变量；提取单木的点云特征变量；使用随机森林方法筛选最优特征变量并进行树种分类。本发明的验证结果表明，与其他使用遥感方法进行树种分类的方法相比，总体精度提升了15%左右；Kappa系数提升了0.13左右。

Description

一种LiDAR波形综合特征的单木识别方法

技术领域

本发明涉及森林资源管理与保护技术领域，具体涉及一种LiDAR波形综合特征的单木识别方法。

背景技术

精确的树种分类对于森林资源调查、动态监测以及生物多样性研究以及模拟指定树种的单木生长有重要意义。同时，这些信息也可以用于森林资源调查、动态监测以及生物多样性研究，从而为小尺度和中等尺度的森林资源规划和集约管理提供精确的数据支撑。常规的森林树种调查方法主要依赖于野外调查及利用大比例尺航片判读等，其精度往往不高，且难于在大区域上实用化推广。激光雷达(LiDAR，Light Detection and Ranging)是通过发射激光束击打物体表面并分析其返回信号的一种主动遥感技术。通过LiDAR可获得高精度的地球表面及地表实体的高度信息，如地形和植被等可达到亚米级的垂直精度。现有研究表明，LiDAR可穿透森林冠层获得其三维结构特征，特别适合植被覆盖度高且森林结构复杂地区森林信息提取。

近年来基于小光斑LiDAR数据进行树种分类研究为：等2009年在《RemoteSensing of Environment》第113卷上发表的“Classifying species of individual trees byintensity and structure features derived from airborne laser scanner data”，该研究在已有的单木位置和冠幅信息基础上，提取了单木的多个高度(如最大高，平均高和高度百分位数等)及统计(如高度分布的峰度和偏度等)特征变量，并据此对挪威寒带森林中的针叶和阔叶树种进行了分类。Kim等2009年在《Remote Sensing ofEnvironment》第113卷上发表的“Tree species differentiation using intensity data derivedfrom leaf-on and leaf-off airborne laser scanner data”，该研究分别从“未落叶”和“落叶”两种状态下提取LiDAR点云中的强度信息，并融合这些信息用于北美温带森林中的针叶和阔叶树种分类。Heinzel等2011年在《International Journal of Applied EarthObservation and Geoinformation》第13卷上发表的“Exploring full-waveform LiDARparameters for tree species classification”，该研究结合地面调查数据，通过从LiDAR波形数据中提取强度，波宽以及目标物返回点数等信息并结合线性可分性(LDA)分析进行变量筛选，对欧洲中部的6个树种进行了遥感分类。然而，以上方法大多适用于树种组成较为简单的森林分类研究，在林分组成和结构复杂的森林中分类精度不高。且仅从单一的角度去挖掘LiDAR数据，即并未将“点云”数据中包含的三维空间信息及“波形”数据中包含的几何与能量信息相结合，未能充分挖掘LiDAR数据的潜力。

近年来，机载小光斑“全波形”LiDAR传感器逐步商用化并开始应用于林业研究中。该技术可获得森林冠层完整的后向散射信号，并记录了更为详细的几何和物理能量信息，从而一定程度上弥补了传统LiDAR技术的弱点。而且，借助特定的波形信号处理及信息提取方法，可从“全波形”数据中挖掘更为丰富的森林冠层描述特征变量，理论上能够更为丰富地反映不同树种的结构特征。同时，全波形LiDAR技术并未摒弃传统技术的优点，而是对其进行了提升，即借助特定方法可以从中同时提取出波形和点云特征变量，两者融合用于后续的遥感分析。故本发明将从全波形LiDAR数据中提取与森林冠层结构相关度高的波形和点云特征变量，并用于亚热带地区的典型树种分类。其创新点和特色如下：1)通过综合特征变量，从冠层三维空间和几何能量维度挖掘全波形LiDAR数据中包含的树木生物物理特性信息；2)借助精度均减系数(MDA)对以上特征变量进行重要性判定，筛选最优分类因子，从而利于机理解释、方法移植及大尺度推广和应用。

发明内容

发明目的：针对现有技术中存在的不足，本发明提出一种LiDAR波形综合特征的单木识别方法，即首先从全波形LiDAR数据中提取点云，然后融合点云和波形特征变量用于的单木识别方法；有效提高树种分类的总体精度，易于推广应用等特点。

技术方案：为了实现上述发明目的，本发明采用的技术方案为：

一种LiDAR波形综合特征的单木识别方法，包括以下步骤：

1)借助机载小光斑全波形LiDAR传感器进行数据采集；传感器记录每束激光脉冲返回的完整波形信息；

2)LiDAR波形数据预处理：

A)噪声水平估计和数据平滑：首先把原始数据转换到频率域，再将频率较高的低值部分作为噪声水平的判断标准；然后选用高斯滤波器进行平滑。

B)高斯拟合分解及波形数据点云化：基于回波数据是多个高斯函数的累加这一假设，对波形数据采用非线性最小二乘法进行拟合；然后通过局部最大峰值检测滤波算法从处理后的波形数据上提取离散点云，每个离散点中记录了返回信号的能量和振幅信息；

C)生成数字地形：首先对从波形数据中提取出离散点云进行分类，然后对末次回波进行Kraus滤波处理用以去除非地面点，最后使用滤波后的末次回波数据并借助自然邻近法插值生成数字地形模型DTM；

3)单木定位和冠幅提取：

A)对地面以上点云进行中值滤波，然后将点云中的高度信息栅格化生成数字表面模型DSM；将DSM减去数字地形模型，得到归一化植被高度CHM；

B)通过局部最大值法确定单木树顶所在位置；

C)单木冠幅的确定：首先以树顶为中心拟合16个半径方向上的冠幅剖面，然后计算到局部最小值的水平距离，最后将这些距离值进行水平方向上的平均从而得到冠幅半径；

4)基于发射能量及传感器与地物的距离信息对LiDAR波形数据进行校正；

5)构建体元框架并进行LiDAR波形的结构化分解

首先利用三维体积单元划分地表以上的三维空间(每个单元的大小依据脉冲采样间隔及光斑大小进行设置)；然后在每个单元内，汇总穿过其中的最大脉冲能量并基于数字地形模型进行高度归一化；最后将每列体元柱信息汇总得到伪垂直波；

6)提取单木的波形特征变量

a)单木分割的树冠范围内，针对每组高斯分解后的每个波形，计算其能量(ER)、振幅(WT)和回波次数(NT)的均值(μ)及标准差(σ)，以此作为单木的特征变量(第一组)；

b)基于伪垂直波提取特征变量(第二组)：HOME(成功应用于估算冠层垂直结构和郁闭度)、波形距离(WD，通常与树高联系紧密)、高度中位数比(该变量对冠层结构变化很敏感)、波峰数量、返回波形能量、冠层外层粗糙度(ROUGH，用于描述植被表面空间结构)、前坡度角(FS，用于描述植被冠层上部的变异性)和VDR(冠层高度和HOME的差值除以冠层高)。以上变量也分别在单木冠幅内计算均值(μ)和标准差(σ)；

7)提取单木的点云特征变量

在单木冠幅内提取了点云特征变量：a.高度百分位数变量组(h₁₀,h₂₅,h₅₀,h₇₅,h₉₀,h₉₅)，即激光返回点的高度分布百分位数(10th,25th,50th,75th,90th,95th)；b.冠层返回点云密度变量组(d₂,d₄,d₆,d₈)，即在各百分位高度等级(20th,40th,60th,80th)以上的激光返回点在所有返回点中所占的百分比；c.最小和平均高度(h_min/h_mean)，即归一化高度的最小(或平均)值；d.高度变异系数(h_cv)，即归一化高度的变异系数(标准差与平均数的比值)；e.覆盖度(CC_2m/CC_mean)，即高于2m(或平均高)的激光返回点所占的百分比；

8)使用随机森林方法筛选最优特征变量并进行树种分类

A)随机森林分类是由很多决策树分类模型{h(X,Θ_k),k＝1,2,...}组成的组合分类模型，且参数集{Θ_k}是独立同分布的随机向量，在给定自变量X下，每个决策树分类模型都由一票投票权来选择最优的分类结果；通过k轮训练，得到一个分类模型序列{h₁(X),h₂(X),…h_k(X)}，再用它们构成一个多分类模型系统，其最终的分类决策函数为：

$H\left(x\right)=\arg \underset{Y}{\max }\underset{i=1}{\overset{k}{\mathrm{\Σ}}}I(h_i\left(x\right)=Y)---\left(4\right)$

其中，H(x)表示组合分类模型，h_i是单个决策树分类模型，Y表示输出变量，I为示性函数；

B)通过随机森林方法生成重要性变量Xm，该变量通过累积每个t节点的权重减少量p(t)Δi(s_t,t)来评估预测量Y，由此平均树N_T结构写为：

$\mathrm{IMP}\left(X_m\right)=\frac{1}{N_T}\underset{T}{\mathrm{\Σ}}\underset{t\∈T:\mathrm{\υ}\left(s_t\right)=X_m}{\mathrm{\Σ}}p\left(t\right)\mathrm{\Δi}(s_t,t)---\left(5\right)$

式中，p(t)是当样本到达t时N_t/N的比值，v(s_t)是s_t的分裂变量；

依据上述理论，精度均减系数(MDA)作为测试重要性的指标，表达了模型引入一个变量后均方根误差的减少量。

C)将随机森林模型中的决策树数量设置为100，每个分叉的变量数为4；将上述波形和点云特征变量作为输入参数代入模型；模型自动剔除引入参数直到分类总体精度趋于稳定；然后再依据MDA指数选择重要性最高的前3个特征变量作为最优特征变量；

D)使用筛选后的最优特征变量和随机森林分类器再次进行分类，并借助总体分类精度(即对角线像元数/总体像元数×100％)和卡帕系数这两个定量化指标对分类精度进行评价；卡帕系数计算公式：

$K_{\mathrm{hat}}=\frac{N\underset{i=1}{\overset{r}{\mathrm{\Σ}}}{x}_{\mathrm{ii}}-\underset{i=1}{\overset{r}{\mathrm{\Σ}}}\left(x_{i+}{x}_{+i}\right)}{N^2-\underset{i=1}{\overset{r}{\mathrm{\Σ}}}\left(x_{i+}{x}_{+i}\right)}---\left(6\right)$

式中，r为总的类别数，x_ii为对角线上的像元数，x_i+和x_+i是列和行的总像元，N是总像元数。

步骤1)中，所述的LiDAR传感器为奥地利的Riegl LMS-Q680i。

步骤1)中，所使用的遥感平台为运-5固定翼飞机。

步骤1)中，所述传感器的采样间隔为1ns。

步骤3)的B)中，通过圆形的搜索框在CHM上“滑行”遍历栅格图层，通过逐个比较搜索框内的高度信息来确定最大值点；搜索框的半径通过高度信息和参数β₀及β₁来确定；β₀及β₁则通过地面实测数据拟合来获取，公式为：

CW(m)＝β₀+β₁×h²      (1)

式中，CW为地面实测冠幅半径，h为树高(m)，β₀和β₁为模型参数。

步骤4)中校正的公式为：

$W_i^c=W_i/W^e---\left(2\right)$

$I_i^c=(I_i\·s_i^k)/(I^e\·s_0^k)---\left(3\right)$

上式中，为校正后返回波内每个高斯波的波宽，为校正后返回波内每个高斯波内每个高斯波的能量强度，W_i为原始返回波内每个高斯波的波宽，W^e为发射波的波宽，I_i为原始返回波的能量强度，I^e为发射波的能量强度，传感器到反射物体的距离，为标称距离，k为变化系数。

有益效果：与现有技术相比，本发明的特色和创新点如下：通过综合LiDAR特征变量，从冠层三维空间和几何能量维度挖掘全波形LiDAR数据中包含的树木生物物理特性信息；并借助精度均减系数(MDA)对以上特征变量进行重要性判定，筛选最优分类因子；利于机理解释、方法移植及大尺度推广和应用。实验验证结果表明，通过本发明对北亚热带天然次生混交林的主要树种进行树种分类，与其他使用遥感方法(注：对比条件为：单一遥感数据源；且在类似复杂森林类型内)进行树种分类的方法相比总体精度提升了15％左右；Kappa系数提升了0.13左右。

附图说明

图1是LiDAR点云数据预处理流程图；a从波形中提取的点云数据；b滤波后去除了地面以上点；c地面点插值得到数字地形模型(DTM)；d DSM减去数字地形模型(DTM)从而得到归一化植被高度(CHM)；

图2是单木冠幅内的点云和波形脉冲可视化图；

图3是部分波形特征变量(第二组)的提取原理示意图。

具体实施方式

下面结合具体实施例对本发明作进一步的说明。

实施例1

一种LiDAR波形综合特征的单木识别方法，以针对一个北亚热带天然次生混交林为主要森林类型的林区内树种分类为例。林区海拔20-261m，面积约1100公顷。主要树种为针叶的马尾松(Pinus massoniana)、杉木(Cunninghamia lanceolata)和湿地松(Pinus elliottii)，以及阔叶的麻栎(Quercus acutissima)、枫香(Liquidambarformosana)和冬青(Ilex chinensis)。林区内根据树种组成、林龄和立地状况等布设了12个方形样地(30×30m)，每个样地内人工判别了单木树种，并实测了胸径、树高和冠幅等森林参数；样地的中心通过差分GPS进行定位，每木的相对位置(即距离样地中心的水平距离和方向角)通过皮尺和森林罗盘仪测定(最后换算为每木的绝对位置坐标)。

1)LiDAR数据获取

借助奥地利Riegl LMS-Q680i的机载小光斑全波形LiDAR传感器进行数据采集。遥感平台为：运-5飞机(飞行高度900m，速度为：55m/s)。获得的LiDAR波形数据参数为：脉冲发射频率400kHz，扫描频率114lines/sec(扫描角±30°)，激光脉冲的光斑半径为45cm。传感器记录了每束激光脉冲返回的完整波形信息，采样间隔为1ns。获得数据的脉冲间距为0.48m，脉冲点密度2.1pulse/m²⁽航向重叠区域的脉冲点密度更高)。

2)LiDAR波形数据预处理

a)噪声水平估计和数据平滑。首先把原始数据转换到频率域，再将频率较高的低值部分作为噪声水平的判断标准。然后选用高斯滤波器进行平滑(内核设置：FWHM＝2.35×σ)，这是由于高斯滤波器在有效平滑数据的同时，还可以最大限度地保持原有曲线的趋势。

b)高斯拟合(分解)及波形数据点云化。基于回波数据是多个高斯函数的累加这一假设，对波形数据采用非线性最小二乘法进行拟合。然后通过局部最大峰值检测滤波算法从处理后的波形数据上提取离散点云(即通过移动窗体逐个判断某点周边的最近4个点的脉冲强度信息，中间值都高于周边的点即为所要提取的点)，每个离散点中记录了返回信号的能量和振幅信息。

c)生成数字地形。LiDAR数据高度归一化的目的是为了得到去除了地形影响的“真实”植被高度，通常采用原始LiDAR数据高度信息减去地形高度得到。因此，精确生成数字地形模型(DTM)是计算归一化植被高度的重要前提。首先对从波形数据中提取出离散点云进行分类，然后对末次回波进行Kraus滤波处理用以去除非地面点(图1.b)(在此基础上进一步采用中值滤波进行平滑，窗体大小3×3)，最后使用滤波后的末次回波数据并借助自然邻近法插值生成数字地形模型(图1.c)。

4)单木定位和冠幅提取

A)对地面以上点云进行中值滤波(窗体大小为3×3)，然后将点云中的高度信息栅格化生成数字表面模型(DSM)。将DSM减去数字地形模型(DTM)从而得到归一化植被高度(CHM)(图1.d)。

B)通过局部最大值法确定单木树顶(即最高点)所在位置。即通过圆形的搜索框在CHM上“滑行”遍历栅格图层，通过逐个比较搜索框内的高度信息来确定最大值点。搜索框的半径通过高度信息和参数β₀及β₁来确定。β₀及β₁则通过地面实测数据拟合来获取。其公式为：

CW(m)＝β₀+β₁×h²      (1)

其中CW为地面实测冠幅半径，h为树高(m)，β₀和β₁为模型参数。

C)单木冠幅则是借助冠幅半径来进行描述。其方法为首先以树顶为中心拟合16个半径方向上的冠幅剖面，然后计算到局部最小值的水平距离，最后将这些距离值进行水平方向上的平均从而得到冠幅半径。

根据所测森林结构参数，选取样地中的典型样本：马尾松(75棵)、杉木(69棵)、湿地松(65棵)、麻栎(65棵)、枫香(57棵)和冬青(45棵)，拟合模型得到β₀为1.781，β₁为0.029(模型的R²为0.63)。单木冠幅则是借助冠幅半径来进行描述。其方法为：首先以树顶为中心拟合16个半径方向上的冠幅剖面，然后计算到局部最小值的水平距离，最后将这些距离值进行水平方向上的平均从而得到冠幅半径。将LiDAR提取单木与地面实测数据进行空间位置对应“链接”(即LiDAR估算和地面实测的树冠中心水平位置在2m范围内则认定为匹配)，得到单木提取正确率为76％。单木树高的RMSE(均方根误差)为0.63m，冠幅的RMSE为0.55m。共提取(430棵)有效树：马尾松(85棵)、杉木(81棵)、湿地松(70棵)、麻栎(72棵)、枫香(73棵)和冬青(50棵)。

5)基于发射能量及(传感器与地物的)距离信息对LiDAR波形数据进行校正：

$W_i^c=W_i/W^e---\left(2\right)$

$I_i^c=(I_i\·s_i^k)/(I^e\·s_0^k)---\left(3\right)$

上式中，为校正后返回波内每个高斯波的波宽，为校正后返回波内每个高斯波内每个高斯波的能量强度，W_i为原始返回波内每个高斯波的波宽，W^e为发射波的波宽，I_i为原始返回波的能量强度，I^e为发射波的能量强度，传感器到反射物体的距离，为标称距离，k为变化系数。由于指定传感器的发射脉冲波宽和强度在一次飞行中基本恒定，故通过在多个条带中随机采集典型样本并做分析，确定W^e设置为发射波的波宽为3.5ns，I^e为200(相对值)，为传感器到反射物体的距离(即飞行高度H减去波形中锚点所在高度h然后除以扫描角θ的cos值)，为标称距离(本申请取平均飞行高度，即900m)，k为变化系数(本申请取经验值：2)。

6)构建体元框架并进行LiDAR波形的结构化分解。

首先利用三维体积单元划分地表以上的三维空间(每个单元的大小依据脉冲采样间隔及光斑大小设置为：0.3×0.3×0.5m)；然后在每个单元内，汇总穿过其中的最大脉冲能量并基于数字地形模型进行高度归一化。最后将每列体元柱信息汇总得到伪垂直波。

7)提取单木的波形特征变量

a)单木分割的树冠范围内，针对每组高斯分解后的每个波形，计算其能量(ER)、振幅(WT)和回波次数(NT)的均值(μ)及标准差(σ)，以此作为单木的特征变量(第一组)。

b)基于伪垂直波提取特征变量(第二组)：HOME(成功应用于估算冠层垂直结构和郁闭度)、波形距离(WD，通常与树高联系紧密)、高度中位数比(该变量对冠层结构变化很敏感)、波峰数量、返回波形能量、冠层外层粗糙度(ROUGH，用于描述植被表面空间结构)、前坡度角(FS，用于描述植被冠层上部的变异性)和VDR(冠层高度和HOME的差值除以冠层高)。以上变量也分别在单木冠幅内计算均值(μ)和标准差(σ)(部分特征变量的提取原理见图2)(计算方法详见表3)。

8)提取单木的点云特征变量

同样在单木冠幅内提取了点云特征变量：a.高度百分位数变量组(h₁₀,h₂₅,h₅₀,h₇₅,h₉₀,h₉₅)，即激光返回点的高度分布百分位数(10th,25th,50th,75th,90th,95th)；b.冠层返回点云密度变量组(d₂,d₄,d₆,d₈)，即在各百分位高度等级(20th,40th,60th,80th)以上的激光返回点在所有返回点中所占的百分比；c.最小和平均高度(h_min/h_mean)，即归一化高度的最小(或平均)值；d.高度变异系数(h_cv)，即归一化高度的变异系数(标准差与平均数的比值)；e.覆盖度(CC_2m/CC_mean)，即高于2m(或平均高)的激光返回点所占的百分比。

9)使用随机森林方法筛选最优特征变量并进行树种分类。

A)随机森林分类是由很多决策树分类模型{h(X,Θ_k),k＝1,2,...}组成的组合分类模型，且参数集{Θ_k}是独立同分布的随机向量，在给定自变量X下，每个决策树分类模型都由一票投票权来选择最优的分类结果。通过k轮训练，得到一个分类模型序列{h₁(X),h₂(X),…h_k(X)},再用它们构成一个多分类模型系统，其最终的分类决策函数为：

$H\left(x\right)=\arg \underset{Y}{\max }\underset{i=1}{\overset{k}{\mathrm{\Σ}}}I(h_i\left(x\right)=Y)---\left(4\right)$

其中，H(x)表示组合分类模型，h_i是单个决策树分类模型，Y表示输出变量，I为示性函数。

B)通过随机森林方法可以生成重要性变量Xm，该变量可以通过累积每个t节点的权重减少量p(t)Δi(s_t,t)来评估预测量Y，由此平均树N_T结构可写为：

$\mathrm{IMP}\left(X_m\right)=\frac{1}{N_T}\underset{T}{\mathrm{\Σ}}\underset{t\∈T:\mathrm{\υ}\left(s_t\right)=X_m}{\mathrm{\Σ}}p\left(t\right)\mathrm{\Δi}(s_t,t)---\left(5\right)$

C)其中，p(t)是当样本到达t时N_t/N的比值，而v(s_t)则是s_t的分裂变量。依据上述理论，精度均减系数(MDA)作为测试重要性的指标，表达了模型引入一个变量后均方根误差的减少量。

D)将随机森林模型中的决策树数量设置为100，每个分叉的变量数为4。将上述波形和点云特征变量作为输入参数代入模型。模型自动剔除引入参数直到分类总体精度趋于稳定。然后再依据MDA指数选择重要性最高的前3个特征变量作为最优特征变量。

E)使用筛选后的最优特征变量和随机森林分类器再次进行分类，并借助总体分类精度(即对角线像元数/总体像元数×100％)和卡帕系数这两个定量化指标对分类精度进行评价。卡帕系数计算公式：

$K_{\mathrm{hat}}=\frac{N\underset{i=1}{\overset{r}{\mathrm{\Σ}}}{x}_{\mathrm{ii}}-\underset{i=1}{\overset{r}{\mathrm{\Σ}}}\left(x_{i+}{x}_{+i}\right)}{N^2-\underset{i=1}{\overset{r}{\mathrm{\Σ}}}\left(x_{i+}{x}_{+i}\right)}---\left(6\right)$

其中r为总的类别数，x_ii为对角线上的像元数，x_i+和x_+i是列和行的总像元，N是总像元数。

使用随机森林方法筛选最优特征变量并进行树种分类。将上述波形和点云特征变量作为输入参数代入模型。模型自动剔除引入参数直到分类总体精度趋于稳定(本次实验运行结果：剩余8个特征变量)。然后再依据MDA指数选择重要性最高的前3个特征变量作为最优特征变量(分别为：波形特征变量HOME的均值，WD的标准差以及点云特征变量h_mean)。使用随机森林分类器进行以单木为对象的树种分类。将以上3个特征变量作为输入参数代入随机森林分类模型(决策树数量：100；每个分叉设置4个变量)，并借助总体分类精度和卡帕系数等定量指标对分类精度进行评价。2个分类等级的混淆矩阵见表1-2。6个树种分类结果(表1)(总体精度＝66.89％；Kappa系数＝0.627)，4个树种分类结果优于6个树种(表2)(总体精度＝76.28％；Kappa系数＝0.684)。

表1 6个树种分类结果验证混淆矩阵

注：像元数已转换为百分比。

表2 4个主要树种分类结果验证混淆矩阵

注：像元数已转换为百分比。

表3波形特征变量描述

Claims (6)

1.一种LiDAR波形综合特征的单木识别方法，其特征在于，包括以下步骤：

1)借助机载小光斑全波形LiDAR传感器进行数据采集；传感器记录每束激光脉冲返回的完整波形信息；

2)LiDAR波形数据预处理：

A)噪声水平估计和数据平滑：首先把原始数据转换到频率域，再将频率较高的低值部分作为噪声水平的判断标准；然后选用高斯滤波器进行平滑。

B)高斯拟合分解及波形数据点云化：基于回波数据是多个高斯函数的累加这一假设，对波形数据采用非线性最小二乘法进行拟合；然后通过局部最大峰值检测滤波算法从处理后的波形数据上提取离散点云，每个离散点中记录了返回信号的能量和振幅信息；

C)生成数字地形：首先对从波形数据中提取出离散点云进行分类，然后对末次回波进行Kraus滤波处理用以去除非地面点，最后使用滤波后的末次回波数据并借助自然邻近法插值生成数字地形模型DTM；

3)单木定位和冠幅提取：

A)对地面以上点云进行中值滤波，然后将点云中的高度信息栅格化生成数字表面模型DSM；将DSM减去数字地形模型，得到归一化植被高度CHM；

B)通过局部最大值法确定单木树顶所在位置；

C)单木冠幅的确定：首先以树顶为中心拟合16个半径方向上的冠幅剖面，然后计算到局部最小值的水平距离，最后将这些距离值进行水平方向上的平均从而得到冠幅半径；

4)基于发射能量及传感器与地物的距离信息对LiDAR波形数据进行校正；

5)构建体元框架并进行LiDAR波形的结构化分解

首先利用三维体积单元划分地表以上的三维空间；然后在每个单元内，汇总穿过其中的最大脉冲能量并基于数字地形模型进行高度归一化；最后将每列体元柱信息汇总得到伪垂直波；

6)提取单木的波形特征变量

a)单木分割的树冠范围内，针对每组高斯分解后的每个波形，计算其能量ER、振幅WT和回波次数NT的均值μ及标准差σ，以此作为单木的第一组特征变量；

b)基于伪垂直波提取第二组特征变量：HOME、波形距离、高度中位数比、波峰数量、返回波形能量、冠层外层粗糙度、前坡度角和VDR；以上变量也分别在单木冠幅内计算均值μ和标准差σ；

7)提取单木的点云特征变量

在单木冠幅内提取了点云特征变量：高度百分位数变量组：h₁₀，h₂₅，h₅₀，h₇₅，h₉₀；冠层返回点云密度变量组：d₂，d₄，d₆，d₈；c.最小和平均高度h_min/h_mean；高度变异系数h_cv；覆盖度CC_2m/CC_mean；

8)使用随机森林方法筛选最优特征变量并进行树种分类

A)随机森林分类是由很多决策树分类模型{h(X,Θ_k),k＝1,2,...}组成的组合分类模型，且参数集{Θ_k}是独立同分布的随机向量，在给定自变量X下，每个决策树分类模型都由一票投票权来选择最优的分类结果；通过k轮训练，得到一个分类模型序列{h₁(X),h₂(X),…h_k(X)}，再用它们构成一个多分类模型系统，其最终的分类决策函数为：

$H\left(x\right)=\arg \underset{Y}{\max }\underset{i=1}{\overset{k}{\mathrm{\Σ}}}I(h_i\left(x\right)=Y)---\left(4\right)$

其中，H(x)表示组合分类模型，h_i是单个决策树分类模型，Y表示输出变量，I为示性函数；

B)通过随机森林方法生成重要性变量Xm，该变量通过累积每个t节点的权重减少量p(t)Δi(s_t,t)来评估预测量Y，由此平均树N_T结构写为：

$\mathrm{IMP}\left(X_m\right)=\frac{1}{N_T}\underset{T}{\mathrm{\Σ}}\underset{t\∈T:\mathrm{\υ}\left(s_t\right)=X_m}{\mathrm{\Σ}}p\left(t\right)\mathrm{\Δi}(s_t,t)---\left(5\right)$

式中，p(t)是当样本到达t时N_t/N的比值，v(s_t)是s_t的分裂变量；

C)将随机森林模型中的决策树数量设置为100，每个分叉的变量数为4；将上述波形和点云特征变量作为输入参数代入模型；模型自动剔除引入参数直到分类总体精度趋于稳定；然后再依据MDA指数选择重要性最高的前3个特征变量作为最优特征变量；

D)使用筛选后的最优特征变量和随机森林分类器再次进行分类，并借助总体分类精度和卡帕系数这两个定量化指标对分类精度进行评价；卡帕系数计算公式：

$K_{\mathrm{hat}}=\frac{N\underset{i=1}{\overset{r}{\mathrm{\Σ}}}{x}_{\mathrm{ii}}-\underset{i=1}{\overset{r}{\mathrm{\Σ}}}\left(x_{i+}{x}_{+i}\right)}{N^2-\underset{i=1}{\overset{r}{\mathrm{\Σ}}}\left(x_{i+}{x}_{+i}\right)}---\left(6\right)$

式中，r为总的类别数，x_ii为对角线上的像元数，x_i+和x_+i是列和行的总像元，N是总像元数。

2.根据权利要求1所述的LiDAR波形综合特征的单木识别方法，其特征在于：步骤1)中，所述的LiDAR传感器为奥地利的Riegl LMS-Q680i。

3.根据权利要求1所述的LiDAR波形综合特征的单木识别方法，其特征在于：步骤1)中，所使用的遥感平台为运-5固定翼飞机。

4.根据权利要求1所述的LiDAR波形综合特征的单木识别方法，其特征在于：步骤1)中，所述传感器的采样间隔为1ns。

5.根据权利要求1所述的LiDAR波形综合特征的单木识别方法，其特征在于：步骤3)的B)中，通过圆形的搜索框在CHM上“滑行”遍历栅格图层，通过逐个比较搜索框内的高度信息来确定最大值点；搜索框的半径通过高度信息和参数β₀及β₁来确定；β₀及β₁则通过地面实测数据拟合来获取，公式为：

CW(m)＝β₀+β₁×h²             (1)

式中，CW为地面实测冠幅半径，h为树高(m)，β₀和β₁为模型参数。

6.根据权利要求1所述的LiDAR波形综合特征的单木识别方法，其特征在于：步骤4)中校正的公式为：

$W_i^c=W_i/W^e---\left(2\right)$

$I_i^c=(I_i\·s_i^k)/(I^e\·s_0^k)---\left(3\right)$

上式中，为校正后返回波内每个高斯波的波宽，为校正后返回波内每个高斯波内每个高斯波的能量强度，W_i为原始返回波内每个高斯波的波宽，W^e为发射波的波宽，I_i为原始返回波的能量强度，I^e为发射波的能量强度，传感器到反射物体的距离，为标称距离，k为变化系数。