李 華，寧化展，丁希旺

（濱州市水利勘測設(shè)計研究院測量隊，山東 濱州256612）

0 引 言

數(shù)字地形模型（DTM）是以離散分布的平面點(diǎn)來模擬連續(xù)分布的地形，是野外地表勘測成果的數(shù)字化展現(xiàn)，廣泛的應(yīng)用于地理信息系統(tǒng)各領(lǐng)域中。Delaunay三角網(wǎng)是DTM的主要實(shí)現(xiàn)形式，用一系列互不交叉、重疊的連接在一起的三角形網(wǎng)來表示地形［1］。

礦產(chǎn)儲量，簡稱儲量，一般是指經(jīng)過地質(zhì)勘探、基建勘探和生產(chǎn)勘探后，經(jīng)勘查證實(shí)存在礦體，探明其空間分布、產(chǎn)狀、形態(tài)、規(guī)模和質(zhì)量，可以開發(fā)利用的原地礦產(chǎn)資源量［2］。根據(jù)各種探礦工程和技術(shù)手段所得到的資料及信息（一般是鉆孔資料），通過一定的計算方法計算礦產(chǎn)的地下埋藏量，這一系列的工作稱之為礦產(chǎn)儲量（資源量）計算。

礦產(chǎn)資源儲量的計算方法有斷面法、算術(shù)平均法、等值線法等，這些算法都是在礦產(chǎn)資源勘探時期估算礦產(chǎn)資源儲量的計算方法，這些算法并不適于計算礦體采空區(qū)的開采量。為此，提出了一種基于礦體區(qū)域的采空區(qū)礦體儲量計算方法，根據(jù)礦體采空區(qū)范圍及回采率計算出礦體資源已開采量和礦體資源剩余儲量。

1 礦體區(qū)域Delaunay三角網(wǎng)的生成

基本Delaunay三角網(wǎng)生成有以下四個步驟：凸殼的生成過程、初始三角網(wǎng)生成過程、離散點(diǎn)插入初始三角網(wǎng)過程、局部優(yōu)化過程。礦產(chǎn)資源的儲量計算是在已知礦產(chǎn)資源分布邊界，以及區(qū)域鉆孔資料數(shù)據(jù)計算所得，根據(jù)鉆孔的稀疏程度還要對離散鉆孔點(diǎn)進(jìn)行插值運(yùn)算。因此論述的Delaunay三角網(wǎng)生成有以下六個步驟：

1）點(diǎn)差值過程：由于鉆孔點(diǎn)不多以及區(qū)域邊界太長，為了使生成的三角網(wǎng)最優(yōu)及礦產(chǎn)量的計算更準(zhǔn)確需要對邊界點(diǎn)及鉆孔點(diǎn)進(jìn)行插值運(yùn)算。

2）凸殼生成過程：凸殼［3］是包含所有離散點(diǎn)的最小多邊形，一般的凸殼構(gòu)造方法只是找出了最少點(diǎn)的多邊形，對于多點(diǎn)在一條凸殼邊上的情況只是找出了兩個端點(diǎn)而中間點(diǎn)卻沒有找出來，這種情況對礦體的建模有影響，為此采用一種夾角與距離最?。?］的凸殼構(gòu)造方法。

3）初始三角網(wǎng)生成過程：利用凸殼點(diǎn)生成初始三角網(wǎng)（如圖1所示）。

4）點(diǎn)插入過程：將點(diǎn)集中的離散點(diǎn)（不包含凸殼點(diǎn)）依次插入［5］，在初始三角形鏈表中找出三角形外接圓包含插入點(diǎn)的三角形（稱為該點(diǎn)的影響三角形），剔除查找出的影響三角形的公共邊，將此插入點(diǎn)同查找到的影響三角形的全部頂點(diǎn)連接起來，完成一個點(diǎn)在Delaunay三角形鏈表中的插入，如圖2所示。

5）局部優(yōu)化過程：局部優(yōu)化算法是為了生成Delaunay三角網(wǎng)。算法的基本含義：對由兩個公共邊組成的四邊形進(jìn)行判斷，如果其中一個三角形的外接圓包含第四個頂點(diǎn)，則這個四邊形的對角線互換，如圖3所示。

