王麗英，夏玉紅，徐 艷，趙元丁

(1. 遼寧工程技術(shù)大學(xué),遼寧 阜新 123000； 2. 北京金景科技有限公司，北京 100094)

機(jī)載激光雷達(dá)(light detection and ranging，LiDAR)技術(shù)可快速獲取具有三維坐標(biāo)和一定屬性(強(qiáng)度等)的海量、不規(guī)則空間分布三維點(diǎn)云[1]，為三維地理空間信息獲取提供了全新的技術(shù)手段。獨(dú)特的優(yōu)勢(shì)使得點(diǎn)云數(shù)據(jù)需求量隨之劇增。然而點(diǎn)云中常存在局部數(shù)據(jù)缺失，產(chǎn)生缺失數(shù)據(jù)的原因包括：條帶不重疊；系統(tǒng)無法記錄數(shù)據(jù)；完全透射或吸收性地物導(dǎo)致無反射回波；遮擋。缺失數(shù)據(jù)在點(diǎn)云中呈現(xiàn)出空洞現(xiàn)象(如圖1所示)，因而也被稱作數(shù)據(jù)空洞。其存在導(dǎo)致相應(yīng)區(qū)域的地物細(xì)節(jié)不能體現(xiàn)，從而影響數(shù)據(jù)的完整性，進(jìn)而對(duì)后續(xù)點(diǎn)云數(shù)據(jù)的處理、分析和應(yīng)用產(chǎn)生消極影響。其中，前兩類缺失數(shù)據(jù)可以通過補(bǔ)飛的方式重新獲取；然而，后兩類缺失數(shù)據(jù)則是無法避免的。因此，需辨識(shí)并填補(bǔ)各類數(shù)據(jù)缺失。

圖1 數(shù)據(jù)空洞示例(不同灰度代表高程不同)

經(jīng)典的數(shù)據(jù)空洞填補(bǔ)方法包括最低點(diǎn)填補(bǔ)、鄰近填補(bǔ)及基于地形的填補(bǔ)。最低點(diǎn)填補(bǔ)[2-6]提取數(shù)據(jù)空洞的外邊界中高程最低值作為數(shù)據(jù)缺失區(qū)域的高程；鄰近填補(bǔ)[7-10]將與待填補(bǔ)點(diǎn)距離最小的數(shù)據(jù)點(diǎn)的高程作為待填補(bǔ)點(diǎn)的高程；基于地形的填充方法[11-15]針對(duì)點(diǎn)云數(shù)據(jù)柵格化過程中出現(xiàn)的空白區(qū)域，對(duì)其邊界進(jìn)行一維形態(tài)學(xué)濾波得到其邊界各點(diǎn)地形高，利用邊界點(diǎn)構(gòu)建不規(guī)則三角網(wǎng)進(jìn)行地形內(nèi)插填充。上述3種方法均利用數(shù)據(jù)空洞的外邊界點(diǎn)的高程信息填補(bǔ)數(shù)據(jù)空洞，未考慮數(shù)據(jù)空洞的影響區(qū)域。只有根據(jù)數(shù)據(jù)空洞邊界輪廓線確定其影響區(qū)域的大小，進(jìn)而用影響區(qū)域內(nèi)的點(diǎn)云構(gòu)造逼近的連續(xù)曲面才能實(shí)現(xiàn)填補(bǔ)數(shù)據(jù)與周圍數(shù)據(jù)的光順連接。因而，影響區(qū)域的確定是數(shù)據(jù)空洞填補(bǔ)的最重要因素。另外，針對(duì)同一數(shù)據(jù)空洞，采用不同的空間插值方法填補(bǔ)空洞的效率及精度均不同。因而選擇最優(yōu)的空間插值方法是數(shù)據(jù)空洞填補(bǔ)的另一重要因素。針對(duì)上述兩個(gè)因素，本文提出一種數(shù)據(jù)空洞辨識(shí)及填補(bǔ)方法。該方法基于空洞輪廓凸包等距擴(kuò)大原理確定數(shù)據(jù)空洞的有效影響區(qū)域；利用影響區(qū)域內(nèi)的數(shù)據(jù)點(diǎn)集，針對(duì)空洞的高程特性，研究不同空間插值算法補(bǔ)償空洞的精度，進(jìn)而選擇最優(yōu)的空間插值方法并補(bǔ)償空洞。

1 缺失數(shù)據(jù)填補(bǔ)原理

