基于LUT的快速3D氣道樹骨架線提取

劉明威， 顧力栩

(上海交通大學 生物醫學工程學院，上海 200240)

基于LUT的快速3D氣道樹骨架線提取

劉明威， 顧力栩

(上海交通大學 生物醫學工程學院，上海 200240)

為提高肺部支氣管骨架線的提取效率，提出并使用一種基于look-up-table(LUT) 的腐蝕細化算法. 分析建立腐蝕模型，并根據該模型優化建立了LUT，以該LUT為依據通過索引查找對原始數據進行快速腐蝕細化，對得到的腐蝕結果進行剪枝處理以得到最終的骨架線. 實驗結果表明：提取過程中LUT的應用從根本上降低了腐蝕細化中判斷的復雜度，將復雜的簡單點判斷問題轉化為LUT中的查詢問題，從而極大地優化了腐蝕細化中關鍵的腐蝕過程. 相比傳統方法，基于LUT的腐蝕細化算法顯著提高了骨架線的提取速度，較傳統細化法提速近22.95倍.

計算機輔助診斷； 氣管樹； 虛擬支氣管鏡；骨架線；查找表

肺癌是威脅人類生命健康的惡性腫瘤疾病，中國目前總體肺癌發病率增長了19.48%，這很可能是加劇的環境污染問題如大規模霧霾引起的[1]. 早期診斷能大幅度提高病患存活率，但目前早期診斷率僅15%[2]. 常見肺癌早期診斷手段包括CT成像與光學內窺鏡. CT成像擁有密度分辨率高、操作簡單等優勢，但傳統CT只能提供單層的斷面信息，不直觀，且在診斷過程中病人需承受大量放射輻射；光學內窺鏡可在零輻射下觀察到管腔內部真實紋理色彩，但操作上要比CT成像復雜得多.

隨著計算機輔助診斷技術的發展，虛擬支氣管鏡VB[3](virtual bronchoscopy)技術在肺癌診斷領域正逐漸得以普遍應用. 該技術利用從CT數據中重建的氣管模型，對肺部疾病診斷有重大意義. 臨床表明VB不僅可以顯示有無氣道異常如支氣管狹窄等，還可以配合導航系統引導經支氣管鏡肺活檢，明顯提高支氣管鏡在肺部外周病灶診斷效率. 相對于傳統內窺鏡而言，VB引導下的導航技術具有以下優點：1)定位準確，診斷效能高；2)可快速引導支氣管鏡，縮短檢查時間，減少病患痛苦；3)可減少支氣管鏡檢查次數；4)可提供自由的虛擬觀測角度，降低支氣管異常的辨別難度； 5)可部分代替X線引導的支氣管鏡活檢，減少輻射；6)可避免盲目穿刺，降低并發癥發生率. 大量文獻證實虛擬內窺鏡能在不失診斷準確率的基礎上，保障診斷的微創性[4-5].

在虛擬內窺鏡實現過程中，快速的骨架線提取一直是其中的難點和要點. 而在VB中，一個良好提取的骨架線不僅可以提供直觀的支氣管拓撲結構，還可以為管腔內瀏覽起到引導作用. 與人體其他管腔如結腸、血管相比，支氣管結構更為復雜[6]，體現在：1)復雜的拓撲結構，擁有多級分支，多通路；2)各級分支粗細不一；3)個體差異性，不同人支氣管雖在級數上有一定相似性，但大小和分支位置卻存在差異. 目前三維骨架提取方式主要有距離變換法、中軸變換法和腐蝕細化算法. 然而由于支氣管的復雜性，這些方式在支氣管骨架線的提取上都略顯性能不足. 隨著當今4DCT[7]數據的逐步普及，傳統方法的運算速度更顯得不盡人意.

本文提出一種基于腐蝕細化的改進三維骨架線提取方法，以克服傳統腐蝕細化方法中復雜度高、耗時長的問題. 本方法充分利用了當今計算機多核以及高內存的優勢，在傳統三維腐蝕過程中引入LUT. 相關實驗證實，本方法能夠顯著提高骨架線的提取速度.

