CT數據一致性條件及其應用綜述

湯少杰， 俞恒永， 牟軒沁

(1.西安郵電大學 自動化學院， 西安 710121; 2.西安機器人智能系統國際科技合作基地， 西安 710121; 3.馬薩諸塞大學 洛厄爾分校電氣與計算機工程系， 馬薩 01854; 4. 西安交通大學 電子與信息工程學院， 西安 710121)

0 引言

當今計算機斷層掃描(Computed Tomography, CT)技術在醫療、工業、安檢等眾多領域有著不可替代的作用。X射線CT技術的發展經歷了長期的過程。1901年，第一個諾貝爾物理學獎因發現X射線被授予德國物理學家R?ntgen(圖1)。這是X射線CT技術的重要物理前提。

圖1 (a)德國物理學家R?ntgen;(b)其夫人左手X光照片

1917年,奧地利數學家Radon提出二維Radon變換及其逆變換[1]，成為現代CT重建理論的基礎。1979年,諾貝爾生理學或醫學獎因CT研發做出貢獻被授予英國工程師Hounsfield[2]和美籍南非裔物理學家Cormack(圖2)。

圖2 (a)英國工程師Hounsfield；(b)美籍南非裔物理學家Cormack

在應用推動下，二維CT之后出現了多排CT，投影幾何從二維平行束發展到等角扇束、等距扇束，掃描軌跡也從直線軌跡、圓形軌跡發展到一般軌跡。適應于多排CT的FDK近似重建算法在1984年由Feldkamp、Davis 、Kress 三人提出[3]。隨著探測器工業的進步，多排CT之后又出現了錐束CT，同時，錐束重建理論也飛速發展?；谌SRadon變換及其逆變換，Tuy和Smith提出了作為錐束CT重建充要條件的Tuy-Smith完備理論及相應重建算法[4-5]。隨后，基于三維Radon變換導函數，Grangeat提出了算法流程更為適用的錐束CT重建算法[6]。1993年，Wang等提出螺旋錐束CT成像方式及相應的近似重建算法[7]，2000年，Kachelrie?等提出了ASSR近似重建算法算法[8]。進入21世紀后，錐束CT重建理論研究迎來了關鍵時期。在此期間，三維FBP類型的Katsevich算法[9]以及三維BPF類型算法[10-12]相繼被提出并實現，從而攻克了長久以來的難題——螺旋(甚至更為復雜掃描軌跡)錐束CT精確重建。

眾所周知，CT成像包含了眾多軟硬件技術，任一部分的不完善都會導致CT重建圖像中表現出偽影[13]，影響CT圖像質量與后續診斷。對上述CT重建算法與偽影校正感興趣的讀者推薦深入研究本文相關參考文獻[19-21]。

CT理論中還包含數據一致性條件方面的研究[22]。數據一致性條件有著非常多樣的表現形式和性質。其中，局部形式包括原函數(如共軛一致性)與微分等式；而全局形式包括積分等式、不變量和不變式。局部形式比較適合于部分數據補全等，而全局形式比較適合于全局參數估計或偽影校正等。這些不同的表現形式，能在一定程度上影響數據一致性條件的利用難度，以及對特定CT應用的適應性。本文重點綜述CT數據一致性條件形式和性質及其在CT成像中的應用。

1 映射算子及成像幾何

CT理論中的各種映射算子[如Radon變換、X-ray變換、發散束(扇束、錐束)變換等]將物體函數映射為非物體函數。由于映射算子固有性質顯示出了一些獨立于物體函數的獨特性質(例如，對于二維平行束正弦圖中，每個投影的積分都是常數)，將這些特性統稱為數據一致性條件[22]。本文映射算子所涉及的CT成像投影幾何如圖3所示。

圖3 各種CT成像投影幾何(a)圓形軌跡平行束；(b)圓形軌跡等角扇束；(c)圓形軌跡等距扇束；(d)圓形軌跡等角錐束；(e)圓形軌跡等距錐束

