基于幾何圖像矩的有限角度錐形束ct圖像重建方法
【技術(shù)領(lǐng)域】
[0001] 本發(fā)明涉及一種醫(yī)學(xué)圖像重建方法，特別涉及了一種CT圖像重建方法。
【背景技術(shù)】
[0002] 作為目前一種常規(guī)有效的臨床醫(yī)學(xué)診斷工具，X射線計(jì)算機(jī)斷層成像技術(shù)（x-ray Computerized Tomography, CT)為臨床醫(yī)生的診斷提供了豐富的人體器官組織信息。但是 由相關(guān)研究表明：一次完整的CT掃描通常伴隨著較高程度的電離輻射，而高劑量電離輻射 可誘發(fā)人體新陳代謝異常乃至癌癥、白血病等疾病。因此，如何在降低X射線使用劑量的同 時(shí)，保證重建圖像質(zhì)量滿足臨床診斷要求成為醫(yī)學(xué)圖像處理領(lǐng)域研究的重點(diǎn)。
[0003] 臨床上減少病患輻射量的重要方法之一就是減小CT掃描范圍，即將探測(cè)器的旋 轉(zhuǎn)角度范圍限制在某個(gè)小于標(biāo)準(zhǔn)的區(qū)間內(nèi)，從而在總體上大幅減少了患者所受X射線輻射 量。雖然限制CT設(shè)備掃描范圍能夠降低患者所受X射線輻射量，但同時(shí)會(huì)造成所獲CT投 影數(shù)據(jù)部分缺失，即獲得的是不完備投影數(shù)據(jù)，使重建CT圖像質(zhì)量明顯下降，以至于無法 滿足臨床診斷的需要。
[0004] 目前針對(duì)有限角度CT圖像重建方法主要是基于統(tǒng)計(jì)模型的迭代重建方法。雖然 該方法在掃描角度有限的情況下仍能獲得較好的重建結(jié)果，但是該方法需要對(duì)目標(biāo)函數(shù)進(jìn) 行上百次的迭代求解，重建CT圖像所需計(jì)算時(shí)間大幅增加，遠(yuǎn)遠(yuǎn)超過經(jīng)典的解析重建方 法，不能滿足臨床CT實(shí)時(shí)顯像的要求。
[0005] 經(jīng)典的解析重建方法雖然重建速度快，但是在掃描角度有限的情況下會(huì)丟失圖像 原有細(xì)節(jié)信息，從而導(dǎo)致重建圖像中出現(xiàn)大量偽影和噪聲，嚴(yán)重影響了對(duì)特征點(diǎn)的分辨。

【發(fā)明內(nèi)容】