1 骨架線基本定義與現有方法綜述

1.1 骨架線基本定義

骨架(Skeleton)是指與原物體具有一致連通性和拓撲結構的細曲線的一種理想表示，就是位于物體內部，且能體現其拓撲特征的簡化圖形[8]. 在三維情況下提取的骨架線需要滿足如下要求[9]：1) 連通性，骨架線應與目標物體擁有相同的連通性；2) 居中性 ，骨架線盡量保持在管腔的中央；3) 單體素性，組成骨架線的體素要求為單元寬度. 對于簡單的柱狀組織，中心線只有唯一的一條，而對于支氣管樹或血管樹這種復雜的管狀結構，中心線會有多條，也被稱作為骨架線.

1.2 現有方法綜述

就目前而言，二維的中心線提取算法已經相對成熟. 相對而言，三維中心線的提取方式只有距離變換法、中軸變換法和腐蝕細化算法.

距離變換法[10]是指根據不同目的，將物體的體素點進行距離標識. 目前主流的方法是基于雙距離場的方法[11]：首先對物體進行邊界距離變換，得到體素到達邊界的邊界距離場DFB (distance from boundary)；其次指定源點，進行源點距離變換，得到源距離場DFS (distance from source)；最后通過邊界距離場中局部最大值的篩選，配合源距離場進行骨架點的連接. 中軸變換法[12]在二維上可以用最大內切圓盤的相關術語定義，而在三維上則表現為最大內切球中心的集合. 其中圓盤中心到達邊界的距離稱為圓盤法線，一個目標內點的最大內切圓在邊界上至少有兩個切點，而每個對稱點有兩個或多個圓盤法線. 腐蝕細化算法[13]的思想則更為樸素，根據物體歐拉特性[14]以及連通性的不變性，給出在刪除后不影響原物體拓撲結構的簡單點，并以此為依據均勻、對稱地從不同方向對目標物體邊緣上的簡單點(Simple Voxel, SV)[15]進行剝離，直至僅剩下單像素寬的體素集合，去除該集合中的偽骨架點后得到最終的骨架線.

以上方法各自均有不同的局限性和一些共同不足. 距離變換法中，在構建雙距離場時不僅因需要對各體素點以迭代的方式進行距離標記而導致高時耗，還因需要與目標物體同樣大小的內存空間儲存距離場而帶來高內存消耗. 同時，依據不同的歐氏距離定義建立的距離場會產生不同的候選中心點，直接影響骨架線結果的一致性[16]，且候選中心點通常是不連續的，連接時會造成骨架線在曲率較大位置的離心偏移，難以保證居中性. 中軸變換法中，由于中軸變換是理想的數學模型描述，主要基于數學的幾何推理證明，從而導致使用機器語言進行實現的難度大、算法復雜度高，常見的模型有地表火模型[17]，距離曲面脊線模型[18]. 以上兩種方法存在共同的缺陷就是易受邊緣噪聲干擾，且難以保證骨架的連續性. 腐蝕細化算法的不足在于需要不斷地對邊緣的簡單點進行判斷，現有的判斷模型復雜度均不夠理想，導致當目標物體體素點過多時，帶來巨大的時耗.

2 改進的腐蝕細化方法

考慮到基于距離場的算法對于邊緣噪聲的敏感性和內存的大量消耗，以及傳統腐蝕細化算法中對邊緣點反復判斷的高復雜度，本文提出一種基于LUT的快速三維腐蝕細化算法. 該方法保留了傳統腐蝕細化算法的穩定性，以及結果中對目標物體拓撲結構和連通性的完好保留性，并通過使用LUT的思想簡化細化過程，避免重復判斷，顯著提高了細化速度. 該方法首先從原始數據中用區域增長的方式分割出支氣管二值圖像，建立腐蝕細化模型并建立LUT，然后查找LUT進行腐蝕細化，最后對細化的結果進行剪枝處理.

