基于地面三維激光掃描的露天礦采剝工程量計(jì)算方法

王 森 何 群 劉善軍 毛亞純

（東北大學(xué)資源與土木工程學(xué)院，遼寧沈陽110819）

采場驗(yàn)收測量是露天礦山測量中的一項(xiàng)主要工作，而采剝工程量計(jì)算則是其中最頻繁、最重要的工作，既為露天礦檢查生產(chǎn)計(jì)劃的執(zhí)行情況提供依據(jù)，又為安排和調(diào)整未來生產(chǎn)計(jì)劃提供不可或缺的資料。目前，露天采場驗(yàn)收測量手段主要使用全站儀對采場斷面進(jìn)行測量［1］，而后利用斷面法計(jì)算采掘礦（巖）塊體積，再與礦（巖）塊體密度相乘，即可獲得采掘量。由于該方法屬于點(diǎn)式測量，因而存在采樣密度低、體積計(jì)算誤差大、測量周期長等不足。近年來，隨著具有精度高、速度快等優(yōu)勢的面式測量手段（包括三維激光掃描、無人機(jī)航攝測量等）［2-3］的興起，大大提高了礦山測量效率。

目前，三維激光掃描系統(tǒng)是較為先進(jìn)的獲取地面空間多目標(biāo)三維數(shù)據(jù)長距離影像的測量技術(shù)，被稱為測繪領(lǐng)域繼GPS技術(shù)之后的又一技術(shù)革命［4］。三維激光掃描獲取的數(shù)據(jù)集合——點(diǎn)云，是空間目標(biāo)實(shí)體的真實(shí)表達(dá)，能夠用于生成高精度DTM。但在大多數(shù)情況下，一次掃描范圍不足以獲得整個場景，需要利用部分重疊的高空間信息點(diǎn)云進(jìn)行拼接［5］；另一方面，在進(jìn)行多期點(diǎn)云變化監(jiān)測、土石方量計(jì)算時，需要將不同坐標(biāo)系統(tǒng)的點(diǎn)云歸算到同一坐標(biāo)系統(tǒng)下［6］，這就需要對點(diǎn)云進(jìn)行精確配準(zhǔn)。因此，點(diǎn)云數(shù)據(jù)精確配準(zhǔn)便成為土石方量計(jì)算的一個關(guān)鍵環(huán)節(jié)。目前，Besl等［7］提出的ICP算法及改進(jìn)ICP算法是應(yīng)用最廣、效果最佳、最為理想的點(diǎn)云精配準(zhǔn)算法。但該算法對初始粗配準(zhǔn)的輸入要求較高，否則，會影響收斂速度，并且可能出現(xiàn)局部收斂現(xiàn)象［8-10］。Zheng等［11］提出了一種利用點(diǎn)云融合影像數(shù)據(jù)的配準(zhǔn)方法，該方法利用相機(jī)和LiDAR掃描儀分別獲取研究區(qū)域的影像和點(diǎn)云2種數(shù)據(jù)，將SIFT算法對影像提取出的特征點(diǎn)用作改進(jìn)ICP算法過程中的對應(yīng)點(diǎn)，以此進(jìn)行配準(zhǔn)迭代。該方法在配準(zhǔn)過程中增加了一種影像數(shù)據(jù)源，并且在ICP精確配準(zhǔn)中僅考慮了特征點(diǎn)配準(zhǔn)，未顧及到全局點(diǎn)的配準(zhǔn)計(jì)算。

本研究在分析上述成果的基礎(chǔ)上，對原有SIFTICP點(diǎn)云配準(zhǔn)算法進(jìn)行改進(jìn)，并利用露天礦三維激光掃描數(shù)據(jù)進(jìn)行算法驗(yàn)證，實(shí)現(xiàn)基于地面三維激光點(diǎn)云數(shù)據(jù)的露天礦采剝工程量精準(zhǔn)計(jì)算。

1 新SIFT-ICP算法原理