圖1 初始三角網(wǎng)

圖2 點(diǎn)插入

圖3 局部有化過程

6）三角形剔除過程：根據(jù)已知邊界條件，將不符合條件的三角形剔除，如圖4所示。若只是按照三角形三個點(diǎn)都在邊界上這個條件來剔除可能會剔除一些合法的三角形，因此還要加上內(nèi)切圓圓心是否在邊界內(nèi)這個條件。算法描述如下：①判斷三角形的三個頂點(diǎn)是否位于邊界上，②若果不是都在邊界上刪除三角形，如果均位于邊界上，求出內(nèi)切圓的圓心，判斷是否位于邊界內(nèi)，若是三角形保留否則三角形剔除［4］，如圖4所示。

圖4 三角形剔除

2 采空區(qū)域Delaunay三角網(wǎng)的生成

礦山儲量動態(tài)管理中，遇到的主要問題是兩個三維實(shí)體：采空區(qū)與礦體［6］。需要管理的內(nèi)容是二者相交部分的礦產(chǎn)量，此礦產(chǎn)量的計算首先要求采空區(qū)的平面區(qū)域與礦體Delaunay三角網(wǎng)區(qū)域的相交Delaunay三角網(wǎng)區(qū)域，如圖5所示，再求出此Delaunay三角網(wǎng)區(qū)域的礦產(chǎn)量。

圖5 采空區(qū)域

根據(jù)采空區(qū)平面二維邊界求其與礦體表面Delaunay三角網(wǎng)相交部分的Delaunay三角網(wǎng)的算法思想如下：

1）如圖5所示，內(nèi)邊界點(diǎn)也要參與Delaunay三角網(wǎng)的構(gòu)建；

2）以內(nèi)邊界的一條邊開始查找此邊影響的三角形（即：與此邊相交的三角形）并依次對查找的三角形編號記錄下來，判斷三角形在邊界外側(cè)的頂點(diǎn)與此三角形跟影響邊的交點(diǎn)組成的是三角形還是四邊形，若是三角形直接連接，若是四邊形根據(jù)Delaunay三角形優(yōu)化原則連接四邊形的一條對角線構(gòu)成兩個三角形，如圖6所示。

影響三角形的查找方法如下：以起始邊的一個端點(diǎn)開始，首先查找以此端點(diǎn)為頂點(diǎn)的三角形并判斷那個三角形與此邊相交求出交點(diǎn)，如此得到第一個影響三角形；找出與第一個三角形影響邊共邊的三角形并求出此三角形與起始邊的交點(diǎn)，如此得到第二個影響三角形，依次循環(huán)直到求出的交點(diǎn)為起始邊的另一端點(diǎn)，這樣就找出了所有影響三角形及它們與影響邊的交點(diǎn)；

圖6 邊界影響三角形

3）依照步驟2）求出其它邊界邊的影響三角形和每個影響三角形與邊界邊的交點(diǎn)，并按步驟2）的方法生成邊界外側(cè)的三角網(wǎng)；

4）邊界點(diǎn)和邊界上的交點(diǎn)與邊界內(nèi)的點(diǎn)重新生成Delaunay三角網(wǎng)，如圖7所示。

圖7 三角網(wǎng)生成

3 基于Delaunay三角網(wǎng)的儲量計算方法

根據(jù)Delaunay三角網(wǎng)模型，直接利用Delaunay三角形進(jìn)行體積計算，體積計算如下：先求出每個Delaunay三角形所構(gòu)成立柱體的體積Vi，再把區(qū)域內(nèi)所有立柱體的體積相加，從而得到整個礦區(qū)內(nèi)礦產(chǎn)體積儲量V，再乘以礦產(chǎn)密度，就是礦體重。因此，地質(zhì)體的體積計算也就是計算每個三棱柱［7］，如圖8所示的體積，然后累積和就是地質(zhì)體的體積。

圖8 礦體線性模型

已知三棱柱上表面三角形三個點(diǎn)A，B，C和下表面三角形三個點(diǎn)D，E，F(xiàn)的（X，Y，Z）坐標(biāo)，求三棱柱的體積算法思想如下：

