陳 馳,楊必勝,田 茂,李健平,鄒響紅,吳唯同,宋易恒

1. 武漢大學(xué)測(cè)繪遙感信息工程國家重點(diǎn)實(shí)驗(yàn)室，湖北 武漢 430079; 2. 武漢大學(xué)時(shí)空數(shù)據(jù)智能獲取技術(shù)與應(yīng)用教育部工程研究中心，湖北 武漢 430079

車載移動(dòng)測(cè)量系統(tǒng)實(shí)現(xiàn)了地基側(cè)面視角激光點(diǎn)云與影像數(shù)據(jù)的覆蓋，可高效采集被測(cè)環(huán)境中高分辨率立面激光雷達(dá)(light detection and ranging,LiDAR)點(diǎn)云與全景影像等數(shù)據(jù)，是一種新型的攝影測(cè)量與遙感手段。車載移動(dòng)測(cè)量系統(tǒng)(mobile mapping system，MMS)提供的高分辨率LiDAR點(diǎn)云、全景影像等多種對(duì)地觀測(cè)數(shù)據(jù)，可應(yīng)用于條帶地形測(cè)圖、高速公路資產(chǎn)管理、影像城市構(gòu)建等領(lǐng)域。多源數(shù)據(jù)融合是車載MMS應(yīng)用發(fā)展的趨勢(shì)。車載MMS多源數(shù)據(jù)配準(zhǔn)則是融合應(yīng)用需要解決的首要問題。

全景影像與車載激光點(diǎn)云是典型的車載MMS測(cè)繪數(shù)據(jù)，其配準(zhǔn)方法主要分為4類：

(1) 基于多傳感器標(biāo)定的半自動(dòng)方法，即直接使用位置與姿態(tài)系統(tǒng)(position and orientation system,POS)輸出方位元素值整體補(bǔ)償傳感器平臺(tái)標(biāo)定參數(shù)的方法進(jìn)行車載激光點(diǎn)云和影像的配準(zhǔn)[1-2]，該類方法通常需要使用標(biāo)定場實(shí)現(xiàn)標(biāo)定。研究表明由于配準(zhǔn)誤差來源的多樣性與硬件同步與幾何標(biāo)定的局限性，使用硬件同步與幾何檢校的方法進(jìn)行數(shù)據(jù)配準(zhǔn)后，配準(zhǔn)差依然存在[3-6]。同時(shí)，由于不同時(shí)段GPS信號(hào)質(zhì)量與IMU漂移各不相同，該類方法亦不能解決此類配準(zhǔn)誤差。

(2) 基于幾何配準(zhǔn)基元匹配的2D(全景影像)-3D(激光點(diǎn)云)自動(dòng)方法。該類方法在LiDAR點(diǎn)云與影像中提取共軛配準(zhǔn)基元匹配對(duì)，通常利用已有的建筑物線框模型作為先驗(yàn)知識(shí)，通過線框輪廓與激光點(diǎn)云輪廓的相似性，進(jìn)行兩者間的配準(zhǔn)[3,7-8]。如何穩(wěn)健實(shí)現(xiàn)共軛配準(zhǔn)基元的自動(dòng)提取與匹配是學(xué)術(shù)界尚未完全解決的研究問題。

(3) 基于互信息最大化匹配的2D-3D自動(dòng)方法。該類方法通常將點(diǎn)云數(shù)據(jù)依據(jù)全景影像成像模型成像，計(jì)算所得的點(diǎn)云強(qiáng)度/深度影像與全景影像的互信息測(cè)度[9]，使用Nelder-Mead方法[10]求解最高互信息解，實(shí)現(xiàn)配準(zhǔn)[4]，但需要預(yù)先對(duì)激光點(diǎn)云強(qiáng)度進(jìn)行標(biāo)定。

(4) 基于多視立體匹配點(diǎn)云與激光點(diǎn)云3D(影像點(diǎn)云)-3D(激光點(diǎn)云)的配準(zhǔn)方法[5]。該類方法利用POS系統(tǒng)的輸出值作為初始配準(zhǔn)，使用非剛性最近點(diǎn)迭代方法實(shí)現(xiàn)二者的配準(zhǔn)，但是迭代最鄰近點(diǎn)(iterative closet point，ICP)方法對(duì)初始轉(zhuǎn)換參數(shù)的近似值要求較高[11]，在POS提供的初始值偏差較大的情況下，難以保證該方法的收斂，算法適應(yīng)性有限。

