基於邊界距離場的血管中心路徑提取方法
2023-04-28 03:53:21
專利名稱：：基於邊界距離場的血管中心路徑提取方法
技術領域：
：本發明屬於生物醫學
技術領域：
，特別涉及血管的中心路徑提取方法，可用於各種管腔組織的中心路徑提取，並應用到虛擬內窺鏡中。
背景技術：
：醫學圖像中的血管精確提取和描述具有廣泛的應用，如梗塞或者狹窄化的血管疾病診斷及不同時段的病人圖像配準等；血管的提取方法是計算機血管輔助診療系統的關鍵。雖然目前存在很多種血管提取方法，但因成像模式、應用領域、人機互動以及其他因素的影響而各不相同，目前尚沒有一種通用的提取方法能夠適用於所有模式的醫學圖像。可視化技術的發展使得人們可以對所獲得的三維數據通過各種重建方法進行顯示以獲取有用的信息。對於血管疾病診斷而言，利用可視化技術對提取後的血管進行三維顯示，不僅能夠對血管外表面進行觀察，還可以對血管內部進行觀察，即血管虛擬內窺鏡。要實現血管虛擬內窺鏡，必須提取血管的中心路徑，即自動漫遊路徑。中心路徑的提取方法一般分為三類人工指定，拓撲細化和距離變換。人工指定方法要求用戶在體數據切片上人工指定路徑中心點，所有的點連起來組成中心路徑。這與人工漫遊一樣，雖然比較費時，但精確度較高，能夠更準確的接近感興趣區域。拓撲細化又稱為剝洋蔥法，是一個較為傳統的方法，主要思想是通過反覆剝除物體最外層體素，且每一步只刪除不會影響模型拓撲結構的點，直到只剩下單體素寬度的骨架為止，將所獲得的骨架作為中心路徑。拓撲細化是一個迭代過程，它能夠保持模型的整體連接性，但計算量十分龐大，而且不能直接連通到任意給定的兩點，即對於不在骨架上的點，需要事先連接到骨架上。距離變換主要基於這樣的事實物體的中軸是物體內部到邊界的最大距離點的集合。這種方法通常是先對二值體數據進行距離變換，然後將變換後的距離數據看作一個有向或無向的加權圖，最後採用Dijkstra最短路徑生成算法得到中心路徑。距離變換法的速度較快，但在保持模型的拓撲結構方面不如拓撲細化法。用於血管漫遊的中心路徑提取多採用距離變換法，原因主要有兩點一是醫學模型通常數據量較大，採用拓撲細化法耗時太多；二是提取中心路徑是為了在漫遊時指導視點的移動，其主要目的在於更好的觀察血管內部結構。目前存在很多基於距離變換的中心路徑提取算法。德國的Marburg大學數學與計算機學院的JanEgger等在CT血管數據中利用距離變換法提取血管的骨架作為虛擬漫遊路徑，取得了較好的效果，見文獻"EggerJ,Mostarki6Z，GroUkopfS,etal.Afastvesselcenterlineextractionalgorithmforcathetersimulation.TwentiethIEEEInternationalSymposiumonComputer國BasedMedicalSystems,Maribor,Slovenia.2007."。Bitter等和Sato等用邊界距離場和Dijkstra最短路徑提取算法提取結腸等數據的中心路徑，見文獻"SatoM,BitterI,BenderMA,etal.TEASAR:Tree-structureextractionalgorithmforaccurateandrobustskeletons.Proc.ofIEEEPacificGraphics2000,HongKong,China.2000:281-289."禾卩"BitterI,SatoM,BenderM,etal.CEASAR:ASmooth,Accurateandrobustcenterlineextractionalgorithm.Proc.ofIEEEVisualization2000，SaltLakeCity,USA.2000:45-52."。上述現有的方法主要存在以下三個缺點(1)邊界距離場變換是針對二值數據而言，因此需要預先提取血管組織。目前的血管提取方法，如匹配濾波器，主動輪廓模型以及水平集等方法雖然能夠很好的將血管提取出來，但速度慢是一個很大的應用瓶頸；(2)以二值數據為目標的邊界距離場計算，對於不同位置的血管組織，其初始的邊界距離值是相同的，不利於快速提取到位於血管中心的路徑點；(3)根據邊界距離場提取的中心路徑無法保證中心路徑上的點處於血管管腔橫截面的重心處，會導致血管的虛擬漫遊過程與人類的觀察習慣產生偏差。
發明內容本發明的目的在於克服上述已有技術的不足，提供一種基於邊界距離場的血管中心路徑高精度快速提取方法，以實現對血管的快速提取和血管中心路徑的高精度提取，使血管虛擬漫遊過程更符合人類的觀察習慣。實現本發明目的的技術方案，包括如下步驟A.對輸入的CTA血管造影體數據進行八叉樹分解，得到每一分塊的信息；B.根據每個分塊的最大值和最小值與雙閾值之間的關係，遍歷存在血管體素的葉子節點提取血管組織；C.根據提取的血管組織中所有點的坐標，計算提取後的血管組織邊界距離場；D.在邊界距離場的基礎上建立最大代價生成樹；E.提取最大代價生成樹的主幹作為血管組織的初始中心路徑；F.利用重心法修正提取的血管組織初始中心路徑，使路徑上的點都處於血管管腔橫截面的重心位置，得到一條符合人類觀察習慣的最終血管中心路徑。上述步驟B所述的根據每個分塊的最大值和最小值與雙閾值之間的關係，遍歷存在血管體素的葉子節點提取血管組織，按如下步驟進行Bl.從深度為0的八叉樹節點開始，判斷其最大值和最小值與雙閾值之間的關係，如果有至少一個閾值在最大值和最小值範圍內，則該節點內存在血管體素，否則忽略該節點；B2.判斷存在血管體素的節點的深度值與設定的八叉樹深度值之間的關係，如果節點的深度值與所設定的八叉樹的深度值相等，則該節點為葉子節點，並遍歷該葉子節點提取血管組織，否則判斷該節點的子節點的最大值和最小值與雙閾值之間的關係；B3.如此循環，直到將所有包含血管體素的葉子節點遍歷完畢為止，完成血管組織的提取。上述步驟C所述的計算提取後的血管組織邊界距離場，按如下步驟計算Cl.根據提取的血管組織中所有點的坐標，利用原CTA血管造影體數據中血管體素的梯度值倒數和拉普拉斯變換值之和作為血管體素的邊界距離初始值其中，p為提取的血管組織中任意一點的坐標，/e為原CTA血管數據中點/7的梯度值，厶為原CTA血管數據中點p的拉普拉斯變換值；C2.根據邊界距離初始值用如下公式計算點/7的邊界距離值/^。,(/0:(P)=min(min(/咖(A)+10),min(/加"(P2)+14),min(/麵(/3)+17》其中，義^(A)為點P的最近鄰鄰域點A的邊界距離初始值，,^(A)為點P的面對角鄰域點p2的邊界距離初始值，,M(p3)為點p的體對角鄰域點A的邊界距離初始值；C3.利用八叉樹遍歷所有包含血管體素的葉子節點，完成提取後的血管組織邊界距離場計算。上述步驟F所述的利用重心法對提取的血管組織初始中心路徑進行修正，按如下步驟進行Fl.對於初始中心路徑上的每一點，取與該點法向量垂直，並且與該點相距為10個體素的點組成的圓形二維圖像Q;F2.利用如下重心法公式求得Q的重心C:C=Wc卜v/、，《,、其中，w0cj)=a(/2D(x,y)—w)，a=—1，w=max(/2D(xj)),/^(x，力為Q中坐標為(x，力處的像素值;F3.連接所有重心點作為最終的血管中心路徑。本發明由於採用了八叉樹對CTA血管造影體數據進行分解，並利用原CTA血管數據中的相關信息計算血管組織的邊界距離場，以及利用重心法修正血管組織初始中心路徑，因此與現有技術相比具有如下優點(1)在血管提取和邊界距離場計算過程中只遍歷包含血管體素的節點，節省了時間；(2)利用原CTA血管數據中的相關信息計算邊界距離場，使處於不同位置的體素邊界距離值相差更為明顯，提高了中心路徑提取的精度；(3)利用重心法修正後的中心路徑更符合人類的觀察習慣，有利於血管疾病的診斷。圖l為本發明的流程圖2為本發明對CTA血管造影體數據進行八叉樹分解的示意圖；圖3是本發明中三類鄰域距離示意圖4為本發明中使用的兩個CTA血管造影體數據三維重建結果圖5為本發明中對頭部動脈瘤數據用不同方法提取的中心路徑與手動提取的中心路徑之間的距離對比示意圖6為本發明中對腦血管數據用不同方法提取的中心路徑與手動提取的中心路徑之間的距離對比示意圖7為利用本發明提取的中心路徑進行血管虛擬內窺鏡的效果圖。具體實施例方式本發明的核心思想是利用八叉樹對輸入的CTA血管造影體數據進行分解，並利用八叉樹中存儲的信息加速血管組織的提取和邊界距離場的計算;利用原CTA血管數據中的相關信息提高邊界距離場計算的精確度；採用重心法對提取的血管初始中心路徑進行修正，提高中心路徑的精度。參照圖1，本發明具體步驟如下步驟l，對輸入的CTA血管造影體數據進行八叉樹分解，得到每一分塊的信息。(la)按照如下方式定義八叉樹的屬性；StructOctreecharLevelX,LevelY，LevelZ;〃節點X,Y，Z方向上的深度，unsignedintXmin，Xmax,Ymin,Ymax,Zmin,Zmax;〃節點的坐木示範圍，boolIsLeaf;〃是否子節點，intMin_value，Max_value;〃節點內的最小最大體數據值，Octree*Children[8];〃節點的八個子節點；(lb)根據設定的八叉樹深度初始化八叉樹，逐層為當前八叉樹節點分配八個子節點；八叉樹中最上層的節點表示整個體數據場，深度為0，如圖2(a)所示，它被剖分成編號為0~7的八個深度為1的子節點，如圖2(b)所示；當深度為2時，每個深度為1的節點被剖分為8個字節點，如圖2(c)所示；以此類推，直到八叉樹節點的成員變量LevelX，LevelY和LevelZ與預先設定的深度值相等時，剖分結束。此時八叉樹節點為葉子節點，IsLeaf為True;剖分過程中同時記錄每一分塊的信息；體數據八叉樹分解的樹型表示如圖2(d)所示。步驟2，根據每個分塊的最大值和最小值與雙閾值之間的關係，遍歷存在血管體素的葉子節點提取血管組織。(2a)從深度為O的八叉樹節點開始，判斷其最大值和最小值與雙閾值之間的關係，如果有至少一個閾值在最大值和最小值範圍內，則該節點內存在血管體素，否則忽略該節點；(2b)判斷存在血管體素的節點的深度值與設定的八叉樹深度值之間的關係，如果節點的深度值與所設定的八叉樹的深度值相等，則該節點為葉子節點，並遍歷該葉子節點提取血管組織，否則判斷該節點的子節點的最大值和最小值與雙閾值之間的關係；(2c)如此循環，直到將所有包含血管體素的葉子節點遍歷完畢為止，完成血管組織的提取。步驟3，計算提取後的血管組織邊界距離場。(3a)根據提取的血管組織中所有點的坐標，利用原CTA血管造影體數據中血管體素的梯度值倒數和拉普拉斯變換值之和作為血管體素的邊界距離初始值其中，p為提取的血管組織中任意一點的坐標，/c為原CTA血管數據中點p的梯度值，/￡為原CTA血管數據中點;j的拉普拉斯變換值；(3b)根據邊界距離初始值用如下公式計算點^的邊界距離值y^,(p):(P)=min(min(/劍")+10),min(/劍(P2)+14),min(/Jtort(p3)+17》其中，為點P的最近鄰鄰域點A的邊界距離初始值，/stort(P2)為點P的面對角鄰域點A的邊界距離初始值，,^(^)為點;7的體對角鄰域點A的邊界距離初始值；三類鄰域距離如圖3所示；(3c)利用八叉樹遍歷所有包含血管體素的葉子節點，完成血管組織的邊界距離場計算。步驟4，構建提取後的血管組織最大代價生成樹。8借鑑Dijkstra動態規劃思想，首先將提取的所有血管體素看成樹的節點，除根節點S外，每個節點都指向一個父節點，從而形成一棵連接所有血管體素的有向樹。在樹的生成過程中利用邊界距離值，定義^指向萬的權值為『AB-DFB(5)，即如果萬為^的父節點，貝U^連通到5的代價為DFB(5)。最大代價生成樹的建立過程如下(4a)從源點S開始，將它所有26-鄰域點均指向S，並以S為父節點，然後把它們的邊界距離值放入一個排序隊列《中；^(4b)從排序隊列《中取出邊界距離值最大的點M作為當前處理的點，讓M的26-鄰域中未被處理過的點均指向M，並以M為父節點，則這些連接邊的權值都是DFB(M)，再把這些點的邊界距離值放入排序隊列《中；(4c)如此循環，直到將所有提取的血管體素都處理完為止。在建立最大代價生成樹過程中，採用快速隊列排序方法以提高效率。首先建立一個動態査找表，表的索引範圍為0255;接著在(4a)中把每個邊界距離值歸到0255，再把該點的體素序號放入對應的查找表中；在(4b)中始終從索引值最大的不為空的査找表中提取體素序號，直到所有的查找表都為空。步驟5，提取最大代價生成樹的主幹作為血管初始中心路徑。最大代價生成樹主幹的提取一般分為兩步(5a)如果用戶沒有指定結束點，則默認源距離值最大的點為路徑終點五；(5b)從終點五開始，根據節點的指向關係直接連通到源點S作為血管初始中心路徑。步驟6，利用重心法對血管初始中心路徑進行修正，使路徑上的點都處於血管管腔橫截面的重心位置，得到一條符合人類觀察習慣的最終血管中心路徑。(6a)對於初始中心路徑上的每一點，取與該點法向量垂直，並且與該點相距為10個體素的點組成的圓形二維圖像Q;(6b)利用如下重心法公式求得Q的重心C:其中，wO,力-a(/2D(x,力—m)，"=一1，w=max(/2D(>:,:v))，/20(;:，力為￡2中坐標為(>;，>0處的像素值；(6c)連接所有重心點作為最終的血管中心路徑。本發明的優點可通過以下實驗進一步說明1.實驗條件本發明的實驗是在德國圖賓根大學的兩組CTA血管造影體數據上進行的。一組是頭部動脈瘤數據，解析度為512體素X512體素X512體素，實際大小為100X100X100mm3，三個方向上的相鄰體素間隔均為0.1953mm，數據量為128MBytes。另一組是腦血管數據，解析度為416體素x512體素xll2體素，實際大小為171mmx210.5mmx46mm，三個方向上的相鄰體素間隔均為0.412mm，數據量為22.7MBytes。頭部動脈瘤數據的三維重建結果如圖4(a)所示，腦血管數據的三維重建結果如圖4(b)所示。實驗中提取了四種路徑，一是多次手動提取的路徑求平均作為一條最優路徑，並將該路徑作為與其他三條由方法所提取的路徑進行比較的基準，稱之為基準路徑；二是基於邊界距離場方法提取的路徑，但未使用梯度值及拉普拉斯變換值作為邊界距離變換初始值，且未使用重心法修正，稱之為二類路徑；三是基於邊界距離場方法提取的路徑，使用梯度值及拉普拉斯變換值作為邊界距離變換初始值，但未使用重心法修正，稱之為三類路徑；四是本發明基於邊界距離場方法提取的路徑，既使用了梯度值及拉普拉斯變換值作為邊界距離變換初始值，也使用了重心法修正，稱之為四類路徑。實驗在這兩組體數據中比較用八叉樹分解的方法進行血管提取所用的時間與用全局遍歷方法進行血管提取所用的時間；比較用八叉樹分解的方法進行邊界距離場計算所用的時間與用全局遍歷方法進行邊界距離場計算所用的時間；比較二類路徑、三類路徑、四類路徑與基準路徑的距離。2.實驗結果實驗中兩種方法在兩組數據中所用時間對比如表1和表2所示。兩組數據中所提取的二類路徑，三類路徑，四類路徑與基準路徑之間的距離對比如圖5和圖6所示。表1頭部動脈瘤數據八叉樹方法與全局遍歷方法血管提取及邊界距離場計算用時比較tableseeoriginaldocumentpage10表2腦血管數據八叉樹方法與全局遍歷方法血管提取及邊界距離場計算用時比較tableseeoriginaldocumentpage10表1中，提取閾值在三種不同情況下，八叉樹大小為16體素xl6體素xl6體素時，血管提取用時與邊界距離場計算用時之和分別為1.808秒，1.719秒，1.625秒，比同閾值下全局遍歷所用時間之和4.227秒，4.242秒，4.204秒分別要少57.2%，59.4%，61.3%。表2中，提取閾值在三種不同情況下，八叉樹大小為64體素x64體素x64體素時，血管提取用時與邊界距離場計算用時之和分別為0.037秒，0.036秒，0.030秒，比同閾值下全局遍歷所用時間之和3.969秒，3.938秒，3.906秒分別要少99%，99%，99.2%。從表1和表2中可以看出，加入八叉樹分解後，血管的提取速度和邊界距離場的計算速度得到了大幅度的提高，節省了時間消耗。實驗中共提取了兩段血管的中心路徑，圖5中的路徑共有85個點構成，圖6中的路徑共有68個點構成。兩圖中，黑色線為二類路徑與基準路徑的距離，藍色線為三類路徑與基準路徑的距離，紅色線為四類路徑與基準路徑的距離；圖的橫坐標為路徑上點的序號，縱坐標為該點與基準路徑上對應點的距離值，單位為體素；圖中*，+，o均表示路徑上的點。圖5中，黑色線的均值為3.69，即二類路徑與基準路徑之間的距離均值為3.69個體素；藍色線的均值為2.85，即三類路徑與基準路徑之間的距離均值為2.85個體素；紅色線的均值為2.21，即四類路徑與基準路徑之間的距離均值為2,21個體素。圖6中，黑色線的均值為3.46,即二類路徑與基準路徑之間的距離均值為3.46個體素；藍色線的均值為2.36，即三類路徑與基準路徑之間的距離均值為2.36個體素；紅色線的均值為1.68，即四類路徑與基準路徑之間的距離均值為1.68個體素。直觀上，在兩圖中，黑色線在藍色線之上，而藍色線又在紅色線之上，即四類路徑最接近基準路徑。因此，利用梯度值倒數及拉普拉斯變換值作為邊界距離變換初始值，且使用重心法修正從最大代價生成樹提取的路徑可使最終路徑最接近血管空腔的真實中心路徑。利用四類路徑對頭部動脈瘤數據進行血管虛擬內窺鏡的效果如圖7(a)所示，利用四類路徑對腦血管數據進行血管虛擬內窺鏡的效果如圖7(b)所示。上述實驗結果中，本發明用八叉樹方法進行血管提取和邊界距離場計算，與目前存在的全局遍歷方法進行血管提取和邊界距離場計算相比，用時要少；利用原CTA血管造影體數據的相關信息得到的距離場，比目前存在的利用二值血管數據得到的距離場精度要高；利用重心法修正後的路徑比直接提取基於邊界距離場的最大代價生成樹主幹得到的中心路徑更貼近人類的觀察習慣。1權利要求1.一種基於邊界距離場的血管中心路徑提取方法，包括如下過程A.對輸入的CTA血管造影體數據進行八叉樹分解，得到每一分塊的信息；B.根據每個分塊的最大值和最小值與雙閾值之間的關係，遍歷存在血管體素的葉子節點提取血管組織；C.根據提取的血管組織中所有點的坐標，計算提取後的血管組織邊界距離場；D.在邊界距離場的基礎上建立最大代價生成樹；E.提取最大代價生成樹的主幹作為血管組織的初始中心路徑；F.利用重心法修正提取的血管組織初始中心路徑，使路徑上的點都處於血管管腔橫截面的重心位置，得到一條符合人類觀察習慣的最終血管中心路徑。2.根據權利要求1所述的方法，其中步驟A所述的對輸入的CTA血管造影體數據進行八叉樹分解，按如下步驟進行Al.按照如下方式定義八叉樹的屬性；StructOctreecharLevelX,LevelY,LevelZ;〃節點X，Y,Z方向上的深度，unsignedintXmin，Xmax，Ymin，Ymax,Zmta,Zmax;〃節點的坐標範圍，boolIsLeaf;〃是否子節點，intMin—value,Max—value;〃節點內的最小最大體數據值，Octree*Children[8];〃節點的八個子節點；A2.根據所設定的八叉樹的深度值LevelX，LevelY，LevelZ對輸入的CTA血管造影體數據進行逐層分解，同時利用八叉樹記錄每一分塊的信息。3.根據權利要求1所述的方法，其中步驟B所述的根據每個分塊的最大值和最小值與雙閾值之間的關係，遍歷存在血管體素的葉子節點提取血管組織，按如下步驟進行Bl.從深度為0的八叉樹節點開始，判斷其最大值和最小值與雙閾值之間的關係，如果有至少一個閾值在最大值和最小值範圍內，則該節點內存在血管體素，否則忽略該節點；B2.判斷存在血管體素的節點的深度值與設定的八叉樹深度值之間的關係，如果節點的深度值與所設定的八叉樹的深度值相等，則該節點為葉子節點，並遍歷該葉子節點提取血管組織，否則判斷該節點的子節點的最大值和最小值與雙閾值之間的關係；B3.如此循環，直到將所有包含血管體素的葉子節點遍歷完畢為止，完成血管組織的提取。4.根據權利要求1所述的方法，其中步驟C所述的計算提取後的血管組織邊界距離場，按如下步驟計算Cl.根據提取的血管組織中所有點的坐標，利用原CTA血管造影體數據中血管體素的梯度值倒數和拉普拉斯變換值之和作為血管體素的邊界距離初始值Xwrt(P)=l//c+//_其中，；為提取的血管組織中任意一點的坐標，/e為原CTA血管數據中點p的梯度值，/,為原CTA血管數據中點p的拉普拉斯變換值；C2.根據邊界距離初始值用如下公式計算點p的邊界距離值O)=min(min(/劍)+10)，min(/咖(p2)+14),min(/加,,(p3)+17》其中，,^(A)為點P的最近鄰鄰域點A的邊界距離初始值，,^(A)為點P的面對角鄰域點a的邊界距離初始值，y;。rt073)為點/的體對角鄰域點^的邊界距離初始值；C3.利用八叉樹遍歷所有包含血管體素的葉子節點，完成提取後的血管組織邊界距離場計算。5.根據權利要求1所述的方法，其中步驟F所述的利用重心法對提取的血管組織初始中心路徑進行修正，按如下步驟進行Fl.對於初始中心路徑上的每一點，取與該點法向量垂直，並且與該點相距為IO個體素的點組成的圓形二維圖像Q;F2.利用如下重心法公式求得Q的重心C:&，義L1>(1，力，JX;c,力乂其中，w(;cj)=c0—w)，《=—1，w-max(/2D(x,力)，/200,力為0中坐標為(x，力處的像素值；F3.連接所有重心點作為最終的血管中心路徑。全文摘要本發明公開了一種基於邊界距離場的血管中心路徑高精度快速提取方法，主要解決在CTA血管造影體數據中血管提取速度慢和中心路徑提取精度不高的問題。該方法包括利用八叉樹對體數據進行分塊，並記錄所有分塊信息；利用每一分塊的最大值和最小值與雙閾值的關係，對血管進行提取和邊界距離場計算；用原CTA數據中血管體素的梯度值倒數和拉普拉斯變換值之和作為邊界距離初始值；建立基於邊界距離場的最大代價生成樹，並提取樹的主幹作為血管的初始中心路徑；用重心法修正初始中心路徑，使路徑上的點都處於血管管腔橫截面的重心處，得到最終的血管中心路徑。本發明具有血管提取速度快，中心路徑提取精度高的優點，適用於各種管腔組織的中心路徑提取。文檔編號G06T7/00GK101551862SQ200910022500公開日2009年10月7日申請日期2009年5月13日優先權日2009年5月13日發明者呂新榮,寧文博,曹豔豔,李晉舟,林查,梁承志,王向華,華鄒,高新波申請人:西安電子科技大學