基于时空重建的森林林木形态监测算法研究

梁 轩，王耀力，常 青，翟银枝

（1.太原理工大学信息与计算机学院，山西晋中 030600；2.山西省晋中市榆次区国有乌金山林场，山西 晋中 030619）

传统林木形态测量技术依靠人力测量，费时费力。通过无人机或地面车辆搭载相机[1-4]，短时间内收集大量图像，为实现自动化林木监测提供数据基础。同时定位与建图（SLAM）和运动推断结构(SFM)方法[5-7]适用于静态场景的重建，但无法解决动态变化的林木重建问题。构建时空地图[8-9]可以提高时空数据关联的准确性，这通常为自动驾驶应用设计，解决不同视角与不同时间的数据关联方法。董靖[10]实现了农作物的时空重建，但前端实时提取SIFT 特征将对数据收集前端造成较大计算负担。文献[11-12]通过计算农作物重建结果的几何特征，得出农作物生长趋势数据，但分析模型根据农作物设计，难以直接应用于树木重建结果分析中。

为此，提出了基于时空重建的森林林木形态监测算法，该算法由三部分组成：第一部分是搭建基于同步定位与建图的视觉数据收集前端，用于估计相机位姿；第二部分是后端建立图像间鲁棒数据关联，构建完整的时空重建；第三部分是设计通过稠密点云预处理，进行树木形态变化监测。

1 V-SLAM数据收集前端

1.1 RGB-D针孔相机模型

V-SLAM 数据收集前端的传感器是RGB-D 深度相机，RGB-D 深度相机基于针孔相机模型。该模型包括四个坐标系：世界坐标系、相机坐标系、成像平面坐标系与像素平面坐标系。假设现实世界空间点的三维坐标为(Xw,Yw,Zw)T，记为Pw。对应相机平面坐标系坐标为(Xc,Yc,Zc)T，记为Pc。对应成像平面坐标系的坐标为(Xi,Yi)T，记为Pi。对应像素平面坐标系的坐标为(u,v)T，记为P。设成像平面到相机光心O的距离（焦距）为f。

世界坐标系与相机坐标系转换如图1 所示。

图1 世界坐标系与相机坐标系转换

由模型推导出，世界坐标系的空间点转换到相机平面坐标系下的公式为：

R为3×3 矩阵表征相机旋转量，t表为3×1 矩阵表征相机平移量，都是系统待计算的数据，假设相机初始位置为空间原点时，旋转与平移矩阵实际上代表相机位姿矩阵T。

相机平面坐标系下空间点坐标转换到像素平面坐标系下齐次坐标的公式为：

K为相机内参矩阵，根据实际情况进行标定得到。图2、3 是相机平面到像素平面的转换过程示意图。dx与dy是成像平面到像素平面横轴与纵轴的放缩倍数，u0与v0是像素平面坐标原点在横轴与纵轴平移量。

图2 相机平面到成像平面转换

图3 成像平面到像素平面转换

1.2 ORB特征提取与PnP匹配

ORB 特征[13]在相机的旋转与图像的放缩情况下都有鲁棒性，且在相机自动增益、自动曝光与明暗变化的情况具有良好不变性。

提取ORB 特征比提取SIFT 特征[14]耗费更少时间，具体见后面的实验。数据收集前端运行过程中，系统实时性与鲁棒性制约相机移动速度，常用方法是借助GPU 加速并行处理，但在数据收集前端中使用GPU 过于奢侈。搭载相机的无人机或者手持嵌入式设备方式下使用GPU 也不便捷。提取ORB 特征可以适应普通CPU 的系统需求，降低前端设备成本，具有更好的实时性。

ORB 特征提取后，采用PnP 匹配模式估计相机位姿变化，通过随机抽样一致算法去除大部分匹配错误的特征。采用的前端相机可以直接输出深度图，由此求出三维空间特征点坐标，克服了单目相机采用PnP 匹配模式时需要初始化的缺点。采用非线性方法来计算PnP 匹配下的相机位姿估计，具体细节见1.3 节。

1.3 相机位姿估计与最小化重投影误差

在相机位姿估计中，采用PnP 匹配模式求解，把相机位姿估计看成一个非线性最小二乘问题。如图4 所示，求解最小二乘问题最优解，实际上就是最小化空间点重投影误差的过程。重投影误差是估计值与测量值的差值，估计值是重投影的像素坐标，测量值是相机测得的像素坐标P，由以下公式表示：

图4 重投影误差原理

将多空间点重投影误差最小二乘项相加求最小值，整体误差取得最小值时，求出的相机位姿就是最优解。经过RANSAC 排除错误匹配对，且优化前相机位姿与真实值比较接近，采用最小二乘法最小化重投影误差不容易陷入局部最优解中。

1.4 跟踪模块与局部建图模块

数据收集前端包括跟踪模块与局部建图模块两个核心模块。前端系统结构如图5 所示。

图5 前端系统结构

跟踪模块在相机当前帧中发现ORB 特征，完成PnP 匹配估计相机位姿。相机当前帧输入前需要预处理，因为采用的是基于机器人操作平台（ROS）的RGB-D 深度相机，RGB 模块与深度模块工作频率通常不同，采用ROS 的多传感器时间软同步机制实现两模块输入在相近时间内的非严格对准。跟踪模块对每帧图像做关键帧准入判定，关键帧主要保存了ORB 特征与相机估计位姿。

