潘 斌,蔡志剛,劉微微,曾玉娟,王 博
(1. 嘉善縣自然資源和規劃局,浙江 嘉興 314100; 2. 浙江省測繪科學技術研究院,杭州 310000; 3.南京航空航天大學, 南京 210010)
近年來激光雷達技術發展迅速,其可以快速獲取空間三維信息,進行大場景的空間探測,廣泛應用于城市分析、無人駕駛、電力巡檢等領域,激光雷達測量結果以點云形式存儲,但點云數據通常密度高、體積大,對數據處理提出挑戰。特別是對大規模無組織點云,構建以特征為索引的點云查詢方式能提高查詢效率,對后續的區域分割和三維重建等處理有重要意義。三維點云中基于點的索引已經有很多研究成果,包括格網索引[1]、k-d樹索引[2]、八叉樹索引[3]和其他混合索引方法[4-7]等。基于線特征的索引卻很少被關注,特別是在直線特征豐富的人造場景。線特征既可以是邊界點的統一表現,又可以是平面的交線集合,對點云的表達和重建工作有著重要意義。進行線特征的索引前,需對三維線段進行檢測和提取,倪歡等[8]將三維線段的檢測工作主要分為基于圖像的間接法,基于點和基于面的直接法。
基于圖像的方法依靠三維點云對應的二維圖像提取特征,或將三維點云轉化二維圖像,在圖像上提取二維線段,再將結果重投影至三維空間。X. Ge 等[9]將點云進行分割,投影至不同平面,從不同的二維平面中獲取特征點。LI H等[10]利用三維點云對應的圖像提取邊緣特征,再與三維點云進行匹配確定邊緣點,通過形態學方法得到最終特征。基于圖像的方法需要考慮投影策略,構建不同平面之間的對應關系,在轉化過程中會造成信息的損耗,對三維點云的幾何信息利用不充分。
基于點的方法直接對邊界點進行檢測,然后通過線段擬合獲得點云的線特征。陳華偉等[11]根據曲率突變點隱含特征線的結論,通過區域聚類提取備選點集,按照一定的搜索原則進行特征線的細化,在機械模型和藝術品模型上取得良好的實驗效果。T. Hackel 等[12]從點的鄰域中提取特征,使用二元分類器預測每個點輪廓分數作為構建候選輪廓的基礎,在候選輪廓中選擇最佳輪廓作為特征。Demantké. J[13]提出了一種多尺度方法PCA計算激光雷達點云上的穩健幾何特征,在不同半徑大小的球面鄰域上計算三維特征,基于局部結構張量的特征值組合描述鄰域的形狀,檢出線性局部特征。基于點的方法直接對特征點進行處理,因此對噪聲的魯棒性較低。
基于面的檢測方法通過法線或曲率特征將點云聚類分割為不同面,通過面面相交或面上邊緣檢測得到特征線。LIN Y B等[14]首先將點云根據不同視點轉換成陰影圖像的集合,在陰影圖像中提取二維線段及其支撐區域,再反投影至三維點云空間得到特征線,但該方法只能提取平面交線。孫憲龍等[15]通過投影直方圖進行點云分割,利用霍夫變換進行直線檢測,對直線進行求交獲得建筑物的墻線特征。張帆等[16]將三維點云投影至球面坐標系,在全景圖上進行邊緣檢測得到包含直線特征點的邊緣點,通過霍夫變換提取直線段。Lin Y等[17]提出基于面的檢測方法提取線特征,通過局部k-means聚類將點云分割成面,再在平面上進行邊界點提取和線段分組。Sheng Q等[18]將點云投影至不同擬合平面,對平面進行網格化和二值化,在二值化網格平面采用差分編碼的方式提取邊界點,對邊界線特征進行差分編碼以便同名直線匹配。基于平面的方法容易面臨確定交叉線段端點的問題,也容易出現線段的錯檢和漏檢。
文中針對城市無組織點云數據,提出一種城市機載激光雷達點云的直線特征提取方法。首先對點云基于區域生長原理進行平面聚類,分割成不同平面,然后將對應點云投影至擬合平面,通過穩健的二維直線提取方法提取出良好的線特征集合,最后通過空間直線索引方法完成對線特征的索引。
文中方法提取城市點云直線特征方法為:基于區域增長的點云平面聚類、基于網格的邊界特征提取、建立線特征空間索引。首先計算法向量和曲率對點云進行排序,利用區域增長算法進行聚類,完成點云平面的分割與擬合;然后將點云投影至擬合平面,進行網格化和二值化,檢測直線特征并進行直線擬合;最終將直線重投影至三維空間,建立空間索引。具體技術路線如圖1所示。

