利用航測數(shù)據(jù)進行礦山生態(tài)修復土石方量計算模型的研究

張 鴿,垢元培

(河北省地質礦產(chǎn)勘查開發(fā)局第二地質大隊(河北省礦山環(huán)境修復治理技術中心),河北 唐山 063000)

我國是礦產(chǎn)資源大國,礦山資源促進了經(jīng)濟的高速發(fā)展,同時礦產(chǎn)資源的開發(fā)利用也造成了嚴重的生態(tài)破壞和環(huán)境污染,成為制約我國經(jīng)濟、社會長遠發(fā)展的重要因素。因此,開展礦山生態(tài)修復顯得尤為重要。在對礦山特別是露天開采礦山進行生態(tài)修復時,要對其形成的坡面進行不同程度的處理,在保證邊坡穩(wěn)定性的前提下,進行生態(tài)環(huán)境修復措施[1]。礦山地質環(huán)境治理和礦山生態(tài)修復工程過程中涉及了大量的土石方量計算工作。土石方量計算結果的精確度不僅影響工程的投資預算和設計方案的優(yōu)化,同時也影響到業(yè)主單位與施工單位的利益[2],以及工程的結算和施工方案的實施程度,是工程量計量的重要環(huán)節(jié)。因此,如何得到準確的土石方量計算成果顯得格外重要。

對土石方量計算的研究大多為計算方法、計算原理或計算步驟,而少有從計算模型的角度進行研究。文獻[3—4]研究了計算方法適用的條件、適用范圍和精度要求等。文獻[5—8]研究了土石方量計算原理、方法與步驟,介紹了幾種主要應用軟件、TIN模型模擬地面特征及DEM數(shù)據(jù)進行土方量計算的基本步驟和注意事項。

隨著無人機低空攝影測量技術的快速發(fā)展,在礦山生態(tài)修復中無人機的利用也得到了空前發(fā)展。無人機攝影測量的優(yōu)點是作業(yè)靈活,可快速獲取高分辨率數(shù)字影像,同時受現(xiàn)場地質地理環(huán)境影響較小,工作中既提高了作業(yè)效率,又節(jié)約了大量的人力物力,大幅度提高了測量人員的安全保障[9]。

本文以無人機航測數(shù)據(jù)中的點云數(shù)據(jù)為基礎,利用散點法的優(yōu)勢,以地形性質線為特征線,進行不規(guī)則三角網(wǎng)模型的建立,實現(xiàn)礦山生態(tài)修復中土石方量計算模型的構建;并通過與DEM數(shù)據(jù)計算的結果進行對比分析完成精度驗證,從而快速、準確地得到土石方量計算成果。

1 土石方量計算數(shù)據(jù)類型

土石方量計算常用軟件主要有南方CASS軟件、飛時達軟件、GIS(ArcGIS)軟件和土石方大師軟件4種。其中,南方CASS軟件和飛時達軟件均以CAD軟件二次開發(fā)為基準的平臺,均可采用三角網(wǎng)模型計算;GIS(ArcGIS)軟件以DEM數(shù)據(jù)進行計算;土石方大師軟件為國內(nèi)創(chuàng)新團隊自主研發(fā)軟件,采用三角網(wǎng)模型計算,可實現(xiàn)三維化。這幾款軟件均可適用于礦山生態(tài)修復中土石方量計算。本文主要研究以無人機航測點云數(shù)據(jù)為基礎,通過三角網(wǎng)法建立模型。

CASS軟件系統(tǒng)是較常用的測繪軟件,在土石方量計算中提供了方格網(wǎng)法、DTM法、等高線法和斷面法等豐富的土方計算方法[10],對不同的工程和地形條件可靈活采用合適的土石方計算方法。而在礦山生態(tài)環(huán)境修復項目中主要采用DTM法建立土石方量計算模型。FastTFT(飛時達)是一款基于CAD平臺二次開發(fā)的專業(yè)土石方計算繪圖軟件[11],針對各種復雜地形情況和工程實際要求,提供了方格網(wǎng)法、三角網(wǎng)法、斷面法、道路斷面法、田塊法、整體估算法6種土石方量計算方法,同時對于土方挖填量的結果可進行分區(qū)域調配優(yōu)化。GIS軟件在土石方量計算中使用最多的為ArcGIS軟件,ArcGIS 軟件具有土石方量計算的工具。其計算原理與三角網(wǎng)法原理類似[12],采用克里金插值法,利用TIN和DEM數(shù)據(jù)通過ArcGIS軟件多樣化的功能模塊,優(yōu)化了傳統(tǒng)方法的空間關系,從而使ArcGIS軟件在土石方計算上具有優(yōu)越性。土石方大師軟件具有多圖層、實景影像結合、三維立體化展示、模型數(shù)據(jù)種類多樣、成果快速生成、智能化操作計算等優(yōu)點,可計算多種三維模型的方量與面積等成果,為土石方量計算提供了豐富的三維模型展示[13]。此軟件為國內(nèi)創(chuàng)新團隊自主研發(fā)土石方量計算專用軟件,模型全面、計算精度高。