2.1 支氣管二值數據的獲取

采用有監控的動態閾值區域生長配合形態學閉操作來獲取組成支氣管，流程及效果如圖1所示. 首先以氣管壁與空氣CT值的差異為依據，通過應用DMQ (dynamic marking QUEUE)的數據結構進行增長，提高了提取的速度. 在頂層CT中選擇氣管入口作為增長原點壓入DMQ，目標點離開DMQ時對其6鄰域進行增長判斷，若符合則壓入DMQ. 在判斷目標點是否歸入增長區域時，采用了18鄰域平均值比較的方式，最大可能增加分割目標點數. 對每一輪新增的點與上一輪增長結果進行比較，以進行溢出監控. 該方法得到一個26連通的初步結果，為平滑邊緣、填充管腔的空隙，在其基礎上使用公式IMG ·M = (IMG ? M) ⊕M進行形態學閉操作，用同一模版M分別對IMG進行膨脹腐蝕操作得到最終的二值分割結果.

2.2 細化方法

細化過程具體分為4步：1)腐蝕模型的建立； 2)基于腐蝕模型建立LUT(提高速度的關鍵步驟)；3)基于LUT進行快速細化； 4)對腐蝕結果進行剪枝處理.

圖1 支氣管分割流程

2.2.1 腐蝕模型的建立

在三維空間V3中，設中心點為c(x,y,z),對于任意點p(i,j,k)有如下的鄰域定義：

Nb6(p) = {p(i,j,k)∈V3| |i-x| + |j-y|+ |k-z| = 1},

Nb26(p) = {p(i,j,k)∈V3| |i-x|<2 ∩|j-y|<2∩|k-z| <2},

Nb18(p) = {p(i,j,k)∈V3| |i-x| + |j-y|+ |k-z| < 3}∩ Nb26(p).

將屬于目標物體的體素定義為前景點Sfg，其余定義為背景點Sbg. 在前景點集合中存在子集Ssub-fg，若其中所有元素均符合以下四種判斷準則[19]，則稱之為簡單點SV(simplevoxel).

a) ?Sfgin Nb26and Sum(Sfgin Nb26)>1；

b) ?Sbgin Nb6；

c) ?Sfgin Nb26and is 26Connected(Sfg)；

d) ?Sbgin Nb6and is 6Connected By Nb18(Sbg).

a、b條件可以直接通過相應鄰域中的前景點、背景點分布獲得，而c、d條件的判斷才是關鍵.

對于條件c，提出一種26鄰域區域增長法，如圖2所示. 根據條件c的描述易證明：若Nb26中的前景點可以彼此26連通，那么在該Nb26中，一定有且只有一個26連通的前景點區域. 因此，取26鄰域中的任意Sfg，在該Nb26中進行26鄰域的前景點區域增長，得到增長結果集合S′. 將其與26鄰域中的前景點構成的集合S進行比較，若兩者點數相同，則符合條件c.

在條件d的判斷上，并不能完全與條件c等同，因為在Nb18中的Sbg中，6連通的數量并不唯一. 充分考慮到這一點，利用Nb18中Sbg的6鄰域局部區域增長法，使用列表套列表的結構來儲存連通的6鄰域數組. 在遍歷結束后，將Nb6中的背景點集合S與連通的6鄰域列表中的Si依次進行比較，得到是否可以刪除.

圖2 條件c判斷

2.2.2 鄰域分布分析與LUT的建立

分析腐蝕模型可以看出，判斷一個前景點是否為簡單點是由a、b、c、d這4個條件的滿足情況決定的. 進一步分析，條件a為非孤立點的判斷，取決于Nb26；條件b為邊界點的判斷，取決于Nb6；條件c為Nb26的背景點的26連通情況，取決于Nb26；條件d為Nb26前景點的18連通情況，同樣取決于Nb26. 因此可得出結論，對于任意體素點p，可根據其Nb26(p)分布情況VC(voxelcombination)得到其簡單性.

