基于GOCAD的巢湖北部區(qū)域的三維地質建模

王祥浩

(安徽理工大學 測繪學院，安徽 淮南 232001)

在地質學當中，地層的空間分布、巖體的構造關系是最為基本的信息，這些幾何信息用以揭示地質結構的發(fā)展歷史，估計區(qū)域的應力場分布。同時也因為這些地質幾何結構與礦產(chǎn)資源的分布、人造工程結構的穩(wěn)定性和地下水循環(huán)有著密切的聯(lián)系，因而它們也在礦產(chǎn)開采、工程建設及環(huán)境分析方面有重要應用，有著巨大的經(jīng)濟及社會價值。

為了有效地展現(xiàn)通過地質調(diào)查手段獲取的巖層構造數(shù)據(jù)，人們傳統(tǒng)上通常利用地質圖作為工具，它的主要方法是將構造面與地面的交線用幾何投影的辦法描繪到水平面上，或者依據(jù)野外調(diào)查數(shù)據(jù)繪制剖面圖，將這些代表著地質結構的空間上連續(xù)信息以二維形式呈現(xiàn)。地質工作者依據(jù)這些圖鑒來推測立體的地質結構。然而隨著地理信息系統(tǒng)和計算機圖形學技術的發(fā)展，出現(xiàn)了對地質結構建立三維模型的方法，該方法在強烈起伏的復雜地質構造中優(yōu)勢明顯，因為地質結構可以被直觀地觀察到，能夠幫助地質工作者更好地理解地質體在地下深處的延展情況，避免做出不符合實際情況的解釋。

從20世紀80年代之后，用計算機進行三維地質建模的這一領域不斷有新的成果涌現(xiàn)，為其在實際應用中奠定了堅實的基礎。三維地質建模的概念最早由加拿大Simon W.Houlding 于1993年[1]提出，但在此之前已經(jīng)有人對三維地質建模的理論進行探討,1991年A.B.Ekoule 解決了非凸輪廓線的三維物體重構問題[2],Raper J F、Thomas R Fihsher、Rongxing Li等人又進行了大量的研究，主要包括空間數(shù)據(jù)的模型與結構、數(shù)據(jù)的三維可視化、三維矢量化地圖的數(shù)據(jù)結構等[3-5]。在數(shù)據(jù)結構方面國內(nèi)學者也有很多成果，李青元等人研究由5組拓撲關系結構組成的3D 矢量模型[6]。陳軍等人綜合三維GIS矢量數(shù)據(jù)結構和基于表面剖分的3DCM空間數(shù)據(jù)模型提出單純性拓撲空間的數(shù)據(jù)模型[7-8]。郭達志等人基于八叉樹理論開發(fā)了礦山信息存儲結構模型[9-10]。龔建雅研究了矢量與柵格集成的三維地質模型[11]。李清泉和李德仁研究了TIN-CSG 混合數(shù)據(jù)模型[12]。侯恩科和吳立新研究了面向地質建模的三維拓撲數(shù)據(jù)模型[13]。在軟件發(fā)展方面，法國Nancy理工大學GOCAD(Geological Object Computer Aided Design)是地質領域公認的主流建模軟件，它考慮地質數(shù)據(jù)來源的多樣性，并提供多種建模方法。

本文力圖以GOCAD軟件為基礎并使用由GIS技術所管理的相關地質數(shù)據(jù)來構建巢湖北部區(qū)域的三維地質模型。

1 研究區(qū)概況及數(shù)據(jù)介紹

1.1 研究區(qū)介紹

本文所選取的建模區(qū)域位于安徽省巢湖市北部山區(qū)地帶，其范圍大致為東經(jīng)117°47′～117°54′,北緯31°36′～31°42′,其東北西三面環(huán)山，山脈走向為北東35°～40°。最高山峰大尖山海拔350 m，一般山峰海拔高100～300 m,最低處獅子口海拔高程僅20 m。巢湖市北部山區(qū)屬于楊子地層區(qū)，以古生界發(fā)育為特點，其中，尤以晚古生界和下、中三疊統(tǒng)發(fā)育最完整，該區(qū)域屬于半湯復式背斜西翼，由3個二級褶皺組成，自東向西為喻府大村向斜、鳳凰山背斜和平頂山向斜。褶皺軸面傾向北西，褶皺軸跡方向為20°～ 30°,樞紐均向SSW傾伏，傾伏角15°～26°。褶皺在平面上表現(xiàn)為一個斜體字母“M”形。褶皺多被斷裂破壞，并有小型巖漿活動[14]。

