(浙江工業(yè)大學 計算機科學與技術學院，浙江 杭州 310023)

3D模型在計算機輔助設計、計算機圖形學和科學可視化等領域中扮演重要角色。作為重要研究對象，3D模型被廣泛應用，包括三維模型特征提取[1]、孔洞修復[2]、模型重建[3]和遙感等。傳統(tǒng)計算機圖形學中的三維形狀表示方法一般表示的是模型表面信息，并不描述模型內部結構和模型拓撲信息[4-5]。而模型的三維表示是“昂貴的”，許多應用需要模型的可替代的緊湊表示[6]。有這樣的一種“線狀”或“棒狀”表示方法，稱之為“骨架表示”或“曲線骨架”[7]。Blum[8]以中軸的形式建立了骨架化的基礎，該骨架由具有較低維度的對稱平面/軸組成。理論上講，Blum的骨架或中軸是使用草火傳播過程定義的。假設物體是干草材質，在物體邊界用火點燃，火焰以均勻的速度在物體內部傳播。骨架被定義為兩條獨立火線相遇點的集合。

模型骨架提取方法可以分為以下幾類：1) 拓撲細化法(Topology thinning)，又稱為模擬燒草模型法。此類算法從邊界開始，反復迭代地逐層剝離離散后的模型形成骨架。絕大部分細化方法依靠簡單點作用于離散空間(體素)，去除之后不會改變形狀的拓撲結構的點(體素)稱之為簡單點[9]。為提供拓撲保護，文獻[10]提出了兩種分別用于提取“表面骨架”和“曲線骨架”的8 次子迭代算法，其中每次子迭代只能去除某種特定類型的邊界點。文獻[11]提出了一種拓撲保護的三維骨架化方法，該方法主要計算體素標記為D6距離模型的面骨架或曲線骨架?？紤]笛卡爾超立方網格，文獻[12]提供了一個數(shù)學框架，其中給出了迭代細化過程去除點的顯式布爾條件。2) 基于距離場的方法。首先，需要估計一個距離圖，然后視中軸為距離圖的脊，通過脊跟蹤算法在此距離圖中抽取得到中軸。物體內部每點的距離值DT被定義為此點離模型邊界點的最短距離[13]。由于骨架對模型整體而言是處于中心位置應具有較大的DT值，因此DT能作為骨架提取的一種重要參考指標。Bitter等[14]通過梯度搜尋方法檢測非均勻梯度的相鄰部分，并標記這些點為候選體素。Bouix等[15]用散度計算選取優(yōu)先函數(shù)并由用戶自定義閾值來進行體素的簡單點移除。3) 基于拓撲與幾何分析的方法，通過構造模型的Voronoi圖和Reeb圖來得到骨架。Voronoi圖將空間劃分成子區(qū)域，Voronoi圖的內部邊和面可用于提取形狀的中軸面(骨架)的近似值。筆者提出了一種基于距離場的骨架提取方法，通過尋找距離域的“脊線” 從而找到骨架。筆者主要的貢獻在于提出了最大內切球算法來對候選體素進行約束，能很好地去除冗余體素點。筆者提供了幾個模型的實驗結果并將之與經典骨架提取方法進行比較，實驗結果表明：筆者方法能得到一個質量較好的三維模型骨架。

1 算法流程

筆者提出了基于距離場的骨架提取方法，通過包圍盒技術確定包裹模型的最小立方體空間。根據所求得的包圍盒將模型體素化，從而將模型劃分為一個個立方塊空間。再計算每塊體素的距離值，找出單位空間內距離值最大的體素作為骨架點的潛在體素。比較過后，所得的候選體素并非全部屬于骨架點所在體素，提出最大內切球逼近算法來對體素做出約束，去除冗余的體素，從而找出最終的骨架點。最后，根據所求得的體素，按照一定規(guī)則有序連接形成骨架。具體算法流程如圖1所示。

圖1 算法流程Fig.1 Pipeline of skeleton extraction method

如圖1所示，筆者方法主要步驟如下：