在二值圖像前提下，任意體素點只有0或1兩種取值. 如圖3所示，若對26鄰域指定排列順序，則可等同于一個26bits的二進制數，進而任意一種VC將對應唯一的26bits二進制數，而相對的任意一個26bits二進制數也能唯一映射為一種VC. 基于以上分析，可聯想到建立LUT輔助腐蝕細化. 在腐蝕細化過程中，將目標點的26鄰域VC轉化為26bits的索引，在LUT中進行查詢判斷. 文獻[20]顯示LUT在二維的細化中得到了很好的應用，值得一提的是，由于三維計算量巨大，空間分布復雜性以及受拘束于過去計算機的內存，三維LUT鮮有人嘗試.

為避免建立LUT過程中的重復判斷，對體素的鄰域的空間分布狀態進行特征分析. 對于任意立方體，給出如下對于相似體素分布SVC(similarvoxelcombination)的定義：

1) 通過任意旋轉變換可以完全重合；

2) 通過左右手坐標系對換可以完全重合，即

圖3 26鄰域轉換

(1)

(2)

(3)

條件1)可解釋為依次以立方體每一面朝向的坐標系方向為觀察方向，而后以該方向為軸向，使用如式(1)～(3)中對應方向的旋轉矩陣對立方體進行圍繞該軸向的三維旋轉，為保證立方體空間位置不變，這里旋轉的角度應該分別取0°、90°、180°、270°，這樣通過旋轉可得4*6 共24種SVC；條件2)可解釋為條件1)基礎上的鏡像變換，因此同樣存在24種SVC. 因此可得出結論，對于各向異性的立方體，存在24+24共48種SVC，應用在二值化的Nb26情況下，粗略估算在全空間2^26中存在2^26/48約為1398101種SVC.

(4)

而在二值化的Nb26情況下，并非所有分布都是各向異性，存在著大量的自反對稱，經過式(4)的計算，其中i表示Nb26中前景點的個數，R(i) 表示Nb26中重復的組合個數，由于二值互補性，因此只需計算13次(0～12個前景點的情況)，外加13個前景點的情況. 針對i取值不同時需要具體列出對稱情況進行計算，但由于其復雜的旋轉對稱性，需要使用計算機輔助邏輯進行計算，最終得到不同前景點個數下SVC如表1所示，SVC一共有1426 144種，并在實際運算遍歷中加以檢驗.

表1 26鄰域中前景體素點與SVC關系表

基于以上分析，本文提出多線程并行配合雙布爾數組的算法建立LUT，數組大小為2^26，下標代表該索引值所唯一對應的VC. 其中一個數組為LUT，是后續腐蝕細化中所依賴的參照，其取值表示該VC下的中心點是否為簡單點；另一個數組為MT(marking table)，其取值用來表示該VC是否已被計算. 將所有2^26個索引等分為多個索引段，平均分配給多個子線程，配合多核CPU進行計算，以達到最高的計算效率.

所有子線程將共享LUT和MT，每個子線程將同時并行遍歷被分配的索引段. 在具體計算索引時，每一次簡單點判斷的成本很高，但應用雙數組算法可有效避免重復計算. 對于一個索引值，首先在MT中查找判斷該索引值是否已被計算，若尚未計算，則將該索引值轉換為其對應的VC，并根據前文的判斷條件和方法得出其是否為簡單點；而后通過正方體的旋轉對稱性，獲取該VC的所有其他47個SVC，并轉換為索引序列IA(index array)；最后在LUT中將IA中的索引值賦予與該VC相同的簡單性，同時在MT中將IA中的這些索引值標記為已判斷. 雙數組的優勢也在此，通過標記MT避免了遍歷到SVC所對應索引值時潛在的重復簡單性判斷，從而大大減少了判斷的次數，整體加速了LUT的建立過程. 為更為直觀，下面給出該流程框架的偽代碼.

函數 1: CreateLUT(創建LUT 的主函數)

輸入: 需要的線程數目 threadNumber.

輸出: 計算好的lutArray.

for eachiinthreadNumberdo

