魏 嘉，唐 杰，武港山，張 楊

（1.中國石油化工股份有限公司石油物探技術(shù)研究院，江蘇南京211103；2.南京大學(xué)，江蘇南京210083）

復(fù)雜油氣藏三維地質(zhì)建模技術(shù)已成功應(yīng)用于油氣藏開發(fā)等方面，然而，隨著油氣藏勘探開發(fā)的不斷深入，對油氣藏地質(zhì)建模工作提出了更高的要求［1］。油氣藏地質(zhì)建模的工作流程分為構(gòu)造建模和屬性建模兩個階段，構(gòu)造建模是復(fù)雜油氣藏地質(zhì)建模的基礎(chǔ)，其技術(shù)的形成比屬性建模技術(shù)晚，但發(fā)展速度明顯快于屬性建模技術(shù)，尤其是隨著計算幾何和三維可視化等計算機(jī)技術(shù)的發(fā)展，促進(jìn)了三維構(gòu)造建模相關(guān)技術(shù)方法的進(jìn)步和發(fā)展［2］。早在1984年，Haldorson提出了油田尺度下用于油藏動態(tài)模擬的隨機(jī)模擬儲層建模方法［3］；Mallet于1989年和1992年先后發(fā)表了兩篇關(guān)于離散光滑插值［4－5］建模方法的文章，標(biāo)志著三維構(gòu)造建模技術(shù)中的地質(zhì)曲面技術(shù)獲得了突破；為了不斷提高三維地質(zhì)模型的描述精度，層位面和斷層面三角網(wǎng)格重建、地質(zhì)曲面空間拓?fù)浞治龊偷刭|(zhì)塊體劃分等技術(shù)得到了長足的發(fā)展［6－11］。復(fù)雜油氣藏三維地質(zhì)建模精度的提高，必然會產(chǎn)生巨大的計算量和數(shù)據(jù)量，這對三維地質(zhì)模型的數(shù)據(jù)管理和三維可視化顯示提出了挑戰(zhàn)；有效地組織和管理精細(xì)地質(zhì)模型的海量數(shù)據(jù)、實時地進(jìn)行三維地質(zhì)模型海量數(shù)據(jù)的三維可視化顯示成為急需解決的技術(shù)難題。

三維可視化技術(shù)改變了利用三維地震資料進(jìn)行油氣藏表征和建模的工作流程，并成為降低油氣勘探風(fēng)險和不確定性的有效技術(shù)手段［12］。由于受到客觀技術(shù)條件的制約（計算機(jī)內(nèi)存和圖形卡顯存總是不能滿足三維地震數(shù)據(jù)體數(shù)據(jù)量的要求），為了能夠?qū)崿F(xiàn)海量地震規(guī)則數(shù)據(jù)體的實時顯示，許多學(xué)者對在三維地震數(shù)據(jù)體繪制中如何減少I／O操作的問題進(jìn)行了研究［13－16］，主要方法包括減少體繪制過程中一次讀入的數(shù)據(jù)量和使用壓縮方法減少需要讀入的總數(shù)據(jù)量；另外，基于數(shù)據(jù)分塊和層次細(xì)節(jié)（LOD）模型技術(shù)進(jìn)行三維地震等規(guī)則數(shù)據(jù)體的多分辨率數(shù)據(jù)組織和管理技術(shù)也得到了廣泛的應(yīng)用［17－20］，通過基于擴(kuò)展四叉樹或八叉樹等數(shù)據(jù)結(jié)構(gòu)，建立分塊多層的多分辨率模型進(jìn)行LOD數(shù)據(jù)建模［21］，可以顯著提高海量地震數(shù)據(jù)體存取和顯示的效率。

三維地質(zhì)模型數(shù)據(jù)主要包括以空間多邊形（三角網(wǎng)格為主）表達(dá)的地質(zhì)層位面和斷層面（三維構(gòu)造框架模型）以及以空間多面體（任意六面體和四面體等）表達(dá)的地層網(wǎng)格（三維地層實體模型的基本單元）；三維地震數(shù)據(jù)體主要采用規(guī)則化立方體進(jìn)行數(shù)據(jù)的描述。對于非規(guī)則化的地質(zhì)曲面數(shù)據(jù)和地層網(wǎng)格數(shù)據(jù)的多分辨率數(shù)據(jù)組織和管理，有兩點基本要求：①不同分辨率的數(shù)據(jù)不能改變層位面、斷層面和地質(zhì)體的空間形態(tài)；②能夠為三維地質(zhì)模型三維可視化實時顯示提供快速的數(shù)據(jù)存取。

1 基于外存的地質(zhì)曲面海量數(shù)據(jù)多分辨率組織和繪制

目前國內(nèi)外已提出的解決大規(guī)模三角網(wǎng)格實時顯示問題的方法可以分為層次細(xì)節(jié)模型和多分辨率模型兩類。層次細(xì)節(jié)模型方法是指預(yù)先在原模型數(shù)據(jù)的基礎(chǔ)上，生成若干不同分辨率的近似模型，在實時繪制中根據(jù)當(dāng)前的視點參數(shù)選用相應(yīng)的近似模型進(jìn)行繪制。這類算法的優(yōu)點是簡單易用、速度快；缺點是近似模型之間不連續(xù)，顯示時會有跳躍感，而層次的增多必然導(dǎo)致整個模型占用的存儲空間增大。多分辨率模型方法是將網(wǎng)格組織成一個基網(wǎng)格加上一系列細(xì)化操作?？梢暬@示時，隨著視點的移動，多分辨率網(wǎng)格將實時自適應(yīng)更新。相對于層次細(xì)節(jié)模型來說，這種方式建立起來的場景更加符合人的視覺特性，模型的細(xì)節(jié)是連續(xù)變化的。但是，多分辨率網(wǎng)格的建立和自適應(yīng)更新要耗費額外的時間，其技術(shù)難點是要在精度和速度之間選擇一個合適的平衡點。