算法流程如圖2所示。首先，依據(jù)空洞幾何特征及其空間分布構(gòu)建數(shù)據(jù)空洞圖，檢測(cè)空洞的存在；其次，對(duì)其進(jìn)行連通區(qū)域標(biāo)記，將數(shù)據(jù)空洞分割為單一數(shù)據(jù)空洞；然后，針對(duì)大面積的單個(gè)數(shù)據(jù)空洞，提取其邊界輪廓并將其等距擴(kuò)大獲得數(shù)據(jù)空洞的有效影響區(qū)域；最后，利用影響區(qū)域內(nèi)點(diǎn)云，借助最優(yōu)空間插值算法完成缺失數(shù)據(jù)填補(bǔ)。

圖2 數(shù)據(jù)空洞填補(bǔ)流程

1.1 數(shù)據(jù)空洞辨識(shí)

用數(shù)據(jù)空洞圖表征數(shù)據(jù)空洞的存在。數(shù)據(jù)空洞圖為一幅二值圖像，可用具有一定分辨率的二維格網(wǎng)(格網(wǎng)邊長(zhǎng)取平均點(diǎn)間距)覆蓋整個(gè)場(chǎng)景目標(biāo)的xy平面，并定位其中空的格網(wǎng)單元獲得，如圖3所示。其中，空洞為黑色像素，像素值為0；非空洞為白色像素，像素值為1。

圖3 數(shù)據(jù)空洞

1.2 數(shù)據(jù)空洞分割

對(duì)數(shù)據(jù)空洞圖中的0值像素進(jìn)行8連通區(qū)域標(biāo)記，將數(shù)據(jù)空洞圖分割并標(biāo)記為單一數(shù)據(jù)空洞。進(jìn)而統(tǒng)計(jì)各個(gè)空洞的面積，若其面積小于給定閾值Ta，則取該空洞的最鄰近點(diǎn)的高程作為該空洞內(nèi)各個(gè)數(shù)據(jù)點(diǎn)的高程；否則，需進(jìn)入下一步確定該空洞的有效影響區(qū)域，并利用影響區(qū)域內(nèi)點(diǎn)云，借助最優(yōu)空間插值算法完成缺失數(shù)據(jù)填補(bǔ)。

1.3 數(shù)據(jù)空洞的有效影響區(qū)域確定

在數(shù)據(jù)空洞的填補(bǔ)處理中，為了保證填補(bǔ)數(shù)據(jù)與周圍點(diǎn)云實(shí)現(xiàn)最佳融合，填補(bǔ)數(shù)據(jù)需利用數(shù)據(jù)空洞的周圍一定區(qū)域內(nèi)點(diǎn)云數(shù)據(jù)，運(yùn)用空間插值理論確定。其中，周圍數(shù)據(jù)的范圍及空間插值方法直接影響空洞填補(bǔ)的效率及精度。本文首先基于數(shù)據(jù)空洞的輪廓向外等距擴(kuò)大確定空洞的有效影響區(qū)域，進(jìn)而用影響區(qū)域內(nèi)的點(diǎn)云構(gòu)造逼近的連續(xù)曲面以實(shí)現(xiàn)填補(bǔ)數(shù)據(jù)與周圍數(shù)據(jù)的最佳融合，詳細(xì)方案如下：

(1) 提取數(shù)據(jù)空洞的外邊界輪廓，搜尋位于其上的數(shù)據(jù)點(diǎn)作為其外邊界輪廓點(diǎn)。

(3) 計(jì)算Bi到平面α的距離方差E和Bi投影到平面α所得邊界輪廓曲線(記作K)上的數(shù)據(jù)點(diǎn)的曲率均值k′。

(1)

式中，di為Bi到α的距離；X0為距離di的均值；m為K中的數(shù)據(jù)點(diǎn)數(shù)量；κi為數(shù)據(jù)點(diǎn)i處的曲線曲率。

(4) 根據(jù)E和k′將K向外等距得到擴(kuò)展的邊界輪廓(記作K′，K′=cK)，向外等距的距離c即最佳影響區(qū)域。其中，c=δE+εk′,δ+ε=1。如圖4所示。

圖4 影響區(qū)域

1.4 最優(yōu)空間插值方案確定及數(shù)據(jù)空洞填補(bǔ)

針對(duì)同一數(shù)據(jù)空洞，利用影響區(qū)域內(nèi)的數(shù)據(jù)點(diǎn)集，研究不同空間插值算法填補(bǔ)空洞的精度，進(jìn)而選擇最優(yōu)的空間插值方法并填補(bǔ)空洞。

2 試驗(yàn)數(shù)據(jù)