[0006] 為了解決上述【背景技術(shù)】提出的技術(shù)問題，本發(fā)明旨在提供基于幾何圖像矩的有限 角度錐形束CT圖像重建方法，能夠在減小掃描范圍---即投影數(shù)據(jù)不完備條件下重建出符 合臨床診斷要求、高質(zhì)量的錐形束CT圖像。
[0007] 為了實(shí)現(xiàn)上述技術(shù)目的，本發(fā)明的技術(shù)方案為：
[0008] 基于幾何圖像矩的有限角度錐形束CT圖像重建方法，包括以下步驟：
[0009] (1)獲取有限角度掃描條件下的不完備的錐形束投影數(shù)據(jù)；
[0010] (2)利用Grangeat公式將步驟（1)獲得的錐形束投影數(shù)據(jù)轉(zhuǎn)化為已知三維Radon 數(shù)據(jù)；
[0011] (3)對(duì)步驟⑵獲得的已知三維Radon數(shù)據(jù)進(jìn)行投影幾何矩變換；
[0012] (4)建立步驟（3)獲得的已知三維Radon數(shù)據(jù)的投影幾何矩變換與幾何圖像矩之 間的關(guān)系式；
[0013] (5)利用步驟⑷建立的關(guān)系式，從步驟⑶獲得的已知三維Radon數(shù)據(jù)的投影幾 何矩變換中計(jì)算出幾何圖像矩；
[0014] (6)建立未知三維Radon數(shù)據(jù)與幾何圖像矩之間的關(guān)系式；
[0015] (7)根據(jù)步驟（6)建立的關(guān)系式，從步驟（5)獲得的幾何圖像矩中估計(jì)出未知三維 Radon數(shù)據(jù)的投影幾何矩變換；
[0016] (8)利用投影幾何矩逆變換從步驟（7)獲得的未知三維Radon數(shù)據(jù)的投影幾何矩 變換中求出未知三維Radon數(shù)據(jù)；
[0017] (9)將步驟（2)得到的已知三維Radon數(shù)據(jù)和步驟（8)得到的未知三維Radon數(shù) 據(jù)進(jìn)行數(shù)據(jù)拼合，獲得補(bǔ)全的三維Radon數(shù)據(jù)；
[0018] (10)通過三維Radon逆變換，根據(jù)步驟（9)獲得的補(bǔ)全的三維Radon數(shù)據(jù)重建出 CT圖像。
[0019] 進(jìn)一步地，在步驟（1)中，所述有限角度掃描條件下是指在[隊(duì)360() -識(shí)]范圍內(nèi)旋 轉(zhuǎn)掃描，其中，
[0020] 進(jìn)一步地，步驟（2)的具體過程如下：
[0021] 首先利用Grangeat公式將錐形束投影數(shù)據(jù)轉(zhuǎn)化為已知三維Radon數(shù)據(jù)關(guān)于P的 一階導(dǎo)數(shù)：
[0023] 式（1)中，/?/'〇>萬）是已知三維Radon數(shù)據(jù)關(guān)于P的一階導(dǎo)數(shù)，cos β = S0/SCD， CD是放射線與探測(cè)器平面的交點(diǎn)，則SC D是放射源S與交點(diǎn)C D之間的距離，SO是放射源S 與坐標(biāo)原點(diǎn)〇之間的距離，P為坐標(biāo)原點(diǎn)〇與特征點(diǎn)C之間的距離；S為從坐標(biāo)原點(diǎn)0指 向特征點(diǎn)c的單位向量；;是對(duì)步驟（1)中獲取的不完備錐形束投影數(shù)據(jù)進(jìn)行 加權(quán)計(jì)算后的投影數(shù)據(jù)值，直線t與線段0CD垂直，且坐標(biāo)原點(diǎn)0到直線t的垂直距離為s， P，q分別為探測(cè)器平面的橫坐標(biāo)軸和縱坐標(biāo)軸，α為線段〇CD與探測(cè)器平面橫坐標(biāo)軸p的 夾角；
[0024] 然后對(duì)已知三維Radon數(shù)據(jù)關(guān)于P的一階導(dǎo)數(shù)/?7\>幻進(jìn)行積分，即獲得已知三 維Radon數(shù)據(jù)，所述三維Radon數(shù)據(jù)定義為：
[0026] 式⑵中，/:(?是三維圖像f在i處的灰度值，萬=(sin 0 cos φ, sin Θ sin φ, cos Θ)， Θ是單位向量S:與坐標(biāo)軸z的夾角。
[0027] 進(jìn)一步地，將S0/SA作為權(quán)值，對(duì)步驟（1)中獲取的不完備錐形束投影數(shù)據(jù)進(jìn)行加 權(quán)計(jì)算后得到投影數(shù)據(jù)值久w/(s(pfiU)，其中，SA是放射源S與直線t上任意點(diǎn)A的距離。
[0028] 進(jìn)一步地，所述步驟（3)中的投影幾何矩變換的定義為：
[0030] 式（3)中，ip(H)為在方向兩上的p階投影幾何矩變換數(shù)據(jù)。
[0031] 進(jìn)一步地，所述步驟（4)中的幾何圖像矩的定義為：
[0033] 式⑷中，f (X，y, z)是三維圖像f在點(diǎn)（X，y, z)處的灰度值。
[0034] 則已知三維Radon數(shù)據(jù)的p階投影幾何矩變換與幾何圖像矩之間的關(guān)系式為：
[0038] 進(jìn)一步地，所述步驟（5)的具體過程如下：
[0039] 首先將式（5)改寫成如式（7)所示的矩陣形式：
[0046] Μ為所用幾何矩的最大階數(shù)；
[0047] 然后對(duì)式（7)采用矩陣除法，獲得幾何圖像矩向量：
[0049] 進(jìn)一步地，所述步驟(6)中的建立幾何圖像矩與未知三維Radon數(shù)據(jù)的投影幾何 矩變換之間的關(guān)系式為：
[0051] 式（12)中，Σμ(?1')足未知Radon數(shù)據(jù)的投影幾何矩變換，g;為未掃描角度范圍內(nèi) 的向量，即
[0053] 進(jìn)一步地，所述步驟（7)的具體過程為：
[0054] 首先將式（12)改寫成如式（13)所示的矩陣形式：

