專利名稱：一種數(shù)字減影血管造影運(yùn)動(dòng)偽影消除方法及其系統(tǒng)的制作方法
技術(shù)領(lǐng)域：
本發(fā)明涉及醫(yī)學(xué)圖像處理技術(shù)領(lǐng)域。特別是一種造影圖像中運(yùn)動(dòng)偽影的消除方法 和偽影消除系統(tǒng)。
背景技術(shù)：
隨著現(xiàn)代科學(xué)技術(shù)的發(fā)展，一些先進(jìn)的技術(shù)手段和計(jì)算機(jī)科學(xué)技術(shù)正在不斷地應(yīng) 用于醫(yī)學(xué)領(lǐng)域。特別是在放射診斷和治療學(xué)科中，計(jì)算機(jī)圖像處理技術(shù)正發(fā)揮著越來(lái)越重 要的作用。
目前，數(shù)字減影血管造影(DigitalSubtractionAngiography，DSA)是血管可視化的關(guān) 鍵技術(shù)，國(guó)內(nèi)最早于1984年引進(jìn)了 DSA設(shè)備，并在全國(guó)迅速推廣。在普通的X光圖像中，由 于血管與其周圍組織的對(duì)比度非常小，若要單純對(duì)血管進(jìn)行細(xì)微觀察就較為困難。為了獲 得更清晰的血管圖像，通常將造影劑(不透光的碘溶液)注入需要診斷的血管中以提高血管 的對(duì)比度。病人在檢查臺(tái)上接受多次X射線的照射，第一次是未注射造影劑時(shí)拍攝的X光 片，稱為蒙片(MaskImage)。經(jīng)導(dǎo)管向病變區(qū)域或檢查區(qū)域血管內(nèi)快速注入有機(jī)碘水造影 齊U，注入血管造影劑之后拍攝的X光片稱為活片或盈片(Livelmage)。從注入造影劑開(kāi)始到 血管內(nèi)造影劑濃度達(dá)到峰值這段時(shí)間內(nèi)，使檢查部位連續(xù)成像。在這一系列圖像中，取蒙片 圖像與造影劑濃度達(dá)到峰值的盈片圖像執(zhí)行減法操作，可以去除非感興趣組織(如肌肉、骨 骼)對(duì)X射線吸收形成的背景圖像，從而得到清晰的血管造影圖像。假設(shè)在血管周圍組織、成像條件完全不變的情況下，將蒙片與活片相減，如果活片 與蒙片的背景(如肌肉、骨骼)完全對(duì)應(yīng)，且有相同的灰度分布，則這些背景將被完全消除， 即理論上可以得到僅包含血管的清晰數(shù)字減影圖像。但在臨床應(yīng)用中，這種理論假設(shè)總是 不成立的。由于蒙片和活片序列在不同時(shí)刻拍攝，在成像的同時(shí)，病人會(huì)存在各種運(yùn)動(dòng)，其 中有些運(yùn)動(dòng)是不可避免的，比如呼吸、吞咽、眼球轉(zhuǎn)動(dòng)、肌肉運(yùn)動(dòng)、臟器蠕動(dòng)、病人的躁動(dòng)不 配合以及在造影劑注人體時(shí)由于病人的自然反應(yīng)而引起的局部運(yùn)動(dòng)等。這些運(yùn)動(dòng)會(huì)使活片 和蒙片之間存在并非由造影劑產(chǎn)生的差異，從而使減影后的圖像出現(xiàn)難以預(yù)料的背景，稱 ^jiE 云力 i^di (MotionArtifact)。在早期DSA中，減少運(yùn)動(dòng)偽影的技術(shù)主要是在如何去避免病人的運(yùn)動(dòng)，以及 改進(jìn)減影方法(如雙能量減影技術(shù))。近年來(lái)的研究中，運(yùn)動(dòng)偽影主要通過(guò)圖像配準(zhǔn) (ImageRegistration)的方法被校正。圖像配準(zhǔn)就是將不同條件下獲取的兩幅或多幅圖像 進(jìn)行匹配、疊加的過(guò)程，實(shí)質(zhì)上是兩個(gè)圖像坐標(biāo)系之間的矩陣變換的計(jì)算過(guò)程，可以為剛性 配準(zhǔn)或非剛性配準(zhǔn)。圖像配準(zhǔn)也可以看作是一個(gè)運(yùn)動(dòng)估計(jì)問(wèn)題，為了獲得兩幅圖像之間的 對(duì)應(yīng)關(guān)系，必須進(jìn)行運(yùn)動(dòng)估計(jì)。運(yùn)動(dòng)估計(jì)就是自動(dòng)計(jì)算圖像中局部運(yùn)動(dòng)或某個(gè)結(jié)構(gòu)位移的 技術(shù)。根據(jù)DSA原理，造影劑濃度達(dá)到峰值的盈片圖像與蒙片圖像相減，從而得到清晰 的血管造影圖像。從注入造影劑開(kāi)始到血管內(nèi)造影劑濃度達(dá)到峰值這段時(shí)間內(nèi)(約有數(shù)秒)，檢查部位連續(xù)成像，將產(chǎn)生一系列血管逐漸清晰的盈片圖像，即在一個(gè)DSA的DICOM文件中包含了幾十張盈片圖像數(shù)據(jù)。在DSA圖像序列中，運(yùn)動(dòng)像素不僅與同一幅盈片圖像中 周圍的像素相關(guān)(空間相關(guān)性)，而且也與前后幾幅盈片圖像中周圍的像素相關(guān)(時(shí)間相關(guān) 性)?；跁r(shí)空信息的分析方法近年來(lái)在數(shù)字視頻處理中已經(jīng)得到了廣泛的應(yīng)用，但是由于 DSA圖像序列中大部分是利用價(jià)值不高的、造影劑濃度沒(méi)有達(dá)到峰值的盈片圖像(包括到達(dá) 峰值前以及部分消退圖像)，所以導(dǎo)致消除運(yùn)動(dòng)偽影的效果不佳。

