基于馬爾科夫隨機場的三維激光雷達路面實時分割

朱 株，劉濟林

（浙江大學 信息與電子工程學系，浙江 杭州310027）

精確感知三維環境是地面自主智能車輛的一項重要工作.為保證安全行駛，自主車輛通常裝有激光雷達和相機組成的傳感器系統，用于探測障礙區域和地面區域.

早期自主機器人平臺一般采用單線激光雷達作為距離傳感器，具有數據量小及響應速度快的優點，但是存在探測范圍有限的缺點.隨著VELODYNE 64線激光雷達的問世，這種三維激光雷達開始替代單線激光雷達，廣泛用于自主車平臺［1-3］.

現有的三維激光雷達地面分割算法主要有如下2類.

1）基于柵格地圖的地面分割算法.這類算法將數據柵格化，投影到一個柵格地圖中，分析每個柵格中的三維點云特性來確定該柵格的障礙物屬性.Kammel等［4］提出了一種基于柵格中點云最大高程差的地面檢測方法，如果柵格最大高程差小于某個閾值，柵格將被標記為地面區域.Himmelsbach等［5］將三維數據以極坐標的方式柵格化，并且用非參數的方法來擬合地平面，分離地面和障礙物.Guo等［6］將柵格圖以無向圖的方式處理，用馬爾科夫隨機場將柵格分為4類，從中選取地面.雖然基于柵格的檢測方法較為穩定，但是該方法的檢測精度取決于柵格大?。ㄍ?0 cm×20 cm），相比于原始三維點云數據精度（0.2 cm），柵格算法精度較??；受限于計算速度，柵格地圖范圍不大，典型的250×150柵格圖相當于50 m×30 m的檢測范圍，相比于原始數據200 m×200 m的檢測范圍，大量數據丟失；由于三維數據分布的非均勻性，大量柵格內并沒有數據，造成了存儲和處理的浪費.

2）基于圖的地面檢測.三維激光雷達數據結構的特點：任意三維點與其在同一掃描線上相鄰角度的點以及前后掃描線相同角度的點組成4鄰域系統，可以將三維點云以無向圖的方式組織起來.Montemerlo等［7］使用三維數據環之間的距離來判斷單個三維點是否屬于地面區域（三維激光雷達地面反射點在二維平面上呈現出環的樣式，障礙物點的環間距要遠小于地的面點間距）.Moosmann等［8］使用局部凹凸性分割不同物體.Douillard等［9］采用Mesh方式構造圖，并且通過計算梯度來確定地面點.

相比于柵格類的分割算法，基于圖的算法精度較高，并且能處理全部雷達數據，但是現有算法均使用單個三維點作為圖的節點，特征單一，容易受到噪聲干擾，在起伏道路環境中魯棒性較差.

本文提出一種新的基于圖的算法.將掃描線χ-y平面上的投影曲線進行分段，構建以線段為節點的馬爾科夫隨機場，通過對線段特征分析建立能量函數，并采用圖分割算法求全局最優解.算法不僅能在100 ms內完成，而且在不同場景下分割效果都優于現有基于圖的分割算法.

1 算法概述

觀察雷達數據在χ-y平面上的投影可以發現，地面數據和障礙物數據在線型上的明顯區別.地面區域如圖1（a）環狀圖案所示，單條掃描線障礙區域放大后如圖1（b）、（c）所示：在同一掃描線上呈斷開狀或是折角狀.據此，本文提出的地面分割算法采用線段特征.

圖1 掃描線在x-y平面上投影的樣式Fig.1 Projection of scan lines onχ-y plane

圖2 地面分割算法流程圖Fig.2 Flow chart of ground segmentation algorithm

圖2為地面分割算法的流程框圖，該流程分為3個步驟：掃描線投影分段；線段端點精確定位；建立無向圖馬爾科夫隨機場并且求解.

首先輸入激光雷達原始數據，并且進行預處理，包括異常距離三維點濾除以及采用GPS進行點云位姿修正.

考慮到在χ-y平面的投影中每條掃描線處于同一直線上的相鄰點屬性相同，算法首先對每條掃描線進行快速分割.由于存在噪聲，分割后相鄰線段端點（角點）可能出現定位錯誤，通過角點檢測，準確選取端點，增加算法的精確度.采用線段作為圖節點，不僅大大減少了節點數量（從每幀大約130 000個點減少到20 000個點左右）加快圖分割的處理速度，而且能獲取多個基于線段的特征，增加了算法的魯棒性.最后通過圖分割算法，獲取全局優化后的分割結果，提升算法的性能.

