孫健東,陳 需,周 宇,王 韜,林鈺淇,李 方,于鑫鑫,劉 鑫
(1.華北科技學院,河北 三河 065201;2.中國礦業大學(北京) 煤炭資源與安全開采國家重點實驗室,北京 100083;3.神華準格爾能源有限責任公司 黑岱溝露天煤礦,內蒙古 鄂爾多斯 010300;4.中國礦業大學(北京) 能源與礦業學院,北京 100083)
拋擲爆破-拉斗鏟倒堆工藝是世界范圍內露天礦山的主流開采工藝,拋擲爆破后爆堆的形態一方面反映了拋擲爆破的效果,可用于指導爆破參數的優化與設計,另一方面又決定了系統的有效作業量,用于指導后續倒堆作業中工作面的設計,因此爆堆形態的采集是拋擲爆破作業中的一項重要工作[1]。然而,當前現場采用的激光掃描技術手段存在程序復雜、架設位置受限、作業存在安全風險、掃描作業滯后等諸多問題,在很大程度上影響了爆堆采集工作的效率及精度,因此有必要探索更為安全、可靠、快捷的測繪手段。
無人機獲取影像數據具有時效性強、成本低、影像分辨率高、靈活性高等優點[2],諸多學者近年來圍繞無人機在露天礦山測繪等領域的應用開展了研究工作。劉光偉等[3]利用無人機對伊敏露天礦采場、煤場的采剝工程量測繪問題進行了研究;方翔等[4]構建了舟山某大型露天采石礦山的三維模型,實現了礦區三維全景漫游;趙紅澤等[5]提出無人機在露天礦地形建模的應用研究;秦秀山[6]等利用無人機航測技術開展了非接觸式現場工程地質調查等等,上述研究工作充分驗證了無人機應用到礦山測繪工作中的可行性。
綜上,利用無人機搭載圖像傳感器進行爆堆空中掃描,可以以大范圍、高精度的方式系統化感知爆破工作面,有效反映爆堆的外觀、位置、高度等屬性,有效提升模型的生產效率。整套技術方法簡單可靠,爆破后30min內即可完成掃描作業,5~10h內即可交付工作面三維模型,且滿足比例尺1∶500測繪精度要求,為爆堆形態的采集提供了全新的技術手段。
我國黑岱溝露天礦采用拋擲爆破-拉斗鏟倒堆工藝剝離煤層以上約35~50m的巖石,拋擲爆破后,約30%的巖石直接拋擲進入排土場,剩余物料由拉斗鏟進行倒堆處理,黑岱溝穿爆隊利用MDL高精度激光掃描儀進行爆堆掃描作業,激光發射器分辨率為1cm,精度達10cm,掃描角分辨率0.01°,精度為0.02°,掃描范圍:水平角度0~360°,垂直角度-45°~80°,掃描距離700m,掃描速度為7200p/h[7],掃描儀重量約9.7kg,再加上三腳架等裝置后,設備系統總重接近15kg。在現場應用過程中,主要存在以下問題:
1)移設不便。設備架設于高臺階工作面或倒土堆頂部,受掃描角度限制,需要多次移設。
2)受環境及氣候影響較大。工作面凹凸不平,系統的找平定位較為繁瑣,遇到雨、雪等惡劣天氣時,作業更加困難。
3)存在安全隱患。爆破前掃描作業時,高臺階坡頂距離煤層底板約70~80m,激光掃描儀架設在高臺階邊緣時較為危險,如架設在安全距離以外則無法避免掃描盲區。
4)掃描作業滯后。一次掃描作業約耗時2h,拋擲爆破作業在下午5點左右進行,受交接班影響及冬季照明條件限制,一般次日上午10點開始掃描作業,然而此時推土機已經處理了部分物料,因此采集的數據不能反應真實拋擲情況(如圖1所示)。
5)當現場工期接續緊張時,往往不進行爆前采空區掃描作業。然而受煤層底板起伏、采空區堆積物料等因素影響,后期采集到的數據無法對上述因素造成的影響進行分析。
鑒于上述原因,研究提出利用無人機傾斜攝影技術開展爆堆的快速掃描作業。
像控點指無人機航測前,利用RTK、全站儀等設備在拍攝目標測區內建立的具有標志性的真實坐標點,可對后期航測成果進行校正,從而提高模型精度[8]。合理的像控點布設方案會提高原始爆堆三維模型的精度,布設時應盡可能在待測區域內均勻布置。但在現場布點工作過程中發現了以下問題:
1)拋擲爆破區域地形復雜,大部分區域人員無法抵達,難以進行像控點布設工作。
2)坑底采空區位置受倒土堆及高臺階遮擋,無法正常接收RTK信號。
3)工程接續緊張,爆破后20min內推土機駛入爆堆作業,因此在爆堆上無法布置像控點。
綜合考慮上述問題,研究提出利用區域網法的密周邊布點方案,所有像控點均為平高點,爆堆區域大約為550m×250m,選擇布設11個像控點,其中6個作為內業加密像控點,5個作為檢驗點。在保證安全的前提下,像控點盡可能均勻分布在爆堆周邊(如圖2所示)。另一方面,從礦方穿爆隊獲取拋擲爆破高臺階坡頂線及各鉆孔的平面及高程坐標,起到輔助加密的作用,極大程度上提高了未來模型的精度,降低了像控點布設工作量。
2.2.1 試驗無人機參數概述
本次外業航測選取大疆Phantom 4Pro2.0專業版航測無人機,其參數見表1。無人機的相機參數見表2。
2.2.2 無人機外業飛行航高設計
相對航高是指無人機測繪時距離拋擲爆破工作面區域的平均高度,計算公式為:
H=(f×GSD)/a
(1)
式中,H為相對航高,m;f為無人機搭載相機焦距,mm;GSD為地面分辨率,cm;a為像元尺寸,mm。
大疆Phantom 4Pro2.0無人機搭載的是非量測數碼相機,采用8.8mm定焦鏡頭,1英寸CMOS約為13.2×8.8mm,分辨率大小為5472×3648,計算可得到具體的像元尺寸為0.00241mm。根據航測程度要求,1∶500測圖比例尺條件下地面分辨率值應小于等于5cm。根據式1可知,1∶500測圖比例尺對應的最大航高限制為182.57m,研究中為了確保地面分辨率值小于等于2.5cm,因此設定最大航高為90m。
2.2.3 像片重疊度設計
像片重疊具體包括航向重疊和旁向重疊,分別指同一航線方向、相鄰航線上的相鄰相片重疊。合理地規劃航向重疊度以及旁向重疊度可以提高航測產品的精度[9],且影響到無人機的飛行作業時間,在航高一定的情況下,像片重疊度越高,無人機飛行作業時間越長,產品精度也會一定程度提高。重疊度計算公式為:
qx%=(px/lx)×100%
(2)
qy%=(py/lx)×100%
(3)
式中,qx為航向重疊度;qy為旁向重疊度;p為相鄰兩張相片航向重疊部分的長度;py為兩條相鄰航帶之間相鄰相片重疊部分的長度;lx為整個像幅影像的長度。
根據《低空數字航空攝影規范》(CH/Z 4005—2010)的要求[10],無人機航拍的像片航向重疊度應為60%~80%,旁向重疊度應設置為15%~60%。考慮拋擲爆破區域高臺階的遮擋作用,結合國內學者相關研究[11],為避免模型幾何結構的黏連,影像重疊度最多可設計為80%~90%。因此研究中設計航向重疊度及旁向重疊度均為80%。
2.2.4 無人機巡航速度設計
最大巡航速度用于限制無人機的飛行速度,防止影響動態模糊,計算方法為:
Vmax=(δmax×GSD)/t
(4)
式中,Vmax為最大巡航速度,m/s;GSD為地面分辨率,cm;t為相機曝光時間,s;δmax為容許像移值。
根據航測經驗,光圈為f=2.8,感光度ISO為200時,陰天、多云、晴天的相機合適的曝光時間分別為1/200s,、1/1000s、1/2000s,根據式4可知,晴天條件下巡航速度為10m/s,電子快門速度設置為1/2000s時,δmax=20%。為控制容許像移值δmax以保證拍攝質量,需要減少曝光時間且提高的感光度ISO,但感光度提高后通常會導致影像質量降低。研究中經過多次嘗試,確定晴天條件下電子快門速度設置為1/8000s、無人機巡航速度為10m/s,該參數設置條件下δmax=5%,可以得到較好的拍攝效果。
爆堆區域大小約為550m×250m,考慮周圍的像控點,實際航測區域大小約為700m×390m,可將爆堆區域周圍的像控點全覆蓋,保證三維模型的精度。現場試驗時正處冬季,無人機低溫條件下續航時間不足25min,因此為安全起見,無人機每次飛行時間控制在20min以內。研究采用單鏡頭無人機傾斜攝影技術,共需要3個架次飛行,實際耗時約70min,布設像控點的時間約為60min,外業總作業時間約為130min,共獲取1432張原始爆堆圖像,部分原始爆堆影像如圖3所示。