[0061] 然后將步驟（5)獲得的幾何圖像矩向量ΨΜ代入式（13)，計(jì)算出未知三維Radon 數(shù)據(jù)的投影幾何矩變換(冗）.
[0062] 進(jìn)一步地，所述步驟（8)中投影幾何矩逆變換的公式為：
[0064] 其中，是未知三維Radon數(shù)據(jù)；
[0065] 所述步驟（10)中三維Radon逆變換的公式為：
[0067] 其中，八幻足重建的CT圖像。
[0068] 采用上述技術(shù)方案帶來的有益效果：
[0069] (1)本發(fā)明能夠在減小掃描范圍一即投影數(shù)據(jù)不完備條件下重建出符合臨床診 斷要求、高質(zhì)量的錐形束CT圖像；
[0070] (2)由于本發(fā)明不涉及迭代運(yùn)算，所以計(jì)算時(shí)間少于統(tǒng)計(jì)迭代方法；
[0071] (3)本發(fā)明可使用較少的錐形束投影數(shù)據(jù)重建出高質(zhì)量的錐形束CT圖像，故本發(fā) 明能在保證重建CT圖像質(zhì)量的前提下，有效地減少患者所受X射線輻射量。
【附圖說明】
[0072] 圖1是本發(fā)明方法的基本流程示意圖；
[0073] 圖2是錐形束投影的CT成像原理示意圖；
[0074] 圖3是步驟⑵中單位向量的幾何關(guān)系示意圖；
[0075] 圖4是步驟（2)中直線t的幾何關(guān)系示意圖；
[0076] 圖5是使用傳統(tǒng)的濾波反投影重建算法得到的重建圖像；
[0077] 圖6是使用本發(fā)明得到的重建圖像。
【具體實(shí)施方式】
[0078] 以下將結(jié)合附圖，對(duì)本發(fā)明的技術(shù)方案進(jìn)行詳細(xì)說明。
[0079] 基于幾何圖像矩的有限角度錐形束CT圖像重建方法，如圖1所示，包括如下步 驟：
[0080] 以三維Shepp-logan頭模型為計(jì)算機(jī)仿真實(shí)驗(yàn)對(duì)象，尺寸為512*512*60,模擬CT 機(jī)的X射線源到旋轉(zhuǎn)中心和探測(cè)器的距離分別為570mm和1040mm。錐形束CT設(shè)備旋轉(zhuǎn)掃 描范圍為[γ，360° -γ]，其中γ =50° ;角度采樣值為1160;每個(gè)采樣角對(duì)應(yīng)672個(gè)探 測(cè)器單元，探測(cè)器單元的大小為1. 407mm。
[0081] 步驟（1)，獲取有限角度掃描條件下的錐形束投影數(shù)據(jù)燈〇(ρ??),〇·如附 圖2、3、4所示，P為坐標(biāo)原點(diǎn)0與特征點(diǎn)C之間的距離，直線t與線段0CD垂直，且 坐標(biāo)原點(diǎn)〇到直線t的垂直距離為S 為從坐標(biāo)原點(diǎn)0指向特征點(diǎn)C的單位向量， 5 = (sin Θ a>Si/?, sin Θ sin φ, cos Θ)。由于CT設(shè)備未按常規(guī)沿軌道在[0°，360° )范圍內(nèi) 旋轉(zhuǎn)掃描一圈，而只是在有限角度[γ，360° -γ]范圍內(nèi)旋轉(zhuǎn)掃描。因此，對(duì)于已獲得的錐 形束投影數(shù)據(jù)，其對(duì)應(yīng)的單位向量i的參數(shù)φ e |_y,36011 -ρμ值得注意的是， ΙΟ^φ)和(360() - φ, 360'3)角度范圍內(nèi)的錐形束投影數(shù)據(jù)未知。
[0082]