優化的梯度最短路徑骨架提取算法

楊晨暉,劉 聰

(廈門大學信息科學與技術學院,福建 廈門 361005)

骨架的概念最早是由Blum在1967年提出的,當時他稱骨架為“中軸”(medial axis)[1],Blum還分別給出了2種骨架典型定義[1-2]．二維圖像的骨架提取方法大體可分為4大類:

1) 基于拓撲和幾何分析的方法;2) 細化的方法;3) 基于距離場的方法;4) 基于廣義勢場的方法．這些方法提取的結果依然不理想,主要表現在骨架的不穩定性和骨架的連通性方面,Bai等[3]引入邊界離散曲線演化(discrete curve evolution)在保證骨架連通性的前提下,剪除對視覺不重要的骨架分枝;Golland等[4]利用物體模型的拓撲先驗知識,提出用固定拓撲骨架(fixed topology skeleton)來減少骨架的冗余分枝,消除骨架的不穩定性．Choi等[5]引入尺度因子,提出距離場內基于連續準則的骨架提取算法,在一定程度上解決了骨架的連通性問題;丁頤等[6]提出一種基于種子生長的方法來連接斷裂的骨架．

本文引入劉俊濤等[7]提出的梯度最短路徑的算法,得到連通的、準確的物體模型骨架,并在此基礎上提出基于梯度的優化方法,大大地減少了程序運行的時間．

1 基于距離變換的梯度最短路徑算法

1.1 圖像梯度

本文利用歐氏距離變換對二值化的人體圖進行距離變換,然后利用一階微分形象化的定義,在行和列方向上各自構建一個簡單的3×3算子(如圖1所示),計算圖像的梯度．

圖1 自定義的行方向和列方向上的算子

Fig.1 Custom operators in the row and column

如圖2(c)可知,距離變換(DT)的梯度圖(|DT|)突出顯示了圖2(b)中灰度值變化比較顯著像素點(圖2(c)中前景點中黑色的部分)．實際這些像素點的集合已經組成了原人體圖的中脊線．由于梯度化的結果也突出模型輪廓邊緣的像素點,很難用閾值來消除非脊點的影響;即使消除這些影響,提取出的骨架也難以保證單像素性和連通性．因此單純依靠梯度圖來獲取模型的骨架不能滿足要求.

圖2 二值化人體、距離變換、梯度和局部極大值點示意圖Fig.2 The binary image of human body, the image of distance transformation, the grads image, the local maximum points image

1.2 關鍵點

設p是物體A內的任意一個像素點(或體素點),表示p點的距離變換值,S是點p鄰域點組成的集合．對任意的像素點(或體素點),都有DT(q)≤DT(p),則稱p為物體A的距離變換局部極大值點[7]．

本文是基于8鄰域搜索得到的距離變換局部極值點．

觀察圖2(b)發現,潛在骨架點的距離變換值應該比鄰域的距離值大．因此通過尋找距離圖中的局部極大值點,可以得到潛在的骨架點．如圖2(d)所示,粗點表示距離變換圖的局部極大值點,而且可以看到,局部極大值點全部落在梯度圖的脊線上．

由于最短路徑算法是按照梯度方向連接每個骨架點,局部極大值點必然會在最長的連接曲線上．利用這一性質,本文用一個極大值點來取代連通的所有局部極大值點,稱這一替代的極大值點為關鍵點．關鍵點是連通的局部極大值點中梯度最小的那一個像素點(體素點)[7]．

易知,關鍵點也是局部極大值點,所以關鍵點都坐落在梯度圖的脊線上,可以認為這些關鍵點是物體骨架的一個子集．之后的工作就是尋求一種方法把這些關鍵點按照梯度的方向連接起來,得到物體的曲線骨架．

1.3 基于梯度的最短路徑

本文采用的最短路徑算法是迪杰斯克拉(Dijkstra)算法．根據該算法的需要,必須根據梯度圖中的像素點來構建鄰接矩陣．

因此可以根據梯度權的定義,來構建圖像的一個無向圖和鄰接矩陣．值得本文關注的是如果圖像的兩個點不相鄰,他們邊的權值是無窮大的．這就要求必須保證新構建的無向圖所有頂點都處在同一個連通的分支,否則只要有孤立點出現,就會造成算法的死循環.

2 基于梯度的優化方法

基于梯度的最短路徑算法關鍵在于構建鄰接矩陣．鄰接矩陣的大小直接決定了算法時間的開銷．希望找到一種方法,在保證最終骨架效果的前提下,減少前景點的數量,降低鄰接矩陣的大小,從而提高算法運行效率．

由1.1中的圖像的梯度可知,對圖像求梯度能突出圖像的邊緣．針對距離變換圖,能突出潛在骨架點．根據梯度的這一特性,可以對距離圖求多個方向的梯度,找出所有潛在的骨架點,以達到減少搜索點的數量．這種方法特別在細節較多的情況時,更能起到好的效果．