1) 模型內部體素化。此過程主要包括包圍盒計算和模型外部體素去除。由于筆者方法的核心思想是基于距離場來找到處于中心的那條“脊線”。因此，模型外部體素可用信息不大。根據需求，對于包圍盒的選取筆者選擇的是AABB盒，主要是其操作簡單且存儲空間小。在判斷體素是否外部體素時，充分利用了模型邊界點的法向量指向模型外部這一信息，根據向量積的正負，成功地去除掉外部體素。

2) 單位空間內最大距離值體素的選取及最大內切球約束。模型內部體素化后，從邊界到中心，體素的距離值是遞增的。因此，位于模型中心的骨架的距離值一定是最大的，所要找的就是這類體素。因此，事先設定一單位空間，然后找出該空間內的最大距離值的體素，比較完所有體素后得到初始篩選的候選體素。然而，這些體素并非全部都屬于骨架。因此，利用最大內切球再對這些候選體素進行約束，對于每個候選體素找出其最大內切球，將此最大內切球所包裹的所有體素再進行距離值比較，選擇其中距離值最大的體素與最大內切球球心距離值(半徑)作比較。若該體素距離值大于該最大內切球半徑，則將此體素作為該局部的骨架點。若該體素距離值小于該最大內切球半徑，則將此最大內切球的球心作為該局部的體素點。

3) 骨架點連接形成骨架。此部分主要是考慮骨架點連接的規(guī)則，根據已求得的骨架點形成骨架。必須保證骨架點連接不能重復、亂序和錯漏。

2 最大內切球逼近

三維模型骨架就是模型M中所有最大內切球球心的集合[6]。筆者提出的最大內切球逼近算法主要思想：從模型內部任一點出發(fā)，計算其與表面最近點之間的矢量并向模型內部延伸來調整點的位置，從而將初始點位置從低距離值往高距離值的場內移動，經多次調整內部點位置便可以找到三維模型的最大內切球及其球心。

設三維網格模型M(x,y,z)={vk,k=1,2,…,N}，其中(x,y,z)是頂點v的坐標，N是頂點數(shù)目，則距離值函數(shù)為

(1)

式中v為距離點p最近的模型表面點。圖2為Human模型局部放大距離場灰度編碼圖，顏色越淺代表距離值越大(模型中心為灰白色，模型邊界為灰黑色)。

圖2 距離場灰度編碼Fig.2 Distance field color coding

定義函數(shù)Extend_Len(p′,p″,l)來調整點的位置。假設給定點p′和點p″，沿矢量p′p″延長距離l得到點p?，則

p?(x?,y?,z?)=Extend_Len(p′,p″,l)

(2)

結合公式(1，2)可得

(3)

最大內切球逼近算法流程如下：

Step1最近點計算。針對模型內部任一點p1，遍歷模型表面點，找到與其距離最近的模型表面點v1，并計算p1的距離值以及矢量v1p1。

Step2反向延伸調整位置，找到距離新點的模型最近點并計算兩點間的矢量。沿矢量v1p1方向向p1的距離值等值線內部延長距離l2，其中l(wèi)2=β×L(0<β<1.0)，L為模型網格平均邊長。由式(2)確定模型內新點p2=Extend_Len(v1,p1,l2)，并找到距離p2最近的模型邊界點v2，計算p2的距離值以及矢量v2p2。

Step3重復執(zhí)行Step 2，計算pk的距離值以及矢量vkpk，直到‖vk+1pk+1‖<‖vkpk‖為止，即當pk+1的距離值小于pk的距離值，此時pk的距離值達到最大。則以pk為球心，‖vkpk‖為半徑可以得到模型最大內切球，其球心為模型骨架點。

圖3為horse模型的最大內切球逼近效果圖。圖3(a)為模型內部初始點位置，圖3(b)為其逼近的最大內切球效果圖，其中β取0.1，耗時52 ms。

圖3 最大內切球逼近Fig.3 Maximum inscribed ball approximation

3 三維模型體素化

體素化是將模型的幾何形式表示轉換成一種最接近該物體的離散表示形式，產生體數(shù)據集。表示模型的空間體素跟表示圖像的二維像素比較相似，只不過從二維的面擴展到三維的立方體單元，而且基于體素的三維模型技術被廣泛應用。筆者算法主要思路就是操作離散化模型，找到屬于骨架的體素點。

3.1 三維模型體素化

