李 峰,崔希民,袁德寶,王 強,吳亞軍
(1.防災科技學院,三河 065201;2.中國礦業(yè)大學(北京)地球科學與測繪工程學院,北京 100083)
利用高分辨率遙感影像自動提取建筑物是攝影測量和計算機視覺領域的熱點之一。3D 建筑物在城市規(guī)劃與管理、緊急事態(tài)處置以及災害救援等方面有著重要的應用。由于缺少高程信息,以及受到建筑物陰影和樹木遮閉的影響,利用航空或衛(wèi)星影像數(shù)據(jù)提取建筑物信息的精度往往不高。最近新興的機載Li-DAR(light detection and ranging)技術為城市建筑物的識別與量測提供了新方法。加載在飛機平臺上的機載LiDAR 系統(tǒng)可以直接獲取地表、建筑物、樹木和汽車等地物對象的大量3D 點云。與航空或衛(wèi)星影像相比,機載LiDAR 點云不受陰影和投影差的影響,因此在城市建筑物提取方面有著較大優(yōu)勢。
從LiDAR 點云中提取建筑物的通常做法是:首先,從非地面點云中分離出地面點云并生成數(shù)字高程模型(digital elevation model,DEM);然后,用全部點云生成的數(shù)字表面模型(digital surface model,DSM)減去DEM 得到包含非地面點云的歸一化的數(shù)字表面模型(normalized DEM,nDSM);最后,從nDSM 中區(qū)分建筑物點和植被點。該過程用到的點云和影像信息包括高程、面積、形狀、表面的平滑程度、色彩和紋理[1-2]。Zhang 等[3]聯(lián)合區(qū)域生長法和平面擬合的方法借助高差從非地面點云中分離出建筑物和樹木;Rottensteiner 等[4]將高差、歸一化差值植被指數(shù)(normalized difference vegetation index,NDVI)和表面粗糙度參數(shù)應用到證據(jù)推理理論中,從而提取建筑物;Khoshelham 等[5]將首次和末次回波DSM 的高差以及從彩紅外影像中提取的NDVI應用到貝葉斯決策函數(shù)中進行建筑物分類;Zhou等[6]通過點云分布的規(guī)則性、水平性、平整度和法線分布等參數(shù),使用支持向量機(support vector machine,SVM)法探測植被點類,根據(jù)區(qū)域生長法確定建筑物面片;Meng 等[7]基于形態(tài)學的方法識別并恢復建筑物區(qū)域,以面積和緊湊性剔除非建筑物區(qū);Salah 等[8]引入灰度共生矩陣的同質(zhì)性和熵來區(qū)別建筑物和植被;于海洋等[9]結(jié)合面向?qū)ο蟮姆椒ê蚐VM 技術提取了地震后復雜環(huán)境下倒塌建筑物的信息;李婷等[10]采用點云的高程、法向量和平面投影密度信息從車載激光點云中分離建筑物點。
目前多數(shù)算法需要結(jié)合高分辨率的多光譜影像才能取得較高的精度,計算方法復雜、效率較低;另外,建筑物的投影差、陰影、影像的空間分辨率、時間差異和配準等造成的融合誤差很難消除,建筑物提取精度受到影響。因此,本文首先直接基于原始機載LiDAR 點云,以漸進式不規(guī)則三角網(wǎng)(TIN)加密法過濾出非地面點云,構建非地面點的nDSM;然后通過自定義的形態(tài)學分割算子斷開樹木與建筑物之間可能的連接,利用區(qū)域生成法獲得建筑物和植被的面域;最后根據(jù)點云的大坡度密度性質(zhì)來快速區(qū)分建筑物的面域。
采用的機載LiDAR 點云來源于德國攝影測量、遙感與地理信息協(xié)會(DGPF)提供的專門用于3D建筑物重建和城市地物分類的ISPRS 測試工程數(shù)據(jù)。試驗數(shù)據(jù)由Leica ALS50 掃描儀獲取,相對航高為500 m,點云平均密度為9.2個/m2,點云包括多回波和反射強度信息。試驗區(qū)域范圍是德國Vaihingen 市,包含建筑物、低植被(草地)、樹木、汽車、道路及人工地面等地類,有輕微的地形起伏,為了提高算法的執(zhí)行效率,本文選擇了3個區(qū)域進行測試。區(qū)域1、區(qū)域2 和區(qū)域3 的面積分別為23 552.19 m2,30 145.28 m2和30 285.94 m2;區(qū)域1 的建筑物以尖頂房為主,大小各異,形狀復雜,植被較少但與建筑物連接緊密;區(qū)域2 的地形起伏較大,有4 座大型平頂樓房,樓房形狀復雜,高度呈階梯狀分布,茂密高大的植被包圍了全部建筑物;區(qū)域3 以尖頂房為主,大小不一,植被繁多且以中低高度為主。3個區(qū)域的數(shù)據(jù)都很復雜,適合檢測建筑物面域提取算法的效果。試驗區(qū)域分布如圖1 所示。
圖1 3個測試區(qū)域的位置分布Fig.1 Positions of three test areas
為了更好地識別建筑物,需要從非地面點云中分離出地面點云。首先,由于激光的多路徑效應、飛鳥、空中漂浮物和激光掃描儀的系統(tǒng)誤差可能造成不必要的高位或低位的粗差點云,所以通過生成LiDAR點云的高程直方圖的方法來確定閾值,剔除低位和高位粗差點云的噪聲點云。然后,采取漸進式TIN 加密法對地面點云進行分類[11]。具體算法:將量測的最大建筑物尺寸作為最大格網(wǎng)的邊長,從每個格網(wǎng)中選擇最低的高程點作為種子點生成初始地表TIN,以待判點到三角面片的迭代距離d 和與結(jié)點的迭代角度α 作為判斷依據(jù),對于滿足一定的最大距離閾值dmax和最大角度閾值αmax的點添加到TIN 中,組成新的TIN;對于超過dmax和αmax的地形突變點(如懸崖和建筑物邊緣),將最近的三角形結(jié)點作為對稱中心,生成一個虛擬的鏡像點,計算該點到最近三角形面片的距離,如果該距離大于dmax,則剔除該點,否則保留為地面點類并添加到TIN 中。程序迭代運行直到?jīng)]有更多的點被添加到TIN 中為止,最后TIN 中的點即為初步分類的地面點類。每次迭代用到的dmax和αmax分別設定為地形高差和地形坡度的中位數(shù)。
漸進式TIN 加密法同時受迭代距離和迭代角度的共同約束,因此會拒絕一些迭代角度較大而迭代距離較小的地面點,并認為它們是非地面點。為了分類出所有的地面點,用初步分類的地面點構建DEM,如果點云與DEM 的高程差小于LiDAR 點云的最大高程誤差0.3 m 時,就將該點云歸類為地面點類,最后將新分類的地面點重新生成DEM。
LiDAR 點云的地面點濾波后,剩余的非地面點類包含草地、灌木、樹木、房屋和汽車等??紤]到城市房屋的高度約3 m,則將與DEM 高度差小于3 m的點簡單地歸類為低植被點,剩余的LiDAR點云主要包括樹木和房屋點類。此時,這些點云中存在著少量的孤立點云,它們組成樹木或屋頂面的幾率很小。為了剔除這些孤立點云,遍歷每個點云周圍其他點云的個數(shù),如果點數(shù)小于4,則認為是孤立點,予以刪除。
此時剩余的未分類點云用來生成二值化格網(wǎng),只要空格網(wǎng)最鄰近的點和落入格網(wǎng)內(nèi)的點為地面點、低植被點或者孤立點中的任意一種,則格網(wǎng)被賦值為0;否則賦值為1。同時,對應于同樣大小的格網(wǎng)生成二維高程數(shù)組,如果格網(wǎng)為非空格網(wǎng),取落入格網(wǎng)內(nèi)的最高點作為格網(wǎng)高程;如果格網(wǎng)為空格網(wǎng),則取距離格網(wǎng)最近點的高程作為格網(wǎng)高程。
建筑物和樹木生成的二值化格網(wǎng)仍然存在大量建筑物和植被粘連在一起的現(xiàn)象,為進一步分離二者,需要判斷非0 中心格網(wǎng)周圍8 鄰域格網(wǎng)的分布情況。圖2(a)展示的是水平、垂直和對角線值分布均為0-1-0 的情形,當它們兩側(cè)的空白格網(wǎng)分別不同時為0 時,則中心格網(wǎng)賦值0。圖2(b)(c)展示的是非0 中心格網(wǎng)周邊分布值為0-1-0-1 的情形,如果其余5個空白格網(wǎng)中有任意一個格網(wǎng)為0,則中心格網(wǎng)賦值0。
圖2 建筑物與植被的分割算子Fig.2 Segmentation operators of building and vegetation
為了保證建筑物與樹木之間最大程度的分離,參照LiDAR 點云之間的高差,采用圖像處理中的區(qū)域生長算法繼續(xù)分割房屋與樹木。首先,搜索值為1 的格網(wǎng)作為種子點,將種子點的8 鄰域范圍內(nèi)點的高程與種子點的高程進行對比,如果小于一定的高差閾值,則將該點包含到種子點區(qū)域,新搜索的點重新作為種子點重復以上搜索及合并的過程,直到?jīng)]有更多的點可以合并為止。高差閾值dh可由房頂最大坡度smax和點云平均間距c 計算得到,有時點云平均間距小于實際間距,所以擴大2 倍以保證得到更準確的高差;計算的高差再與LiDAR 點云的高程誤差dh0比較,取其中較大者作為高差閾值,即
式中dh0一般在0.15~0.3 m 之間。
在建筑物和植被面域被分割出來后,需要判斷各個面域的歸屬類別。統(tǒng)計各個面域的面積,如果小于預先設置的最小面積閾值,就刪除此塊面域。另外,引入LiDAR 點云大坡度密度數(shù)來區(qū)分建筑物面域和植被面域。因為規(guī)則建筑物的大坡度值多數(shù)分布在建筑物的邊緣處,而植被點云的高程差異相對較大,它的大坡度值比建筑物來更容易向中心擴展;因此,LiDAR 點云的大坡度密度數(shù)可以用來分辨二者。首先,依據(jù)格網(wǎng)高程計算LiDAR 點云的坡度值g,即
式中:d z/d x 表示坡度在水平方向上的變化率;d z/d y表示坡度在垂直方向上的變化率。
然后,統(tǒng)計各個面域內(nèi)坡度值大于61°的坡度數(shù),用坡度數(shù)除以此面域的面積可得大坡度密度數(shù)ρs。對于ρs小于密度閾值的面域,就識別為建筑物面域;否則認為是植被面域。經(jīng)試驗獲知,ρs處于0.5~2.5個/m2的范圍內(nèi)。
識別出的建筑物面域中仍然存在由采集過程產(chǎn)生的少量數(shù)據(jù)空白,為了填充這些小空白區(qū)域,需要預先量測出最大數(shù)據(jù)空白的面積,以這個填充面積除以單位格網(wǎng)的面積,取整后可得數(shù)學形態(tài)學閉運算所需要的窗口(結(jié)構元素)B 的大小。數(shù)學形態(tài)學閉運算具有填充比窗口小的缺口或孔洞,連通小間隔的間斷以及平滑邊界的功能,閉運算還可以部分彌補建筑物面域的缺失。閉運算是用結(jié)構元素B對初始提取的建筑物面域A 進行膨脹,緊接著再進行一次腐蝕的結(jié)果,其公式為
式中:A 為提取建筑物面域的集合;B 為結(jié)構元素;符號·表示閉運算;⊕表示膨脹運算;Θ 表示腐蝕運算。
測試算法時所用的3個區(qū)域參數(shù)和測試結(jié)果分別見表1 和圖3。
表1 建筑物面域提取參數(shù)Tab.1 Parameters for extracting building regions
圖3 建筑物面域提取結(jié)果Fig.3 Results of extracting building regions
為了評價建筑物提取的精度,使用以下評價因子來估計提取結(jié)果精度[12],即
式中:T 表示提取數(shù)據(jù)與參考數(shù)據(jù)相匹配的部分;F 表示未與參考數(shù)據(jù)匹配的提取數(shù)據(jù);N 表示未與提取數(shù)據(jù)匹配的參考數(shù)據(jù);C 表示正確提取的結(jié)果占參考數(shù)據(jù)的比率;Q 表示正確提取的結(jié)果占全部數(shù)據(jù)的比率;B 表示未匹配的提取結(jié)果占正確提取結(jié)果的比率;M 表示未匹配的參考數(shù)據(jù)占正確提取結(jié)果的比率。
結(jié)合點云的高程和反射強度信息,在高密度的點云上人工繪制參考建筑物區(qū)域的線性輪廓,再用本文算法自動提取的建筑物面域共同計算以上質(zhì)量評價因子,具體結(jié)果見表2。
表2 3個區(qū)域建筑物面域提取結(jié)果的精度Tab.2 Precision of extracting building regions from three test areas (%)
從表2 看出,3個區(qū)域的建筑物提取精度結(jié)果中,區(qū)域1 的Q 最好,B 最高,而M 最小;區(qū)域2 的Q 最差,B 一般,但M 最大;區(qū)域3 的Q 一般,B 最小,而M 居中。綜合來看,區(qū)域1 的提取效果最好,區(qū)域2 的提取效果最差,區(qū)域3 的提取效果一般。區(qū)域1 均存在多余提取和丟失建筑物的情況,區(qū)域2 和區(qū)域3 的多余提取較少而丟失率較大。結(jié)合圖3(b)(c)可以看出,主要是因為區(qū)域2 和區(qū)域3 中均存在未探測出的建筑物。經(jīng)查看原始點云和影像數(shù)據(jù)后發(fā)現(xiàn),有2 棟未探測出的建筑物的屋頂為高吸收率的材料(如瀝青),所以導致屋頂點云數(shù)量異常稀少,算法無法進行識別。假設不存在這2 棟高吸收率的房屋,則區(qū)域2 和區(qū)域3 的M 將會很小,質(zhì)量將有大的提升。區(qū)域1 的多余探測率較高,是因為其中大部分房屋被茂密的植被緊密包圍,造成算法無法準確區(qū)分植被與建筑物。
在有茂密植被覆蓋區(qū)域的復雜城市環(huán)境下,本文嘗試探測并提取了高度在3 m 以上的建筑物面域,選取的3個典型區(qū)域的提取精度都達到了91%以上,可以滿足自動識別建筑物的要求,但不能識別出覆蓋有高吸收率材料的屋頂。此外,本文未考慮高度在3 m 以下的點云少量小建筑物存在的情形。因此,在未來的工作中,可以結(jié)合高分辨率的多光譜影像來彌補點云缺失所造成的建筑物丟失的問題,測試算法對大面積城市區(qū)域的適應性,繼續(xù)研究3 m以下建筑物的識別算法,以求達到更準確地提取城市建筑物面域的目標。
志謝:感謝德國DGPF 為本文提供了試驗數(shù)據(jù)支持。
[1]Awrangjeb M,Zhang C,F(xiàn)raser C S.Improved building detection using texture information[C]//International Archives of Photogrammetry,Remote Sensing and Spatial Information Sciences,2011,38(3/W22):143-148.
[2]成曉倩,樊良新,趙紅強.基于圖像分割技術的城區(qū)機載Li-DAR 數(shù)據(jù)濾波方法[J].國土資源遙感,2012,24(3):29-32.Cheng X Q,F(xiàn)an L X,Zhao H Q.Filtering of airborne LiDAR data for cityscapes based on segmentation[J].Remote Sensing for Land and Resources,2012,24(3):29-32.
[3]Zhang K,Yan J,Chen S C.Automatic construction of building footprints from airborne LiDAR data[J].IEEE Transactions on Geoscience and Remote Sensing,2006,44(9):2523-2533.
[4]Rottensteiner F,Trinder J,Clode S,et al.Building detection by fusion of airborne laser scanner data and multi-spectral images:Performance evaluation and sensitivity analysis[J].ISPRS Journal of Photogrammetry and Remote Sensing,2007,62(2):135-149.
[5]Khoshelham K,Nedkov S,Nardinocchi C.A comparison of bayesian and evidence- based fusion methods for automated building detection in aerial data[C]//International Archives of the Photogrammetry,Remote Sensing and Spatial Information Sciences,2008,37(B7):1183-1188.
[6]Zhou Q Y,Neumann U.Fast and extensible building modeling from airborne LiDAR data[C]//Proceedings of the 16th ACM Sigspatial International Conference on Advances in Geographic Information Systems.New York,2008-11-05(07).
[7]Meng X,Wang L,Currit N.Morphology-based building detection from airborne LiDAR data[J].Photogrammetric Engineering and Remote Sensing,2009,75(4):437-442.
[8]Salah M,Trinder J,Shaker A.Evaluation of the self- organizingmap classifier for building detection from LiDAR data and multispectral aerial images[J].Journal of Spatial Science,2009,54(2):1-20.
[9]于海洋,程 鋼,張育民,等.基于LiDAR 和航空影像的地震災害倒塌建筑物信息提?。跩].國土資源遙感,2011,23(3):77-81.Yu H Y,Cheng G,Zhang Y M,et al.The detection of earthquake- caused collapsed building information from LiDAR data and aerophotograph[J].Remote Sensing for Land and Resources,2011,23(3):77-81.
[10]李 婷,詹慶明,喻 亮.基于地物特征提取的車載激光點云數(shù)據(jù)分類方法[J].國土資源遙感,2012,24(1):17-21.Li T,Zhan Q M,Yu L.A classification method for mobile laser scanning data based on object feature extraction[J].Remote Sensing for Land and Resources,2012,24(1):17-21.
[11]Axelsson P E.DEM generation from laser scanner data using adaptive TIN models[C]//International Archives of Photogrammetry,Remote Sensing and Spatial Information Sciences,2000,33(B4):110-117.
[12]Hermosilla T,Ruiz L A,Recio J A,et al.Evaluation of automatic building detection approaches combining high resolution images and LiDAR data[J].Remote Sensing,2011,3(6):1188-1210.