畢業(yè)設(shè)計(jì)--基于matlab實(shí)現(xiàn)的ct重建算法仿真比較研究

1、<p>　　學(xué) 士 學(xué) 位 論 文</p><p>　　THESIS OF BACHELOR</p><p>　?。?009—2013年）</p><p>　　題 目 基于Matlab實(shí)現(xiàn)的CT重建算法仿真比較研究 </p><p>

2、;　　學(xué) 院： 信息工程 系 電子信息工程 </p><p>　　專業(yè)班級(jí)： </p><p>　　學(xué)生姓名： 學(xué)號(hào)： </p><p>　　指導(dǎo)教師：

3、 職稱： </p><p>　　起訖日期： 2013.3.18 —— 2013.6.7 </p><p>　　基于Matlab實(shí)現(xiàn)的CT重建算法仿真比較研究</p><p>　　專 業(yè)： 學(xué) 號(hào)： </p><p>　　學(xué)生

4、姓名： 指導(dǎo)教師： </p><p><b>　　摘要</b></p><p>　　計(jì)算機(jī)斷層成像（即CT（Computerized Tomography））技術(shù)是一項(xiàng)近年來快速發(fā)展多學(xué)科交叉的先進(jìn)技術(shù)，除了醫(yī)學(xué)應(yīng)用，還應(yīng)用在工業(yè)無損檢測(cè)、射電天文學(xué)、精密儀器反演等多個(gè)重要領(lǐng)域。CT圖像重建是CT技術(shù)的核心，重建算法

5、的優(yōu)劣直接關(guān)系到對(duì)檢測(cè)結(jié)果判斷的準(zhǔn)確性。醫(yī)學(xué)CT普遍采用解析成像技術(shù)，這主要因?yàn)獒t(yī)用CT可以采集到完整的投影數(shù)據(jù)。然而在工業(yè)應(yīng)用中，一般得不到完整的投影數(shù)據(jù)，因此要利用不完整的投影數(shù)據(jù)得到較清晰的圖像，必須采用迭代重建法。</p><p>　　針對(duì)以上情況，本課題作了以下兩方面的分析比較。</p><p>　　1.能獲取完整投影數(shù)據(jù)時(shí)，用計(jì)算機(jī)模擬了解析法的圖像重建。發(fā)現(xiàn)圖像重建的速度快，

6、且重建的圖像與原圖的視覺相似度接近1。</p><p>　　2.不能獲取完整投影數(shù)據(jù)時(shí)，用計(jì)算機(jī)模擬了迭代算法的圖像重建。先用ART迭代算法引入松弛參數(shù)，分析了松弛參數(shù)重建圖像的影響。接著用SIRT算法重建圖像，這種算法所得的圖像比ART算法較好。最后在ART算法基礎(chǔ)上采用全差分對(duì)重建圖像進(jìn)行優(yōu)化，圖像質(zhì)量有了明顯提高。</p><p>　　關(guān)鍵字：圖像重建算法；解析重建；迭代重建<

7、/p><p>　　Based on Matlab realize the CT reconstruction algorithm simulation comparative study</p><p><b>　　Abstract</b></p><p>　　Computed tomography (CT (Computerized Tomogr

8、aphy)) technology is a rapid development in recent years, multidisciplinary advanced technology, in addition to medical applications, but also used in industrial nondestructive testing, radio astronomy, precision instrum

9、ents and other important areas of the inversion. CT image reconstruction CT technology is the core of the reconstruction algorithm is directly related to the accuracy of the test results determine. Medical CT imaging com

10、monly used analy</p><p>　　For the above, this project made the following two aspects of the analysis and comparison.</p><p>　　1. can get the full projection data, using a computer simulation of

11、the analytical method of image reconstruction. Found fast image reconstruction and the reconstructed image with the original visual similarity is close to 1.</p><p>　　2.can not get the full projection data,

12、using a computer simulation of the iterative algorithm for image reconstruction. First with the introduction of ART relaxation parameter iterative algorithm analyzes the impact of relaxation parameters of the reconstruct

13、ed image. Followed by SIRT algorithm for reconstruction of images, the resulting image of this algorithm is better than the ART algorithm. Finally ART algorithm based on a fully differential optimized for image reconstru

14、ction, image quality</p><p>　　Keywords:Image reconstruction algorithm；analytic reconstruction；</p><p>　　iterative reconstruction</p><p><b>　　目錄</b></p><p>&

15、lt;b>　　摘 要I</b></p><p>　　AbstractII</p><p><b>　　第一章 緒論1</b></p><p><b>　　1.1 引言1</b></p><p>　　1.2 CT的分類2</p><p>　　1.3

16、CT的組成結(jié)構(gòu)4</p><p>　　1.4 論文的組織結(jié)構(gòu)4</p><p>　　第二章 CT圖像重建的基本原理6</p><p>　　2.1 所需數(shù)據(jù)的掃描獲取6</p><p>　　2.2 圖像重建算法6</p><p>　　第三章 CT圖像解析重建算法9</p><p>　　

17、3.1 反投影重建法（BP）9</p><p>　　3.2 濾波反投影圖像重建法（FBP）11</p><p>　　3.3 Radon變換與傅立葉變換12</p><p>　　3.3.1 Radon變換12</p><p>　　3.3.2 傅立葉變換的定義13</p><p>　　3.4 中心切片定理14&

18、lt;/p><p>　　3.5 平行投影重建算法17</p><p>　　3.5.1 算法原理17</p><p>　　3.5.2 濾波（卷積）反投影18</p><p>　　3.5.3 Radon反變換重建算法18</p><p>　　3.6 解析重建算法的實(shí)現(xiàn)和對(duì)比分析19</p><p&

19、gt;　　第四章 迭代重建算法22</p><p>　　4.1 迭代重建算法的思想原理22</p><p>　　4.1.1 迭代重建的數(shù)學(xué)模型22</p><p>　　4.2 代數(shù)重建法（ART）25</p><p>　　4.3 聯(lián)立迭代重建法（SIRT）26</p><p>　　4.4 全差分迭代算法（TV

20、）27</p><p>　　4.5 迭代重建算法的實(shí)現(xiàn)和對(duì)比分析27</p><p><b>　　第五章 結(jié)論37</b></p><p><b>　　參考文獻(xiàn)38</b></p><p><b>　　附錄40</b></p><p><

21、b>　　謝辭60</b></p><p><b>　　第一章 緒論</b></p><p><b>　　1.1 引言</b></p><p>　　自從X射線發(fā)現(xiàn)后，醫(yī)學(xué)上就開始用它來探測(cè)人體疾病。但是，由于人體內(nèi)有些器官對(duì)X線的吸收差別極小，因此X射線對(duì)那些前后重疊的組織的病變就難以發(fā)現(xiàn)。于是，美國與英

22、國的科學(xué)家開始了尋找一種新的東西來彌補(bǔ)用X線技術(shù)檢查人體病變的不足。</p><p>　　CT(電子計(jì)算機(jī)體層成像)是70年代初放射診斷的項(xiàng)重突破，CT不是X線攝影，而是用X線對(duì)人體掃描（圖1.1所示），取得信息,經(jīng)電子計(jì)算機(jī)處理而獲得的重建圖像它能使傳統(tǒng)的X線檢查難以顯示的器官及其病變顯示成像，而且圖像逼真，解剖關(guān)系明確，從而擴(kuò)大了人體的檢查范圍，大大提高了病變的早期檢出率和診斷準(zhǔn)確率。這種檢查簡(jiǎn)便、安全、無痛