Lindstrom等［22］在1996年首先提出了針對數(shù)字高程（DEM）模型的實時動態(tài)多分辨率地形繪制算法，并不斷改進(jìn)了相關(guān)算法［23－24］。該算法基于一種四叉樹的數(shù)據(jù)結(jié)構(gòu)，將大范圍的DEM數(shù)據(jù)用四叉樹思想遞歸地分割成小的地形塊，再根據(jù)視覺范圍來控制地形塊的分裂和合并。該算法處理的三角網(wǎng)格是規(guī)則的DEM模型，需要在矩形參數(shù)域上等間距規(guī)則采樣。另外，Hoppe等［25－26］提出了針對非規(guī)則三角形網(wǎng)格（TIN）模型的漸進(jìn)網(wǎng)格算法，該算法首先將網(wǎng)格簡化，并記錄其簡化過程，進(jìn)而將簡化的數(shù)據(jù)組織成一個層次結(jié)構(gòu)，并定義了分裂或合并操作的有效性規(guī)則，從而將三角網(wǎng)格組織成一個基網(wǎng)格和一系列頂點關(guān)聯(lián)操作，方便了顯示不同的分辨率。

1.1 地質(zhì)曲面海量三角網(wǎng)格的分塊

由于對三維地質(zhì)模型中地質(zhì)曲面建模精度要求的不斷提高，在三維層位面和斷層面建模過程中，需要處理上百萬、甚至于上千萬的三角網(wǎng)格，而計算機(jī)的內(nèi)存總是有限的。我們對Hoppe等［25－26］提出的算法進(jìn)行了改進(jìn)，在創(chuàng)建地質(zhì)曲面海量三角網(wǎng)格多分辨率數(shù)據(jù)時，提出并實現(xiàn)了基于外存的海量三角網(wǎng)格多分辨率數(shù)據(jù)創(chuàng)建算法，該算法采用二叉樹的數(shù)據(jù)結(jié)構(gòu)對三角網(wǎng)格模型進(jìn)行劃分，建立層次結(jié)構(gòu)，依次對每塊數(shù)據(jù)進(jìn)行多分辨率模型創(chuàng)建，并將中間結(jié)果輸出到臨時性的數(shù)據(jù)文件，最后再反向創(chuàng)建多分辨率網(wǎng)格的數(shù)據(jù)文件。

數(shù)據(jù)分塊的方法有八叉樹劃分、四叉樹劃分和二叉樹劃分等。八叉樹劃分和四叉樹劃分等在用于劃分規(guī)則數(shù)據(jù)體（例如圖像和規(guī)則體數(shù)據(jù)）時非常有效，因為這些數(shù)據(jù)的拓?fù)湫畔ⅲ礃颖局g的相互位置關(guān)系是隱含在樣本的索引值上的。這樣相鄰兩塊數(shù)據(jù)之間分辨率容易統(tǒng)一。而三角網(wǎng)格的拓?fù)湫畔⑹秋@式表達(dá)的，即三角網(wǎng)格不僅有其幾何信息（各頂點的空間位置），還要顯式地記錄拓?fù)湫畔ⅲǜ黜旤c之間的連接性信息）。各數(shù)據(jù)塊分別處理，相鄰兩塊數(shù)據(jù)之間的拓?fù)湫畔⒑茈y做到一致，輕則造成邊界處的數(shù)據(jù)裂縫，重則導(dǎo)致多分辨率模型崩潰。因此，為減少相鄰數(shù)據(jù)塊之間的同步壓力，在分析了三角網(wǎng)格二維流形的特征以后，我們選擇二叉樹作為分塊的數(shù)據(jù)結(jié)構(gòu)（圖1）。其基本算法流程為：①掃描三角網(wǎng)格文件，計算三角網(wǎng)格的頂點數(shù)、三角網(wǎng)格的個數(shù)以及三維構(gòu)造框架模型的包圍盒大小；②根據(jù)包圍盒的大小以及設(shè)定的每塊數(shù)據(jù)的最大頂點數(shù)，沿包圍盒最長方向計算二叉樹劃分順序和劃分層數(shù)；③在二叉樹劃分順序和劃分層數(shù)計算的基礎(chǔ)上，構(gòu)建二叉樹的層次父子結(jié)構(gòu)；④根據(jù)三角網(wǎng)格頂點的空間位置和三角網(wǎng)格邊簡化的相關(guān)性，將其歸屬到相應(yīng)層次的葉節(jié)點上；⑤建立三角網(wǎng)格頂點與二叉樹不同層次葉節(jié)點的對應(yīng)關(guān)系，并設(shè)置三角網(wǎng)格3個頂點的合并關(guān)系；⑥輸出三角網(wǎng)格數(shù)據(jù)分塊的劃分結(jié)果。

圖1 三角網(wǎng)格頂點的層次結(jié)構(gòu)

1.2 基于拓?fù)湟恢碌亩喾直媛誓Ｐ偷膭?chuàng)建