露天采場常用的采剝驗(yàn)收量計(jì)算方法有斷面法、規(guī)則方格網(wǎng)法、等高線法、DTM（Digital terrain model）三角網(wǎng)法等［1］。DTM由Miller教授于1956年提出，能夠反映地表真實(shí)特征，用于各種線路設(shè)計(jì)，各種工程面積、體積、坡度計(jì)算，任意兩點(diǎn)間可視性判斷以及任意斷面圖繪制等［12］。與其他幾種方法對比，DTM能夠全面真實(shí)地描述地形信息，差值計(jì)算采剝工作量的精度較高［13-14］，能夠滿足礦山驗(yàn)收測量的要求。DTM三角網(wǎng)法計(jì)算露天礦采掘工程量，是根據(jù)實(shí)測數(shù)據(jù)生成三角網(wǎng)，利用三角網(wǎng)計(jì)算擬定設(shè)計(jì)高程下的三棱柱體積，對指定區(qū)域內(nèi)每個三棱柱體積進(jìn)行求和得到整個區(qū)域在擬定高程下的體積總和，不同時間序列、同一區(qū)域、同等設(shè)計(jì)高程下的體積差值即為該區(qū)域的采掘工程量［15］。在整個過程中，點(diǎn)云精確配準(zhǔn)直接影響了三棱柱體積的計(jì)算精度，進(jìn)而影響了整個區(qū)域工程量的計(jì)算精度。

為構(gòu)建研究區(qū)域完整的三維模型，需要將多個視點(diǎn)下獲取的點(diǎn)云進(jìn)行重新定位，生成一個統(tǒng)一坐標(biāo)系下的三維數(shù)據(jù)點(diǎn)集，即點(diǎn)云配準(zhǔn)［16］。點(diǎn)云配準(zhǔn)分為粗配準(zhǔn)和精配準(zhǔn)，粗配準(zhǔn)是對2個坐標(biāo)系統(tǒng)相差較大的點(diǎn)云進(jìn)行配準(zhǔn)，粗配準(zhǔn)后的點(diǎn)云能夠達(dá)到坐標(biāo)系統(tǒng)大體一致。ICP精配準(zhǔn)算法是在粗配準(zhǔn)的基礎(chǔ)上進(jìn)行全局優(yōu)化配準(zhǔn)，通過全局點(diǎn)不斷迭代尋找對應(yīng)點(diǎn)的關(guān)系，計(jì)算旋轉(zhuǎn)矩陣R及平移系數(shù)T，使得所有對應(yīng)點(diǎn)之間歐式距離的平方和達(dá)到最小，以完成配準(zhǔn)。ICP算法精確配準(zhǔn)的初始狀態(tài)如果輸入不佳，即會影響后續(xù)迭代解算，導(dǎo)致迭代收斂速度慢、迭代不收斂或無法收斂到正確結(jié)果，易導(dǎo)致出現(xiàn)局部最優(yōu)解，使得配準(zhǔn)誤差偏大。

為優(yōu)化配準(zhǔn)算法，本研究首先設(shè)計(jì)了一種無靶標(biāo)點(diǎn)的雙線性內(nèi)插SIFT特征點(diǎn)粗配準(zhǔn)方法，然后將該方法的配準(zhǔn)結(jié)果作為輸入數(shù)據(jù)，再進(jìn)行ICP配準(zhǔn)，從而將SIFT算法與ICP算法的優(yōu)勢相結(jié)合，實(shí)現(xiàn)精確配準(zhǔn)。由于該方法面向單一數(shù)據(jù)源（點(diǎn)云），而原有的SIFT-ICP算法面向2種數(shù)據(jù)源（可見光影像、點(diǎn)云數(shù)據(jù)）［11］，兩者具有明顯區(qū)別，因此該方法可稱為新SIFT-ICP算法。該方法首先利用SIFT算法在實(shí)體模型上提取出特征點(diǎn)集，然后在特征點(diǎn)集上利用雙線性內(nèi)插方法求解點(diǎn)云點(diǎn)集對，以此進(jìn)行點(diǎn)云粗配準(zhǔn)。

SIFT算法由 Lowe［17］提出，該算法設(shè)計(jì)的圖像特征描述子能夠?qū)τ趫D像平移、旋轉(zhuǎn)、縮放等特性保持不變。該算法的主要實(shí)現(xiàn)步驟為：①使用高斯函數(shù)構(gòu)建不同尺度空間；②在高斯差分（Difference of Gaussian，DOG）尺度空間進(jìn)行局部極值點(diǎn)檢測；③為極值點(diǎn)分配128維方向參數(shù)；④生成特征點(diǎn)描述。