為解決序列全景影像與激光點(diǎn)云的穩(wěn)健配準(zhǔn)問題，本文研究了一種2D(全景影像)-3D(激光點(diǎn)云)幾何基元對(duì)配準(zhǔn)與3D(影像點(diǎn)云)-3D(激光點(diǎn)云)配準(zhǔn)相結(jié)合的自動(dòng)化配準(zhǔn)方法：首先使用提取到的配準(zhǔn)基元對(duì)進(jìn)行2D-3D粗配準(zhǔn)，進(jìn)而使用全景影像虛擬影像生成的多視立體匹配密集影像點(diǎn)云與激光點(diǎn)云進(jìn)行3D-3D精化配準(zhǔn)結(jié)果，消除配準(zhǔn)基元提取與匹配誤差對(duì)配準(zhǔn)結(jié)果的影響，從而保證算法穩(wěn)健高精度收斂。

1 配準(zhǔn)方法

本文提出的車載MMS激光點(diǎn)云與序列全景影像的自動(dòng)配準(zhǔn)流程如圖1所示。該方法主要分為配準(zhǔn)基元提取、粗配準(zhǔn)、精配準(zhǔn)3個(gè)步驟：

(1) 配準(zhǔn)基元提取。選取天際線(建筑物上邊線)作為配準(zhǔn)基元，在LiDAR點(diǎn)云中采用層次化城市場景目標(biāo)提取方法[12]提取建筑物立面并規(guī)則化生成點(diǎn)云天際線矢量，在序列全景影像中，通過虛擬成像方法，將球面展開成像模型的全景影像轉(zhuǎn)換為框幅式影像，并在該虛擬影像上結(jié)合GrabCut算法[13]與弧到點(diǎn)距離累計(jì)(chord-to-point distance accumulation，CPDA)角點(diǎn)檢測(cè)算法[14]生成天際線矢量。

(2) 基于配準(zhǔn)基元圖匹配的2D-3D粗配準(zhǔn)。依據(jù)提取到的天際線矢量，構(gòu)建配準(zhǔn)基元圖，并通過配準(zhǔn)基元圖編輯距離最小化準(zhǔn)則進(jìn)行匹配，組成共軛配準(zhǔn)基元對(duì)，自生成的共軛配準(zhǔn)基元對(duì)，解算2D-3D粗配準(zhǔn)模型，通過投票方法，獲得攝影和測(cè)量坐標(biāo)系與LiDAR參考坐標(biāo)系之間的初始轉(zhuǎn)換關(guān)系。

(3) 基于迭代最鄰近點(diǎn)方法的3D-3D精確配準(zhǔn)。自序列全景影像虛擬成像生成序列虛擬影像，使用運(yùn)動(dòng)結(jié)構(gòu)恢復(fù)(structure from motion，SfM)方法恢復(fù)序列虛擬影像在攝影測(cè)量坐標(biāo)系中的外方位元素，并通過多視立體匹配的方法由序列影像生成多視立體匹配(multi-view stereo，MVS)影像密集點(diǎn)云，將粗配準(zhǔn)解算獲得空間坐標(biāo)轉(zhuǎn)換關(guān)系作為初始值，使用ICP算法變種，實(shí)現(xiàn)全景影像MVS密集點(diǎn)云與LiDAR點(diǎn)云之間的最優(yōu)配準(zhǔn)，獲得序列影像的精確配準(zhǔn)參數(shù)。

圖1 車載MMS激光點(diǎn)云與序列全景影像自動(dòng)配準(zhǔn)方法Fig.1 Flow chart of automatic registration of vehicles-borne mobile mapping system panoramic images sequence and laser point cloud

1.1 配準(zhǔn)基元提取

基于特征的地基視角影像與激光點(diǎn)云數(shù)據(jù)配準(zhǔn)方法中，多采用線段[15-17]、建筑物外框[18]、天際線[19]以及Sift描述子[20-21]作為特征基元。城市場景中有豐富的建筑物立面，且其提取方法有大量的成熟研究[22]。天際線指建筑物頂部邊線，建筑物頂部邊線與天空相接，反差較大，通常由于位置較高，不易被環(huán)境地物遮擋，出現(xiàn)頻率與建筑物立面出現(xiàn)頻率一致，較為適合作為配準(zhǔn)基元特征。因此，本文采用天際線作為LiDAR點(diǎn)云與全景影像之間的配準(zhǔn)基元。

