沉積地層系統(tǒng)三維實(shí)體模型構(gòu)建方法

朱良峰 ，孫建中，張成娟

（1. 華東師范大學(xué) 地理信息科學(xué)教育部重點(diǎn)實(shí)驗(yàn)室，上海 200062；2. 中國科學(xué)院 上海高等研究院，上海 201203；3. 上海市城市化生態(tài)過程與生態(tài)恢復(fù)重點(diǎn)實(shí)驗(yàn)室，上海 200062）

1 引 言

沉積地層是一種常見的層狀地質(zhì)實(shí)體，構(gòu)建沉積地層環(huán)境下的三維地質(zhì)模型對(duì)實(shí)際的地質(zhì)分析與工程設(shè)計(jì)具有很大的應(yīng)用價(jià)值[1]。近10年來，國內(nèi)外在沉積地層系統(tǒng)三維模擬方面開展了卓有成效的研究，相關(guān)的研究進(jìn)展主要體現(xiàn)在建模數(shù)據(jù)標(biāo)準(zhǔn)化、三維實(shí)體模型表示方法、三維模型重構(gòu)技術(shù)等3個(gè)方面。目前，沉積地層系統(tǒng)建模數(shù)據(jù)標(biāo)準(zhǔn)化、三維實(shí)體模型表示方法這兩方面所涉及的問題已基本解決，并在學(xué)術(shù)界和工程界形成共識(shí)；而在三維地質(zhì)模型重構(gòu)技術(shù)方面，雖然提出了一些建模算法，但當(dāng)前仍缺乏簡單通用、高效自動(dòng)的建模方法，亟需結(jié)合各種復(fù)雜的地質(zhì)背景條件，對(duì)已有的構(gòu)模方法進(jìn)行優(yōu)化、改進(jìn)及驗(yàn)證工作。

工程鉆探法是一種傳統(tǒng)的地下空間信息獲取方式，鉆孔資料具有簡單、直觀、準(zhǔn)確、詳細(xì)的優(yōu)點(diǎn)?；阢@孔數(shù)據(jù)的建模方法一直是沉積地層系統(tǒng)三維重構(gòu)技術(shù)的主流，常用的有地層骨架構(gòu)模法[2－3]、地層層位法[4]、鉆孔-層面模型法[5]、地層縱向分塊處理法[6]等。經(jīng)過實(shí)際工程驗(yàn)證發(fā)現(xiàn)：現(xiàn)有的建模方法及實(shí)現(xiàn)算法都存在或多或少的缺陷，其中一個(gè)突出的問題就是這些方法缺乏合理完善的缺失地層處理機(jī)制，導(dǎo)致計(jì)算機(jī)建模結(jié)果與實(shí)際地層分布情況差別較大，無法真實(shí)地反映缺失地層的空間展布特征。缺失地層引起的地質(zhì)界面不連續(xù)會(huì)增加地層層面插值擬合的復(fù)雜度。雖然有研究者提出可以通過添加虛擬鉆孔[7]或地質(zhì)剖面[4,7]的方法將缺失地層界面控制信息追加到建模流程之中，但在具體實(shí)現(xiàn)的過程中，手工添加虛擬鉆孔或地質(zhì)剖面具有很大的主觀性，同時(shí)也是一項(xiàng)極為繁瑣的工作，不適合于用計(jì)算機(jī)進(jìn)行自動(dòng)或半自動(dòng)的建模。

為解決這一難題，在綜合前人成果的基礎(chǔ)上，本文將對(duì)基于鉆孔數(shù)據(jù)的沉積地層系統(tǒng)三維實(shí)體模型生成技術(shù)展開更為深入的研究。本文將首先總結(jié)歸納沉積地層系統(tǒng)的幾何特征；然后基于多層DEM構(gòu)模原理，提出一種新穎的三維地層模型生成技術(shù)“鉆孔-層面-實(shí)體法”，給出具體的實(shí)現(xiàn)過程與關(guān)鍵算法；最后通過上海城建工程中的一個(gè)研究實(shí)例來驗(yàn)證該技術(shù)的有效性。

2 沉積地層系統(tǒng)的幾何特征

沉積地層系統(tǒng)中的各個(gè)地層單元未經(jīng)變形或輕微變形，沒有被斷層或節(jié)理切割，形態(tài)相對(duì)簡單規(guī)則。按照地層的完整程度與空間分布特征，沉積地層系統(tǒng)中的地層可分為兩類[6]：一類是完整地層，它在研究區(qū)內(nèi)連續(xù)分布，其頂、底界面是完整的曲面；另一類是缺失地層，它在研究區(qū)內(nèi)呈不連續(xù)分布，其頂、底界面分別由帶有“空洞”的連續(xù)曲面組成，或是由多個(gè)不連續(xù)的曲面拼合而成。