2 掃描線投影分割

由于地面不平整以及激光雷達本身具有的檢測噪聲，掃描線在χ-y平面上投影含有噪聲.采用基于最大模糊線段的方法對掃描線投影進行快速分割.

模糊線段（blurred segment，BS）定義為介于直線l1：aχ+by=μ以及直線l2：aχ+by=μ+ω之間的有限離散二維點集，如圖3所示.ν為BS的寬度：

圖3 模糊線段Fig.3 Blurred segment

一個二維點集模糊線段的參數可以通過計算該點集凸包絡的支持線來獲?。?0］，計算復雜度與點的數目成線性關系.分割算法流程如下：

1）將單條掃描線的映射點順序輸入，如果當前點和前一個點的距離大于閾值T p-d，則將前一個點作為當前線段的終止端點，當前點作為新一條線段的起始端點.相鄰2個點的屬性都標記為分離點.

2）若當前點和前一個點的距離小于閾值Tp-d，則進行模糊線段參數計算.如果模糊線段的寬度大于閾值Tv，則將前一個點作為當前線段的終止端點，當前點作為新一條線段的起始端點.相鄰2個點的屬性都標記為連接點.

三維激光雷達在不同距離下有不同的分辨率，Tp-d和Tv均為距離相關的自適應參數：

式中：D為當前點在χoy平面的投影到原點的距離，μ1、μ2為比例因子.由于 VELODYNE HDL64E雷達在不同距離，掃描線在地面上的點間隔不同，對此閾值使用比例因子與距離的乘積.考慮到在10 Hz運行下，激光雷達角分辨率最大為0.18°，2個連續掃描點間距理論上為D×0.18π／180.作為閾值將該值適當放大，取3倍值，計算中μ1=0.009 42.式（3）用于確定分割線段的最大寬度.該閾值不僅受到距離的影響，也受到不同地形環境的影響.本文使用經驗值，平坦地形下μ2=0.003，在顛簸環境下μ2=0.005.

分割結果見圖4，相鄰的線段用2種不同的灰度表示.

3 角點檢測

圖4 模糊線段分割結果Fig.4 Segmentation by blurred segment

模糊線段有一定的寬度，當連接點很多時不能準確定位在局部曲率最大的角點上.可以通過角點檢測算法重新定位連接點.

實空間R上的算子A被定義為以下形式：

二維曲線Γ定義為參數化方程a（s），核函數φ定義為高斯函數，其中u為均值，t為方差：

將a（u）在u=0附近進行泰勒展開：

可得：

積分算子在χ處的值在曲線法方向nχ上正比于曲線曲率κχ，選擇局部曲率最大的點作為角點.角點定位前后分割結果如圖5（a）、（b）所示.

圖5 角點定位前后不同的分割結果Fig.5 Segment results before and after corner detection

4 建立無向圖馬爾科夫隨機場

馬爾科夫隨機場框架建立了全局標記優化模型，常被用于解決場景分類問題［11］.定義G= （V，E）為無向圖，節點v i∈V；邊 （vi，v j）∈E對應于相鄰的節點對.每條邊（vi，vj）∈E被賦以非負權重W（v i，v j），用于描述相鄰節點vi及v j間的相似度.

標記推測算法能為每個節點找到概率最大的標記，例如圖分割［12］.定義L為一個有限的標記集合，f為標記賦值函數，f（v i）為每個節點賦以標記.假設標記在區域內變化緩慢，而在區域邊緣變化劇烈.標記賦值的準確度由能量方程來衡量：

式（8）中等號右側第一項是數據項，是對節點本身的約束；第二項是平滑項，是對能量方程解的平滑性約束［13］.

本文算法利用激光雷達原始數據結構，構建無向圖.把每條線段當作一個節點，所有包含與當前節點線段相鄰三維點的線段，與該節點線段組成一個相鄰集合.定義與節點相鄰的同一掃描線線段為側鄰線段（SN），定義與節點相鄰的前一掃描線線段為上鄰線段（UN），定義與節點相鄰的后一掃描線線段為下鄰線段（DN）.根據定義，每條線段只有左右2條側鄰線段，可以有多條上下鄰線段.

根據線段內點數依據閾值Tp-n將線段分為2類：長線段（Sl）和短線段（Ss）.由于地面較為平坦，長線段大多屬于高概率地面區域（Agh），只有在以下情況長線段被定義為高概率的障礙區域（Aoh）：

1）長線段端點標記為分離點，且該點距雷達原點的距離與側鄰線段分離點距原點距離的差大于閾值Td.