23、苦、無創(chuàng)傷、無危險(xiǎn)，它促進(jìn)了醫(yī)學(xué)影像診斷學(xué)的發(fā)展，CT的研制成功被譽(yù)為自倫琴發(fā)現(xiàn)X射線以后，放射診斷學(xué)上最重要的成就。發(fā)明者亨斯菲爾德和科馬克共同獲得了1979年的諾貝爾獎(jiǎng)金。CT最初只用于頭部檢查，1974年又出現(xiàn)了全身的CT。在短短10余年間，CT已遍及全球，從第一代發(fā)展到第五代。我國各大城市醫(yī)院所使用的CT屬第三代。全身CT可以作頭、胸、腹、骨盆的橫斷掃描，也可作甲狀腺、脊柱、關(guān)節(jié)和軟組織及五官等小部位的區(qū)域掃描。CT最適于查明占

24、位性病變?nèi)缒[瘤、囊腫、增大的淋巴結(jié)、血腫、膿腫和肉芽腫的大小、形態(tài)、目和侵犯范圍，它可以決定某些器官癌腫的分期和是否能進(jìn)行手術(shù)切除。在某些情況下，CT還能區(qū)別病變的病理特性如實(shí)性、囊性、血管性、炎性、鈣性、脂肪等。</p><p>　　CT成像技術(shù)的基本原理是根據(jù)人體中不同組織對(duì)X射線的衰減作用不同，用X線束對(duì)人體某部一定厚度的層面進(jìn)行掃描，由探測(cè)器接收透過該層面的X線，轉(zhuǎn)變?yōu)榭梢姽夂?，由光電轉(zhuǎn)換變?yōu)殡娦盘?hào)，再經(jīng)

25、模擬/數(shù)字轉(zhuǎn)換器（analog/digital converter）轉(zhuǎn)為數(shù)字，輸入計(jì)算機(jī)處理。圖像形成的處理有如對(duì)選定層面分成若干個(gè)體積相同的長(zhǎng)方體，稱之為體素（voxel）。掃描所得信息經(jīng)計(jì)算而獲得每個(gè)體素的X線衰減系數(shù)或吸收系數(shù)，再排列成矩陣，即數(shù)字矩陣（digital matrix），數(shù)字矩陣可存貯于磁盤或光盤中。經(jīng)數(shù)字/模擬轉(zhuǎn)換器（digital/analog converter）把數(shù)字矩陣中的每個(gè)數(shù)字轉(zhuǎn)為由黑到白不等灰度的小方

26、塊，即像素（pixel），并按矩陣排列，即構(gòu)成CT圖像。所以，CT圖像是重建圖像。每個(gè)體素的X線吸收系數(shù)可以通過不同的數(shù)學(xué)方法算出。</p><p><b>　　如圖1.1所示。</b></p><p>　　圖1.1 計(jì)算機(jī)斷層成像示意圖</p><p><b>　　1.2 CT的分類</b></p><

27、;p>　　在投影重建的原理和方法的框架下，根據(jù)投影方式不同產(chǎn)生了不同種類的CT成像術(shù)。</p><p><b>　　⑴透射斷層成像</b></p><p>　　透射斷層成像（Transmission Computed Tomography,TCT，簡(jiǎn)稱CT）系統(tǒng)中，從發(fā)射源射出的X射線穿透物體到達(dá)接收器。射線在通過物體時(shí)被物體吸收一部分，余下部分被接收器接收。由

28、于物體各部分對(duì)射線的吸收不同，所以接收器獲得的射線強(qiáng)度實(shí)際上反映了物體各部分對(duì)射線的吸收情況。由于透射成像術(shù)中采用X射線，所以也叫X-CT。接收器接收到的模擬信號(hào)經(jīng)模數(shù)轉(zhuǎn)換器轉(zhuǎn)換成數(shù)字信號(hào)后，把代表著不同角度下的投影數(shù)據(jù)送給計(jì)算機(jī)，它再運(yùn)用復(fù)雜的數(shù)學(xué)方法重建出物體截面的圖像。在醫(yī)學(xué)上它主要針對(duì)骨質(zhì)成像。X—CT成像原理如圖1.2所示：</p><p>　　圖1.2 X射線CT成像原理</p><

29、;p><b>　?、瓢l(fā)射斷層成像</b></p><p>　　在發(fā)射斷層成像（Emission Computed Tomography,ECT）系統(tǒng)中，發(fā)射源在物體內(nèi)部，一般是將具有放射性的離子注入物體內(nèi)部，從物體外檢測(cè)放射出來的量。通過這種方法可以了解離子在物體內(nèi)運(yùn)動(dòng)情況的分布，從而可以檢測(cè)到與生理相關(guān)的狀況?，F(xiàn)在通常用的ECT主要有兩種：①正電子成像PET（Positron Emi

30、ssion Tomography）；②單光子成像SPECT（Single Positron Emission Computed Tomography）。它們的區(qū)別是，采用單光子輻射器時(shí)，檢測(cè)器所接收的光子計(jì)數(shù)便是投影測(cè)量數(shù)據(jù)；而采用正電子輻射器時(shí)，由于在正電子輻射點(diǎn)的毫米級(jí)范圍內(nèi)，一個(gè)正電子會(huì)俘獲一個(gè)負(fù)電子，并產(chǎn)生沿相反方向運(yùn)動(dòng)的兩個(gè)中子，因而測(cè)量數(shù)據(jù)是在一個(gè)很小時(shí)間間隔內(nèi)由兩個(gè)性質(zhì)相反的檢測(cè)器檢測(cè)到的光子對(duì)的數(shù)目。ECT主要針對(duì)代謝過

31、程的成像。</p><p><b>　?、谴殴舱癯上?lt;/b></p><p>　　磁共振成像（Magnetic Resonance Imaging,MRI）在早期稱為核磁共振（Nuclear Magnetic Resonance,NMR）。它的發(fā)明者斯坦福大學(xué)的F.Bloch和哈佛大學(xué)的E.M.Purcell為此獲得了1952年度的諾貝爾物理獎(jiǎng)。它的工作原理簡(jiǎn)介如下：

32、氫核以及其它具有奇數(shù)個(gè)質(zhì)子和中子的核，包含具有一定磁動(dòng)量和旋量的質(zhì)子。如果把它們放在磁場(chǎng)中，它們會(huì)像陀螺在地球重力場(chǎng)中一樣在磁場(chǎng)中進(jìn)動(dòng)。一般情況下質(zhì)子在磁場(chǎng)中時(shí)隨意排列著的，當(dāng)一個(gè)適當(dāng)強(qiáng)度和頻率的共振場(chǎng)信號(hào)作用于物體時(shí)，質(zhì)子吸收能量并轉(zhuǎn)向與磁場(chǎng)相交的朝向。如果此時(shí)將共振場(chǎng)信號(hào)除去，質(zhì)子吸收的能量釋放并被接收器檢測(cè)到。根據(jù)檢測(cè)到的信號(hào)就可以確定質(zhì)子的密度。通過控制所用的共振場(chǎng)信號(hào)和磁場(chǎng)強(qiáng)度，可每次檢測(cè)到沿透過目標(biāo)中一條線的信號(hào)。換句話說，

33、檢測(cè)到的信號(hào)是MRI信號(hào)沿直線的積分。所以檢測(cè)目標(biāo)的工作成了投影重建的問題。由這些信號(hào)可以得到極好的器官解剖圖，在細(xì)節(jié)上遠(yuǎn)遠(yuǎn)超過由X射線CT產(chǎn)生的圖像。由于它能夠很好的保留豐富的軟組織圖像，所以在腦功能研究上，MRI勝過ECT。</p><p><b>　?、瘸晹鄬映上?lt;/b></p><p>　　超聲斷層成像（Ultrasound CT,UCT）與上述的CT不同的