土石方量計算數(shù)據(jù)的類型主要有兩大類,分別為散點法(點云法)和數(shù)字高程模型(DEM法),兩種方法的特征見表1。

表1 土石方量計算數(shù)據(jù)類型特征

無人機點云數(shù)據(jù)是通過無人機航測手段,獲取的具有位置信息和顏色信息的三維坐標數(shù)據(jù)集合。日常計算中通常對點云數(shù)據(jù)進行稀釋后使用,在測量工作中主要采用三角網(wǎng)法。

數(shù)字高程模型(DEM)是以數(shù)字形式按一定的結構組織在一起的,能夠表示實際地形特征、空間分布的一種實體地面模型,也是地形地理位置和地形起伏的數(shù)字描述[14],可用于土石方量的計算。DEM數(shù)據(jù)與無人機點云數(shù)據(jù)如圖1所示。

圖1 同一區(qū)域DEM數(shù)據(jù)與點云數(shù)據(jù)對照

2 計算模型建立方法的研究

無人機航測得到的點云一般用顏色表示高程,以反映地形。但在日常工作中,均需要根據(jù)精度進行點云抽稀,而抽稀后的點云數(shù)據(jù)為“規(guī)則無規(guī)律”的三維數(shù)據(jù),無法真實反映復雜地形,具有局限性、容易造成數(shù)據(jù)處理錯誤、費時費力、受人為因素影響大、計算誤差大等弊端。

應用點云數(shù)據(jù)計算土石方量時需要使用地性線。地性線也稱地形性質線[15],是地貌形態(tài)的骨架線,是描述地貌形態(tài)時的控制線,這條線在實際地形中是不存在的,是在內(nèi)業(yè)數(shù)據(jù)制圖中得到的派生數(shù)據(jù)。

2.1 土石方計算模型的建立

采用Delaunay三角剖分與地性線邊界虛擬點進行高程點云的有限元分析(finite element analysis,FEA)。應用此方法經(jīng)三角剖分后含有最大化最小角特征及最重要的空圓特性,可生成最接近于規(guī)則化的三角形的三角網(wǎng),同時具有點云中任意的4個點不共圓的唯一特性,這種方法明顯優(yōu)于常用的四叉樹和八叉樹方法及網(wǎng)格前沿算法。而地形性質線的使用有效彌補了航測點云的“規(guī)則無規(guī)律”性。將Delaunay方法與地形性質線相結合,彌補了因復雜地形造成的計算誤差、計算錯誤等問題,尤其是在礦山生態(tài)修復項目等高山地地形。

空圓特性是指在一個三角形的外接圓范圍內(nèi),包括外接圓范圍線上,除了三角形的3個頂點,不再包含點云數(shù)據(jù)中的任何其他點。其中,四點共圓屬于一個特例,空圓特性如圖2所示。

圖2 空圓特性示意

最大化最小角特性在本文中主要指,經(jīng)抽稀后的無人機點云數(shù)據(jù)所組成的三角形中兩個相鄰的所組成凸四邊形的對角線。最大化最小角特性如圖3所示。

圖3 最大化最小角特性示意

由于必須使用到地形性質線,應先將地形性質線按抽稀后的點云密度進行同間距或小于點云密度的間距插點,也可以稱為地形性質線上的分割點。由于此分割點沒有高程信息,需先根據(jù)地形線分割開的周圍點云數(shù)據(jù)進行高程賦值。而三角網(wǎng)的建立采用“插值法”,即由外向內(nèi)的方法。首先建立一個覆蓋全部點云的三角形;然后將點云數(shù)據(jù)進行逐一插入計算,從而逐一找出符合空圓特性和最大化最小角特性的三角形;最后建立起不規(guī)則三角網(wǎng)數(shù)據(jù)模型。假設插入一點后的外切圓上的3個點坐標為A(x1,y1)、B(x2,y2)、C(x3,y3),則其外切圓的圓心(x0,y0)和半徑(r)計算公式分別為

(1)

(2)

(3)

2.2 模型體積計算

在三角網(wǎng)數(shù)據(jù)中逐一選出三角形,其頂點坐標為A(x1,y1,h3)、B(x2,y2,h3)、C(x3,y3,h3),體積計算公式為

