【翻译论文】Multispectral LiDAR Data for Land Cover Classification of Urban Areas（2017）

翻译论文

翻译论文

DOI: 10.3390/s17050958

文章目录

• 摘要：

  • 关键词：
  • 一、引言：
  • 二、多光谱激光雷达系统的演进历程：
  • 三、分类方法：

  • • 3.1基于图像的分类方法：

    • 3.2基于点的分类方法：

    • • 3.2.1多波长激光雷达点合并技术：
      • 3.2.2地面滤波技术：
      • 3.2.3NDVI计算公式：
      • 3.2.4数据聚类方法：

    • 3.3辐射矫正

• 四、研究区域和数据集

• 五、土地覆盖分类结果

• • 5.1 基于图像的分类结果分析

    • 5.2 基于点的分类结果探讨
    • 5.3 辐射校正对分类精度的影响分析

  • 六、 讨论

  • 结论

  • 致谢

  • 参考文献

摘要：

机载光探测和测距（LiDAR）系统主要基于单一波长运行，其核心功能是评估物体表面反射光谱的范围和强度。近期，多光谱激光雷达传感器的出现显著拓展了这一技术的应用范围，这些传感器能够采集不同波长的光谱数据。通过记录物体光谱反射率的多样性特征，本研究旨在探索多光谱LiDAR数据在土地覆盖分类中的有效性。研究采用了两种不同的分析方法：基于图像的分类技术和基于点的分类技术。在基于图像的方法中，LiDAR点云被转化为强度和高度图像，随后应用最大似然分类器进行识别。而在基于点的分类中，经过地面滤波处理后，计算归一化差异植被指数（NDVI）以辅助分类。研究数据集来自加拿大安大略省奥沙瓦的城市区域，分为四个类别：建筑、植被、道路和草地。实验结果表明，图像分类和三维点分类的准确率分别为89.9%和92.7%。此外，通过辐射校正模型处理强度数据，有效补偿了系统失真和地形高度变化导致的衰减效应。重复分类过程后，整体准确度未见显著提升。

关键词：

多光谱激光雷达；土地覆盖；地面过滤； NDVI；辐射校正

一、介绍

近年来，关于机载LiDAR技术在土地覆盖分类研究中应用高度和强度数据，已有大量研究进行了深入探讨 [1-4]。已有研究将激光雷达衍生高度表面与多光谱航空/卫星图像结合，并采用DSM或归一化数字表面模型(nDSM)格式。另有研究将多光谱航空/卫星图像与LiDAR高度和强度数据相结合 [7-9]。由于过去大多数研究将LiDAR强度或高度数据转换为二维图像，因此如强度[2-4,7,8]、多重返回[2,7]、DSM、数字地形模型(DTM)[2,4]和nDSM[5,6,8]等典型LiDAR图像被创建。在这些研究中，应用了基于像素的传统监督分类技术，如最大似然分类器[5,9]、基于规则的分类方法[2,3,8]和高斯混合模型[7]。而通过使用面向对象的分类技术[4,6]，则能够更好地考虑不同对象的空间相干性，从而有效避免基于像素分类结果中的噪声。

Brennan 和 Webster [2] 利用基于规则的分类方法，将由LiDAR数据生成的五个波段进行分割和归类，以确定土地覆盖类别。五个波段包括：地物表面模型（DSM）、数字地形模型（DTM）、反射强度、多次反射和归一化高度。分割过程分为四个层次，首先通过对图像像素应用平均强度、强度标准差、平均DSM、平均归一化高度和平均多重回报阈值进行处理，从而将图像分割成对象。第一层次将图像对象分为水体、低地和高地物体三类，随后将升高物体进一步细分为树木和建筑物，低地物体则被划分为低植被、道路和潮间带等子类。这种分类过程持续到第四层次，最终将图像对象划分为十类。经过验证，十类和七类分类的总体准确率分别为94%和98%。然而，该研究主要依赖于激光雷达数据的高度信息。此外，将分类规则中使用的阈值应用于对象的平均值（例如，平均强度）存在局限性，尤其是在建筑物边缘和地面返回在归一化高度范围内形成一个图像对象的情况下。此外，一些茂密的针叶树由于无法被激光束穿透而呈现出单回波，这些树被错误地归类为建筑物，因为其分离主要依赖于多重回报带。

在文献[3]中，系统性分析了三种不同地点的八种LiDAR衍生特征的敏感性。首先，基于LiDAR回波的五个波段进行了图像对象分割，包括裸土、第一回波、最后回波、高度和强度特征。其次，除了平均强度和紧密度外，还计算了基于高度的六个特征，如均值、标准差、同质性、对比度、熵和相关性。最后，通过决策树将图像对象分类为五个土地覆盖类别，分类准确率达到90%以上。值得注意的是，分类过程主要依赖于平均高度、高度标准差和平均强度特征，而其他特征的引入并未显著提升分类效果。这可能是因为测试区域的景观较为简单（如无密集的树木和建筑干扰）。此外，在分割过程中，强度波段的权重显著低于其他波段（仅为0.1），而其他波段的权重设为1。这是因为在数据采集过程中，人工调整了强度的记录权重。结果表明，沿不同航线的强度数据存在显著差异，这表明强度数据在本研究中未能得到充分应用。

Antonarakis等人。 [4] 首先采用基于监督学习的图像分析方法，将地形划分为九类。他们首次提取了8个LiDAR衍生波段，包括冠层表面模型、地形模型、植被高度模型（VHM）、强度模型、强度差异模型、偏度模型、峰度模型和冠层百分比模型。随后，基于三个市区数据集，应用决策树进行分类。该方法在三个数据集上实现了超过93%的整体准确率。然而，该方法未能考虑到道路和建筑物这两个基本类别，在其中一个调查区域中存在部分建筑物。VHM通过从树冠表面模型（第一个返回值生成）中减去数字地形模型（最后一次返回值生成）来计算。值得注意的是，激光雷达接收器中的噪声会干扰VHM的计算，导致一些最后的返回值比第一个脉冲返回具有更高的高度。此外，LiDAR点的三角不规则网络插值引入了图像构建的误差，且高海拔区域的观测数据未能准确反映水面状况。