車載LiDAR點(diǎn)云數(shù)據(jù)中天際線提取，等價(jià)于建筑物立面的提取。本文采用層次化城市場景目標(biāo)提取方法[12]，進(jìn)行車載LiDAR點(diǎn)云數(shù)據(jù)中建筑物立面的提取。對(duì)于提取到的建筑物立面點(diǎn)云塊，使用RMBR算法[23]提取其輪廓并進(jìn)行規(guī)則化，將提取到的外部多邊形作為候選配準(zhǔn)基元。全景影像成像模型為球狀(見圖2)。球面為不可展開平面，在進(jìn)行球面展開生成二維光譜信息存儲(chǔ)矩陣(全景影像)的過程中，產(chǎn)生了較大的變形。全景影像中的地物幾何結(jié)構(gòu)受此影響存在較大的變形扭曲。本文提出一種全景影像虛擬成像天際線提取方法，該方法首先使用虛擬成像方法將全景影像轉(zhuǎn)化為框幅式虛擬影像，繼而進(jìn)行直方圖閾值分割確定大致天際線位置，繼而使用GrabCut分割方法[13]實(shí)現(xiàn)天際線位置的優(yōu)化，最后利用CPDA角點(diǎn)檢測(cè)算法生成天際線矢量。令物方某點(diǎn)在虛擬成像平面上的像空間坐標(biāo)系坐標(biāo)為[xv,yv,-f]T，在全景球模型上的成像點(diǎn)為[xp,yp,zp]T，對(duì)應(yīng)展開模型上的像素位置為(r,c)，則全景球上的像點(diǎn)與虛擬成像像點(diǎn)之間的共線關(guān)系可以表示為式(1)

(1)

式中，h、w分別為全景影像高與寬；r、c分別為當(dāng)前像素在全景影像上的行列號(hào)。結(jié)合預(yù)定義的虛擬相機(jī)內(nèi)方位元素、像素大小、虛擬面陣大小，通過對(duì)式(1)的反解，即可得出任意焦距、方位的普通鏡頭虛擬成像影像。天空由于反射強(qiáng)度較高，與圖像主體的城市環(huán)境有較大的灰度差異，且與建筑物相接處有明顯的邊緣特征(圖3(a)、(b))。首先，對(duì)圖像直方圖進(jìn)行統(tǒng)計(jì)，選擇灰度中間值(128)與灰度高峰值的中間值作為分割閾值(圖3(c));然后，使用該分割閾值對(duì)圖像進(jìn)行分割，將高于閾值的區(qū)域作為候選天空像素區(qū)域，低于閾值的部分作為候選地物區(qū)域(圖3(d))?；谥狈綀D閾值分割即可確定天際線初始位置，繼而使用直方圖閾值分割天際線初始位置作為GrabCut算法[13]的初始分割。優(yōu)化地物前景與天空背景的分割結(jié)果，提高天際線提取質(zhì)量具體步驟如下：①依據(jù)直方圖閾值分割結(jié)果與天空背景位于圖像上部且連續(xù)分布的先驗(yàn)知識(shí)，指定GrabCut初始分割，將圖像定義為如圖3(e)的3個(gè)區(qū)域即，天空區(qū)域(頂部)、可能是天空區(qū)域(中部)可能是地物區(qū)域(底部)；②將此初始分割作代入到GrabCut算法中進(jìn)行迭代求解直到算法收斂(圖3(f))；③依據(jù)分割獲得的二值影像，進(jìn)行輪廓提取，獲得天際線。然后采用CPDA方法[14,24]實(shí)現(xiàn)輪廓線上的角點(diǎn)檢測(cè)，生成天際線矢量數(shù)據(jù)(圖3(g))。在全景影像與激光點(diǎn)云中提取到的天際線矢量將用于構(gòu)建配準(zhǔn)基元圖，繼而參入2D全景影像-3D激光點(diǎn)云粗配準(zhǔn)模型參數(shù)的解算。

圖2 全景影像虛擬成像方法Fig.2 Virtual imaging of panoramic image

1.2 基于配準(zhǔn)基元圖匹配的2D-3D粗配準(zhǔn)

采用2D-3D共線配準(zhǔn)模型描述全景影像與LiDAR點(diǎn)云之間的幾何關(guān)系。令虛擬影像上的天際線角點(diǎn)配準(zhǔn)基元的像點(diǎn)為[xv,yv,-f]T，其對(duì)應(yīng)的全景影像的像點(diǎn)坐標(biāo)為[xp,yp,zp]T，對(duì)應(yīng)的建筑物角點(diǎn)在LiDAR數(shù)據(jù)中的坐標(biāo)為(xw,yw,zw)，則兩角點(diǎn)之間的共線關(guān)系可表述為式(2)

(2)

