基于形態學方法的牙頜模型自動分割算法

陳 虎 馬珂楠 胡畢富

1(北京大學口腔醫院口腔醫學數字化研究中心 北京 100081)

2(北京航空航天大學機械工程及自動化學院飛行器制造工程系 北京 100191)

0 引 言

在口腔醫學領域中,CAD/CAM技術的引入正逐漸改變傳統口腔治療方式,而牙科CAD系統作為口腔數字化醫療中的重要環節,對其設計策略和算法的研究都具有重要的意義[1]?？谇会t學包括口腔修復、正畸、種植以及頜面修復等主要方面,每種不同治療方式需要構建的CAD系統差異性很大,但很多情況下都需要從牙頜的數字化模型分離開牙齒和牙齦部分,實現對牙頜模型的分割,例如口腔正畸中的虛擬排牙等。受不同牙頜的形狀差異性,以及牙齒畸形、牙列缺損等因素影響,目前的牙科CAD軟件大多是以繁瑣、低效的手動交互方式完成分割的,對牙頜模型分割的自動化程度較低。

牙頜模型的分割屬于數字幾何處理領域中的三維網格分割問題,近年來三維網格分割已成為一個研究熱點[2],但已有的網格模型分割方法不能實現對任意形狀的模型的自動分割。目前針對具有復雜形狀的牙頜模型,已提出了較多的分割方法。文獻[3-4]從全景圖像中構建了兩幅距離圖像實現了牙齒的自動分割,避免了直接處理復雜的三維網格數據問題,但也使得這種分割方式不適合于處理三維牙頜模型的牙科CAD系統。目前比較常見的方法是通過分析牙頜模型的離散曲率實現對牙齒和牙齦的分離。Kronfeld等[5]利用牙齒和牙齦邊界附近頂點的曲率值較大的特點,構建了活動輪廓線,之后通過能量法局部優化邊界線實現牙頜分割。但該方法不易處理活動輪廓線出現不連續的牙頜模型。Kumar等[6]通過設置曲率閾值選取出高曲率值的頂點獲得齦緣線輪廓,在分離開牙齦和牙齒后,將齦緣線上的尖點作為相鄰牙的分割邊界點分離開相鄰牙齒。該方法也難以處理牙齒畸形等情形,比較依賴于標準理想的牙頜模型。為了提高分割的質量,在曲率獲取齦緣線的基礎上,利用形態學操作獲取更精細的牙齒邊界線。Zhao等[7]在加入較多交互的方式下,利用形態學操作實現了牙頜模型的分割。Yuan等[8]通過形態學操作分離開牙齒和牙齦曲面,并實現了對相鄰牙齒的分離以及對牙齒曲面的修補,該方法同樣需要較多的手動交互操作。Wu等[9]在形態學操作分割牙頜模型的基礎上,實現了少量交互下對牙齒牙齦以及相鄰牙齒的分割。除了利用牙頜網格模型的曲率特性實現分割外,還有其他反映網格曲面特征的信息可用于牙頜模型的分割。Zou等[10]通過計算網格上的調和場來反映網格的幾何細節,再手動選取具有約束作用的特征點,實現牙頜模型的分割,在之后的改進方法中[11],通過對牙頜模型的多層切分實現了自動選擇特征點,從而完成對牙頜模型的自動分割,但該方法在處理相鄰牙齒貼靠緊密的情況時難以將牙齒分割開。Hwang等[12]也基于調和場進行牙頜模型的分割,并采用狄利克雷邊界條件,凸分割等后處理方法,在具有噪聲和網格碎片的臨床實際的牙頜模型上取得了一定的效果。Fan等[13]采用離散Morse理論分析曲面結構,把分割問題轉化成從曲面上冗余特征線族中提取分割邊界線,再通過優化方法從底層分割中聚合產生子部分,獲得了較好的分割效果。Kim等[14]通過對牙頜模型進行截面投影,隨后應用區域增長算法實現牙齒的自動分割。He等[15]則基于曲率信息和局部分布密度信息初步劃分牙齒-牙齦邊界和牙齒間邊界,然后通過區域增長算法實現牙齒分割。Pupykina等[16]利用可視化工具包進行交互,基于用戶選擇的少量關鍵特征點,通過模糊綜合評判的方法確定牙齒的邊界。近年來也出現了許多基于機器學習的方法。Xu等[17]建立了牙齒分割的深度卷積神經網絡,將牙頜網格面片的特征作為輸入,經過兩層神經網絡,第一層用于區分牙齒和牙齦區域,第二層用于區分不同牙齒區域,最終通過模糊聚類對結果進行優化,取得了較好的分割效果。Lian等[18]建立了MeshSegNet,實現了牙頜模型多尺度的特征提取和學習,采用圖約束(Graph-Constrained)學習模塊進行局部特征提取,并通過密集融合(dense fusion)策略進行更高層次的特征學習,最后采用圖割算法對分類結果進行優化調整。Sun等[19]基于圖卷積神經網絡實現了牙齒和牙齦的同時分割。Zhang等[20]將牙頜三維模型映射到二維的調和參數空間,轉換為圖像信息,從而利用卷積神經網絡進行牙齒分割,隨后將二維圖像上的分割結果映射回三維實體模型中,并通過改進的模糊聚類和分割方法完成邊界的調整和優化。這些基于機器學習的方法雖然在一些數據集上取得了較高的準確率,但是受訓練數據規模、質量等影響較大,在處理數據質量較差的臨床實際模型時可能無法達到理想的效果。另一方面,由于訓練數據通常需要由專業人員進行大規模和高質量的標注,成本高昂,訓練和測試時對機器運算能力也有非常高的要求。

