一種基于反距離加權方法的層狀三維地質界面擬合算法

李 楠，雷玲玲，肖克炎，

(1．中國地質科學院礦產資源研究所，北京 100037;2．長江大學地球科學學院，湖北 荊州434023)

一種基于反距離加權方法的層狀三維地質界面擬合算法

李 楠1，雷玲玲2，肖克炎1，2

(1．中國地質科學院礦產資源研究所，北京 100037;2．長江大學地球科學學院，湖北 荊州434023)

三維空間中地質曲面模擬是礦產勘查、大比例尺成礦預測等地學領域的重要研究內容之一。與地球化學、地球物理數據的曲面擬合相比，地質界面擬合更加復雜且具有一般性。從功能需求角度出發，地質界面擬合需要處理層狀地層、不整合地層以及存在構造的地層等多種情況;從算法實現角度講，該擬合過程不同于物化探數據通過某一屬性值進行曲面擬合的情況，因為在地質界面模擬過程中，不存在這樣的特征值進行曲面擬合，因此不能直接使用物化探數據等值面提取的思想模擬地質接觸界面。針對上述問題，提出了一種基于反距離加權法移動立方格算法的層狀地層曲面擬合算法。該方法模型適合絕大多數層狀地學三維曲面模擬。

地質界面;空間插值;移動立方格算法;等值面;隱式地質建模

0 引言

在地質勘探過程中，勘探人員通過鉆孔、槽探、坑道、填圖和編錄等技術手段可以獲取地質體和地質界面的三維控制點信息。在計算機領域，這些控制點被稱為空間散亂點或離散點數據。由離散點集出發，重構該點集所屬的曲面，稱為曲面重建。筆者所討論的地質界面重建屬于該類問題。曲面重建，即構造出對象幾何模型，是當今三維數據場建模的核心問題之一。曲面重建問題抽象為：給定三維空間曲面S上的一組采樣點P，通過數學模型找到一個曲面 S'，使得 S'合理地逼近 S(Hoppe et al，1992)。其主要包括兩大類算法：一類是參數曲面模型法，另一類是曲面提取法。

參數曲面模型法主要包括逼近樣條曲面和插值樣條曲面兩種類型。其中，插值樣條通過離散點，而逼近樣條不過離散點。該方法的優點是曲面的獲得不是由用戶給定某一特征值得到的，而是由某一數學模型擬合得到的。但無論上述哪種方法，生成的地質界面均為單值面，不能很好地滿足地質實際需要。文獻對上述方法做了詳細分析(毛先成，2006)。

曲面提取方法是地學曲面重構的主要研究內容。其核心思想是按照給定的閾值從體數據中抽取空間曲面，然后再按照傳統方法進行繪制。這類方法主要包括輪廓線法(Keppel，1975)，等值面法(Lorensen et al，1987)，幾何變換(Boissonnat，1984)等方法。提取等值面法由于其操作簡單、精確度較高等原因一直是曲面重建研究的重要內容。這種方法的核心思想是用戶給定一組值，算法根據這組值提取對應的等值面。

地質體接觸界面提取與其他地學數據等值面提取的主要區別是接觸界面提取的對象往往是一組地質認識信息而不是某一特征值，如，地球化學、地球物理等數據的分析值等。

綜上所述，在上述兩種曲面擬合算法思想的基礎上，提出基于反距離加權方法和移動立方格(Marching Cubes，MC)算法的層狀三維地質界面擬合算法。該算法的核心思想：首先由離散點進行空間插值得到三維柵格模型;其次建立增量函數f(Δ z)=zIDW－z，求出最優擬合曲面值 f(Δz)=ziso;最后，基于Marching Cubes算法和ziso追蹤得到等值面S，即為層狀地質界面。

1 算法思想

1．1 算法基本概念

介紹一些相關的概念及術語定義以便后文論述使用。

(1)曲面擬合(Edward，2003)。所謂曲面擬合，就是根據實際試驗測試數據，求取函數f(x，y，z)與變量x，y，z之間的解析式，使其所確定的曲面通過或近似通過所有的實驗測試點。也就是說，是所有實驗測試點擬合的曲面M'能夠合理逼近控制點所在曲面M。