2 數據一致性條件及其應用

2.1 共軛一致性條件

共軛一致性條件(Conjugate Consistency Condition，CCC)[23～26]一般應用于二維投影幾何中，也可以近似推廣到三維投影幾何。

2.1.1 共軛一致性條件

在二維投影幾何中，共軛一致性條件特指[23-26]：

<圓形軌跡平行束>

p2(s,θ)=p2(-s,θ+π)

(1)

<圓形軌跡等角扇束>

(2)

<圓形軌跡等距扇束>

(3)

式中，θ為圓形軌跡平行束投影角度；β為圓形軌跡扇束投影角度；γ為圓形軌跡等角扇束扇角；s為圓形軌跡平行束的探測器像素在探測器坐標系上的空間位置；u為圓形軌跡等距扇束的探測器像素在探測器坐標系上的空間位置；D為源點到探測器中心點距離。

2.1.2 特點

共軛一致性條件的物理意義：顛倒任意一條投影路徑兩端的X射線源點與探測器像素點位置，理想情況下所得到投影數據嚴格相等。

共軛一致性條件的特點：① 原理簡單易掌握；② 在投影域具有局部性；③ 投影原函數形式，不需要積分或微分；④ 可方便地從二維擴展到三維。

2.1.3 應用

Tang等基于二維共軛一致性條件，構造并證明“共軛投影誤差加權和等于零”，并可將其應用于運動校正[23]。2005年與2006年，Tang等將三維共軛一致性條件分別用于圓形錐束與螺旋錐束CT重建算法設計中，通過三維加權方案來有效實現錐束偽影校正[24-25]。2018年，Li等嘗試將三維加權方案應用到錐束CT重建迭代算法設計中，取得了一定效果[26]。

圖4 基于共軛一致性條件的錐束偽影校正[24](a)校正前圖像；(b) 校正后圖像(經原作者同意)

2.2 微分限制條件

微分限制條件(Differential Constraint Condition，DCC)[24]應用于二維平行束投影幾何中。

2.2.1 微分限制條件

在二維平行束投影幾何中，微分限制條件特指[27]：

<圓形軌跡平行束>

(4)

2.2.2 特點

2.2.3 應用

2011年Tang等基于二維平行束微分關系[31]：

<平行束>

(5)

得到了微分限制條件[式(4)]，并基于兩步Hilbert變換[28]將其應用于二維梯度圖像重建[27]。

圖5 基于兩步Hilbert變換與微分限制條件的二維梯度圖像重建[27](a)(b)無截斷；(c)(d)有截斷(經原作者同意)

2.3 Helgason-Ludwig一致性條件

Helgason-Ludwig一致性條件(Helgason-Ludwig Consistency Condition，HLCC)[29-31]應用于二維平行束投影幾何中。

2.3.1 Helgason-Ludwig一致性條件

在二維平行束投影幾何中，Helgason-Ludwig一致性條件特指[31]：

<圓形軌跡平行束>

(6)

式中，左側為二維平行束投影矩，右側為二維圖像矩作為系數的cosθ與sinθ的m階多項式。

2.3.2 特點

Helgason-Ludwig一致性條件的物理意義：反映了二維平行束投影矩與二維圖像矩之間的等式關系。

Helgason-Ludwig一致性條件的特點：① 投影矩為cosθ與sinθ的m階多項式；② 在投影域具有全局性；③ 投影積分形式；④ 速降函數空間的充要條件。

2.3.3 應用

2000年，Basu等基于Helgason-Ludwig一致性條件構造并證明了一個確定投影數據的所對應變量(s,θ)的條件[32]，該方法可用于MRI與三維電鏡動態成像。Helgason-Ludwig一致性條件還可用于CT投影數據截斷補全或FOV擴大[33-36]，任一角度投影數據估計[37-38]，投影數據任意損失補全[39]，以及基于運動估計的動態成像[40-43]、射束硬化偽影校正[44-47]、散射校正[48]。