現有的牙科CAD系統在處理牙頜模型的分割時,都是以手動交互的方式分出牙齒再進行其他處理,頻繁地交互操作難以提高CAD系統的使用效果。但受到牙頜的形狀差異性的影響,以及牙齒排列畸形、牙列缺損等口腔病癥的出現,使得牙頜模型的自動分割一直沒有好的解決方案,還沒有通用的方法,已有方法都或多或少存在交互和準確性問題。

本文在現有形態學方法處理牙頜模型分割的基礎上[21],提出了一種自動閉合牙頜模型輪廓特征線的方法。牙頜三角網格模型在經過形態學處理得到單個頂點寬度的特征線后,有時存在不連續的情形。首先通過計算特征頂點的復雜度[8,21]判斷出特征線的端點,之后通過距離約束需要配對連接的端點以及無效的端點,再逐層分別搜索端點周圍的非特征頂點并標記,直到兩配對的端點都搜索到同一非特征頂點,以該頂點作為兩個配對端點的中點,根據最短路徑原則分別與兩端點相連,從而使得特征線封閉。該方法解決了原來需要手動處理特征線不封閉的問題,實現了對牙頜模型的自動分割。

1 網格曲面上的形態學操作

1.1 三角網格的描述

對于三角網格模型M通?？杀硎緸镸={V,T},其中V是網格頂點的集合,T是網格上三角面片面的集合。

(1)

三角網格拓撲信息的建立主要考慮網格M的頂點V、邊E和三角面片T的信息以及相互之間的遍歷和訪問。

牙頜模型曲面的特征提取主要可歸結為三角網格分割問題,特征線一般是分割的邊界線。實現對牙頜模型的分割關鍵在于齦緣線的提取,在牙頜模型中,齦緣線及相鄰牙齒之間的特征線附近的頂點平均曲率值比其他區域偏大,可以通過平均曲率值初步區分出齦緣線附近的特征區域。為了計算網格頂點上的離散曲率,對于每個頂點v∈V的主曲率可分別表示為κ1(v)和κ2(v),從而平均曲率表示為：