創(chuàng)建層位面或斷層面多分辨率網(wǎng)格模型時，首先由葉子節(jié)點開始，依次對每個節(jié)點的網(wǎng)格頂點數(shù)量進(jìn)行簡化。在簡化過程中，不僅需要保證層位面和斷層面的空間形態(tài)基本保持不變，同時也要保證相鄰節(jié)點中網(wǎng)格模型的拓?fù)湫畔⒁恢?。為此，我們提出了邊界三角網(wǎng)格整形處理方法，其中關(guān)鍵環(huán)節(jié)包括三角網(wǎng)格頂點簡化的誤差計算、三角網(wǎng)格形狀的評估因子計算和數(shù)據(jù)分塊的邊界處理。

1.2.1 三角網(wǎng)格頂點簡化的誤差計算

在簡化三角形頂點的過程中，首先在每個頂點計算一個誤差系數(shù)矩陣，當(dāng)該頂點被合并到另一個頂點，也就是簡化三角網(wǎng)格的一條邊時，該系數(shù)矩陣也被加到那個頂點的系數(shù)矩陣上。當(dāng)三角網(wǎng)格的一條邊簡化后，三角網(wǎng)格的形狀在該條邊的鄰接范圍內(nèi)發(fā)生了變化，其新頂點到原來的鄰接三角網(wǎng)格的距離的平方和即為簡化誤差，可以用公式（1）計算：

其中，v表示空間一點；EQ（v）表示空間一點v到某平面的距離誤差；D（v）表示空間一點到某個平面的距離，

A是一個3×3的系數(shù)矩陣；b為系數(shù)向量；c為系數(shù)標(biāo)量。

1.2.2 三角網(wǎng)格形狀的評估因子計算

在三角網(wǎng)格剖分計算過程中，出現(xiàn)內(nèi)角非常小的狹長三角形，對于后期的數(shù)據(jù)處理和圖形繪制都是不利的，現(xiàn)有的大部分三角網(wǎng)格計算方法要求三角網(wǎng)格應(yīng)盡量趨于正三角形。為了避免計算結(jié)果出現(xiàn)不可接受的誤差，我們引入Gueziec［27］提出的三角網(wǎng)格形狀的評估因子：

式中：li為邊長；ω為三角片的面積。

γ值為0～1。γ越小，說明三角網(wǎng)格的狹長度越大，越不利于三角網(wǎng)格的處理；γ越大，說明三角片越接近于等邊三角形。為防止狹長三角片的出現(xiàn)，再計算刪除邊后所有目標(biāo)點的一階鄰接三角片形狀評估因子，如果存在小于指定閾值的值，表示有三角網(wǎng)格狹長，則在其計算出的簡化誤差上加上一個值，推后其在網(wǎng)格簡化隊列中的位置。

1.2.3 三角網(wǎng)格數(shù)據(jù)分塊的邊界處理

對葉節(jié)點處理完畢后，根據(jù)二叉樹層次結(jié)構(gòu)，依次將相鄰葉節(jié)點剩余的頂點和三角網(wǎng)格以及相鄰處的三角網(wǎng)格合并起來（圖2，紅圈處表示分塊處理的塊內(nèi)三角網(wǎng)格，需要進(jìn)行合并處理），再進(jìn)行簡化處理。

圖2 在分塊邊界上合并相鄰的節(jié)點

通過計算誤差系數(shù)矩陣，保證了地質(zhì)曲面空間形態(tài)的相似性；通過計算三角網(wǎng)格形態(tài)評價因子，避免了狹長三角形的出現(xiàn)，可以提高數(shù)據(jù)處理的穩(wěn)定性；通過綜合應(yīng)用誤差系數(shù)矩陣計算、三角網(wǎng)格形態(tài)評價因子計算和邊界三角網(wǎng)格的特殊處理方法，實現(xiàn)了地質(zhì)曲面多分辨率數(shù)據(jù)組織的幾何一致和拓?fù)湟恢拢▓D3：圖2是分塊處理的地質(zhì)曲面形態(tài)，與實際地質(zhì)曲面的拓?fù)潢P(guān)系存在一定的差異；進(jìn)行三角網(wǎng)格的整形優(yōu)化處理后得到圖3，保證了地質(zhì)曲面的幾何一致和拓?fù)湟恢拢?/p>

1.3 大規(guī)模三角網(wǎng)格多分辨率模型的實時繪制

由于地質(zhì)曲面多分辨率模型中的點、三角網(wǎng)格以及網(wǎng)格頂點的分裂和合并操作都是以連續(xù)的方式存放的，因此可以利用內(nèi)存映射機(jī)制將存在磁盤上的數(shù)據(jù)映射進(jìn)內(nèi)存，借助操作系統(tǒng)的頁面換入／換出機(jī)制實現(xiàn)數(shù)據(jù)調(diào)度。當(dāng)多分辨率數(shù)據(jù)量較大時，不能一次性映射進(jìn)內(nèi)存，需要分批將磁盤文件映射進(jìn)內(nèi)存。