發(fā)明內(nèi)容
本發(fā)明的目的是提供一種數(shù)字減影血管造影運(yùn)動(dòng)偽影消除方法及其系統(tǒng)，要解決 數(shù)字減影血管造影圖像中的運(yùn)動(dòng)偽影的消除技術(shù)問(wèn)題，并解決盈片圖像與蒙片圖像精細(xì)配 準(zhǔn)的技術(shù)問(wèn)題。為實(shí)現(xiàn)上述目的，本發(fā)明采用如下技術(shù)方案一種數(shù)字減影血管造影運(yùn)動(dòng)偽影消 除方法，其特征在于步驟如下步驟1、讀??；將數(shù)字減影血管造影儀器與計(jì)算機(jī)連接，在計(jì) 算機(jī)中安裝操作系統(tǒng)和圖像處理系統(tǒng)，用圖像處理軟件讀入一組DSA圖像，并將DSA圖像進(jìn) 行預(yù)處理。步驟2、選點(diǎn)；采用邊緣檢測(cè)算法找到每幅DSA圖像中血管影像的邊緣和邊緣點(diǎn)， 計(jì)算出邊緣點(diǎn)的梯度值，根據(jù)閾值來(lái)判斷該點(diǎn)是否作為備選控制點(diǎn)，并選取至少10個(gè)備選 控制點(diǎn)。步驟3、構(gòu)建DSA時(shí)空體；將該組DSA 二維圖像以時(shí)間T為軸按拍攝順序間隔排列 成平行六面體空間。步驟4、時(shí)空切片；用切片指令對(duì)DSA時(shí)空體做垂直方向的時(shí)空切片。步驟5、連軌跡；用圖像處理軟件將各時(shí)空切片對(duì)應(yīng)的各控制點(diǎn)連線，得到各個(gè)控 制點(diǎn)的運(yùn)動(dòng)軌跡。步驟6，分析DSA像素的時(shí)空運(yùn)動(dòng)特性；選擇運(yùn)動(dòng)軌跡變化最大的10個(gè)反映DSA圖 像背景結(jié)構(gòu)運(yùn)動(dòng)規(guī)律的控制點(diǎn)。步驟7、三角剖分；用三角剖分的方法將所有的控制點(diǎn)作為三角形的頂點(diǎn)，把圖像 剖分成不均勻的三角形網(wǎng)格。步驟8、仿射變換(affinetransfomation)；采用仿射變換法進(jìn)行初步幾何校正。步驟9、時(shí)空卷繞；用Warping非線性變換法參考前一幅圖像中控制點(diǎn)的變形位 置，逐幅確定圖像非控制點(diǎn)的位移。步驟10、優(yōu)化；采用Powell優(yōu)化算法來(lái)優(yōu)化相似度準(zhǔn)則函數(shù)求極值的過(guò)程。步驟11、配準(zhǔn)；選擇三次樣條插值算法進(jìn)行亞像素級(jí)的配準(zhǔn)。步驟12、灰度校正；逐幅對(duì)圖像的灰度失真進(jìn)行回朔校正。步驟13，對(duì)數(shù)減影；用對(duì)數(shù)減影法進(jìn)行減影，獲得消除偽影的減影圖像。所述邊緣檢測(cè)算法是Carmy算子或高斯-拉普拉斯算子 (Gauss—Laplacianoperator)ο一種數(shù)字減影血管造影運(yùn)動(dòng)偽影消除系統(tǒng)，包括DSA設(shè)備、實(shí)現(xiàn)運(yùn)動(dòng)偽影消除的 DSA軟件工作站、DSA膠片打印機(jī)、圖像顯示系統(tǒng)、外部存儲(chǔ)設(shè)備、網(wǎng)絡(luò)傳輸系統(tǒng)，它的DSA軟 件工作站包括以下部分輸入接口 用于接收DSA設(shè)備輸出的DICOM圖像。
存儲(chǔ)器用于存儲(chǔ)操作系統(tǒng)、圖像處理軟件、圖像分析數(shù)據(jù)和圖像資料； 數(shù)字減影部分采用本發(fā)明的算法進(jìn)行數(shù)字減影。DSA軟件工作站用于讀入DSA圖像，并將DSA圖像進(jìn)行預(yù)處理；采用邊緣檢測(cè)算 法找到每幅DSA圖像中血管影像的邊緣和邊緣點(diǎn)，計(jì)算出邊緣點(diǎn)的梯度值，根據(jù)閾值來(lái)判 斷該點(diǎn)是否作為備選控制點(diǎn)，并選取至少10個(gè)備選控制點(diǎn)；用于構(gòu)建DSA時(shí)空體；將該組 DSA 二維圖像以時(shí)間T為軸按拍攝順序間隔排列成平行六面體空間；用于進(jìn)行時(shí)空切片；用 切片指令對(duì)DSA時(shí)空體做垂直方向的時(shí)空切片；用于將各時(shí)空切片對(duì)應(yīng)的各控制點(diǎn)連線， 得到各個(gè)控制點(diǎn)的運(yùn)動(dòng)軌跡；用于分析DSA像素的時(shí)空運(yùn)動(dòng)特性；選擇運(yùn)動(dòng)軌跡變化最大 的10個(gè)反映DSA圖像背景結(jié)構(gòu)運(yùn)動(dòng)規(guī)律的控制點(diǎn)；用三角剖分的方法將所有的控制點(diǎn)作為 三角形的頂點(diǎn)，把圖像剖分成不均勻的三角形網(wǎng)格；用于采用仿射變換法進(jìn)行初步幾何校 正；用于作時(shí)空卷繞；用Warping非線性變換法參考前一幅圖像中控制點(diǎn)的變形位置，逐幅 確定圖像非控制點(diǎn)的位移；采用Powell優(yōu)化算法來(lái)優(yōu)化相似度準(zhǔn)則函數(shù)求極值的過(guò)程；選 擇三次樣條插值算法進(jìn)行亞像素級(jí)的配準(zhǔn)；用于逐幅對(duì)圖像的灰度失真進(jìn)行回朔校正；用 對(duì)數(shù)減影法進(jìn)行減影，獲得消除偽影的減影圖像。膠片打印控制部分控制DSA圖像的膠片打印。圖像處理部分對(duì)DICOM圖像進(jìn)行適當(dāng)?shù)奶幚怼＞W(wǎng)絡(luò)傳輸控制對(duì)DSA圖像在PACS系統(tǒng)中的網(wǎng)絡(luò)傳輸進(jìn)行控制。報(bào)告系統(tǒng)自動(dòng)生成DSA圖文檢查報(bào)告。顯示器用于顯示DSA圖像。輸出接口 用于接打印機(jī)。還包括用于與遠(yuǎn)程計(jì)算機(jī)或互聯(lián)網(wǎng)連接的通訊接口。與現(xiàn)有技術(shù)相比本發(fā)明具有以下特點(diǎn)和有益效果這種數(shù)字減影血管造影運(yùn)動(dòng) 偽影消除方法，采用時(shí)空體(Spatio-TemporalVolume)結(jié)構(gòu)來(lái)描述DSA圖像序列；在采 用常規(guī)方法(如Carmy算子)得到的備選控制點(diǎn)的基礎(chǔ)上，對(duì)DSA時(shí)空體進(jìn)行時(shí)空切片 (Spatio-TemporalSlice)得到的各個(gè)控制點(diǎn)的運(yùn)動(dòng)軌跡，選擇運(yùn)動(dòng)軌跡變化最大的若干個(gè) 反映DSA圖像背景結(jié)構(gòu)運(yùn)動(dòng)規(guī)律的控制點(diǎn)；控制點(diǎn)不是僅僅分布在在蒙片或減影盈片上， 而是貫穿在所有的蒙盈片序列上；然后基于這些控制點(diǎn)進(jìn)行圖像配準(zhǔn)和數(shù)字減影以消除運(yùn) 動(dòng)偽影。我們選擇的控制點(diǎn)不是僅僅分布在在蒙片或減影盈片上，而是貫穿在所有的蒙盈 片序列上，更便于分析DSA圖像的時(shí)空運(yùn)動(dòng)特性，也為下一步進(jìn)行時(shí)空卷繞(Warping)打下 基礎(chǔ)。在DSA圖像序列中的每個(gè)圖像上，以所有的控制點(diǎn)作為三角形的頂點(diǎn)，把圖像剖分成 不均勻的Delaunay三角形網(wǎng)格。在采用常規(guī)方法得到的備選控制點(diǎn)的基礎(chǔ)上，在DSA圖像序列中的每個(gè)圖像上， 以所有的控制點(diǎn)作為三角形的頂點(diǎn)，把圖像剖分成不均勻的Delaimay三角形網(wǎng)格；從第 一幅盈片圖像開(kāi)始，結(jié)合前一幅圖像中控制點(diǎn)的位置來(lái)確定本幅圖像的Warping參數(shù)， Warping采用多項(xiàng)式Warping方式；以此類推，直到完成對(duì)最后一幅圖像的Warping，即以時(shí) 空卷繞的方式逐幅圖像地對(duì)盈片圖像進(jìn)行非線性校正，逐步達(dá)到與蒙片的精細(xì)配準(zhǔn)；然后 進(jìn)行數(shù)字減影以消除運(yùn)動(dòng)偽影。本發(fā)明克服了傳統(tǒng)不能充分利用DSA圖像序列中大部分造影劑濃度沒(méi)有達(dá)到峰值的盈片圖像的缺點(diǎn)，可以進(jìn)行更加精細(xì)的圖像配準(zhǔn)。由于采用的控制點(diǎn)不是僅僅分布在在蒙片或減影盈片上，而是貫穿在所有的蒙 盈片序列上，更便于分析DSA圖像的時(shí)空運(yùn)動(dòng)特性，也為下一步進(jìn)行時(shí)空卷繞打下基礎(chǔ)。在 DSA圖像序列中的每個(gè)圖像上，以所有的控制點(diǎn)作為三角形的頂點(diǎn)，把圖像剖分成不均勻的 Delaunay三角形網(wǎng)格，本發(fā)明具有充分利用時(shí)空信息，通過(guò)時(shí)空卷繞處理，逐幅圖像地對(duì)盈 片圖像進(jìn)行非線性校正，逐步提高與蒙片的精細(xì)配準(zhǔn)，數(shù)字減影后有效地去除運(yùn)動(dòng)偽影的 技術(shù)進(jìn)步。本發(fā)明的應(yīng)用能夠有效地去除運(yùn)動(dòng)偽影，改善DSA圖像質(zhì)量，進(jìn)而提高醫(yī)生的診 斷準(zhǔn)確性，提高工作效率。