subThreadsi=start_new_thread(MultiThreadCalculate(lutArray,mtArray,rangei) )

endfor

wait_sub_threads(subThreads)

write_to_bin_file(lutArray)

函數2:MultiThreadCalCulate(用于計算指定區段的子函數)

輸入: 公共使用的lutArray，用于標記是否計算過的mtArray，需要遍歷的區段range.

輸出: 計算指定區段range后的lutArray， 標記指定區段range后的mtArray

foreachiinRange

ifnotmtArrayithen

isSimple=isSimplePoint(i)

similarIndexArray=getSimilarIndex(i)

foreachjinsimilarIndexArray

lutArrayj=isSimple

mtArrayj=true

endfor

endif

endfor

2.2.3 基于LUT進行迭代腐蝕細化

在建立好LUT的基礎上，傳統的腐蝕細化過程便得以優化，其中對于SV的判斷可以直接依據LUT的結果，而不是重新對前文4個定義條件的判斷. 在進行細化前先將圖像中的前景點存入臨時數組PointListfg，以避免重復遍歷，提高效率. 而后進行如下操作：

1)在PointListfg中根據Nb26體素分布情況尋找邊緣點，得到PointListborder.

2)遍歷PointListborder中的點，將Nb26體素分布轉換為26bits索引值，根據LUT的判斷值進行簡單點的判斷刪除.

3)重復步驟1)、2)直至PointListborder中再也無法找到簡單點為止.

2.2.4 后續剪枝處理

實驗證實，通過上述腐蝕細化算法所得結果會在一些較粗的主氣管上出現一些額外的毛刺，這無疑會對后續的內窺鏡路徑的規劃帶來不必要的分枝，需要進行剪枝處理.

根據毛刺分支長度短的特點，設計了一種尺寸自適應的剪枝算法. 記當前的前景點集合為S，首先根據骨架線的體素數設置適當的剪枝閾值Threshold，設閾值百分比為th%，這里th取經驗值1.7，因此經計算得到 Threshold=S×th%. 接著遍歷細化后的骨架線，得到所有端點PointListend，以PointListend中的點為起點進行遍歷，直至分支點，若其長度小于Threshold則將該支刪除.

3 實驗結果

3.1 實驗環境與實驗數據

實驗數據采用從上海肺科醫院獲得的5組肺CT圖像序列，每組序列均包含有最少265層到最多468層CT斷層圖像，每層圖像的分辨率均是512*512 且灰度值在-2048到4196區間內. 用于實驗的PCCPU為Intel(R)Core(TM)i7-3610QM，RAM為6GB，OS為Win7 64-bit. 算法的實現基于VTKv6.1.0開源庫，配合QT在VS2008IDE下進行開發和實驗. 3.2 LUT的建立

為獲取最適合本實驗環境的線程數量，以得到最快的LUT建立速度，以線程數作為輸入參數進行實驗，得到如圖4所示的結果.

圖4 LUT建立線程時間函數

可見當線程數少于5個時，增加線程可顯著減少用時，5個線程的處理過程表現最為優秀，當線程多于5個時，用時無明顯減少.

3.3 骨架線提取結果與比較

為比較基于LUT的方法與傳統腐蝕細化算法的效率，采用如下的對比試驗，使用相同的數據集進行重復試驗100次，表2列出了實驗結果.

在腐蝕細化階段，基于LUT的平均時間為482.74 ms，而基于傳統腐蝕細化算法的平均時間為11 077.7 ms，可見本文方法提高了22.95倍速度.

3.4 骨架線的提取效果

圖5為本文算法提取骨架線的結果與基于距離變換算法[21]的結果，其中圖5(a)、(b)、(c)三組結果中，上方為本文算法的結果，下方為距離變換算法的結果. 通過對比可以看出，基于距離變換的方法會在末端發生斷裂，不連續，在同一支氣管段上可能產生多分支，不方便路徑規劃，且操作中需要手工選定端點，既不方便也無法保證連續性；而本文方法在支氣管骨架線的提取上既可以保證骨架線在細末端的連續性，不發生斷裂，也無需手動，具備可重復性. 圖5中的三組數據顯示，該算法很好地保持了支氣管拓撲結果，同時也保障骨架線的中心性、連續性、一致性和可重復性，為內窺鏡的漫游提供可能.