圖1 技術路線圖
點云的幾何特征信息能表現物體的外輪廓,點云法向量和曲率是點云數據幾何特征的代表。首先基于PCA算法對點云法向量進行估計,再利用基于曲率排序的區域增長算法對點云進行聚類并分割,對分割后的區域進行平面擬合,最后對相鄰平面進行合并判斷,可以準確獲取點云平面。
1.1.1 點云法向量和曲率的計算
曲面上一點的法向量與該點的切平面垂直,點云的法向量估計一般利用鄰域點構建局部曲面,通過計算局部曲面的法向量得到結果。

(1)
最小二乘原理擬合局部平面求解法向量,可等價為求解式(2)的最小值:
(2)
其中,n為局部平面的法向量,它對應平面在空間中三個正交基投影方差最小的向量。
P的最小值求解可以轉化為計算協方差矩陣的特征值,最小特征值對應的特征向量即為法向量。
在求解前對所有鄰域點進行去重心化以簡便運算,新的鄰域點集表示為{pj|(xj,yj,zj),j=1,2,…,N},其中:
(3)
協方差矩陣為:
(4)

(5)
為保證法向量的方向全局一致,可以計算兩個相鄰點云的法向量內積,若內積為正,則表示同向,內積為負,表示反向,將鄰域點方向調整為相反。
曲率是指曲線或曲面上點的曲率,在三維點云中,某點鄰域范圍曲率變化越大,則越不光滑,選擇這樣的點作為邊界點的候選點。
曲率計算式為:
(6)
式中:cu為當前點的曲率,λ0,λ1,λ2為協方差矩陣A奇異值分解求得的特征值。
1.1.2 區域增長
一般而言,區域增長方法是根據同一物體區域內像素的相似性質來聚集像素點的方法,從初始區域(種子點)開始,將相鄰的具有相似性質的像素或其他區域歸并到目前的區域中從而逐步增長區域,直至沒有可以增長的點為止。在點云中,利用區域增長方法,選取恰當的種子點和相似判別依據,可以提取出目標點云中的所有平面。當點云曲率越小,該點鄰域變化越小,區域越平坦,符合平面特征,從平坦區域開始增長可以有效減少區域總數,因此在曲率小的點云中選取種子點。文中相似判別規則如下:


(7)
若Dplane和Dpo同時小于設定閾值,則認定點pj屬于當前區域。
區域增長聚類的步驟如下:
1)計算所有點曲率,并按照曲率從小到大排序,曲率小于一定閾值的點作為種子點。
2)選擇曲率最小點作為首要種子點,利用K鄰近算法搜索鄰域點,如果鄰域點曲率小于等于設定閾值,則作為新的種子點加入種子點集。若大于設定閾值,對其進行判斷:①種子點與該鄰域點的法向量內積結果大于設定閾值;②該鄰域點的Dplane和Dpo,同時小于閾值D。滿足上述兩個判斷條件的鄰域點歸并到當前種子點的區域。
3)設置區域最小容納點數和最大容納點數。對于符合判斷條件的點,將其從原始點云剔除,歸入種子點所在的點集。剩余點云中生成區域的點數小于最小容納數時停止增長,得到區域分割結果。
經過區域增長后,可得到點云平面初步分割結果,但由于鄰近點及閾值選擇存在一定誤差,屬于同一平面的點云可能被誤分割為多個平面,因此需要進行相鄰平面的合并判斷。當相鄰平面的平面夾角小于規定閾值T,進行合并處理。
文中采用平面投影法檢測邊緣特征,首先將區域點云投影到擬合平面上,進行網格化和二值化,然后在二值化圖像上提取邊界特征。