(4)

(5)

V=Sh

(6)

式中,S表示三角形面積;h表示高程;h0表示參考高程;V表示體積。

根據(jù)h的正負可以計算出挖方量和填方量,挖方量與填方量的和為總體積,最后將所有體積相加得到總方量。

3 礦山生態(tài)修復中土石方工程類型

根據(jù)工作實際和礦山生態(tài)修復工程的需要,可以將礦山生態(tài)修復中土石方工程按照礦山類型、礦山地形、破壞程度、修復方法、工作方法、施工流程等分為山體平行挖方或填方工程、依坡挖方或填方工程及山體挖削臺階工程。以上3種工程類型基本涵蓋了礦山生態(tài)修復的不全類型,可真實、直觀反映礦山修復后的地形地貌。

3.1 山體平行挖方或填方工程

山體平行挖方或填方工程指對獨立山體進行挖方或其上平臺進行填方,其工程部位不依附于周圍山體邊坡或地形,施工后對周圍工程不產(chǎn)生影響,如圖4所示。

圖4 山體平行挖方前、后對照

3.2 依坡挖削方或填方工程

依坡挖削方或填方工程指進行工程部位一側或四周依附于邊坡或山體,當進行土石方量計算模型設計時需考慮邊坡的影響。利用地形性質線劃分出邊坡與平臺范圍,形成依坡挖方或填方工程,施工后會增加或減少邊坡面積,如圖5所示。

圖5 依坡挖削方前、后對照

3.3 山體挖削臺階工程

山體挖削臺階工程指對開采過后的礦山破碎邊坡進行削坡降段,在模型建立過程中需使用地形性質線將挖削方形成的臺階與邊坡進行分割,最終形成逐層臺階與邊坡,如圖6所示。

圖6 山體挖削臺階前、后對照

4 礦山生態(tài)修復中土石方量計算模型

利用航測數(shù)據(jù)進行礦山生態(tài)修復土石方量計算模型的研究,點云數(shù)據(jù)是基礎,而地形性質線(如平臺邊界線、坎頂線、坎底線、邊坡邊界線、臺階與邊坡處分界線等)對土石方計算起到了非常關鍵的作用。

依據(jù)礦山生態(tài)修復中土石方工程類型和點云數(shù)據(jù)計算的特點,可將土石方量計算模型分成3類,分別為山體平行挖方或填方模型、依坡挖方或填方模型以及山體挖削臺階模型。模型建立流程如圖7所示。

圖7 模型建立流程

以唐山市某礦山生態(tài)修復治理項目中土石方量計算數(shù)據(jù)為例。本文項目航飛任務采用大疆Phantom 4 RTK無人機進行航測數(shù)據(jù)采集,采用ContextCapture軟件進行空三數(shù)據(jù)處理,采用Cloud Compare或Global Mapper軟件將點云數(shù)據(jù)處理成5 m間距,采用土石方大師軟件進行點云數(shù)據(jù)土石方三維建模與方量計算,采用ArcGIS軟件進行DEM數(shù)據(jù)的土石方量計算,最終通過對比分析得到點云數(shù)據(jù)計算的精度分析。

4.1 山體平行挖方或填方模型

山體平行挖方或填方模型主要適用于區(qū)域性整體性挖削方工程或區(qū)域性平臺進行填方工程。此模型適用性最為普遍、通用性最強,由于此模型不依附任何地形,因此外部范圍線即可作為地形性質線。

山體平行挖方或填方模型的建立方法是以挖填區(qū)域外部范圍為地形性質線,以外圈地形性質線和點云為基礎建立三角網(wǎng)模型,如圖8所示,再根據(jù)三角網(wǎng)法計算兩期間土石方量。

圖8 山體平行挖方或填方模型

通過對實例數(shù)據(jù)的計算,使用點云法,采用山體平行挖方或填方模型所計算的土石方挖方量為25 288.1 m3。

使用ArcGIS軟件對相同區(qū)域DEM數(shù)據(jù)進行土石方量計算作為“山體平行挖方或填方模型”計算方法的精度驗證,如圖9所示,計算結果為25 383.1 m3。

圖9 同區(qū)域DEM數(shù)據(jù)模型

4.2 依坡挖方或填方模型

依坡挖方或填方模型主要適用于依附大型山體挖削方工程、基坑填方或依附山體的平臺堆坡工程。此模型選擇繪制地形性質線最為重要,點云間距越大地形性質線對計算的方量精度的影響越大。

依坡挖方或填方模型的建立方法是以挖填區(qū)域外部范圍和坡底(頂)為地形性質線,以外圈地形性質線、坡底(頂)地形性質線和點云為基礎建立三角網(wǎng)模型,如圖10所示,再根據(jù)三角網(wǎng)法計算兩期間土石方。