式中，Rv為虛擬成像過程中的相機(jī)旋轉(zhuǎn)矩陣；Rf為旋轉(zhuǎn)平移解耦計(jì)算中的旋轉(zhuǎn)矩陣；Tf為相機(jī)平移參數(shù)，λ為比例參數(shù)；h、w分別為全景影像高與寬；r、c分別為當(dāng)前像素在全景影像上的行列號(hào)。將全景影像直接地理定向數(shù)據(jù)經(jīng)導(dǎo)航坐標(biāo)系到攝影測(cè)量坐標(biāo)系轉(zhuǎn)換得到全景影像外方位元素(Rf,Tf)初始值[25]。使用該值將提取到的LiDAR天際線角點(diǎn)反投影到全景影像面上，繼而提出一種最小圖像編輯距離(graph edit distance，GED)圖匹配算法，實(shí)現(xiàn)配準(zhǔn)基元匹配。對(duì)于車載MMS街景數(shù)據(jù)，選擇提取到的n個(gè)配準(zhǔn)基元角點(diǎn)組成圖節(jié)點(diǎn)集E。選取圖頂點(diǎn)集合E中的任意3個(gè)頂點(diǎn)組合為三角形，作為圖G的生成核，圖的邊連接規(guī)則定義為：對(duì)于核三角形頂點(diǎn)，直接連接構(gòu)成完全圖；對(duì)于非核三角形頂點(diǎn)，則計(jì)算該頂點(diǎn)到核三角形3個(gè)頂點(diǎn)的距離，將最短長度的邊作為連接邊。使用GED度量影像配準(zhǔn)基元圖Gimage=(Eimage,Vimage)與激光點(diǎn)云配準(zhǔn)基元圖Glas=(Elas,Vlas)在局部核三角形匹配對(duì)定義的圖轉(zhuǎn)換下的全局圖邊相似度，即將最優(yōu)圖匹配的問題轉(zhuǎn)化為最小圖編輯距離查找問題(圖4)。Gimage、Glas與之間的GED定義為式(3)

GED(Glas,Gimage,T)=

(3)

式中，cost(op,T)為將Gimage經(jīng)過圖編輯之后完全匹配Glas所需要的K步操作(添加、刪除、替換)所對(duì)應(yīng)的圖編輯代價(jià)函數(shù)，本文設(shè)置該代價(jià)函數(shù)為旋轉(zhuǎn)圖邊差；length(·)為圖邊長度累加函數(shù)；T為局部核三角形匹配對(duì)定義的圖轉(zhuǎn)換。首先，使用KD樹[26]對(duì)所有的核三角形內(nèi)角值三維點(diǎn)坐標(biāo)數(shù)據(jù)構(gòu)建索引，提高匹配效率；其次，遍歷核三角形匹配集計(jì)算各個(gè)核三角形匹配對(duì)應(yīng)圖轉(zhuǎn)換的GED，并進(jìn)行排序，將最小GED對(duì)應(yīng)的T作為最佳匹配，繼而采用旋轉(zhuǎn)平移解耦合計(jì)算方法[11]，計(jì)算基于2D-3D幾何特征對(duì)的粗配準(zhǔn)參數(shù)。

由于幾何配準(zhǔn)基元提取與匹配存在誤差，導(dǎo)致粗配準(zhǔn)結(jié)果精度與穩(wěn)健性受限。為解決該問題，本文將粗配準(zhǔn)空間坐標(biāo)轉(zhuǎn)換關(guān)系作為初始值，使用ICP算法變種，最優(yōu)配準(zhǔn)全景影像密集點(diǎn)云與LiDAR點(diǎn)云，間接獲得序列影像的精確配準(zhǔn)參數(shù)，實(shí)現(xiàn)數(shù)據(jù)穩(wěn)健高精度配準(zhǔn)。

圖3 全景影像虛擬影像幾何配準(zhǔn)基元(天際線)提取Fig.3 Geometric primitive (skyline) extraction on virtual images generated by panoramic images

圖4 基于圖匹配的幾何配準(zhǔn)基元對(duì)生成方法Fig.4 Conjugate geometric primitives generated by GED minimization

1.3 基于迭代最鄰近點(diǎn)方法的3D-3D精確配準(zhǔn)

將全景影像密集點(diǎn)云構(gòu)建的問題轉(zhuǎn)換為傳統(tǒng)短焦針孔相機(jī)成像數(shù)據(jù)的密集點(diǎn)云構(gòu)建的問題。采用增量光束法平差SfM算法實(shí)現(xiàn)虛擬相機(jī)外參數(shù)標(biāo)定[27]，然后使用Daisy算法[28]生成MVS影像點(diǎn)云，實(shí)現(xiàn)物方的重構(gòu)。繼而最小化MVS點(diǎn)云與LiDAR點(diǎn)云之間的距離實(shí)現(xiàn)配準(zhǔn)參數(shù)精化，其配準(zhǔn)幾何模型可表示為式(4)

Mlas=λRMMVS+T

(4)