根據(jù)不同的地質(zhì)成因，本文將沉積地層系統(tǒng)中的缺失地層分為3類。第1類是由非沉積引起的缺失地層，即地層的缺失區(qū)域根本就沒有發(fā)生過該地層的沉積，它一般是由地球內(nèi)動(dòng)力地質(zhì)作用引起的。第2類是由侵蝕引起的缺失地層，即該地層沉積后，在上覆地層沉積前就遭受到了剝蝕，它一般是由地表外動(dòng)力地質(zhì)作用（如河流侵蝕、風(fēng)化剝蝕）引起的。非沉積引起的缺失地層的缺失邊界受控于下伏地層的頂面。侵蝕引起的缺失地層的缺失邊界受控于上覆地層的底面。在實(shí)際的地層沉積過程中，還可能形成復(fù)合類別的缺失地層（本文稱其為第3類缺失地層）。復(fù)合缺失地層的缺失邊界既受控于上覆地層的底面，也受控于下伏地層的頂面。

基于大量的實(shí)例分析，可總結(jié)歸納出沉積地層系統(tǒng)所具備的5個(gè)幾何特征：

①成層性。沉積地層系統(tǒng)可以按照一定的標(biāo)準(zhǔn)劃分成若干個(gè)“層”，屬于同一“層”的地質(zhì)單元具有大體一致的沉積時(shí)代和物理力學(xué)指標(biāo)，在某種程度上可以看作是同一類型的物質(zhì)。

②順序性。在給定的研究場區(qū)，可按照一定的標(biāo)準(zhǔn)和規(guī)則整理出一個(gè)標(biāo)準(zhǔn)的地層沉積順序。

③連續(xù)性。沉積地層系統(tǒng)中的完整地層在研究區(qū)域內(nèi)是連續(xù)分布的，缺失地層在局部區(qū)域會(huì)出現(xiàn)不連續(xù)，但可將缺失區(qū)域視為頂、底界面重合（零厚度）的連續(xù)地層進(jìn)行特殊處理；因此，沉積地層系統(tǒng)中的所有地層都可視為連續(xù)層狀地質(zhì)體。

④封閉性。各個(gè)地層單元均由頂面、底面和側(cè)面封閉而成，頂面和底面在水平面上的投影完全重合[8]；在三維建模時(shí)，只需要重構(gòu)頂面和底面，側(cè)面可自動(dòng)生成。

⑤單值性。任一平面點(diǎn)P(x, y)在某一地層頂面或底面上對(duì)應(yīng)的高程值惟一。

沉積地層系統(tǒng)三維建模的結(jié)果被稱作三維地層模型，也常被稱作“蛋糕”模型（“l(fā)ayer-cake”model）[1]。但實(shí)際的地層單元空間展布特征遠(yuǎn)比蛋糕中的分層復(fù)雜。由于缺失地層的存在，地層面與地層面之間會(huì)發(fā)生交切關(guān)系，進(jìn)而導(dǎo)致地層在空間分布上的不連續(xù)，各個(gè)地層單元的幾何形態(tài)趨于復(fù)雜?；阢@孔數(shù)據(jù)，快速構(gòu)建地質(zhì)合理、高精度的三維地層模型并不是一件容易的工作。在實(shí)際的三維地層構(gòu)模算法設(shè)計(jì)與軟件開發(fā)中，必須充分考慮并體現(xiàn)由缺失地層引起的地層尖滅、侵入、透鏡體等特殊地質(zhì)現(xiàn)象對(duì)三維地層模型的影響[8－10]。

3 構(gòu)模原理

沉積地層系統(tǒng)的一個(gè)重要幾何特征在于其各個(gè)地層及其控制界面按照一定的順序有規(guī)律的堆疊在一起[1]。地層頂、底界面的重構(gòu)是沉積地層系統(tǒng)三維模擬過程中最基礎(chǔ)、最核心的工作。多層 DEM構(gòu)模技術(shù)是近年來被廣泛使用的構(gòu)建沉積地層系統(tǒng)三維實(shí)體模型的方法，它具有輸入數(shù)據(jù)簡單、建模流程清晰、模型生成速度快、自動(dòng)化程度較高等顯著優(yōu)點(diǎn)，目前已有多種實(shí)現(xiàn)算法，如地層層位法[4]、鉆孔-層面模型法[5]和地層骨架構(gòu)模法[2－3]。經(jīng)過實(shí)踐驗(yàn)證發(fā)現(xiàn)，地層層位法[4]和鉆孔-層面模型法[5]僅適用于沖積環(huán)境下非沉積作用引起的缺失地層（第1類缺失地層），無法重構(gòu)侵蝕作用引起的地層缺失。基于多層DEM構(gòu)模技術(shù)改進(jìn)的地層骨架構(gòu)模法[3]考慮了河流侵蝕作用對(duì)地層的切割影響，在建模的過程中對(duì)地層層面進(jìn)行了兩次插值與高程調(diào)整，能夠重構(gòu)侵蝕作用引起的缺失地層（第2類缺失地層），但難以推廣應(yīng)用于第1類和第3類缺失地層。

