基于八叉樹和混合搜索樹的地質曲面快速求交方法

孫黎明 魏迎奇 蔡紅 嚴俊 宋建正 喬蕓蕓

摘要：為處理地質界面之間的空間相交關系，提出一種新的針對三角地質曲面的快速求交方法。該方法融合優化八叉樹法和OBB搜索樹方法，可以更快速準確地剔除遠離交線的其他三角形。求交剩余的三角形得到交線，應用三角網局部重構和網格優化算法修正交線附近的三角網，最終分割交線兩側的地質曲面，完成2個地質曲面的離散化求交過程。與AABB、OBB和空間分解法相比，該方法在大數據量三角曲面求交中效率優勢明顯，可以快速準確處理地質模型構建和分析中的曲面求交問題，為三維地質模型自動化構建的實現提供有效支撐。

關鍵詞：地質建模；地質曲面；曲面求交；八叉樹；OBB搜索樹；三角網重構；模型切割；區域離散化

中圖分類號：P221.1

文獻標志碼：B

0 引 言

三角形的地質曲面是構成三維地質模型的基本幾何元素，在大區域復雜地質模型的自動化構建中，通常需要處理地質界面之間的空間相交關系，這是三維地質建模應用研究的熱點問題之一。[1-6]在建模過程中，不整合地層邊界計算、斷層面與地層的交線和任意切割面與地層的交線等[7-8]應用廣泛。三維地質模型的分析應用，例如剖面制作、巷道模擬和地下工程等[9]，需要準確快速得到相交區域的空間位置、拓撲信息和內部構造，以對模型區域有更準確的理解，同時也需要進行多次曲面求交和計算分析。[10]隨著地質建模數據量越來越龐大，尤其是地震數據重構斷層和地層時，三維地質建模需要一種快速穩定的三角曲面求交算法。地質曲面的求交過程根據2個已有的三角曲面集合，得到相交曲面之間的交線，并重構2個三角曲面和交線兩側的三角形，完成曲面的離散化。如果循環判斷求交，那么所有的三角形都要參與求交計算，算法的效率為O（n2）。

實際進行求交的三角形對直接影響三角曲面求交速度，因而算法的核心是快速準確地找到與交線相關的三角形。目前，曲面求交和物體碰撞檢測的方法主要有空間分解法（八叉樹、BSP樹和KD樹）和層次包圍盒法（AABB、OBB、球體包圍盒和k-DOP）等。地質曲面求交通常直接采用物體碰撞檢測的方法，LINDENBECK等[11]開發TRICUT軟件，用Rapid庫解決三角網求交問題，應用廣泛。CGAL計算幾何算法庫改進了AABB樹求交方法 [12]，基于OBB的地質曲面切割方法也廣泛應用并取得良好效果[13]，但在數據量較大時，構建層次搜索樹十分耗時。YANG等[8]提出一種基于混合包圍盒的針對真三維模型的切割方法，可解決多種幾何體一體切割的問題。ELSHEIKH等[14]提出一種基于追蹤法的可靠地質曲面求交方法，需要利用三角網拓撲追蹤曲面交線，但尋找起始追蹤點較困難。

本文提出一種改進的地質曲面快速求交方法，用結構化的八叉樹算法與改進的OBB混合搜索樹相結合的方法快速求取三角曲面交線，并重構交線所影響的三角形，然后將算法應用于地質建模分析的典型問題中，包括曲面切割地質模型，斷層面求交和地質區域離散化等。為將復雜地質曲面求交問題簡化，將地層求交分解為任意兩兩三角曲面求交。

1 理論方法和求交步驟