其他研究者致力于利用 LiDAR 衍生的高程表面进行土地覆盖分类，其中多光谱图像中的nDSM（数字高程模型）被广泛应用。黄等人。 [5] 通过融合高分辨率 RGB 航拍图像和近红外波段图像，结合 LiDAR 衍生的 nDSM，采用基于像素的最大似然分类方法，实现了对建筑物、树木、道路和草地四个土地覆盖类别的识别。该方法的总体分类准确率达到 88.3%。进一步引入基于知识的分类和校正系统，分类准确率提升至 93.9%。该技术通过设定一组阈值，基于高度、高差、平滑度、各向异性平滑度、强度、NDVI、转换植被指数、面积和形状等因素进行检测，从而实现四种土地类别的分类。陈等人。 [6] 将 LiDAR 衍生的 nDSM 与 Quick-Bird 图像相结合，用于地形分类。首先，从 Quick-Bird 图像中提取归一化差异水指数和 NDVI 两个波段，与 nDSM 数据进行融合。其次，采用分层面向对象的分类方法，包括图像分割和对象阈值应用。分级分类过程在九个土地类别上的总体准确率达到 89.4%。然而，该方法未能有效区分道路与空地，因为两者在光谱和高程特征上具有相似性。值得注意的是，上述两项研究均未将 LiDAR 强度数据纳入分析框架。为了提高分类精度，研究者将光学图像重新采样至较粗分辨率，以匹配与 LiDAR 生成的混合像素一致的分辨率。这种处理方式导致部分像素同时具有多种土地覆盖特性，从而引入了分类错误。

其他研究使用带有LiDAR数据（高度和强度）的多光谱图像来利用多种陆地物体的光谱范围（如可见光和近红外(NIR)）的反射率变化。查拉尼亚等人。[7]基于航拍图像在可见光谱中测量的高度、强度、高度变化、多次回波和亮度带生成四个LiDAR波段。然后，利用高斯混合模型对四个类别（道路、草地、建筑物和树木）的训练数据进行建模。通过期望最大化算法估计模型参数和后验概率。随后，利用这些参数对测试数据集进行分类，分类的总体准确率达到85%。结果表明，高度变化在分类中起着关键作用，其中通过排除高度带获得了最差的分类结果。此外，通过排除航拍图像，整体精度有所下降。通过引入多重回报的方法，对道路和建筑物的分类性能进行了改进。然而，由于草地斑块的错误分类，导致其他地形覆盖物的分类精度下降。

哈特菲尔德等人。[8] 将1米分辨率的多光谱航空影像与LiDAR数据相结合。从LiDAR数据中提取了两个波段，即强度和nDSM，同时从多光谱航拍图像中获得了NDVI。分类和回归树对这些波段组合进行了性能测试。经过测试，LiDAR nDSM、多光谱图像和NDVI的组合在八个土地覆盖类别中表现出最佳的总体准确度（89.2%）。航拍图像中的阴影特征对分类结果具有显著影响。此外，由于LiDAR强度数据需要按照[8]中的方法进行校准，因此裸地和草本（草）类之间出现了误分类现象。辛格等人。[9] 将Landsat Thematic Mapper™图像与LiDAR衍生成的三个波段相结合，包括强度、冠层高度模型和nDSM。通过最大似然分类器，将土地覆盖划分为六类。考虑到不同分辨率的Landsat TM图像（1米、5米、10米、15米和30米分辨率），研究者进行了多组波段组合的测试。将1米分辨率的Landsat TM图像与三个LiDAR波段结合的分类组合取得了最高准确度（85%）。分类结果受到两个主要因素的影响：第一，在构建2D LiDAR图像时，LiDAR数据间的差距导致了理解上的偏差；第二，LiDAR强度数据未标准化到标准范围，这影响了分类效果。

近期研究表明，激光雷达强度数据的辐射校正对土地覆盖分类精度具有显著影响。辐射校正的目标是将记录的强度数据转换为陆地物体的光谱反射率。基于雷达距离方程，研究了飞行高度、距离、入射角、传感器孔径大小和大气衰减等系统和环境参数，以校正激光雷达强度数据[10,11]。采用经辐射校正的LiDAR强度数据后，城市区域的整体分类精度分别提升至7.4% [12]、9.4-12.8% [11]和3.8-16.5% [13]。

二、多光谱激光雷达系统的历史发展

在过去的几年里，针对多光谱LiDAR系统，人们进行了大量探索。实验室开发的多光谱LiDAR系统已成功收集了不同波长的数据 [14-16]。分析了从地面激光扫描(TLS)平台获取的多光谱LiDAR数据，以提取生物物理和/或生化植被参数 [17-21]。多光谱机载LiDAR系统已报道了多项探索，它整合了多种机载LiDAR系统，并结合同一研究区域的多架次飞行任务 [22-24]。

开发了基于实验室的多光谱LiDAR系统，获取了波长为531、550、660和780 nm [14]，以及556、670、700和780 nm [15]的数据，以评估森林冠层的三维结构。石等人。 [16] 研发了一种校准技术，用于校准基于实验室的多光谱LiDAR系统在556、670、700和780 nm波长下的反向散射强度测量。该方法综合考虑了入射角和表面粗糙度的影响。随后，研究者定义并探讨了多种植被指数，以提升分类的精确度。

其他调查采用了 TLS 平台以收集多光谱 LiDAR 数据。例如，[18] 开发了一个双波长全波形 TLS 平台，其工作波长为 NIR（1063 nm）和中红外（MIR，1545 nm）。该平台通过记录林冠返回的全波形数据，以测量其三维结构。芬兰大地测量研究所开发了一种高光谱 LiDAR (HSL) 系统，可传输 400-2500 nm 的连续光谱数据。在 500-980 nm 范围内设置了七个波段进行户外实验，以便根据其光谱响应区分人造目标和植被。随后，[21] 设计了一种便携式地面全波形 TLS 系统，用于收集内华达山脉国家森林的数据。该系统通过比较叶子在 1548 nm 处的吸收强度，成功将叶子与木质材料区分开来。

在多光谱机载LiDAR领域，Briese等学者[22]开发了一种实用的多波长机载LiDAR数据辐射校准工作流程。该方法基于全波形观测（包括范围、幅度和回波宽度）、飞行轨迹以及原位参考目标。本研究的数据集基于三个月内完成的三个飞行任务，这些任务均遵循相同的飞行计划。所采用的传感器包括三个RIEGL型号：VQ-820-G（532 nm）、VQ-580（1064 nm）和LMS-Q680i（1550 nm）。这些传感器在本研究中分别用于不同的任务。在辐射校准过程中，参考目标的原位测量采用了不同传感器的配合，这些测量是在零入射角的干燥条件下完成的。此外，本研究涵盖了从8月到12月不同季节的数据采集。由于各个飞行任务的地面条件存在差异，这导致校准后的强度值出现变化。例如，校准后的1064 nm波长反射率显示出比其他传感器更高的值。