現(xiàn)有的多層DEM實(shí)現(xiàn)算法對(duì)沉積順序清晰、呈水平或近于水平的完整地層非常有效，但對(duì)地質(zhì)構(gòu)造活動(dòng)強(qiáng)烈、包含多種缺失地層的復(fù)雜場區(qū)則建模效果較差，模型無法達(dá)到應(yīng)有的精度，建模結(jié)果與實(shí)際情況出入較大。筆者認(rèn)為，之所以會(huì)出現(xiàn)這種情況，主要是因?yàn)楝F(xiàn)有的構(gòu)模算法沒有全面考慮缺失地層形成的地質(zhì)背景條件，無法將控制缺失地層空間展布特征的地質(zhì)原理轉(zhuǎn)化為計(jì)算機(jī)能夠識(shí)別、編程的構(gòu)模規(guī)則。針對(duì)這種缺陷，需要結(jié)合缺失地層形成的具體地質(zhì)背景條件，對(duì)現(xiàn)有的構(gòu)模算法進(jìn)行優(yōu)化、改進(jìn)，以提高建模的精度與效率。

基于多層DEM構(gòu)模原理，本文提出一種新穎的三維地層模型生成技術(shù)“鉆孔-層面-實(shí)體法”。這種技術(shù)克服了先前技術(shù)方法的局限，既能構(gòu)建非沉積引起的缺失地層，又能重構(gòu)侵蝕引起的缺失地層，還能處理非沉積和侵蝕疊加形成的復(fù)合缺失地層，以及多個(gè)缺失地層組合在一起的情況。該技術(shù)的適應(yīng)性強(qiáng)，準(zhǔn)確性高，能夠重構(gòu)包含任意數(shù)目缺失地層的復(fù)雜沉積地層系統(tǒng)的精細(xì)三維實(shí)體模型。按照該技術(shù)構(gòu)建三維地層模型共有11個(gè)關(guān)鍵環(huán)節(jié)。下面詳細(xì)介紹建模過程中的各個(gè)具體步驟。

4 建模步驟

4.1 界定建模場區(qū)，提取相關(guān)鉆孔信息

對(duì)于一個(gè)特定的研究場區(qū)，可能會(huì)有數(shù)目眾多的鉆孔。建模前，需要基于一個(gè)規(guī)范的數(shù)據(jù)庫格式，將全部鉆孔信息存儲(chǔ)到關(guān)系數(shù)據(jù)庫中。建模過程的第1步，就是從鉆孔數(shù)據(jù)庫中提取待建模場區(qū)的鉆孔信息，包括鉆孔類型、平面位置以及具體的地層分層信息。這些鉆孔數(shù)據(jù)將作為后續(xù)的地層屬性判別和地層頂、底界面插值擬合的初始數(shù)據(jù)。

4.2 建模場區(qū)總體地層編號(hào)與完整性判定

首先，將鉆孔所揭示的地層按照沉積順序進(jìn)行編號(hào)，生成一個(gè)涵蓋建模場區(qū)全部地層的“地層層序表”。進(jìn)行地層編號(hào)的規(guī)則是“從老到新，逐層遞增”，即：最老的地層編號(hào)為1，記為S1；其他地層按照沉積順序（或形成時(shí)代）逐層加 1，遞增編號(hào)，記為 S2、S3… Sn。

然后，根據(jù)鉆孔對(duì)地層的揭示程度，將鉆孔分為完全鉆孔和不完全鉆孔兩類。將建模場區(qū)最新地層（即位于最頂部的地層）和最老地層（即位于最底部的地層）都探測到的鉆孔為完全鉆孔；否則為不完全鉆孔。

接著，根據(jù)鉆孔所揭示的地層信息，將建模場區(qū)的地層分為完整地層和缺失地層兩類。在建模場區(qū)內(nèi)，全部完全鉆孔和可能探及此層的不完全鉆孔都揭示到的地層為完整地層；否則為缺失地層。

需要指出的是，為了處理方便，本文假定整個(gè)建模場區(qū)的最頂層和最底層地層一定是完整地層。如果建模場區(qū)地層不滿足這一條件，將在最頂層上方或最底層下方自動(dòng)添加虛擬的完整地層，然后在完成“生成三維實(shí)體模型”（見第4.10節(jié)）操作后，自動(dòng)剝除添加的虛擬地層。

4.3 鉆孔數(shù)據(jù)離散化

鉆孔分層信息是對(duì)上下相鄰的地層的分界面的描述，能夠揭示地層在鉆孔位置處的豎向分布狀況。建模過程的第3步，主要工作是將鉆孔所揭示的地層分界面離散為三維空間中的一系列散點(diǎn)，并使用二分拓?fù)浣Y(jié)構(gòu)[6,9]來記錄鉆孔所蘊(yùn)含的地層分層信息。在二分拓?fù)浣Y(jié)構(gòu)中，每一個(gè)散點(diǎn)均是兩個(gè)地層的分界面上的控制點(diǎn)，使用該分界面的上層和下層地層編號(hào)來描述該散點(diǎn)的二分拓?fù)鋵傩浴Ｃ恳粋€(gè)散點(diǎn)的數(shù)據(jù)結(jié)構(gòu)如下：

