朱 征,包騰飛,3,鄭東健,龔 健,胡雨菡,談家誠
(1.河海大學水文水資源與水利工程科學國家重點實驗室,江蘇 南京 210098;2.河海大學水利水電學院,江蘇 南京 210098; 3.三峽大學水利與環(huán)境學院,湖北 宜昌 443002)
我國西南地區(qū)滑坡、泥石流、崩塌等地質災害多發(fā),運動土石沖入河道就會截斷河流,形成堰塞湖。僅2008—2017年,有記錄的滑坡堰塞壩超過100座[1]。該地區(qū)歷史上曾發(fā)生過多起大型堰塞堵江事件,如1786年大渡河[2]、1933年疊溪[3]、1967年唐古棟[4]、2000年易貢[5-6]、2008年唐家山[7]、2010年舟曲[8]、2014年紅石巖[9]。堰塞湖形成迅速,牽涉范圍廣,具有很高的危險性。2018年10和11月發(fā)生的白格堰塞事件,位于四川省白玉縣則巴村和西藏自治區(qū)江達縣白格村交界處(31.08°N, 98.71°E)。白格堰塞湖是同一地點連續(xù)兩次大規(guī)模高位滑坡阻斷金沙江干流形成的堰塞湖[10],造成上游多個鄉(xiāng)被淹,下游多座橋梁被沖毀,四川、西藏及云南多地緊急轉移安置數(shù)萬人,引起了社會的廣泛關注。
堰塞體往往結構松散、穩(wěn)定性差,極易受外界環(huán)境影響發(fā)生形態(tài)變化,存在隨時潰決的風險。堰塞體一旦失穩(wěn),其蓄積的大量水體將沖擊沿線數(shù)百公里,造成的破壞難以估量。作為了解災情和科學決策的重要依據(jù),快速、準確地獲取現(xiàn)場信息必不可少。然而,受災地區(qū)一般伴有通信與電力中斷、交通癱瘓、基礎設施遭受重大破壞等緊急情況,堰塞災害的事發(fā)地也存在交通不便、環(huán)境原始、人員物資輸送困難等問題。這些因素制約了常規(guī)測繪手段的開展。
無人機傾斜攝影作為一種新興的測繪技術,在災害救援[11]、建筑[12]、城市規(guī)劃[13]、測繪[14]、考古[15]、水利規(guī)劃[16]等領域已經有許多應用。Kirnbauer等[17]利用傾斜攝影量化研究了高山雪線的變化。Doneus[18]使用航空影像進行傾斜攝影成圖,從而還原了Teurnia古城遺址。Mergili等[19]以Piz Cengalo-Bondo滑坡事件的傾斜攝影重建成果為基礎,對該災害事件進行了后分析。Tonkin等[20]使用地面控制點提升了無人機傾斜攝影重建精度,并在一處野外邊坡開展了驗證試驗。Harwin等[21]利用無人機采集海邊懸崖的圖像,并以重建點云實現(xiàn)了亞分米級侵蝕監(jiān)測。Eltner等[22]在野外搭設多目相機系統(tǒng),監(jiān)測了雨水沖刷引起的土壤流失。在滑坡災害研究方面,帥向華等[23]將傾斜攝影技術用于在云南魯?shù)榈卣鹨l(fā)的滑坡測繪;許建華等[24]利用無人機傾斜攝影對地震裂縫、滑坡等幾何信息進行量測,從而為烈度評估提供參考;許強等[25-26]將傾斜攝影用于白格堰塞堵江事件的災害調查;王曉剛等[27]利用旋翼無人機傾斜攝影技術在茂縣高位滑坡災害進行了應急測繪。陳誠等[28]采用無人機進行了表面流場測量。
傾斜攝影在眾多領域的運用中表現(xiàn)出該技術具有測繪效率高、人力資源占用較少、成果呈現(xiàn)方式多元化等優(yōu)點,具備用于堰塞現(xiàn)場的非接觸式應急測繪的潛力。
本文介紹了傾斜攝影的建模原理,闡述了無人機傾斜攝影的步驟要點,并在白格堰塞現(xiàn)場進行了無人機飛行作業(yè)。重建了白格堰塞現(xiàn)場的三維數(shù)字模型,并生成了數(shù)字高程模型(DEM)、數(shù)字地表模型(DSM)、數(shù)字正射圖像(DOM)。無人機傾斜攝影三維重建技術在白格堰塞現(xiàn)場測繪方面表現(xiàn)出了良好的適用性。
傾斜攝影是運用無人機云臺攝像機,以一定傾斜角度獲取地面物體的圖像。通過設定拍攝間隔,保證圖像間存在重疊。利用圖像間的重疊部分進行空中三角測量計算,還原相機姿態(tài)并獲得稀疏點云。再以稀疏點云為基礎進行稠密重建得到稠密點云,從而實現(xiàn)傾斜攝影三維重建。
空中三角測量(aerial triangulation)包含4個部分:圖像特征檢測、圖像特征點匹配、點云初值計算、誤差優(yōu)化。
圖像特征檢測指的是在圖像中找出具有一定代表性的像素點,并且用描述算子對該點進行描述。在眾多算法中,SIFT(scale-invariant feature transform)算法[29]得到廣泛應用。該方法對旋轉、尺度縮放、亮度變化保持不變性,同時對視角變化、仿射變換、噪聲也保持一定程度的穩(wěn)定性。
特征點匹配將圖像間的具有一定相似度的特征點進行配對。如圖1所示,傾斜攝影圖像間存在重疊,重疊區(qū)域內存在大量的相似特征點。相似度計算主要有兩類方法,一類是按候選點描述算子間的空間距離按閾值過濾,另一類則是采用聚類算法尋找同類點。由于匹配時不可避免地存在一些錯誤匹配,這需要應用雙重檢驗機制來提高準確性,如隨機采樣一致性(RNSAC)[30]算法。
圖1 基于重疊圖像的特征檢測原理
點云初值計算是利用特征點間的配對關系進行的特征點三維坐標計算的過程。根據(jù)經典的小孔成像模型,像素點、實物點、焦點之間滿足式(1)的關系。在此成像模型下,實物點坐標和像素點坐標之間呈線性關系。當以多視角拍攝物體時,相機焦點、實物點、像素點間存在一定關系,即外極幾何約束。如圖2所示,兩張圖像拍攝于不同位置,對應的焦點位置為F1、F2,其中,像素點A1、A2指向同一實物點A,那么在極線平面F1AF2上存在關系式(2)。
圖2 基于兩視圖的實物點坐標計算原理
式中:K為內部參數(shù)矩陣;f為相機焦距;u、v為相機主點位置;R、t分別為相機成像坐標系與世界坐標間的旋轉矩陣和平移矩陣;X、Y、Z為實物點坐標;x、y為像素點坐標。
(2)
式中:A1、A2分別為像素點A1、A2的坐標(x1、y1)和(x2、y2);K1、K2分別為兩相機的內部參數(shù)矩陣;E為本質矩陣。
式(2)中僅本質矩陣E中存有未知量,應用8組以上的特征點對即可求解E全部未知量。再通過SVD分解得到相機間相對姿態(tài),以此代入式(1)可計算出每個特征點的三維坐標。匹配點組數(shù)越多,該過程的準確度越高。
實際成像過程與理論存在差異,如鏡頭畸變、投影模型非線性、特征點的檢測與匹配錯誤。這些由各處積累的誤差會降低實物點三維坐標計算的精準度。目前,得到廣泛運用的思想是最小化整體重投影誤差平方(式(3))。重投影誤差平方和最小化的實現(xiàn),普遍采用的是光束平差法(bundle adjustment)。這是一種啟發(fā)式的阻尼高斯牛頓法,在幾何視覺中廣泛使用,具有較高的求解速度和穩(wěn)健性。
(3)
式中:P為實物點向像素點的投影矩陣;Xi為實物點坐標;xi為像素點坐標。
空中三角測量通過對圖像集特征點的三維投影計算,重建了特征點對應的實物點云。然而,特征點數(shù)量相對全部像素點而言是稀疏的,因此空中三角測量得到的點云也是天然稀疏的,包含信息豐富度較低。
多視影像密集匹配(multiple view stereo)以稀疏點云為基礎,以一定方法尋找具有一致性的像素點加密點云,從而增強點云的可利用度。點云擴散法是較為常用的稠密重建方法,主要思想是以初始點云為基礎按照鄰近關系逐步填充未重建像素。在此過程中會產生部分飛離點,以深度值與鄰近區(qū)域一致程度可以剔除這些飛離點。
經過多年的發(fā)展,傾斜攝影及三維重建技術的理論水平和應用能力隨著研究的深入不斷提升。攝像器材、數(shù)字成像技術、三維重建算法和計算機性能的發(fā)展使得三維重建更加容易實現(xiàn),其成圖精度更高、速度更快、成圖規(guī)模更大。
原始數(shù)據(jù)指傾斜攝影作業(yè)采集到的圖像數(shù)據(jù)、POS(position and orientation system)數(shù)據(jù)(位置、姿態(tài))和相機參數(shù)。圖像由機載相機拍攝,POS數(shù)據(jù)由機載傳感器記錄,相機參數(shù)(徑向畸變、成像中心點等)由廠商標定。
圖像數(shù)據(jù)使用無人機機載相機獲取,采集的圖像應當清晰、重疊度高、視角豐富。傾斜攝影作業(yè)前應擬定好飛行線路、攝像重疊度、拍攝數(shù)量等參數(shù),選擇光照充足、少云正光的時段進行飛行作業(yè),避免云朵移動造成被攝物體陰陽變化及斜光照射產生的大背光區(qū)。
采集到圖像后,有必要對圖像進行篩查。一方面,檢查采集成果是否達到計劃的數(shù)量、范圍,重點關注地形起伏、細節(jié)復雜、光影交錯的區(qū)域,若不滿足要求,應重飛補拍;另一方面,篩查圖像質量,防止低質量圖像增大重建誤差,應舍棄失焦、雜物遮擋的圖像,有雜點、移動物的圖像可作蒙蓋處理。
圖像集經過預處理后,首先進行空中三角測量計算??罩腥菧y量通過圖像間的像素特征匹配生成大量連接點,利用連接點的匹配信息計算獲得連接點的三維位置(稀疏點云)。POS數(shù)據(jù)和圖像畸變的誤差會隨迭代得到修正。然后,采用多視影像密集匹配技術,利用像素網(wǎng)格匹配對稀疏點云進行加密計算,生成稠密點云。由于稠密點云數(shù)據(jù)量極大,通常采用切塊降低計算機的壓力。切塊也有利于計算機群的平行處理。隨后,稠密點云經過優(yōu)化處理降低數(shù)據(jù)冗余,通過三角網(wǎng)格化最終獲得不規(guī)則三角網(wǎng)模型。利用不規(guī)則三角網(wǎng)模型與像素間的匹配關系,將原始圖像紋理貼附到不規(guī)則三角網(wǎng)模型上,可得到帶紋理的三維真實感模型。
傾斜攝影三維重建的基本成果是三維網(wǎng)格模型(3D mesh),三維網(wǎng)格模型由不規(guī)則三角網(wǎng)構成。對三維網(wǎng)格模型進行光滑處理可生成不同精度的數(shù)字高程模型(DEM);三維網(wǎng)格模型經材質附著可得到數(shù)字地表模型(DSM);數(shù)字地表模型經正射投影可生成數(shù)字正射圖像(DOM)。
這些成果可以被運用到廣泛的領域,數(shù)字高程模型可用于地理信息系統(tǒng)(GIS)、正射投影影像可用于測繪成圖(Mapping)、數(shù)字地表模型可用于建筑信息模型(BIM)。
金沙江是長江上的重要支流,位于青藏高原東緣。金沙江落差大、流速快,兩岸邊坡陡峻,周圍地形溝壑交錯,山頂與谷底的落差可達1 km。頻繁的地質活動和強烈的風化作用使沿線的滑坡、崩塌事件非常多,由此引發(fā)的堰塞災害是沿線居民生產生活安全的重大威脅。圖3為衛(wèi)星記錄的白格堰塞的初期征兆、堰塞規(guī)模、潰后狀態(tài)等光學影像信息。
圖3 白格堰塞事件全階段衛(wèi)星光學圖像
2018年10月11日,第1次滑坡發(fā)生,滑坡體積2 400萬m2、堆積高度85 m;1天后,堰塞湖蓄水達2.9億m3,堰塞體自潰留下左岸殘積體[31-32];11月3日,第1次滑坡區(qū)域原位發(fā)生第2次滑坡,下滑土石350萬m2,沖壓第1次滑坡殘積體形成135 m高的堰塞體,數(shù)天內水位快速上漲,搶險人員采取人工開挖導流槽的方式主動降低壩體高度;11月12日堰塞湖蓄水達5.24億m3,水位達到人工導流槽高程2 956.40 m并開始泄流,3天后金沙江水位基本穩(wěn)定,險情解除。這起連續(xù)滑坡堰塞事件造成上游江達縣波羅鄉(xiāng)、白玉縣金沙鄉(xiāng)等先后被淹;導致下游318國道金沙江大橋損毀嚴重,云南迪慶藏族自治州、麗江市部分道路中斷、學校被淹、農田被沖毀。
白格堰塞事件發(fā)生突然、滑坡體量大、周邊環(huán)境原始,造成的經濟損失巨大、威脅人口眾多,具備大型堰塞災害的典型特征,對研究無人機傾斜攝影技術在堰塞現(xiàn)場的適用性有很高的參考價值。
采用無人機對堰塞體展開了傾斜攝影作業(yè),使用的無人機型號是DJI Phantom 4 RTK,該設備質量1.3 kg、單次續(xù)航能力約30 min,綜合GPS、北斗和伽利略全球衛(wèi)星導航系統(tǒng)(global navigation satellite system,GNSS),機載相機型號為DJI FC6310R(1 英寸 CMOS;5 472像素×3 648像素分辨率;3軸云臺穩(wěn)定系統(tǒng))。
如圖4所示,無人機飛行兩個架次分別重點采集右岸泄流槽和左岸殘積體數(shù)據(jù)。為保證三維數(shù)據(jù)生產成果質量,航線采用大重疊率:航向重疊率為80%,旁向重疊率為70%。這兩架次飛行共獲取1 497幅JPEG格式圖像(29.9 G像素)及每張圖像的位置信息、姿態(tài)信息。如圖5所示(圖中箭頭指向為水流方向),圖像涵蓋了堰塞區(qū)的重點部位,如泄流槽右岸邊坡、河谷左岸殘積體、河谷右岸滑坡面等。
圖4 白格堰塞體無人機傾斜攝影作業(yè)
圖5 白格堰塞區(qū)現(xiàn)場無人機航拍影像
三維重建需要圖像數(shù)據(jù)、POS數(shù)據(jù)(位置、姿態(tài))、相機參數(shù)。完成現(xiàn)場采集后,導出無人機內置GNSS中記錄的位置數(shù)據(jù)、云臺穩(wěn)定器中記錄的姿態(tài)數(shù)據(jù)以及相機拍攝的圖像數(shù)據(jù)。
無人機GNSS記錄的位置數(shù)據(jù)存在一定誤差,這部分誤差可在空中三角測量計算過程中由光束平差法迭代優(yōu)化,調整后的位置數(shù)據(jù)精度得到提高。姿態(tài)信息(俯仰角、橫滾角、偏航角)由云臺穩(wěn)定器記錄,云臺穩(wěn)定器的硬件性能保證了角度抖動量低于0.02°,具備高精度。相機參數(shù)在無人機出廠時經過校正,其徑向畸變、橫向畸變、焦距、成像中心點等參數(shù)較準確。
圖像經過初步處理后,進行空中三角測量計算。以GNSS定位最大誤差為調整限值對位置數(shù)據(jù)采取平差調整策略,姿態(tài)數(shù)據(jù)采用云臺穩(wěn)定器記錄數(shù)據(jù)原值,相機參數(shù)采用了出廠校正值。如圖6所示,空中三角測量計算總計生成442 947個連接點(平均每張圖像包含1 223個連接點),覆蓋面積為417萬m2,像素點的地面分辨率平均值為5.79 cm,重建均方根誤差(RMS)為0.62個像素點。
圖6 白格堰塞區(qū)傾斜攝影三維重建
完成空中三角測量計算后,進行三維場景重建。白格堰塞區(qū)的場景范圍較大,以堰塞體為建模中心(31.083 8°N,98.715°E)設定建模區(qū)域(1 796 m×2 322 m×573 m)。三維場景重建對計算機內存的需求較高(310 GB),因此采用自適應切塊方法。自適應切塊是以稀疏點云密度分度預估稠密點云的密度分布,將高密區(qū)域分成小瓦片,從而保證單個瓦片的內存需求在硬件容量內。三維重建成果如圖7所示。其中數(shù)字高程模型描述的是堰塞區(qū)地面高程信息,可用于堰塞體覆蓋面積、體積、堆積高度、坡度等特征的計算,作為現(xiàn)場作業(yè)規(guī)劃的直接依據(jù);數(shù)字地表模型是在數(shù)字高程模型的基礎上,進一步涵蓋了除地面以外的其他地表信息的高程,豐富的信息用于判別堰塞體的結構分布、材料類型、顆粒級配等,可為堰塞體穩(wěn)定評估提供參考。
圖7 白格堰塞區(qū)傾斜攝影建模三維成果
本節(jié)三維場景重建處理所采用的工具為三維圖像建模軟件ContextCapture和PhotoScan。
無人機傾斜攝影獲取現(xiàn)場圖像數(shù)據(jù)、POS數(shù)據(jù),經過空中三角測量計算、三維場景重建,獲得堰塞區(qū)三維網(wǎng)格模型。以此為基礎可生成白格堰塞區(qū)的數(shù)字高程模型(DEM)、數(shù)字地表模型(DSM)、數(shù)字正射圖像(DOM)等形式的成果如圖8所示。數(shù)字正射圖像是利用數(shù)字高程模型對無人機圖像進行投影差改正、鑲嵌、剪裁生成的數(shù)字正射影像數(shù)據(jù)集,具有精度高、信息豐富、直觀真實等特點,可作為減災規(guī)劃決策的參考依據(jù)(圖8(a))。
圖8 白格堰塞區(qū)傾斜攝影建模平面成果
結合傾斜攝影重建成果,在數(shù)字三維模型中進行坐標讀取、距離量測、坡比計算等工作,獲得了堰塞體形態(tài)、結構等方面的信息,為救災減災提供服務。
人工泄流槽初始尺寸為:坎底高程2 952.52 m、長220 m、頂寬42 m、底寬3 m、最大深度15 m[32]。經測量,沖刷擴大后泄流槽左岸邊坡角度為35°~38°、右岸邊坡角度為32°~36°,與初始邊坡角度35°基本一致;槽底高程約為2 865.50 m,較初始降低85 m;槽底寬度約為80 m,約為初始底寬的26倍;槽底長度約720 m,約為初始底長的3.5倍;槽頂寬度約為300 m,約為初始頂寬的7倍。
堆積體最高點位于左岸邊坡(高程2 974.93 m),距泄流槽左坡頂距離約160 m。以最大縱剖面為基準,距上游坡腳水平距離約230 m(23°~25°),距下游坡腳水平距離約500 m(11°~14°)。
a. 堰塞區(qū)具有范圍寬、地表起伏大等特點,無人機分批次變高度飛行可覆蓋堰塞區(qū)的重點部位,如泄流槽邊坡、殘積體、滑坡面等。同時,采用大重疊率的航線規(guī)劃可保證無人機影像間的有效匹配,從而保證三維重建的高質量。
b. 白格堰塞區(qū)三維重建計算表明,大量的高分辨率影像數(shù)據(jù)會產生較大的計算壓力,尤其占用計算機內存資源。以堰塞體中心設定建模區(qū)域,采用切塊分批計算策略可有效緩解計算壓力。
c. 測算白格堰塞區(qū)傾斜攝影重建模型得出:殘積體最高點高程2 974.93 m,上游平均坡度為24°、下游平均坡度13°。人工泄流槽經沖刷擴大后,其側坡坡度與初始坡度基本一致,均約為35°,其槽底下切至2 865.50 m高程、寬底擴至80 m,槽底長度約720 m。