布里斯等人。 [23] 对上述三个 RIEGL 传感器获取的多波长机载 LiDAR 数据进行了校准。 LiDAR 数据通过两次飞行任务（各配备两个传感器的飞机）在短时间内（仅四天）捕获，以确保研究地点在所有波长下的反射特性更为稳定。在 532 nm 处的校准强度数据非常暗，而 1064 nm 处的则更为明亮。此外，RIEGL VQ-820-G 相对于其他两个底拖传感器呈现出不同视角的LiDAR数据。一般来说，由于时间表面变化、大气条件和水分含量的影响，单个飞行任务的表面条件存在显著差异[23]。

然后，将LiDAR点转换为分辨率1米的光谱图像，并结合后续处理步骤进行综合分析。通过分析Optech和RIEGL传感器的数据，提取了三个关键特征，包括幅度（强度）、回波宽度和表面高度。随后，利用监督分类算法和支持向量机方法，将地形划分为六类，涵盖土壤、低植被、道路和砾石、高植被、建筑物屋顶和水体。通过测试不同特征组合，模型的总体准确率在84.3%至97.4%之间得到了验证。将3D LiDAR点转换为2D光谱图像会影响冠层下物体的冠层反射率信息，其中冠层无法与林下植被和土壤实现有效区分。本研究仅考虑了每个完整波形的第一个返回信号进行处理。然而，树木、建筑物屋顶或低矮植被等地物覆盖情况可能涉及多个返回信号。此外，在对比RIEGL和Optech的振幅信息时，建筑物屋顶未能与土壤或低矮植被完全分离。一个可能的原因是，强度数据来源于不同任务执行时间的不一致测量。由于天气条件和地表状态随时间变化，同一物体的强度值也会随之变化。因此，地表高度和回波宽度被确定为主导的特征，而幅度信息则作为补充数据[24]。

2014年，加拿大安大略省沃恩市的Teledyne Optech成功开发出世界上首款商用机载多光谱激光雷达传感器，该设备被命名为“Optech Titan”。该创新技术不仅支持白天和夜晚的多光谱主动数据采集，还推动了激光雷达在信息提取方面的广泛应用。该传感器具备三种工作波长，在三个独立通道中实现了不同视角的点云采集：其中C1通道采用1550 nm的红外线（MIR）进行3.5°前视扫描，C2通道则在0°最低点使用1064 nm的红外线进行扫描，而C3通道则采用532 nm的绿色光进行7°前方扫描。具体参数可在表1[25]中找到详细说明。

相比于仅使用单色波长 LiDAR 数据，采用多光谱 LiDAR 数据在三个不同波长上的数据融合能够显著提升数据的可靠性和准确性。已有研究表明，使用 Optech Titan 设备收集的多光谱 LiDAR 数据在土地覆盖分类方面表现出了良好的效果。研究者如威奇曼等人（参考文献[26]）对不同光谱模式进行了深入研究，发现强度值在土地覆盖分类中具有显著的区分能力。通过整合 LiDAR 强度和高度数据生成栅格图像，并结合先进的图像分类技术（参考文献[27]）实现精准的分类。

在作者之前开展的一项研究工作中，研究团队将最大似然分类器应用于单一强度图像、组合三强度图像以及组合三强度图像与DSM [28]的结合。实验结果表明，当采用组合三强度图像进行分类时，将地形划分为六类的总体准确率为65.5%。与单一强度图像相比，该方法的分类精度提升了17%。此外，结合DSM后，分类系统的整体准确率提升至72.5%。在方法验证过程中，研究团队从三个通道中提取了强度值，并计算了相应的光谱指数。在市区范围内进行验证测试，光谱指数能够将低覆盖植被与高覆盖植被与建成区域区分开来，总体准确率达到96% [29]。

与基于单独的多光谱图像或高分辨率航空/卫星LiDAR高程数据(DSM)[30]的研究成果相比，本研究中LiDAR高程和强度数据的结合使用显著改善了分类效果。尽管LiDAR系统能够获取高密度且精确的三维点云数据，但以往大多数研究将三维点云转换为二维强度和/或高程图像，以便应用图像分类技术。然而，这种转换过程会导致三维信息（即z分量）的丢失，从而可能导致不完整或不准确的分类结果。本研究旨在展示基于多光谱机载LiDAR数据的土地覆盖分类能力。具体而言，本研究的主要目标包括：（1）探索现有图像分类技术在多光谱LiDAR数据分类中的应用；（2）开发一种能够融合多波长激光雷达数据的方法；（3）开发一种基于三维多光谱LiDAR点的自动土地覆盖分类方法；（4）评估激光雷达数据辐射校正对土地覆盖分类结果的影响。本研究的结构安排如下：首先介绍本研究的方法论，包括基于图像和点的分类技术；随后在第3节详细阐述方法的理论基础；第4节介绍研究区域和使用的数据集；第5节展示土地覆盖分类结果；第6节对结果进行讨论和分析；最后在第7节总结主要研究发现和结论。

三、方法

多光谱LiDAR数据被用于将土地覆盖分类为四个类别：建筑物、树木、道路和草地。该系统开发了两种独立的分类方法，分别是基于图像的分类技术和基于点的分类技术。分类过程的流程图如图1所示。基于图像的分类技术是通过从LiDAR数据的高度和强度生成多个波段来实现的。具体而言，利用从三个不同的光谱波长收集的强度数据生成了三个强度图像。此外，通过利用高度数据生成数字高程模型（DSM）。将这三个强度图像与DSM相结合后，应用最大似然分类器进行分类。基于点的分类技术则直接作用于三维点云数据。首先，将来自LiDAR三个通道的点进行组合，并对每个LiDAR点分配相应的强度值。随后，通过地面滤波技术将非地面点与地面点区分开来。接着，基于三个强度值计算三个光谱指数，将非地面点分类为建筑物和树木，将地面点分类为道路和草地。这两种分类方法分别应用于原始LiDAR强度数据和经过辐射校正后的LiDAR强度数据。为了验证分类结果，研究者采用了与原始LiDAR数据同时获取的航拍图像进行对比测试。测试过程中，采用了总体准确度、生产者准确度、用户准确度以及Kappa统计量作为评估指标。