初始時，前端頂點保存基網(wǎng)格頂點。每一幀中，根據(jù)視點參數(shù)判斷活躍頂點是否需要合并或者分裂。當(dāng)活躍頂點在可視范圍內(nèi)且誤差大于指定閾值時拆分該點，當(dāng)前端頂點不在可視范圍或父節(jié)點誤差小于閾值時合并該點。最后采用映射緩沖區(qū)對象的方法來加速繪制速度。緩沖區(qū)對象是一種數(shù)據(jù)存儲倉庫，用于在本地視頻內(nèi)存中存儲數(shù)據(jù)。在頂點數(shù)組繪制模式下，將活躍三角形直接生成到顯存可以顯著提升繪制速率，繪制結(jié)果如圖4所示。圖4中紅框表示實際的可視區(qū)域，該區(qū)域內(nèi)采用高分辨率的三角網(wǎng)格數(shù)據(jù)進(jìn)行顯示，該區(qū)域以外的部分采用較低分辨率的三角網(wǎng)格數(shù)據(jù)。通過應(yīng)用上述技術(shù)基本達(dá)到了海量三角網(wǎng)格多分辨率實時繪制的要求（圖5）。

2 海量三維地層網(wǎng)格數(shù)據(jù)的多分辨率數(shù)據(jù)組織和可視化

地質(zhì)曲面組合構(gòu)建的是三維地層框架模型，但是，僅僅有層位面和斷層面，對于三維地質(zhì)建模是不夠的，還需要對地層的內(nèi)部實體進(jìn)行建模。地層實體建模的方法是基于已經(jīng)建立的地層框架模型，對模型中的地層按照沉積約束條件進(jìn)行地層網(wǎng)格剖分，進(jìn)而建立以三維地層網(wǎng)格為單元的地層實體模型［2］。由于復(fù)雜油氣藏的三維地質(zhì)建模的精度要求越來越高，三維地層網(wǎng)格的剖分粒度越來越細(xì)，一個三維地質(zhì)模型的地層網(wǎng)格數(shù)據(jù)量往往達(dá)到了數(shù)百GB級別，目前的三維可視化顯示技術(shù)很難滿足實時性的要求。為此，我們提出了一種基于數(shù)據(jù)分塊的層次細(xì)節(jié)模型的三維地層網(wǎng)格數(shù)據(jù)組織和管理技術(shù)。

2.1 三維地層網(wǎng)格數(shù)據(jù)的數(shù)據(jù)分塊

空間數(shù)據(jù)分塊是一種基于外存數(shù)據(jù)調(diào)度策略的技術(shù)，該技術(shù)充分利用地質(zhì)模型的空間分布性質(zhì)將網(wǎng)格數(shù)據(jù)在空間上劃分為若干子塊，每個子塊可獨立進(jìn)行處理。這樣對快速地進(jìn)行視錐判斷和圖形裁減非常有利，而且可以采用多線程的方法對多個數(shù)據(jù)子塊并行化地進(jìn)行處理，從而加快數(shù)據(jù)處理的速度，減少響應(yīng)等待時間，改善三維地質(zhì)建模過程的交互性能。

在基于數(shù)據(jù)分塊的調(diào)度系統(tǒng)中，數(shù)據(jù)塊的大小是直接影響系統(tǒng)效率的重要因素，過小的分塊導(dǎo)致過多的I／O中斷，降低磁盤的讀取效率；過大的分塊不僅會導(dǎo)致讀入過多的無用數(shù)據(jù)，而且較長的單次讀取時間會造成整個系統(tǒng)的響應(yīng)緩慢。因此，必須慎重選擇數(shù)據(jù)分塊的大小。考慮到操作系統(tǒng)進(jìn)行頁面調(diào)度和磁盤緩存的因素，為減少缺頁中斷的次數(shù)，我們將數(shù)據(jù)分塊設(shè)置為32×32×32，這樣每個數(shù)據(jù)塊大小為8B×32×103＝256KB，并根據(jù)分塊尺度的信息建立相關(guān)模型數(shù)據(jù)的索引，從而提高數(shù)據(jù)的存取效率。

2.2 面向地層網(wǎng)格數(shù)據(jù)的LOD模型

傳統(tǒng)的三維可視化體繪制方法并不區(qū)分場景中對象的狀態(tài)，所有位于視錐范圍之內(nèi)的對象都會被以最高精度進(jìn)行繪制，一旦需要繪制的對象數(shù)據(jù)量急劇增大，這種繪制全部細(xì)節(jié)的方法便會力不從心。事實上，由于屏幕分辨率有限，無限制地提高繪制精度對于最終顯示效果的提高作用微乎其微。當(dāng)對象在場景中距離視平面較遠(yuǎn)時，可以采用較低分辨率的數(shù)據(jù)進(jìn)行繪制，而并不會明顯降低顯示的質(zhì)量。在處理大規(guī)模地層網(wǎng)格數(shù)據(jù)方面，我們采用了LOD模型進(jìn)行數(shù)據(jù)的組織和管理。

目前面向大規(guī)模數(shù)據(jù)實時繪制的LOD技術(shù)主要分為靜態(tài)LOD方法和動態(tài)LOD方法［28］。靜態(tài)LOD方法以離線方式預(yù)先對數(shù)據(jù)模型進(jìn)行處理，生成一系列不同細(xì)節(jié)層次的網(wǎng)格模型，繪制時根據(jù)顯示參數(shù)和硬件性能選擇一個合適的層次進(jìn)行繪制。靜態(tài)LOD方法將耗費資源的數(shù)據(jù)生成過程作為預(yù)處理，因此實時性較好，但在視點實時變換及細(xì)節(jié)層次突變的處理上不如動態(tài)LOD方法。動態(tài)LOD方法則是在繪制時才生成相應(yīng)的細(xì)節(jié)層次進(jìn)行渲染，因此需要對模型作實時簡化計算，對硬件性能要求較高［29］?？紤]到大規(guī)模數(shù)據(jù)的實時處理效率，我們?yōu)楹Ａ康貙泳W(wǎng)格數(shù)據(jù)設(shè)計了一種靜態(tài)LOD模型。