檢測2個三角曲面的相交部分，核心是快速找到相交的三角形。無論曲面三角形數量多大，實際相交的三角形都很少，直接構建層次搜索樹效率很低，因此快速剔除大部分不相交的三角形，只留下交線附近區域的三角形，然后再構建層次搜索樹，可大大提升求交速度。采用結構化的八叉樹數據結構高效剔除求交曲面中的大多數不相交三角形，構建基于改進OBB的搜索樹，快速準確地求取交線，改進曲面求交方法流程見圖1。采用靜態結構化數據結構[15-16]表達1個三角曲面，即用1個大于三角網的六面體等分后表達地質界面和2個曲面相交的三角形，建立特殊的層次二叉搜索樹，根節點根據三角形的AABB表示，葉節點則用OBB代替，最后對三角形求交得到交線。

1.1 優化八叉樹的粗略求交

空間八叉樹是二維平面中四叉樹在三維空間的擴展，是用均勻細分空間描述三維場景中幾何對象的一種樹狀數據結構，其每個父節點用一個正方體表示，每個父節點都有8個子節點，是將父節點均勻細分為8個小正方體得到的。八叉樹通常用于三維空間中海量物體的快速查詢和空間索引，目前廣泛應用于大數據量的三維GIS，包括三維智慧城市、三維數字管線、三維BIM可視化和3D大場景游戲等。

八叉樹的特點是可以用多層樹快速索引每個節點，但對于海量數據離散點構建的地質曲面和細分到單一節點的八叉樹，每次計算的迭代數量龐大，且八叉樹的層級非常多，計算效率更低，每次循環計算耗時很長。本文提出一種非對稱的空間八叉樹模型，用于快速查找2個三角曲面交線所在的三角形。與傳統空間八叉樹模型不同，該模型不是每層都按照8n細分，而是達到一定的邊界條件就采用結合OBB的方法進行求交。設計八叉樹的邊界條件，計算八叉樹停止的條件，以省略八叉樹細分到葉節點的過程，提高求交的計算效率。優化八叉樹求交過程見圖2。

地質曲面所在的三維笛卡爾坐標系定義為R，定義三角網中的點VR、三角形TR、地層三角網SR。任意點P（x， y， z）對應立方體ER，立方體的邊長為三角網包圍盒在x、y和z這3個方向上長度的最大值l。依次將每個父節點分為8個子節點，葉節點滿足與同行列子節點內的曲面點個數接近條件，并且當個數最小時八叉樹就不再細分。與交線不相關的不相交三角形見圖3，剔除遠離交線后的三角網見圖4。

利用八叉樹結構計算1個子節點內的三角形，可以很快找到距離另一個曲面最近的三角形集合（見圖3）。顯然，該方法得到的結果仍然包括許多不相交的三角形，如果單用此算法求交，效率不高，因此本文結合改進OBB層次搜索樹的方法進一步剔除不相交的三角形。

1.2 OBB混合搜索樹的精確求交

OBB是1個緊密的包圍盒，是包含該物體且相對于坐標軸任意方向的最小六面體，可提高檢測效率。2個相交地質曲面的OBB包圍盒見圖5。

構建OBB搜索樹和計算OBB比計算AABB更繁瑣耗時。根據幾何坐標計算物體的OBB，因為三角形OBB需要多2個向量計算，所以大數據量時直接構建OBB搜索樹復雜低效。GOTTSCHALK等[17]提出一種有效的實現方法，采用AABB構建中間所有層次的根節點，將第1.1節中計算的六面體包圍盒直接作為AABB，而葉節點采用OBB方法，OBB混合搜索樹見圖6。

搜索樹采用二叉樹結構存

儲，最終葉節點采用OBB求交，能最大限度減少精確判斷相交的三角形對。綜合利用2種包圍盒的優點，可減少計算量，加快搜索樹的構建速度。[18]

根據可能相交的三角形構建層次搜索樹的簡要步驟如下：

（1）選擇根節點。將所有的三角形最大包圍盒投影到3個坐標軸，找到3個方向中的最長軸，選擇其中間值為最大值與最小值之和的1/2。

（2）構建搜索樹。根據平均三角形AABB的大小，選擇樹的深度層次。樹的深度選擇是影響搜索求交速度的重要因素。經過第1.1節剔除工作后，求交的三角形只占曲面三角形總數的小部分，層次不需太深。