根據實際效果和時間開銷的考慮,本文只從0°、45°、90°和135° 4個方向求距離圖的梯度．然后結合4個方向的梯度圖,取每個梯度圖對應位置灰度值最大的值作為當前像素點灰度值,形成新的梯度圖．

對比圖3和圖2(d)容易發現,圖3包括了所有潛在的骨架點和局部極大值點(連同無法用閾值消除的輪廓邊界點),而且中間的脊線更加的“粗”,這樣有利于最短路徑算法的搜索．觀察圖3,容易發現,在脊線的周圍有很多“毛刺”,這不利于鄰接矩陣的構建．通過研究每個像素的灰度值發現,周圍的“毛刺”的灰度值與中間脊線的灰度值相差較大．因此,本文利用閾值來消除大部分的“毛刺”．

本文圖像的梯度是利用卷積實現的,而卷積是基于浮點操作的．因此閾值的選取也應該根據浮點的結果來選?。驗楸疚闹皇侨コ糠窒袼攸c,閾值越小,原圖的拓撲結構保存得就越完整．通過實驗的統計,閾值的取值范圍為ε∈[0,8)為宜．圖4是閾值為5的結果圖．

雖然通過閾值的方法消除了大部分的“毛刺”點,但還有部分的離散點分布在脊線的位置．本文通過對閾值后的結果圖進行取輪廓操作,脊線即為輪廓最大的那些像素點所組成的閉合集合．

雖然中間的脊線能正確的提取,但所包含的像素點還是相對較多．優化的目的是在保證拓撲結構的前提下,盡量減少最短路徑搜索的像素點的數量．通過研究發現,拓撲細化的方法能夠滿足要求,圖5為用輪廓法消除“毛刺”的結果圖．簡單的細化過程類似“燒草”模型[1],雖然細化的操作不能保證得到的骨架接近中軸,但本文通過梯度化的處理等操作后,細化能保證得到的初步骨架與梯度化距離圖中的脊線位置重合．同時,為了保證關鍵點在細化過程中不丟失,本文在細化過程中,關鍵點被定義是端節點,不會被刪除．實驗證明,關鍵點都會落在細化的初步骨架上．

由于本文提出的優化方法在整個過程始終保持著連通性,能滿足構建最短路徑算法中鄰接矩陣的要求．

圖3 綜合4方向梯度圖形成的新的梯度圖

Fig.3 New grads image which is conbined the grads in four directions

圖4 閾值為5的結果圖

Fig.4 Image of the result when the threshold is 5

圖5 用輪廓法消除“毛刺”點的結果圖

Fig.5 The result of eliminating burr by the contour method

3 算法描述

為了獲得物體的曲線骨架,本文從二值圖像的物體模型出發,計算其二值圖的距離變換,然后梯度化距離圖,尋找潛在的骨架點;之后用基于梯度的優化方法構建最短路徑算法的鄰接矩陣;最后用基于梯度的最短路徑算法,連接所有的關鍵點,得到最后的骨架．算法步驟如下:

1) 計算輸入二值圖像I的距離變換DT,得到距離變換圖SDT;

2) 求SDT中所有距離變換局部極大值點所組成的集合C={c1,c2,…},并求SDT中距離變換值最大的點s∈C;

3) 應用圖3所定義的算子計算距離圖SDT的梯度DT,得到梯度圖S;

4) 根據集合C={c1,c2,…}和梯度圖S尋找關鍵點集合K={k1,k2,…}?C;

5) 應用基于梯度的優化方法,得到點集G;

6) 根據點集G構建鄰接矩陣,以s為起始點,應用基于梯度的最短路徑算法,找到s到K的所有路徑集合Ps→K={p1,p2,…},其中pi為s到關鍵點ki的最短路徑;

7) 以Ps→K為基礎構建一幅圖像,該圖像即為I的骨架圖．

4 實驗結果與分析

本文是在Pentium(R) Dual-Core,E6500 @2.93 G,內存2.0 G,操作系統為Window XP的PC機上,利用OpenCV2.3.1庫和Visual Studio 2008實現算法．選取MPEG-7 CE Shape-1 Part B、Kimia99和Kimia216[8]中部分有代表性的二值圖像,如棍形物體、多邊物體、非剛性物體等,應用本文的優化算法．同時對比分析了未優化前的劉俊濤等[7]人提出基于梯度的路徑算法(簡稱優化前算法)．

4.1 誤 差

圖6 優化前后梯度誤差和骨架點數目關系圖(ε=3)Fig.6 The connection between the number of skeleton point and gradient error before and after optimization

4.2 時 間