經(jīng)過鉆孔數(shù)據(jù)離散化操作，可得到一個(gè)散點(diǎn)數(shù)據(jù)集P，集合了地層頂、底界面高程插值擬合的樣本數(shù)據(jù)。另外，需要指出，不僅鉆孔數(shù)據(jù)可離散為散點(diǎn)，其他可用的地質(zhì)數(shù)據(jù)（如地質(zhì)剖面圖、地層等厚度線）都可離散化為三維空間中的散點(diǎn)，并使用二分拓?fù)浣Y(jié)構(gòu)來記錄地層分層信息，然后添加到散點(diǎn)數(shù)據(jù)集P中作為地層控制界面插值擬合的樣本點(diǎn)。

4.4 定義骨架三角網(wǎng)

沉積地層系統(tǒng)中各個(gè)地層的頂、底界面在水平面上的投影完全重合，因此，可定義一個(gè)骨架三角網(wǎng)作為擬合各個(gè)地層頂、底界面空間幾何形態(tài)的基準(zhǔn)三角網(wǎng)。所謂骨架三角網(wǎng)，是指以鉆孔孔口平面坐標(biāo)x、y為基準(zhǔn)，在建模場區(qū)外邊界的約束下，采用標(biāo)準(zhǔn)的三角化算法構(gòu)網(wǎng)并細(xì)分加密的一個(gè)不規(guī)則三角網(wǎng)（TIN）[4－5]。骨架三角網(wǎng)不僅明確界定了待構(gòu)建的三維地層模型的外邊界，還隱式定義了各個(gè)地層頂、底界面上下嚴(yán)格一致的幾何關(guān)系和拓?fù)潢P(guān)系。使用骨架三角網(wǎng)能減少后續(xù)的地層層面插值擬合、層面交切、高程調(diào)整、體元生成等處理步驟的復(fù)雜度，增強(qiáng)后續(xù)步驟實(shí)現(xiàn)算法的穩(wěn)健性。

4.5 插值擬合地層頂、底界面

沉積地層系統(tǒng)中的每一個(gè)地層都由頂面三角網(wǎng)和底面三角網(wǎng)封閉而成（側(cè)面可自動(dòng)生成），這兩個(gè)三角網(wǎng)都是基于骨架三角網(wǎng)生成的，因此，所有地層的頂、底界面都具有上下嚴(yán)格一致的幾何關(guān)系和拓?fù)潢P(guān)系。在建模過程的第5步，從散點(diǎn)數(shù)據(jù)集P中提取各個(gè)地層頂、底界面的樣本點(diǎn)坐標(biāo)，利用這些散點(diǎn)插值擬合相應(yīng)界面骨架三角網(wǎng)頂點(diǎn)的初始高程。

在擬合編號(hào)為 i的地層 Si的頂面三角網(wǎng) Fi_Top時(shí)，從散點(diǎn)數(shù)據(jù)集P中提取m_ iBelowID為i的散點(diǎn)作為插值擬合的樣本點(diǎn)；在擬合Si的底面三角網(wǎng)Fi_Bottom時(shí)，使用m_iAboveID為i的散點(diǎn)作為插值擬合的樣本點(diǎn)。在進(jìn)行插值時(shí)，可使用距離反比加權(quán)法、最鄰近點(diǎn)法、自然鄰近點(diǎn)法、徑向基函數(shù)法、Kriging法等常用的插值算法。這些算法概念簡單，方便快捷，既支持內(nèi)插，又支持外推，實(shí)際應(yīng)用效果良好。

4.6 判定缺失地層類別

這一步的工作是根據(jù)插值擬合出的缺失地層頂、底界面初始高程來自動(dòng)推算缺失地層的類別。對(duì)于任一缺失地層Si，其頂面三角網(wǎng)記為Fi_Top，底面三角網(wǎng)記為Fi_Bottom。如果Si在鉆孔Bj處缺失，則鉆孔Bj被稱作地層Si的缺失鉆孔，缺失地層的上覆地層 Si+1與下伏地層 Si-1的分界點(diǎn)（簡稱為缺失點(diǎn)）記為 P0，P0的高程記為 Z0。Si的頂面三角網(wǎng)Fi_Top在鉆孔 Bj處的頂點(diǎn)記為 PTop，底面三角網(wǎng)Fi_Bottom在鉆孔Bj處的頂點(diǎn)記為PBottom。經(jīng)過前一步插值擬合，可求得 PTop、PBottom的初始高程，分別記為 ZTop、ZBottom。使用式（1）計(jì)算 PTop、PBottom與P0之間的距離d1、d2：

由于 ZTop≥ZBottom，所以 d1≥d2。如果 d1=0、d2=0，即PTop、PBottom與缺失點(diǎn)P0重合，我們認(rèn)為地層Si恰巧在鉆孔Bj位置處尖滅，可作為0厚度的完整地層來對(duì)待；如果 d1≤0、d2＜0，則地層 Si在鉆孔 Bj處非沉積缺失（圖 1（a））。如果 d1＞0、d2≥0，則地層 Si在鉆孔 Bj處侵蝕缺失（圖 1（b））。如果d1＞0、d2＜0，則地層 Si在鉆孔Bj處復(fù)合缺失（圖1（c））。