(2)地理學第一定律(Tobler，1970)。任何事物之間均相關，而距離較近的事物總是比距離較遠的事物影響更大。這里的距離不僅包含空間上的，也包含時間上的距離。

(3)移動立方格算法(Lorensen et al，1987)。移動立方格算法本質是從一個三維數據場中提取出一個等值面，所以也稱為“等值面”算法。其算法思想是在規模大小為Nx×Ny×Nz的體數據場中，由每8個頂點組成1個體元，每個頂點都有1個數值，體元每條棱上的數值都是線性變化的。等值面的生成就是根據兩相鄰頂點的增量函數值生成等值點，然后將每個體元中的所有等值點按組合方式連接起來，生成指定值的等值面的過程。在標量場曲面的增量函數建立完成后，可以使用上述方法追蹤地質體曲面。

(4)K-近鄰查詢。給定查詢點和正整數K，從數據集中距離查詢點最近的K個數據，當K=1時，為最近領域查詢。

1．2 算法主要思想

筆者提出基于反距離加權方法和移動立方格算法的層狀三維地質界面擬合算法。反距離加權方法假設未知點的值受較近控制點的影響比較遠控制點影響更大。影響程度用點之間距離乘方的倒數表示。這種特性使得該方法可以應用于模擬連續曲面，最經典的例子是構建DTM模型。在三維空間中，盡管空間插值的結果是規則網格，但其數學本質不會發生變化，即較近控制點的影響比較遠控制點影響更大。因此，如果待插值為高程值z，則空間插值的結果表示由采樣點擬合的曲面M。此時就可以使用移動立方格算法提取曲面M'，使之能夠合理逼近控制點所在曲面M。

1．3 算法流程

如1．2節所述，總結筆者提出的算法總體流程如下。

步驟1：基于反距離加權(IDW)方法對地層采樣點的高程值z進行空間插值。

步驟2：建立插值結果與待插點實際高程值z'的增量函數f(Δz)。

步驟3：選擇增量函數f(Δz)的某個值，作為曲面提取值。

步驟4：應用移動立方格(MC)法，提取f(Δz)的某個值處的等值面。

步驟5：計算曲面光照、設置顏色、重建曲面拓撲結構等。

綜上所述，根據算法的主要思想和總體流程，總結算法流程如圖1所示。

圖1 算法總體流程圖

2 算法實現

本算法主要包含兩個主要內容：(1)快速三維空間插值算法;(2)移動立方格算法。

在三維空間中，隨著點數量的增加，K-近鄰查詢非常耗時。例如，在三維空間中對12 000個點進行反距離加權空間插值，設待插點達到1 000萬，如果使用線性方法進行方法其時間復雜度將達到0(12 000×10 000 000)。可以看出，上述方法的時間復雜度是用戶無法接受的。

筆者提出使用八叉樹和R-Tree樹兩種空間索引結構，完成三維散亂點云的K-近鄰查詢，提高三維空間插值算法的效率。

2．1 算法1：三維離散點快速空間插值算法

輸入：離散點集合 P{(x，y，z)|x，y，z為任意實數}

輸出：柵格體S

步驟1：建立三維點云P的八叉樹空間索引。步驟2：建立八叉樹空間索引的三維R-Tree樹索引。

步驟3：通過矩陣計算獲得待插點的搜索范圍。

步驟4：應用反距離加權方法對地層分界點的高程值z進行空間插值，得到S。

使用移動立方格算法，提取地質曲面。

在移動立方格方法中，由于每個體元存在8個頂點，每個頂點可能有正和負兩種狀態(等值面經過0狀態)，因此，每個體元按照其頂點的正、負分布，共有28=256種不同的組合狀態。盡管判斷等值面將與哪些體元相交在理論上很容易理解，但是，根據這256種狀態求出每個體元的等值面是個相當繁瑣且容易出錯的過程。實際工作中，根據互補對稱性和旋轉對稱性可以將這256種狀態簡化為15種情況(蔣先剛，2003)(圖2)。但是，簡化的模型是一個存在二義性的模型，在文中采用漸近判別法(Nielson et al，1991)來消除二義性。

圖2 模型簡化的15種情況

2．2 算法2：提取地質分界曲面

輸入：柵格體S，待提取等值面的值ziso

輸出：三角網Surface

步驟1：查找包含值ziso的立方格。

