土壤孔隙結構三維重建軟件系統的設計與實現

呂菲

土壤結構尤其是土壤孔隙結構特征不僅決定了許多關鍵的土壤物理性質(水分持留、傳導能力等)，而且也控制著發生在土壤中的諸多重要生態過程(如土壤水分運動、養分遷移等)。快速獲取原狀孔隙結構的定量信息是從根本上認識土壤結構與功能間關系及有關過程內在機制的必要前提，也是土壤學工作者長期關注的熱點和難點問題之一。

土壤微形態學研究中通常采用的切片方法不僅制備過程費時費力，而且根據二維的切片信息推斷三維的土壤結構容易出現偏差甚至錯誤。近20多年來，CT(Computerized Tomography)、核磁共振(NMR)等無損、顯微探測技術的長足進步使得高精度數字圖像的快速獲取成為可能，其在醫學、工程技術等領域的成功應用、推廣為我們提供了可供借鑒的新思路。在土壤學研究領域，國內近年也逐漸有學者開展了此方面的研究，如馮杰和郝振純(2002年)、賀秀斌等(2005年)、呂菲等(2007年)等。但就目前國內的研究現狀來看，由于缺少通用的或有針對性的數字圖像分析、三維重建及可視化軟件工具，已經大大制約了研究工作的深入、全面開展。針對這一問題，本文嘗試開發基于計算機圖形學的、集圖像數據管理、數字圖像分析處理及三維可視化為一體的土壤孔隙結構分析與三維重建軟件系統，以期為當前有關領域的研究提供一個功能強大、簡單易用的操作平臺。

1 總體設計思路

由于實際土壤在孔隙形態學上的復雜性，為了便于三維網絡模型的建立與求解，孔隙體和孔喉通常被簡化成簡單的幾何形狀:孔隙體多以球形或立方體表示，孔喉則以圓柱形或其他具有代表性的非環形管道表示。雖然這樣形成的網絡模型是三維的，但是在視覺上與實際土壤內部結構還有較大差距。土壤孔隙結構三維重建系統是將一系列二維切片信息合成為三維立體影像的軟件，可以直接顯示出一個生動的虛擬實體，更直觀、更準確的觀察樣本土壤的內部結構，有助于了解大孔隙對土壤有效水力學性質的影響。

軟件開發環境包括:采用Delphi7.0作為軟件包的開發工具，OpenGL作為三維圖形顯示的開發工具。硬件環境:Intel奔騰4 CPU 2.4 GHz，512 M內存，GeForce4顯卡。土壤孔隙結構三維重建系統主要由主模塊、讀取圖片模塊、圖片觀察模塊、三維重構模塊以及數據存儲模塊組成。主模塊是系統的主體，負責完成系統的所有操作和各個模塊的激活;讀取圖片模塊將由X射線斷層攝影術(CT)、核磁共振(NMR)、掃描電鏡等無損技術獲得的二維連續切片灰度圖像數據讀取到系統中;在圖片觀察模塊中可以對這些二維切片進行縮放觀察、開運算、二值化處理等;系統在三維重構模塊中將處理后的二維連續切片按照設定的灰度值范圍以及閾值進行三維重建，重建出的三維立體影像可以實現旋轉、拖動、平移、光照控制等功能;整個過程中用到的數據都統一在數據儲存模塊中進行存放或讀取。整個系統的操作流程見圖1。

2 三維重建算法

三維重構是本軟件系統的核心功能，其算法設計將直接影響系統的運行效率、計算精度和可視化效果。用于三維重建的算法和技術有多種，如移動立方體算法(Marching Cubes，簡稱為MC算法)、剖分立方體法(Dividing Cubes)、射線投射法(Ray Casting)、剪切—曲變法(Shear-Wrap)等(張季和王宜杰，2006年)。其中，由Lorensen和Cline(1987年)提出的MC算法是發展較為成熟且得到最廣泛應用的一種表面重建算法。其基本思路是通過在三維網格頂點上的采樣構建出平面的組合來近似表達空間曲面，即將空間分割為一系列的立方體，然后將每一立方體在三維數據場中移動，并對立方體的頂點與等值曲面進行判斷比較，用合適數量的多邊形取代立方體，所有產生的多邊形將形成許多個等值表面，這些表面就是對輸入數據的近似表達(蔣先剛，2006年)。雖然MC算法的內存開銷較大，但其在編程上容易實現，且得到的三維圖像品質較高。本文中采用的MC算法的程序框圖如圖2所示。

算法首先建立單個立方體與等值曲面相交狀態的索引值表。對于一個已知的等值面的值，立方體每個定點處的值或大于或小于或等于該值。若頂點數據值大于或等于等值面的值，則該頂點位于等值面之上，賦值為1;否則該頂點位于等值面之下，賦值為0。等值面只與那些兩端點取不同值的邊相交。因為立方體有8個頂點，每個頂點有兩種可能的狀態值，故立方體與等值面相交的情況最多有256種(28)。將各頂點的狀態值拼成一個字節進行編碼即可構成一張索引值表(見圖3)，對應256種可能的相交情況。如索引值為二進制數00001000，指出等值面與立方體的e3，e4，e11三條邊相交。頂點v4高于等值面，其余頂點都低于等值面。