除了地層Si恰巧在缺失鉆孔位置處尖滅的情況之外，如果地層Si在全部缺失鉆孔位置處均為非沉積缺失，認(rèn)為該地層屬于第1類缺失地層。如果地層Si在全部缺失鉆孔位置處均為侵蝕缺失，則可認(rèn)為該地層屬于第2類缺失地層。如果地層Si在某些缺失鉆孔位置處為復(fù)合缺失，或者在部分缺失鉆孔位置處呈非沉積缺失、在部分缺失鉆孔位置處侵蝕缺失，則認(rèn)為該地層屬于第3類缺失地層。

圖1 根據(jù)地層頂、底界面初始高程判定缺失地層類別Fig.1 Determining the types of missing strata in terms of the initial elevations of the top and bottom surfaces

4.7 生成缺失地層處理次序

經(jīng)過地層頂、底界面插值擬合生成了各地層頂、底界面三角網(wǎng)頂點(diǎn)的初始高程。對(duì)于完整地層，其頂面和底面三角網(wǎng)不會(huì)交叉，也不會(huì)與其他完整地層的頂、底界面三角網(wǎng)相交叉。但對(duì)于缺失地層，其頂面或底面三角網(wǎng)則可能會(huì)與其上覆地層或下伏地層的頂、底界面三角網(wǎng)交叉。需要根據(jù)缺失地層的類別，分別對(duì)缺失地層及其控制地層的頂、底界面進(jìn)行交切處理與高程調(diào)整。

如果一個(gè)地層系統(tǒng)中存在多個(gè)缺失地層，在進(jìn)行地層界面交切處理與高程調(diào)整前，需要按照地層的切錯(cuò)關(guān)系排出相應(yīng)的優(yōu)先級(jí)，生成缺失地層的處理次序?！般@孔-層面-實(shí)體法”按照下面的規(guī)則自動(dòng)進(jìn)行缺失地層排序工作：

①第3類缺失地層可視為第1類與第2類缺失地層疊加合成的結(jié)果。因此，在一個(gè)地層系統(tǒng)中，如果有第3類缺失地層，需要將第3類缺失地層分解為“第1類缺失地層+第2類缺失地層”，即將第3類缺失地層重復(fù)記錄，它既屬于第1類缺失地層，又屬于第2類缺失地層。

②第1類缺失地層的優(yōu)先級(jí)高于第2類缺失地層。在進(jìn)行缺失地層交切處理與高程調(diào)整時(shí)，先處理第1類缺失地層，再處理第2類缺失地層。

③如果有多個(gè)第 1類缺失地層，按照“從老到新，逐個(gè)處理”的規(guī)則，先處理下層的缺失地層，再處理其上的缺失地層。

④如果有多個(gè)第 2類缺失地層，按照“從新到老，逐個(gè)處理”的規(guī)則，先處理上層的缺失地層，再處理其下的缺失地層。

經(jīng)過排序，生成一個(gè)只含第1類缺失地層和第2類缺失地層的缺失地層處理隊(duì)列。

4.8 缺失地層控制界面交切處理與高程調(diào)整

這一步的主要工作是按照前一步生成的缺失地層處理次序，逐個(gè)對(duì)缺失地層控制界面進(jìn)行交切處理與高程調(diào)整。對(duì)于第1類和第2類缺失地層，因其受控于不同的地質(zhì)界面，使用的處理方法也有差異，下面分別描述。

4.8.1 第1類缺失地層

如果地層Si屬于第1類缺失地層，則Si的缺失邊界受控于下伏地層Si-1的頂面。如圖2（a）所示，經(jīng)過4.5節(jié)的插值，擬合出Si的頂、底界面及其下伏地層 Si-1的頂面、上覆地層 Si+1的底面的初始形態(tài)。其中，Si的頂面三角網(wǎng)記為Fi_Top，底面三角網(wǎng)記為 Fi_Bottom，Si-1的頂面三角網(wǎng)記為 F(i-1)_Top，Si+1的底面三角網(wǎng)記為F(i+1)_Bottom。

首先，處理Si的頂面三角網(wǎng)Fi_Top。將Fi_Top與F(i-1)_Top進(jìn)行相交運(yùn)算，求出這兩個(gè)三角網(wǎng)的所有交點(diǎn)，作為新頂點(diǎn)插入到骨架三角網(wǎng)中，同時(shí)修改Fi_Top、F(i-1)_Top以及其他全部地層頂、底界面三角網(wǎng)的頂點(diǎn)和拓?fù)潢P(guān)系。正常情況下，F(xiàn)i_Top中頂點(diǎn)的高程應(yīng)不低于F(i-1)_Top上對(duì)應(yīng)的頂點(diǎn)。但在Si缺失的區(qū)域，會(huì)出現(xiàn)Fi_Top頂點(diǎn)的高程低于F(i-1)_Top上對(duì)應(yīng)頂點(diǎn)的高程的情況。這時(shí)，需要逐個(gè)抬高 Fi_Top頂點(diǎn)的高程值，將其設(shè)定為與 F(i-1)_Top上對(duì)應(yīng)頂點(diǎn)的高程值相等（如圖 2（a）所示）。修改后的 Fi_Top記為F'i_Top。