雖然優化前后提取的骨架效果一致,但在時間的消耗上相差很大．本文通過提取馬、手掌、錘子、人體、5角圖案和鉛筆6種典型物體的曲線骨架,對比說明了本文優化的方法在時間開銷上的優勢．表1列出了6種物體運用優化前和優化后算法的時間開銷．其中前4種中閾值ε=4,后2種ε=6,保證優化前后算法提取的結果是一致的．

表1 優化前和優化后時間開銷對比Tab.1 Time consuming compared before and after optimization

顯然,本文提出的優化方法大大地減少了基于梯度的最短路徑算法的輪廓內前景點數量,減少了程序運行的時間．

4.3 閾 值

優化方法中,閾值ε的大小決定著算法時間的開銷．本文在不影響提取效果的基礎上,分析了4種典型剛性與非剛性、規則與不規則物體的時間開銷與ε的關系．圖7表明,當閾值ε達到某個臨界值時,繼續增大其值,將不會明顯地減少程序的運行時間．

圖8 加入邊界噪聲前后效果圖對比Fig.8 The result of the image before and after joing the counter nosiy

圖7 程序運行時間與閾值ε關系圖Fig.7 The connection between the running time of program and the threshold

4.4 邊界噪聲

本文利用距離變換和圖像的梯度來獲取物體的骨架,由于其固有的特性,在一定程度上克服了邊界噪聲對骨架的影響．圖8顯示了不同類型物體加入噪聲后,本文的優化算法提取的骨架結果圖．可以發現,加入邊界噪聲前后的骨架圖在拓撲結構上沒有發現大的變化,僅在形狀上有細微的差別,如線狀的平滑度．

4.5 算法復雜度分析

本文假設圖像大小為m×n,鄰域為8鄰域．下面主要從時間方面來分析本文算法的復雜度．

根據第3節所述的過程,算法時間的消耗主要表現在:1) 距離變換的計算;2) 圖像梯度的計算;3) 圖像輪廓的尋找;4) 細化;5) 最短路徑的尋找．

距離變換是基于歐氏的變換,算法的時間復雜度控制在O(m×n)．若采用基于邊界跟蹤的算法,雖然不能改變總體的復雜度,但大大地減少了計算的時間開銷．

圖像的梯度是利用圖像的卷積計算得到的結果,給定一個k×k的卷積核(即本文所提到的算子),卷積計算的時間開銷O((k×k)×(m×n)),其中k一般取3或5等奇數常數．因此,計算圖像梯度的時間復雜度為O(m×n)．

本文利用openCV提供的findContours函數來獲取二值圖像的輪廓．而findContours是根據suzuki85實現的,時間復雜度不會超過O(m×n)．

細化是掃描圖像的每一個像素,通過8鄰域的關系判定每個像素是否為簡單點．本文只迭代的做一次細化,因此細化的時間的復雜度為O(1×8×(m×n))=O(m×n)．

搜索兩個關鍵點之間的基于梯度的最短距離,引用Dijkstra算法,其時間復雜度為O(N2)(其中N為第3節優化后的輪廓內前景點數量)．

5 結 論

本文詳細介紹了基于梯度最短路徑算法相關定義、原理、算法流程和復雜度分析．針對算法中鄰接矩陣過大,導致算法時間開銷過大的問題,提出了基于梯度的優化方法．并且采用閾值去噪和輪廓細化的方法進一步提高算法的效率,選取了MPEG-7和Kimia中部分有代表性的圖片,分析對比了優化前后算法的時間開銷和骨架效果．結果證明:在很好的得到單像素連通骨架的同時,本文的算法效率得到很大的提高.最后還研究了優化方法中閾值ε的大小對算法時間開銷的影響．

[1] Blum H.A transformation for extracting new descriptors of shape[M].Boston:MIT Press,1967：362-380.

[2] Blum H.Biological shape and visual science(part I)[J].Journal of Theoretical Biology,1973,38(2):205-287.

[3] Bai X,Latecki L J,Liu W.Skeleton pruning by contour partitioning with discrete curve evolution[J].IEEE transactions on Pattern Analysis and Machine Intelligence,2007,29(3):449-462.

[4] Golland P,Eric W,Grimson L.Fixed topology skeletons[C]∥Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition.Hilton Head Island,South Carolina,USA:IEEE Computer Society Press,2000:10-17.

[5] Choi W P,Lam K M,Siu W C.Extraction of the euclidean skeleton based on a connectivity criterion [J].Pattern Recognition,2003,36(3):721-729.

[6] 丁頤,劉文予,鄭宇化.基于距離變換的多尺度連通骨架算法[J].紅外與毫米波學報,2004,24(4):281-285.

[7] 劉俊濤,劉文予,吳彩華,等.一種提取物體線形骨架的新方法[J].計算機學報,2008,34(6):617-622.

[8] MPEG MPEG7 CE Shape-1 Part B[EB/OL].2004.http:∥www.imageprocessingplace.com/root_files_V3/image_databases.htm.