表2 實驗結果對比表

(a) 512*512*265 (b) 512*512*429 (c) 512*512*438

4 實驗討論

4.1 基于LUT算法的優越性

本文同時還與基于距離場[21]算法得出的速度結果進行比較，由于使用不同的數據集，因此定義指標ts(thining speed)來表示算法的效率作為比較依據，ts表示毫秒內處理的立方點數, 單位為103voxels·ms-1. 經計算，文獻[21]的速度為220.16 ts，而本文的速度為45 253.39 ts，提高了205倍. 結合本文與以上兩種優化過的傳統骨架線提取方式[21-22]的比較，可以充分看出基于LUT的算法能在保證骨架線拓撲結構完整的基礎上，顯著提高運行速度.

總的來說，基于LUT的算法在本質上優化了傳統細化過程中的腐蝕細化環節，從而在整體上提高了骨架線提取的速度. 不僅如此，LUT的另一大優勢是將LUT的建立與LUT的參照分離，這樣在簡單點的定義發生變化時，可以直接更新LUT的數據，而不需要更改后續整個程序.

4.2 參數與復雜度分析

基于LUT的腐蝕細化算法中，主要的輸入參數體現于兩點：一是LUT建立過程中線程數的選擇，二是應用LUT進行腐蝕細化時，不同輸入物體的體素數量.

在LUT建立過程中，線程數是重要的輸入參數，并影響最終LUT建成的時耗. 實驗結果表明5個線程為最優，分析其原因為：當線程池大小為CPU數+1時，可達到最優利用率，即當執行密集型任務時既使偶爾因頁缺失故障或者其他原因導致暫停時，此額外線程也可確保CPU的時鐘周期不被浪費. 實驗PC 為四核CPU，因此5個線程可達到最優化. 算法的空間復雜度為2*2^26*1 byte = 128 MB，用于緩存雙數組中的LUT和MT，在LUT建成結束后，得到的結果為一個64MB的二進制文件，可以直接作為查找腐蝕的依據. 在時間上，對腐蝕細化的判斷條件c的判斷由于應用了26鄰域內區域增長法，只用一重循環，是一個線性過程，將該步驟的時間復雜度降低為O(n)；而對于條件d的判斷，在區域增長上的平均時間復雜度依舊為O(n)，而在列表比較上的平均時間復雜度為O(n2)，兩步之間為串行，因而整體的復雜度依舊為O(n2)，但較傳統方法相比依舊有所提升.

在腐蝕細化的過程中，由于LUT已經確立，因此主要的輸入參數為物體的體素數. 對3.3中的實驗結果進行線性擬合，得到t=0.007 2v-152.91，擬合優度R2為0.988 5，其中t表示腐蝕細化的時耗，單位為ms，v表示物體的體素數. 因此可得，腐蝕細化的時耗正比于輸入物體的體素數. 算法的空間復雜度取決于物體的體素數，會隨著CT采集時分辨率的不同而有所變化，另外需考慮LUT占用的64MB；在時間復雜度上，查詢LUT接近于O(1)，對體素點進行擦除的復雜度也為O(1)，因此整體復雜度仍為O(1).

4.3 基于LUT算法的不足

與基于其他理論的骨架線提取方式如距離變換法相比，基于LUT的腐蝕細化法不足在于最終得到的骨架線是組成骨架線的離散點，缺少點與點之間拓撲關系的描述，為后續骨架線拓撲結構分析帶來了一定的困難. 同時剪枝過程也有一定的局限性，無法精確區分細末較短的支氣管與偽分支. 而在LUT的使用中依舊可以作適當優化，如在查詢表的尺寸上，若是對內存有苛刻的要求，可將LUT壓縮為8MB，以bit的形式來儲存每一個索引的索引值. 由于對bit的檢索需要額外的比較，其效率可能會有一定下降，這是一個棘手的問題.