圖3 部分原始爆堆影像
在外業數據采集后,為了獲取拋擲爆破爆堆的原始三維模型,需要進行多視影像聯合平差,充分考慮傾斜攝影圖像之間的幾何變形和遮擋關系。結合無人機POS系統記錄的影像外方位元素,采取由粗到精的金字塔匹配策略,在每級影像上進行同名點自動匹配和自由網光束法平差[12],關鍵步驟如下:
SIFT指尺度不變特征變換,該算法用于解決影像中存在的光照影響、目標遮擋、噪聲等問題,當進行目標物的旋轉、縮放等操作時可保證其不會改變,算法關鍵步驟如下:
3.1.1 構建高斯尺度空間
對圖像進行不同程度的模糊處理是構建高斯尺度空間的關鍵,利用不同的高斯核可以得到不同程度的爆堆模糊圖像,通過爆堆圖像與高斯函數二者進行卷積計算實現爆堆模糊圖像的處理。具體如下:
L(x,y,σ)=G(x,y,σ)*I(x,y)
(5)
式中,*表示卷積運算;(x,y)為圖像空間坐標;I(x,y)為圖像所在位置(x,y)的像素值;σ為尺度空間因子;G(x,y,σ)為尺度可變高斯函數,其中:
2個尺度圖像之差為:
D(x,y,σ)=[G(x,y,kσ)-G(x,y,σ)]*
I(x,y)=L(x,y,kσ)-L(x,y,σ)
(7)
3.1.2 DoG空間極值檢測
該步驟的目的是得到高斯尺度空間的最大或者最小點,比較圖像中標記像素點與其周圍的8個點的大小,若8個點都處于同一尺度空間內,之后再與它上下兩個緊挨著的尺度空間的18個點進行比值大小,最終得到極值點。
3.1.3 關鍵點方向分配
根據SIFT算法的定義可知,要想實現算子尺度不變特征變換,就必須先得到關鍵點的具體的方向參數,方向參數要根據這個點的空間方向特點來計算[13],具體流程包括:
計算各點梯度的模:
m(x,y)=