(2)

在牙頜模型的特征線附近都是凹下去的“谷線”區域,這些區域頂點的平均曲率值都是偏小的負值,通常直接設置曲率閾值ε來獲得初始的特征點集F：

F={v∈V|κH(v)≤ε}

(3)

根據曲率閾值初步篩選出的特征點集有較多的噪聲點,而且能夠表示齦緣線和相鄰牙齒分割線的特征點也存在缺陷,需要做后續處理,才能獲得適用于牙頜模型分割的只有一個頂點寬度的特征線。

1.2 形態學操作

形態學方法可以處理初始特征點集并獲得只有一個頂點寬度可用于模型分割的特征線。R?ssl等[21]最早將形態學操作應用于三角網格曲面的特征線提取,并將用于圖像處理的形態學操作轉化到處理三角網格曲面特征區域上,其中使用的幾種算子主要包括：膨脹、腐蝕、開閉操作、骨架化以及裁剪操作等。該方法主要在于利用頂點的鄰接關系來表示形態學中的結構元,為此定義了任意頂點vi的鄰接關系nhd,對于頂點vi的1-ring鄰接關系定義為：

nhd{i}={i}∪{j|?edge(vi,vj)}

(4)

鄰接關系可以逐次擴展為n-ring鄰接關系nhdn：

nhdn+1=nhd(nhdn{i})n>1

(5)

1.2.1膨脹和腐蝕操作

膨脹操作可以將比較接近但不連續的特征區域合并在一起,并且特征區域向周圍擴展,比原來的特征區域連續性更好。對特征區域F以nhdn為結構元的膨脹可表達為：

dilaten(F)={j|?i∈F：j∈nhdn{i}}

(6)

膨脹操作是把特征區域F中頂點v的n-ring鄰接點全部都加入到特征集F中。

腐蝕操作的作用和膨脹操作相反,可以將特征區域細化,達到簡化區域,使得特征結構更清晰的目的。對特征區域F以nhdn為結構元的腐蝕可表達為：

eroden(F)={j|nhdn{j}?F}

(7)

膨脹和腐蝕操作通常是組合來應用的,網格上的特征區域需要經過一系列的膨脹或腐蝕操作,兩者不同的組合會出現不同的操作效果,其中包括常用的開操作、閉操作等,這些操作會使得特征區域在保持原有基本形狀,獲得較清晰的骨架特征。

1.2.2骨架化和裁剪操作

骨架化操作可以把網格上的特征區域細化為單個頂點寬度的特征線,骨架化操作的實質是在判斷應去除特征區域中的哪些頂點。特征區域中的頂點根據其復雜度的不同可以分為三類[13]。對于頂點vi的復雜度CP(i)定義為：

(8)

式中：Cnei(i)表示頂點vi依次按順時針的1-ring鄰接頂點{v1,v2,…,vm},如果vk∈Cnei(i)且vk∈F,就令Cnei(i)k=1,否則Cnei(i)k=0。當CP(i)≥4時,頂點vi就是復雜點(Complex Vertex),構成的點集記作C;當CP(i)=0時是中心點(Center Vertex),構成的點集記作⊙;其他的特征點為外側點(Disk Vertex),其點集記作○,通常CP(i)=2。

骨架化操作就是去除外側點并保留復雜點的過程,可用數學式表示為：

(9)

骨架化操作產生短的無關的分支可以通過裁剪操作消除。令S?F表示獲得的骨架的特征點集合,裁剪操作為：

(10)

2 牙頜模型自動分割算法

2.1 算法概述

本文的牙頜模型分割算法中使用的牙頜模型是通過三維掃描儀根據牙頜的石膏模型獲取的點云數據生成的,牙頜數字化模型是一張三角網格曲面。本文的算法首先使用形態學方法生成齦緣線以及相鄰牙交線這兩種能用于牙齒分割的特征線,再通過端點的距離約束用最短路徑原則連接特征線斷開的地方,得到封閉特征線用于牙頜模型的分割[22]。如圖1所示,算法的主要步驟有：