3.1基于图像的分类技术

基于LiDAR的点云数据处理流程首先生成栅格图像。通过不同波长的强度测量值，生成三个独立的强度图像。点的高程信息则用于生成高度图像（即数字高程模型，DSM）。为了确保每个像素包含足够的点，选择1米的像素大小。通过计算每个像素的平均强度值或高程，并将其分配给该像素位置。随后，采用3×3的移动平均窗口填充图像中的空缺区域[31]。将不同波长的强度图像叠加，即合并强度带，并在生成的强度带中包含DSM信息。根据研究区域的航拍图像，选择四个具有代表性的区域作为训练样本，提取每个类别的光谱特征。不同类别的样本具有不同的光谱特征，例如，建筑类由不同颜色屋顶的样本组成，如白色、灰色和红色屋顶。此外，植被类别包含不同植被覆盖程度的样本。最后，将两个波段的叠加图像应用最大似然分类器进行监督分类。该分类器考虑了每个像素/点属于特定类别的概率，并考虑了类别间的变异性。

3.2基于点的分类技术

基于点的分类技术分为四个阶段，这些阶段被直接应用于3D点云。首先，将不同通道的点云进行组合，估计每个点的三个强度值。其次，利用地面过滤技术，基于高程属性将非地面点与地面点区分开来。第三，根据不同通道记录的强度值，计算非地面点和地面点的NDVI值。第四，使用 Jenks自然间断优化方法定义阈值，随后将LiDAR点聚类到不同的类别中。有关此技术的更多细节将在后续章节中详细说明。

3.2.1 多波长激光雷达点合并

由于新的多光谱传感器获取不同波长的 LiDAR 数据，因此会针对相同的覆盖区域收集点云，但具有与不同波长相关的不同强度值。然而，合并这些点云并预测每个单点在所有波长的强度值，使可用数据更加密集和可靠。虽然 Optech Titan 同时在三个波长下工作，但它在三个通道中以不同角度获取 LiDAR 点。因此，不同通道中同一对象的采集点可能不会在同一位置完全重合。 3D 空间连接技术可以为合并所有通道的点提供一种可能的解决方案，其中一个通道的点的强度值分配给另一个通道的最近点 [26]。但是，这种技术可能会导致点之间的不正确匹配，如图 2 所示，并通过以下场景进行解释。情况 (1) 表示来自通道 C2 和 C3 的完美点匹配。在情况 (2) 中，来自 C2 的一个点可以与来自 C3 的两个不同点匹配两次，因为该点是两个点的最近邻居。情况 (3) 显示了来自 C2 的两个可能的相邻点，它们到 C3 的一个点的距离相同。情况 (4) 表示在预定义半径的球体内没有从 C2 到从 C3 的点的相邻点。因此，不能将每个点的强度值与 th 相同

在文献[29]中，该文提出了一种整合来自三个通道的LiDAR数据的新方法。研究者将来自每个通道的LiDAR数据划分为单元格尺寸为1 m的网格。每个单元格内所有点的平均强度值被分配到该单元格的中心位置。随后，研究者通过网格单元的平均强度值构建了三个光谱指数的网格。接下来，采用基于点位置的光谱指数网格的双线性插值技术，将三个光谱指数的值内插到每个LiDAR点中，其中相邻单元中心被用于计算插值过程。最后，具有三个光谱指数的LiDAR点被用于实现土地/水体识别和土地覆盖分类任务[29]。整个研究过程总结为图3。

在进行土地/水域区分时，这种方法是可以接受的，但当需要将地形分类为不同的土地覆盖类型时，就会导致错误分类。这种现象主要归因于两个方面。首先，利用单元格内各点的平均强度值。这些点可能属于或不属于同一个土地覆盖类型，因此，一个单元格可能代表多个土地覆盖类型。其次，通过双线性插值方法获取3D点的光谱值。因此，具有多个返回的点的光谱值会被错误地分配相同的值，例如来自树枝和树下裸土的点会被分配相同的光谱值。因此，为了正确预测一个点的强度值，需要通过从相邻通道的点计算中值来实现。本研究工作的点合并过程可以描述为

令pi、pj、ph分别代表C1、C2、C3中的点，其中i、j、h分别取值1, 2, 3, . . . , nC1、nC2、nC3，而nC1、nC2、nC3分别表示在C1、C2、C3中采集的LiDAR点的总数。每个通道中的LiDAR点首先通过K-d树数据结构进行组织，以便实现多维范围搜索的有效应用。分别对应C2和C3的pi的相邻点NC2pi和NC3pi在预设搜索半径®的球体内被获取，具体过程如下所述：

r值作为1米长，符合两个条件：第一条是有足够的点数量，第二条是不包含来自不同特征的点。NC2pi和NC3pi点在后，按照它们的强度值按升序排列。C2和C3中pi的强度值IC2pi和IC3pi分别计算得到。

中值强度值被用来消除任何强度数据噪声。当无法找到相邻点时，强度值将被赋值为零值。等式 (1)-(4) 用于 C2 中的任意点 pj，从而获得相邻点 NC1pj 和 NC3pj，以及分别来自 C1 和 C3 的强度值 IC1pj 和 IC3pj。同样的方法也适用于 C3 中的任何点 ph，其相邻点 NC1ph 和 NC2ph 以及强度值 IC1ph 和 IC2ph 则分别来自 C1 和 C2。LiDAR 点被组合，重复的点 (nd) 通过 MATLAB 函数“unique”被去除，而唯一的 xyz LiDAR 点则被检测并考虑用于分类过程；这样总点数（N）=nC1+nC2+nC3-nd，每个 LiDAR 点具有六个属性：x、y、z、IC1、IC2 和 IC3。

3.2.2.地面滤波

地面滤波器旨在通过图4所示的决策规则将地面点与其他非地面点区分开来。作为地面滤波器的第一步，将基于偏度的平衡统计分析应用于点的高程。自然测量的数据导致严格遵循正态分布模型[32]。因此，假设在LiDAR数据中收集的地面点严格遵循正态分布模型，而其他物体点（非地面点）可能会干扰这一分布[33,34]。通过从LiDAR数据中去除这些物体点，即可获得地面点。高阶矩特征，如偏度，能够有效地描述LiDAR点云的空间分布特征。偏度(Sk)定义为：