圖2 三維點投影為二維點示意圖
(8)
即可得到空間點投影到擬合平面上的投影坐標。
將平面所屬點全部投影至擬合平面后,建立網格對離散點進行像素化和二值化。獲取投影點的橫縱坐標最大值和最小值,假設其最大值、最小值分別為xmax,xmin,ymax,ymin。設定單位網格邊長為s,網格的行數R和列數C為:
(9)
其中,y=[x]為取整函數。
(10)
如圖3所示,將投影點分配完成后,部分網格無投影點存在。由于邊界存在于包含落點的網格中,可對網格構建二值圖像掩膜,將有投影點和無投影點的網格進行更清晰的區分:有投影點的網格賦值1;無投影點的網格值賦0,網格轉換成僅由像素值為0和1表示的二值圖像。

圖3 網格劃分示意圖
為了增強邊界以及避免提取一些無必要的內輪廓,對生成的二值圖像進行形態學的閉運算操作,利用結構元素的原點遍歷二值圖像所有像素點,其中,若結構元素中至少有一個1與原圖重合,則膨脹結果賦中該像素點賦值為1,否則為0;若重合區域中結構元素的1均能對應于原圖像中的1,則賦值結果中該像素點賦值為1,否則為0。
在二值化網格平面中進行直線檢測與擬合,獲取線段的起始點和結束點坐標,將點坐標重投影至三維空間后,得到能夠表達線段的起始點坐標(xs,ys,zs)和結束點坐標(xe,ye,ze),將所有線段存入線段集中,便于索引與查詢。
先建立點云的固定網格索引,根據查詢內容找到符合條件的線段集合,在網格索引中找到起始點所在網格,以起始點所在網格開始進行索引,如圖4所示。假設網格邊長為la,線段在x,y,z3個方向上的斜率分別為kx,ky,kz,構建新的搜索點(xi,yi,zi),即:

圖4 線特征索引示意圖
(11)
xi=xs+kxla,
yi=ys+kyla,
zi=zs+kzla.
(12)
相應地下一個搜索點在此基礎上進行計算:
xi+1=xi+kxla,
yi+1=yi+kyla,
zi+1=zi+kzla.
(13)
對搜索點進行判別:①判斷搜索點坐標是否大于結束點,大于時,直接搜索結束點所在網格,記錄網格節點位置,結束搜索;②當小于結束點,判斷搜索點是否仍在當前網格中,若是則直接計算下一搜索點替代當前搜索點,繼續進行判斷;③若已超出當前網格,記錄新網格的節點位置,繼續構建下一搜索點進行判斷。
通過以上流程,可求出線段經過的所有網格,完成線特征的空間索引。
為驗證線特征索引方法,文中實驗選取浙江省嘉興市嘉善縣地區的機載激光雷達點云進行實驗,該點云區域中包含多種建筑類別,框選其中典型區域進行實驗,點云具體信息如表1所示。實驗環境為IntelCorei3-6100CPU和4 GB內存計算機。首先進行線特征的初次提取,得到線段集,通過線特征索引查詢線段鄰域點,組成新的點云再次提取線特征,直至滿足迭代條件。并與平面相交提取直線特征的方法進行對比。

表1 機載激光雷達點云信息
點云影像如圖5(a)、5(c)所示,按高程進行灰度賦色,點云缺失處一般為水體或建筑物陰影。區域1為居民區,建筑低矮且分布密集;區域2為高樓,建筑高且分布較為稀疏。利用本文方法進行線特征提取,整體結果如圖5(b)、5(d)所示。

圖5 機載雷達點云影像和線特征提取結果
為直觀表達線特征提取的結果,將點云與線特征提取結果融合顯示如圖6所示。區域1范圍較大,建筑密集,選取局部區域進行展示,如圖6(a)所示;圖6(b)展示了區域2提取結果,區域2中主要為高層建筑,墻體線較為密集,選取高程不同的三處建筑物進行局部放大,展示結果見圖6(c)。