1.2 數(shù)據(jù)介紹

本文收集與建模有關的地形數(shù)據(jù)與地質數(shù)據(jù)。地形數(shù)據(jù)是以20 m為等高距的1∶5萬等高線地形圖，并利用GIS軟件將其矢量化后儲存于點線圖層中。地質數(shù)據(jù)是該地區(qū)的1∶5萬地質圖，圖上主要表明地下地層出露于地表的情況(稱為露頭數(shù)據(jù))和斷層線在地表的分布情況，他們在圖上分別以面狀與線狀展示，將其存儲在GIS軟件中。這兩類數(shù)據(jù)都是立體的，但被投影于水平面上，因此，它們都以二維的形式來呈現(xiàn)。與露頭數(shù)據(jù)相連是其產(chǎn)狀數(shù)據(jù)，包括走向與傾角，利用產(chǎn)狀數(shù)據(jù)可以估計地層或斷層面在地下的延展情況。

2 建模流程

2.1 地表形態(tài)的模擬

地表形態(tài)可以看作是地下地質結構出露于地表并被外營力所影響而形成的部分，因此，地表形態(tài)也是地下地質結構形態(tài)的重要組成部分。以現(xiàn)有的1∶5萬等高線地形圖為基礎模擬地表形態(tài)要對其做插值處理。常用的插值方法有很多，針對地質建模特點，GOCAD開發(fā)了專門的對應算法即離散光滑插值算法(Discrete Smooth Interpolation，DSI)。它適用于構建復雜模型和處理模型表面不連續(xù)的情況[15]。具體實現(xiàn)過程：首先將已經(jīng)在GIS中矢量化的等高線數(shù)據(jù)轉換為DXF格式，再將其導入GOCAD中保存在Objects下的Curve類中。在PointSet選項中使用New-From PointSet Curve or Surface命令將等高線離散成點文件。再根據(jù) Wizard-Surface Creation-From Data命令下的導向將這些離散點文件用DSI插值方法生成光滑曲面。其結果如圖1所示。

圖1 由等高線建立地表模型

2.2 地層形態(tài)的模擬

地層形態(tài)是對成層巖石幾何特點的描述，是地質構造的重要部分，用以下方法對其進行近似模擬。

2.2.1 露頭數(shù)據(jù)的分割及產(chǎn)狀的確定

模擬的第一步需要對已知數(shù)據(jù)進行進一步細化。如前所述，巢湖北部出露于地表的地質構造呈M型分布，其走向與傾角數(shù)據(jù)隨地點的不同而發(fā)生明顯的變化，在之后的建模過程中需要將產(chǎn)狀一致的露頭數(shù)據(jù)劃分為一個單元進行處理，為此，需要將已經(jīng)矢量化的數(shù)據(jù)根據(jù)區(qū)域內(nèi)產(chǎn)狀的不同而劃分為不同的面狀單元。在已有的1∶5萬地質圖中產(chǎn)狀數(shù)據(jù)已大致按不同區(qū)域加以標注，可以據(jù)此劃分數(shù)據(jù)。但是在某些區(qū)域，比如連接兩片區(qū)域的轉折處，地圖上沒有對此標注產(chǎn)狀數(shù)據(jù)，而這些地區(qū)在后期模型幾何形態(tài)的構建上又十分重要，因此，有必要利用地圖上已有信息去間接推斷此處的產(chǎn)狀。

根據(jù)原有地質圖上的等高線信息，可以選用三點法去求其產(chǎn)狀(見圖2)。其步驟如下：1)在所求產(chǎn)狀的層面上取3個不在同1條直線上相距不遠的點ABC，并連接最高點A與最低點C。2)根據(jù)點B的高程在線段AC上使用比例尺確定一點D，使其高程等于點B高程并連接DB，可以看出線段DB就是所求層面的走向線。3)過點C在圖上做1條與DB相平行的直線，則其為高程與點C相同的另一條走向線。4)在線段DB上任取一點E并設其在水平面上的投影點為O，過點E作DB的垂線交過點C的走向線為F，則EF在水平面上的投影OF為巖層的傾向線，使用量角器讀出其方位角后即可確定該層面傾向。5)由圖可見，角EFO的度數(shù)即是該層面的傾角，線段EO的高度即是DC兩點的高程差，OF的長度可根據(jù)比例尺在地圖上加以量算而得。將EO除以DC所得數(shù)值便是傾角的正切值，由此可得其度數(shù)。在獲取傾向與傾角后該區(qū)域的產(chǎn)狀便隨之確定。