其中 N 代表 LiDAR 点的总数量，Zi 表示高程数据，其中 i ∈ {1, 2, . . . , N}。S 代表数据的标准差，µ 代表高程的算术平均值，根据等式 (6) 和 (7) 进行计算。

首先，将点云的高程按升序排列。随后，通过等式（5）计算所有点的偏度。如果计算得到的偏度大于零，则移除具有最高高程的点并归类为非地面点。剩余的点用于计算偏度，并重复上述步骤，直至点云的偏度平衡（Sk = 0）。完成偏斜平衡后，剩余的点被分类为潜在的地面点，并假设其在特定的斜率范围内。输出分离是基于每个LiDAR点与其相邻点的斜率变化进行细化的。为了进一步区分非地面点，应用阈值（S_thrd）将具有显著高斜率的点标记为非地面点。最后，将剩余的地面点划分为网格结构，以去除高程异常值。对于每个网格区域，计算其最小高程（Z_min），并应用高程阈值（E_thrd）进行筛选。

3.2.3 NDVI计算

NDVI 值的计算方法与 [35] 中定义的强度数据类似，如下所示：

其中，MIR、NIR和G分别代表MIR、NIR和绿色波段的强度值。NDVI值落在区间-1至1之间。然而，当一个点在两个通道中的强度值均为零时，NDVI将无法计算。此时，该点将被标记为未分类点。

3.2.4 数据聚类

Jenks自然间断优化方法用于确定阈值范围（NDVI_thrd），以便根据NDVI值对LiDAR点进行聚类[36]。该优化方法旨在最小化类内方差并最大化类间方差。给定NDVI值范围为[a,…,b]，其中-1 ≤ a < b ≤ 1，且阈值（NDVI_thrd）位于该范围内。阈值（NDVI_thrd）被识别为将地面点归为道路和草地，而非地面点归为建筑物和树木。通过最大化类间平方差，其中M为NDVI值的平均值，M1和M2分别为第一类和第二类的平均值。首先计算M。然后，将点分为范围为[a,…,NDVI_thrd]和[NDVI_thrd,…,b]。计算平均值M1和M2。最后，由式（11）得到最优阈值（NDVI_thrd）。

3.3辐射矫正

其中，Pt代表发射的激光能量，Dr为孔径直径，R为测量范围，βt是激光束宽度，ηsys是系统效率，ηatm是大气衰减系数。激光横截面积σ由投影目标区域A、扫描角度θ和被照射表面的光谱反射率ρ共同决定。在本研究中，假设其他参数不变，采用基于雷达距离方程的辐射校正模型，通过保持R、A和θ不变，将Pr值转换为ρ值，以消除系统相关失真。入射角θ的定义基于两种情况：一种是入射激光束与物体表面法线之间的角度 [11]；另一种是考虑表面斜率所决定的入射角与扫描角的组合 [11]。有关辐射校正模型的详细信息，读者可参考 [11,13]。

四、研究区域和数据集

研究区位于加拿大安大略省奥沙瓦市。于2014年9月3日的飞行任务中，通过Optech Titan多光谱LiDAR传感器采集单个条带的LiDAR数据点。该设备在1075米高程、±20度扫描角度、200千赫兹/通道的脉冲重复频率下，通过三个通道采集了LiDAR点频率和40赫兹的扫描频率。每个通道的平均点密度为3.6每平方米，点间距约为0.5米。所获取的数据包括轨迹位置信息，以及每个通道的LAS文件格式中带有多个返回（最多4个返回）的三维点云。这些LAS数据文件包含每个LiDAR点的xyz坐标、原始强度值、扫描角度和GPS时间信息。

LiDAR 条带中的一个子区域被裁剪为550米乘以380米的区域，用于测试目的。研究区域涵盖了地面上多种类型的地物特征，包括建筑设施、道路、停车场、灌木丛、树木以及草地覆盖的空地。在测试过程中，来自C1、C2和C3三个通道的点数分别为712,171、763,507和595,387个。不同通道记录的点数变化反映了不同波长下陆地物体相互作用的特性，例如植被的绿色度。同时获取的航拍图像与LiDAR数据实现了地理对齐，用于验证土地覆盖分类的准确性，如图5所示。

因为3D参考点无法获取，所以采用了多边形来提取各类别（如建筑、树木、道路和草地）的参考点。在各类中绘制的多边形内部的所有点，包括道路、草地以及建筑物，都被标记为同一类别，而树类的顶层被用作参考点。各类的总得分为45618分，具体分布情况可见表2。

五、土地覆盖分类结果

研究区域划分为四种不同的土地覆盖类型，包括建筑物、树木（包括树木和灌木）、道路（沥青面、停车场和裸土）以及草地（包括绿/干草和湿地）。通过在预设多边形内部设置多个采样点来进行精度评估。首先，围绕研究区域的中心位置进行数字化处理，以减少分类混淆的可能性。在采用图像分类技术时，可能会出现由混合像素导致的分类模糊。随后，研究区域内的所有采样点均通过地理参考航拍图像进行标注。接着，构建混淆矩阵并计算包括总体准确率、生产者准确率和用户准确率在内的准确性指标，同时计算Kappa统计量来评估分类性能。

5.1 基于图像的分类结果

LiDAR点用于根据强度值（即分别来自通道1、2和3的C1、C2和C3）创建三个光栅图像，以及来自高度数据（即DSM）的光栅图像，空间分辨率为1 m，如图6所示。将三个强度图像叠加，形成组合强度带（CIB）和基于DSM的光栅图像（CIBs_DSM）。图7a展示了CIB的伪彩色合成结果，其中C1以红色显示，C2以绿色，C3以蓝色；图7b展示了CIBs_DSM的伪彩色合成结果，其中DSM以红色显示，C3以绿色，C2以蓝色。在识别不同类别的训练特征后，将最大似然分类器应用于CIB和CIBs_DSM这两个波段组合。图7c、d展示了来自未带DSM和带DSM的组合强度带的分类结果。表3和表4分别展示了仅使用CIB和使用CIBs_DSM时的混淆矩阵、总体准确度和Kappa统计量。其中，CIB的总体准确度和Kappa值分别为77.3%和0.675，而CIBs_DSM的总体准确度和Kappa值分别为89.9%和0.855。