八叉樹結(jié)構(gòu)是由四叉樹結(jié)構(gòu)推廣到三維空間而形成的一種三維數(shù)據(jù)結(jié)構(gòu)［30－31］，它能夠?qū)θS空間形成有效劃分和高效的分層次管理，因此適合于組織空間數(shù)據(jù)的LOD模型。空間八叉樹通常采用自頂向下的方法構(gòu)建，節(jié)點之間采用指針的方式表達(dá)派生關(guān)系，優(yōu)點是查詢效率較高，缺點是大量的指針數(shù)據(jù)占據(jù)了本不充裕的存儲空間。我們針對地層網(wǎng)格存儲方式具有的規(guī)則數(shù)據(jù)體的特殊性，采用自下向上的順序構(gòu)建地層網(wǎng)格的八叉樹LOD模型。不同于傳統(tǒng)的八叉樹結(jié)構(gòu)，節(jié)點之間的指向關(guān)系不需要用指針進(jìn)行存儲，而是通過節(jié)點索引得到，節(jié)省了存儲大量指針數(shù)據(jù)的空間，更有利于海量數(shù)據(jù)模型的多分辨率處理。

在構(gòu)建層次細(xì)節(jié)模型時，首先依據(jù)數(shù)據(jù)分塊的大小對原始地層的規(guī)則數(shù)據(jù)體在X，Y和Z 3個維度上進(jìn)行空間劃分，形成最高級分辨率的數(shù)據(jù)塊陣列，并作為葉節(jié)點存入八叉樹的最低層。然后根據(jù)八叉樹的特點，在X，Y和Z 3個維度上分別對緊鄰的數(shù)據(jù)塊進(jìn)行二合一操作，這樣空間關(guān)系上緊鄰的8個數(shù)據(jù)塊將被合并為一個數(shù)據(jù)塊，形成次高級分辨率的數(shù)據(jù)塊陣列，并存入八叉樹的上一層。如果在某個維度上的數(shù)據(jù)塊個數(shù)為奇數(shù)，則增補(bǔ)一個空白數(shù)據(jù)塊進(jìn)行填充，以保證八叉樹結(jié)構(gòu)的完整性。這個過程順著八叉樹以自底向上的順序迭代進(jìn)行，直到八叉樹某個層次上的數(shù)據(jù)塊個數(shù)達(dá)到閾值為止。

2.3 地層網(wǎng)格LOD模型的數(shù)據(jù)流水線

地層網(wǎng)格多分辨率LOD模型的構(gòu)建過程中，控制點的生成需要大量的浮點計算，CPU資源的開銷很大；將數(shù)據(jù)塊寫入到輸出文件則需要與慢速的磁盤I／O系統(tǒng)協(xié)同工作，經(jīng)常造成CPU等待，資源閑置。這兩個階段在整個數(shù)據(jù)流程中耗時最長，同時系統(tǒng)資源卻沒有得到充分利用。為了充分發(fā)揮硬件性能，提高程序的運行效率，我們將LOD模型的構(gòu)建過程進(jìn)行4級流水線處理，使用4個線程分別負(fù)責(zé)流水線中的一級（圖6）。

圖6 地層網(wǎng)格數(shù)據(jù)控制流圖示

2.4 基于外存的海量地層數(shù)據(jù)調(diào)度與可視化技術(shù)

目前基于外存的數(shù)據(jù)調(diào)度和繪制方面的研究主要集中在地形數(shù)據(jù)的顯示方面，特別是數(shù)字高程模型的三維可視化顯示［32－34］，而針對三維地質(zhì)模型的地層網(wǎng)格數(shù)據(jù)的研究則相對較少。我們設(shè)計了海量地層數(shù)據(jù)的調(diào)度與可視化框架（圖7），以外部磁盤上的海量數(shù)據(jù)文件組織與管理為基礎(chǔ)，兼顧磁盤讀寫效率與用戶交互性能，實現(xiàn)了層狀實體模型的三維可視化顯示和線、道、折線方向的剖面數(shù)據(jù)快速提取與顯示。

圖7 基于外存的數(shù)據(jù)調(diào)度與可視化框架

2.4.1 多分辨率數(shù)據(jù)調(diào)度策略

海量的地層網(wǎng)格數(shù)據(jù)經(jīng)過多分辨率組織與存儲技術(shù)處理后，形成了不同精細(xì)程度的幾個層次細(xì)節(jié)模型，每一層次的模型又在空間上以數(shù)據(jù)分塊的形式進(jìn)行組織并存儲。繪制時可根據(jù)當(dāng)前的顯示參數(shù)和所用機(jī)器的硬件性能選擇一個合適的細(xì)節(jié)層次，同時也不需要將該層次上所有數(shù)據(jù)塊全部調(diào)入主存，只需根據(jù)當(dāng)前的顯示參數(shù)實時調(diào)入正在或?qū)⒁L制的部分?jǐn)?shù)據(jù)塊。我們采用了一種基于優(yōu)先級隊列的自頂向下的多分辨率數(shù)據(jù)調(diào)度策略：優(yōu)先級隊列采用最大堆實現(xiàn)，存放當(dāng)前進(jìn)入調(diào)度系統(tǒng)的數(shù)據(jù)塊索引，按照優(yōu)先權(quán)函數(shù)值進(jìn)行排列，越靠近隊列頭的數(shù)據(jù)塊就越應(yīng)該在調(diào)度的下次操作中被細(xì)化。每當(dāng)有新的數(shù)據(jù)塊索引從優(yōu)先級隊的尾部入隊，或具有最大優(yōu)先權(quán)的數(shù)據(jù)塊索引從優(yōu)先級隊的頭部出隊，優(yōu)先級隊列都要進(jìn)行自上向下或自下向上的堆調(diào)整以保證正確的優(yōu)先權(quán)排序。