圖6 線特征和點云結合顯示結果
根據實驗結果,兩個區域內均提取出分布均勻的建筑物線特征。區域1為居民區,建筑低矮且密集分布,房屋類型較為一致。屋頂線特征提取結果很好,大部分樓棟的屋頂線結構完整,且保留了坡形屋頂的完整結構。但是由于居民區建筑密集,受限于掃描角度和屋頂遮蔽,機載雷達采集的點云在垂直面上的分布較為稀疏,無法形成完整的墻面點云,在二值化網格平面上提取到的邊界特征有限,因此墻體線提取效果弱于屋頂線。區域2以高樓建筑為主,建筑線特征提取完整,墻體線分布均勻,周邊的低矮建筑屋頂線特征也完整保留。但是由于屋頂結構較為復雜,短線段特征較多,后處理過程過濾異常線段時將部分短線特征被刪除,從而導致部分建筑屋頂線特征缺失。
為驗證文中方法效率,與平面相交法進行對比實驗,在兩個區域內隨機選擇一個小區域進行精確率和召回率的比較。其中精確率表示被正確提取的線段數占被提取出的線段總數的比例;召回率表示被正確提取的線段數占測試區域所有線段數的比例。表2展示了平面相交方法與文中方法在兩個區域中提取到的線段數、耗時和小區域的精確率和召回率。圖7展示了兩種方法的對比效果。

圖7 與平面相交法的對比效果

表2 與平面相交法的對比結果
根據實驗結果,文中方法耗時比平面相交法長,在區域1和區域2的耗時分別是平面相交法的1.21倍和1.18倍。這是因為平面相交法僅需通過判斷平面是否相鄰以及求相鄰平面的交線方程提取線段,而文中方法需要在獲取平面后進行二值化網格化、直線提取、線特征索引和迭代等過程。
但文中方法提取到的線段數遠高于平面相交法獲取的線段數,其中區域1文中方法提取到的線段數是平面相交法的5.44倍,區域2文中方法提取到的線段數是平面相交法的5.58倍。這是由于平面相交法僅能提取兩平面的交線,無法求解屋頂邊線方程,并且在面關系復雜的局部可能有失效情況出現,比如區域1墻面點云缺失,影響交線方程求解,導致平面相交法無法提取到該處的線段特征,區域2中高層建筑設計復雜,屋頂線輪廓曲折且部分墻面缺失,平面相交法無法取得好的效果。而文中方法通過平面分割投影提取線段,除了邊界直線,還能提取建筑的交面線和邊緣線等線特征,例如區域2的屋頂線特征,可以通過提取屋頂平面的邊界輪廓擬合線段特征。因此文中方法與平面相交法雖然在精確率上差距較小,但在召回率上文中方法明顯優于平面相交法。計算線段提取效率在區域1,平面相交法每秒提取7.157條線段,文中方法每秒提取32.141條線段;在區域2,平面相交法每秒提取4.328條線段,文中方法每秒提取20.437條線段。文中方法在線段的提取效率上遠高于平面相交法。
為定量評價本文算法的線特征提取效果,選擇直線擬合精度作為評價指標,直線擬合精度是擬合點到擬合直線段的平均距離誤差,其中擬合點是直線檢測中參與線段擬合的點集。根據建立的線特征索引結構,隨機選擇500條直線參與計算,記錄所有直線中的最大誤差、最小誤差與平均誤差。為確保結果穩定性,對每個區域隨機選擇相同數量的直線進行兩次驗證實驗,實驗結果如表3所示。區域1與區域2的點平均距離分別為0.123 3 m和0.131 9 m,提取的直線最大誤差與點平均距離相近,平均誤差小于點平均距離的一半,對比之下可以看到文中算法具有較高的線特征提取精度。

表3 直線提取精度對比 m
與基于面的平面相交法進行對比,文中提取建筑物的線段特征具有較大優勢。與基于點的比較有代表性的PCA方法進行對比分析,對比效果如圖8所示,PCA方法能夠提取建筑物傾斜直線上的特征點,但提取效果不佳,且需要人工選取搜索范圍參數閾值,適用性一般。相比于PCA方法,文中所提算法對建筑物特征點具有良好的提取效果,提取結果中噪聲較少,可以較好地提取相同方向的相鄰直線上的特征點,精確區分出鄰近直線段。

圖8 與PCA方法的對比效果
文中提出了一種城市機載激光雷達點云直線特征提取方法。該方法基于區域增長的點云平面分割和穩健的二維線段檢測,并通過三維線段索引方法進行迭代處理雜亂線段,最終完成對城市點云的直線特征的提取。與傳統的基于點的方法不同,文中方法可以更直觀地展示城市結構特征,方便后續建立點、線、面三種幾何表征的關系。但是文中方法主要基于直線特征進行,如何增強地面道路和城市綠化等不規則線段的提取效果和如何進行曲線索引都是今后需要進行探索的內容。