圖2 “三點法”測量產(chǎn)狀

2.2.2 露頭數(shù)據(jù)投影至地表

劃分好的露頭數(shù)據(jù)在計算機中依然是以平面加以儲存的，要想使其在三維空間中表現(xiàn)真實的形態(tài)，必須對其賦予所有的高程數(shù)據(jù)。因為露頭數(shù)據(jù)直接出露于地表，其高程值就是地表高程值，這意味著只需要將露頭數(shù)據(jù)投影到地形表面即可獲取其實際形態(tài)。

為了對其投影，首先將露頭數(shù)據(jù)導入至GOCAD中，此時，它是以平面多邊形的形式展現(xiàn)的，在GOCAD中作為曲線(Curve)對象儲存，這種形式不能表現(xiàn)平面經(jīng)過投影而被拉伸的形態(tài)，要對其進一步加工，具體來說應使原始的多邊形轉換為不規(guī)則三角網(wǎng)(TIN)的形式來儲存面型數(shù)據(jù)。根據(jù)GOCAD軟件Application-Wizard-Surface Creation-From Date and Borders with bad Zs命令下的導向提示即可生成不規(guī)則三角網(wǎng)。此時即可對用三角網(wǎng)所表現(xiàn)的平面露頭數(shù)據(jù)進行投影操作，在菜單欄下選擇Property-By vertical projection，并在由此彈出的對話框內(nèi)的Vertical Property Client空框中輸入要被投影的對象，在Vertical property server from空框中輸入投影至的對象，分別在空框New Property name和Property內(nèi)輸入字母Z(這里代表高程)，表示將投影至的對象的高程值賦給被投影的對象。執(zhí)行完此命令后計算機將參照投影至的對象來調(diào)整被投影對象三角網(wǎng)內(nèi)各個頂點的相對位置使其獲得高程值。圖3顯示了將一小塊露頭數(shù)據(jù)投影至地表的效果，左圖為投影前的數(shù)據(jù)，右圖為投影之后的結果。

圖3 將露頭數(shù)據(jù)投影至地表

2.2.3 地下延伸體的計算構建

在完成對露頭數(shù)據(jù)的投影后，出露于地表的地層其幾何形態(tài)已經(jīng)確定，但是埋于地下的部分卻無法直接測定。因為地層是一個連續(xù)體，其地表部分的產(chǎn)狀在小范圍內(nèi)一般不會與地下部分的產(chǎn)狀有明顯的不一致，因此，可以采用間接的方式用地表巖層的產(chǎn)狀來代替它在地下延伸部分的產(chǎn)狀，從而獲取在一定深度范圍內(nèi)巖層的空間形態(tài)。

假定在一定深度范圍內(nèi)，地下巖層產(chǎn)狀與地表一致，同時考慮到地層在空間中一般是密實的連續(xù)體，那么可以據(jù)此采用構建柱狀三維體模型的方式來模擬此范圍內(nèi)的地層形態(tài)。這是因為假設巖層在該范圍內(nèi)產(chǎn)狀一致，從而它在空間中的傾斜方向也相同，而柱體當中的母線之間相互平行，都在空間中的同一個方向傾斜，因此，柱體的這一幾何特點可以很好地表達巖層的產(chǎn)狀。

在柱體形體模型中，它的上底面可直接用露頭數(shù)據(jù)代替，而其下底面可通過平移上底面獲得。為了使它表現(xiàn)真實的地層形態(tài)，需要根據(jù)露頭數(shù)據(jù)的產(chǎn)狀來確定下底面相對于上底面的平移量。為此，先建立坐標系(見圖4)，其中X軸正向為東，Y軸正向為北，Z軸正向，假設圖中的平面A為原始的露頭數(shù)據(jù)，將其平移至A1，設想在A中任意取一點C，它也隨A的平移過程落在A1上，可稱其為C1,即C在A1上的對應點，過點C做垂直于地面的垂線，同時過點C1向這條垂線作與地面相平行的射線，使兩線交于點B，則角B為直角。由此圖可知，若平移的方位與地層的產(chǎn)狀一致，那么由平移所形成的直角三角形CC1B中,角CC1B即為該地層產(chǎn)狀中的傾角，C1B即為產(chǎn)狀的傾向線，CB為深度。由上述的幾何關系可以給出點C1相對于點C在直角坐標系下坐標的改變量，設dx，dy分別代表其在X方向與Y方向的坐標增量，h代表深度,α與β分別代表傾角與傾向,則它們之間的關系可表示為