34、地方是超聲射線在人體組織中不一定按直線傳播，尤其是在軟組織和硬組織的分界面會(huì)產(chǎn)生折射效應(yīng)。因此超聲成像主要用于人體中沒有硬組織（如骨骼）部分的檢測(cè)，特別是婦女乳房癌的檢測(cè)。</p><p><b>　?、煞瓷鋽鄬映上?lt;/b></p><p>　　反射斷層成像（Reflection CT,RCT）也是利用投影重建的原理工作的。常見的一個(gè)例子是雷達(dá)系統(tǒng)，其中的雷達(dá)圖是物體

35、反射的回波所產(chǎn)生的。例如在前視雷達(dá)（FLR）中，雷達(dá)發(fā)射器從空中向地面發(fā)射無線電波；雷達(dá)接收器在特定角度所接收到的回波強(qiáng)度是地面發(fā)射量的一個(gè)掃描段的積分。</p><p>　　通過發(fā)展硬件設(shè)備和改進(jìn)軟件算法，可以不斷降低噪聲和提高圖像的清晰度。CT技術(shù)在臨床的許多領(lǐng)域得到愈來愈廣泛的應(yīng)用。應(yīng)該指出，上面介紹的各類CT的功能是相互補(bǔ)充而不是相互替代的。</p><p>　　1.3 CT的組成

36、結(jié)構(gòu)</p><p>　　以上幾種成像技術(shù)盡管使用的物理系統(tǒng)各不相同，但所有系統(tǒng)所測(cè)量到的數(shù)據(jù)均具有物體某種感興趣的物理特性在空間分布的積分形式。成像技術(shù)廣義上講大致可分為三個(gè)部分：信息的收集、重建以及處理和傳遞。</p><p>　　⑴數(shù)據(jù)采集，即獲得攜帶有物體內(nèi)部特征信息的數(shù)據(jù)。所采用的物理測(cè)量系統(tǒng)包括放射源、檢測(cè)器及配套的機(jī)械、電路設(shè)備。為進(jìn)行計(jì)算機(jī)處理，系統(tǒng)還需將測(cè)量得到的模擬量轉(zhuǎn)

37、化成數(shù)字量。</p><p>　?、茍D像重建，即利用各種反演算法從所測(cè)數(shù)據(jù)重建出物體內(nèi)部某種物理量的分布圖。這一步是成像技術(shù)的核心。</p><p>　?、菆D像后處理，包括圖像去噪、銳化，特征提取以及圖像的存儲(chǔ)、壓縮和傳輸?shù)取?lt;/p><p>　　1.4 論文的組織結(jié)構(gòu)</p><p>　　論文主要研究輻射成像技術(shù)中的透射斷層成像理論。講述了

38、CT重建的基本原理</p><p>　　第一章緒論部分，首先介紹了計(jì)算機(jī)斷層成像術(shù)（CT）的發(fā)展歷史，及其醫(yī)學(xué)意義根據(jù)射線產(chǎn)生的物理方式對(duì)CT具體的分類。之后詳細(xì)介紹了透射斷層成像的研究現(xiàn)狀，進(jìn)而對(duì)于目前投影重建中存在的問題進(jìn)行了闡述，從而引出了論文的工作目標(biāo)和內(nèi)容。</p><p>　　第二章詳細(xì)介紹了CT圖像重建系統(tǒng)的基本原理，根據(jù)Beer定理對(duì)斷層進(jìn)行體素劃分，通過聯(lián)立方程組求解方程

39、組的方式重建圖像。但當(dāng)斷層被劃分為更小更多的體素時(shí)，即使采用今天最先進(jìn)的計(jì)算機(jī)技術(shù)解這一組方程組也不是件容易的事。所以這個(gè)方法在后來實(shí)際的圖像重建中已不再使用，從而引出了現(xiàn)在CT圖像重建中常用的方法解析法。</p><p>　　第三章詳細(xì)介紹了現(xiàn)在CT常用的圖像重建算法解析法的算法原理，研究進(jìn)展。先講述了直接反投影重建算法，但由于直接反投影法存在圖像偽跡，從而提出了濾波反投影重建算法。在運(yùn)算的過程中根據(jù)卷積定理從

40、而提出了卷積反投影重建算法，最后對(duì)上述三種解析算法重建的圖像進(jìn)行了對(duì)比分析。</p><p>　　第四章迭代重建算法部分，上述解析算法是目前醫(yī)用CT中最常用的重建算法。但它并不適用于所有投影重建場(chǎng)合。在實(shí)際應(yīng)用中，有時(shí)無法測(cè)得大量的投影，有時(shí)投影不是均勻的分布在180度或360度范圍內(nèi)。例如，為了避免物體（如臟器）運(yùn)動(dòng)，或?yàn)榱藴p小劑量，投影數(shù)據(jù)采集不足；又如某物體在一個(gè)方向上尺寸特別長(zhǎng)，則在這一方向上的投影數(shù)據(jù)無

41、法或較難測(cè)得。再如用成像法勘探地球資源時(shí)，用孔穴法采集的數(shù)據(jù)，數(shù)量不足又不均勻等上述濾波反投影法無能為力，于是提出了圖像重建的迭代重建算法。本章介紹了迭代算法的成像模型及其成像原理。最后演示了一下迭代算法重建圖像的效果。</p><p>　　第五章總結(jié)部分，本章對(duì)論文中講述的CT重建算法進(jìn)行了認(rèn)真的總結(jié)和比較。探討了尚需進(jìn)一步解決的問題及今后有待努力的方向。</p><p>　　第二章 C

42、T圖像重建的基本原理</p><p>　　2.1 所需數(shù)據(jù)的掃描獲取</p><p>　　假設(shè)采用單射線源平移—旋轉(zhuǎn)方式來獲取圖像重建所需的數(shù)據(jù)，該方法需要一個(gè)射線發(fā)生器和一個(gè)射線探測(cè)器。用X射線發(fā)生器發(fā)出單色細(xì)束X射線，對(duì)成像對(duì)象進(jìn)行投射，如圖2.1所示。設(shè)X線源發(fā)出的X射線強(qiáng)度是，探測(cè)器接收到的X射線強(qiáng)度是一次掃射完成后將X射線與探測(cè)器平移一個(gè)小的步長(zhǎng)，再次對(duì)成像對(duì)象進(jìn)行投射，探測(cè)器又

43、接收到一個(gè)新的透射強(qiáng)度，在該投射角度上反復(fù)小步長(zhǎng)平移投射并記錄檢測(cè)的數(shù)據(jù)后，將X射線源與探測(cè)器繞被測(cè)對(duì)象旋轉(zhuǎn)一個(gè)小的角度后在進(jìn)行平行投射并平移，如此重復(fù)多次，直到旋轉(zhuǎn)180度為止。取得這些不同投射位置的透射強(qiáng)度后，就可以用這些數(shù)據(jù)，采用一定的算法重建出X射線源環(huán)繞的斷層的圖像。</p><p>　　圖2.1 單射線源平移旋轉(zhuǎn)掃描示意圖</p><p>　　2.2 圖像重建算法</p&