(8)
計算各點的梯度方向:
式中,m(x,y)是(x,y)所處梯度的模值;θ(x,y)為(x,y)所處梯度方向。然后利用直方圖統計關鍵點臨域內像素對應的梯度方向和幅值,其中直方圖的峰值就是關鍵點的主方向。
3.1.4 關鍵點特征描述
關鍵點對于之后圖像之間的相互匹配至關重要,需要利用向量對關鍵點進行表示,稱為描述子。θ角度為關鍵點的主方向,將關鍵點所在坐標軸旋轉θ角度,得到的新坐標為:

(x,y∈[-radius,radius])
(10)
旋轉后以主方向為中心取8*8的窗口,計算各像素點的梯度方向和幅值,利用高斯窗口進行加權運算后,在降維后得到的小塊上繪制8個方向的梯度直方圖并計算各方向的累加值,即可得到關鍵描述子,描述子生成效果如圖4所示。

圖4 關鍵點描述子生成效果
3.1.5 關鍵點匹配
通過上述步驟可得關鍵點的特征方向以及特征方向的尺寸大小,即最終所需要的關鍵點,由于采集到爆堆影像之間具有較高的重疊度(80%),每張爆堆圖片可理解為由若干個關鍵點所組成,因此可對所有最終確定的關鍵點進行互相匹配,建立關鍵點之間的匹配關系,進而建立爆堆圖像之間的匹配關系。
3.1.6 匹配結果
研究在Matlab平臺下編程實現了上述步驟,以臨近高臺階南部的兩幅爆堆圖像為例,對兩幅爆堆圖像分別進行特征點提取,第一幅圖像提取出4502個特征點,如圖5(a)所示,第二幅圖像提取出4606個特征點,如圖5(b)所示。圖5中圓圈半徑大小代表特征點尺度大小;橫線表示特征點方向。最后對兩幅爆堆圖像進行特征點匹配,共匹配到753對特征點,如圖5(c)所示。