(1)

因為面是由點構成的，由此公式再計算出柱體中上底面的任意一個點在平移至下底面后坐標的改變量后，下底面的位置也就自然確定了。

圖4 延伸體計算

計算出體模型的上下兩底面后還要考慮斷層面對它的影響。斷層面的構建可由已知產(chǎn)狀的斷層線向下延伸獲得，利用式(1)計算地表的斷層線，依產(chǎn)狀向下平移一定深度后所處的位置，得到它后再將它與地表斷層線連接成閉合曲線便可生成斷層面的形態(tài)。地層如果與斷層面相交，應進行面相交計算切除多余地層，沒有相交的地方應對地層面進行延伸，直到與斷層面相交。

巖層在空間中是連續(xù)分布的密閉整體，通過對體模型上下底面位置的計算和與斷層面的關系處理，已經(jīng)確定了巖層所占據(jù)的空間范圍，接下來應對兩底面之間的空間進行填充，使其成為閉合的空間區(qū)域。

在GOCAD中生成體對象有兩種方法，分別是SGrid(網(wǎng)格模型對象)和Solid(實體模型對象)，其中，實體模型對象是通過點集、封閉曲面或單個面及厚度來創(chuàng)建，它們均需要非常詳細的數(shù)據(jù)來支持，且生成的模型因不能實現(xiàn)和其他程序交互計算所以在實際應用中有局限性。網(wǎng)格對象模型是用六面體格子組成的集合，在確定地層界面及網(wǎng)格的密度后軟件將使用六面體格子“填充”層面間的地層，且可以通過調(diào)整六面體單元的大小來控制擬合的細致程度?；谝陨系膶Ρ龋P者選用網(wǎng)格模型對象來構建三維體模型。

具體實現(xiàn)過程中應在Workflows選項下建立3D reservoir Grid Builder任務，并根據(jù)該任務的工作流程導向提示完成模型建立，主要步驟如下：1)選擇參考面Reference Horizon，并設置網(wǎng)格方向，這里參考面即是地層的頂面，后期的六面體單元將以此為基準面加以構建。2)Areal Gridding Specifications，這一步主要是在已經(jīng)選中的參考面上繪制平面網(wǎng)格(Gridding)決定其空間朝向。3)Create Straight Pillars，通過此步驟完成在豎直方向上對格網(wǎng)六面體指向的確定，在這個過程中需要指定另外一個參考面，即是之前計算出來的地層模型的下底面，軟件可以自動計算上下兩底面的連接線。4)Build Stratigraphic Grid，此步驟最終完成六面體網(wǎng)格模型的建立，它通過計算每個格網(wǎng)單元的中心和邊角位置來確定格網(wǎng)模型的形態(tài)，同時可以控制網(wǎng)格單元的密度來確定模型的精細程度。

圖5 部分地質模型

巢湖北部部分三維地質模型見圖5，圖中代表的是巢湖北部五通組地層的部分模型。

3 結束語

GOCAD是一款功能豐富的地質建模軟件，可通過離散光滑插值算法建立表面模型，使其更加貼近真實的地質體。對于創(chuàng)建相似的層面可使用投影方法直接將原有層面的屬性賦予被新創(chuàng)建的層面，增加了創(chuàng)建工作的靈活性。在創(chuàng)建空間體模型時，GOCAD提供了SGrid(網(wǎng)格模型對象)的方法，可依據(jù)層面的特點自動給出格網(wǎng)的形態(tài)，也可手動調(diào)整，保證模型的精確性。

本文通過收集巢湖北部區(qū)域的地形數(shù)據(jù)和地質數(shù)據(jù)，使用GOCAD軟件，創(chuàng)建此地區(qū)的地質模型，直觀、清晰地展現(xiàn)該地質結構的空間關系。

在建模過程中使用三點法確定小片區(qū)域內(nèi)的產(chǎn)狀，加強了模型的精度，同時建立了將產(chǎn)狀數(shù)據(jù)直接轉換為在直角坐標系中的偏移量公式，化簡建模過程，為三維地質建模進行有參考價值的探索。