基于幾種包圍盒技術理論，在選擇包圍體類型時充分考慮各種包圍體類型的優(yōu)缺點以及筆者方法的需求性，最終選擇軸對齊包圍盒(AABB)，主要是因為其構造比較簡單，存儲空間小。

設三維網格模型M(x,y,z)={vk,k=1,2,…,N}，其中(x,y,z)是頂點v坐標，N是頂點數(shù)目。遍歷這N個頂點，找到Xmin,Xmax,Ymin,Ymax,Zmin,Zmax這6 個值，從而確定包圍盒的8 個頂點。依次為：P1(Xmin,Ymax,Zmin)，P2(Xmin,Ymin,Zmin)，P3(Xmin,Ymin,Zmax)，P4(Xmin,Ymax,Zmax)，P5(Xmax,Ymax,Zmin)，P6(Xmax,Ymin,Zmin)，P7(Xmax,Ymin,Zmax)，P8(Xmax,Ymax,Zmax)。圖4為Human模型AABB包圍盒，其中圖4(a)為Human模型初始圖，圖4(b)是其包圍盒生成圖。

圖4 包圍盒生成Fig.4 Bounding box generation

3.2 模型內部體素提取

離散化操作中，選取精度為d(2 個單位的平均網格邊長度)來將包圍盒所包裹的空間劃分成一個個體素塊。根據得到的包圍盒信息，全部體素數(shù)量應為Num=Nx×Ny×Nz。其中，Nx=(Xmax-Xmin)/d，Ny=(Ymax-Ymin)/d，Nz=(Zmax-Zmin)/d。

這些體素主要分為3 部分：模型外部體素、模型邊界體素和模型內部體素。實驗中，實驗外部體素通過判定條件剔除，模型邊界體素也不可能是骨架所在的體素，真正需要的是模型內部體素。圖5(a)提供了人體模型體素化后(已剔除外部體素)的效果圖，為更清楚觀察體素化后的內部情況，圖5(b)提供了身軀的局部放大圖。其中，精度為2 個單位的平均網格邊長度。為清楚觀察每塊體素，放置半徑為0.5 個平均網格邊長度單位的小球在每塊體素的中心。外部體素判定是通過向量乘積判定的，如圖6所示。對于任一點p(x,y)，計算其與模型表面最近的點，可以利用模型表面點法向量始終指向模型外部這一信息來判定。圖6(a)中，模型內部點p(x,y)，p1(x1,y1)是與之最近的模型表面點，由法向量指向外部可以得到：pp1·n1>0，其中n1是p1(x1,y1)的法向量。圖6(b)中，對于模型外部點p(x,y)，pp2·n2<0，n2是p2(x2,y2)的法向量。

圖5 模型體素化Fig.5 Model voxelization

圖6 點的判定Fig.6 Point determination

4 三維模型骨架提取

三維模型骨架提取主要包括對體素化后的模型進行距離值局部最大值求解以及最大內切球約束。距離值局部最大值求解是用一個小的立方塊空間去遍歷整個模型內部，不能存在錯漏，找出該空間內的距離值最大的體素。該篩選后的體素包含處于骨架上的體素點以及需要去除的冗余體素點。最大內切球約束是對篩選后的體素進行最大內切球逼近，找出最大內切球所包圍的距離值最大的體素，再將該最大值與該最大內切球球心距離值(半徑)作比較。若該體素的距離值較大，則將該體素點作為骨架點。若最大內切球半徑大，則將最大內切球球心作為骨架點。經過1～2 次的最大內切球約束，能較好地去出冗余體素。最后再按照一定規(guī)則將所求得的骨架點進行有序連接形成骨架。

4.1 距離值局部最大求解

對于體素化模型，只有模型中心位置的體素屬于骨架點，本部分主要考慮如何剔除大量不屬于骨架位置的體素。實驗過程中，考慮使用邊長為l的立方塊遍歷模型內部空間，找出每個空間內距離值最大的體素塊，從而完成模型體素快的初步剔除。為提升效率，對每一體素塊的數(shù)據結構定義為

Voxel=(x,y,z,distance,flag)

式中：x,y,z為體素的位置信息；distance為體素的距離值信息；flag為體素的標志位，表示體素是否被比較，體素標志位設置初始值flag=0。為將筆者方法闡述清楚，先在二維情況下對筆者方法作個引入，如圖7所示。