圖5 爆堆圖像特征點檢測及匹配示意圖
光束法區域網空中三角測量是以一張影像組成的一束光線作為平差單元,以中心投影的共線性方程作為平差的基礎方程,對各個光束旋轉和平移使得公共光線的最優交會。最后將整個觀測區域轉化至控制點坐標系統,確定加密點的地面坐標及像片的外方位元素[14]。
其中,以中心投影的共線性方程為[15](詳細推導過程見文獻[15]):
vx=a11ΔXs+a12ΔYs+a13ΔZs+a14Δφ+
a15Δω+a16Δκ+a11ΔX+a12ΔY+a13ΔZ-lx
vy=a21ΔXs+a22ΔYs+a23ΔZs+a24Δφ+
a25Δω+a26Δκ+a21ΔX+a22ΔY+a23ΔZ-ly
(11)
式(11)用矩陣表示為:
式中,
t=[ΔXs,ΔYs,ΔZs,Δφs,Δωs,Δκs]
那么此時所對應的法方程為:
如果對空三的成果進行精度檢驗,可利用外業布設好的檢查點進行計算,檢查點的平面中誤差、高程中誤差均根據式(14)求解。
式中,m為檢查點中誤差;Δ為檢查點測量坐標值和解算值之差;n為檢查點的個數。
研究利用ContextCapture軟件進行爆堆圖像的空三加密。得到了每張原始爆堆影像的精確的外方位元素,如圖6所示。

圖6 爆堆影像的空三加密示意圖
SIFT算法能夠通過影像中特征點的提取實現影像精準匹配,但匹配形成的特征點較為稀疏。例如3.1節中兩幅爆堆圖片僅有753對特征點,因此還需要利用多視影像密集匹配得到同名點的空間坐標信息形成密集的三維點云。這一過程中利用的算法是PMVS(Patch-based Multi-view Stereo),算法主要步驟為:初始特征匹配,生成初始稀疏的匹配點,形成下一步匹配傳播的種子點;擴散,由稀疏的種子點擴展得到密集的點云[16,17];過濾,刪選面片,確保結果的準確性。
3.3.1 初始種子面片生成
首先,利用Harris算子和DOG算子提取各幅影像的特征點,通過多幅影像間特征點的匹配構建稀疏面片集,并選出參考影像與待匹配影像。以各幅影像作為參考影像,依次與待匹配影像進行匹配,并存儲于圖像塊之中。然后將每一對匹配得到的特征點對進行前方交會計算地面點的三維坐標式(15),形成種子面片。
式中,x、y、z分別代表同名點坐標和物方點深度信息;X、Y、Z為物方點的三維坐標;P為投影矩陣。
3.3.2 面片擴散
基于上一步得到的種子面片在其周圍生成待擴散的面片,其三維坐標和法向量由種子面片初始化,根據待擴散面片在影像上的投影位置灰度信息的相似度來優化其具體的三維坐標和法向量。如果面片p想要擴散到其它區域生成新的面片,必須滿足式(16):
C(p)={Ci(x′,y′)|p∈Q(x,y),
|x-x′|+|y-y′|=1}
(16)
式中,C(p)為圖像塊領域集合;C(x′,y′)為圖像塊;Q(x,y)為圖像塊中所有面片集合。
3.3.3 面片過濾
重建過程中可能會生成誤差較大的面片,因此需要通過對以下幾類面片進行過濾以保證結果的準確性[18]:
1)處于同一圖像塊中且平均相關系數差異較大的面片。
2)根據光照一致性約束,計算格網中各面片深度信息并統計獲取可視圖像數量,若小于某閾值,則予以剔除。
3)對于每個面片p,統計準可視圖像中屬于p所在圖像塊及相鄰圖像塊中的面片總數,若p的臨近面片占面片總數的比例低于1/4,則濾掉該點。
基于上述三維建模關鍵技術,使用Matlab、Context Capture等平臺進行爆堆影像的三維重建。建模采用高斯-克里格3度帶投影,中央子午線為111°,平面坐標系采用國家2000坐標系,高程采用1985國家高程基準。具體成果如圖7所示。

