步驟（2)，利用Grangeat公式將已知錐形束投影數(shù)據(jù)轉(zhuǎn)化為已知三維Radon數(shù)據(jù) 關(guān)于p的一階導(dǎo)數(shù)。Grangeat公式為：
[0084] 其中，幻足已知三維Radon數(shù)據(jù)關(guān)于P的一階導(dǎo)數(shù)；如附圖2、3、4所示，C〇si3 = S0/SCD，S0是放射源S與坐標(biāo)原點(diǎn)0之間的距離，CD是放射線與探測器平面的交 點(diǎn)，則SCD是放射源S與交點(diǎn)C D之間的距離；P為坐標(biāo)原點(diǎn)0與特征點(diǎn)C之間的距離；舜為 從坐標(biāo)原點(diǎn)0指向特征點(diǎn)C的單位向量；;^/(S(p幻,?)是對步驟（1)中獲取的不完備錐形 束投影數(shù)據(jù)進(jìn)行加權(quán)計算后的投影數(shù)據(jù)值；所用的權(quán)值可優(yōu)選為SO/SA，SA是 放射源S與直線t上點(diǎn)的距離；直線t與線段0CD垂直，且坐標(biāo)原點(diǎn)0到直線t的垂直距離 為s，p，q分別為探測器平面的橫坐標(biāo)軸和縱坐標(biāo)軸，α為線段0CD與探測器平面橫坐標(biāo)軸 P的夾角；
[0085] 然后對已知三維Radon數(shù)據(jù)關(guān)于P的一階導(dǎo)數(shù)/?'/'(ρ??)進(jìn)行積分即可獲得已知三 維Radon數(shù)據(jù)&瓦)〇:
[0086] 步驟（3)，計算已知三維Radon數(shù)據(jù)財(兩)的投影幾何矩變換公式為：
[0088] 其中，為在方向S上的P階投影幾何矩變換。
[0089] 通過式（2)計算階數(shù)P分別為[0，1，2，···，Μ]的一系列投影幾何矩變換，Μ為所使 用的幾何矩的最大階數(shù)。在本實(shí)例中Μ =20。然后，將這些投影幾何矩變換值按一定順序 組成向量Φμ<Χ):
[0091] 步驟（4)，針對已掃描角度范圍，計算系數(shù)矩陣

[0094] 步驟（5)，已知三維Radon數(shù)據(jù)的投影幾何矩變換與幾何圖像矩之間關(guān)系的矩陣 形式為：
[0096] 則利用公式（5)可得幾何圖像矩向量為：
[0098] 幾何圖像矩向量ΨΜ由一系列階數(shù)不超過Μ的幾何圖像矩G_組成，如下所示：
[0101] 步驟（6)，針對未掃描角度范圍，計算系數(shù)矩陣￡1?〇〇:
[0103] 其中，_為未掃描角度范圍內(nèi)的向量，BP