然后，處理Si的底面三角網(wǎng)Fi_Bottom。從概念上來講，F(xiàn)i_Bottom與F(i-1)_Top的相應(yīng)頂點(diǎn)高程理應(yīng)完全一致。但由于Si是缺失地層，在4.5節(jié)分別插值求取這兩個(gè)三角網(wǎng)頂點(diǎn)高程時(shí)，使用的樣本點(diǎn)并不一致，這導(dǎo)致插值擬合出的對(duì)應(yīng)頂點(diǎn)高程也不完全相同。這需要通過界面交切處理與高程調(diào)整來實(shí)現(xiàn)兩個(gè)界面的一致。將Fi_Bottom與F(i-1)_Top進(jìn)行相交運(yùn)算，求出這兩個(gè)三角網(wǎng)的所有交點(diǎn)，作為新頂點(diǎn)插入到骨架三角網(wǎng)中，同時(shí)修改 Fi_Bottom、F(i-1)_Top以及其他全部地層頂、底界面三角網(wǎng)的頂點(diǎn)和拓?fù)潢P(guān)系；接著，逐個(gè)比較Fi_Bottom與F(i-1)_Top的對(duì)應(yīng)頂點(diǎn)高程是否相同：如果Fi_Bottom頂點(diǎn)的高程低于（或高于）F(i-1)_Top上對(duì)應(yīng)頂點(diǎn)的高程，則抬高（或壓低）Fi_Bottom頂點(diǎn)的高程值（如圖 2（a）所示），將其設(shè)定為與F(i-1)_Top對(duì)應(yīng)頂點(diǎn)的高程值相等。修改后的 Fi_Bottom記為 F’i_Bottom。

交切處理與高程調(diào)整后的缺失地層Si的頂、底界面如圖2（b）所示。

圖2 第1類缺失地層控制界面交切處理、高程調(diào)整與一致性處理Fig.2 Intersecting surfaces, adjusting elevations and keep consistency for type 1 missing stratum

4.8.2 第2類缺失地層

如果地層Si屬于第2類缺失地層，則Si的缺失邊界受控于上覆地層Si+1的底面。如圖3（a）所示，經(jīng)過4.5節(jié)的插值，擬合出Si的頂、底界面及其下伏地層 Si-1的頂面、上覆地層 Si+1的底面的初始形態(tài)。其中，Si的頂面三角網(wǎng)記為Fi_Top，底面三角網(wǎng)記為 Fi_Bottom，Si-1的頂面三角網(wǎng)記為 F(i-1)_Top，Si+1的底面三角網(wǎng)記為F(i+1)_Bottom。

首先，處理Si的底面三角網(wǎng)Fi_Bottom。將Fi_Bottom與 F(i+1)_Bottom進(jìn)行相交運(yùn)算，求出這兩個(gè)三角網(wǎng)的所有交點(diǎn)，作為新頂點(diǎn)插入到骨架三角網(wǎng)中，同時(shí)修改Fi_Bottom、F(i+1)_Bottom以及其他全部地層頂、底界面三角網(wǎng)的頂點(diǎn)和拓?fù)潢P(guān)系。正常情況下，F(xiàn)i_Bottom中頂點(diǎn)的高程應(yīng)不高于F(i+1)_Bottom上對(duì)應(yīng)的頂點(diǎn)。但在Si缺失的區(qū)域，會(huì)出現(xiàn)Fi_Bottom頂點(diǎn)的高程高于 F(i+1)_Bottom上對(duì)應(yīng)頂點(diǎn)的高程的情況。這時(shí)，需要逐個(gè)壓低 Fi_Bottom頂點(diǎn)的高程值，將其設(shè)定為與F(i+1)_Bottom上對(duì)應(yīng)頂點(diǎn)的高程值相等（如圖3（a）所示）。修改后的Fi_Bottom記為F’i_Bottom。

然后，處理Si的頂面三角網(wǎng)Fi_Top。將Fi_Top與F(i+1)_Bottom進(jìn)行相交運(yùn)算，求出這兩個(gè)三角網(wǎng)的所有交點(diǎn)，作為新頂點(diǎn)插入到骨架三角網(wǎng)中，同時(shí)修改Fi_Top、F(i+1)_Bottom以及其他全部地層頂、底界面三角網(wǎng)的頂點(diǎn)和拓?fù)潢P(guān)系；接著，逐個(gè)比較 Fi_Top與F(i+1)_Bottom的頂點(diǎn)高程是否相同：如果Fi_Top頂點(diǎn)的高程低于（或高于）F(i+1)_Bottom上對(duì)應(yīng)頂點(diǎn)的高程，則抬高（或壓低）Fi_Top頂點(diǎn)的高程值（如圖3（a）所示），將其設(shè)定為與F(i+1)_Bottom對(duì)應(yīng)頂點(diǎn)的高程值相等。修改后的Fi_Top記為F'i_Top。