圖7 骨架提取原理圖Fig.7 Skeleton extraction schematic

如圖7(a)所示，對一已體素化的形狀，圖7(a)中用點來表示體素，該點位于體素中心。圖7(a)中所有點初始設置為灰色，對應著三維情況下標志位初始值flag=0。圖7(a)中方框對應著三維情況下的立方塊，所要比較的就是方框內所包裹的所有體素的距離值，距離值最大的體素標記為深黑色，相應地三維空間內的flag=1，其余的比較過的體素標記為淡灰色，對應的flag=-1。因此，對于所有體素點進行局部距離值最大求解時，會對標志位flag先進行判定，如此可大大提高效率。圖7(b)為局部距離值最大求解的結果，所選中的初始骨架點已標記為深黑色。

圖8為人體模型距離值局部最大求解得到的候選體素,采用的立方塊邊長為l=3.25×d。從圖8可以看出：距離值局部最大求解已經從人體模型中剔除了大量的體素塊。但骨架點主要集中于人體模型中心位置，因此需要再次對這些體素塊進行選擇。主要考慮使用最大內切球約束的方法去挑選屬于骨架的體素。

圖8 候選體素Fig.8 Candidate voxel

4.2 最大內切球約束

此部分是對距離值局部最大求解得到的候選體素篩選做詳細說明。如圖7(b)所示，并不是所有深黑色點均為骨架點，需要進行選擇，筆者主要是通過最大內切球的方法進行選擇。如圖7(c)所示，對深黑色體素點進行最大內切圓逼近，找出該最大內切圓內所包裹的體素點，并將其中距離值最大的體素點的距離值與最大內切圓圓心距離值(半徑)進行比較。若體素的距離值相比最大內切圓半徑大，則選擇該體素為骨架點。反之，則令該最大內切圓圓心作為骨架點。對應三維中，是對所有標志位flag=1的體素進行最大內切球逼近，對于該球所包裹的體素點進行比較。選中為骨架點的體素點，令其標志位flag=2。對未選中的體素點，令其標志位flag=-1。最終標志位flag=2的體素即為所選中的骨架點，如圖7(d)所示。

圖9為對圖8中候選體素進行最大內切球約束效果圖。圖9(a)是經過1 次最大內切球約束的結果，但仍然有冗余點的存在，需要進一步的選擇，如人體模型的胸腔部分。圖9(b)為2 次約束的效果圖。因為最大內切球逼近能保證所找的骨架點較好地居于骨架正中心，因此圖9(b)中所選擇的骨架點有較好的中心性。

圖9 最大內切球約束Fig.9 Maximum inscribed ball constraint

4.3 骨架生成

本部分主要介紹利用圖9所選的骨架點進行骨架生成。如上文所述，所得到的骨架點標志位flag=2。基于距離優(yōu)先的思想，骨架生成過程如下：

1) 設定2 存儲序列A1,A2，用A1存儲所有flag=2的骨架點。

2) 從A1中任意選擇并剔除一個骨架點，令其flag=3然后置于A2。

3) 遍歷A2，計算每個骨架點與A1的所有骨架點之間的距離長度并進行排序，選擇最小的那對骨架點。將該對骨架點連接，并將此flag=2的骨架點從A1剔除且令flag=3，然后置于A2。

4) 多次重復步驟3)，直至A1為空。

在實驗過程中對步驟3)的連接進行了限制，要求選擇的一對骨架點連線不能超出模型外部(通過判斷連線四分點是否都處于模型內部)，否則舍棄此種連接，選擇其他相離最近的一對骨架點。

需要說明的是，經過距離值局部最大求解和最大內切球約束后，不會出現(xiàn)骨架點距離較近的情況。對于分支較好的模型，當分支骨架點flag=3，也能較好地參與步驟3)的選擇連接，因為步驟4)的停止條件是A1為空，即不存在flag=2的骨架點，因此不會出現(xiàn)漏連的情況。對于拓撲結構復雜的模型，處理的效果還不足，將來需要對連接規(guī)則進一步的優(yōu)化。

5 實驗結果與分析