[0104] 步驟（7)，未知Radon數(shù)據(jù)的投影幾何矩變換與幾何圖像矩之間關(guān)系的矩陣形式 為：
[0106] 因?yàn)椴襟E（5)計算出的幾何圖像矩ΨΜ與拜、皮均無關(guān)，所以可將步驟（5)獲得的 幾何圖像矩ΨΜ以及步驟（6)獲得的系數(shù)矩陣:!QM(rY)代入公式（9)，即可計算出未知Radon 數(shù)據(jù)的投影幾何矩變換$_〇〇, 是一系列階數(shù)不大于Μ的投影幾何矩組成 的向量（〇彡P(guān)彡M)。
[0107] 步驟（8)，投影幾何矩逆變換公式如下：
[0109] 則利用公式（10)和步驟（7)獲得的Zp〇〇即可求得未知Radon數(shù)據(jù)
[0110] 步驟（9)，將估計出的未知Radon數(shù)據(jù)按照掃描角度順序和已知Radon數(shù)據(jù)合并， 即可獲得補(bǔ)全的Radon數(shù)據(jù)。
[0111] 步驟（10)，三維Radon逆變換公式如下：
[0113] 利用三維Radon逆變換即可從補(bǔ)全的Radon數(shù)據(jù)中重建出CT圖像/(幻
[0114] 為了對比本發(fā)明所示方法的效果，使用傳統(tǒng)的濾波反投影重建算法對不完備錐形 束投影數(shù)據(jù)進(jìn)行重建，得到的重建圖像如圖5所示；而使用本發(fā)明方法獲得的重建圖像如 圖6所示。比較圖5和圖6的重建結(jié)果，可以看出本發(fā)明方法所得結(jié)果中偽影相對較少，而 且頭模型的輪廓更加完整和清晰。因此可以認(rèn)為本發(fā)明方法能夠在部分投影數(shù)據(jù)缺失的條 件重建出偽影更少、圖像質(zhì)量更好的結(jié)果。
[0115] 以上實(shí)施例僅為說明本發(fā)明的技術(shù)思想，不能以此限定本發(fā)明的保護(hù)范圍，凡是 按照本發(fā)明提出的技術(shù)思想，在技術(shù)方案基礎(chǔ)上所做的任何改動，均落入本發(fā)明的保護(hù)范 圍之內(nèi)。
【主權(quán)項】
1. 基于幾何圖像矩的有限角度錐形束CT圖像重建方法，其特征在于，包括以下步驟： (1) 獲取有限角度掃描條件下的不完備的錐形束投影數(shù)據(jù)； (2) 利用Grangeat公式將步驟（1)獲得的錐形束投影數(shù)據(jù)轉(zhuǎn)化為已知三維Radon數(shù) 據(jù)； (3) 對步驟（2)獲得的已知三維Radon數(shù)據(jù)進(jìn)行投影幾何矩變換； (4) 建立步驟（3)獲得的已知三維Radon數(shù)據(jù)的投影幾何矩變換與幾何圖像矩之間的 關(guān)系式； (5) 利用步驟（4)建立的關(guān)系式，從步驟（3)獲得的已知三維Radon數(shù)據(jù)的投影幾何矩 變換中計算出幾何圖像矩； (6) 建立未知三維Radon數(shù)據(jù)與幾何圖像矩之間的關(guān)系式； (7) 根據(jù)步驟（6)建立的關(guān)系式，從步驟（5)獲得的幾何圖像矩中估計出未知三維 Radon數(shù)據(jù)的投影幾何矩變換； (8) 利用投影幾何矩逆變換從步驟（7)獲得的未知三維Radon數(shù)據(jù)的投影幾何矩變換 中求出未知三維Radon數(shù)據(jù)； (9) 將步驟⑵得到的已知三維Radon數(shù)據(jù)和步驟⑶得到的未知三維Radon數(shù)據(jù)進(jìn) 行數(shù)據(jù)拼合，獲得補(bǔ)全的三維Radon數(shù)據(jù)； (10) 通過三維Radon逆變換，根據(jù)步驟（9)獲得的補(bǔ)全的三維Radon數(shù)據(jù)重建出CT圖 像。2. 如權(quán)利要求1所述基于幾何圖像矩的有限角度錐形束CT圖像重建方法，其特征在 于：在步驟（1)中，所述有限角度掃描條件下是指在[<p,范圍內(nèi)旋轉(zhuǎn)掃描，其中， 0 ()<φ<360°，3. 如權(quán)利要求2所述基于幾何圖像矩的有限角度錐形束CT圖像重建方法，其特征在 于，步驟（2)的具體過程如下： 首先利用Grangeat公式將錐形束投影數(shù)據(jù)轉(zhuǎn)化為已知三維Radon數(shù)據(jù)關(guān)于P的一階 導(dǎo)數(shù)：式（1)中，是已知三維Radon數(shù)據(jù)關(guān)于P的一階導(dǎo)數(shù)，cos β = S0/SCD，Cd是 放射線與探測器平面的交點(diǎn)，則SCd是放射源S與交點(diǎn)C D之間的距離，SO是放射源S與坐 標(biāo)原點(diǎn)0之間的距離，P為坐標(biāo)原點(diǎn)0與特征點(diǎn)C之間的距離，S為從坐標(biāo)原點(diǎn)0指向特征 點(diǎn)C的單位向量，是對步驟（1)中獲取的不完備錐形束投影數(shù)據(jù)進(jìn)行加權(quán)計 算后的投影數(shù)據(jù)值，直線t與線段OCd垂直，且坐標(biāo)原點(diǎn)0到直線t的垂直距離為s，p，q分 別為探測器平面的橫坐標(biāo)軸和縱坐標(biāo)軸，α為線段OCd與探測器平面橫坐標(biāo)軸p的夾角； 然后對已知三維Radon數(shù)據(jù)關(guān)于P的一階導(dǎo)數(shù)進(jìn)行積分，即獲得已知三維 Radon數(shù)據(jù)，所述三維Radon數(shù)據(jù)定義為：式⑵中，/(幻是三維圖像f在添處的灰度值，Θ 是單位向量海__與坐標(biāo)軸Z的夾角。4. 如權(quán)利要求3所述基于幾何圖像矩的有限角度錐形束CT圖像重建方法，其特征在 于：將SO/SA作為權(quán)值，對步驟（1)中獲取的不完備錐形束投影數(shù)據(jù)進(jìn)行加權(quán)計算后得到投 影數(shù)據(jù)值其中，SA是放射源S與直線t上任意點(diǎn)A的距離。5. 如權(quán)利要求3所述基于幾何圖像矩的有限角度錐形束CT圖像重建方法，其特征在 于，所述步驟（3)中的投影幾何矩變換的定義為：式（3)中，Lp 〇i)為在方向i上的p階投影幾何矩變換數(shù)據(jù)。6. 如權(quán)利要求5所述基于幾何圖像矩的有限角度錐形束CT圖像重建方法，其特征在 于，所述步驟（4)中的幾何圖像矩的定義為：式（4)中，f (X，y, z)是三維圖像f在點(diǎn)（X，y, z)處的灰度值。 則已知三維Radon數(shù)據(jù)的p階投影幾何矩變換與幾何圖像矩之間的關(guān)系式為：7. 如權(quán)利要求6所述基于幾何圖像矩的有限角度錐形束CT圖像重建方法，其特征在 于，所述步驟（5)的具體過程如下： 首先將式（5)改寫成如式（7)所示的矩陣形式：M為所用幾何矩的最大階數(shù)； 然后對式（7)采用矩陣除法，獲得幾何圖像矩向量：8. 如權(quán)利要求6或7所述基于幾何圖像矩的有限角度錐形束CT圖像重建方法，其特征 在于，所述步驟(6)中的建立幾何圖像矩與未知三維Radon數(shù)據(jù)的投影幾何矩變換之間的 關(guān)系式為：式（12)中，IpOO是未知Radon數(shù)據(jù)的投影幾何矩變換，為未掃描角度范圍內(nèi)的向 量，即9. 如權(quán)利要求8所述基于幾何圖像矩的有限角度錐形束CT圖像重建方法，其特征在 于，所述步驟（7)的具體過程為： 首先將式（12)改寫成如式（13)所示的矩陣形式：然后將步驟（5)獲得的幾何圖像矩向量ΦΜ代入式（13)，計算出未知三維Radon數(shù)據(jù) 的投影幾何矩變換10. 如權(quán)利要求9所述基于幾何圖像矩的有限角度錐形束CT圖像重建方法，其特征在 于，所述步驟（8)中投影幾何矩逆變換的公式為：其中，是未知三維Radon數(shù)據(jù)； 所述步驟（10)中三維Radon逆變換的公式為：I. 其中，八幻是重建的CT圖像。
【專利摘要】本發(fā)明公開了基于幾何圖像矩的有限角度錐形束CT圖像重建方法，步驟依次為，將獲得的錐形束投影數(shù)據(jù)轉(zhuǎn)化為已知Radon數(shù)據(jù)、獲得已知Radon數(shù)據(jù)的投影幾何矩變換、根據(jù)已知Radon數(shù)據(jù)的投影幾何矩變換中計算出幾何圖像矩、根據(jù)幾何圖像矩估計出未知Radon數(shù)據(jù)的投影幾何矩變換、經(jīng)逆變換求出未知Radon數(shù)據(jù)、將步已知Radon數(shù)據(jù)和未知Radon數(shù)據(jù)進(jìn)行數(shù)據(jù)拼合，獲得補(bǔ)全的三維Radon數(shù)據(jù)并重建CT圖像。本發(fā)明能夠在減小掃描范圍的條件下重建出符合臨床診斷要求、高質(zhì)量的錐形束CT圖像。
【IPC分類】G06T15/06, G06T11/00
【公開號】CN105354868
【申請?zhí)枴緾N201510644674
【發(fā)明人】戴修斌, 劉天亮, 胡棟, 晏善成, 石丹丹, 鄧黃建
【申請人】南京郵電大學(xué)
【公開日】2016年2月24日
【申請日】2015年10月8日