5.2 基于点的分类结果

该技术主要包含地面滤波和NDVI计算两个关键步骤。首先，这些计算分别用于高程和强度属性的处理。随后，处理结果被划分为不同的类别。地面滤波过程从倾斜平衡机制入手，旨在将地面点与非地面点区分开来。在此基础上，通过基于坡度的变化机制对地面点进行细化处理。具体而言，研究了每个点与其周围点的斜率特征；如果某点的斜率超过阈值S_thrd（设定为10°），则该点将被归类为非地面点。值得注意的是，部分海拔较高的点并未被归类为非地面点。因此，处理后的输出结果被划分为单元格大小为25m的网格结构。对于每个网格单元，如果其高程值高于最小高程值3m（设定为E_thrd=3m），则该单元将被归类为非地面点。地面点主要包括道路和草类，而非地面点则包括建筑物和树木类。随后，根据公式（8）至（10）计算NDVI值。表5展示了通过Jenks_break优化方法对非地面和地面点分别进行NDVI计算，以实现将建筑物与树木和道路与草类分开的效果。植被特征（即树木或草类）在NIR、MIR和绿色波长处表现出较强的反射特性。因此，在计算NDVI值时，植被点的NDVI值将显著高于建成区的NDVI值。具体而言，对于某个特定点，如果其满足NDVINIR-MIR、NDVINIR-G或NDVIMIR-G均低于NDVI_thrd，则该点的覆盖类型应归类为建筑物或道路类；若至少有一个NDVI指标超过NDVI_thrd，则应归类为植被类。图8展示了基于三个NDVI指标的三维分类点云图。表6至表8则提供了不同分类情况下的混淆矩阵、总体准确度和Kappa统计量。

在两个主要来源的影响下，该区域的分类结果存在明显偏差。在强度预测期间，对于任何给定的点，若在其搜索半径内未找到相邻点，则将该点的强度值设为零。该点的NDVI值将被标记为未分类点，而非数字。表6-8中，错误分类情况以灰色标注，其中道路和草地被错误分类为建筑物或树木，反之亦然。

5.3 辐射校正对分类精度的影响

LiDAR强度数据经过系统范围和扫描角度的衰减校正处理。通过等式（12）和（13）计算每个LiDAR点的强度ρ值。随后，基于三个通道构建校正后的强度图像，并重复应用基于图像的分类技术进行分类。当仅使用LiDAR强度数据进行辐射校正时，整体精度提升幅度约为1.7%；而当引入数字高程模型（DSM）时，精度提升幅度显著降低，仅为0.6%。值得注意的是，建筑类物体的分类准确率显著提升，约为10.8%；而树木类物体的分类准确率则提高了约3.1%。在对辐射校正后的LiDAR点进行相同分类流程处理时，未观察到整体准确性的显著提升。这可能与LiDAR数据仅在有限的扫描角度下获取有关。

六、 讨论

总的来看，分类结果证实了利用多光谱激光雷达数据将地形划分为四类的能力，即建筑、树木、道路和草地。两种分类技术分别达到了89.9%和92.7%的总体分类准确率。此前的研究在相同的四个土地覆盖类别中实现了85-89.5%的总体分类准确度。他们在研究中结合了多光谱航空/卫星图像和从LiDAR数据[5,6]以及基于LiDAR高度和强度数据[7-9]得出的nDSM，而本研究仅依赖于LiDAR数据。多光谱LiDAR数据的可用性消除了对多光谱航空/卫星图像进行分类的需要。图9总结了两种分类技术所达到的总体分类准确率。

研究表明，融合LiDAR强度和高度数据具有重要意义。实验证明，在将DSM与强度图像融合后，分类精度显著提升，较之前提升了12%以上。这一现象在图7c中可明显观察到，研究区域中部分像素（以红色矩形标注）被误分类为树木（即高植被）。此外，许多原本属于路面的像素被误分类为建筑物，而这些像素应被正确识别为草（即低植被）或道路。在DSM融合后，如图7d所示（以黑色矩形标注），分类精度得到了显著提升。进一步分析表明，当考虑DSM时，生产者和用户的分类准确度分别得到了显著提升。具体而言，生产者对建筑物、树木和道路类别的准确度分别提升了29.7%、9.1%和13.3%。此外，用户对树木和道路分类的准确率分别提升了15.4%和28.8%。

点分类技术在土地覆盖分类中产生了多种整体精度。采用 NDVINIR-MIR、NDVINIR-G 和 NDVIMIR-G 时，其总体准确度分别为 77.8%、92.7% 和 88.0%。如前所述，未分类点和地面过滤是导致分类错误的两个主要原因。未分类点的误差范围为 0.0至2.3%，而地面过滤的误差范围为 0.0至4.4%。

在对单一类别进行关注时，约35.6%的树点被省略（生产者准确率为64.4%），并且在使用NDVINIR-MIR时这些点被错误地归类为建筑物。这种遗漏导致建筑物类别的错误分类率约为37%（用户准确率为63%）。这种现象主要由植被的水分含量引起，其中干燥的植被在C1 [38]中表现出较高的强度值。结果，NDVINIR-MIR产生的值较低，导致树点被错误地归类为建筑物。同样地，约19.7%的草地点被省略（生产者准确率为80.3%），主要归类为道路。此外，分别使用NDVINIR-MIR和NDVIMIR-G时，约31.1%和14.5%的道路点被错误分类为草。

此外，使用 NDVIMIR-G 时，约 20.7%（生产者准确度达 79.3%）的建筑物点被错误分类为树木。这些屋顶材料在 C1 中表现出极高的强度值，导致 NDVIMIR-G 被引入，从而引发分类错误。结果，这种遗漏导致树类被错误分类的比例约为 14.4%，用户准确率降至 85.6%。

