梁 軒，王耀力，常 青，翟銀枝

（1.太原理工大學(xué)信息與計算機學(xué)院，山西晉中 030600；2.山西省晉中市榆次區(qū)國有烏金山林場，山西 晉中 030619）

傳統(tǒng)林木形態(tài)測量技術(shù)依靠人力測量，費時費力。通過無人機或地面車輛搭載相機[1-4]，短時間內(nèi)收集大量圖像，為實現(xiàn)自動化林木監(jiān)測提供數(shù)據(jù)基礎(chǔ)。同時定位與建圖（SLAM）和運動推斷結(jié)構(gòu)(SFM)方法[5-7]適用于靜態(tài)場景的重建，但無法解決動態(tài)變化的林木重建問題。構(gòu)建時空地圖[8-9]可以提高時空數(shù)據(jù)關(guān)聯(lián)的準確性，這通常為自動駕駛應(yīng)用設(shè)計，解決不同視角與不同時間的數(shù)據(jù)關(guān)聯(lián)方法。董靖[10]實現(xiàn)了農(nóng)作物的時空重建，但前端實時提取SIFT 特征將對數(shù)據(jù)收集前端造成較大計算負擔(dān)。文獻[11-12]通過計算農(nóng)作物重建結(jié)果的幾何特征，得出農(nóng)作物生長趨勢數(shù)據(jù)，但分析模型根據(jù)農(nóng)作物設(shè)計，難以直接應(yīng)用于樹木重建結(jié)果分析中。

為此，提出了基于時空重建的森林林木形態(tài)監(jiān)測算法，該算法由三部分組成：第一部分是搭建基于同步定位與建圖的視覺數(shù)據(jù)收集前端，用于估計相機位姿；第二部分是后端建立圖像間魯棒數(shù)據(jù)關(guān)聯(lián)，構(gòu)建完整的時空重建；第三部分是設(shè)計通過稠密點云預(yù)處理，進行樹木形態(tài)變化監(jiān)測。

1 V-SLAM數(shù)據(jù)收集前端

1.1 RGB-D針孔相機模型

V-SLAM 數(shù)據(jù)收集前端的傳感器是RGB-D 深度相機，RGB-D 深度相機基于針孔相機模型。該模型包括四個坐標(biāo)系：世界坐標(biāo)系、相機坐標(biāo)系、成像平面坐標(biāo)系與像素平面坐標(biāo)系。假設(shè)現(xiàn)實世界空間點的三維坐標(biāo)為(Xw,Yw,Zw)T，記為Pw。對應(yīng)相機平面坐標(biāo)系坐標(biāo)為(Xc,Yc,Zc)T，記為Pc。對應(yīng)成像平面坐標(biāo)系的坐標(biāo)為(Xi,Yi)T，記為Pi。對應(yīng)像素平面坐標(biāo)系的坐標(biāo)為(u,v)T，記為P。設(shè)成像平面到相機光心O的距離（焦距）為f。

世界坐標(biāo)系與相機坐標(biāo)系轉(zhuǎn)換如圖1 所示。

圖1 世界坐標(biāo)系與相機坐標(biāo)系轉(zhuǎn)換

由模型推導(dǎo)出，世界坐標(biāo)系的空間點轉(zhuǎn)換到相機平面坐標(biāo)系下的公式為：

R為3×3 矩陣表征相機旋轉(zhuǎn)量，t表為3×1 矩陣表征相機平移量，都是系統(tǒng)待計算的數(shù)據(jù)，假設(shè)相機初始位置為空間原點時，旋轉(zhuǎn)與平移矩陣實際上代表相機位姿矩陣T。

相機平面坐標(biāo)系下空間點坐標(biāo)轉(zhuǎn)換到像素平面坐標(biāo)系下齊次坐標(biāo)的公式為：

K為相機內(nèi)參矩陣，根據(jù)實際情況進行標(biāo)定得到。圖2、3 是相機平面到像素平面的轉(zhuǎn)換過程示意圖。dx與dy是成像平面到像素平面橫軸與縱軸的放縮倍數(shù)，u0與v0是像素平面坐標(biāo)原點在橫軸與縱軸平移量。

圖2 相機平面到成像平面轉(zhuǎn)換

圖3 成像平面到像素平面轉(zhuǎn)換

1.2 ORB特征提取與PnP匹配

ORB 特征[13]在相機的旋轉(zhuǎn)與圖像的放縮情況下都有魯棒性，且在相機自動增益、自動曝光與明暗變化的情況具有良好不變性。

提取ORB 特征比提取SIFT 特征[14]耗費更少時間，具體見后面的實驗。數(shù)據(jù)收集前端運行過程中，系統(tǒng)實時性與魯棒性制約相機移動速度，常用方法是借助GPU 加速并行處理，但在數(shù)據(jù)收集前端中使用GPU 過于奢侈。搭載相機的無人機或者手持嵌入式設(shè)備方式下使用GPU 也不便捷。提取ORB 特征可以適應(yīng)普通CPU 的系統(tǒng)需求，降低前端設(shè)備成本，具有更好的實時性。