2.4.2 緩存管理方案

基于外存的數(shù)據(jù)組織技術(shù)要求數(shù)據(jù)調(diào)度模塊不斷地訪問外存設(shè)備獲取新的數(shù)據(jù)來更新當(dāng)前顯示的舊數(shù)據(jù)。從用戶瀏覽地質(zhì)體數(shù)據(jù)的方式來看，當(dāng)前可見的數(shù)據(jù)塊在不久的將來依然可見的可能性是相當(dāng)高的（幀－幀相關(guān)性），也就是說對外存數(shù)據(jù)的訪問具有時間局部性。因此，為了減少磁盤I／O交互次數(shù)，有效利用之前已經(jīng)讀入內(nèi)存的數(shù)據(jù)，我們設(shè)計了一個針對地層LOD模型的緩存管理方案，將內(nèi)、外存數(shù)據(jù)交互這個瓶頸對整個多分辨率數(shù)據(jù)調(diào)度與可視化系統(tǒng)性能帶來的不利影響降低到最小。

2.4.3 利用頂點緩沖的大規(guī)模繪制

經(jīng)過調(diào)度系統(tǒng)輸出的數(shù)據(jù)塊隨后進(jìn)入裝配與繪制過程。我們將需要繪制的數(shù)據(jù)塊中的控制點在空間上進(jìn)行拓?fù)潢P(guān)系檢查，以四邊形網(wǎng)格的形式作為線、道、折線方向提取的剖面及其屬性。為了更有效地將網(wǎng)格數(shù)據(jù)提交到渲染引擎，我們將剖面數(shù)據(jù)組織成頂點、索引、屬性3個數(shù)據(jù)流的形式，利用OpenGL提供的頂點數(shù)組機(jī)制和頂點緩沖對象的擴(kuò)展技術(shù)，減少OpenGL函數(shù)調(diào)用的次數(shù)，提高網(wǎng)格數(shù)據(jù)的傳輸效率，減少了每幀渲染時間，提高了幀率，保證數(shù)據(jù)顯示的實時性（圖8）。

圖8是對一個20GB數(shù)據(jù)量的三維地層網(wǎng)格數(shù)據(jù)進(jìn)行三維可視化測試的結(jié)果，從圖8a到圖8d的變化可以看出：在不適用緩存的情況下，幀率極不穩(wěn)定，經(jīng)常出現(xiàn)嚴(yán)重的掉幀，顯示過程中頻繁出現(xiàn)畫面停滯，嚴(yán)重影響了軟件交互的性能。隨著緩沖區(qū)中緩存數(shù)據(jù)塊的增加，幀率的穩(wěn)定性有較大幅度提高，當(dāng)緩存數(shù)據(jù)塊達(dá)到1 000個以上，雖仍有部分輕微的掉幀情況出現(xiàn)，但數(shù)據(jù)繪制的過程自然流暢，用戶交互操作反映較為靈敏，基本滿足了三維地質(zhì)模型的實時可視化要求。

圖8 不同緩存方案的幀率曲線

3 實際數(shù)據(jù)的算法驗證

利用本文提出的方法對某實際地震構(gòu)造解釋數(shù)據(jù)進(jìn)行了層位面和斷層面的建模，并在地層框架模型的基礎(chǔ)上建立了地層實體模型，通過三維構(gòu)造建模的全過程對以上所述相關(guān)算法進(jìn)行了系統(tǒng)驗證。

試驗工區(qū)的三維地質(zhì)建模包括了3個層位面和若干個斷層面，每個層面的離散數(shù)據(jù)點數(shù)為130×104～150×104。圖9給出了其中一個層位面（共有133×104個數(shù)據(jù)點）的解釋數(shù)據(jù)在不同分辨率時的建模效果。從圖9a到圖9d可見，不同分辨率的三角網(wǎng)格數(shù)據(jù)空間形態(tài)和拓?fù)潢P(guān)系一致，并且在硬件條件相同的情況下，計算效率優(yōu)于某商業(yè)軟件得到的結(jié)果。

圖9 相同層位不同分辨率三維顯示結(jié)果

圖10給出了試驗工區(qū)地層實體模型不同分辨率數(shù)據(jù)三維可視化顯示結(jié)果。該模型地層網(wǎng)格的精度為0.29m×0.29m；地層網(wǎng)格的層間距為0.29m；數(shù)據(jù)塊的數(shù)量為84 105塊；原始數(shù)據(jù)量為22 047.6MB。為了三維可視化顯示的流暢性，建立了7級LOD多分辨率數(shù)據(jù)，多分辨率數(shù)據(jù)的實際數(shù)據(jù)量為42.02GB，大約是原始數(shù)據(jù)量的1.91倍，生成7級LOD模型共耗時1 440s。