（3）遞歸分割。同時將每個三角形和點的包圍盒自頂向下遍歷并存放指針到每一層合適的位置節點，直到單一三角形為止，向下一直到葉節點都是三角形，葉節點用OBB代替AABB。

（4）求交。自頂而下判斷包圍盒是否相交，如果遞歸到葉節點，根據三角形的3個點構造OBB，采用SAT方法求交。[19-20]

1.3 交線區域網格重構

對得到的所有三角形對逐一求交，即可得到最終的交線。一對空間三角形相交的結果有3種情況：點、線段和相交三角形。地質曲面求交中最常見的是線段，線段首尾相連形成交線。交線可以有多條，也可能是封閉環。目前，空間三角形求交的算法主要分為標量判別法和矢量判別法[21]，本文選擇適于地質曲面求交的MLLER[22]三角形求交法。

加入新的交線后，將交線點分別加入到交線兩側的三角網中，旋轉變換（不是投影）到合適的二維平面。[23]采用局部三角形修復方法，若交點在三角形的邊上，則在此處添加新點，將這條邊分割；若交點在三角形內部，則在內部添加新點，連接頂點構成新三角形，盡量滿足Delaunay法則，加入新點后的三角網重構見圖7。

新的三角網交線區域必須進行網格質量優化。若新加入的點形成狹長三角形和拓撲錯誤三角形，在地質模型空間分析、體元劃分和數值評估中都無法滿足要求，則需要進行三角網質量優化，局部修改

重構。根據交線點所在的三角形的網格拓撲找到周圍三角形[24-25]，優化網格質量，包括邊界交換、合并三角形、點重構和插入新點。[26]

將局部修復好的三角形重新變換回三維空間，將原來的曲面以交線為邊界分為2個曲面，判斷其他三角形頂點相對交線的位置，若分別屬于不同的曲面，則求取2個曲面的交線。多個地質曲面求交見圖8。

2 實驗與應用

2.1 算法性能測試對比

在一臺型號為Intel i3-2330M、CPU為2.1 GHz、內存為6 GB的計算機上測試算法效率，并與其他2種方法進行對比，見圖9。靜態結構化數據搜索的速度與分割的包圍盒的個數有關，與三角形個數關系不大，時間復雜度為O（2n），保證用于構建OBB混合搜索樹的三角形個數隨求交三角形個數的增加而平穩增長。該改進算法在數據量大的曲面求交中效率優勢明顯。

2.2 改進算法在三維地質建模中的應用

由于地質構造的復雜性，地質曲面之間的空間關系非常復雜，因此在建模過程中，地層與地層或斷層相交時，通常需要求取交線并離散化，最常用的方法是任意剖面切割法。三維曲面求交方法為基礎算法，可以很好地應用到地質自動化建模和分析中。

在三維地質建模中，構建含斷層的地質模型是研究熱點。根據地質剖面圖確定斷點、斷距和傾角等基本信息，構建初始斷層面，得到現有斷層面與地層面的求交線，切割上、下2個地層面，重構地層和斷層面邊界。保證模型和數據的一致性，檢查拓撲關系的正確性是核心。含簡單斷層的地層模型見圖10。

復雜的地質構造有許多界面相交，分割線通過一系列線段表示，不同區域之間用多個線段隔開。為自動化模型重構和模型數值分析[27]，需要分離交線分割的不同區域，手動方法繁瑣且不可靠。為實現自動分離不同區域，嘗試改進前文提出的求交算法以處理簡單區域劃分的問題。先確定交線類型，再根據交線追蹤區域邊界，形成內部環。遇到邊界就停止追蹤，形成多個相鄰多邊形封閉區域。根據三角形和交線的相對位置，劃分所有三角形到其所屬的區域內，完成地質區域離散，見圖11。

3 結束語