精配準(zhǔn)包括精配準(zhǔn)模型初值穩(wěn)健估計(jì)與迭代最鄰近點(diǎn)模型參數(shù)優(yōu)化。令虛擬影像上一像點(diǎn)為m，其在MVS攝影測(cè)量坐標(biāo)系CMVS中對(duì)應(yīng)點(diǎn)為MMVS，在激光點(diǎn)云參考系Cw中對(duì)應(yīng)點(diǎn)為Mlas則其之間的共線性關(guān)系可表達(dá)為式(5)和式(6)

sdecouplem=A[Rdecouple|tdecouple]Mlas

(5)

sbundlem=A[Rbundle|tbundle]MMVS

(6)

式中，[Rbundle|tbundle]，[Rdecouple|tdecouple]分別為SfM與粗配準(zhǔn)過程中計(jì)算獲得的影像外參數(shù)矩陣，A為已知相機(jī)內(nèi)參數(shù)矩陣。聯(lián)立式(5)和式(6)可得式(7)

(7)

令Pdecouple(Xi,Yi,Zi)t，Pbundle(xi,yi,zi)t，為相機(jī)在Cw與CMVS中的坐標(biāo)重心化后的相機(jī)位子。則其比例與平移參數(shù)可通過式(8)、式(9)計(jì)算

(8)

T=Pdecouple-λRPbundle

(9)

對(duì)于多張配準(zhǔn)影像，采用投票聚類的方法實(shí)現(xiàn)參數(shù)的穩(wěn)健估計(jì)[29]，實(shí)現(xiàn)轉(zhuǎn)換參數(shù)的解算。采用剛體空間相似變換模型作為MVS與LiDAR點(diǎn)云之間的配準(zhǔn)模型。使用點(diǎn)到面距離差[30]作為點(diǎn)對(duì)誤差，則誤差方程可表示為式(10)

(10)

式中，mi、Mi分別為MVS與LiDAR點(diǎn)；wi為點(diǎn)對(duì)權(quán)；ηi為Mi點(diǎn)法向量；(λ3D-3D,R3D-3D,T3D-3D)表示需要求解的最小MVS與LiDAR點(diǎn)云點(diǎn)到面距離的空間相似變換。將粗配準(zhǔn)解算獲得的轉(zhuǎn)換參數(shù)作為ICP算法初始值，對(duì)式(10)進(jìn)行迭代解算，避免由于POS數(shù)據(jù)質(zhì)量、幾何基元配準(zhǔn)提取與匹配錯(cuò)誤導(dǎo)致ICP迭代錯(cuò)誤收斂。為避免MVS點(diǎn)云外點(diǎn)對(duì)ICP算法的影響，采用相對(duì)運(yùn)動(dòng)閾值算法[31]對(duì)限制點(diǎn)對(duì)搜索，保證穩(wěn)健性收斂。記最優(yōu)攝影測(cè)量與LiDAR參考系之間的坐標(biāo)轉(zhuǎn)換為(λ3D-3D,R3D-3D,T3D-3D)，則修正后的虛擬影像外方位元素(Rv-cam,Tv-cam)可表示為式(11)

(11)

全景影像在虛擬成像過程中與虛擬影像之間存在已知的空間旋轉(zhuǎn)，即虛擬成像旋轉(zhuǎn)矩陣Rv，虛擬影像與全景影像共成像中心，則精化后的全景影像外方位元素(Rpano,Tpano)可表示為式(12)

(12)

2 配準(zhǔn)試驗(yàn)與分析

為驗(yàn)證本文方法的有效性，采用兩種車載MMS系統(tǒng)采集的兩組LiDAR點(diǎn)云與序列全景影像數(shù)據(jù)集對(duì)提出的車載MMS點(diǎn)云與全景影像自動(dòng)配準(zhǔn)方法進(jìn)行驗(yàn)證。LiDAR點(diǎn)密度、全景影像分辨率等相關(guān)數(shù)據(jù)參數(shù)見表1。數(shù)據(jù)集1與數(shù)據(jù)集2中場景分別為典型城區(qū)道路與居民區(qū)，均為車載MMS數(shù)據(jù)采集中常見的場景，具有較好的代表性。數(shù)據(jù)集1、2中包含的激光點(diǎn)云數(shù)據(jù)與樣例全景影像見圖5。

表1 試驗(yàn)數(shù)據(jù)集描述

圖5 試驗(yàn)數(shù)據(jù)集Fig.5 Experiment datasets

2.1 配準(zhǔn)結(jié)果