2）長線段端點標記為連接點，且它與側鄰線段夾角與直角的差小于閾值Ta.

3）長線段與上鄰線段的梯度大于閾值Tg.梯度定義為線段包含三維點平均高度差與線段在χ-y平面上的平均距離的比值.

短線段大多屬于高概率障礙區域，只有在以下情況長線段被定義為高概率的地面區域：

1）短線段端點標記為連接點，與側鄰定義為高概率地面區域的長線段夾角與直角的差值大于閾值Ta.

2）短線段與上鄰定義為高概率地面區域的長線段的梯度小于閾值Tg.

本文采用線段點數閾值的經驗值Tp-n=6.Td、Tg均為障礙物閾值，障礙物會造成同一掃描線上線段的分離以及相鄰掃描線上線段間產生梯度.這2個閾值的取值與地形或者障礙物區分敏感度有關.實驗中，城市道路的Td=40，Tg=tan 12°（0.212）；鄉村顛簸道路的Td=60，Tg=tan 25°（0.466）.Ta主要用于區分城市道路中的馬路牙子和地面，實驗中Ta=30.

根據以上準則，可得：

η是符合判斷準則下，標記為地面的先驗概率值，實驗中η=0.8.

為保證局部區域標記的一致性，定義平滑項：

φ（Δh）是根據相鄰線段在三維空間中的平均高度差Δh來度量相似性的函數，k為經驗系數（實驗中k=0.01）：

最后使用圖分割的方法（Graph-Cut）優化能量方程，得出標記結果.

5 實驗過程及其結果

實驗平臺如圖6所示.三維激光雷達安裝在車頂部，除少部分盲區外，能較為全面地獲取車體周圍的三維環境信息.

圖6 三維激光雷達與自主車Fig.6 3D LIDAR and autonomous landing vehicle

分別在城市平坦路面場景和鄉村起伏道路場景采集數據，將本文算法結果與局部凹凸性算法［8］、Mesh梯度算法［9］以及人工標記的結果進行對比.圖7中第1行是城市平坦路面場景；第2行是鄉村起伏道路場景；圖7（a）、（f）是對2種場景的相機圖片；圖7（b）、（g）是三維激光雷達數據人工標記的結果，在本文中該圖為灰度圖，實際中該圖的其中地面區域標記為綠色，障礙區域標記為紅色，黑色區域為未知區域；圖7（c）、（h）為局部凹凸性算法的結果；圖7（d）、（i）為 Mesh梯度算法結果；圖7（e）、（j）為本文算法的結果.

在城市平坦路面場景下，Mesh梯度算法（圖7（d））結果和本文算法（圖7（e））的結果都較好，接近人工標記的結果（圖7（a））.局部凹凸性算法（圖7（c））結果較差，以局部平面的法向量方向作為判斷準則，用單個三維點的鄰接4點擬合平面.因此，在構建圖時，刪除了長度大于某個閾值的邊，以防止遠距離點用于平面擬合，造成障礙物邊緣判斷錯誤.由于在遠距離，激光雷達相鄰點本身相距很遠，導致刪除長邊后的圖出現大量未知區域，影響了檢測效果.

在鄉村起伏道路場景下，地面起伏變化大，使得基于單個三維點梯度特征的Mesh梯度算法（圖7（i））檢測效果下降，出現了大量地面區域漏檢點；路邊大量低矮植物造成了地面區域誤檢點.局部凹凸性算法（圖7（h））不僅在遠距離存在未知區域；在近距離，由于地面不平整，使得局部平面的法向量變化很大，導致地面區域漏檢，路邊大量低矮植物又造成局部平面的法向量變化緩慢，產生地面區域誤檢點.本文算法（圖7（j））效果優于Mesh梯度算法和局部凹凸性算法：由于本算法使用線段的梯度特征，減少噪聲對單點梯度的影響，線段夾角特征和分離點距離特征的判斷增加了地面區域判斷的可靠性，最后使用馬爾科夫隨機場進行全局優化使得結果具有更好的魯棒性.

圖7 2種場景的分割結果Fig.7 Segmentation results of two kinds of scenes

式中：NTP為地面點被正確標記的數目，NFP為障礙點被標記為地面點的數目，NFN為地面點被標記為障礙點的數目，NTN為障礙點正確標記的數目.

城市場景下的算法準確率如圖8所示，橫坐標為連續幀的編號，縱坐標為準確率.Mesh梯度算法和本文算法都有較高的準確率.