44、gt;<p>　　CT與普通的X線攝影術(shù)之間有著非常重要的區(qū)別。在CT技術(shù)中，組織對(duì)X線的局部衰減特性被用于離散成像；而在通常的X線攝影術(shù)中，這種衰減信息則重疊在X線底片上。組織對(duì)X線的這種局部衰減特性是X線與被檢測(cè)物體之間的若干相互作用過程的產(chǎn)物。例如光電吸收過程和康普頓散射過程。這些過程中的每一種都有它自己的發(fā)生幾率。幾率也是輻射能量的函數(shù),因?yàn)閺腦線管產(chǎn)生的X線由全部能譜組成，具有“線衰減系數(shù)”的某種組織的衰減性質(zhì)是

45、一個(gè)復(fù)雜的函數(shù)，隨著輻射情況的變化，它可以表現(xiàn)出不同的數(shù)值。</p><p>　　對(duì)于能量為的單能射線，假設(shè)其投射的物體是均勻且為同一介質(zhì)，被投射物體的線性衰減系數(shù)為，厚度為，則在入射線強(qiáng)度為，出射X射線強(qiáng)度為時(shí)有</p><p>　　即（2.1）</p><p>　　這就是著名的Lambert-Beer定律。</p><p>　　

46、實(shí)際被投射的物體是非均勻的，如果在X射線束通過的路徑上，介質(zhì)是非均勻的，則可以將X射線穿過的介質(zhì)沿掃描路徑人為劃分為大小相等的個(gè)小方塊，每個(gè)小塊的厚度為，塊內(nèi)視為均勻介質(zhì)。這樣，每個(gè)小塊有一個(gè)均勻的線性衰減系數(shù)，如圖所示，設(shè)各個(gè)方塊的線性衰減系數(shù)分別為。X射線通過每一個(gè)方塊的衰減為，射線通過第二個(gè)方塊的衰減為，通過第個(gè)方塊的衰減為。將上述各式代入整理，有：</p><p>　　圖2.2 X射線通過非均勻介質(zhì)時(shí)的體

47、素劃分</p><p><b>　?。?.2）</b></p><p>　　即（2.3）</p><p>　　式中，是通過射線強(qiáng)度檢測(cè)得到的檢測(cè)值，是一個(gè)已知量，稱為投影。這樣，上式中，等號(hào)右邊可視為一個(gè)常數(shù)，等待求解的是等式左邊的每一個(gè)方塊的衰減系數(shù)。由于未知量太多，這個(gè)方程沒有唯一解?？紤]到圍繞受檢體的一次掃描可以得到

48、一個(gè)有關(guān)衰減系數(shù)的方程，如果將整個(gè)感興趣斷層都劃分為大小相等的小方塊，每個(gè)小方塊由于一個(gè)固定的衰減系數(shù)，這樣，上述單束射線環(huán)繞掃描得到的每個(gè)穿透強(qiáng)度值都可以形成一個(gè)關(guān)于衰減系數(shù)的方程。只要掃描線不重復(fù)且足夠多，這些衰減系數(shù)的確定值肯定可以通過解線性方程組的方法解出唯一實(shí)數(shù)解。把每個(gè)方塊的衰減系數(shù)值用灰度表示，就可以重建出以衰減系數(shù)為特征的斷層圖像。這種圖像重建方法稱為方程組法。</p><p>　　值得說明的是，

49、當(dāng)上圖中時(shí)，投影就是</p><p><b>　?。?.4）</b></p><p>　　式中，衰減系數(shù)是隨路徑連續(xù)變化的函數(shù)。上式說明，入射X射線強(qiáng)度與出射X射線強(qiáng)度之比經(jīng)對(duì)數(shù)運(yùn)算后，表示沿射線路徑上衰減系數(shù)的線積分，而投影與射線穿越介質(zhì)的路徑長(zhǎng)度成正比。</p><p>　　然而，一般二維斷層圖像至少應(yīng)劃分出160×160個(gè)像素，

50、如要用方程組法重建這個(gè)斷層圖像，就要25600個(gè)獨(dú)立的方程聯(lián)立求解，這是一個(gè)相當(dāng)費(fèi)時(shí)的任務(wù)。盡管1967年世界首臺(tái)CT試驗(yàn)機(jī)采用的就是這個(gè)圖像重建方法，但當(dāng)斷層被劃為更小更多的體素時(shí)，即使采用今天最先進(jìn)的計(jì)算機(jī)技術(shù)，解一組這樣的方程組也不是件容易的事。所以，這種方法在后來的實(shí)際應(yīng)用中已不再使用，而是采用下一章將要介紹的反投影解析法。</p><p>　　第三章 CT圖像解析重建算法</p><

51、p>　　通過上一章已經(jīng)知道CT圖像重建的基本原理，但由于使用最簡(jiǎn)單直接的方程組法不易求得斷層的各個(gè)像素的衰減系數(shù)。這一章介紹現(xiàn)在CT最常用的圖像重建算法——反投影解析法。主要介紹了由投影重建圖像的系統(tǒng)模型，引進(jìn)Radon變換和傅立葉變換，它們對(duì)應(yīng)于X 射線CT的二維重建問題；且給出了Radon變換在圖像重建中的具體形式，即截面函數(shù)沿著直線的線積分等于其 Radon變換，而圖像重建過程即是將截面函數(shù)沿許多不同角度下直線的線積分所產(chǎn)生

52、的投影數(shù)據(jù)進(jìn)行逆Radon變換從而得到截面函數(shù)。另外由逆Radon變換的推導(dǎo)給出了兩種常用的圖像解析重建法：濾波反投影法和卷積反投影法。</p><p>　　3.1 反投影重建法（BP）</p><p>　　反投影是一個(gè)求和的過程，它把所有角度上的數(shù)據(jù)都加在一起。它把投影域中的數(shù)據(jù)沿著原來的投影路徑“涂抹”回去，但不改變數(shù)據(jù)的值。</p><p>　　直接反投影重建

53、是指，斷層平面中某像素的線性衰減系數(shù)可由平面內(nèi)所有經(jīng)過該像素射線的投影和得出，而整幅圖像重建可以看做為所有方向上掃描投影的累加。通俗地講，所謂直接反投影，就是在對(duì)某體面一個(gè)方向的掃描完成后，以得到的投影為灰度值沿著掃描路徑經(jīng)過的體素回抹到體素對(duì)應(yīng)的像素上。改變方向后的多次掃描形成多次回抹，同一像素上多次回抹的灰度累加即完成圖像重建。</p><p>　　圖3.1 斷層的像素和射線模型</p><

54、;p>　　圖3.2 反投影重建示例</p><p>　　如圖3.2所示，其中（a）為原圖像像素，（b）為反投影重建后圖像，（c）為（b）中像素值除以投影線數(shù)后得到的平均值。有上述例子可以看出原圖中像素值為零的像素，經(jīng)反投影重建后不再為零。即有偽跡。</p><p>　　直接反投影圖像產(chǎn)生星狀偽跡的原因在于，該方法是把取自有限物體空間的投影均勻地回抹（反投影）到了射線所及的無限空間的各

55、個(gè)像素上，包括原來像素值為0的點(diǎn)。</p><p>　　圖3.3 反投影重建的等效系統(tǒng)</p><p>　　設(shè)為處于處的一個(gè)點(diǎn)源，即二維斷面上的一個(gè)沖擊函數(shù)，這時(shí)的輸出就是系統(tǒng)的沖擊響應(yīng)。</p><p>　　已知（3.1）</p><p>　　可見，反投影重建過程后得到的圖像不是原來的圖像，系統(tǒng)存在星狀偽跡。</p