交切處理與高程調(diào)整后的地層Si的頂、底界面如圖3（b）所示。

4.9 缺失地層與鄰接地層界面一致性處理

這一步的主要工作是按照4.7節(jié)生成的缺失地層處理次序，逐個(gè)處理缺失地層與鄰接地層的公共界面，以保證上下鄰接地層公共界面的一致性，即：所有地層（無論是完整地層，還是缺失地層）的頂面三角網(wǎng)與其上覆地層的底面三角網(wǎng)完全一致，所有地層的底面三角網(wǎng)與其下伏地層的頂面三角網(wǎng)完全一致。

如果地層Si屬于第1類缺失地層，則逐個(gè)修改上覆地層 Si+1底面三角網(wǎng) F(i+1)_Bottom的頂點(diǎn)高程，使其與 F'i_Top完全相同。修改后的 F(i+1)_Bottom記為F'(i+1)_Bottom。如果地層Si屬于第2類缺失地層，則逐個(gè)修改下伏地層Si-1頂面三角網(wǎng)F(i-1)_Top的頂點(diǎn)高程，使其與 F'i_Bottom完全相同。修改后的 F(i-1)_Top記為 F'(i-1)_Top。一致性處理后的地層頂、底界面分別如圖2（c）（第1類缺失地層）、圖3（c）（第2類缺失地層）所示。

圖3 第2類缺失地層控制界面交切處理、高程調(diào)整與一致性處理Fig.3 Intersecting surfaces, adjusting elevations and keep consistency for type 2 missing stratum

4.10 生成三維實(shí)體模型

經(jīng)過前面各步的處理，完成了各個(gè)地層頂、底界面三角網(wǎng)的重構(gòu)工作，生成了由三角面片組成的地層頂、底界面。這一步的工作是：從下到上逐個(gè)處理各個(gè)地層，以其頂面三角網(wǎng)和底面三角網(wǎng)為空間基準(zhǔn)，擠壓各個(gè)地層頂面三角網(wǎng)和底面三角網(wǎng)[4]，形成塊體模型；也可使用三維空間離散化技術(shù)對(duì)頂面三角網(wǎng)和底面三角網(wǎng)圍成的空間進(jìn)行剖分，生成用三棱柱體元充填的三維實(shí)體模型。

4.11 三維可視化分析

最后，對(duì)生成的地層實(shí)體模型進(jìn)行常規(guī)的可視化展示與專業(yè)的三維分析，包括三維觀察、交互定位、模型剖切、隧洞開挖、體積量算等?；谌庵w元表示的三維實(shí)體模型非常適合于空間查詢與空間分析，能夠直觀、準(zhǔn)確、全面地展示三維地層模型的內(nèi)部信息。

5 研究實(shí)例

研究場區(qū)位于上海市普陀區(qū)蘇州河畔，面積約800 m×600 m。如圖4（a）所示，研究區(qū)內(nèi)6個(gè)鉆孔揭示了4個(gè)地層（從下到上依次簡記為S1～S4），其中地層S1、S4是完整地層，S2、S3是缺失地層。圖4（b）～4（f）是使用這6個(gè)鉆孔數(shù)據(jù)自動(dòng)生成的各地層的控制界面以及三維實(shí)體模型。其中，圖4（b）展示了地層分界面的空間分布，圖4（c）、4（d）是從不同視角查看的三維實(shí)體模型，圖4（e）是過鉆孔C2、C6、C5、C4生成的剖面，圖4（f）是對(duì)三維實(shí)體模型進(jìn)行剖切操作后生成的柵狀圖。

從圖4中可以看出：地層S2、S3在某些鉆孔位置處（如C2）為非沉積缺失，在某些鉆孔位置處（如C4）為侵蝕缺失，地層的頂、底界面受非沉積作用和侵蝕作用共同控制，因此，屬第3類缺失地層；使用“鉆孔-層面-實(shí)體法”自動(dòng)重構(gòu)出來的地層頂、底界面形態(tài)自然，缺失地層的控制邊界也很合理。為了定量評(píng)估這個(gè)模型的精度，筆者使用實(shí)際工程建設(shè)過程中新增的鉆孔和基坑開挖數(shù)據(jù)進(jìn)行模型檢驗(yàn)。結(jié)果發(fā)現(xiàn)，在本文所構(gòu)建的三維地層模型中，完整地層頂、底界面高程與實(shí)測數(shù)據(jù)的誤差在5 cm以內(nèi)，缺失地層頂、底界面高程與實(shí)測數(shù)據(jù)的誤差小于8 cm?？傮w來看，模型的精度較高，可直接應(yīng)用于工程實(shí)際。這個(gè)實(shí)例直觀地展示了“鉆孔-層面-實(shí)體法”處理復(fù)合缺失地層的有效性和合理性，說明其具有自動(dòng)重構(gòu)復(fù)雜地質(zhì)條件下形成的沉積地層系統(tǒng)的能力。