ICP算法是對SIFT粗配準(zhǔn)結(jié)果進(jìn)行基于最小二乘運(yùn)算的最優(yōu)匹配方法［18］，通過重復(fù)進(jìn)行“確定對應(yīng)點(diǎn)集、計(jì)算最優(yōu)剛體變換過程”等步驟，直至某個表示正確匹配的收斂準(zhǔn)則得到滿足為止。該算法主要實(shí)現(xiàn)步驟為：①求解待配準(zhǔn)點(diǎn)云在參考點(diǎn)云中的最近點(diǎn)，記為對應(yīng)點(diǎn)集，計(jì)算最近點(diǎn)集距離平方和；②根據(jù)對應(yīng)點(diǎn)集計(jì)算轉(zhuǎn)換參數(shù)，由轉(zhuǎn)換參數(shù)對待配準(zhǔn)點(diǎn)云進(jìn)行坐標(biāo)轉(zhuǎn)換；③轉(zhuǎn)換后的點(diǎn)云與參考點(diǎn)再次求解最近點(diǎn)集，計(jì)算最近點(diǎn)集之間的距離平方和，若2次距離平方和的差值小于既定閾值，則停止迭代，否則，重復(fù)①～③步。可見，SIFT算法能夠?qū)崿F(xiàn)特征點(diǎn)檢測，提取對應(yīng)點(diǎn)關(guān)系，ICP算法能夠完成最優(yōu)匹配。

綜合上述分析，新SIFT-ICP算法的主要實(shí)現(xiàn)步驟為：

（1）采用Delaunay三角剖分算法構(gòu)建露天采場DTM，并對模型進(jìn)行柵格化處理。

（2）在柵格化的模型上采用SIFT算法檢測像平面直角坐標(biāo)系下的特征點(diǎn)集。

（3）對于在像平面直角坐標(biāo)系中的任意一點(diǎn)（x，y）（圖1），采用雙線性內(nèi)插方法進(jìn)行轉(zhuǎn)換，將像平面直角坐標(biāo)系x-y上的特征點(diǎn)（x，y）轉(zhuǎn)換至點(diǎn)云對應(yīng)的大地空間直角坐標(biāo)系X-Y-Z在XY面上的投影，即XY平面直角坐標(biāo)系下的坐標(biāo)（X，Y），

式中，x1、y1為DTM 柵格化后的像平面長與寬，（X1，Y1）與（X2，Y2）分別為像平面坐標(biāo)（圖1）。根據(jù)X-Y平面直角坐標(biāo)系下的坐標(biāo)，在研究區(qū)域中的DEM（圖1）上搜索確定出Z值，便可得到點(diǎn)云在大地空間直角坐標(biāo)系中的特征點(diǎn)集。

（4）根據(jù)特征點(diǎn)集，對原始點(diǎn)云數(shù)據(jù)進(jìn)行粗配準(zhǔn)。

（5）由粗配準(zhǔn)后的點(diǎn)云數(shù)據(jù)進(jìn)行ICP點(diǎn)云數(shù)據(jù)精配準(zhǔn)。

2 算法試驗(yàn)

2.1 數(shù)據(jù)獲取

本研究試驗(yàn)區(qū)域位于遼寧省鞍山市鞍千礦礦區(qū)，該礦始建于2005年，于2006年8月正式投產(chǎn)。礦區(qū)由3個采區(qū)組成，分別為許東溝、啞巴嶺、西大背，3個采區(qū)均已形成固定幫坡。試驗(yàn)區(qū)域啞巴嶺采場最低露天開采標(biāo)高為-132 m。試驗(yàn)采用Riegl VZ-4000三維激光掃描儀，最大測量距離為4 000 m，最近測量距離為5 m，垂直視場掃描范圍為100°（+30°～-30°），水平視場掃描范圍為0～360°，角度分辨率優(yōu)于0.000 5°，測點(diǎn)精度為15 mm。采用該型掃描儀對露天采場進(jìn)行精細(xì)掃描，分別獲取2期三維激光掃描數(shù)據(jù)，數(shù)據(jù)采集時間與采場驗(yàn)收測量同步（圖2）。

2.2 SIFT-ICP算法應(yīng)用

2.2.1 SIFT算法實(shí)施

對于完成預(yù)處理后的數(shù)據(jù)進(jìn)行Delaunay三角剖分建模，并對柵格化后的模型采用SIFT算法進(jìn)行特征點(diǎn)集提取，而后進(jìn)行同名點(diǎn)匹配，匹配結(jié)果如圖3所示。

由于鞍千礦啞巴嶺采場經(jīng)過多年開采，外圍已經(jīng)形成固定邊幫，模型匹配得到的點(diǎn)集主要分布于礦坑邊幫區(qū)域。由于SIFT算法采用128維特征描述的局限性，匹配點(diǎn)會存在一些錯誤點(diǎn)，有必要刪除。錯誤點(diǎn)識別主要考慮兩個方面：

（1）局部區(qū)域是否滿足平行關(guān)系，由于邊幫匹配得到的特征點(diǎn)沒有發(fā)生變形，在局部區(qū)域僅存在點(diǎn)云坐標(biāo)系引起的系統(tǒng)誤差，沒有匹配時的離散隨機(jī)誤差，匹配點(diǎn)之間的連線應(yīng)該滿足平行關(guān)系，本研究將該類點(diǎn)判定為有效特征點(diǎn)。