本節(jié)算法已經在Microsoft visual studio 2005開發(fā)環(huán)境下實現(xiàn)，利用OSG進行模型的渲染。程序的運行環(huán)境為IntelCore2 Duo E7500，2.93 GHz CPU，2 GB內存。為證明筆者方法的實用性，本部分展示了幾個模型的骨架提取效果圖以及與其他骨架提取方法比較。

圖10是從圖9骨架點所生成的人體模型骨架，頂點數(shù)12 500，面片數(shù)目25 000。耗時7.6 s。從圖10可以看出：筆者方法所提取的骨架基本上能獲取到人體模型的拓撲結構，且能較好地居于模型正中心，因為所有被選取的骨架點都經過了最大內切球的約束，能很好地保證骨架的中心性。

圖10 人體模型骨架Fig.10 Human body model skeleton

圖11(a)是Torso模型初始圖，頂點數(shù)目99 447，面片198 384。圖11(b)是文獻[16]提取效果圖，圖11(c)是筆者方法提取效果圖。從圖11可以看出：筆者方法提取的結果圖與文獻[16]提取的結果圖在效果上相近，都能大致反映原始模型的拓撲結構。文獻[16]提取的結果圖在骨架平滑性上處理的不夠，但相比筆者方法，其提取的骨架更居于模型中心。提取Torso模型耗時44 s，主要是因為在體素化模型并去除模型外部體素過程中耗費了大量時間。距離值局部最大求解過程中，筆者方法采用的立方塊邊長l=2.5×d個平均網格邊長度，體素化精度為d=8 個平均網格邊長度。

圖11 骨架對比Fig.11 Skeleton contrast

圖12為筆者方法對horse模型提取效果圖，頂點數(shù)目19 851，面片數(shù)目39 698。圖12(a)為初始模型，圖12(b，c)分別為Wu等[17]和Lien等[18]提取的骨架效果圖，圖12(d)為筆者方法效果圖。采用的立方塊邊長l=4×d個平均網格邊長度，體素化精度為d=2 個平均網格邊長度，消耗時間為7.8 s。從圖12可以看出：提取的效果與文獻[17]和文獻[18]相近，筆者方法提取的骨架基本上能獲取到模型的拓撲結構。但筆者方法提取的骨架更接近模型中心位置，比如horse模型腿部位置，而文獻[17]和文獻[18]做的相對差一些。

圖12 骨架對比Fig.12 Skeleton contrast

圖13為Liu等[16]、Tagliasacchi等[19]和筆者算法對Homer模型的骨架提取效果圖。Homer模型頂點數(shù)目5 103，面片數(shù)目10 202，筆者采用的立方塊邊長l=3.0×d個平均網格邊長度，體素化精度為d=1.6 個平均網格邊長度，耗時2.8 s。從圖13可以看出：這3 種方法提取的模型骨架都能較好地體現(xiàn)輸入模型的拓撲結構，相比Tagliasacchi等[19]的效果，筆者方法提取的骨架更為完整，圖13(c)模型腳部骨架有一定的缺失。Liu等[16]提取的骨架雖然結構完整，但其平滑性不如筆者方法。但Liu等[16]、Tagliasacchi等[19]的魯棒性較強，能處理許多復雜模型，筆者在這一點上做的還不夠，將來需要改進。

圖13 骨架對比Fig.13 Skeleton contrast

6 結 論

提出了一種基于距離場的全自動骨架生成方法。首先計算輸入的三維模型的包圍盒，并根據包圍盒尺寸將該空間離散化，分散成一個個立方塊空間；之后，通過計算剔除模型外部體素，保留模型內部體素并計算體素的距離信息，從而完成模型的體素化；然后，對存儲著距離信息的模型內部體素進行距離值局部最大求解，完成對體素的初步選?。唤又?，利用最大內切球逼近算法對保留下來的候選體素進行最大內切球約束得到骨架點；最后，根據選擇的骨架點，按照一定的規(guī)則進行連接形成骨架。筆者方法最大的創(chuàng)新之處在于用最大內切球逼近算法約束候選體素，能有效地去除冗余體素并對候選體素的位置進行一定的矯正。由于在模型體素化過程中要耗費大量時間，對于頂點數(shù)較多的模型消耗的時間會相對長些，如圖11所示的Torso模型。未來方向將考慮如何提升體素化效率以及處理更多復雜模型。