受實驗數據源限制，未進行大規模驗證. 同時受實驗室硬件條件限制，缺少在其他配置相當的計算機上進行重復驗證，也沒有獲取該算法在更差運行環境中的表現. 更加充足的數據源和豐富的實驗平臺將能夠對該算法有一個更全面的評估，這也是未來工作的目標.

5 結 語

基于LUT的腐蝕細化算法在傳統的三維腐蝕細化算法中引入LUT的思想，通過鄰域區域增長法優化了簡單點腐蝕模型的判斷，并通過二值化26鄰域旋轉對稱性的分析進行快速的LUT建立，在腐蝕細化中使用LUT取代傳統算法中對簡單點的復雜重復判斷過程，從而顯著提高了骨架線的提取速度. 實驗結果表明：本文方法在速度上是傳統細化法方法[22]的22.95倍，傳統距離場方法[21]的近220倍，可以完好地保留肺部氣管樹的拓撲結構，具備易操作性、連續性和可重復性，為內窺鏡瀏覽導航提供可能.

[1] CHEN W, ZHENG R, ZENG H, et al. Annual report on status of cancer in China, 2011[J]. Chinese Journal of Cancer Research, 2015, 27(1): 2.

[2] 韋春暉. 肺癌早期診斷進展[J]. 臨床肺科雜志, 2010, 15(8):1136-1138.

WEI Chunhui. Progress in early diagnosis of lung cancer[J]. Journal of Clinical Pulmonary Medicine,2010,15(8):1136-1138.

[3] JUNG S Y, PAE S Y, CHUNG S M, et al. Three-dimensional CT with virtual bronchoscopy: a useful modality for bronchial foreign bodies in pediatric patients[J]. European Archives of Oto-Rhino-Laryngology, 2012, 269(1): 223-228.

[4] WERNER H, SANTOS J R L D, FONTES R, et al. Virtual bronchoscopy for evaluating cervical tumors of the fetus[J]. Ultrasound in Obstetrics & Gynecology, 2013, 41(1): 90-94.

[5] DE WEVER W, VANDECAVEYE V, LANCIOTTI S, et al. Multidetector CT-generated virtual bronchoscopy: an illustrated review of the potential clinical indications[J]. European Respiratory Journal, 2004, 23(5): 776-782.

[6] DENIZ A, HOFFMAN E A, GEOFFREY M L, et al. Segmentation and analysis of the human airway tree from three-dimensional X-ray CT images.[J]. IEEE Transactions on Medical Imaging, 2003, 22(8):940-50.

[7] COLE A J, O'Hare J M, McMahon S J, et al. Investigating the potential impact of four-dimensional computed tomography (4DCT) on toxicity, outcomes and dose escalation for radical lung cancer radiotherapy[J]. Clinical Oncology, 2014, 26(3): 142-150.

[8] 陳剛, 呂煊, 王志成, 等. 肺 CT 圖像的血管骨架化方法[J]. 計算機科學, 2013, 40(5): 274-278.

CHEN Gang, LV Xuan, WANG Zhicheng, et al. Vessel Skeletonization Method for Lung CT Images[J].Computer Science, 2013, 40(5): 274-278.

[9] ZHOU Y, TOGA A W. Efficient skeletonization of volumetric objects[J]. IEEE Transactions on Visualization & Computer Graphics, 1999, 5(3):196-209.

[10]WAN M, DACHILLE F, KAUFMAN A. Distance-field based skeletons for virtual navigation[C]//Proceedings of the Conference on Visualization'01. IEEE Computer Society, 2001: 239-246.

[11]WANG S, WU J, WEI M, et al. Robust curve skeleton extraction for vascular structures[J]. Graphical Models, 2012, 74(4): 109-120.

[12]潘鵬, 賀三維, 吳艷蘭, 等. 曲邊多邊形中軸提取的新方法[J]. 測繪學報, 2012, 41(2): 278-283.