尽管先前的研究显示，辐射校正和归一化处理可能在一定程度上改善了强度均匀性[11-13]，但在本研究中并未观察到这一现象。这种差异可能源于所使用的多光谱LiDAR数据在扫描角度方面存在显著差异：本研究采用的是窄扫描角度，而其他研究则采用了宽扫描角度。这种扫描方式的差异导致了能量的明显损失，尤其是在接近扫描线边缘的区域。此外，与目标相关的参数，包括目标到传感器的距离、目标尺寸、激光入射角以及目标表面的照度，对LiDAR强度数据的影响可能与以往研究不完全相同。另一个重要因素是研究区域内与数据采集相关的环境参数，如气溶胶、瑞利散射、气溶胶吸收以及大气衰减，这些参数保持不变，从而对LiDAR强度数据的影响与先前的研究保持一致。最后，Optech Titan传感器的三个通道通过其内部传递函数精确调节，输出的信号强度呈现良好的线性特性。值得注意的是，由于信号传输系统的完善，这种情况下并未出现能量损失。因此，在类似的实验条件下，可能不需要总是进行辐射校正处理。

结论

本研究探讨了多光谱LiDAR数据在城市地区土地覆盖分类问题中的应用。多光谱数据由Optech Titan传感器捕获，该传感器在1550 nm、1064 nm和532 nm三个波长范围内运行。通过两种分类方法将多光谱LiDAR数据划分为建筑、树木、道路和草地四种类别。第一种方法是基于图像的分类，其中将LiDAR强度和高度数据转换为图像数据。构建了两个波段组合：一个是三强度图像，另一个是三强度图像与DSM的结合。对该方法进行评估时发现，仅依赖LiDAR强度数据，实现了约77.32%的总体分类准确率。当结合DSM数据时，分类准确率提升至89.89%。

第二种方法是基于点分类的，通过将三个通道中的3D LiDAR点进行整合，并将每个LiDAR点的三个强度值作为预处理步骤进行计算。随后，通过偏斜平衡和斜率变化的分析，对LiDAR数据进行地面滤波处理，将数据划分为地面点和非地面点。接着，计算了NDVI值，并运用 Jenks break优化方法确定阈值，将地面点分类为道路和草地，非地面点分类为建筑物和树木。采用NDVINIR-MIR、NDVINIR-G和NDVIMIR-G方法时，该技术分别获得了77.83%、92.70%和88.02%的总体准确率。此外，基于雷达距离方程的物理模型对强度数据进行了辐射校正处理，校正过程中考虑了系统参数和地形因素的影响。值得注意的是，在对特定数据集进行辐射校正后，对土地覆盖分类结果的影响并不显著。整个研究工作凸显了使用多维LiDAR数据（包含强度和高度信息）相较于单维LiDAR数据（仅包含强度或高度信息）的优势。

致谢

这项研究工作获得了加拿大自然科学与工程研究委员会 (NSERC) (RGPIN-2015-03960) 和安大略延龄草奖学金 (OTS) 的专项拨款支持。作者还需感谢 Paul LaRocque 和 Teledyne Optech 提供来自世界上第一个多光谱 LiDAR 系统 Optech Titan 的 LiDAR 数据。

参考文献