為解決三維地質建模過程和模型分析中的三角曲面求交問題，提出基于優化八叉樹的搜索算法和OBB混合搜索樹相結合的曲面求交方法。通過結構化搜索快速剔除遠離交線的不相交三角形，結合OBB混合層次搜索樹，快速求取三角曲面的交線，并將該方法應用到地質建模和分析中，驗證算法的可靠性和實用性。采用優化的八叉樹快速剔除大量不相交三角形，可提升求交算法的效率，除三角曲面外還可以擴展到其他類型曲面。將該算法應用到地質分析中，解決參數曲面切割地質模型、斷層構建和區域離散化問題。將該算法作為基礎算法工具應用到三維地質模型的構建和分析中，為快速自動化建模分析提供一種新的方法。

改進的曲面求交方法可以較好地處理地質曲面求交的相關問題，未來可以嘗試采用GPU并行計算的方法提升運算效率。在解決地質曲面求交相關應用問題上，工作流程自動化和減少人工干預是未來努力的主要方向。

參考文獻：

[1] 鄧飛， 王瑞， 王美平， 等. 復雜三維地層建模及快速射線追蹤的研究與實現[J]. 大慶石油地質與開發， 2007， 26（1）： 113-118. DOI： 10.3969/j.issn.1000-3754.2007.01.028.

[2] LI N， BAGAS L， LINDSAY M， et al. An irregular triangle mesh buffer analysis method for boundary representation geological object in three-dimension[J]. Earth Science Informatics， 2017， 10（2）： 149-167.

[3] FU Q， WU Z， WANG X， et al. An algorithm for finding intersection between ball B-spline curves[J]. Journal of Computational and Applied Mathematics， 2018， 327： 260-273.

[4] WELLMANN J F， de la VARGA M， MURDIE R E， et al. Uncertainty estimation for a geological model of Sandstone greenstone belt， Western Australia-insights from integrated geological and geophysical inversion in a Bayesian inference framework[J]. Geological Society， London， Special Publications， 2017， 453（1）： 41. DOI： 10.1144/SP453.12.

[5] CALCAGNO P， COURRIOUX G， LOPEZ S， et al. How geological architecture helps 3D modelling[C]// Proceedings of 4th Meeting of European 3D GeoModelling Community. Orleans， 2018.

[6] JIANG Q P， JIE T， YUAN C F. Fast triangle mesh surface intersection algorithm based on uniform grid[J]. Computer Engineering， 2008， 34（21）： 172-174.

[7] BISTACCHI A， MASSIRONI M， dal PIAZ G V， et al. 3D fold and fault reconstruction with an uncertainty model： An example from an Alpine tunnel case study[J]. Computers & Geosciences， 2008， 34（4）： 351-372.

[8] YANG Y， LIU X J， ZHOU B， et al. Research of 3D GIS section for true three-dimensional geo-spatial model[C]// Proceedings of 2010 International Conference on Multimedia Technology. Ningbo： IEEE， 2010： 1-4. DOI： 10.1109/ICMULT.2010.5631185.

[9] 孫黎明， 李青元， 譚海， 等. 基于Netgen的層狀地質體四面體網格劃分方法[J]. 計算機輔助工程， 2013， 22（3）： 47-52. DOI： 10.3969/j.issn.1006-0871.2013.03.010.

[10] FRANK T， TERTOIS A L， MALLET J L. 3D-reconstruction of complex geological interfaces from irregularly distributed and noisy point data[J]. Computers & Geosciences， 2007， 33（7）： 932-943. DOI： 10.1016/j.cageo.2006.11.014.

[11] LINDENBECK C H， EBERT H D， ULMER H， et al. TRICUT： A program to clip triangle meshes using rapid and triangle libraries and visualization toolkit[J]. Computers & Geosciences， 2002， 28（7）： 841-850. DOI： 10.1016/S0098-3004（01）00110-8.