建立相交狀態的索引值表后對每個立方體內的等值面進行三角化。通過以下三個條件:1)繞軸旋轉一定角度;2)圖形相對軸面鏡面對稱;3)頂點的位置狀態取反，將256種位置狀態減少到不同構的15種(見圖4)。這樣就更加容易建立起三角面集合近似表達物體表面。

在三角化過程中物體表面與立方體邊的交點可以用線性插值的方法得到，精確的插值點取決于插值過程中用到的量值(對應二維圖像中的閾值)。例如，P1和P2是鄰接邊的兩個端點，V1和V2是這兩個點對應的采樣值(對應二維圖像中的灰度值)，物體表面與邊的交點P的坐標可由以下公式計算得到:

用一定規格的立方體分割空間，確定網格點的坐標及其對應的采樣值(對應二維圖像中的灰度值)。按順序選擇一個立方體，根據其8個頂點的不同標識獲得索引值并查表獲得當前立方體各相交邊所在，依次在各交邊上找到交點坐標，根據索引值和交點坐標形成當前立方體內三角形片的構成，選擇下一個立方體直到所有立方體都找到與之對應的三角形片構成。這樣構成的三角形片整體就是三維表面的近似表示。

從二維圖像提取出三維空間信息后，需要采用三維圖形開發工具將模型數據進行可視化。本文中的重建算法和三維可視化通過Borland Delphi來編程實現，其中使用了與Delphi7.0兼容的OpenGL 1.1。

實現的系統功能主要包括:1)繪制三維物體:通過讀取三角形數據鏈表繪制需要的三維物體;2)查看物體:建立視點，指定觀察角度、方向和觀察的范圍，實現對象的旋轉、縮放;3)指定顏色模式和光照:對繪制的三維物體加以著色、添加材質、光照，使其逼真于真實環境下的物體特性。

3 實例分析

供試土壤取自黃淮海沖積平原的中國科學院河南封丘生態實驗站內，土壤質地為壤質粘土(國際制)。按常規分析方法測定土壤容重為1.54 g/cm3，孔隙度為0.42。

CT掃描樣本用高20 cm、內徑10 cm的PVC管采集，采樣深度35 cm～70 cm。樣本運至中國科學院蘭州寒區旱區研究所進行測定，采用Siemens公司生產的SOMATOM-PLUS螺旋式CT機，綜合圖像處理系統為FLS-7。最大掃描分辨率為62.5μm，相鄰兩層切片間距為1mm。選取樣本中形狀保留完好的2 cm厚土柱進行掃描，共獲得連續切片數字圖像20張。

首先打開土壤三維重構系統并讀取切片圖像集，進入二維切片圖像操作界面，選取原狀土柱CT掃描獲得的8張連續圖片，在界面中可設置三維重構區域的灰度值邊界以及二值化閾值。閾值的選擇可根據實際土壤孔隙度按照單閾值法確定。

生成的土壤孔隙結構三維重建圖像如圖5所示，網格體數目設置為128×128×7，圖5中黑色實體部分為土壤孔隙，圖5中顯示出較大孔隙(＞750μm)的位置也與試驗土樣相一致。在此界面中可以對三維重構模型進行旋轉、平移或選取特定部分顯示等功能，更好的了解土壤的內部微觀結構。

4 結語

土壤孔隙結構三維重建系統可以根據CT，NMR等獲得的二維序列灰度圖像，通過計算機圖像處理與可視化技術，構建出三維土壤孔隙的虛擬實體，逼真、直觀、準確的顯示出土壤微觀結構，為土壤微觀模型的構建提供了一個可視化的軟件平臺。

在今后的研究中，需在以下幾個方面對該系統作進一步完善和改進:1)土壤孔隙結構三維重構系統僅生成了土壤孔隙的虛擬實體，缺乏對二維圖片以及三維重構模型分析計算的功能，不能得到實際土壤有效水力學性質方面的結果。2)對每個立方體進行三角形構造過程中，并沒有考慮其周圍立方體的情況，造成部分計算的冗余，應該進一步改進重建算法提高重建速度與精度。3)數據存儲功能尚不完善，一次重建后得到的圖形數據應當以標準格式儲存起來，以便下次直接調用。儲存格式可以與其他三維繪圖工具軟件(AutoCAD，3DS等)相一致，方便相互調用。

［1］ W.E.Lorensen，H.E.Cline.Marching Cubes：A High Resolution 3D Surface Construction Algorithm［J］.Computer Graphics，1987，21(4)：163-169.

［2］ 馮 杰，郝振純.CT掃描確定土壤大孔隙分布［J］.水科學進展，2002，13(5)：611-617.

［3］ 呂 菲，劉建立，何 娟.利用CT數字圖像和網絡模型預測近飽和土壤水力學性質［J］.農業工程學報，2008，24(5)：10-14.

［4］ 曾 箏，董芳華，陳 曉，等.利用MATLAB實現CT斷層圖像的三維重建［J］.CT理論與應用研究，2004，13(2)：24-29.

［5］ 蔣先剛.基于Delphi的數字圖像處理工程軟件設計［M］.北京：中國水利水電出版社，2006：291.

［6］ 張 季，王宜杰.醫學圖像三維重建方法的比較研究［J］.醫學信息，2006，19(5)：948-950.