1.Yan, W.Y .; Shaker, A.; El-Ashmawy , N. Urban land cover classification using airborne LiDAR data: A review.Remote Sens. Environ. 2015, 158, 295–310. [CrossRef]
2. Brennan, R.; W ebster, T.L. Object-oriented land cover classification of lidar-derived surfaces. Can. J. Remote Sens.2006, 32, 162–172. [CrossRef]
3. Im, J.; Jensen, J.R.; Hodgson, M.E. Object-based land cover classification using high-posting-density LiDARdata. GISci. Remote Sens. 2008, 45, 209–228. [CrossRef]
4. Antonarakis, A.S.; Richards, K.S.; Brasington, J. Object-based land cover classification using airborne LiDAR.Remote Sens. Environ. 2008, 112, 2988–2998. [CrossRef]
5. Huang, M.-J.; Shyue, S.-W.; Lee, L.-H.; Kao, C.-C. A knowledge-based approach to urban feature classificationusing aerial imagery with Lidar data. Photogramm. Eng. Remote Sens. 2008, 74, 1473–1485. [CrossRef]
6. Chen, Y .; Su, W.; Li, J.; Sun, Z. Hierarchical object oriented classification using very high resolution imagery and LIDAR data over urban areas. Adv. Space Res. 2009, 43, 1101–1110. [CrossRef]
7. Charaniya, A.P .; Manduchi, R.; Lodha, S.K. Supervised parametric classification of aerial LiDAR data.In Proceedings of the IEEE Computer Society Conference on Computer Vision and Pattern Recognition Workshops, Washington, DC, USA, 27 June–2 July 2004.
8. Hartfield, K.A.; Landau, K.I.; van Leeuwen, W.J.D. Fusion of high resolution aerial Multispectral and LiDAR data: Land cover in the context of urban mosquito habitat. Remote Sens. 2011, 3, 2364–2383. [CrossRef]
9. Singh, K.K.; V ogler, J.B.; Shoemaker, D.A.; Meentemeyer, R.K. LiDAR-Landsat data fusion for large-area assessment of urban land cover: Balancing spatial resolution, data volume and mapping accuracy . ISPRS J.Photogramm. Remote Sens. 2012, 74, 110–121. [CrossRef]
10. Höfle, B.; Pfeifer, N. Correction of laser scanning intensity data: Data and model-driven approaches. ISPRS J.Photogramm. Remote Sens. 2007, 62, 415–433. [CrossRef]
11. Yan, W.Y .; Shaker, A.; Habib, A.; Kersting, A.P . Improving classification accuracy of airborne LiDAR intensity data by geometric calibration and radiometric correction. ISPRS J. Photogramm. Remote Sens. 2012, 67, 35–44.[CrossRef]
12. Habib, A.F.; Kersting, A.P .; Shaker, A.; Yan, W.-Y . Geometric calibration and radiometric correction of LiDAR data and their impact on the quality of derived products. Sensors 2011, 11, 9069–9097. [CrossRef] [PubMed]
13. Yan, W.Y .; Shaker, A. Radiometric correction and normalization of airborne LiDAR intensity data for improving land-cover classification. IEEE T rans. Geosci. Remote Sens. 2014, 52, 7658–7673.
14. Woodhouse, I.H.; Nichol, C.; Sinclair, P .; Jack, J.; Morsdorf, F.; Malthus, T.J.; Patenaude, G. A Multispectral
Canopy LiDAR Demonstrator Project. IEEE Geosci. Remote Sens. 2011, 8, 839–843. [CrossRef]
15 Wei, G.; Shalei, S.; Bo, Z.; Shuo, S.; Faquan, L.; Xuewu, C. Multi-wavelength canopy LiDAR for remote sensing of vegetation: Design and system performance. ISPRS J. Photogramm. Remote Sens. 2012, 69, 1–9.[CrossRef]
16.Shi, S.; Song, S.; Gong, W.; Du, L.; Zhu, B.; Huang, X. Improving Backscatter Intensity Calibration for Multispectral LiDAR. IEEE Geosci. Remote Sens. 2015, 12, 1421–1425. [CrossRef]
17.Wallace, A.; Nichol, C.; Woodhouse, I. Recovery of Forest Canopy Parameters by Inversion of Multispectral LiDAR Data. Remote Sens. 2012, 4, 509–531. [CrossRef]
18.Danson, F.M.; Gaulton, R.; Armitage, R.P .; Disney , M.; Gunawan, O.; Lewis, P .; Pearson, G.; Ramirez, A.F.Developing a dual-wavelength full-waveform terrestrial laser scanner to characterize forest canopy structure.Agric. For. Meteorol. 2014, 198–199, 7–14. [CrossRef]
19.Hakala, T.; Suomalainen, J.; Kaasalainen, S.; Chen, Y . Full waveform hyperspectral LiDAR for terrestrial laser scanning. Opt. Express 2012, 20, 7119. [CrossRef] [PubMed]
20.Puttonen, E.; Hakala, T.; Nevalainen, O.; Kaasalainen, S.; Krooks, A.; Karjalainen, M.; Anttila, K. Artificial target detection with a hyperspectral LiDAR over 26-h measurement. Opt. Eng. 2015, 54, 013105. [CrossRef]
21.Douglas, E.S.; Martel, J.; Li, Z.; Howe, G.; Hewawasam, K.; Marshall, R.A.; Schaaf, C.L.; Cook, T.A.;Newnham, G.J.; Strahler, A.; et al. Finding leaves in the forest: The dual-wavelength Echidna Lidar.IEEE Geosci. Remote Sens. 2015, 12, 776–780. [CrossRef]
22.Briese, C.; Pfennigbauer, M.; Lehner, H.; Ullrich, A.; Wagner, W.; Pfeifer, N. Radiometric calibration of multi-wavelength airborne laser scanning data. In Proceedings of the XXII ISPRS Congress, Melbourne,Australia, 25 August–1 September 2012.
23.Briese, C.; Pfennigbauer, M.; Ullrich, A.; Doneus, M. Multi-wavelength airborne laser scanning for archaeological prospection. In Proceedings of the XXIV International CIPA Symposium, Strasbourg, France,2–6 September 2013.
24.Wang, C.-K.; Tseng, Y .-H.; Chu, H.-J. Airborne dual-wavelength LiDAR data for classifying land cover. Remote Sens. 2014, 6, 700–715. [CrossRef]
25.Teledyne Optech—Titan Brochure and Specifications, 2015. Optech Titan Multispectral LiDAR System—High Precision Environmental Mapping. Available online: http://www.teledyneoptech.com/wp-content/uploads/Titan-Specsheet-150515-WEB.pdf (accessed on 30 August 2016).
26.Wichmann, V .; Bremer, M.; Lindenberger, J.; Rutzinger, M.; Georges, C.; Petrini-Monteferri, F. Evaluating the potential of multispectral airborne LiDAR for topographic mapping and land cover classification.In Proceedings of the ISPRS Geospatial Week, La Grande Motte, France, 28 September–3 October 2015.
27.Bakuła, K.; Kupidura, P .; Jełowicki, Ł. Testing of land cover classification from multispectral airborne laser scanning data. In Proceedings of the XXIII ISPRS Congress, Prague, Czech Republic, 12–19 July 2016.
28.Morsy , S.; Shaker, A.; El-Rabbany , A. Potential use of multispectral airborne LiDAR data in land cover classification. In Proceedings of the 37th Asian conference on Remote Sensing, Colombo, Sri Lanka,17–21 October 2016.
29.Morsy , S.; Shaker, A.; El-Rabbany , A.; LaRocque, P .E. Airborne multispectral LiDAR data for land-cover classification and land/water mapping using different spectral indexes. In Proceedings of the XXIII ISPRS Congress, Prague, Czech Republic, 12–19 July 2016.
30.Zhou, W. An object-based approach for urban land cover classification: Integrating LiDAR height and intensity data. IEEE Geosci. Remote Sens. 2013, 10, 928–931. [CrossRef]
31.Reuter, H.I.; Nelson, A.; Jarvis, A. An evaluation of void-filling interpolation methods for SRTM data. Int. J.Geogr. Inf. Sci. 2007, 21, 983–1008. [CrossRef]
32.Duda, R.O.; Hart, P .E.; Stork, D.G. Pattern Classification, 2nd ed.; Wiley: New York, NY , USA, 2001; p. 55.
33.Bartels, M.; Wei, H.; Mason, D.C. DTM generation from LIDAR data using skewness balancing.In Proceedings of the 18th IEEE International Conference on Pattern Recognition (ICPR), Hong Kong,
China, 20–24 August 2006.
34.Yan, W.Y .; Morsy , S.; Shaker, A.; Tulloch, M. Automatic extraction of highway light poles and towers from
mobile LiDAR data. Opt. Laser T echnol. 2016, 77, 162–168. [CrossRef]
35.Rouse, J., Jr.; Haas, R.H.; Schell, J.A.; Deering, D.W. Monitoring vegetation systems in the Great Plains with ERTS. In Third Earth Resources T echnology Satellite-1 Symposium NASA SP-351; Freden, S.C., Mercanti, E.P .,Becker, M.A., Eds.; NASA: Washington, DC, USA, 1974; V olume I, p. 309.
36.Jenks, G.F. The data model concept in statistical mapping. In International Y earbook of Cartography;Frenzel, K., Ed.; Kirschbaum V erlag: Bonn, Germany , 1967; V olume 7, pp. 186–190.
37.Jelalian, A.V . Laser radar systems. In Proceedings of the Electronics and Aerospace Systems Conference,Arlington, V A, USA, 29 September–1 October 1980.
38.Gao, B.C. NDWI—A normalized difference water index for remote sensing of vegetation liquid water fromspace. Remote Sens. Environ. 1996, 58, 257–266. [CrossRef]

Morsy S, Shaker A, El-Rabbany A. multi-spectral LiDAR data in the context of urban land cover classification[J]. Sensors, 2017, 17(5): 958.