專利名稱：一種對擴散張量磁共振圖像進行恢復處理的方法
技術領域：
本發(fā)明涉及醫(yī)療診斷圖像處理技術領域，特別是涉及一種對擴散張量磁共振圖像進行恢復處理的方法。

背景技術：
擴散張量磁共振成像(diffusion tensor imaging，DTI)是新近發(fā)展的一種成像技術，它使得無創(chuàng)傷地研究人類組織結構(例如進行人腦神經(jīng)纖維)和進行多發(fā)性硬化癥以及腦中風等疾病的診斷成為可能。臨床中，由于病人體內熱噪聲等的影響使得DTI圖像受到了賴斯噪聲的污染。為了便于進行張量計算和纖維追蹤等后處理，首先需要對圖像進行恢復處理。
目前，DTI圖像恢復的主流為基于偏微分方程(partial differentialequation，PDE)的方法。現(xiàn)在，采用PDE對DTI圖像恢復均基于歐幾里得梯度?；跉W幾里得梯度的圖像恢復方法雖然具有平滑效果好且處理速度快的優(yōu)點但是卻沒有仿射不變梯度所具有的優(yōu)良特性仿射不變性。
仿射不變性是指在經(jīng)歷了像平移、旋轉或縮放這樣的仿射變換后，其性質(如數(shù)值、算法、收斂等等)仍然保持不變的特性。
現(xiàn)有技術在處理圖像時，為了保留仿射不變的優(yōu)良特性，其處理對象為標量圖像。使用該方法進行DTI向量圖像恢復，需要把向量各元素分別進行處理；該方法的缺陷是由于向量各元素獨立進行演化、因而無法保留元素間的相關性，從而導致原始圖像的很多有用信息丟失。


發(fā)明內容
為了消除或減輕賴斯噪聲的污染、對圖像進行恢復處理，保留仿射不變的優(yōu)良特性，克服現(xiàn)有技術在處理圖像時丟失很多有用信息的缺陷，本發(fā)明提出了以下技術方案。
1.一種對擴散張量磁共振圖像進行恢復處理的方法，使用計算機，包括以下步驟 第一步，獲取研究對象的圖像信息； 所述的獲取，是在磁共振成像系統(tǒng)上進行； 所述的圖像信息，包括未施加梯度脈沖獲得的DTI圖像信息、施加加權梯度脈沖獲得的DTI圖像信息即每層擴散加權的DWI圖像信息； 所述的每層擴散加權DWI圖像信息，是在每層中以6個以上的不共線方向對研究對象進行加權成像從而獲得該層的6組以上擴散加權DWI圖像信息； 第二步，把每層DWI圖像信息轉化為向量圖像信息，其表示式為I＝(I1(x，y)，...，Im(x，y))T； 其中In，n＝1，2，...，m表示第n加權脈沖加權圖像，(x，y)為二維圖像中的橫、縱坐標，T表示矩陣轉置； 第三步，按自然層數(shù)為序或者按另行編排的先后層數(shù)為序進行排列，并令最前層為第1層、最后層為第P層； 令層數(shù)計數(shù)器為1； 設置步長為dt，設置在每一層的循環(huán)次數(shù)為Q； 第四步，令循環(huán)次數(shù)計數(shù)器為1； 第五步，對當前層數(shù)、當前循環(huán)次數(shù)的H和J進行計算； 所述的H和J是兩個相互獨立的仿射不變算子； 上式中，fjk＝[f1jk，...，fmjk]T j，k＝x，y T表示矩陣轉置；當m＝1時，fjk是標量； 對任意尺度的向量，J由下式構造 由上式可得， 第六步，計算仿射不變梯度Gaff和歐幾里得曲率κ； κ＝J |▽f|3 上式中，fi(x，y，0)＝Ii(x，y)(i＝1，...，m) “▽”表示梯度； 第七步，確定向量圖像仿射不變各向異性平滑的數(shù)值 fit＝Gaffκ1/3|▽fi|； 第八步，根據(jù)以下公式計算得到恢復后的圖像信息并存儲該信息； f(x，y，t)＝f(x，y，t-1)+dt×ft-1(t＝1，...，N)； 上式中，f(x，y，t)＝(f1(x，y，t)，...，fm(x，y，t))T，ft＝(f1t，...，fmt)T； 第九步，檢查循環(huán)次數(shù)計數(shù)器的計數(shù)值是否為Q，如果為Q則跳過第十步、直接進行第十一步，如果不為Q則進行下一步； 第十步，當?shù)诰挪綑z查循環(huán)次數(shù)計數(shù)器的計數(shù)值不為Q時，循環(huán)次數(shù)計數(shù)器加1后返回第五步； 第十一步，檢查層數(shù)計數(shù)器的計數(shù)值是否為P，如果為P則跳過第十二步、直接進行第十三步，如果不為P則進行下一步； 第十二步，當?shù)谑徊綑z查層數(shù)計數(shù)器的計數(shù)值不為P時，層數(shù)計數(shù)器加1后返回第四步； 第十三步，存儲圖像信息；進入結束。
2.所述的設置在每一層的循環(huán)次數(shù)為20。
本發(fā)明的有益效果是在消除或減少賴斯噪聲的污染、對圖像進行恢復處理中，不僅平滑效果好，保留了仿射不變的優(yōu)良特性，還具有突出的優(yōu)點即保留更多的有用信息，從而使圖像更清晰、更明了，展示的圖像細節(jié)也更多。
此外，本發(fā)明方法的處理對象可以是任意維數(shù)的向量圖像。