圖6 基于Helgason-Ludwig一致性條件的射束硬化偽影校正[47](a)單能譜重建圖像；(b)寬能譜重建圖像；(c)校正后圖像(經原作者同意)

2.4 扇束數據一致性條件

扇束數據一致性條件(Fan-beam Data Consistency Condition，FDCC)[47]應用于二維扇束投影幾何中。

2.4.1 扇束數據一致性條件

在二維等角扇束投影幾何中，2005年，Chen等證明的扇束數據一致性條件特指[49]：

<圓形軌跡等角扇束>

(7)

該式本質上是投影域一維有限Hilbert逆變換[50-51]，其中

(8)

(9)

β′=π+β+γ+sin-1(p′/R)

(10)

(11)

2.4.2 特點

扇束數據一致性條件的物理意義：投影數據可由共軛側投影數據的一維有限Hilbert逆變換得到。

扇束數據一致性條件的特點：① 比共軛一致性條件更復雜，可看作不同的插值方法；② 投影積分形式；③ 存在奇異積分。

2.4.3 應用

扇束數據一致性條件可用于缺失數據修補恢復[49]，也可用于運動偽影校正[52]。

2.5 扇束數據充要條件

扇束數據充要條件[53]應用于直線軌跡二維扇束投影幾何中。

2.5.1 扇束數據充要條件

在二維等角扇束投影幾何中，2013年Clackdoyle證明的扇束數據充要條件特指[53]：

<直線軌跡等角扇束>

(12)

式中，τ對應于源點在直線軌跡上的位置。該充要條件也有等距扇束形式，此處不再贅述。

2.5.2 特點

扇束數據充要條件的特點：① 投影積分形式；② 不存在奇異積分；③ 相比通過變量代換推導的扇束投影數據Helgason-Ludwig一致性條件，扇束數據充要條件形式上更類似于平行束Helgason-Ludwig一致性條件，因此使用將更方便。

2.5.3 應用

扇束數據充要條件可能應用于丟失數據補全。

2.6 截斷數據一致性條件

截斷數據一致性條件(Truncated Data Consistency Condition)[54-55]可應用于一般軌跡二維平行束或扇束投影幾何中。

2.6.1 截斷數據一致性條件

在二維平行束或扇束投影幾何中，2015年Clackdoyle等證明的圓形軌跡截斷數據一致性條件特指[54]：

<圓形軌跡等角扇束>

(13)

(14)

(15)

(16)

2.6.2 特點

截斷數據一致性條件的物理意義：Bn(x1)為x1的n階多項式，直線軌跡扇束數據一致性條件在圓形軌跡上的應用。

截斷數據一致性條件的特點：① 適用于特定截斷數據情形；② 投影積分形式；③ 避免了奇異積分。

2.6.3 應用

截斷數據一致性條件可用于運動偽影校正[54-55]。

2.7 積分不變量

積分不變量(Integral Invariants)[56]可應用于發散束(扇束、錐束)投影幾何中。

2.7.1 積分不變量

在發散束幾何中，2006年Wei證明了大量的積分不變量[56]，因為篇幅原因僅列舉如下直線軌跡等角扇束與等角錐束各一例：

<直線軌跡等角扇束>

(17)

式中，τ對應于源點在垂直坐標軸上的位置；b為任意實數，但要求投影不發生截斷。

<直線軌跡等角錐束>

(18)

2.7.2 特點

積分不變量的物理意義：I2與IIV分別為源點在垂直坐標軸上的位置τ不變量。

積分不變量的特點：① 與運動對稱群理論相關；② 投影積分形式；③ 特定形式的積分不變量需奇異積分。

圖7 積分不變量的奇異性[58](a)非奇異情形;(b)奇異積分情形(經原作者同意)