采用剪切自國(guó)際攝影測(cè)量與遙感協(xié)會(huì)(International Society for Photogrammetry and Remote Sensing，ISPRS)第三工作組(http:∥www.itc.nl/isprswg III-3/filtertest/)提供的包含不同的數(shù)據(jù)空洞類型的實(shí)測(cè)點(diǎn)云數(shù)據(jù)，如圖5所示。其中，F(xiàn)Site5如圖1所示。數(shù)據(jù)由Optech ALTM系統(tǒng)獲取，圖5、圖1中各個(gè)數(shù)據(jù)的激光腳點(diǎn)平均點(diǎn)密度分別約為0.75、0.14個(gè)/m2，回波次數(shù)為2。

圖5 試驗(yàn)數(shù)據(jù)

3 試 驗(yàn)

本文算法在CPU Core(TM)i5-2400 3.10 GHz、內(nèi)存4 GB、Windows 7旗艦版系統(tǒng)上使用Matlab 7.11.0平臺(tái)編程實(shí)現(xiàn)。

3.1 影響區(qū)域的有效性

試驗(yàn)旨在驗(yàn)證提出的數(shù)據(jù)空洞影響區(qū)域的確定方案的有效性。從CSite4中剪切一塊點(diǎn)云數(shù)據(jù)，如圖5(c)中區(qū)域所示。上述剪切數(shù)據(jù)可用作數(shù)據(jù)空洞的標(biāo)準(zhǔn)數(shù)據(jù)。試驗(yàn)針對(duì)剪切后區(qū)域(即數(shù)據(jù)空洞)。首先利用本文提出的方案確定數(shù)據(jù)空洞的影響區(qū)域；然后利用某一空間插值算法(反距離插值)填補(bǔ)空洞內(nèi)數(shù)據(jù)點(diǎn)(填補(bǔ)的數(shù)據(jù)點(diǎn)的平面坐標(biāo)和標(biāo)準(zhǔn)數(shù)據(jù)中數(shù)據(jù)點(diǎn)的平面坐標(biāo)一一對(duì)應(yīng)，但各個(gè)數(shù)據(jù)點(diǎn)的高程不同)；最后對(duì)比填補(bǔ)結(jié)果和標(biāo)準(zhǔn)數(shù)據(jù)，計(jì)算空洞填補(bǔ)結(jié)果的中誤差(見式(1))，見表1中第1列所示。同時(shí)，為了驗(yàn)證上述影響區(qū)域方案的有效性，將影響區(qū)域等距擴(kuò)大1倍、2倍，或縮小50%、75%，得到不同的影響區(qū)域(如圖6所示)，并基于上述影響區(qū)域結(jié)果計(jì)算空洞填補(bǔ)結(jié)果中誤差，見表1中第2～5列。

(2)

圖6 不同的影響區(qū)域及其內(nèi)的點(diǎn)云

本文確定的影響區(qū)域影響區(qū)域等距縮小50%影響區(qū)域等距縮小75%影響區(qū)域等距擴(kuò)大1倍影響區(qū)域等距擴(kuò)大2倍0.0641.7310.9780.1450.132

由表1可知，本文的影響區(qū)域確定方案對(duì)應(yīng)的空洞填補(bǔ)結(jié)果的中誤差僅為0.064。在相同條件下，對(duì)比其他影響區(qū)域，本文的影響區(qū)域確定方案對(duì)應(yīng)的空洞填補(bǔ)結(jié)果的中誤差最小，從而驗(yàn)證本文提出的數(shù)據(jù)空洞的影響區(qū)域方案的有效性。

3.2 不同數(shù)據(jù)空洞類型的最優(yōu)空間插值方法確定

對(duì)于不同的數(shù)據(jù)空洞類型，基于數(shù)據(jù)空洞的影響區(qū)域內(nèi)的數(shù)據(jù)點(diǎn)，分別采用經(jīng)典的克里金插值、線性插值、最臨近插值和反距離插值等空間插值方法對(duì)數(shù)據(jù)空洞進(jìn)行填補(bǔ)，研究不同空間插值算法填補(bǔ)的數(shù)據(jù)空洞的精度，進(jìn)而選擇最優(yōu)的空間插值方法填補(bǔ)空洞。其中，填補(bǔ)結(jié)果的精度評(píng)價(jià)采用了交叉驗(yàn)證的方法(由于沒有數(shù)據(jù)空洞區(qū)域的標(biāo)準(zhǔn)數(shù)據(jù)，無法直接對(duì)比)，即首先填補(bǔ)數(shù)據(jù)空洞，然后利用填補(bǔ)數(shù)據(jù)反向插值影響區(qū)域內(nèi)的激光點(diǎn)的高程，進(jìn)而對(duì)比插值結(jié)果和影響區(qū)域內(nèi)已知數(shù)據(jù)，并采用中誤差指標(biāo)進(jìn)行定量評(píng)價(jià)。