從圖10中可以看出，本文提出的地層網(wǎng)格多分辨率組織技術(shù)，不僅保證了地層網(wǎng)格組成的地質(zhì)體在不同分辨率的情況下幾何形態(tài)和拓?fù)潢P(guān)系一致，同樣，也保證了地層網(wǎng)格的屬性（圖10中的顏色顯示）在不同分辨率數(shù)據(jù)中，空間變化規(guī)律一致。本文提出的地層網(wǎng)格多分辨率數(shù)據(jù)組織方法不僅能夠支撐三維可視化顯示，而且也為三維地質(zhì)模型的交互操作提供了便利（圖11）。

4 結(jié)束語

1）針對海量地質(zhì)三角網(wǎng)格曲面的實時顯示問題，本文提出了一種基于二次誤差和邊簡化操作的三角網(wǎng)格多分辨率模型的數(shù)據(jù)組織及其創(chuàng)建算法。該算法首先計算模型中所有邊的簡化代價，并從小到大進(jìn)行排序；然后依次簡化模型中三角網(wǎng)格的邊，進(jìn)行向下采樣；最后在恢復(fù)分辨率的同時，創(chuàng)建該三角網(wǎng)格的多分辨率模型。與其它算法相比，本文提出的算法可以處理海量三角網(wǎng)格，不受內(nèi)存等系統(tǒng)資源的限制。

2）針對海量地層網(wǎng)格數(shù)據(jù)的特點，本文提出了一種基于空間八叉樹結(jié)構(gòu)分塊處理非規(guī)則化地層實體模型的多分辨率數(shù)據(jù)組織方法，該方法綜合運用數(shù)據(jù)分塊和層次細(xì)節(jié)模型等技術(shù)手段，為海量地層網(wǎng)格數(shù)據(jù)生成適合于多分辨率應(yīng)用的LOD模型，并以此建立了地層網(wǎng)格數(shù)據(jù)的多層級索引，實現(xiàn)了海量地層網(wǎng)格數(shù)據(jù)的高效存取。

3）以海量地層網(wǎng)格數(shù)據(jù)組織與管理為基礎(chǔ)，本文提出了一種針對海量地層網(wǎng)格數(shù)據(jù)的調(diào)度與可視化框架，充分運用基于外存技術(shù)的各種手段，實現(xiàn)了包含地層屬性的層狀實體模型的線、道、折線方向的剖面數(shù)據(jù)快速提取與顯示，同時有效提升了大規(guī)模數(shù)據(jù)調(diào)度時的磁盤讀寫效率，改善了大規(guī)模數(shù)據(jù)繪制時的交互操作性能。

4）三維地質(zhì)模型海量數(shù)據(jù)的組織管理和實時三維可視化是三維地質(zhì)建模的關(guān)鍵技術(shù)，本文從用戶體驗出發(fā)，解決了海量三角網(wǎng)格數(shù)據(jù)、海量地層網(wǎng)格數(shù)據(jù)等三維地質(zhì)模型基本數(shù)據(jù)的多分辨率數(shù)據(jù)組織與三維可視化實時顯示的關(guān)鍵問題。上述技術(shù)與其它三維建模技術(shù)（地質(zhì)曲面重建、地質(zhì)曲面拓?fù)浞治?、地層網(wǎng)格剖分和相控屬性建模等技術(shù)）是研發(fā)面向復(fù)雜油氣藏的三維地質(zhì)建模軟件的技術(shù)基礎(chǔ)，這些技術(shù)的發(fā)展將進(jìn)一步促進(jìn)地震綜合解釋和地質(zhì)建模的一體化。

［1］魏嘉.地質(zhì)建模技術(shù)［J］.勘探地球物理進(jìn)展，2007，30（1）：1－6 Wei J.Review of geologic model building techniques［J］.Progress in Exploration Geophysics，2007，30（1）：1－6

［2］魏嘉，唐杰，岳承琪，等.三維地質(zhì)構(gòu)造建模技術(shù)研究［J］.石油物探，2008，47（4）：319－327 Wei J，Tang J，Yue C Q，et al.Study of 3－D geological structure model building［J］.Geophysical Prospecting for Petroleum，2008，47（4）：319－327

［3］Haldorson H，Lake L.A new approach to shale management in filed scale simulation models［J］.The Journal of Petroleum Technology，1984，24（4）：1－4

［4］Mallet J L.Discrete smooth interpolation［J］.ACM Transactions on Graphics，1989，8（2）：121－144

［5］Mallet J L.Discrete smooth interpolation in geometric modeling［J］.Computer Aided Design，1992，24（4）：177－191

［6］Euler N.A new tool to seal 3Dearth model：a cut with constraint［J］.Expanded Abstracts of 68thAnnual Internat SEG Mtg，1998，321－324

［7］楊欽，徐永安，陳其明，等.三維約束Delaunay三角化的研究［J］.計算機(jī)輔助設(shè)計與圖形學(xué)學(xué)報，2000，12（8）：590－594 Yang Q，Xu Y A，Chen Q M，et al.Research on 3D constrained delaunay triangulation［J］.Journal of Computer Aided Design and Computer Graphics，2000，12（8）：590－594

［8］Thore P，Shtuka A，Lecour M，et al.Structural uncertainties：determination，management，and applications［J］.Geophysics，2002，67（3）：840－852