2.7.3 應用

積分不變量可用于運動偽影校正[56]，成像系統幾何校正[57]，也可用于稀疏角度投影的射束硬化校正[58]。根據運動的相對性，當源點不動物體沿著垂直方向運動時，I2與IIV積分不變量無法感知到物體運動。因此，對運動偽影進行校正時，需要清晰理解所采用的積分不變量對不同運動的感知能力。

圖8 基于積分不變量的射束硬化偽影校正[58](a)校正前圖像；(b)校正后圖像注：稀疏角度投影數為64，采用環形模式不導致奇異積分(經原作者同意)

2.8 錐束數據一致性條件

錐束數據一致性條件[59]顧名思義可應用于錐束投影幾何中。

2.8.1 錐束數據一致性條件

在三維錐束投影幾何中，2016年Clackdoyle等[59]證明了一族錐束數據一致性條件：

<圓形軌跡等距錐束>

(19)

(20)

(21)

2.8.2 特點

錐束數據一致性條件的物理意義：Mn(β)為cosβ與sinβ的n階多項式[59]。

錐束數據一致性條件的特點：① 投影積分形式；② 存在奇異積分；③ 相比通過變量代換推導的扇束投影數據Helgason-Ludwig一致性條件，錐束數據一致性條件形式上更類似于平行束Helgason-Ludwig一致性條件，因此使用將更方便。

2.8.3 應用

錐束數據一致性條件可應用于丟失數據補全，也可用于投影非線性或全局偽影等整體缺陷因素的校正，例如，射束硬化和散射校正或全局運動參數估計等。

2.9 約翰方程

約翰方程(John’s Equations)[60]可應用于三維錐束投影幾何中。

2.9.1 約翰方程

在三維錐束投影幾何中，1938年John[60]證明了一組有四個獨立變量的超雙曲方程：

<三維錐束>

i,j=1,2,3

(22)

2.9.2 特點

約翰方程的物理意義：三維錐束投影必須滿足四個獨立變量的超雙曲方程組。

約翰方程的特點：① 投影微分形式；② 可用傅里葉變換也可用PDE進行數值計算；③ 推導到不同投影幾何與掃描軌跡較困難易發生錯誤。

2.9.3 應用

約翰方程可應用于三維PET投影數據重排[61]、基于三維傅里葉變換的錐束投影數據補全[62]、基于PDE的投影數據補全[63]、心臟運動偽影校正[64]、螺旋錐束CT投影數據重排[65]、變螺距錐束CT重建[66]、散射校正[67-68]。另外，2005年Sidky等針對X射線源點位于二維圓柱面上的三鞍線情形，提出一個不同的錐束投影數據一致性條件[69]。與約翰方程稍有不同，并猜想可用于缺失錐束投影數據補全，并推斷如猜想成立，則三鞍線X射線源軌跡下整個三維ROI中的圖像函數可重建。

圖9 基于約翰方程的散射校正[67](即PC-VI算法)的原理性示意圖(經原作者同意)

圖10 基于約翰方程的散射校正[68]其中Uncorrected為未校正結果，SI為空域插值結果，PC-VI為[67]結果，JECC為[68]結果(經原作者同意)

2.10 頻域雙楔形零能量區域性質

頻域雙楔形零能量性質[22]可應用于平行束、扇束投影幾何中。

2.10.1 頻域雙楔形零能量性質

2010年Mazin等將二維平行束投影幾何的頻域雙楔形零能量性質推廣到二維扇束投影幾何，該性質作為二維投影數據必須滿足的必要條件也可看作一種獨特的數據一致性條件[22]。

2.10.2 特點

頻域雙楔形零能量性質的物理意義：ROI內的圖像函數能量落在頻域雙楔形零能量區域外，ROI外的圖像函數能量落在全頻域。

