血管的三維重建VIP

血管的三維重建問(wèn)題 摘要 越來(lái)越多的在醫(yī)學(xué)領(lǐng)域及生物領(lǐng)域需要大量的切片工作 來(lái)研 究切片表面的組織形態(tài) 因而如何能將這些切片在利用計(jì)算機(jī)進(jìn)行 復(fù)原 顯得很重要 在對(duì)管道的半徑的求解中 我們采用的是枚舉法 即逐個(gè)將 100 切片的最大內(nèi)切圓的半徑都求出來(lái) 最后取其平均值 即為我 們所求的管道的半徑 在半徑求出的同時(shí)把最大內(nèi)切圓的圓心也求 出來(lái) 要求得切片的最大內(nèi)切圓 首先對(duì)切片圖像 0 1 矩陣轉(zhuǎn)化 然后利用 edge 函數(shù)找到圖形的輪廓 利用 bwmorph 函數(shù)求得圖形 的輪廓 然后利用 find 函數(shù)求出輪廓與骨架上的坐標(biāo) 利用兩點(diǎn)間 的距離公式 我們就能求出骨架上每個(gè)點(diǎn)到輪廓的距離 在所有點(diǎn) 到輪廓的最短距離中的那個(gè)最大值 就是我們要求的最大內(nèi)切圓半 徑 此圓心即為中軸線上的點(diǎn) 在中軸線與中軸線在 XY YZ ZX 平面的投影的計(jì)算中 我們主要 采用多項(xiàng)式擬合 polyfit 的方法 從得出擬合曲線可以看出與實(shí)際圖 形很接近 最后我們擬合出了血管的三維形態(tài) 關(guān)鍵詞 最大內(nèi)切圓 半徑 中軸線 多項(xiàng)式擬合 一 問(wèn)題重述 如今在醫(yī)學(xué)領(lǐng)域 生物學(xué)領(lǐng)域等都需要了解生物組織 器官 等的斷面形態(tài) 以便于能夠準(zhǔn)確的研究該樣品 得出具有說(shuō)服力的 研究成果 例如 將樣本染色后切成厚約 1 m 的切片 在顯微鏡下 觀察該橫斷面的組織形態(tài)結(jié)構(gòu) 如果用切片機(jī)連續(xù)不斷地將樣本切 成數(shù)十 成百的平行切片 可依次逐片觀察 這時(shí)候我們?nèi)绻覀?需要知道血管的三維形態(tài) 根據(jù)拍照并采樣得到的平行切片數(shù)字圖 象 運(yùn)用計(jì)算機(jī)就可重建組織 器官等的三維形態(tài) 假設(shè)某些血管 可視為一類特殊的管道 該管道的表面是由球心沿著某一曲線 稱 為中軸線 的球滾動(dòng)包絡(luò)而成 本文給出了某管道的相繼 100 張平行切片圖象 記錄了管道與 切片的交 圖象文件名依次為 0 bmp 1 bmp 99 bmp 格式 均為 BMP 寬 高均為 512 個(gè)象素 pixel 我們要做的就是根據(jù) 這些切片計(jì)算管道的中軸線與半徑 給出具體的算法 并繪制中軸 線在 XY YZ ZX 平面的投影圖 二 問(wèn)題分析 本題是一個(gè)求最優(yōu)解的問(wèn)題 通過(guò)具體的算法得出最接近實(shí)際 血管的中軸線與半徑 并能較準(zhǔn)確的繪制出中軸線在 XY YZ ZX 平 面的投影圖 得出直觀的圖形 需要處理大量的圖片 得出大量的 數(shù)據(jù)信息 所以為處理該問(wèn)題的程序中需要大量的循環(huán) 所以程序 的完成與運(yùn)行皆很復(fù)雜 本題管道為由球心沿著某一曲線的球滾動(dòng)包絡(luò)而成 且管道中 軸線與每張切片被假設(shè)有且只有一個(gè)交點(diǎn) 則每張切片一定包含球 的最大截面 即過(guò)球心的圓面 則每張切片的最大內(nèi)切圓就是過(guò)球 心的圓面 由于本題中取坐標(biāo)系的 Z 軸垂直于切片 第 1 張切片為平面 Z 0 第 100 張切片為平面 Z 99 Z z 切片圖象中象素的坐標(biāo)依它們 在文件中出現(xiàn)的前后次序?yàn)?256 256 z 256 255 z 256 255 z 255 256 z 255 255 z 255 255 z 255 256 z 255 255 z 255 255 z 則每個(gè)切片的最大內(nèi)切圓的圓心的 z 坐標(biāo)已給 再加上 x y 則中心 軸就是每個(gè)切片的最大內(nèi)切圓的圓心 x y z 連成的光滑曲線 而管 道的半徑就是這個(gè)球的半徑 即最大內(nèi)切圓的半徑 進(jìn)而再求出中 軸線在 XY YZ ZX 平面的投影圖 那么這道題的問(wèn)題最后就轉(zhuǎn)化為 求每張切片的最大內(nèi)切圓的圓心與半徑問(wèn)題 三 模型假設(shè) 1 管道的表面是由球心沿著某一曲線 稱為中軸線 的球滾動(dòng)包絡(luò) 而成 2 管道中軸線與每張切片有且只有一個(gè)交點(diǎn) 3 球半徑固定 4 切片間距以及圖象象素的尺寸均為 1 四 符號(hào)說(shuō)明 一張切片內(nèi)骨架點(diǎn) i 到邊緣點(diǎn) 的距離 第 i 張圖的最大內(nèi)切圓半徑 半徑的平均值 五 模型的建立與求解 問(wèn)題一 求管道的半徑與每張切片的最大內(nèi)切圓的圓心 要求得管道的半徑 首先需要求出每個(gè)切片的最大內(nèi)切圓 為了得出比較精確的結(jié)果 我們采用枚舉法 對(duì)與 100 張切片 考 慮到求得的最大內(nèi)切圓半徑由于誤差的存在肯定會(huì)有所不同 我們 將所有切片的最大內(nèi)切圓半徑求出后 然后取其平均值 為得到每張圖片的最大內(nèi)切圓半徑 首先先將二值圖片轉(zhuǎn)換 成 0 1 矩陣 矩陣橫縱坐標(biāo)對(duì)應(yīng)原圖像的直角坐標(biāo)系位置 1 代表 黑像素 0 代表白像素 為了后面 find 函數(shù)尋找矩陣中為 1 的坐 標(biāo)并記錄 通過(guò)相關(guān)函數(shù)找出輪廓 edge 和骨架 bwmorph 記錄 輪廓 q i j 和骨架 p i j 的坐標(biāo) 利用兩點(diǎn)間的公式去求骨架上的點(diǎn)到邊界上所有點(diǎn)的最小距離即骨 架上此點(diǎn)的內(nèi)切圓的半徑 1 骨架上所有點(diǎn)的最大內(nèi)切圓即為該切片的最大內(nèi)切圓 2 這里 i 代表第 i 張切片 程序見(jiàn)附錄 1 其半徑即為最大內(nèi)切圓半徑 通過(guò)循環(huán)找出最大半徑并記錄該 半徑所對(duì)的坐標(biāo)點(diǎn) 基于此我們能求出所有切片的最大內(nèi)切圓的半 徑與圓心 數(shù)據(jù)見(jiàn)附錄 5 通過(guò)求得的 100 張切片的最大內(nèi)切圓的半徑我們?nèi)∑淦骄档?出更加精確的管道的半徑 問(wèn)題二 求中心軸與中心軸在 XY YZ ZX 平面的投影 x y z 為切片的最大內(nèi)切圓的圓心坐標(biāo) 1 利用問(wèn)題一中求得數(shù)據(jù) x y z 然后根據(jù)多項(xiàng)式擬合 polyfit 得出 x y z 的擬合方程為 中軸線在 Z Y 平面的投影擬合多項(xiàng)式 中軸線在 X Y 平面的投影擬合多項(xiàng)式 中軸線在 Z X 平面的投影擬合多項(xiàng)式 程序見(jiàn)附錄 2 得到的擬合曲線如下 圖 1 為中軸線在 XY 平面的投影 曲線代表多項(xiàng)式擬合的結(jié)果 折線代表真實(shí)數(shù)據(jù)的圖形 圖 2 為中軸線在 ZY 平面的投影 曲線代表多項(xiàng)式擬合的結(jié)果 折線代表真實(shí)數(shù)據(jù)的圖形 圖 為中軸線在 ZX 平面的投影圖 曲線代表多項(xiàng)式擬合的結(jié)果 折線代表真實(shí)數(shù)據(jù)的圖 形 從以上擬合圖形中我們可以看出 通過(guò)多項(xiàng)式擬合的投影圖與 通過(guò)坐標(biāo)點(diǎn)得出的投影圖基本相同 所以可以近似認(rèn)為擬合的投影 圖即為我們所要求得的中軸線在 XY YZ ZX 平面的投影 2 同樣利用多項(xiàng)式擬合 polyfit 得出管道的中軸線在三維坐標(biāo)中的 擬合曲線如下 程序見(jiàn)附錄 3 圖 為中軸線在三維坐標(biāo)的立體圖 3 利用球坐標(biāo)以及以上求出的最大半徑 以及擬合時(shí)在 z 時(shí)的坐標(biāo) x1 y1 利用 plot3 程序如附錄四 圖 血管的三維重建 六 模型評(píng)價(jià)及總結(jié) 解決本題的首要前提是要理解像素的概念 在解此題的最初由 于對(duì)圖片的像素沒(méi)有明確的理解而浪費(fèi)了很長(zhǎng)時(shí)間 在我們的模型 中還有尚有不完美的地方 在求切片的最大內(nèi)切圓半徑與圓心的時(shí) 候 沒(méi)有編寫出能夠一下就將所有切片的數(shù)據(jù)都求出的程序 而是 僅僅編寫出一個(gè)圖片的程序 然后逐一求出 在這個(gè)問(wèn)題上花費(fèi)的 時(shí)間長(zhǎng) 且過(guò)程繁瑣 但是在進(jìn)行曲線擬合時(shí) 由于運(yùn)用 知道了相關(guān)函數(shù)的用法 運(yùn)用 polyfit 進(jìn)行擬合 很大程度上節(jié)省 了時(shí)間 擬合出的相對(duì)來(lái)說(shuō)比較滿意 立體切片的三維重建在醫(yī)學(xué) 與生物學(xué)等領(lǐng)域的作用很大 因?yàn)檫@些領(lǐng)域往往會(huì)研究一些物體切 片的組織結(jié)構(gòu) 所以此問(wèn)題的實(shí)現(xiàn)能夠解決很多現(xiàn)實(shí)問(wèn)題 為一些 研究領(lǐng)域提供很大的幫助 七 參考文獻(xiàn) 1 董辰輝 MATLAB2008 全程指南 電子工業(yè)出版社 2 張紅云等 圖像骨架的提取的應(yīng)用 2010 年第四期 3 4 附錄 1 jieguo zeros 100 4 J0 imread E 99 BMP for i 1 512 for j 1 512 j0 i j 1 J0 i j 轉(zhuǎn)化為黑色為 1 白色為 0 為了后面 find 函數(shù)尋找矩陣中 為 1 的坐標(biāo)并記錄 end end lk edge j0 sobel gj bwmorph j0 skel inf x0 y0 v0 find lk 找到輪廓灰色區(qū)域 a0 b0 c0 find gj 找到骨架灰色區(qū)域 m length a0 骨架灰色區(qū)域個(gè)數(shù) n length x0 輪廓灰色區(qū)域個(gè)數(shù) jl zeros m n 建立 0 矩陣為求內(nèi)切圓半徑 cf zeros m 2 建立兩列 0 矩陣為存放中心點(diǎn)坐標(biāo) for i 1 m for j 1 n p1 a0 i q1 b0 i p2 x0 j q2 y0 j 骨架 輪廓坐標(biāo)賦值 jl i j sqrt p1 p2 2 q1 q2 2 end zx zxxh min jl i 骨架上一點(diǎn)到輪廓的最短距離 即骨架上各個(gè)點(diǎn)的內(nèi)切圓的半徑 cf i 1 zx cf i 2 zxxh end zd zdxh max cf 1 找到其中最大的半徑 并把半徑和圓心坐標(biāo)存儲(chǔ) x a0 zdxh 256 與題目所給坐標(biāo)軸對(duì)應(yīng) y b0 zdxh 256 與題目所給坐標(biāo)軸對(duì)應(yīng) jieguo k 1 1 x x 軸坐標(biāo) jieguo k 1 2 y jieguo k 1 3 k 1 jieguo k 1 4 zd 2 polyfit Z X 7 plot Z X 0 10 99 polyval ans 0 10 99 p s polyfit Z X 7 polyfit Z Y 7 plot Z Y 0 10 99 polyval ans 0 10 99 p s polyfit Z Y 7 plot X Y 160 10 188 polyval ans 160 10 188 3 format long px polyfit z x 7 x1 polyval px z 取擬合后再 z 點(diǎn)時(shí)的點(diǎn) py polyfit z y 7 y1 polyval py z figure 1 plot3 x1 y1 z r t linspace 0 pi 25 存儲(chǔ)平均分 pi 為 25 份存儲(chǔ) p linspace 0 2 pi 25 存儲(chǔ)平均分 2pi 為 25 份存儲(chǔ) theta phi meshgrid t p 存儲(chǔ)角度 for i 1 100 ceen 1 x1 ceen 2 y1 ceen 3 z x 29 49288 sin theta sin phi ceen i 1 y 29 49288 sin theta cos phi ceen i 2 z1 29 49288 cos theta ceen i 3 hold on 使圖像不被覆蓋 plot3 x y z1 axis equal hold off end 5 切片號(hào)半徑x 坐標(biāo)y 坐標(biāo)z 坐標(biāo)切片號(hào)半徑x 坐標(biāo)y 坐標(biāo)z 坐標(biāo) 129 069 160005130 414 11411750 228 284 160215230 414 11411751 329 160125330 414 11311852 429 069 160135430 414 11211953 529 069 160245530 11112054 629 069 160255629 732 11112055 729 160165729 698 11112056 829 017 160475829 698 11112057 929 160185929 53 8114258 1028 862 160196029 547 5115659 1128 862 1607106129 547 5115660 1228 862 1608116229 614 3116261 1328 862 1609126329 614 3116262 1429 017 16010136429 614 3116263 1529 017 16012146529 614 3516164 1629 017 16013156629 614 3516165 1729 017 16014166729 428 2616366 1829 017 16016176829 411 3516167 1929 017 16017186929 275 2616368 2029 017 16018197029 4284616369 2129 017 16019207129 6144616370 2229 017 16020217229 6144616371 2329 017 16021227329 6144616372 2429 017 16022237429 6146515873 2529 069 16021247529 7326815774 2629 069 16021257629 7326515875 2729 069 16021267729 5478115276 2829 155538115277 2929 275538115278 3029 27569813511779 3129 42869813611680 3229 61469813611681 3329 61469813711582 3429 61469813811483 3529 61469813811484 3629 73269813911385 3729 73269813911386 3829 73269813911387 3929 73269814011288 4029 61469814011289 4129 5475317267