下面結(jié)合附圖對(duì)本發(fā)明做進(jìn)一步詳細(xì)的說(shuō)明。圖1 減影圖像中出現(xiàn)的偽影示意圖。圖2 :DSA蒙片與盈片的部分序列示意圖。圖3 采用邊緣檢測(cè)算法選取控制點(diǎn)示意圖。圖4 時(shí)空卷繞的示意圖。圖 5 時(shí)空體(Spatio-temporalVolume)的結(jié)構(gòu)示意圖。圖6 =DSA圖像的Delaunay三角剖分示意圖。圖7 圖像卷繞的示意圖。圖8 基于時(shí)空切片的DSA運(yùn)動(dòng)偽影消除方法的流程圖。圖9 本發(fā)明的偽影消除與DSA設(shè)備減影效果對(duì)比圖。圖10 仿射變換示意圖。圖11 本發(fā)明的系統(tǒng)示意圖。
具體實(shí)施例方式消除偽影的主要手段就是要找到背景部分像素的運(yùn)動(dòng)規(guī)律，實(shí)現(xiàn)精確的圖像配準(zhǔn) 并減影。雖然DSA圖像序列中大部分是造影劑濃度沒(méi)有達(dá)到峰值的盈片圖像，如圖1所示。 但這些圖像中只是血管部分沒(méi)有被完全突出出來(lái)(變化中)，其它部分如骨骼、肌肉的灰度 值并沒(méi)有明顯變化(部分像素的位置會(huì)發(fā)生變化)。我們完全可以利用DSA圖像序列中的時(shí) 空相關(guān)信息，研究位于背景部分的控制點(diǎn)(ControlPoints)的運(yùn)動(dòng)規(guī)律，再通過(guò)圖像卷繞實(shí) 現(xiàn)圖像配準(zhǔn)。參見(jiàn)圖8所示，本發(fā)明的具體實(shí)施步驟如下DSA圖像資料由河南大學(xué)附屬醫(yī)院淮 河醫(yī)院(三級(jí)甲等醫(yī)院)提供，采用的DSA設(shè)備型號(hào)為美國(guó)GE公司LCVPlusDSA系統(tǒng)。步驟1、讀入DSA圖像序列(DIC0M格式)并進(jìn)行預(yù)處理。預(yù)處理包括灰度拉伸等 簡(jiǎn)單的處理，如果原始圖像質(zhì)量比較好也可不做預(yù)處理。如圖2所示。步驟2、采用Carmy算子找到圖像的邊緣點(diǎn)。由于沒(méi)必要把所有用Carmy算子檢測(cè)出的點(diǎn)都作為控制點(diǎn)，所以在Carmy邊緣點(diǎn)計(jì)算出該點(diǎn)的梯度值后，根據(jù)閾值來(lái)判斷該點(diǎn) 是否作為備選控制點(diǎn)。本實(shí)例選取約100個(gè)備選控制點(diǎn)。如圖3所示。Canny邊緣檢測(cè)算子是JohnF. Canny于1986年開(kāi)發(fā)出來(lái)的一個(gè)多級(jí)邊緣檢測(cè)算法，Carmy邊緣檢測(cè)算法stepl 用高斯濾波器平滑圖象；step2:用一階偏導(dǎo)的有限差分來(lái)計(jì)算梯度的幅值和方向；step3:對(duì)梯度幅值進(jìn)行非極大值抑制；step4:用雙閾值算法檢 測(cè)和連接邊緣。為了滿足這些要求Carmy使用了變分法，這是一種尋找滿足特定功能的函 數(shù)的方法。最優(yōu)檢測(cè)使用四個(gè)指數(shù)函數(shù)項(xiàng)的和表示，但是它非常近似于高斯函數(shù)的一階導(dǎo) 數(shù)。當(dāng)噪聲點(diǎn)對(duì)邊緣檢測(cè)有較大的影響，效果更好的邊緣檢測(cè)器是高斯_拉普拉斯 (LoG)算子。它把高斯平滑濾波器和拉普拉斯銳化濾波器結(jié)合起來(lái)，先平滑掉噪聲，再進(jìn)行 邊緣檢測(cè)，所以效果更好。步驟3、利用DSA圖像序列構(gòu)建DSA時(shí)空體，將該組DSA 二維圖像以時(shí)間T為軸按 拍攝順序間隔排列成平行六面體空間。步驟4、對(duì)時(shí)空體做若干個(gè)垂直方向的時(shí)空切片(Spatio-Temporal-Slice)。DSA 時(shí)空體中運(yùn)動(dòng)點(diǎn)目標(biāo)下一時(shí)刻的運(yùn)動(dòng)方向可歸納為與時(shí)空切片平行；與時(shí)空切片垂直； 與時(shí)空切片成夾角三種可能性。如圖5所示。當(dāng)目標(biāo)控制點(diǎn)沿時(shí)空切片方向運(yùn)動(dòng)時(shí)，控制點(diǎn)將被多個(gè)時(shí)空切片依次捕獲，呈現(xiàn) 在時(shí)空切片圖像中為一條斜線。當(dāng)目標(biāo)控制點(diǎn)沿水平方向運(yùn)動(dòng)時(shí)，如果控制點(diǎn)運(yùn)動(dòng)速度小 于時(shí)空切片間隔，則控制點(diǎn)將被同一個(gè)時(shí)空切片多次捕獲，在時(shí)空切片圖中表現(xiàn)為一段有 限長(zhǎng)線段。當(dāng)目標(biāo)控制點(diǎn)運(yùn)動(dòng)方向與時(shí)空切片成夾角時(shí)，它可以分解成平行方向和垂直方 向，它的運(yùn)動(dòng)軌跡為上述兩個(gè)方向軌跡的合成。步驟5、將各時(shí)空切片對(duì)應(yīng)的各控制點(diǎn)連線，得到各個(gè)控制點(diǎn)的運(yùn)動(dòng)軌跡。步驟6，選擇運(yùn)動(dòng)軌跡變化最大的10個(gè)反映DSA圖像背景結(jié)構(gòu)運(yùn)動(dòng)規(guī)律的控制點(diǎn)； 本發(fā)明選擇的控制點(diǎn)不是僅僅分布在在蒙片或減影盈片上，而是貫穿在所有的蒙盈片序列 上，更好地記錄了 DSA圖像的時(shí)空運(yùn)動(dòng)特性。如圖7所示。步驟7、用三角剖分的方法(DelaunayTriangulationMethod)。將所有的控制點(diǎn)作 為三角形的頂點(diǎn)，把圖像剖分成不均勻的Delaimay三角形網(wǎng)格；如圖6所示。在DSA圖像序 列中的每個(gè)圖像上，以所有的控制點(diǎn)作為三角形的頂點(diǎn)，把圖像剖分成不均勻的Delaimay 三角形網(wǎng)格。三角剖分是代數(shù)拓?fù)鋵W(xué)里最基本的研究方法，是極為重要的一項(xiàng)預(yù)處理技術(shù)。 通過(guò)Delaimay三角剖分，可以將復(fù)雜的多邊形轉(zhuǎn)化相對(duì)比較簡(jiǎn)單的三角形集合。在DSA偽 影區(qū)域會(huì)出現(xiàn)大量的三角形，這些區(qū)域正是要重點(diǎn)控制的區(qū)域。步驟8、剖分完成后，采用仿射變換進(jìn)行初步幾何校正。如圖10所示。仿射變 換是空間直角坐標(biāo)變換的一種，它是一種二維坐標(biāo)到二維坐標(biāo)之間的線性變換，保持二 維圖形的“平直線”和“平行性”，其可以通過(guò)一系列的原子變換的復(fù)合來(lái)實(shí)現(xiàn)，包括平移 (Translation)、縮放(Scale) JlR(Flip)J^R(Rotation)和剪切(Shear)。幾何上，兩 個(gè)向量空間之間的一個(gè)仿射變換或者仿射映射由一個(gè)線性變換接上一個(gè)平移組成。在有限 維的情況，每個(gè)仿射變換可以由一個(gè)矩陣A和一個(gè)向量b給出，它可以寫(xiě)作A和一個(gè)附加的 列b。一個(gè)仿射變換對(duì)應(yīng)于一個(gè)矩陣和一個(gè)向量的乘法，而仿射變換的復(fù)合對(duì)應(yīng)于普通的矩 陣乘法，只要加入一個(gè)額外的行到矩陣的底下，這一行全部是O除了最右邊是一個(gè)1，而列 向量的底下要加上一個(gè)1。步驟9、用Warping非線性變換法參考前一幅圖像中控制點(diǎn)的變形位置，逐幅確定 圖像非控制點(diǎn)的位移；如圖4所示。從第一幅盈片圖像開(kāi)始，結(jié)合前一幅圖像中控制點(diǎn)的位置來(lái)確定本幅圖像的Warping參數(shù)。Warping采用多項(xiàng)式Warping方式，多項(xiàng)式的項(xiàng)數(shù)與控 制點(diǎn)數(shù)相同。本實(shí)例中控制點(diǎn)數(shù)為10，多項(xiàng)式如下
<formula>formula see original document page 9</formula>
解線性方程組，得系數(shù)后矩陣求逆，然后按控制點(diǎn)逐點(diǎn)Warping。以此類推，直到完成 對(duì)最后一幅圖像的Warping。由于綜合考慮了 DSA時(shí)空體的整體和局部時(shí)空信息，以時(shí)空 Warping的方式能逐幅圖像地對(duì)盈片圖像進(jìn)行非線性校正，逐步達(dá)到與蒙片的精細(xì)配準(zhǔn)。步驟10、本實(shí)例中配準(zhǔn)中所選取相似性測(cè)度為互信息(Mutuallnformation)準(zhǔn) 貝U。本發(fā)明利用時(shí)空信息進(jìn)行高精度DSA配準(zhǔn)，計(jì)算量較大，所以采用Powell優(yōu)化算法來(lái) 優(yōu)化相似度準(zhǔn)則函數(shù)求極值的過(guò)程，使得有效地克服互信息函數(shù)的局部極值，大大地提高 了配準(zhǔn)精度，達(dá)到亞像素級(jí)。Powell優(yōu)化算法的配準(zhǔn)速度快、精度高，魯棒性較好。步驟11、在DSA圖像配準(zhǔn)中，即使亞像素(或子像素)的錯(cuò)誤配準(zhǔn)也會(huì)在減影圖像 上產(chǎn)生大量的偽影。為配合精確的時(shí)空卷繞，選擇較為復(fù)雜的三次樣條插值算法進(jìn)行亞像 素級(jí)的配準(zhǔn)。三次樣條插值(簡(jiǎn)稱Spline插值)是通過(guò)一系列形值點(diǎn)的一條光滑曲線，數(shù)學(xué) 上通過(guò)求解三彎矩方程組得出曲線函數(shù)組的過(guò)程。步驟12、逐幅對(duì)圖像的灰度失真進(jìn)行回朔校正；采用對(duì)數(shù)減影方法進(jìn)行減影，即 在進(jìn)行減影運(yùn)算前先將盈片和蒙片所有點(diǎn)的灰度值取反，分別做對(duì)數(shù)運(yùn)算，獲得的值作為 當(dāng)前點(diǎn)的灰度值，然后盈片減蒙片獲得減影結(jié)果圖像。步驟13，對(duì)數(shù)減影；用對(duì)數(shù)減影法進(jìn)行減影，獲得消除偽影的減影圖像。參見(jiàn)圖9 所示，左側(cè)的兩幅圖像是DSA設(shè)備減影的效果，右側(cè)的兩幅圖像是采用本發(fā)明方法的進(jìn)行 DSA運(yùn)動(dòng)偽影的消除后的效果。參見(jiàn)圖11所示，一種數(shù)字減影血管造影運(yùn)動(dòng)偽影消除系統(tǒng)，包括DSA設(shè)備、實(shí)現(xiàn)運(yùn) 動(dòng)偽影消除的DSA軟件工作站、DSA膠片打印機(jī)、圖像顯示系統(tǒng)、外部存儲(chǔ)設(shè)備、網(wǎng)絡(luò)傳輸系 統(tǒng)，它的DSA軟件工作站包括以下部分輸入接口 用于接收DSA設(shè)備輸出的DICOM圖像。存儲(chǔ)器用于存儲(chǔ)操作系統(tǒng)、圖像處理軟件、圖像分析數(shù)據(jù)和圖像資料。數(shù)字減影部分采用本發(fā)明的方法進(jìn)行數(shù)字減影。中央處理單元用于讀入DSA圖像，并將DSA圖像進(jìn)行預(yù)處理；采用邊緣檢測(cè)算法 找到每幅DSA圖像中血管影像的邊緣和邊緣點(diǎn)，計(jì)算出邊緣點(diǎn)的梯度值，根據(jù)閾值來(lái)判斷 該點(diǎn)是否作為備選控制點(diǎn)，并選取至少10個(gè)備選控制點(diǎn)；用于構(gòu)建DSA時(shí)空體；將該組DSA 二維圖像以時(shí)間T為軸按拍攝順序間隔排列成平行六面體空間；用于進(jìn)行時(shí)空切片；用切 片指令對(duì)DSA時(shí)空體做垂直方向的時(shí)空切片；用于將各時(shí)空切片對(duì)應(yīng)的各控制點(diǎn)連線，得 到各個(gè)控制點(diǎn)的運(yùn)動(dòng)軌跡；用于分析DSA像素的時(shí)空運(yùn)動(dòng)特性；選擇運(yùn)動(dòng)軌跡變化最大的 10個(gè)反映DSA圖像背景結(jié)構(gòu)運(yùn)動(dòng)規(guī)律的控制點(diǎn)；用三角剖分的方法將所有的控制點(diǎn)作為三 角形的頂點(diǎn)，把圖像剖分成不均勻的三角形網(wǎng)格；用于采用仿射變換法進(jìn)行初步幾何校正；用于作時(shí)空卷繞；用Warping非線性變換法參考前一幅圖像中控制點(diǎn)的變形位置，逐幅確定圖像非控制點(diǎn)的位移；采用Powell優(yōu)化算法來(lái)優(yōu)化相似度準(zhǔn)則函數(shù)求極值的過(guò)程；選擇 三次樣條插值算法進(jìn)行亞像素級(jí)的配準(zhǔn)；用于逐幅對(duì)圖像的灰度失真進(jìn)行回朔校正；用對(duì) 數(shù)減影法進(jìn)行減影，獲得消除偽影的減影圖像。膠片打印控制部分控制DSA圖像的膠片打印。圖像處理部分對(duì)DICOM圖像進(jìn)行適當(dāng)?shù)奶幚?。網(wǎng)絡(luò)傳輸控制對(duì)DSA圖像在PACS系統(tǒng)中的網(wǎng)絡(luò)傳輸進(jìn)行控制。報(bào)告系統(tǒng)自動(dòng)生成DSA圖文檢查報(bào)告。顯示器用于顯示DSA圖像。輸出接口 用于接打印機(jī)。還包括用于與遠(yuǎn)程計(jì)算機(jī)或互聯(lián)網(wǎng)連接的通訊接口。參見(jiàn)圖4所示，圖像卷繞是一種非線性的方法，是通過(guò)指定一系列控制點(diǎn)的位 移來(lái)定義空間變換的圖像變形處理。非控制點(diǎn)的位移根據(jù)控制點(diǎn)進(jìn)行插值來(lái)確定。圖 像卷繞一般是在兩維圖像上進(jìn)行，但是在DSA時(shí)空體上的控制點(diǎn)具有空間和時(shí)間上的相 關(guān)性，在對(duì)某一盈片圖像進(jìn)行圖像卷繞時(shí)，如果我們還考慮到其它圖像對(duì)本次圖像卷繞 的影響，就可以使得圖像卷繞的效果更加精細(xì)。因此本發(fā)明中把卷繞推廣為時(shí)空卷繞 (Spatio-Temporalffarping),即在對(duì)一幅盈片進(jìn)行圖像卷繞時(shí)，要結(jié)合前一幅圖像中相關(guān) 點(diǎn)的位置關(guān)系。參見(jiàn)圖5所示，本發(fā)明采用了時(shí)空體(Spatio-TemporalVolume)結(jié)構(gòu)來(lái)描述DSA圖 像序列。時(shí)空體定義了一個(gè)包括χ，γ和T (時(shí)間維度)三維坐標(biāo)軸的平行六面體空間，把空 間信息(XY坐標(biāo))和時(shí)間信息(Τ坐標(biāo))整合在一個(gè)數(shù)據(jù)結(jié)構(gòu)中。采用時(shí)空體的好處是分析 圖像序列特征時(shí)，把整個(gè)圖像序列看作一個(gè)整體，而不是離散的幀，可以充分保證空間和時(shí) 間的連續(xù)性，能更好地表現(xiàn)圖像序列中結(jié)構(gòu)對(duì)象的運(yùn)動(dòng)特征。在DSA圖像時(shí)空體的基礎(chǔ)上，本發(fā)明進(jìn)一步對(duì)時(shí)空體進(jìn)行時(shí)空切片 (Spatio-TemporalSlice)來(lái)分析DSA背景部分像素的運(yùn)動(dòng)規(guī)律。時(shí)空切片是由圖像序列中 相同位置的像素條帶按時(shí)序組成的一幅圖像，或者說(shuō)是從圖像序列的同一個(gè)位置提取出的 一些相關(guān)像素組合而成的一幅二維圖像。再形象一些來(lái)說(shuō)，假如有一個(gè)(X，Y, T)三維DSA時(shí) 空體，其中(X，Y)是圖像維，T是時(shí)間維。我們沿平行于T軸的方向?qū)SA時(shí)空體進(jìn)行切割， 這樣將得到一個(gè)剖面，稱為時(shí)空切片，時(shí)空切片可以看作是包括時(shí)間維的一幅二維圖像。時(shí) 空切片平行于T軸，但相對(duì)于X、Y軸的位置，即切割的方向、位置，可以根據(jù)具體情況來(lái)確 定。引入時(shí)空切片可以方便地獲取DSA圖像序列中背景部分運(yùn)動(dòng)的動(dòng)態(tài)信息，并能減少計(jì) 算復(fù)雜性。為了精確配準(zhǔn)圖像，需要選擇控制點(diǎn)，這一部分采用常規(guī)方法。DSA圖 像序列中的運(yùn)動(dòng)目標(biāo)是在空間域和時(shí)間域上做連續(xù)的運(yùn)動(dòng)，因此它在DSA時(shí)空體 (Spatio-Temporal-Volume)上的軌跡點(diǎn)具有空間和時(shí)間上的相關(guān)性，本發(fā)明利用時(shí)空切片 對(duì)用常規(guī)方法得到的控制點(diǎn)進(jìn)行運(yùn)動(dòng)估計(jì)，即分析運(yùn)動(dòng)像素的時(shí)空相關(guān)特征。通過(guò)對(duì)控制 點(diǎn)運(yùn)動(dòng)軌跡的研究，選擇一批最有代表性的控制點(diǎn)，既可以保證配準(zhǔn)精度，也提高了 DSA圖 像配準(zhǔn)的魯棒性。參見(jiàn)圖6所示，本發(fā)明提出的基于時(shí)空切片的控制點(diǎn)的選擇方法與其它方法的一個(gè)顯著的區(qū)別是，我們選擇的控制點(diǎn)不是僅僅分布在在蒙片或減影盈片上，而是貫穿在所有的蒙盈片序列上，更便于分析DSA圖像的時(shí)空運(yùn)動(dòng)特性，也為下一步進(jìn)行時(shí)空卷繞 (Warping)打下基礎(chǔ)。在DSA圖像序列中的每個(gè)圖像上，以所有的控制點(diǎn)作為三角形的頂 點(diǎn)，把圖像剖分成不均勻的Delaimay三角形網(wǎng)格。 參見(jiàn)圖7所示，通過(guò)時(shí)空卷繞處理，逐幅圖像地對(duì)盈片圖像進(jìn)行非線性校正，逐步 達(dá)到與蒙片的精細(xì)配準(zhǔn)。由于充分利用時(shí)空信息，提高了配準(zhǔn)精度，數(shù)字減影后可以有效地 去除運(yùn)動(dòng)偽影。
權(quán)利要求
一種數(shù)字減影血管造影運(yùn)動(dòng)偽影消除方法，其特征在于步驟如下步驟1、讀??；將數(shù)字減影血管造影儀器與計(jì)算機(jī)連接，在計(jì)算機(jī)中安裝操作系統(tǒng)和圖像處理系統(tǒng)，用圖像處理軟件讀入一組DSA圖像，并將DSA圖像進(jìn)行預(yù)處理；步驟2、選點(diǎn)；采用邊緣檢測(cè)算法找到每幅DSA圖像中血管影像的邊緣和邊緣點(diǎn)，計(jì)算出邊緣點(diǎn)的梯度值，根據(jù)閾值來(lái)判斷該點(diǎn)是否作為備選控制點(diǎn)，并選取至少10個(gè)備選控制點(diǎn)；步驟3、構(gòu)建DSA時(shí)空體；將該組DSA二維圖像以時(shí)間T為軸按拍攝順序間隔排列成平行六面體空間；步驟4、時(shí)空切片；用切片指令對(duì)DSA時(shí)空體做垂直方向的時(shí)空切片；步驟5、連軌跡；用圖像處理軟件將各時(shí)空切片對(duì)應(yīng)的各控制點(diǎn)連線，得到各個(gè)控制點(diǎn)的運(yùn)動(dòng)軌跡；步驟6，分析DSA像素的時(shí)空運(yùn)動(dòng)特性；選擇運(yùn)動(dòng)軌跡變化最大的10個(gè)反映DSA圖像背景結(jié)構(gòu)運(yùn)動(dòng)規(guī)律的控制點(diǎn)；步驟7、三角剖分；用三角剖分的方法將所有的控制點(diǎn)作為三角形的頂點(diǎn)，把圖像剖分成不均勻的三角形網(wǎng)格；步驟8、仿射變換；采用仿射變換法進(jìn)行初步幾何校正；步驟9、時(shí)空卷繞；用Warping非線性變換法參考前一幅圖像中控制點(diǎn)的變形位置，逐幅確定圖像非控制點(diǎn)的位移；步驟10、優(yōu)化；采用Powell優(yōu)化算法來(lái)優(yōu)化相似度準(zhǔn)則函數(shù)求極值的過(guò)程；步驟11、配準(zhǔn)；選擇三次樣條插值算法進(jìn)行亞像素級(jí)的配準(zhǔn)；步驟12、灰度校正；逐幅對(duì)圖像的灰度失真進(jìn)行回朔校正；步驟13，對(duì)數(shù)減影；用對(duì)數(shù)減影法進(jìn)行減影，獲得消除偽影的減影圖像。
2.根據(jù)權(quán)利要求1所述的一種數(shù)字減影血管造影運(yùn)動(dòng)偽影消除方法，其特征在于所 述邊緣檢測(cè)算法是Carmy算子或高斯_拉普拉斯算子。
3.一種數(shù)字減影血管造影運(yùn)動(dòng)偽影消除系統(tǒng)，包括DSA設(shè)備、實(shí)現(xiàn)運(yùn)動(dòng)偽影消除的DSA 軟件工作站、DSA膠片打印機(jī)、圖像顯示系統(tǒng)、外部存儲(chǔ)設(shè)備、網(wǎng)絡(luò)傳輸系統(tǒng)，其特征在于 它的DSA軟件工作站包括以下部分輸入接口 用于接收DSA設(shè)備輸出的DICOM圖像；存儲(chǔ)器用于存儲(chǔ)操作系統(tǒng)、圖像處理軟件、圖像分析數(shù)據(jù)和圖像資料；數(shù)字減影部分采用本發(fā)明的方法進(jìn)行數(shù)字減影；DSA軟件工作站用于讀入DSA圖像，并將DSA圖像進(jìn)行預(yù)處理；采用邊緣檢測(cè)算法找 到每幅DSA圖像中血管影像的邊緣和邊緣點(diǎn)，計(jì)算出邊緣點(diǎn)的梯度值，根據(jù)閾值來(lái)判斷該 點(diǎn)是否作為備選控制點(diǎn)，并選取至少10個(gè)備選控制點(diǎn)；用于構(gòu)建DSA時(shí)空體；將該組DSA 二維圖像以時(shí)間T為軸按拍攝順序間隔排列成平行六面體空間；用于進(jìn)行時(shí)空切片；用切 片指令對(duì)DSA時(shí)空體做垂直方向的時(shí)空切片；用于將各時(shí)空切片對(duì)應(yīng)的各控制點(diǎn)連線，得 到各個(gè)控制點(diǎn)的運(yùn)動(dòng)軌跡；用于分析DSA像素的時(shí)空運(yùn)動(dòng)特性；選擇運(yùn)動(dòng)軌跡變化最大的 10個(gè)反映DSA圖像背景結(jié)構(gòu)運(yùn)動(dòng)規(guī)律的控制點(diǎn)；用三角剖分的方法將所有的控制點(diǎn)作為三 角形的頂點(diǎn)，把圖像剖分成不均勻的三角形網(wǎng)格；用于采用仿射變換法進(jìn)行初步幾何校正； 用于作時(shí)空卷繞；用Warping非線性變換法參考前一幅圖像中控制點(diǎn)的變形位置，逐幅確定圖像非控制點(diǎn)的位移；采用Powell優(yōu)化算法來(lái)優(yōu)化相似度準(zhǔn)則函數(shù)求極值的過(guò)程；選擇 三次樣條插值算法進(jìn)行亞像素級(jí)的配準(zhǔn)；用于逐幅對(duì)圖像的灰度失真進(jìn)行回朔校正；用對(duì) 數(shù)減影法進(jìn)行減影，獲得消除偽影的減影圖像； 膠片打印控制部分控制DSA圖像的膠片打??； 圖像處理部分對(duì)DICOM圖像進(jìn)行適當(dāng)?shù)奶幚恚?網(wǎng)絡(luò)傳輸控制對(duì)DSA圖像在PACS系統(tǒng)中的網(wǎng)絡(luò)傳輸進(jìn)行控制； 報(bào)告系統(tǒng)自動(dòng)生成DSA圖文檢查報(bào)告； 顯示器用于顯示DSA圖像； 輸出接口 用于接打印機(jī)。
4.根據(jù)權(quán)利要求1所述的一種數(shù)字減影血管造影運(yùn)動(dòng)偽影消除系統(tǒng)，其特征在于還 包括用于與遠(yuǎn)程計(jì)算機(jī)或互聯(lián)網(wǎng)連接的通訊接口。
全文摘要
一種數(shù)字減影血管造影運(yùn)動(dòng)偽影消除方法，步驟如下步驟1、讀??；步驟2、選點(diǎn)；步驟3、構(gòu)建DSA時(shí)空體；步驟4、時(shí)空切片；步驟5、連軌跡；步驟6，分析DSA像素的時(shí)空運(yùn)動(dòng)特性；步驟7、三角剖分；步驟8、仿射變換；步驟9、時(shí)空卷繞；步驟10、優(yōu)化；步驟11、配準(zhǔn)；步驟12、灰度校正；步驟13，對(duì)數(shù)減影。本發(fā)明屬于圖像處理技術(shù)。采用時(shí)空分析的方法，使DSA圖像配準(zhǔn)更加精確，有效去除運(yùn)動(dòng)偽影，得到清晰的血管造影圖像。進(jìn)而提高醫(yī)生的診斷準(zhǔn)確性，提高工作效率。
文檔編號(hào)G06T5/00GK101822545SQ201010168648
公開(kāi)日2010年9月8日 申請(qǐng)日期2010年5月11日 優(yōu)先權(quán)日2010年5月11日
發(fā)明者張帆, 張新紅 申請(qǐng)人:河南大學(xué)