頻域雙楔形零能量性質的特點：① 在頻域工作可利用FFT高效率性質；② 可用于數據壓縮；③ 可用于投影噪聲能量估計。

2.10.3 應用

2010年Mazin等將頻域雙楔形零能量性質應用于提高迭代重建算法效率的工作中[22]。2013年Bai等將頻域雙楔形零能量性質應用于投影噪聲能量估計，并將其用于統計迭代重建的正則化參數選擇[70]。

圖11 基于頻域雙楔形零能量性質的投影噪聲能量估計[70](a)錐束投影數據；(b)無噪數據傅里葉空間；(c)含噪數據傅里葉空間(經原作者同意)

3 討論

數據一致性條件有著非常多樣的表現形式，如：

(1)平行束扇束錐束投影幾何；

(2)直線圓形螺旋一般掃描軌跡；

(3)局部全局運算形式；

(4)原函數微分積分運算形式；

(5)積分運算為奇異非奇異；

(6)不變量多項式頻域特定性質。

局部形式包括原函數(如共軛一致性)與微分運算形式；而全局形式包括積分運算形式。局部形式允許投影域截斷，而全局形式不允許投影域截斷。非奇異形式一般數值計算穩定，而奇異形式一般數值計算不穩定。這些不同的表現形式，能在一定程度上影響數據一致性條件的利用難度，以及對特定CT應用的適應性。

在Gowers主編的《普林斯頓數學指南(第一卷)》[71]的I.4數學研究的一般目的(第82頁)這篇文章中，作者對“不變式”的性質進行了很好的論述：對于一個不變式，時常尋求它的兩種主要性質，而這兩種性質又時常是向兩個相反方向起作用的。其一是要它盡可能的細，意思是只要兩個對象不等價，不變式就不同。其二是要能夠實際確定何時不變式不同。一個不變式哪怕是很細，如果無法算出來，那就沒有多大用處。所以，最強有力的不變式大概是哪些既能夠計算出來，又不太容易計算出來的不變式。然而，有時證明對象不等價是很困難的，以至于不變式盡管能部分時間有用，也認為是有用而且有趣的。以上論述也完全適用于本文討論的數據一致性條件，因此，應用時也理應同樣謹慎。

CT成像過程中由于被成像物體、系統幾何或掃描過程缺陷[22]，會導致CT重建圖像中表現出偽影[13,15]，進而影響CT圖像質量與后續診斷，因此，需要進行相應偽影校正。CT圖像偽影多種多樣，如部分體積偽影、射束硬化偽影[44-47,58]、運動偽影[23,40-43,52,54-56,64]、采樣偽影[37-39,62-63,69]、電子偽影、探測器余暉偽影、金屬偽影[45-46]、截斷偽影[33-36,49,70]、散射偽影[48,67-68]、風車偽影[14-18]、錐束偽影[24-26]、階梯偽影等。有必要對CT成像偽影產生機制開展深入研究，評估各種數據一致性條件對偽影校正的適應性。數據一致性條件也被用來改進重建算法的計算效率[22,61,65-66]?；跀祿恢滦詶l件的方法可廣泛應用于CT、PET、SPECT和MRI等多種成像模態。

另外，2012年Tang等對發散束CT各種數據一致性條件之間的邏輯關系采用易被科研人員理解的數學語言進行了梳理[72]，結合本文，相信對從事相關研究工作的科研人員會有一定幫助。2011年Wei 等對發散束投影與Radon變換之間提出了一系列的關系[73]，可涵蓋以前相當多的研究成果。都有著很好的參考價值。

4 結論

本文綜述了各種CT數據一致性條件，及其不同表現形式與應用場合。對數據一致性條件的綜合研究有助于從不同側面理解現代CT成像理論的內涵?；跀祿恢滦詶l件設計新偽影校正或圖像重建方法并評估其應用性能，將成為現代CT成像理論與應用研究的重要組成部分，值得業界對其進一步開展深入研究。