鄉村起伏道路場景下的算法準確率如圖9所示，本文算法相比于Mesh梯度算法和局部凹凸性算法，不僅準確率高于這2種算法，而且穩定性也要好于這2種算法.在起伏道路環境中，地面梯度變化劇烈，基于單點局部特征的算法穩定性差不適用于這類場景.本文算法通過多特征全局優化，獲得了較高的穩定性，更適合鄉村起伏道路場景.

本文算法使用C++實現，在Intel酷睿i5雙核2.8 GHz CPU、2 G內存的自主車平臺上運行時間

圖8 城市場景下算法準確率Fig.8 Accuracy of algorithms under urban scene

圖9 鄉村起伏道路場景下算法準確率Fig.9 Accuracy of algorithms under rural bumpy scene

在2種場景下運行100幀數據序列，通過正確率（R）對不同算法的結果與人工標記的真值進行量化比較.為每幀61 ms，其中原始激光雷達距離數據轉換成三維點云25 ms，掃描線分割15 ms，馬爾科夫隨機場建立和求解21 ms，符合實時運行的要求.

6 結 語

本文提出了一種適用于多類場景的實時三維激光雷達地面分割算法.算法首先使用最大模糊線段方法對三維掃描先在平面上的投影曲線進行分割，然后采用角點檢測算法重新精確定位線段端點，最后建立無向圖馬爾科夫隨機場并求解標記地面區域.在城市平坦路面和鄉村起伏路面下的實驗結果表明：本文算法相比于其他算法具有較高的魯棒性，能在起伏路面場景下穩定可靠地為自主車輛提供路面區域信息.

（

）:

［1］GLENNIE C.Calibration and kinematic analysis of the Velodyne HDL-64E S2 Lidar sensor［J］.Photogrammetric Engineering and Remote Sensing，2012，78（4）：339- 347.

［2］HASELICH M，BING R，PAULUS D.Calibration of multiple cameras to a 3D laser range finder［C］∥International Conference on Emerging Signal Processing Applications（ESPA）.Las Vegas：IEEE，2012：25- 28.

［3］楊飛，朱株，龔小謹，等.基于三維激光雷達的動態障礙實時檢測與跟蹤［J］.浙江大學學報：工學版，2012，46（9）：1565- 1571.YANG Fei，Zhu Zhu，GONG Xiao-jin，et al.Real-time dynamic obstacle detection and tracking using 3D Lidar［J］.Journal of Zhejiang University：Engineering Science，2012，46（9）：1565- 1571.

［4］KAMMEL S，PITZER B.Lidar-based lane marker detection and mapping［C］∥Intelligent Vehicles Symposium.Eindhoven：IEEE，2008：1137- 1142.

［5］HIMMELSBACH M，HUNDELSHAUSEN F，WUENSCHE H.Fast segmentation of 3D point clouds for ground vehicles［C］∥Intelligent Vehicles Symposium（IV）.San Diego：IEEE，2010：560- 565.

［6］GUO C，SATO W，HAN L，et al.Graph-based 2D road representation of 3D point clouds for intelligent vehicles［C］∥Intelligent Vehicles Symposium（IV）.Detroit：IEEE，2011：715- 721.

［7］MONTEMERLO M，BECHER J，BHAT S，et al.Junior：The stanford entry in the urban challenge［J］.Journal of Field Robotics，2008，25（9）：569- 597.

［8］MOOSMANN F，PINK O，STILLER C.Segmentation of 3D lidar data in non-flat urban environments using a local convexity criterion［C］∥Intelligent Vehicles Symposium.Xi’an：IEEE，2009：215- 220.

［9］DOUILLARD B，UNDERWOOD J，KUNTZ N，et al.On the segmentation of 3D LIDAR point clouds［C］∥International Conference on Robotics and Automation（ICRA）.Shanghai：IEEE，2011：2798- 2805.

［10］DEBLED RENNESSON I，FESCHET F，ROUYERDEGLI J.Optimal blurred segments decomposition in linear time［C］∥ Discrete Geometry for Computer Imagery.Berlin Heidelberg：Springer，2005：371- 382.

［11］BESAG J.Spatial interaction and the statistical analysis of lattice systems［J］.Journal of the Royal Statistical Society.Series B（Methodological），1974，36（2）：192- 236.

［12］BOYKOV Y，VEKSLER O，ZABIH R.Fast approximate energy minimization via graph cuts［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，2001，23（11）：1222- 1239.

［13］BOYKOV Y，KOLMOGOROV V.An experimental comparison of min-cut／max-flow algorithms for energy minimization in vision［J］.IEEE Transactions on Pattern Analysis and Machine Intelligenc e，2004，26（9）：1124- 1137.