圖4 復(fù)合缺失地層構(gòu)模實(shí)例Fig.4 Example of missing strata induced by combined action

6 結(jié) 語

沉積地層系統(tǒng)三維重構(gòu)的難點(diǎn)主要體現(xiàn)在對(duì)缺失地層地質(zhì)成因的判定以及缺失邊界的界定上。本文提出的“鉆孔-層面-實(shí)體法”克服了先前技術(shù)方法的局限，能充分挖掘鉆孔數(shù)據(jù)所揭示、隱含的地質(zhì)背景信息，自動(dòng)判斷缺失地層的地質(zhì)成因并界定缺失邊界，真實(shí)地反映缺失地層及其控制界面的空間展布特征與實(shí)際交切關(guān)系。與先前的技術(shù)方法相比，“鉆孔-層面-實(shí)體法”具有算法穩(wěn)健、自動(dòng)性高、建模結(jié)果合理、模型精度高等顯著優(yōu)勢。另外，“鉆孔-層面-實(shí)體法”在實(shí)現(xiàn)的過程中，巧妙地將控制地層空間幾何形態(tài)的地質(zhì)原理轉(zhuǎn)化為計(jì)算機(jī)能夠識(shí)別、編程的構(gòu)模規(guī)則，這對(duì)進(jìn)一步研究復(fù)雜地質(zhì)體三維構(gòu)模技術(shù)具有重要的借鑒意義。

[1]TURNER A K. Challenges and trends for geological modelling and visualization[J]. Bulletin of Engineering Geology and the Environment, 2006, 65(2): 109－127.

[2]賀懷建, 白世偉, 趙新華, 等. 三維地層模型中地層劃分的探討[J]. 巖土力學(xué), 2002, 23(5): 637－639.HE Huai-jian, BAI Shi-wei, ZHAO Xin-hua, et al.Discussion on strata partition in three-dimension strata model[J]. Rock and Soil Mechanics, 2002, 23(5): 637－639.

[3]朱良峰, 潘信. 河流侵淤作用下三維地層模型的構(gòu)建[J].巖土力學(xué), 2005, 26(增刊1): 65－68.ZHU Liang-feng, PAN Xin. Reconstruction of 3D stratigraphic model for fluvial erosion and aggrading action[J]. Rock and Soil Mechanics, 2005, 26(Supp.1):65－68.

[4]LEMON A M, JONES N L. Building solid models from boreholes and user-defined cross-sections[J]. Computers& Geosciences, 2003, 29(5): 547－555.

[5]朱良峰, 吳信才, 劉修國, 等. 基于鉆孔數(shù)據(jù)的三維地層模型的構(gòu)建[J]. 地理與地理信息科學(xué), 2004, 20(3):26－30.ZHU Liang-feng, WU Xin-cai, LIU Xiu-guo, et al.Reconstruction of 3D strata model based on borehole data[J]. Geography and Geo-Information Science, 2004,20(3): 26－30.

[6]張芳, 朱合華, 寧民霞. 適于海量數(shù)據(jù)的三維地層建模方法[J]. 巖石力學(xué)與工程學(xué)報(bào), 2006, 25(增刊1): 3305－3310.ZHANG Fang, ZHU He-hua, NING Min-xia. Modeling method of 3D strata suitable for massive data[J]. Chinese Journal of Rock Mechanics and Engineering, 2006,25(Supp.1): 3305－3310.

[7]朱良峰, 吳信才, 潘信. 三維地層模型誤差修正機(jī)制及其實(shí)現(xiàn)技術(shù)[J]. 巖土力學(xué), 2006, 27(2): 268－271.ZHU Liang-feng, WU Xin-cai, PAN Xin. Mechanism and implementation of error correction for 3D strata model[J].Rock and Soil Mechanics, 2006, 27(2): 268－271.

[8]徐能雄, 何滿潮. 層狀巖體三維構(gòu)模方法與空間數(shù)據(jù)模型[J]. 中國礦業(yè)大學(xué)學(xué)報(bào), 2004, 33(1): 103－108.XU Neng-xiong, HE Man-chao. 3D modeling methods and spatial data model of layered rock-mass[J]. Journal of China University of Mining & Technology, 2004,33(1): 103－108.

[9]朱合華, 吳江斌. 基于Delaunay構(gòu)網(wǎng)的地層2D，2.5D建模[J]. 巖石力學(xué)與工程學(xué)報(bào), 2005, 24(22): 4073－4079.ZHU He-hua, WU Jiang-bin. 2D and 2.5D modeling of strata based on Delaunay triangulation[J]. Chinese Journal of Rock Mechanics and Engineering, 2005,24(22): 4073－4079.

[10]趙宏堅(jiān), 周翠英. 基于實(shí)體建模的三維地層構(gòu)造[J]. 巖土力學(xué), 2010, 31(4): 1257－1263.ZHAO Hong-jian, ZHOU Cui-ying. Solid-oriented 3D-stratum construction[J]. Rock and Soil Mechanics,2010, 31(4): 1257－1263.