局部建图模块接收跟踪模块中的关键帧，建立局部稀疏地图。稀疏地图由能在多关键帧都能观测到的ORB 特征产生。通过维护全局地图容器，该容器包含了最基本的地图信息、关键帧信息及其连接信息，局部建图模块可以与之交互，将局部稀疏地图通过必要的整合操作构成全局稀疏地图，并对全局地图中的关键帧做准出判定，减少需要维护的关键帧。

2 时空重建后端

2.1 SFM实现稠密点云重建

SFM[5-7]与SLAM 类似，前者更关注于特征点估计以完整恢复重建结构，后者更关注相机位姿估计以获取更准确的定位信息。前端准确估计了相机位姿，后端将结合当前帧RGB-D 图像，完成稠密点云重建。

2.2 基于词袋的回环检测

SFM 是全局式重建，但前端收集数据却是增量式的，这导致了累积误差。采用基于词袋的回环检测[15]来降低累积误差。基于关键帧的词袋由每帧ORB 特征生成的描述向量经过k-means 聚类构成，通过树状数据结构查询匹配，最后由TF-IDF进行相似度计算确认回环。确认回环后，使用2.3节的逆向投影图形变化有界搜索方法，测出相机位姿作为先验，重新完成时空数据关联，解决从不同时间与不同视角观察同一场景时图像数据关联的关键问题。

2.3 逆向投影图形变化有界搜索

为建立三维重建之间的时空关联，采取逆向投影有界搜索算法在不同时间不同视角测得的关键帧图像之间建立数据关联。其原理公式如下，设a、b为两个时刻：

式（4）-（5）表示，对于空间同一点，已知深度因子、相机内参矩阵以及两时刻相机位姿矩阵时，两时刻对应像素坐标是确定的。为减少误差影响，在该确定像素坐标(ub,vb)上设定一个有界范围δ重新搜索该匹配像素点(u,v)，以建立鲁棒数据关联。

2.4 图优化解决全局光束法平差

后端会把输入的相机位姿与稀疏地图采用基于图优化理论[15]的G2o 库进行一次全局的光束法平差，其原理遵循1.3 节，最小化重投影误差仅用于局部稀疏地图建立，此处是全局稀疏地图优化。如图6 所示，虚线代表相机的运动轨迹，三角形代表单次位姿估计，表征相机运动模型。点虚线代表相机对ORB特征点对应空间点的观测（如Pw2可以在b、c、d时刻的相机位姿中观测到），观测模型对应1.3 节提到的估计值。

图6 图优化模型

2.5 后端系统结构与后期处理模块

如图7 所示，后端将前端输入数据进行优化后开始稠密重建，通过时空重建模块完成稠密点云时空数据关联，最后输出带有时间信息的稠密点云到点云后期处理模块。

图7 后端系统结构与后期处理模块

3 点云后期处理模块

3.1 基于半径的异常点滤波

时空重建后端输出的时空稠密点云带有许多重建噪声，通过基于半径的异常点滤波[17]去除大部分重建噪声。指定内点搜索半径与相邻最小点数，遍历稠密点云中所有点，点云某点搜索半径内点数少于相邻最小点数，则认为该点为噪声点予以排除。

3.2 基于RGB的重建噪声滤波

经过异常点滤波之后的点云仍存在部分重建噪声，来源是位姿估计误差下重建点的偏移，以远点代替了真正需要重建的点，远点往往重建出来带有特殊的RGB 信息，据此将可疑点滤除。

4 实验与结果分析

4.1 设备配置与数据集

实验测试计算机采用Intel CORE i5 10thCPU，Ubuntu 18.04系统。相机为Intel RealsensetmD455 深度相机。实验流程是使用ROS 发布订阅话题以获取相机输入，采用基于图优化的G2o 库，以及PCL 库做点云处理。采集室外实际植物数据以及TUM 数据集构建测试数据集。

4.2 ORB、SIFT特征提取对比实验

TUM 数据集中挑取任意同一帧图像分别提取SIFT 特征与ORB 特征。提取出3 431 个SIFT 特征，耗时0.57 s，平均每个特征耗时为166.13 μs。提取出500 个ORB 特征，用时9 978 μs，平均每个特征耗时为19.956 μs。提取速度提高接近一个数量级，这仅是提取单帧图像得到的提升，当图像规模高达几百帧时，将节省更多时间。ORB 特征在同一帧图像中数量上比SIFT 特征要少，这在特征匹配步骤中能节省时间。

4.3 后期处理点云实验

采用自建的室外树木数据集进行三维重建后进行基于半径的异常点滤波与基于RGB 的重建噪声滤波，如图8 所示，上图为滤波前，下图为滤波后，可以看到大部分环境中的噪声点都被滤除，保留下的点云更加紧凑与清晰。

图8 基于半径的异常点滤波与基于RGB的重建噪声滤波

4.4 树木形态监测结果

采用自建室外数据集进行时空重建结果变化检测，结果如图9 所示。低矮树植高度变化明显，三次不同时间的重建能在同一场景建立起数据关联，得到高度结果如表1 所示。当收集长期树木数据集时，也可以进行对应的形态监测，结果如图10 所示。

表1 低矮植株重建测量

图9 低矮树植数据关联重建结果（上：3月30日；中：4月13日；下：4月27日）

图10 树木数据关联重建结果（上：4月13日；下：4月27日）

5 结束语

该研究具有三方面优势：改进了数据收集前端，跟踪匹配ORB 特征而不是SIFT 特征，大大减少跟踪时的计算量；在后端实现了重建结果的准确时空数据关联；在重建基础上尽量还原树木信息，测量结果误差在5%以内。后续工作可以通过树木形态监测结果估算树木凋落可燃物量。