[12] FABRI A， PION S. CGAL-computational geometry algorithms library[C]// Proceedings of 17th ACM SIGSPATIAL International Symposium on Advances in Geographic Information Systems. Washington， 2009. DOI： 10.1145/1653771.1653865.

[13] 黃松柏， 徐華. 基于動態OBB層次結構的曲面相交算法[J]. 計算機應用研究， 2011， 28（8）： 3181-3184. DOI： 10.3969/j.issn.1001-3695.2011.08.106.

[14] ELSHEIKH A H， ELSHEIKH M. A reliable triangular mesh intersection algorithm and its application in geological modelling[J]. Engineering with Computers， 2014， 30（1）： 143-157. DOI： 10.1007/s00366-012-0297-3.

[15] SAMET H. Foundations of multidimensional and metric data structures[M]. San Francisco： Morgan Kaufmann， 2006.

[16] DASSI F， PEROTTO S， FORMAGGIA L， et al. Efficient geometric reconstruction of complex geological structures[J]. Mathematics and Computers in Simulation， 2014， 106： 163-184. DOI： 10.1016/j.matcom.2014.01.005.

[17] GOTTSCHALK S， LIN M C， MANOCHA D. OBB tree： A hierarchical structure for rapid interference detection[C]// Proceedings of 23rd Annual Conference on Computer Graphics and Interactive Techniques. New York， 1996： 171-180. DOI： 10.1145/237170.237244.

[18] FANG Z G， XU J， JIANG J X， et al. Efficient collision detection using bounding volume hierarchies of OBB-AABBs and its application[C]//Proceedings of 2010 International Conference on Computer Design and Applications. Qinhuangdao： IEEE， 2010. DOI： 10.1109/ICCDA.2010.5541315.

[19] CHANG J W， KIM M S. Efficient triangle-triangle intersection test for OBB-based collision detection[J]. Computers & Graphics， 2009， 33（3）： 235-240. DOI： 10.1016/j.cag.2009.03.009.

[20] VAN DEN BERGEN G. Efficient collision detection of complex deformable models using AABB trees[J]. Journal of Graphics Tools， 2012，2（4）： 37-41. DOI： 10.1080/10867651.1997.10487480.

[21] 鄒益勝， 丁國富， 何邕， 等. 空間三角形快速相交檢測算法[J]. 計算機應用研究， 2008， 25（10）： 2906-2910. DOI： 10.3969/j.issn.1001-3695.2008.10.007.

[22] MLLER T. A fast triangle-triangle intersection test[J]. Journal of Graphics Tools， 1997， 2（2）： 25-30. DOI： 10.1080/10867651.1997.10487472.

[23] MCLAURIN D， MARCUM D， REMOTIGUE M， et al. Repairing unstructured triangular mesh intersections[J]. International Journal for Numerical Methods in Engineering， 2013， 93（3）： 266-275. DOI： 10.1002/nme.4385.

[24] HOPPE H， DEROSE T， DUCHAMP T， et al. Mesh optimization[C]// Proceedings of 20th Annual Conference on Computer Graphics and Interactive Techniques. New York， 1993. DOI： 10.1145/166117.166119.

[25] ATTENE M， FALCIDIENO B. Remesh： An interactive environment to edit and repair triangle meshes[C]// Proceedings of IEEE International Conference on Shape Modeling and Applications. Matsushima： IEEE， 2006. DOI： 10.1109/SMI.2006.29.

[26] WANG D， HASSAN O， MORGAN K， et al. Enhanced remeshing from STL files with applications to surface grid generation[J]. Communications in Numerical Methods in Engineering， 2007， 23（2）： 227-239. DOI： 10.1002/cnm.894.

[27] 季順迎， 趙金鳳， 狄少丞， 等. 面向環境力學的離散元分析軟件研發和工程應用[J]. 計算機輔助工程， 2014， 23（1）： 69-75. DOI： 10.13340/j.cae.2014.01.014.

（編輯 付宇靚）