56、><p>　　根據(jù)信號(hào)與系統(tǒng)理論，對(duì)于圖3.3有：</p><p><b>　?。?.2）</b></p><p>　　式中，兩個(gè)星號(hào)表示二維卷積。</p><p>　　如果想要去除星狀偽跡，可在系統(tǒng)的輸出端再加上一個(gè)濾波器，設(shè)其時(shí)域函數(shù)為，傳遞函數(shù)為，為了使該加了濾波器的系統(tǒng)輸出的圖像等于系統(tǒng)輸入的原圖，也就是要求<

57、/p><p><b>　　（3.4）</b></p><p><b>　　即</b></p><p><b>　?。?.5）</b></p><p>　　對(duì)該式取二維傅立葉變換，得到：</p><p><b>　?。?.6）</b>&l

58、t;/p><p><b>　　即</b></p><p><b>　　（3.7）</b></p><p>　　這是一個(gè)二維濾波器。即使加上某種近似，這種濾波器也不容易實(shí)現(xiàn)。如果的值取無窮，則是一個(gè)不可能實(shí)現(xiàn)的理想濾波器。為此提出了濾波反投影圖像重建法。</p><p>　　3.2 濾波反投影圖像重建法（

59、FBP）</p><p>　　為了消除直接反投影法產(chǎn)生的星狀偽跡，提出了濾波反投影重建的方法，這種方法是在直接反投影重建方法的基礎(chǔ)上增加一個(gè)濾波器。增加濾波器的思路有兩種，如圖3.4所示。</p><p>　　圖3.4 兩種不同的去除星狀偽跡反投影圖像重建思路</p><p>　　所謂濾波反投影重建法是把獲得的投影函數(shù)先做卷積處理，即人為設(shè)計(jì)一種一維濾波器，用它對(duì)所

60、得的投影函數(shù)進(jìn)行改造，然后把這些改造過的投影函數(shù)進(jìn)行反投影和累加等處理，就可以達(dá)到在一定程度上消除星狀偽跡的目的。</p><p>　　在介紹濾波反投影重建法之前先補(bǔ)充濾波反投影重建法要用到的數(shù)學(xué)知識(shí)和投影定理——中心切片定理。</p><p>　　3.3 Radon變換與傅立葉變換</p><p>　　為了有效和快速地對(duì)圖像進(jìn)行處理，常需要將原定義在圖像空間的圖像

61、以某種形式轉(zhuǎn)換到另外一些空間，并利用在這些空間的特有性質(zhì)方便地進(jìn)行一定的加工，最后再轉(zhuǎn)換回圖像空間以得到所需的效果。通常</p><p>　　正變換：圖像空間到其他空間</p><p>　　反變換：其他空間到圖像空間</p><p>　　3.3.1 Radon變換</p><p>　　Radon變換1917年由奧地利數(shù)學(xué)家Radon提出，是C

62、T圖像重建思想的數(shù)學(xué)表達(dá)。Radon從理論上證明了某物理量的二維分布函數(shù)可由該函數(shù)在其定義域內(nèi)的所有方向上的線積分完全確定。Radon變換的意義在于，只要知道了一個(gè)未知二維分布函數(shù)的所有方向上的線積分，那么就能夠求得該二維分布函數(shù)。所謂CT斷層的圖像重建，就是求取能夠反應(yīng)斷層內(nèi)部結(jié)構(gòu)和組成的某物理量的二維分布。斷層掃描獲得的投影實(shí)際上可視為被測(cè)物體斷層內(nèi)部結(jié)構(gòu)和組成的不同方向上的線積分，所以Radon變換的正變換和逆變換正好對(duì)應(yīng)了CT技

63、術(shù)的射線投影和圖像二維分布函數(shù)的重建過程。</p><p>　　對(duì)的Radon變換定義為沿由和定義的直線的線積分。若設(shè)表示物體的一個(gè)二維斷層分布，通過Radon變換的到沿不同方向直線的線積分相當(dāng)于得到物體不同方向的投影。</p><p>　　圖3.5 定義Radon變換的坐標(biāo)系統(tǒng)</p><p><b>　　上述線積分可寫為：</b></

64、p><p><b>　?。?.8）</b></p><p>　　如果借助Delta函數(shù)，上述線積分還可寫為：</p><p><b>　　（3.9）</b></p><p>　　上式通常稱為的Radon變換。</p><p>　　由上述可知，關(guān)于某直線的Radon變換就是沿該直線

65、上的一維投影。</p><p>　　實(shí)際上當(dāng)將也視作可變參數(shù)時(shí)，就是一個(gè)二元函數(shù)，此時(shí)Radon變換是空間域到變換域的映射。域中的每個(gè)點(diǎn)對(duì)應(yīng)于空間域中的一條直線，這里實(shí)際上就是《圖像分析》中直線的Hough變換。</p><p>　　3.3.2 傅立葉變換的定義</p><p><b>　　（3.10）</b></p><p

66、><b>　　其中</b></p><p>　?。▽?shí)自變量的復(fù)值函數(shù)）（3.11）</p><p>　　稱為的Fourier變換，記為。</p><p>　　稱為的Fourier逆變換，記為。</p><p>　　3.4 中心切片定理</p><p>　　中心切片定理的含義是平行投影的一

67、維傅立葉變換等同于原始物體的二維傅立葉變換的一個(gè)切片。將已知投影數(shù)據(jù)通過一個(gè)簡(jiǎn)單的二維傅立葉反變換可以得到物體截面的一個(gè)評(píng)估。</p><p>　　圖3.6 中心切片定理示意圖</p><p>　　通過上圖可知：某斷層（或它對(duì)應(yīng)的圖像）在視角為時(shí)得到的平行投影（函數(shù)）的一維傅立葉變換，等于二維傅立葉變換過原點(diǎn)的一個(gè)垂直切片，且切片與軸相交成角。</p><p>&l

68、t;b>　　數(shù)學(xué)表達(dá)如下：</b></p><p>　　首先定義代表截面的函數(shù)的二維傅里葉變換：</p><p><b>　?。?.12）</b></p><p>　　類似的定義某個(gè)角度下的一條投影的傅里葉變換：</p><p><b>　　（3.13）</b></p>

69、<p>　　舉個(gè)最簡(jiǎn)單的例子：沿著=0 的直線方向，物體投影的傅里葉變換等于頻域里面= 0的情形</p><p><b>　?。?.13）</b></p><p>　　因?yàn)橄辔灰蜃硬辉僖蕾囉?，在此可以將積分分成兩部分，</p><p><b>　　（3.14）</b></p><p>

70、　　從平行投影的定義看，可以將上式括號(hào)中的部分看成是沿著常量的積分</p><p><b>　?。?.15）</b></p><p>　　將其代入(3.14)得到</p><p><b>　　（3.16）</b></p><p>　　等式右邊表示的是投影的一維傅里葉變換；因此可以得到垂直方向投影和投

71、影函數(shù)的二維變換的關(guān)系式：</p><p><b>　?。?.17）</b></p><p>　　這是傅里葉切片定理最簡(jiǎn)單的形式。很明顯物體與坐標(biāo)系統(tǒng)間的位置方向是獨(dú)立的。例如，在圖3.6中，如果,坐標(biāo)系統(tǒng)為,坐標(biāo)系統(tǒng)旋轉(zhuǎn)，式(3.16)定義的投影的傅里葉變換等于物體沿著一旋轉(zhuǎn)了角的直線的二維傅里葉變換。從而引出了成像系統(tǒng)中著名的中心（傅里葉）切片定理。</p&