將全景影像配準(zhǔn)后成像中心位置一定范圍內(nèi)(50 m)的激光點(diǎn)云按全景影像的球狀成像模型，反投影到全景影像上，形成疊加點(diǎn)云顯示的全景影像，從而定性判斷數(shù)據(jù)配準(zhǔn)精度。圖6為一組數(shù)據(jù)集1、2全景影像疊合激光點(diǎn)云顯示圖，其中圖6(a)、(b)各行分別為原始全景影像，使用直接地理定向數(shù)據(jù)、粗配準(zhǔn)影像外方位元素、精配準(zhǔn)影像外方位元素疊加點(diǎn)云數(shù)據(jù)配準(zhǔn)結(jié)果。由于配準(zhǔn)基元提取過程中的不完整以及定位不精確等因素，導(dǎo)致完成粗配準(zhǔn)后配準(zhǔn)誤差依然存在(圖6(a)、(b)第3行)，但相較使用直接地理定向數(shù)據(jù)進(jìn)行配準(zhǔn)的結(jié)果(圖6(a)、(b)第2行)，精度有所提升。

2.2 配準(zhǔn)誤差評(píng)定

為定量分析本章提出配準(zhǔn)算法配準(zhǔn)數(shù)據(jù)集1、2全景影像與LiDAR點(diǎn)云配準(zhǔn)獲得的配準(zhǔn)結(jié)果精度，人工在配準(zhǔn)后的全景影像與激光全景影像上選擇匹配對(duì)，并計(jì)算之間的像素差，作為配準(zhǔn)精度指標(biāo)。激光全景影像即是指以當(dāng)前全景影像位置為成像中心，切割出該成像中心一定范圍(50 m)內(nèi)的點(diǎn)云投影到全景球狀成像模型并展開為二維平面影像。在數(shù)據(jù)集1、2中各自隨機(jī)抽取20張配準(zhǔn)后的全景影像，生成每一張全景影像對(duì)應(yīng)的激光全景影像。在每對(duì)全景與激光全景影像上人工選取10個(gè)同名點(diǎn)，計(jì)算同名點(diǎn)之間的像素差的絕對(duì)值，并取標(biāo)記點(diǎn)對(duì)誤差平均值作為該幀全景影像配準(zhǔn)精度測(cè)度。表2為數(shù)據(jù)集1、2中隨機(jī)抽樣影像圖人工標(biāo)記點(diǎn)對(duì)誤差統(tǒng)計(jì)表。由表可知，在初始配準(zhǔn)誤差較大的條件下，經(jīng)由粗到精的配準(zhǔn)，配準(zhǔn)誤差逐漸減小，配準(zhǔn)后平均誤差在1.5像素左右。圖7為數(shù)據(jù)集1、2配準(zhǔn)過程中，抽樣的20張精度評(píng)定影像的配準(zhǔn)誤差柱狀圖，可見，使用兩步法配準(zhǔn)方法進(jìn)行配準(zhǔn)后，配準(zhǔn)誤差大幅下降，并且?guī)`差標(biāo)準(zhǔn)差亦達(dá)到較低水平，說明配準(zhǔn)結(jié)果精度較高，且具有較好的精度一致性。

表2 數(shù)據(jù)集1、2隨機(jī)人工標(biāo)記點(diǎn)對(duì)投影誤差

LiDAR點(diǎn)云數(shù)據(jù)與序列全景影像數(shù)據(jù)配準(zhǔn)完成后，對(duì)LiDAR點(diǎn)云數(shù)據(jù)與序列全景影像進(jìn)行融合。圖8為LiDAR點(diǎn)云數(shù)據(jù)與序列全景影像數(shù)據(jù)配準(zhǔn)完成后生成的彩色激光點(diǎn)云，由其房屋邊緣細(xì)節(jié)可知，數(shù)據(jù)配準(zhǔn)精度良好。

圖6 由粗到精的激光點(diǎn)云和全景影像自動(dòng)配準(zhǔn)Fig.6 Coarse to fine registration of panoramic images and laser point cloud

圖7 隨機(jī)抽樣影像配準(zhǔn)精度柱狀圖Fig.7 Residual projection error distribution of the random labeled points

圖8 建筑物立面彩色點(diǎn)云細(xì)節(jié)Fig.8 Color point cloud details on the facade of buildings

3 結(jié) 論

本文提出一種車載移動(dòng)測(cè)量系統(tǒng)(MMS)激光點(diǎn)云與序列全景影像自動(dòng)配準(zhǔn)方法。該方法首先使用幾何基元對(duì)實(shí)現(xiàn)2D全景影像到3D點(diǎn)云的粗配準(zhǔn)，繼而通過最小化密集匹配全景影像點(diǎn)云與激光點(diǎn)云3D-3D點(diǎn)間距實(shí)現(xiàn)精配準(zhǔn)。該方法彌補(bǔ)了3D-3D配準(zhǔn)方法需要較為準(zhǔn)確初始值方能收斂的缺陷，同時(shí)打破了幾何基元對(duì)配準(zhǔn)方法對(duì)基元正確提取與匹配的依賴。經(jīng)定性與定量配準(zhǔn)精度分析，本文提出的自動(dòng)配準(zhǔn)方法可以實(shí)現(xiàn)車載MMS激光點(diǎn)云與序列全景影像的像素級(jí)自動(dòng)配準(zhǔn)，配準(zhǔn)成果可應(yīng)用于車載MMS激光點(diǎn)云/影像數(shù)據(jù)融合應(yīng)用。由于本文方法將建筑物上邊線即天際線作為配準(zhǔn)基元，在高速公路等無建筑物環(huán)境區(qū)域不適用，下一步研究將拓展幾何基元提取與匹配方法，致力于研制一種全地形環(huán)境適用的車載MMS點(diǎn)云與全景影像配準(zhǔn)方法。

[1] 侯艷芳, 葉澤田, 楊勇. 基于POS數(shù)據(jù)的車載面陣CCD影像與激光點(diǎn)云融合處理研究[J]. 遙感信息, 2011(4): 76-79.

HOU Yanfang, YE Zetian, YANG Yong, Research on Integrated Processing of Vehicle-borne Array CCD Images and Laser Point Cloud Based on POS Data[J]. Remote Sensing Information, 2011(4): 76-79.

[2] 吳勝浩, 鐘若飛. 基于移動(dòng)平臺(tái)的激光點(diǎn)云與數(shù)字影像融合方法[J]. 首都師范大學(xué)學(xué)報(bào)(自然科學(xué)版), 2011, 32(4): 57-61.

WU Shenghao, ZHONG Ruofei. The 3D Visualization’s Approach Based on Vehicle-borne Laser Scanning[J]. Journal of Capital Normal University (Natural Science Edition), 2011, 32(4): 57-61.

[4] WANG Ruisheng, FERRIE F P, MACFARLANE J. Automatic Registration of Mobile LiDAR and Spherical Panoramas[C]∥Proceedings of 2012 IEEE Computer Society Conference on Computer Vision and Pattern Recognition Workshops (CVPRW). Providence, RI: IEEE, 2012: 33-40.

[5] SWART A, BROERE J, VELTKAMP R, et al. Refined Non-rigid Registration of A Panoramic Image Sequence to a LiDAR Point Cloud[M]∥STILLA U, ROTTENSTEINER F, MAYER H, et al. Photogrammetric Image Analysis. Berlin: Springer, 2011: 73-84.

[6] CHEN C, YANG B, DENG C. Registration of the Panoramic Image Sequence onto Mobile Laser Scanning Point Cloud Using Linear Features[C]∥The International Symposium on Mobile Mapping Technology (MMT2013). Taiwan, China: 2013: 1-9.

[7] TANEJA A, BALLAN L, POLLEFEYS M. Geometric Change Detection in Urban Environments Using Images[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2015, 37(11): 2193-2206.

[8] HOFMANN S, EGGERT D, BRENNER C. Skyline Matching Based Camera Orientation from Images and Mobile Mapping Point Clouds[J]. ISPRS Annals of Photogrammetry, Remote Sensing and Spatial Information Sciences, 2014, II-5: 181-188.

[9] VIOLA P, WELLS III W M. Alignment by Maximization of Mutual Information[J]. International Journal of Computer Vision, 1997, 24(2): 137-154.

[10] NELDER J A, MEAD R. A Simplex Method for Function Minimization[J]. The Computer Journal, 1965, 7(4): 308-313.

[11] RUSINKIEWICZ S, LEVOY M. Efficient Variants of the ICP Algorithm[C]∥Proceedings of the 3rd International Conference on 3-D Digital Imaging and Modeling. Quebec City, Canada: IEEE, 2001: 145-152.

[12] YANG Bisheng, DONG Zhen, ZHAO Gang, et al. Hierarchical Extraction of Urban Objects from Mobile Laser Scanning Data[J]. ISPRS Journal of Photogrammetry and Remote Sensing, 2015, 99: 45-57.

[13] ROTHER C, KOLMOGOROV V, BLAKE A. “GrabCut”: Interactive Foreground Extraction Using Iterated Graph Cuts[C]∥Proceedings of ACM SIGGRAPH 2004. New York, NY: ACM, 2004: 309-314.

[14] AWRANGJEB M, LU Guojun, FRASER C S. Performance Comparisons of Contour-based Corner Detectors[J]. IEEE Transactions on Image Processing, 2012, 21(9): 4167-4179.

[15] CHRISTY S, HORAUD R. Iterative Pose Computation from Line Correspondences[J]. Computer Vision and Image Understanding, 1999, 73(1): 137-144.

[16] 鄧非, 張祖勛, 張劍清. 利用激光掃描和數(shù)碼相機(jī)進(jìn)行古建筑三維重建研究[J]. 測(cè)繪科學(xué), 2007, 32(2): 29-30.

DENG Fei, ZHANG Zuxun, ZHANG Jianqing. 3D Reconstruction of Old Architecture by Laser Scanner and Digital Camera[J]. Science of Surveying and Mapping, 2007, 32(2): 29-30.

[17] 王晏民, 胡春梅. 一種地面激光雷達(dá)點(diǎn)云與紋理影像穩(wěn)健配準(zhǔn)方法[J]. 測(cè)繪學(xué)報(bào), 2012, 41(2): 266-272.

WANG Yanmin, HU Chunmei. A Robust Registration Method for Terrestrial LiDAR Point Clouds and Texture Image[J]. Acta Geodaetica et Cartographica Sinica, 2012, 41(2): 266-272.

[18] LIU L, STAMOS I. Automatic 3D to 2D Registration for the Photorealistic Rendering of Urban Scenes[C]∥Proceedings of IEEE Computer Society Conference on Computer Vision and Pattern Recognition. San Diego: IEEE, 2005: 137-143.

[19] RAMALINGAM S, BOUAZIZ S, STURM P, et al. SKYLINE2GPS: Localization in Urban Canyons Using Omni-skylines[C]∥Proceedings of 2010 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS). Taipei, Taiwan, China: IEEE, 2010: 3816-3823.

[20] SATTLER T, LEIBE B, KOBBELT L. Fast Image-based Localization Using Direct 2D-to-3D Matching[C]∥Proceedings of 2011 IEEE International Conference on Computer Vision (ICCV). Barcelona, Spain: IEEE, 2011: 667-674.

[21] KNOPP J, SIVIC J, PAJDLA T. Avoiding Confusing Features in Place Recognition[M]∥DANIILIDIS K, MARAGOS P, PARAGIOS N. Computer Vision-ECCV 2010. Berlin: Springer, 2010: 748-761.

[22] LI Zhuqiang, ZHANG Liqiang, ZHONG Ruofei, et al. Classification of Urban Point Clouds: A Robust Supervised Approach with Automatically Generating Training Data[J]. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 2017, 10(3): 1207-1220.

[23] KWAK E. Automatic 3D Building Model Generation by Integrating LiDAR and Aerial Images Using a Hybrid Approach[D]. Calgary: University of Calgary, 2013.

[24] AWRANGJEB M, LU Guojun, FRASER C S, et al. A Fast Corner Detector Based on the Chord-to-point Distance Accumulation Technique[C]∥Proceedings of 2009 Digital Image Computing: Techniques and Applications. Melbourne, VIC, Australia: IEEE, 2009: 519-525.

[25] 袁修孝, 張雪萍, 付建紅. 高斯-克呂格投影坐標(biāo)系下POS角元素的轉(zhuǎn)換方法[J]. 測(cè)繪學(xué)報(bào), 2011, 40(3): 338-344.

YUAN Xiuxiao, ZHANG Xueping, FU Jianhong. Transformation of Angular Elements Obtained Via A Position and Orientation System in Gauss-Kruger Projection Coordinate System[J]. Acta Geodaetica et Cartographica Sinica, 2011, 40(3): 338-344.

[26] ZHOU Kun, HOU Qinming, WANG Rui, et al. Real-Time KD-Tree Construction on Graphics Hardware[J]. ACM Transactions on Graphics (TOG), 2008, 27(5): 126.

[27] SNAVELY N, SEITZ S M, SZELISKI R. Modeling the World from Internet Photo Collections[J]. International Journal of Computer Vision, 2008, 80(2): 189-210.

[28] TOLA E, LEPETIT V, FUA P. DAISY: An Efficient Dense Descriptor Applied to Wide-baseline Stereo[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2010, 32(5): 815-830.

[29] YANG Bisheng, CHEN Chi. Automatic Registration of UAV-borne Sequent Images and LiDAR Data[J]. ISPRS Journal of Photogrammetry and Remote Sensing, 2015, 101: 262-274.

[30] CHEN Y, MEDIONI G. Object Modeling by Registration of Multiple Range Images[C]∥Proceedings of 1991 IEEE International Conference on Robotics and Automation. Sacramento, CA: IEEE, 1991: 2724-2729.

[31] POMERLEAU F, COLAS F, FERLAND F, et al. Relative Motion Threshold for Rejection in ICP Registration[M]∥HOWARD A, IAGNEMMA K, KELLY A. Field and Service Robotics. Berlin: Springer, 2010: 229-238.