1）找出上表面三角形點(diǎn)Y值最小的點(diǎn)A和下表面三角形點(diǎn)Y值最大的點(diǎn)E，并以這兩點(diǎn)為基準(zhǔn)做橫切面把三棱柱分成一個直三棱柱和兩個四棱錐體。

2）根據(jù)A，B，C三點(diǎn)的坐標(biāo)值可求出BB1、CC1、C1B1、AC1、AB1的長度，進(jìn)而求出四邊形CC1B1B的面積S.在三角形AB1C1中根據(jù)求出的三邊長可求出A到邊B1C1的高h(yuǎn)，即四棱錐的高。上表面四棱錐的體積v＝S＊h／3.

3）A，E兩點(diǎn)的高程差即是直三棱柱的H，根據(jù)第二步求出的三角形AB1C1三個邊長可求出直三棱柱的橫截面面積S，則直三棱柱的體積V＝S＊H.

4）按照步驟②的方法求出下表面四棱錐的體積，三個體積累計求和即為此三棱柱的體積，多個三棱柱體積的累計即得到地質(zhì)體的體積。

圖9 礦體面狀模型

圖9、10示出了以滕州曹莊煤礦三層煤數(shù)據(jù)及采空區(qū)數(shù)據(jù)生成的煤層模型和采空區(qū)模型。選中煤層或任一采空區(qū)都可計算出其體積、表面積及投影面積。體積及面積信息在窗體狀態(tài)欄左下角顯示，如圖11所示。

圖10所示淺色區(qū)域?yàn)檫x中的礦體采空區(qū)區(qū)域，其體積為2 894 621.0m3與ARCGIS計算的體積2 888 418.6m3相似率達(dá)0.997 8，實(shí)踐證明此算法是正確有效的。

圖10 礦體面狀模型

圖11 狀態(tài)欄

4 結(jié) 論

根據(jù)鉆孔的坐標(biāo)及屬性數(shù)據(jù)庫，利用java3D可視化技術(shù)，通過對Delaunay三角形生成過程中點(diǎn)集進(jìn)行封閉性判斷，生成任意地區(qū)Delaunay三角網(wǎng)，結(jié)合鉆孔屬性數(shù)據(jù)庫，對任意區(qū)域礦產(chǎn)儲量進(jìn)行計算，為計算任意區(qū)域礦產(chǎn)儲量與可視化提供了一種新思路，該方法不僅克服算術(shù)平均法誤差大及其它方法難以自動實(shí)現(xiàn)的缺點(diǎn)，而且為礦山管理提供可視化查詢與分析功能。

［1］ 袁 翰，李偉波，陳婷婷.對構(gòu)建Delaunay三角網(wǎng)中凸殼算法的研究與改進(jìn)［J］.計算機(jī)工程，2007，33（7）：70－72.

［2］ 蔣 銳，宋煥斌，朱杰勇.基于三維柵格數(shù)據(jù)的露天礦山儲量動態(tài)監(jiān)測研究與應(yīng)用［J］.礦產(chǎn)與地質(zhì)，2009，23（5）：469－472.

［3］ 劉永和，王潤懷，齊永安.一種非凸包邊界約束不規(guī)則三角網(wǎng)生成算法［J］.測繪科學(xué)，2008，33（3）：79－81.

［4］ 寧化展，田茂義.帶有邊界條件的Delaunay三角網(wǎng)生成算法的研究與實(shí)現(xiàn)［J］.全球定位系統(tǒng)，2010，35（4）：39－42.

［5］ 徐道柱，劉海硯.大量約束邊條件下Delaunay三角網(wǎng)的快速生成［J］.測繪工程，2007，16（3）：6－10.

［6］ 張尹軍.礦產(chǎn)資源儲量動態(tài)監(jiān)測的問題及建議［J］.商品與質(zhì)量，2010，（7）：74－75.

［7］ 張渭軍，孔金玲，翁曉鵬.基于鉆孔數(shù)據(jù)Delaunay三角網(wǎng)的礦產(chǎn)儲量計算［J］.測繪科學(xué)，2010，35（4）：57－58.