72、gt;<p>　　即某圖像在視角時(shí)的平行投影的傅里葉變換給出了圖像二維傅里葉變換在與軸夾角為的一個(gè)切片，此切片通過原點(diǎn)。換句話說，的傅里葉變換在數(shù)值上對(duì)應(yīng)于沿著圖中所示直線BB的值。</p><p>　　通過將坐標(biāo)系統(tǒng)旋轉(zhuǎn)得到坐標(biāo)系統(tǒng)后，傅里葉切片定理可以得到更廣泛的應(yīng)用。在系統(tǒng)中沿著常量為的直線進(jìn)行的投影可表示成：</p><p><b>　　（3.18）<

73、/b></p><p><b>　　其傅里葉變換</b></p><p><b>　　（3.19）</b></p><p>　　將以上兩式合并，可以發(fā)現(xiàn)</p><p><b>　?。?.20）</b></p><p>　　借助旋轉(zhuǎn)坐標(biāo)定義可以得到

74、關(guān)于坐標(biāo)系統(tǒng)的變換形式 </p><p><b>　?。?.21）</b></p><p>　　等式右邊表示一個(gè)空間域的傅里葉變換，則</p><p><b>　?。?.22）</b></p><p>　　這個(gè)等式是直接射線斷層成像的基礎(chǔ)，從而也證明了中心（傅里葉）切片定理。</p>

75、<p>　　上面的結(jié)果顯示：從取圖像分別在與軸成角度的直線上的一系列投影，并對(duì)這些投影逐個(gè)進(jìn)行傅里葉變換，可以得到圖像的二維傅里葉變換在與軸成角的直線上的值。理論上，當(dāng)投影線束無限多時(shí)，就能求得在頻域上相應(yīng)線束的無限多條直線的值，從而可以利用傅里葉逆變換得到原圖。傅里葉逆變換可以采用直角坐標(biāo)系公式也可以采用極坐標(biāo)公式：</p><p><b>　?。?.23）</b></p&

76、gt;<p><b>　?。?.24）</b></p><p><b>　　其中</b></p><p><b>　?。?.25） </b></p><p>　　但從多個(gè)一維傅立葉函數(shù)擬合二維傅立葉函數(shù)需要足夠精確的傅立葉空間插值。由于理想的插值方法還在探索中，這阻礙了基于中心切片定理的

77、二維傅立葉逆變換圖像重建法的發(fā)展。</p><p>　　3.5 平行投影重建算法</p><p>　　3.5.1 算法原理</p><p>　　傅里葉切片定理僅僅提供了一種斷層成像的簡(jiǎn)單概念模型，實(shí)際應(yīng)用中需要一種不同的算法。</p><p>　　濾波反投影算法是目前廣泛應(yīng)用于所有直線透射斷層成像的算法。此算法源于中心切片定理，不僅極其準(zhǔn)確而

78、且保證了快速執(zhí)行的穩(wěn)定性。它是通過轉(zhuǎn)換極坐標(biāo)中傅里葉逆變換和重新確定積分限來實(shí)現(xiàn)的。</p><p>　　首先介紹平行束投影的反投影算法。先回顧一下傅里葉反變換的公式：</p><p><b>　?。?.26）</b></p><p>　　將頻域內(nèi)的直角坐標(biāo)系換成極坐標(biāo)系，有如下替換：</p><p><b>

79、　?。?.27）</b></p><p><b>　　改變微分項(xiàng)</b></p><p><b>　?。?.28）</b></p><p>　　則傅里葉反變換的極坐標(biāo)形式：</p><p><b>　?。?.29）</b></p><p>　

80、　此積分可以分成兩部分，從和,</p><p><b>　?。?.30)</b></p><p>　　利用特性，上式可以寫成,</p><p><b>　?。?.31）</b></p><p>　　此處。利用(3.21)式用視角處的投影的傅里葉變換代替二維傅里葉變換，可以得到：</p>

81、<p><b>　?。?.32）</b></p><p>　　上式是濾波反投影重建的圖像表達(dá)式。對(duì)上式進(jìn)行不同的數(shù)學(xué)變形，將導(dǎo)致不同的物理解釋，相應(yīng)于不同的重建算法。由此得到的重建算法有濾波（卷積）反投影法(F（C）BP)和Radon反變換重建法。</p><p>　　3.5.2 濾波（卷積）反投影</p><p>　　有卷積定理：

82、函數(shù)卷積的傅立葉變換是函數(shù)傅立葉變換的乘積，即：一個(gè)域中的卷積相當(dāng)于另一個(gè)域中的乘積。例如：時(shí)域中的卷積等于頻域中的乘積。</p><p>　　設(shè)在某一旋轉(zhuǎn)角時(shí)，采的投影，濾波函數(shù)為，則濾波后的投影為</p><p><b>　?。?.33）</b></p><p>　　由于數(shù)據(jù)在計(jì)算機(jī)實(shí)現(xiàn)時(shí)是離散的，故應(yīng)進(jìn)行離散卷積。</p>

83、<p>　　反投影是濾波反投影算法的另一個(gè)核心，其實(shí)就是前面介紹的直接反投影算法，只不過是在反投影前將投影數(shù)據(jù)進(jìn)行了卷積濾波預(yù)處理。由此可得極坐標(biāo)下卷積反投影重建圖像如下：</p><p><b>　?。?.34）</b></p><p>　　3.5.3 Radon反變換重建算法</p><p>　　在CT成像原理中已經(jīng)介紹過，對(duì)圖像

84、某一方向的投影表示沿射線路徑方向上衰減系數(shù)的積分，而對(duì)圖像進(jìn)行某一視角上的Radon變換求得的的就是圖像該方向上的投影，所以對(duì)投影進(jìn)行Radon反變換就能實(shí)現(xiàn)圖像重建。對(duì)濾波后的圖像投影數(shù)據(jù)進(jìn)行反Radon變換得到重建圖像的方法就是Radon反變換重建法。</p><p>　　3.6 解析重建算法的實(shí)現(xiàn)和對(duì)比分析</p><p>　　針對(duì)本章介紹的幾種重建算法，文本用matlab進(jìn)行了模擬

85、重建，下面采用的是調(diào)用系統(tǒng)函數(shù)phantom生成的標(biāo)準(zhǔn)Sheep Logon頭模型。如圖3.7。</p><p>　　圖3.7 Sheep Logon頭模型</p><p>　　用iradon函數(shù)用不同濾波器進(jìn)行濾波反投影重建的仿真圖像如圖3.8所示：</p><p>　　圖3.8 對(duì)標(biāo)準(zhǔn)頭模型解析重建的圖像對(duì)比</p><p>　　下面是對(duì)

86、幾組圖像處理中常用的圖片直接讀入再進(jìn)行解析法模擬重建的效果圖，重建過程中采用了各種不同的濾波器。</p><p>　　圖3.9 采用各種不同的濾波器進(jìn)行濾波反投影重建效果對(duì)比</p><p>　　在圖像重建中，一般用視覺相似度SSIM對(duì)圖像重建的效果進(jìn)行評(píng)估。視覺相似度越接近1說明重建圖像與原圖越相進(jìn)。采用不同的濾波器進(jìn)行濾波反投影重建的圖像與原始圖像相比的視覺相似度及重建所用的時(shí)間如表3

87、.1所示。</p><p>　　表3.1采用不同濾波器濾波反投影重建圖像與原圖的視覺相似度對(duì)比</p><p>　　通過圖3.9各重建圖像效果的對(duì)比和表3.1中視覺相似度的對(duì)比可知，直接反投影重建由于沒有濾波，存在偽跡，圖像不清晰，與原圖的視覺相似度很小。采用了不同的濾波器，相對(duì)于沒有濾波的情況，視覺相似度都很接近1，說明與原圖比較接近，圖像質(zhì)量有了明顯的提高。</p>&l