圖7 模型構建成果
利用5個檢查點對原始爆堆三維模型進行精度驗證,分別計算模型的平面中誤差、高程中誤差,見表3、表4。

表3 實測檢查點坐標 m

表4 檢查點殘差表 m
經計算平面中誤差m=±0.129m,高程中誤差mz=±0.263m。根據《三維地理信息模型數據產品規范》,比例尺1∶500測繪成果I級精度要求平面位置中誤差小于等于0.3m,高程中誤差差小于等于0.5m,無人機測繪成果滿足上述要求。
通過實踐與研究可以看出,與礦方傳統激光掃描技術相比,無人機傾斜攝影技術具有以下優點:
1)測繪作業效率高。原先激光掃描工作需要第二天進行,第三日完成模型。應用無人機技術在爆破后即可開展測繪,且在第二天早調會議前形成成果。
2)設備成本低。礦方傳統的三維激光掃描儀價格約為40萬元(不含相關軟件費用),而研究中采用的大疆Phantom4 Pro2.0專業版航測無人機僅為3萬元。
3)安全系數高。拋擲爆破臺階拋擲前高達70~80m,屬于高危險區域,傳統測繪方法存在安全隱患,使用無人機測繪技術可極大保障人員及設備的安全。
4)作業人員少,工作強度低。傳統的掃描方式需4人協作,而無人機體積小、重量輕、易攜帶,僅需要2人甚至1人即可完成工作,且操作過程簡單,輸入參數后可實現一鍵作業。
5)獲取的數據更加準確、完整。原先激光掃描方式存在時間滯后性,掃描中無法避免掃描盲區,而無人機測繪技術可以及時、無死角地獲取爆堆真實形態。
6)為評價與設計提供翔實資料。無人機技術可以得到爆破前后工作面三維模型,一方面為爆破效果評價提供精準的數據,另一方面為后續拉斗鏟倒堆作業工作面的設計提供了依據。
1)論文針對基于無人機傾斜攝影的拋擲爆破爆堆形態測量方法展開研究,提出了無人機像控點布設方案,充分結合生產現場特點,利用臺階坡頂線坐標、鉆孔坐標等實現了像控點加密,提高了模型精度、降低了外業工作量。
2)經過試驗總結得出,該場景下應用大疆消費機無人機,為保證地面分辨率值小于等于2.5cm,最大航高應為90m、航向重疊度及旁向重疊度為80%、晴天條件下電子快門速度設置為1/8000s、無人機巡航速度為10m/s。
3)對無人機內業數據處理過程進行了詳細的分析,具體包括圖像匹配算法、空三測量原理、影像密集匹配流程等。
4)經過詳細的理論分析及實踐驗證,研究形成了系統化的無人機現場采集作業方案,模型誤差完全滿足比例尺1∶500測繪成果I級精度要求。可直接供類似應用場景下露天礦山的測繪工作參考與借鑒。
5)另外受經費所限,研究中應用的無人機為單攝像頭且不搭載RTK模塊,外業測繪總時間為130min。若采用更為先進的五鏡頭傾斜攝影無人機并搭載RTK模塊,即可實現單架次飛行、免像控作業,將整體測繪時間縮短至20min以內。