ORB 特征提取后，采用PnP 匹配模式估計相機位姿變化，通過隨機抽樣一致算法去除大部分匹配錯誤的特征。采用的前端相機可以直接輸出深度圖，由此求出三維空間特征點坐標(biāo)，克服了單目相機采用PnP 匹配模式時需要初始化的缺點。采用非線性方法來計算PnP 匹配下的相機位姿估計，具體細節(jié)見1.3 節(jié)。

1.3 相機位姿估計與最小化重投影誤差

在相機位姿估計中，采用PnP 匹配模式求解，把相機位姿估計看成一個非線性最小二乘問題。如圖4 所示，求解最小二乘問題最優(yōu)解，實際上就是最小化空間點重投影誤差的過程。重投影誤差是估計值與測量值的差值，估計值是重投影的像素坐標(biāo)，測量值是相機測得的像素坐標(biāo)P，由以下公式表示：

圖4 重投影誤差原理

將多空間點重投影誤差最小二乘項相加求最小值，整體誤差取得最小值時，求出的相機位姿就是最優(yōu)解。經(jīng)過RANSAC 排除錯誤匹配對，且優(yōu)化前相機位姿與真實值比較接近，采用最小二乘法最小化重投影誤差不容易陷入局部最優(yōu)解中。

1.4 跟蹤模塊與局部建圖模塊

數(shù)據(jù)收集前端包括跟蹤模塊與局部建圖模塊兩個核心模塊。前端系統(tǒng)結(jié)構(gòu)如圖5 所示。

圖5 前端系統(tǒng)結(jié)構(gòu)

跟蹤模塊在相機當(dāng)前幀中發(fā)現(xiàn)ORB 特征，完成PnP 匹配估計相機位姿。相機當(dāng)前幀輸入前需要預(yù)處理，因為采用的是基于機器人操作平臺（ROS）的RGB-D 深度相機，RGB 模塊與深度模塊工作頻率通常不同，采用ROS 的多傳感器時間軟同步機制實現(xiàn)兩模塊輸入在相近時間內(nèi)的非嚴格對準。跟蹤模塊對每幀圖像做關(guān)鍵幀準入判定，關(guān)鍵幀主要保存了ORB 特征與相機估計位姿。

局部建圖模塊接收跟蹤模塊中的關(guān)鍵幀，建立局部稀疏地圖。稀疏地圖由能在多關(guān)鍵幀都能觀測到的ORB 特征產(chǎn)生。通過維護全局地圖容器，該容器包含了最基本的地圖信息、關(guān)鍵幀信息及其連接信息，局部建圖模塊可以與之交互，將局部稀疏地圖通過必要的整合操作構(gòu)成全局稀疏地圖，并對全局地圖中的關(guān)鍵幀做準出判定，減少需要維護的關(guān)鍵幀。

2 時空重建后端

2.1 SFM實現(xiàn)稠密點云重建

SFM[5-7]與SLAM 類似，前者更關(guān)注于特征點估計以完整恢復(fù)重建結(jié)構(gòu)，后者更關(guān)注相機位姿估計以獲取更準確的定位信息。前端準確估計了相機位姿，后端將結(jié)合當(dāng)前幀RGB-D 圖像，完成稠密點云重建。

2.2 基于詞袋的回環(huán)檢測

SFM 是全局式重建，但前端收集數(shù)據(jù)卻是增量式的，這導(dǎo)致了累積誤差。采用基于詞袋的回環(huán)檢測[15]來降低累積誤差?；陉P(guān)鍵幀的詞袋由每幀ORB 特征生成的描述向量經(jīng)過k-means 聚類構(gòu)成，通過樹狀數(shù)據(jù)結(jié)構(gòu)查詢匹配，最后由TF-IDF進行相似度計算確認回環(huán)。確認回環(huán)后，使用2.3節(jié)的逆向投影圖形變化有界搜索方法，測出相機位姿作為先驗，重新完成時空數(shù)據(jù)關(guān)聯(lián)，解決從不同時間與不同視角觀察同一場景時圖像數(shù)據(jù)關(guān)聯(lián)的關(guān)鍵問題。

2.3 逆向投影圖形變化有界搜索

為建立三維重建之間的時空關(guān)聯(lián)，采取逆向投影有界搜索算法在不同時間不同視角測得的關(guān)鍵幀圖像之間建立數(shù)據(jù)關(guān)聯(lián)。其原理公式如下，設(shè)a、b為兩個時刻：

式（4）-（5）表示，對于空間同一點，已知深度因子、相機內(nèi)參矩陣以及兩時刻相機位姿矩陣時，兩時刻對應(yīng)像素坐標(biāo)是確定的。為減少誤差影響，在該確定像素坐標(biāo)(ub,vb)上設(shè)定一個有界范圍δ重新搜索該匹配像素點(u,v)，以建立魯棒數(shù)據(jù)關(guān)聯(lián)。

2.4 圖優(yōu)化解決全局光束法平差