圖1是采集到的原始第一梯度加權DWI圖像； 圖2是對圖1的圖像信息采用標量仿射不變梯度擴散方法進行濾波處理后的圖像； 圖3是對圖1的圖像信息采用本方明方法進行恢復處理后的圖像； 圖4是根據(jù)采集到的原始人腦靠近胼胝體附近一個體元進行纖維追蹤而獲得的白質纖維圖像； 圖5是對圖4采用的體元進行恢復處理后再進行纖維追蹤獲得的白質纖維圖像； 圖6是本發(fā)明方法的流程圖。
下面結合附圖和具體實施方式
對本發(fā)明作進一步的詳細說明。

具體實施例方式 實施本發(fā)明方法，需要依據(jù)本發(fā)明所述的原理等等編寫程序軟件，計算機運行這些程序軟件，能自動進行相關的圖像恢復處理操作。
上述程序軟件可以添加到磁共振成像設備的計算機上，直接由磁共振成像設備上的計算機進行圖像恢復處理工作，并從計算機自帶的顯示屏上觀察圖像效果。
上述程序軟件也可以編寫在其它帶顯示屏的計算機上，然后將需要作恢復處理的圖像信息，從磁共振成像設備的計算機上送至進行處理的計算機上，所述計算機會在程序軟件的安排下進行相關的圖像恢復處理操作，并可在顯示屏上觀察圖像效果。如果采用本方法，則計算機硬件一般應滿足型號奔騰4以上、內存512Mb以上、硬盤20G以上。
以下對本發(fā)明進行相關的介紹、解釋和進一步的補充說明。
為了保持PDE方法的優(yōu)點和仿射不變梯度擴散的優(yōu)良特性又要兼顧DTI圖像的向量特性，本發(fā)明提出了一種新的圖像恢復方法。
在本發(fā)明方法中，處理對象為向量整體，故向量間的聯(lián)系在演化過程中仍然保留，即避免了丟失有用信息。
DTI圖像是通過在傳統(tǒng)磁共振成像系統(tǒng)上施加加權梯度脈沖而獲得的加權圖像。通常需要在N(N大于等于6)個不共線方向對組織的每層進行加權成像從而獲得該層的N個擴散加權圖像(DWI)。
進行擴散加權成像(或者擴散張量成像)的目的是根據(jù)DWI來進行張量場的表面擴張系數(shù)(張量跡)等參數(shù)的計算，以便由這些參數(shù)來了解組織特性或進行疾病診斷等，更特別地，依據(jù)DWI可以進行生物體的纖維追蹤。
受賴斯噪聲影響，DWI圖像質量會下降。由含噪DWI圖像計算所得參數(shù)會偏離實際值，而由其追蹤獲得的白質纖維分辨率降低，纖維長度變短，纖維數(shù)量變少且平滑性變差纖維分叉現(xiàn)象較嚴重。
對圖像進行恢復處理也叫圖像去噪、圖像平滑或者圖像恢復。
對DWI圖像進行恢復即消除賴斯噪聲對DWI圖像的影響，需要把每層N個DWI圖像作為整體來處理。N個圖像對應坐標的像素點組成一個N維向量，故圖像整體為向量圖像。所謂向量圖像也就是元素為向量的圖像。
采用偏微分方程對圖像進行恢復相當于使圖像灰度值進行擴散，該擴散的最快速度即梯度如果具有仿射不變性則稱該擴散為仿射不變擴散。
本發(fā)明方法的步驟如下 第一步，獲取研究對象的圖像信息；所述的獲取，是在磁共振成像系統(tǒng)上進行；所述的圖像信息，包括未施加梯度脈沖獲得的DTI圖像信息、施加加權梯度脈沖獲得的DTI圖像信息即每層擴散加權的DWI圖像信息；所述的每層擴散加權DWI圖像信息，是在每層中以6個以上的不共線方向對研究對象進行加權成像從而獲得該層的6組以上擴散加權DWI圖像信息。
第二步，把每層DWI圖像信息轉化為向量圖像信息，其表示式為I＝(I1(x，y)，....，Im(x，y))T； 其中In，n＝1，2，...，m表示第n加權脈沖加權圖像，(x，y)為二維圖像中的橫、縱坐標，T表示矩陣轉置。
第三步，按自然層數(shù)為序或者按另行編排的先后層數(shù)為序進行排列，并令最前層為第1層、最后層為第P層；令層數(shù)計數(shù)器為1；設置步長為dt，設置在每一層的循環(huán)次數(shù)為Q。
為了有序的進行圖像恢復工作，在本步驟中，需要對圖像的層數(shù)進行排列、編號，排列可以有二種方法進行，第一種是按自然層數(shù)為序，第二種是按另行編排的先后層數(shù)為序。
設置步長為dt、設置在每一層的循環(huán)次數(shù)為Q；該二者的具體數(shù)值，可以根據(jù)不同型號的磁共振成像系統(tǒng)、被檢查物體的性質、研究者希望達到的效果，進行設定并調整，并歸納整理出帶有規(guī)律性的數(shù)據(jù)，以利于今后快速、優(yōu)質地開展圖像恢復工作。
第四步，令循環(huán)次數(shù)計數(shù)器為1。
第五步，對當前層數(shù)、當前循環(huán)次數(shù)的H和J進行計算； 所述的H和J是兩個相互獨立的仿射不變算子； 上式中，fjk＝[f1jk，...，fmjk]T j，k＝x，y T表示矩陣轉置；當m＝1時，fjk是標量； 由上式可得 (上式如果沒有錯誤的話，建議保留) 對任意尺度的向量，J由下式構造 由上式可得， 第六步，計算仿射不變梯度Gaff和歐幾里得曲率κ； κ＝J|▽f|3 上式中，fi(x，y，0)＝Ii(x，y)(i＝1，...，m) “▽”表示梯度。
第七步，確定向量圖像各分量的仿射不變各向異性平滑的數(shù)值 fit＝Gaffκ1/3|▽fi|，公式中，i表示第“i”分量。
第八步，根據(jù)以下公式計算得到恢復后的圖像信息并存儲該信息； f(x，y，t)＝f(x，y，t-1)+dt×ft-1(t＝1，...，N)； 上式中，f(x，y，t)＝(f1(x，y，t)，...，fm(x，y，t))T，ft＝(f1t，...，fmt)T表示t時刻向量圖像的灰度值，ft＝(f1t，...，fmt)T表示向量圖像的仿射不變各向異性平滑的數(shù)值。
由上式可知圖像t時刻的灰度值等于圖像t-1時刻的灰度值加上步長和向量圖像t-1時刻的仿射不變各向異性平滑數(shù)值的乘積。
fit＝Gaffκ1/3|▽fi|。
第九步，檢查循環(huán)次數(shù)計數(shù)器的計數(shù)值是否為Q，如果為Q則跳過第十步、直接進行第十一步；如果不為Q則進行下一步。
在本步驟中，循環(huán)次數(shù)計數(shù)器的計數(shù)值為Q，說明當前層的循環(huán)次數(shù)已滿，已經(jīng)達到設定的循環(huán)次數(shù)。
第十步，第九步檢查循環(huán)次數(shù)計數(shù)器的計數(shù)值不為Q，循環(huán)次數(shù)計數(shù)器加1后返回第五步。
在本步驟中，循環(huán)次數(shù)計數(shù)器的計數(shù)值不為Q，說明當前層需要繼續(xù)進行循環(huán)工作。
第十一步，檢查層數(shù)計數(shù)器的計數(shù)值是否為P，如果為P則跳過第十二步、直接進行第十三步；如果不為P則進行下一步。
在本步驟中，是檢查所有的層數(shù)是否進行過恢復處理，并按照二種不同的結論、進行二種不同的后續(xù)工作。
第十二步，第十一步檢查層數(shù)計數(shù)器的計數(shù)值不為P，層數(shù)計數(shù)器加1后返回第四步。
在本步驟中，層數(shù)計數(shù)器的計數(shù)值不為P，說明還有其他層數(shù)的圖像需要繼續(xù)進行恢復處理。
第十三步，存儲圖像信息；進入結束。
在上述的圖像恢復處理的過程中，迭代步數(shù)，即設置在每一層的循環(huán)次數(shù)Q，是根據(jù)圖像質量來確定的，如果圖像越不平滑，質量越差則迭代步數(shù)取值越大。否則，如果圖像質量高則迭代步數(shù)取值較小。
采用上述的方法步驟對向量圖像進行恢復，通過對向量圖像定義仿射不變因子，使得本平滑模型既具有仿射不變梯度各向異性擴散仿射的優(yōu)點又能在平滑過程中保持了向量圖像各分量之間的相關信息。實驗證明，相對于標量平滑方法，本發(fā)明的恢復方法獲得的圖像更加平滑，并且，恢復圖像追蹤到的人腦神經(jīng)纖維更加平滑且其數(shù)量更多，長度更長。
實施例一 在本實施例中，使用GE公司的Signa 1.5T磁共振系統(tǒng)，先采集人腦DTI圖像，然后進行恢復處理。
圖像是對某健康成人腦部用DT-MRI成像獲得。
在本實施例中，我們對人腦進行了6個加權梯度的加權成像，成像體元為256×256×45，即分45層進行成像，每層圖像的尺寸為256mm×256mm。
圖1是采集到的原始第一梯度加權DWI圖像；圖2是對圖1的圖像信息采用標量仿射不變梯度擴散方法進行濾波處理后的圖像；圖3是對圖1的圖像信息采用本方明方法進行恢復處理后的圖像。
進行三圖比較，在圖3可以看出，本發(fā)明方法可以去除原始圖像的大部分噪聲，因而圖像的清晰度和對比度都得到了很大的提高。相對于標量恢復方法(恢復結果示于圖2)，使用本發(fā)明方法恢復性能更佳。
圖4是根據(jù)采集到的原始人腦靠近胼胝體附近一個體元進行纖維追蹤而獲得的白質纖維圖像；圖5是對圖4采用的體元進行恢復處理后再進行纖維追蹤獲得的白質纖維圖像。
從上述體元中提取了一個64×64×20子體元，并對恢復前后的體元采用纖維追蹤軟件MedINRIA進行了纖維追蹤。由含噪體元追蹤到的纖維示于圖4，根據(jù)本方法恢復圖像追蹤到的人腦白質纖維如圖5所示。從圖5可以看出，本發(fā)明方法濾波后的圖像追蹤到的纖維數(shù)量更多，長度更長，且改善了纖維的交叉和粘連現(xiàn)象。
實驗中，步長dt一般根據(jù)經(jīng)驗取為0.2，考慮到圖像恢復速度和恢復質量的矛盾關系，一般在保證圖像質量的情況下進行盡可能少的循環(huán)。本實施例子中把循環(huán)次數(shù)取為20。因為，大量實驗發(fā)現(xiàn)循環(huán)20次和循環(huán)50次的恢復結果相差不大。
本實施例的步驟如下 第一步，在磁共振成像系統(tǒng)上，獲取研究對象的DWI圖像信息；并在Matlab環(huán)境下讀取采集到的DWI圖像。每層中以6個不共線方向對研究對象進行加權成像、從而獲得該層的6組擴散加權DWI圖像信息。
第二步，把每層DWI圖像信息轉化為向量圖像信息，其表示式為I＝(I1(x，y)，...，Im(x，y))T；其中In，n＝1，2，...，m表示第n加權脈沖加權圖像，(x，y)為二維圖像中的橫、縱坐標，T表示矩陣轉置。
第三步，按自然層數(shù)為序進行排列，圖像信息共45層，令最前層為第1層、最后層為第45層； 令層數(shù)計數(shù)器為1；設置步長dt為0.2、設置在每一層的循環(huán)次數(shù)為20。
第四步，令循環(huán)次數(shù)計數(shù)器為1。
第五步，對當前層數(shù)、當前循環(huán)次數(shù)的H和J進行計算。
第六步，計算仿射不變梯度Gaff和歐幾里得曲率κ。
第七步，確定向量圖像仿射不變各向異性平滑的數(shù)值 fit＝Gaffκ1/3|▽fi| 第八步，根據(jù)以下公式計算得到恢復后的圖像信息并存儲該信息； f(x，y，t)＝f(x，y，t-1)+dt×ft-1(t＝1，...，N)； 上式中，f(x，y，t)＝(f1(x，y，t)，...，fm(x，y，t))T，ft＝(f1t，...，fmt)T。
第九步，檢查循環(huán)次數(shù)計數(shù)器的計數(shù)值是否為Q，如果為Q則跳過第十步、直接進行第十一步；如果不為Q則進行下一步。
第十步，第九步檢查循環(huán)次數(shù)計數(shù)器的計數(shù)值不為Q，循環(huán)次數(shù)計數(shù)器加1后返回第五步。
第十一步，檢查層數(shù)計數(shù)器的計數(shù)值是否為P，如果為P則跳過第十二步、直接進行第十三步；如果不為Q則進行下一步。
第十二步，第十一步檢查層數(shù)計數(shù)器的計數(shù)值不為P，層數(shù)計數(shù)器加1后返回第四步。
第十三步，存儲圖像信息；進入結束。
權利要求
1.一種對擴散張量磁共振圖像進行恢復處理的方法，其特征是，使用計算機，包括以下步驟
第一步，獲取研究對象的圖像信息；
所述的獲取，是在磁共振成像系統(tǒng)上進行；
所述的圖像信息，包括未施加梯度脈沖獲得的DTI圖像信息、施加加權梯度脈沖獲得的DTI圖像信息即每層擴散加權的DWI圖像信息；
所述的每層擴散加權DWI圖像信息，是在每層中以6個以上的不共線方向對研究對象進行加權成像從而獲得該層的6組以上擴散加權DWI圖像信息；
第二步，把每層DWI圖像信息轉化為向量圖像信息，其表示式為I＝(I1(x，y)，...，Im(x，y))T；
其中In，n＝1，2，...，m表示第n加權脈沖加權圖像，(x，y)為二維圖像中的橫、縱坐標，T表示矩陣轉置；
第三步，按自然層數(shù)為序或者按另行編排的先后層數(shù)為序進行排列，并令最前層為第1層、最后層為第P層；
令層數(shù)計數(shù)器為1；
設置步長為dt，設置在每一層的循環(huán)次數(shù)為Q；
第四步，令循環(huán)次數(shù)計數(shù)器為1；
第五步，對當前層數(shù)、當前循環(huán)次數(shù)的H和J進行計算；
所述的H和J是兩個相互獨立的仿射不變算子；
上式中，fjk＝[f1jk，...，fmjk]T j，k＝x，y
T表示矩陣轉置；當m＝1時，fjk是標量；
對任意尺度的向量，J由下式構造
由上式可得，
第六步，計算仿射不變梯度Gaff和歐幾里得曲率κ；
κ＝J |▽f|3
上式中，fi(x，y，0)＝Ii(x，y)(-1，...，m)
“▽”表示梯度；
第七步，確定向量圖像仿射不變各向異性平滑的數(shù)值
fit＝Gaffκ1/3|▽fi|；
第八步，根據(jù)以下公式計算得到恢復后的圖像信息并存儲該信息；
f(x，y，t)＝f(x，y，t-1)+dt×ft-1(t＝1，...，N)；
上式中，f(x，y，t)＝(f1(x，y，t)，...，fm(x，y，t))T，ft＝(f1t，...，fmt)T；
第九步，檢查循環(huán)次數(shù)計數(shù)器的計數(shù)值是否為Q，如果為Q則跳過第十步、直接進行第十一步，如果不為Q則進行下一步；
第十步，當?shù)诰挪綑z查循環(huán)次數(shù)計數(shù)器的計數(shù)值不為Q時，循環(huán)次數(shù)計數(shù)器加1后返回第五步；
第十一步，檢查層數(shù)計數(shù)器的計數(shù)值是否為P，如果為P則跳過第十二步、直接進行第十三步，如果不為P則進行下一步；
第十二步，當?shù)谑徊綑z查層數(shù)計數(shù)器的計數(shù)值不為P時，層數(shù)計數(shù)器加1后返回第四步；
第十三步，存儲圖像信息；進入結束。
2.根據(jù)權利要求1所述的一種對擴散張量磁共振圖像進行恢復處理的方法，其特征是，所述的設置在每一層的循環(huán)次數(shù)為20。
全文摘要
本發(fā)明公開了一種對擴散張量磁共振圖像進行恢復處理的方法，涉及醫(yī)療診斷圖像處理技術領域。為克服現(xiàn)有技術丟失很多有用信息的缺陷，提出技術方案第一步，獲取圖像信息；第二步，轉化為向量圖像信息；第三步，設令圖層編號、層數(shù)計數(shù)器、步長、循環(huán)次數(shù)等；第四步，令循環(huán)計數(shù)器為1；第五步，對當前層數(shù)、當前循環(huán)次數(shù)的H和J進行計算；第六步，計算仿射不變梯度Gaff和歐幾里得曲率κ；第七步，確定向量圖像仿射不變各向異性平滑的數(shù)值；第八步，計算得到恢復后的圖像信息；第九步至第十三步是對循環(huán)次數(shù)和處理層數(shù)進行計算和控制，并決定是否結束。有益效果是對圖像進行恢復處理中，保留更多的有用信息，使圖像更清晰、明了，細節(jié)更多。
文檔編號G01R33/54GK101236241SQ20081003396
公開日2008年8月6日 申請日期2008年2月28日 優(yōu)先權日2008年2月28日
發(fā)明者張相芬, 田蔚風, 宏 葉, 李大志 申請人:上海師范大學