步驟2：按照圖2的方式生成三角面片。

步驟3：計算每一個三角形的法向量。

步驟4：結束。

3 實驗與分析

使用西藏某礦區三維地層采樣數據作為實驗數據，該礦區由鉆孔獲得的地層分界界面采樣數據是416個。實驗環境、硬件開發環境：CPU主頻為2.16G Hz和內存為1G的PC機;軟件開發環境：MS Windows XP操作系統和Visual C++7．0。算法效果如圖3、圖4、圖 5、圖 6、圖7 所示。

4 結論

在已有曲面擬合算法思想的基礎上，筆者提出了基于反距離加權方法和移動立方格算法的層狀三維地質曲面擬合算法。本算法的優點在于以數學原理與地理學第一定律為基礎，實現層狀地層的快速、精確擬合。該方法模型能夠適合絕大多數層狀地學三維曲面模擬。不足之處是不能解決不整合、構造等復雜情況，這也是本算法下一步研究的方向。

圖3 地層接觸界面采樣點

圖4 柵格化結果

圖5 地層分界界面擬合結果

圖6 地層分界界面擬合與剖切結果對比

圖7 地質界面與化探異常以及礦體復合顯示

毛先成．2006．三維數字礦床與隱伏礦體立體定量預測研究［D］．長沙：中南大學．

蔣先剛．2003．三維數據場重構與顯示工程軟件設計［M］．北京：萬水出版社．

BOISSONNAT J D．1984．Geometric structures for three-dimensional shape represent ation［J］．ACM Trans Graphics，3(4)：266－286．

EDWARD ANGLE．2003．Interactive Computer Graphics：A Top-Down Approach Using OpenGL［M］．Beijing：Higher Education Press．

HOPPE H，DEROSE T，DUCHAMP T，et al．1992．Surface reconstruction from unorganized points［J］．Computer Graphics，26(2)：71－78

KEPPEL E．1975．Approximating complex surfaces by triangulation of contour lines［J］．IBM Journal of Research and Development，19(1)：2 －11．

LORENSEN W E，CLINE H E．1987．Marching Cubes：A high Resolution 3D Surface Construction Algorithm ［J］．Computer Graphics，21(4)：163 －169．

NIELSON G M，HAMANN B．1991．The asymptotic decider：resolving the ambiguity in marching cubes［C］//Proceedings of Visulization'91，San Diego：CA，83－91．

TOBLER W．1970．A computer movie simulating urban growth in the Detroit region［J］．Economic Geography，46(2)：234－240．

Stratiform 3D geological interface fitting algorithm based on inverse distance weighted method

LI Nan1，LEI Ling-ling2，XIAO Ke-yan1，2

(1．Institute of Mineral Resources，Chinese Academy of Geological Sciences，Beijing 100037，China;2．School of Geosciences，Yangtze University，Jingzhou 434023，Hubei)

Geological curved surface modeling was an important research content in fields of mineral exploration and large-scale mineral assessment．Geological curved surface modeling was more complex，and of generality compared with that of other geological data like geochemical data or geophysical data．Geological curved surface modeling needed to solve many complex situations such as unconformity or faults．In response to these issues，the authors proposed an algorithm of fitting layered geological surface that was based on two algorithms of inverse distance weighted method(IDW)and marching cubes(MC)．The algorithms were characteristic of speedy，accurate fitting for stratified layers on the foundations of mathematical principles and the first theorem of geography，suited most 3D curved surface modeling of the layered earth sciences．

Geological interface;Spatial interpolation;Marching cubes algorithm;Iso-surface;Implicit geological modeling

TP391

A

1674－3636(2012)03－0291－05

10．3969/j．issn．1674-3636．2012．03．291

2012－06－06;編輯：陸李萍

國家高技術研究發展計劃(863)“礦產資源評價數字礦床模型及可視化技術研究”(2006AA06Z114)、國家科技支撐計劃“西部優勢礦產資源潛力評價技術及應用研究”(2006BAB01A01)、礦產勘查三維預測評價信息平臺開發及應用示范研究(1212010012013)、中央國家機關基本業務費基金項目等基金資助

李楠(1980— )，男，副研究員，博士，主要從事地圖制圖與地理信息工程、大比例尺礦產資源評價研究，E-mail：superln1980@yahoo．com．cn