（2）對于分布于礦坑內(nèi)部的特征點(diǎn)，由于進(jìn)行了采掘工作，匹配得到的特征點(diǎn)之間的連線與周圍邊幫難以構(gòu)成平行關(guān)系，故判定為奇異點(diǎn)；對于圖中還有部分點(diǎn)已經(jīng)匹配到模型外圍，可視其為無效特征點(diǎn)。本研究經(jīng)過錯誤點(diǎn)剔除后保留了礦坑外圍的97個有效特征點(diǎn)建立對應(yīng)點(diǎn)集。

2.2.2 粗配準(zhǔn)

依據(jù)SIFT算法匹配以及篩選后得到的特征點(diǎn)集，采用雙線性內(nèi)插方法，計(jì)算出對應(yīng)點(diǎn)集在大地空間直角坐標(biāo)系下的點(diǎn)，根據(jù)對應(yīng)點(diǎn)集對原始點(diǎn)云數(shù)據(jù)進(jìn)行粗配準(zhǔn)。以初期點(diǎn)云數(shù)據(jù)為參考，對末期點(diǎn)云數(shù)據(jù)進(jìn)行配準(zhǔn)，結(jié)果如圖4所示，其中淡色區(qū)域?qū)?yīng)為第一期點(diǎn)云數(shù)據(jù)，深色區(qū)域?qū)?yīng)為第二期點(diǎn)云數(shù)據(jù)。對于經(jīng)過粗配準(zhǔn)后的點(diǎn)云疊加結(jié)果（圖4（b）），由于特征點(diǎn)集主要分布于右下邊幫，在礦坑右下邊幫中直觀上能夠看出深淺相間的點(diǎn)，粗配準(zhǔn)效果較好。

2.2.3 ICP精配準(zhǔn)

SIFT算法提取出的特征點(diǎn)數(shù)量有限，難以實(shí)現(xiàn)礦坑整個區(qū)域的精確配準(zhǔn)。礦山生產(chǎn)過程中，由于采掘作業(yè)、道路鋪設(shè)、礦車碾壓等都會對礦區(qū)形態(tài)造成影響，使得不同時期點(diǎn)云之間的形態(tài)相差較大。為此，在SIFT算法粗配準(zhǔn)的基礎(chǔ)上，采用了一種基于穩(wěn)定區(qū)域的全局匹配方法。該方法通過選取時間序列段內(nèi)穩(wěn)定不變的區(qū)域，在該區(qū)域內(nèi)應(yīng)用ICP精配準(zhǔn)算法求解轉(zhuǎn)換參數(shù)，再根據(jù)轉(zhuǎn)換參數(shù)對全局區(qū)域進(jìn)行轉(zhuǎn)換配準(zhǔn)。

考慮到礦坑內(nèi)部存在非穩(wěn)定點(diǎn)，難以使用ICP算法進(jìn)行精配準(zhǔn)，因此從點(diǎn)云中裁剪出邊幫穩(wěn)定區(qū)域，穩(wěn)定區(qū)域根據(jù)粗配準(zhǔn)結(jié)果的差值進(jìn)行DEM圈定。對粗配準(zhǔn)后的區(qū)域進(jìn)行ICP算法精配準(zhǔn)，得到邊幫穩(wěn)定區(qū)域的配準(zhǔn)結(jié)果，如圖5（a）所示，其中淡色區(qū)域?qū)?yīng)為第一期點(diǎn)云數(shù)據(jù)，深色區(qū)域?qū)?yīng)為第二期點(diǎn)云數(shù)據(jù)。在此基礎(chǔ)上，將ICP算法得到的最終轉(zhuǎn)換參數(shù)應(yīng)用到整個礦區(qū)點(diǎn)云，得到整個區(qū)域的配準(zhǔn)結(jié)果，如圖5（b）所示。

將點(diǎn)云數(shù)據(jù)精配準(zhǔn)結(jié)果與粗配準(zhǔn)結(jié)果進(jìn)行了對比分析，結(jié)果見表1。由表1可知：采用SIFT算法提取模型的特征點(diǎn)，用于配準(zhǔn)的特征點(diǎn)有97個，該類點(diǎn)主要分布于邊幫區(qū)域，并且在空間分布上也相對均勻，能夠滿足粗配準(zhǔn)要求；利用穩(wěn)定的邊幫區(qū)域進(jìn)行精確配準(zhǔn)，2期配準(zhǔn)點(diǎn)個數(shù)為847 874個，配準(zhǔn)點(diǎn)數(shù)大大增加；2種方法的均方根誤差RMSE相差1.169 105 m。

