吴 昊, 张 晓 萌
(长江水利委员会水文局 长江中游水文水资源勘测局,湖北 武汉 430012)
激光扫描测量(Light Detection and Ranging,LiDAR)是一种非接触式高时空高分辨率高效率的三维对地观测技术[1],可获取物体表面密集点云,在全球变化监测、智慧城市[2]、资源调查、环境监测[3]、基础测绘及全球制图[4-7]等国家重大需求和地球系统科学研究中起到十分重要的作用。
刘世振等[8]对典型河段的不同植被覆盖类型进行精度测试,证实了机载激光技术在山区型河道地形测绘中具有较好的适用性。刘少聪等[9]结合金沙江下游梯级水电站水文泥沙观测项目,经精度分析与评定,机载三维激光扫描技术比测平面精度在0.42 m 以内,高程精度均在0.32 m以内。为避免后处理PPK激光点云解算的滞后性,进一步分析机载激光雷达技术在河道测绘中的适用性,本文采用实时RTK激光点云解算。
本试验区位于汉江仙桃以上至碾盘山河段,该段河道洲滩、潜洲较多,传统测量方法施展河道地形测量效率低且安全隐患多。为提高数据获取效率、减轻外业工作强度,本文将LiDAR技术应用于河道地形和断面测量,可弥补摄影测量航片、遥感影像在精度上的不足,且在空间尺度上不受限制。此外,对部分区域仍采用传统RTK测量,以检测LiDAR测量精度。
机载LiDAR系统是一种可昼夜测量的主动式全数字化设备,能够穿透茂密森林的树冠来采集高密度、高精度的LiDAR点云数据。原始数据主要包含激光脚点三维坐标信息、激光反射强度信息及激光被多次反射的回波次数[5]。LiDAR系统由机载GNSS系统、惯性导航系统、激光雷达扫描仪、中心控制单元、飞行器载体、多光谱相机及其他辅助设备构成(见图1)。下面分别对几个关键部分进行介绍。
LiDAR系统搭载在移动系统上,对空气和植被(空中激光)甚至水(测深激光)发出脉冲激光。搭载激光雷达的载体沿着运动方向移动,扫描仪在垂直运动方向上扫描。雷达扫描仪接收返回的激光脉冲,记录测量距离和角度。不同的扫描速度会影响LiDAR系统测得的点数和回波数。若LiDAR波束呈现圆锥扫描的同时飞机持续运动,则地面扫描线会呈现一系列重叠的椭圆形。
全球导航卫星系统(GNSS)提供有关传感器位置(纬度、经度、高度)的准确地理信息。惯性测量单元(IMU)通过加速度计和陀螺仪记录在此位置传感器的精确方位(俯仰角、横滚角、偏航角)。然后转换处理这两个设备记录的数据,形成静态数据。
机载平台搭载微机系统,处理各种设备获取及发送的信息。机载平台控制中心完成数据采集,实现各组成设备精确同步,同时记录存储采集的大量点云数据。除正常搭载的IMU和GNSS设备外,机载平台还可携带相机采集影像。GNSS基站与机载GNSS配合判断平台航行轨迹是否与预设航线重合。微机解压处理GNSS和IMU数据,计算外方位元素,根据激光雷达定位原理,生成点云并实时输出[5]。
原始数据包括每次脉冲的测距值、瞬时扫描角、机载GNSS观测值、千寻CORS站GNSS观测值、惯导IMU姿态值、同步航空影像。首先对机载和千寻CORS站GNSS数据进行同步差分处理,随后解算出每个激光回波脚点的三维坐标值,得到原始激光点云成果。由于数据采集过程中不可避免存在各类误差,因此需对原始点云进行相应系统误差检校和平差处理。检校平差后的激光点,需进行滤波分类[7]。点云滤波后,还可针对已分类的非地面点,精细分类为建筑物、植被等多种类型,也可实现地物的单体提取及建模。通过点云数据预处理、点云滤波、点云分类,可提供平面和高程精度较为准确的点云成果,供后续河道地形测量和断面测量成图使用。在直接提取数据成果的基础上,对点云进行栅格化内插,引入航空影像进行定向处理,可获得DEM和DOM成果。通过人工提取地物可以进一步得到DLG成果[8]。基于机载激光扫描系统的地形测量任务流程如图2所示。
试验区范围为仙桃市胡杨镇蔡滩村至仙桃市城区沿江大道中段轮渡码头。总面积约为8.7 km2,陆上面积约为6.4 km2。长度约为7.5 km,宽度约为1.2 km。任务为1∶5 000测图比例尺的全江地形测量与固定断面测量,测区地形较为丰富,人工建筑物包括堤防、居民地、道路、码头、电力线等。自然地貌包括林地、草地、农用地、河流等。试验区内分布4个汉江固定断面,编号分别为HX138、HX139、HX140、HX141,如图3所示[9]。
根据项目需求,主要投入的仪器设备如表1所列。
表1 主要仪器设备Tab.1 Main equipment
2.4.1机载激光雷达测量
该项目LiDAR施测时间为2021年7月28日、7月29日。测时天气晴朗、风力3~4级,RTK信号良好,利用大疆经纬M300 RTK无人机搭载禅思L1激光雷达扫描系统平台共飞行10个架次。航摄飞行参数见表2,航摄完成后,现场立即对数据进行预处理和质量检查,以确保数据可用。
表2 航摄飞行参数Tab.2 Aerial photography flight parameters
2.4.2激光雷达数据处理软件介绍
(1) 大疆智图。利用大疆智图进行点云数据处理,包含POS解算、点云与可见光数据融合、标准las点云格式输出、作业报告输出,实现点云数据处理一键式操作。选择点云密度为100%的点云进行处理,设置输出坐标系,选择点云精度优化,通过平差优化不同时刻扫描的点云,使得点云整体精度更高,尤其是能大幅降低相邻条带重叠部分的点云厚度,设置点云有效距离参数,滤除激光雷达采集的超过有效距离的点云。
(2) Terrasolid软件。LiDAR点云数据处理包括数据预处理和数据后处理两个关键步骤。采用由芬兰Terrasolid公司研发的Terrasolid系列软件进行数据处理,功能包括激光雷达点云传感器校准、点云数据校正、表面建模、激光雷达和正射影像生成等[10]。
利用Terrasolid软件的TerraScan模块剔除低位粗差和高位粗差。软件滤波算法是Axelsson于2000年提出的一种基于不规则三角网(TIN)[11]的渐进加密算法。该算法关键在于滤波阈值的选取,主要滤波参数为地形坡度角、迭代角、迭代距离,需根据经验选择合适的参考值,可有效保留地面断裂特征,处理曲面不连续的情况,适用于地物地形密集分布的城区[12]。
2.4.3点云数据生成与导出
利用大疆智图软件激光雷达点云处理模块,对点云数据进行预处理,生成las数据格式保存并导出。通过大疆智图直接浏览预处理结果,并可以通过RGB、反射率、高度、回波次数等多种方式进行数据展示。截取试验区一小块区域如图4所示。裸露地表激光反射强度要明显高于植被覆盖地区;当激光脉冲打到地面、道路等表面平滑地物表面时,激光一般为单次回波,激光脚点排列整齐,分布较为规律;当激光脉冲打到树木、建筑物等具有垂直结构的地物表面时,会发生多次回波现象,每个激光脉冲不再只有一个返回值,而是拥有多个返回值,激光脚点的分布较为散乱,不再规则排列,且由于具有高度的地物对激光的遮挡,会造成面向激光照射方向的点云与一部分地面点云在垂直方向上部分重合,而背向激光发射方向的区域由于地物的遮挡会出现空白区域。
2.4.4点云分类
机载LiDAR原始数据预处理后生成的海量点云包括地面点、植被点、建筑物点、水域点、其他地物点及噪声点,且均在同一数据层。点云分类最关键的一个步骤为地面点提取,主要通过反复迭代建立地面三角网模型分离地面点。算法包括两个阶段:第一,基于分离的低点在一定范围内搜寻初始地面点,建立初始模型并在这些点间构建临时的TIN模型;第二,基于原来初始模型,抬高模型,如果某些点和TIN之间的角度和距离符合预设参数,将增添这些新的点置于原来构建的TIN模型中,每个新增的点使得模型越来越接近真实地面[13]。地面点分类效果如图5所示。
2.4.5点云高程点提取
将RTK实测数据和分类后的点云数据叠加导入到EPS软件中,如图6所示。
EPS软件三维测图模块有手动提取最近点的功能,根据RTK实测点已知平面位置,手动提取与RTK实测点最近距离的点云数据,如图7所示。其中粉色点为RTK实测数据,红色点为选取的与RTK点最近的点云点。本次试验包含RTK实测点281个,与RTK实测点距离最近的281个点云数据。
2.4.6传统RTK测量
传统RTK测量采用千寻网络CORS,平面坐标系统为1954年北京坐标系,投影方式为高斯正形投影,高程系统采用1985国家高程基准。利用Trimble R10 RTK分别于树林、草地、裸露地表3种不同植被类型上,采集了281个检查点,并完成了4个汉江固定断面的施测,如表3所列。
表3 传统测量方式统计Tab.3 Statistics of traditional measurement methods
以RTK测量数据作为真值,统计LiDAR测量点云数据中裸露地表、草地、树林3种地形点的平面坐标位置误差ΔS,如图8~10所示;统计的平面误差最大值、最小值、平均值及中误差见表4,统计的平面误差区间分布情况见表5。
表4 地形点平面位置误差统计Tab.4 Point plane position error of topographic m
表5 地形点平面位置误差百分比区间统计Tab.5 Statistics of percentage interval of plane position error of topographical point
从图8~10和表4~5可知:裸露地表平面误差ΔS最大值为0.120 m,最小值为0.003 m,中误差为±0.049 m,其中97%的点位平面误差均在0.1 m以内;草地平面误差ΔS最大值为0.672 m,最小值为0.033 m,中误差为±0.271 m,其中73.7%的点位平面误差均在0.3 m以内;树林平面误差ΔS最大值为0.782 m,最小值为0.027 m,中误差为±0.392 m,其中63.8% 的点位平面误差均在0.4 m以内。所有检查点的平面中误差均满足1∶5 000成图比例尺下平原地区地形图平面中误差为±2.5 m的精度要求。
统计裸露地表、草地、树林3种地形点采用LiDAR测量和RTK测量的高程误差ΔH,如图11~13所示;统计的高程误差最大值、最小值、平均值及中误差见表6,统计的高程误差区间分布情况见表7。
表6 地形点高程误差统计Tab.6 Elevation error of topographic point m
表7 地形点高程误差百分比区间统计Tab.7 Statistics of percentage interval of elevation error of topographical point
从图11~13和表6~7可知:裸露地表高程误差ΔH最大值为0.210 m,最小值为-0.083 m,中误差为±0.075 m,其中99%的点位高程误差均在±0.2 m以内;草地高程误差ΔH最大值为0.326 m,最小值为-0.282 m,中误差为±0.158 m,其中75.8%的点位高程误差均在±0.2 m以内。其中草地高程误差ΔH大于0占比为93.9%,根据现场实际测量情况,试验区草地长势茂盛,在地表覆盖较厚区域,RTK测量时测量杆底部尖锐洞穿草地,故RTK测点为草地底部高程,即真实地面高程,而激光雷达穿透力有限,并不能够完全穿透草地,测量的大多数为草地表面顶部高程,所以本试验区中激光雷达测点高程一般比RTK测点高程高。树林高程误差ΔH最大值为0.329 m,最小值为-0.219 m,中误差为±0.103 m,其中94%的点位高程误差均在±0.2 m以内,而且所有检查点的高程中误差均满足1∶5 000成图比例尺下平原地区地形图高程中误差为±0.25 m的精度要求[14]。
利用试验河段机载三维激光点云数据进行断面截取,与RTK实测固定断面数据进行面积差比较。断面布设如图1所示,利用EPS点云模块在二维窗口上直接提取计划线上的点云数据,详见图14。
综合各地形地貌共选取4个断面进行比测,各断面面积差如表8所列,显示4个断面的面积差均不大于2%,满足SL 257-2017《水道观测规范》要求。
表8 断面面积差统计Tab.8 Statistics of section area difference %
(1) 经精度分析与评定,机载三维激光扫描测量与传统RTK测量比测,裸露地表、草地和树林平面中误差分别在±0.049,±0.271,±0.392 m以内,优于SL 257-2017《水道观测规范》点位平面位置允许中误差±2.5 m的要求;裸露地表、草地和树林高程中误差分别在±0.075,±0.158,±0.103 m以内,优于SL 257-2017《水道观测规范》陆上高程注记点允许中误差±0.25 m的要求;平面和高程均满足河道地形测量有关规范要求,满足长江中游流域平原地区1∶5 000成图比例尺地形图的相关技术要求。该技术可以应用于长江中游流域地形和断面数据获取。
(2) 对比平原地区机载三维激光扫描测量与传统RTK测量可知,裸露地表精度优于树林和草地精度,其中植被越密集的区域机载激光雷达测量误差越大。