后端會把輸入的相機位姿與稀疏地圖采用基于圖優(yōu)化理論[15]的G2o 庫進行一次全局的光束法平差，其原理遵循1.3 節(jié)，最小化重投影誤差僅用于局部稀疏地圖建立，此處是全局稀疏地圖優(yōu)化。如圖6 所示，虛線代表相機的運動軌跡，三角形代表單次位姿估計，表征相機運動模型。點虛線代表相機對ORB特征點對應(yīng)空間點的觀測（如Pw2可以在b、c、d時刻的相機位姿中觀測到），觀測模型對應(yīng)1.3 節(jié)提到的估計值。

圖6 圖優(yōu)化模型

2.5 后端系統(tǒng)結(jié)構(gòu)與后期處理模塊

如圖7 所示，后端將前端輸入數(shù)據(jù)進行優(yōu)化后開始稠密重建，通過時空重建模塊完成稠密點云時空數(shù)據(jù)關(guān)聯(lián)，最后輸出帶有時間信息的稠密點云到點云后期處理模塊。

圖7 后端系統(tǒng)結(jié)構(gòu)與后期處理模塊

3 點云后期處理模塊

3.1 基于半徑的異常點濾波

時空重建后端輸出的時空稠密點云帶有許多重建噪聲，通過基于半徑的異常點濾波[17]去除大部分重建噪聲。指定內(nèi)點搜索半徑與相鄰最小點數(shù)，遍歷稠密點云中所有點，點云某點搜索半徑內(nèi)點數(shù)少于相鄰最小點數(shù)，則認為該點為噪聲點予以排除。

3.2 基于RGB的重建噪聲濾波

經(jīng)過異常點濾波之后的點云仍存在部分重建噪聲，來源是位姿估計誤差下重建點的偏移，以遠點代替了真正需要重建的點，遠點往往重建出來帶有特殊的RGB 信息，據(jù)此將可疑點濾除。

4 實驗與結(jié)果分析

4.1 設(shè)備配置與數(shù)據(jù)集

實驗測試計算機采用Intel CORE i5 10thCPU，Ubuntu 18.04系統(tǒng)。相機為Intel RealsensetmD455 深度相機。實驗流程是使用ROS 發(fā)布訂閱話題以獲取相機輸入，采用基于圖優(yōu)化的G2o 庫，以及PCL 庫做點云處理。采集室外實際植物數(shù)據(jù)以及TUM 數(shù)據(jù)集構(gòu)建測試數(shù)據(jù)集。

4.2 ORB、SIFT特征提取對比實驗

TUM 數(shù)據(jù)集中挑取任意同一幀圖像分別提取SIFT 特征與ORB 特征。提取出3 431 個SIFT 特征，耗時0.57 s，平均每個特征耗時為166.13 μs。提取出500 個ORB 特征，用時9 978 μs，平均每個特征耗時為19.956 μs。提取速度提高接近一個數(shù)量級，這僅是提取單幀圖像得到的提升，當(dāng)圖像規(guī)模高達幾百幀時，將節(jié)省更多時間。ORB 特征在同一幀圖像中數(shù)量上比SIFT 特征要少，這在特征匹配步驟中能節(jié)省時間。

4.3 后期處理點云實驗

采用自建的室外樹木數(shù)據(jù)集進行三維重建后進行基于半徑的異常點濾波與基于RGB 的重建噪聲濾波，如圖8 所示，上圖為濾波前，下圖為濾波后，可以看到大部分環(huán)境中的噪聲點都被濾除，保留下的點云更加緊湊與清晰。

圖8 基于半徑的異常點濾波與基于RGB的重建噪聲濾波

4.4 樹木形態(tài)監(jiān)測結(jié)果

采用自建室外數(shù)據(jù)集進行時空重建結(jié)果變化檢測，結(jié)果如圖9 所示。低矮樹植高度變化明顯，三次不同時間的重建能在同一場景建立起數(shù)據(jù)關(guān)聯(lián)，得到高度結(jié)果如表1 所示。當(dāng)收集長期樹木數(shù)據(jù)集時，也可以進行對應(yīng)的形態(tài)監(jiān)測，結(jié)果如圖10 所示。

表1 低矮植株重建測量

圖9 低矮樹植數(shù)據(jù)關(guān)聯(lián)重建結(jié)果（上：3月30日；中：4月13日；下：4月27日）

圖10 樹木數(shù)據(jù)關(guān)聯(lián)重建結(jié)果（上：4月13日；下：4月27日）

5 結(jié)束語

該研究具有三方面優(yōu)勢：改進了數(shù)據(jù)收集前端，跟蹤匹配ORB 特征而不是SIFT 特征，大大減少跟蹤時的計算量；在后端實現(xiàn)了重建結(jié)果的準確時空數(shù)據(jù)關(guān)聯(lián)；在重建基礎(chǔ)上盡量還原樹木信息，測量結(jié)果誤差在5%以內(nèi)。后續(xù)工作可以通過樹木形態(tài)監(jiān)測結(jié)果估算樹木凋落可燃物量。