2.3 采掘邊界線圈定

在同一露天礦中，一般是多個采掘區(qū)域同時進(jìn)行采掘，不同的采掘區(qū)域一般是在不同的臺階上。在進(jìn)行采掘工作量計(jì)算過程中，需要分別計(jì)算不同采區(qū)的工程量。為此，需要首先圈定出不同采區(qū)的范圍。根據(jù)地面三維激光掃描獲取的點(diǎn)云數(shù)據(jù)進(jìn)行采區(qū)范圍圈定，具體步驟為：①由配準(zhǔn)后得到的2期三維激光點(diǎn)云數(shù)據(jù)，分別利用插值方法得到2期DEM，即：DEM1、DEM2；②根據(jù) 2期 DEM，計(jì)算出DEM之間的差值，ΔDEM=DEM2-DEM1；③由ΔDEM內(nèi)插出指定高差的等值線，選擇高差為0 m的等值線。

根據(jù)精配準(zhǔn)后的點(diǎn)云數(shù)據(jù)，內(nèi)插出了2期DEM，如圖6所示。

對上述2期DEM計(jì)算出ΔDEM，而后圈定了高差為0 m的等值線圖（圖7）。由圖7可知：由于受到道路鋪設(shè)、礦車碾壓等因素的影響，導(dǎo)致提取得到的等高線與實(shí)際采場范圍略有偏差，故對等高線進(jìn)行了修剪，進(jìn)一步圈定了采場的大致范圍，圖7中圈定的深色區(qū)域可視為采場的主要采剝區(qū)域。

2.4 采掘量計(jì)算

采用Delaunay三角剖分算法構(gòu)建露天采場DTM，圖8為2期DTM疊加效果，其中圈定的不規(guī)則深色區(qū)域?qū)?yīng)為第1期DTM，其余區(qū)域?qū)?yīng)為第2期DTM。通過該圖能夠清晰看到采掘范圍，與圖7中等值線范圍大體一致。此外，在該圖上雖然能夠看到有一些星星點(diǎn)點(diǎn)的變化點(diǎn)，并且該類點(diǎn)對可視化以及等值線提取會有一定的影響，但在圈定采掘邊界線后，對于采場大范圍工程計(jì)算來說影響較小。

根據(jù)采剝區(qū)域的不同巖礦類型，將整個區(qū)域劃分為6個區(qū)域，如圖9所示。圖中各區(qū)域的采剝工程量計(jì)算結(jié)果見表2。由表2可知：基于點(diǎn)云數(shù)據(jù)計(jì)算出的結(jié)果與礦山實(shí)際記錄的采礦量偏差較小，根據(jù)《有色金屬礦山生產(chǎn)技術(shù)規(guī)程》，當(dāng)?shù)V巖驗(yàn)收量在20萬t以上時，驗(yàn)收量允許誤差不大于1%，該方法計(jì)算結(jié)果的相對誤差為0.72%，能夠滿足礦山計(jì)算要求。此外，本研究采用三維激光掃描方法獲取點(diǎn)云數(shù)據(jù)，有助于大幅提升露天礦山驗(yàn)收測量效率。

3 結(jié)語

針對ICP精配準(zhǔn)對初始輸入點(diǎn)云要求較高，以及原有SIFT-ICP算法需要使用影像數(shù)據(jù)作為SIFT算法粗配準(zhǔn)數(shù)據(jù)源的問題，將原有SIFT-ICP算法進(jìn)行了改進(jìn)，提出了新的SIFT-ICP算法。該算法僅使用點(diǎn)云一種數(shù)據(jù)源，將SIFT算法的快速粗配準(zhǔn)與ICP全局精配準(zhǔn)優(yōu)勢相結(jié)合，實(shí)現(xiàn)了點(diǎn)云數(shù)據(jù)的快速和精確配準(zhǔn)，能夠在不需要設(shè)定靶標(biāo)點(diǎn)、不需要影像輔助數(shù)據(jù)輸入的情況下，僅使用三維激光掃描數(shù)據(jù)完成快速精確配準(zhǔn)?；谛碌腟IFT-ICP算法并利用地面三維激光點(diǎn)云數(shù)據(jù)進(jìn)行了實(shí)際露天礦山驗(yàn)證試驗(yàn)，結(jié)果表明，該算法能夠快速精確完成配準(zhǔn)，計(jì)算的采剝工作量相對誤差為0.72%，能夠滿足礦山采場驗(yàn)收測量要求。