88、t;p>　　以上介紹的是平行束投影重建算法。CT掃描方式有平行束和扇束掃描兩種掃描方式。扇束投影重建算法大致可分為兩類：一類是重排算法，即把一個(gè)視圖中采集到的扇形數(shù)據(jù)重新組合成平行的射線投影數(shù)據(jù)，再用上面所介紹的平行束算法重建。這種方法對(duì)扇束投影的視角與每一扇束投影各射線間的增角有一定的約束條件。另一種算法是扇束投影直接重建算法，這種算法不必先把數(shù)據(jù)重排，而是根據(jù)扇束投影數(shù)據(jù)自身的特點(diǎn)，發(fā)展出直接用扇形束投影數(shù)據(jù)重建圖像的方法。&

89、lt;/p><p>　　扇束掃描自身的算法有兩種，區(qū)別在于抽樣的方式上，一種是等角抽樣，即射線是按等角分布的，將檢測(cè)器等角分布在以放射源為中心的圓弧或直線上；另一種是等間距抽樣，檢測(cè)器單元在直線上作等間距據(jù)分布。</p><p><b>　　這里不作詳細(xì)介紹。</b></p><p>　　第四章 迭代重建算法</p><p>

90、;　　4.1 迭代重建算法的思想原理</p><p>　　迭代重建算法的概念與前面講述的解析法最大的區(qū)別在于前者一開始就把連續(xù)的圖像離散化。在解析算法中，假設(shè)圖像是連續(xù)的，每個(gè)像素只是一個(gè)點(diǎn)，這些離散的點(diǎn)是以圖像顯示為目的的，這些點(diǎn)的選擇可以是隨意的，與圖像重建無關(guān)。但在迭代算法中，每個(gè)像素是個(gè)小面積。這些小面積在計(jì)算當(dāng)前圖像的投影數(shù)據(jù)時(shí)要用到。像素模型對(duì)重建圖像的質(zhì)量好壞影響很大。圖像重建的解析算法與圖像重建的

91、迭代算法之間的另一個(gè)區(qū)別是解析算法著力對(duì)一個(gè)積分方程求解，而迭代算法著力對(duì)一個(gè)線性方程組求解。</p><p>　　迭代的主導(dǎo)思想是，假設(shè)斷層截面由一個(gè)未知的矩陣組成，然后由測(cè)量投影數(shù)據(jù)建立一組未知向量的代數(shù)方程式，通過方程組求解未知圖像向量。</p><p>　　此方法雖然比變換重建法易于理解，但對(duì)于圖像要求較高的應(yīng)用來說，這一方法缺少精度和速度。然而在不可能獲取大量投影數(shù)據(jù)或者投影不均

92、勻分布于180或360之間情況下，迭代法是一個(gè)非常有效的方法。</p><p>　　4.1.1 迭代重建的數(shù)學(xué)模型</p><p>　　圖4.1 斷層截面與投影的幾何布置</p><p>　　如圖4.1所示，圖中將截面離散化，即在圖像上疊加一個(gè)方格網(wǎng)。每一個(gè)小方格內(nèi)像素值是一常量，令表示第個(gè)單元的常量值，N表示像素的總數(shù)。對(duì)于迭代法而言射線的定義稍有不同。具體地講，

93、此時(shí)穿過平面的射線具有一定的寬度，為了說明這一點(diǎn)在圖中用陰影加重了第根射線，且每一根都是如此，絕大多數(shù)情況下射線寬度大約等于圖像像素單元的寬度。此時(shí)線積分被稱為射線和。</p><p>　　如圖4.1所示，投影數(shù)據(jù)也采用單索引表示，設(shè)是圖中由第 i 根射線測(cè)量的射線和。則的關(guān)系可以表示如下，</p><p><b>　　（4.1）</b></p><

94、;p>　　式中為射線總數(shù)（或投影總數(shù)），而為加權(quán)因子，其代表第個(gè)像素對(duì)第根射線積分的貢獻(xiàn)。其值等于第根射線穿過第個(gè)像素所占部分的面積，如圖中所示。需要注意，大多數(shù)均為0，對(duì)于任意一條射線和，只有很少的像素單元有貢獻(xiàn)。如果值很小的話，可以用傳統(tǒng)的矩陣?yán)碚摲椒▽?duì)式(4.1)進(jìn)行求逆。實(shí)際中對(duì)于128 ×128 的圖像陣，=16384，如果要重建256×256的圖像，值將更大。同時(shí)往往具有與同樣量級(jí)的數(shù)。假設(shè)等于16

95、384的情況，</p><p>　　式(4.1)中的矩陣的大小等于16384 ×16384，而這將使任何直接求解矩陣的逆變得無法實(shí)現(xiàn)。當(dāng)然，當(dāng)測(cè)量數(shù)據(jù)混有噪聲或時(shí)，即使值很小，直接矩陣求逆也是不允許的，這種情況可能不得不使用最小二乘法來求解。</p><p>　　對(duì)于為很大值的情況，有一些很吸引人的迭代方法來求解式(4.1)。這些方法是基于“投影法”的，最早被Kaczmarz提出

96、，此后由Tanabe進(jìn)一步闡明。為了說明這些方法中的計(jì)算步驟，可以把(4.1)式寫成如下形式，</p><p><b>　?。?.2）</b></p><p>　　令。因此，由表示的一幅圖像可以看成是維空間中的一個(gè)點(diǎn)。在這一空間中，上面方程組中每一個(gè)方程式表示一個(gè)超平面。當(dāng)這些方程組的唯一解存在的話，這些超平面的交點(diǎn)為一個(gè)點(diǎn)，這就是(4.2)式方程組的解。這一概念進(jìn)一

97、步用圖4.2說明，為理解方便，這里僅考慮有兩個(gè)像素的情況，它們滿足下面的方程，</p><p><b>　?。?.3）</b></p><p>　　圖4.2 Kaczmarz交替投影法求解線性方程組的示意圖</p><p>　　找尋圖4.2中解的計(jì)算過程如下：先從一個(gè)初始圖像開始，把這一初始向量投影到第一根射線，把投影點(diǎn)再投到第二根射線上，然后

98、將投影點(diǎn)再投回到第一根射線上，依此類推。若唯一解存在，迭代終究會(huì)收斂到那一點(diǎn)。</p><p>　　對(duì)于這一方法的計(jì)算機(jī)實(shí)現(xiàn)而言，首先在可行解中估計(jì)一個(gè)初始解。這一估計(jì)解是維空間的一個(gè)點(diǎn)，以表示。一般情況下可以設(shè)置。將這一初始點(diǎn)投影到式(4.2)中第一個(gè)等式所表示的超平面得到，圖4.2所示為兩個(gè)變量的情況。再將投影到式(4.3)中第二個(gè)等式所表示的超平面得到，依次下去。當(dāng)投影到(4.2)式中的第個(gè)等式所表示的超平

99、面時(shí)，得到的過程可以用下面的數(shù)學(xué)方式表示</p><p><b>　　（4.4）</b></p><p>　　式中，而與其自身的內(nèi)積。</p><p>　　如前述，迭代重建的計(jì)算過程如下：給定一個(gè)初始解，將其投影到方程(4.2)表示的某個(gè)超平面上，得到另一個(gè)解，再將其投向下一個(gè)超平面上，依此循環(huán)，當(dāng)所有方程所表示的超平面都被投影過后，得到。在下