［9］王秀閩，劉洪，孟小紅，等.用于塊體劃分的最小回路法及模型試驗［J］.石油物探，2007，46（3）：249－253 Wang X M，Liu H，Meng X H，et al.Model test and study on minimum circuits for block partition［J］.Geophysical Prospecting for Petroleum，2007，46（3）：249－253

［10］Lorensen W E.Marching cubes：a high resolution 3D surface construction algorithm［J］.Computer Graphics，1987，21（4）：163－169

［11］袁國棟，施寅，勵曉建.一種新的抽取等值面的四面體分解方法［J］.計算機(jī)工程與科學(xué)，2002，20（6）：48－51 Yuan G D，Shi Y，Li X J.A new tetrahedral decomposition scheme for extracting isosurfaces［J］.Computer Engineering ＆ Science，2002，20（6）：48－51

［12］Chambers H，Brown A L.3－D visualization continues to advance integrated interpretation environment［J］.First Break，2003，21（5）：31－34

［13］Hodgson P，Cunnell C，Krueger A，et al.Cluster visualization rises to the seismic processing challenge［J］.First Break，2005，23（3）：33－35

［14］Gray G.Observations of azimuthal anisotropy in prestack seismic data［J］.Expanded Abstracts of CSPG ＆ CSEG Convention，2007，373－377

［15］Molnar S，Cox M，Ellsworth D，et al.A sorting classification of parallel rendering［J］.IEEE Computer Graphics and Applications，1994，14（4）：23－32

［16］Varadhan G，Manocha D.Out－of－core rendering of massive geometric environments［J］.Expanded Abstracts of Proceedings of IEEE Visualization Mtg，2002，69－76

［17］Imma B，Isabel N，Roberto S.Multiresolution volume visualization with a texture－based octree［J］.The Visual Computer，2001，17：185－197

［18］Guthe S，Wand M，Gonser J，et al.Interactive rendering of large volume data sets［J］.Expanded Abstracts of Proceedings of IEEE Visualization Mtg，2002，104－115

［19］Schneider J，Westerm R.Compression domain volume rendering［J］.Expanded Abstracts of Proceedings of the IEEE Visualization Mtg，2003，293－300

［20］魏嘉，唐杰，武港山，等.地震疊前數(shù)據(jù)三維可視化技術(shù)探討［J］.勘探地球物理進(jìn)展，2009，32（1）：12－17 Wei J，Tang J，Wu G，et al.On 3－D visualization of prestack seismic data［J］.Progress in Exploration Geophysics，2009，32（1）：12－17

［21］魏嘉，唐杰，武港山，等.三維地震數(shù)據(jù)體多分辨率數(shù)據(jù)組織與管理技術(shù)研究［J］.石油物探，2010，49（3）：240－244 Wei J，Tang J，Wu G S，et al.Multi－resolution organization and management of massive 3Dseismic data volume［J］.Progress in Exploration Geophysics，2010，49（3）：240－244

［22］Lindstrom P，Koller D，Ribarsky W，et a1.Real－time continuous level of detail rendering of height fields［J］.Computer Graphics，1996，30（1）：109－118

［23］Lindstrom P，Pascucci V.Terrain simplification sim－plified：ageneral framework for view－dependent outof－core visualization［J］.IEEE Transactions on Visualization and Computer Graphics，2002，8（3）：239－254

［24］Lindstrom P，Isenburg M.Lossless compression of hexahedral meshes［C］∥Expanded Abstracts of IEEE Data Compression Conference.New York：IEEE Publishing，2008：192－201

［25］Hoppe H.Progressive meshes［C］∥Expanded Abstracts of Computer Graphics Mtg.New York：Association for Computing Machinery，1996：99－108

［26］Hoppe H.View－dependent refinement of progressive meshes［C］∥Expanded Abstracts of Computer Graphics Mtg.New York：Association for Computing Machinery，1997：189－198

［27］Gueziec A.Surface simplification with variable tolerance［C］∥Expanded Abstracts of 2ndAnnual Internat Symposiumon Medical Robotics and Computer Assisted Surgery Mtg.Wiley：Orthopaedic Research Society，1995：132－139

［28］Duchaineauy M，Wolinsky M.Roaming terrain：realtime optimally adapting meshs［C］∥Expanded Abstracts of IEEE Visualization Mtg.New York：IEEE Publishing，1997：81－88

［29］Rottger S，Heidrich W，Slussallek P，et al.Rea1－time generation of continuous levels of detail for height fields［J］.Expanded Abstracts of Proceedings of the 6thInternational Conference in Central Europe on Computer Graphics and Visualization，1998，315－322

［30］Samet H.The design and analysis of spatial data structures［M］.New Jersey：Addison－Wesley，1990：223－235

［31］Boada I.An octree－based multiresolution hybrid framework［J］.Future Generation Computer Systems，2004，20（8）：1275－1284

［32］Vitter J S.External memory algorithms and data structures：dealing with massive data［J］.ACM Computing Surveys，2001，33（2）：209－271

［33］Ulrich T.Chunked LOD：rendering massive terrains using chunked level of detail control［C］∥Computer Graphics Proceedings Annual Conference Series Expanded Abstracts.Texas，USA：ACM，2002：202－209

［34］Vo H T，Callahan S P，Smith N，et al.iRun：Interactive rendering of large unstructured grids［C］∥Eurographics Symposium on Parallel Graphics and Visualization Expanded Abstracts.Geneva：The Eurographics Association，2007：93－100