圖1 算法概述

Step1初始化特征點集。計算每個頂點v的平均曲率值κH(v),對于給定的曲閾值ε,若滿足κH(v)≤ε,則把該點v記為特征點加入到集合F。

Step2刪除噪點。將集合F中特征點按是否相連分類得到m個特征區域{F1,F2,…,Fm},記num(Fi)為特征區域Fi中頂點個數,刪除num(Fi)/F<εd的特征區域,得到新的特征點集F,εd是設置的刪除比例。

Step3形態學處理。對特征點集依次做膨脹F=dilate(F),腐蝕F=erode(F),骨架化F=skeletonize(F),修剪F=prune(F)四種操作,得到單個頂點寬度的特征線F。

Step4識別特征線上的端點。對于頂點vi∈F,若vi的復雜度CP(vi)≤2,則把vi放入連接端點的集合ConnectP。

Step5對集合ConnectP中的點配對,計算連接中點MidP,根據最短路徑連接配對的端點,得到新的特征點集F。

Step6根據封閉的特征線F分割牙頜模型。

下面具體介紹Step5中的方法。

2.2 自動閉合特征線

特征線的端點集合ConnectP中的頂點既有用于連接間斷處的端點,也有一些未被完全刪除的分支上的末端點,會造成干擾,如圖1(h)所示。需要解決如何從ConnectP中找出需要相互連接的頂點,可以通過設置一個距離范圍δmax來約束端點的連接,避免錯誤連接的出現。

令δ(v)表示頂點u到v的最短距離：

(11)

式中：Pmin(u,v)表示頂點u到v的最短路徑。查找相互連接的端點的過程是不斷計算每個頂點的vi∈ConnectP與其n-ring鄰域內頂點u在δmax范圍內的最短距離δ(v),對n依次取n=1,2,…,n,每計算一次都要判斷頂點u是否在其他端點vj也被訪問,如果訪問過,則把頂點u存入連接中點的集合MidP中,并將頂點u與端點vi、vj根據最短路徑相連,這個過程如圖1(i)、(j)所示。對ConnectP中的端點同時計算,查找中間點,如果是特征線上的分支端,通常在δmax范圍內不會找到與之配對的端點,在特殊情況下如果存在,可以在形成封閉輪廓后,對網格模型做分割時剔除這些小的區域,不影響分割效果。

3 實驗結果

本研究中用于調參的牙頜模型如圖2所示,這些模型是由激光掃描儀對牙頜的石膏模型掃描獲取的曲面信息生成的三角網格曲面,圖2(a)、(b)、(c)中模型的頂點數量分別為1.53×105、1.32×105、1.56×105。程序中的平均曲率閾值取ε=0.66,在刪除噪聲點過程中的刪除比例取εd=0.01,在閉合特征線過程中的距離范圍取δmax=3.5 mm。牙齒與牙齦的分割以及相鄰牙齒之間的分割都比較準確完整,表明自動的牙頜模型分割算法具有良好的分割效果。

圖2 本文算法自動分割效果

本研究收集了另外5種常見模型(圖3)共49顆牙齒進行測試,本算法成功實現了46顆牙齒的自動分割,成功率達到93.88%。

圖3 本文算法在5種常見模型上的測試

4 結 語

本文改進了基于形態學方法的牙頜模型分割算法,解決了依靠交互方式選取牙頜模型特征線端點的問題,實現了牙頜模型的自動分割。實驗表明,該算法的分割效果較好,適用于不同牙頜模型的自動分割。

相關的后續工作包括：研究其他特殊的牙病患者模型對算法的影響,進一步提高算法的適用性;建立對牙頜模型分割效果的評估方法,判斷分割后的模型能否用于患者的治療。