100、一次迭代中，作為初始解，再依次投影到(4.2)式中表示的各個(gè)超平面上，獲得，如此重復(fù)。Tanabe證明，若方程(4.2)有唯一解，則有</p><p><b>　　（4.5）</b></p><p>　　當(dāng)式(4.2)描繪的M個(gè)超平面彼此正交，且原方程有唯一解，經(jīng)一次迭代循環(huán)便能收斂。若M個(gè)超平面彼此成很小的夾角，則解的收斂速度較慢。很明顯，超平面間的夾角很大程度上影

101、響了求解的收斂速度。由此可知，為了快速得到收斂解，有兩種途徑可以解決：⑴采用Gram-Schmidt過程超將平面正交化，如Ramakrishnan等提出了一種兩兩正交的方案；⑵仔細(xì)選擇投影超平面的順序，以增大兩次投影超平面間的夾角。</p><p>　　當(dāng)M>N，式(4.2)則不會(huì)收斂到唯一解，而在超平面交點(diǎn)附近產(chǎn)生振蕩效應(yīng)；當(dāng)M<N，式(4.2)式無唯一解，可能會(huì)有無限多個(gè)解。</p>

102、<p>　　在(4.4)式基礎(chǔ)上發(fā)展了一系列其他的迭代重建法，很多都是其的近似表達(dá)。為了在后面部分更好的說明，不妨將式(4.4)重新整理，</p><p><b>　　（4.5）</b></p><p><b>　　式中</b></p><p><b>　?。?.6）</b></p&

103、gt;<p>　　式(4.5)說明將第次迭代解投影到第個(gè)超平上時(shí)，第個(gè)像素的灰度值（其當(dāng)前值為）可以通過差值?得到，</p><p>　　?（4.7）</p><p>　　4.2 代數(shù)重建法ART</p><p>　　上一節(jié)介紹了迭代算法的基本原理，我們知道迭代法的要旨是：從一個(gè)假設(shè)的初始圖像出發(fā)，采用迭代的方法，將根據(jù)人為設(shè)定并經(jīng)理論計(jì)算

104、得到的投影值同實(shí)驗(yàn)測(cè)得的投影值比較，不斷進(jìn)行逼近，按照某種最優(yōu)化準(zhǔn)則尋找最優(yōu)解。這種方法的優(yōu)點(diǎn)是更容易把與空間幾何相關(guān)的、與測(cè)量值條件相關(guān)的各種校正因子包括進(jìn)去。例如容易進(jìn)行</p><p>　　空間分辯的不均勻性校正；</p><p><b>　　散射衰減的校正；</b></p><p>　　物體幾何形狀約束和平滑性約束等操作。</p&

105、gt;<p>　　現(xiàn)在講述有上述原理提出的迭代算法。首先介紹ART方法。本文用三種算法實(shí)現(xiàn)了ART迭代重建。</p><p>　　4.2.1 Kaczmarz 方法</p><p>　　Kaczmarz算法是ART類算法中最經(jīng)典的一種算法。這種方法是一種所謂的行的操作方法，因?yàn)槊看蔚家皰哌^”矩陣A中所有的行。這種方法引入的松弛參數(shù)。</p><p&g

106、t;　　它的實(shí)現(xiàn)公式為：（4.8）</p><p>　　其中，松弛參數(shù)=1，為射線總條數(shù)。</p><p><b>　　其基本步驟為：</b></p><p>　?、旁O(shè)定任一初始圖像向量，一般情況下??；</p><p>　　⑵迭代步驟：式(4.8)，射線號(hào)循環(huán)選擇。</p><p>　　Tanab

107、e[Tanabe1971]證明了，若測(cè)量方程(4.3)是相合的，上述迭代算法收斂到的一個(gè)解。否則，Eggermont及Censor引進(jìn)了松弛系列，目的是有效的減小噪聲。松弛系數(shù)通常是一個(gè)關(guān)于迭代次數(shù)的一個(gè)函數(shù)，隨著的增大而減小，且滿足。由此引出了下面兩種ART方法——Randomized Kaczmarz方法（隨機(jī)Kaczmarz方法）和symmetric Kaczmarz方法（對(duì)稱Kaczmarz方法）。通過改變從此系數(shù)和迭代次數(shù)來改

108、善重建圖像的質(zhì)量。跟Kaczmarz相比，不同之處在于symmetric Kaczmarz方法的而Kaczmarz方法的只取而Randomized Kaczmarz方法（隨機(jī)Kaczmarz方法）隨機(jī)選擇，但選中的概率正比于。</p><p>　　4.3 聯(lián)立迭代重建法SIRT</p><p>　　對(duì)于(4.4)式方程組，還可以采用數(shù)值計(jì)算中的Jacobi迭代法來求解，與ART迭代法不同，

109、Jacobi迭代法采用的是并行迭代，即僅當(dāng)所有投影數(shù)據(jù)都計(jì)算完以后，才對(duì)圖像值進(jìn)行更新。具體來說，采用(4.7)式計(jì)算?，但的值并不立即改變，在經(jīng)歷完所有的等式計(jì)算后，才修正的值，由此可以看出這種修正是每個(gè)像素對(duì)應(yīng)于所有迭代變化的平均值。這種迭代重建算法被稱為同步迭代重建法SIRT。具體步驟：</p><p>　?、旁O(shè)定任一初始圖像向量，一般情況下?。?lt;/p><p><b>　

110、　⑵迭代步驟：</b></p><p><b>　　(4.11)</b></p><p>　　此算法以減慢收斂速度換取比ART較好的圖像質(zhì)量。</p><p>　　本論題模擬圖像重建時(shí)應(yīng)用的類似的方法有聯(lián)合代數(shù)重建法SART（Simultaneous Algebraic Reconstruction Technique）、平均分量法

111、（Component Averaging method）、Cimmino投影法（Cimmino’s Projection method）、對(duì)角松弛征正交投影法（Diagonally Relaxed Orthogonal Projections）和經(jīng)典Landweber法（The Classical Landweber method）?？梢詤⒖枷嚓P(guān)文獻(xiàn)，在此不一一贅述。</p><p>　　4.4 全差分迭代算法（

112、TV）</p><p>　　在上面ART的算法基礎(chǔ)上，再加一些約束條件形成的另一種算法—TV（total variation)算法。待重建的圖像表示為一個(gè)含有維的向量</p><p><b>　　它的運(yùn)算公式如下：</b></p><p><b>　　（4.12）</b></p><p>　　其中表

113、示二維圖像寬和高；×。</p><p>　　4.5 迭代重建算法的實(shí)現(xiàn)和對(duì)比分析</p><p>　　針對(duì)4.2和4.3介紹的幾種方法，本文進(jìn)行了計(jì)算機(jī)模擬重建，采用的是自定義的Sheep Logan頭模型，大小為128×128，如圖4.3。</p><p>　　圖4.3 Sheep—Logan頭模型</p><p>　　

114、利用旋轉(zhuǎn)角度為0：5：180，平行射線束為75的平行掃描投影數(shù)據(jù)，進(jìn)行Kaczmarz算法、Symmetric Kaczmarz算法、Randomized Kaczmarz算法的重建圖像如下。</p><p>　　圖4.4 ART三種算法分別迭代10次的重建效果</p><p>　　對(duì)于迭代法，迭代次數(shù)，掃描的角度和每個(gè)角度下射線的條數(shù)都會(huì)影響圖像重建的效果。下面將從重建圖像的視覺效果上用