圖10 依坡挖方或填方模型

通過對實例數(shù)據(jù)的計算,使用點云法,采用依坡挖方或填方模型所計算的土石方挖方量為31 978.7 m3。

使用ArcGIS軟件對相同區(qū)域DEM數(shù)據(jù)進行土石方量計算作為“依坡挖方或填方模型”計算方法的精度驗證,如圖11所示,計算結果為32 026.4 m3。

圖11 同區(qū)域DEM數(shù)據(jù)模型

4.3 山體挖削臺階模型

山體挖削臺階模型主要適用于依附山體的挖削方后形成逐層臺階的生態(tài)修復工程,也是土石方量計算模型中最復雜的一種。同時在礦山生態(tài)修復中的作用也較重要,通過點云數(shù)據(jù)與地形性質線的有效繪制,可以實現(xiàn)對臺階模型的精確建模。此模型與“山體平行挖方或填方模型”“依坡挖方或填方模型”的相同點為外部范圍線均為地形性質線。

山體挖削臺階模型的建立方法是以挖填區(qū)域外部范圍和臺階內(nèi)側、臺階外側為地形性質線,以外圈地形性質線、臺階內(nèi)測地形性質線、臺階外側地形性質線和點云為基礎建立三角網(wǎng)模型,如圖12所示,再根據(jù)三角網(wǎng)法計算兩期間土石方量。

圖12 山體挖削臺階模型

通過對實例數(shù)據(jù)的計算,使用點云法,采用山體挖削臺階模型所計算的土石方挖方量為28 290.4 m3。

使用ArcGIS軟件對相同區(qū)域DEM數(shù)據(jù)進行土石方量計算作為“山體挖削臺階模型”計算方法的精度驗證,如圖13所示,計算結果為28 065.5 m3。

圖13 同區(qū)域DEM數(shù)據(jù)模型

5 模型精度分析

通過建立山體平行挖方或填方、依坡挖方或填方和山體挖削臺階3種模型進行土石方量計算。建立模型的點云數(shù)據(jù)為采用Cloud Compare或Global Mapper軟件抽稀后所形成的5 m點云數(shù)據(jù)。無人機航測無法直接得到DEM數(shù)據(jù),因此當?shù)乇碇脖缓徒嬛镙^少時,可使用DSM代替DEM。精度分析使用的DEM數(shù)據(jù)為航測后ContextCapture軟件形成的DSM數(shù)據(jù),而非使用抽稀后的點云數(shù)據(jù)通過克里金法或三維分析法建立TIN數(shù)據(jù)后形成數(shù)字高程模型。通過實例數(shù)據(jù)分析,對比點云法與數(shù)字高程模型法計算的成果及誤差,見表2。

表2 礦山生態(tài)修復中土石方量計算模型精度分析 m3

表中兩種方法方量差占比計算公式為(數(shù)字高程模型法計算方量-點云法計算方量)÷點云法計算方量×100%。通過對計算結果進行分析可知,使用點云法通過對不同地形條件下模型的準確建立,與數(shù)字高程模型法進行比對分析后,山體平行挖方或填方模型的精度為0.38%,依坡挖方或填方模型的精度為0.15%,山體挖削臺階模型的精度為0.79%。對比精度可得,點云法計算方量與數(shù)字高程模型法計算方量之差的占比均小于1%,點云法計算的土石方量與數(shù)字高程模型法計算的結果精度一致,而非精度低。而點云法相比于數(shù)字高程模型法具有操作簡單、適用軟件多、要求專業(yè)技術水平較低等優(yōu)勢。因此,利用航測數(shù)據(jù)進行礦山生態(tài)修復土石方量計算模型的研究,更利于該方法的推廣。

6 結 語

礦山生態(tài)修復過程中,土石方量的計算是一項基礎性工作,也為生態(tài)修復后續(xù)工作提供了依據(jù)[16]。本文通過研究無人機航測點云數(shù)據(jù)特點,結合地形性質線,利用三角網(wǎng)建立適用于不同類型的土石方量計算模型。依據(jù)礦山生態(tài)修復中土石方工程類型和點云數(shù)據(jù)計算的特點,建立山體平行挖方或填方、依坡挖方或填方和山體挖削臺階3種模型。通過項目實例證明,在礦山生態(tài)修復中,針對不同地形條件建立合理的土石方量計算模型,采用點云數(shù)據(jù)計算的精度與采用DEM數(shù)據(jù)計算的精度相近,可快速準確地應用到土石方量計算中,有較大的推廣與應用價值。