3.2.1 無反射回波類型的空洞的最優(yōu)填補(bǔ)方案確定

圖7 有效影響區(qū)域及其內(nèi)點(diǎn)云頂視圖

克里金插值反距離插值最鄰近插值線性插值0.0640.0720.0080.213

由表2可知，從中誤差指標(biāo)來看，最鄰近插值對(duì)應(yīng)的填補(bǔ)結(jié)果中誤差最小，其為無反射回波類型的數(shù)據(jù)空洞填補(bǔ)的最優(yōu)空間插值方案?；谧钹徑逯档目斩刺钛a(bǔ)結(jié)果如圖8所示。

3.2.2 目標(biāo)遮擋類型的空洞的最優(yōu)填補(bǔ)方案確定

試驗(yàn)數(shù)據(jù)如圖5(b)所示。由于目標(biāo)遮擋類的數(shù)據(jù)空洞的影響區(qū)域內(nèi)必然包含部分高大目標(biāo)點(diǎn)云，而基于影響區(qū)域內(nèi)高大目標(biāo)點(diǎn)云填補(bǔ)數(shù)據(jù)空洞必然造成較大誤差。因此，可首先統(tǒng)計(jì)影響區(qū)域內(nèi)點(diǎn)云的高程特性，若出現(xiàn)兩個(gè)正態(tài)分布，則僅取影響區(qū)域內(nèi)高程值小于谷底高程的點(diǎn)云填補(bǔ)數(shù)據(jù)空洞?；诓煌目臻g插值算法影響區(qū)域內(nèi)的點(diǎn)云填補(bǔ)數(shù)據(jù)空洞，并采用交叉驗(yàn)證的方法進(jìn)行精度評(píng)價(jià)，見表3。

圖8 最鄰近插值的填補(bǔ)結(jié)果

克里金插值反距離插值最鄰近插值線性插值0.0980.1310.0410.439

由表3可知，從中誤差指標(biāo)來看，最鄰近插值對(duì)應(yīng)的填補(bǔ)結(jié)果中誤差最小，其為目標(biāo)遮擋類型的數(shù)據(jù)空洞填補(bǔ)的最優(yōu)空間插值方案。

3.2.3 航帶不重疊類型的空洞的最優(yōu)填補(bǔ)方案確定

表4 不同空間插值算法對(duì)應(yīng)的填補(bǔ)結(jié)果中誤差

由表4可知，從平均絕對(duì)誤差指標(biāo)來看，克里金插值對(duì)應(yīng)的填補(bǔ)結(jié)果中誤差最小。上述結(jié)果表明：克里金插值方法為目標(biāo)遮擋類型的數(shù)據(jù)空洞填補(bǔ)的最優(yōu)空間插值方案。

4 結(jié) 語

本文提出了一種缺失數(shù)據(jù)填補(bǔ)方法。該方法利用二值數(shù)據(jù)空洞圖檢測(cè)缺失數(shù)據(jù)的存在，進(jìn)而利用連通區(qū)域標(biāo)記算法將數(shù)據(jù)空洞圖分割為單個(gè)數(shù)據(jù)空洞，然后對(duì)大面積的數(shù)據(jù)空洞的外邊界輪廓進(jìn)行等距擴(kuò)大獲得數(shù)據(jù)空洞的有效影響區(qū)域，最后研究不同的缺失數(shù)據(jù)類型的最優(yōu)空間插值方案，并基于影響區(qū)域內(nèi)的點(diǎn)云，利用最優(yōu)空間插值方案對(duì)不同類型的缺失數(shù)據(jù)進(jìn)行填補(bǔ)。創(chuàng)新工作主要圍繞數(shù)據(jù)空洞的有效影響區(qū)域確定和不同數(shù)據(jù)空洞類型的最優(yōu)空間插值方法兩個(gè)方面。試驗(yàn)基于ISPRS提供的包含了不同類型的空洞的點(diǎn)云數(shù)據(jù)測(cè)試了影響區(qū)域的有效性、不同類型的缺失數(shù)據(jù)的最優(yōu)空間插值方案，從而驗(yàn)證了提出方法的有效性。本文研究的不足之處在于：對(duì)于不同的缺失數(shù)據(jù)類型的最優(yōu)空間插值方案的研究?jī)H采用了4種經(jīng)典的空間插值算法，而未考慮其他空間插值算法；填補(bǔ)結(jié)果不但與數(shù)據(jù)空洞類型有關(guān)，還與數(shù)據(jù)空洞的范圍大小密切相關(guān)，但本文未考慮空洞范圍對(duì)填補(bǔ)結(jié)果的影響。在以后的研究中，應(yīng)考慮上述兩方面及填補(bǔ)效率的問題。