PAN Peng, HE Sanwei, WU Yanlan, et al. A new method for extracting curved-polygon medial axis[J]. Acta Geodaetica et Cartographica Sinica, 2012, 41(2): 278-283.

[13]CHEN W, SUI L, XU Z, et al. Improved Zhang-Suen thinning algorithm in binary line drawing applications[C]// International Conference on Systems and Informatic. Yantai:IEEE, 2012: 1947-1950.

[14]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.

[15]MICHEL C, GILLES B. New Characterizations of Simple Points in 2D, 3D, and 4D Discrete Spaces[J]. IEEE Transactions on Pattern Analysis & Machine Intelligence, 2008, 31(4):637-648.

[16]GENG H, YANG J, TAN W, et al. Fast 3D skeleton extraction of airways and applications to virtual bronchoscopy[C]// The 26th Chinese Control and Decision Conference (2014 CCDC). Changsha:IEEE, 2014: 3879-3884.

[17]LIU L, CHAMBERS E W, LETSCHER D, et al. Extended grassfire transform on medial axes of 2D shapes[J]. Computer-Aided Design, 2011, 43(11): 1496-1505.

[18]GOLDBERGERA L, AMARAL L A N, GLASS L, et al. Physiobank, physiotoolkit, and physionet components of a new research resource for complex physiologic signals[J]. Circulation, 2000, 101(23): e215-e220.

[19]PALAGYI K. A 3D fully parallel surface-thinning algorithm[J]. Theoretical Computer Science, 2008, 406(1-2):119-135.

[20]楊威, 郭科, 魏義坤. 一種有效的基于八鄰域查表的指紋圖像細化算法[J]. 四川理工學院學報(自然科學版), 2008, 21(2): 61-63.

YANG Wei, GUO Ke, WEI Yikun. An effective thinning algorithm for finger print based on 8-neighbors LUT[J]. Journal of Sichuan University of Science and Engineering(Natural Science Edition), 2008, 21(2): 61-63.

[21]JIANG G, GU L. An automatic and fast centerline extraction algorithm for virtual colonoscopy[C]// 27th Annual International Conference of the Engineering in Medicine and Biology Society. Shanghai:IEEE, 2005: 5149-5152.

[22]陳磊, 王勝軍, 鄭全錄, 等. 基于 CT 圖像的三維拓撲細化算法及其在心臟 CAD 中的應用[J]. 計算機應用, 2007, 27(B06): 406-410.

CHEN Lei, WANG Shengjun, ZHENG Quanlu, et al. A 3D topological thinning algorithm based on CT image and its application in cardiac CAD, Computer Application, 2007, 27(B06): 406-410.

(編輯 王小唯， 苗秀芝)

A fast LUT-based airway skeleton extraction algorithm for virtual bronchoscopy

LIU Mingwei, GU Lixu

(School of Biomedical Engineering, Shanghai Jiao Tong University, Shanghai 200240, China)

In order to improve the efficiency of lung airway skeleton extraction, this paper introduces a new look-up-table (LUT) based thinning algorithm. This new approach consists of three major steps: the analysis and creation of the thinning model and the establishment of LUT based on above thinning model result. Thinning process by index-searching uses the LUT. Branch cutting based on the result of thinning obtains the final result. The usage of LUT index-searching during the thinning step transforms simple point judgments into LUT index search and has significantly improved the performance of the whole algorithm. Experimental results demonstrate that the new skeleton algorithm is 22.95 times faster than the existing thinning algorithm.

computer-aided diagnosis; lung airway; virtual bronchoscopy; skeleton; LUT

10.11918/j.issn.0367-6234.201601026

2016-01-01

國家自然科學基金(61271318)

劉明威(1990—)，男，碩士; 顧力栩(1964—)，男，教授，博士生導師

顧力翔，gulixu@sjtu.edu.cn.

文獻標志碼: A 文章編號: 0367-6234(2017)05-0134-07