基于三維運動軌跡的高分辨顯微CT環狀偽影去除方法

李靜，鄭良，陳璐杰，羅守華

東南大學 生物科學與醫學工程學院，江蘇 南京 210000

引言

高分辨顯微CT可達到微米級別分辨率，能對骨骼、牙齒、生物標本等離體標本進行高分辨的三維成像[1]。為了實現系統的高空間分辨率，高分辨顯微CT系統大多采用了X線耦合可見光系統高分辨探測器[2]，該探測器的核心部件是將X射線轉化為可見光的超薄閃爍體，當系統的成像分辨率達到2 μm時，閃爍體的厚度需要小于10 μm才可在景深范圍內確保清晰成像。如此薄的閃爍體，其研制過程中不可避免地會引入損壞或響應不一致的像素點，在圓軌跡重建下，最終的重建圖像會出現環狀偽影[3]。環狀偽影在重建圖中的表現為以旋轉中心為圓心的同心圓（圖1a），在正弦圖中為豎直方向的直線(圖1b)，嚴重的環狀偽影會影響圖像的后續處理及量化分析。

圖1 環狀偽影在高分辨顯微CT重建圖及正弦圖中的表現注：a. 牙簽重建圖像；b. 相應的正弦圖。

目前的環狀偽影校正方法主要分為三類：第一類是硬件校正方法[4]，其中最為常見的方法是通過移動探測器破壞環狀偽影的形成機制，使得探測器像素壞點反投影后的誤差分散到重建圖像各處，使其無法形成環狀；第二類是投影域校正方法[5-6]，Münch等[7]提出對投影數據進行離散傅里葉變換或小波變換，Liu等[8]也提出可利用低通濾波器濾除水平方向的高頻分量來抑制環狀偽影，Sun等[9]通過迭代算法對X射線波譜進行修正從而去除環狀偽影等非線性偽影；第三類是重建域校正方法，常用的算法是使用極坐標變換法[10-11]將直角坐標系中的環狀偽影映射為極坐標系中的線狀偽影，再使用一系列線狀偽影去除算法對其進行檢測以及校正，再將圖像變換到直角坐標系[12-13]。

這些環狀偽影抑制方法雖然都能一定程度上能夠去除環狀偽影，但是也相應地存在一些不足。投影域算法中的低通濾波器法會過濾掉部分有效信號，其次迭代類型算法[14]，大都存在運行時間過長、不適用于實際等問題；重建域校正算法里的極坐標變換法由于多次使用插值操作，造成圖像分辨率下降；而移動探測器的硬件校正方法圖像質量好、重建時間快，不過對硬件精度要求高。圖2所示為本實驗室自研高分辨顯微CT系統，其具有光學放大裝置以及X線耦合可見光成像系統（虛線框所示），其中X線耦合可見光成像系統內包含了易引入環狀偽影的閃爍晶體。因此移動探測器成本高且非常容易影響圖像質量，移動探測器的環狀偽影硬件校正方法在高分辨顯微CT系統并不適用。

圖2 高分辨顯微CT實體圖

針對以上方法的不足并結合高分辨顯微CT系統的結構特點，本文提出了一種基于三維運動軌跡的高分辨顯微CT環狀偽影去除方法，該方法主要包括兩個部分：三維運動軌跡的位移矩陣測量以及結合位移矩陣的FDK重建。位移矩陣測量需要對圖像進行配準，我們采用了Insight Segmentation and Registration Toolkit （ITK）平臺。ITK是一個開源、跨平臺的醫學圖像處理類庫，主要用于醫學圖像的預處理、分割及配準。基于ITK的配準算法簡化了醫學圖像配準程序的開發工作，為醫學圖像引導臨床診斷的研究應用提供參考，是一種成熟的圖像配準平臺。為了增強該方法的實用性，本文中選用了FDK算法進行重建。FDK算法是由Feldkamp、Davis和Kress提出的一種基于圓形軌道掃描的錐束CT近似解析重建算法，具有重建速度快、可以局部重建的優點，是一種應用廣泛的實用型算法。本文提出的環狀偽影去除方法對環狀偽影去除效果好、重建速度快，適合于高分辨顯微CT系統。

1 算法

1.1 三維運動軌跡去環原理

傳統的高分辨顯微CT中待測物體的運動軌跡為轉臺形成的圓軌跡，重建過程中環狀偽影的形成機制如圖3所示[15]。若單個探測器像素存在響應不一致，此誤差值在反投影后會平均分布在探測器像元與射線源焦點的連線上（如圖3中線L所示）。高分辨顯微CT成像過程中探測器與射線源相對靜止，轉臺圍繞著旋轉中心旋轉采圖，這意味著每一個角度的探測器像元與射線源焦點的連線到旋轉中心的距離始終保持不變。反投影后的誤差數據經圓周累加后會出現圖3的所示的情況，偽影較為稀疏的區域會被探測器其他像元測量得到的正確數據修正，而較為密集的區域則在重建圖像中切出一個以旋轉中心為原點、d為半徑的環狀偽影。

圖3 成環機制

本文中加入了可移動的三維平移臺并結合旋轉臺，在空間中形成全新的三維運動軌跡從而破壞圓軌跡的成環機制。在平移臺移動過程中，旋轉臺、源以及探測器的位置均沒有發生改變。圖4a～c是平移臺移動不同位置的重建圖像，移動位置分別為左上，中間，右上。圖中可以看出重建圖像的旋轉中心以及環狀偽影的成像位置沒有改變，變化的是物體在重建圖中的位置。圖4d指的是將a～c三張圖像針對成像物體進行配準后疊加的結果，圖4e是將成像物體去除之后環狀偽影誤差線的疊加結果，可以看出誤差線被分散到了圖像的各處，環狀偽影的成像機制被破壞了。偽影數據被打亂并分布到重建圖像各處。對于某個重建數據點而言在特定角度下獲得的是偽影數據，但是在其他角度下得到的是真實數據，該重建數據點的偽影數據被真實數據平均，因此不會出現某個數據點的所有角度都是偽影數據而導致此點的重建值誤差異常突出的情況。由上述可知通過平移臺移動多位置并完整拍攝n圈圖像，將他們配準后疊加可以得到去環的效果。但是此種方法掃描時間長，且每一次重建都要進行配準操作，導致重建過程過于復雜而缺乏足夠的實用性。

圖4 移動平移臺去環原理

圖5a中將CT的圓周掃描軌跡等角度間隔劃分，紅點和黃點表示只有旋轉臺旋轉待測物體時產生的圓周軌跡。圖5b中是使用三維平移臺以及旋轉臺產生的三維運動軌跡。黑點表示待測物體只有旋轉，灰點表示待測物體平移之后旋轉，圖5b中可以看出物體平移之后旋轉中心沒有變化但是旋轉半徑發生改變，因此形成了圖中所示的待測物體質心的三維運動軌跡。高分辨顯微CT拍攝的一圈投影圖中屬于不同平移臺位置的投影圖間錯排列，我們將平移臺在拍攝過程中需要在不同位置之間不斷抖動的采圖方式稱為單圈抖動拍圖，本文中使用了單圈抖動的采圖方式來去除環狀偽影。

圖5 單圈抖動拍圖注：a. 圓軌跡；b. 三維運動軌跡。

基于三維運動軌跡的高分辨顯微CT環狀偽影去除方法具體流程為：首先使用重復定位精度優秀的三維平移臺移動到多個固定位置拍攝模體，通過三維圖像配準算法進行位移矩陣測量。獲得位移矩陣之后，使用單圈抖動來拍攝待測物體。最后，利用本文提出的基于位移矩陣的FDK算法進行重建，達到去除環狀偽影的效果，流程圖見圖6。

1.2 三維亞像素剛性配準算法實現

本文中三維圖像的亞像素配準基于ITK實現[16]，配準的過程主要包括四個模塊：變換模型、插值方式、相似度計算、優化方法選擇。由于本文中使用了重復定位精度優秀（0.05% F.S）的三維平移臺進行平移臺移動，平移過程中待拍攝物體不存在非剛性變換以及旋轉，因此配準只需要考慮三維剛性平移變換。此外插值方法是三線性插值，相似度計算為均方誤差方程，優化方法選取正則化步長梯度下降法。基于ITK的配準實現偽代碼如下所示：① 設置Transform、Interpolator（插值方式）、metric（相似度計算）、optimizer（優化方法）類型；② 設置fixedImage、MovingImage；③ 設置迭代次數S、最大迭代步長、最小迭代步長等迭代參數；④ 初始化參數；⑤ 重復迭代；⑥ For i = 1 to S do；⑦ 將MovingImage結合相應的位移矩陣、插值方式進行轉換得到TransferImage；⑧ 計算TransferImage與fixedImage之間的metric；⑨ 對位移矩陣進行迭代優化；⑩ End For。

圖6 基于三維運動軌跡的環狀偽影去除方法流程圖

1.3 基于位移矩陣的FDK重建算法

本文中基于位移矩陣的FDK反投影使用的是像素驅動方法，將每個重建像素坐標投影到投影圖上，獲得對應投影圖坐標下的值后累加到該重建像素值上。由于在單圈抖動拍攝的一圈投影圖中每一個角度對應的平移臺位置不同，因此需要在反投影時將每個位置位移矩陣考慮到FDK算法內，使得實際物體上同一個點在不同平移臺位置時的投影值都能累加到同一個CT值上。如圖7所示，以源到探測器平面的垂線與探測器平面的交點為坐標原點，源與坐標原點的連線為y軸，旋轉軸為z軸。α為投影線與y軸之間的夾角，焦距為D。

確定源、探測器、旋轉中心后，通過旋轉物體來獲得不同角度的投影值，基于位移矩陣的FDK算法主要有三步：

第一步，對物體旋轉角為β的投影數據進行預處理，即投影值乘以cosα：

第二步，對不同投影角度的投影數據進行水平方向上的一維濾波，g(x)為一維斜坡濾波器：

第三步，最后沿X射線方向進行三維反投影。我們定義旋轉以及平移的操作為M得到：

得到某點坐標之后，我們需要將經過M變換之后的點（xM，yM，zM）做正投影到探測器上，得到此點的投影點坐標對應的上的值，將這個操作定義為B：

重建的CT值是各角度下通過該體素的所有角度投影的貢獻之和，式(5)中N為投影角度數，因此我們得到：

圖7 錐束的坐標系

2 實驗及結果分析

基于三維運動軌跡的環狀偽影去除算法的驗證實驗設計，本文采用一組仿真數據與一組離體小鼠肱骨真實數據驗證算法的有效性。本文中的小鼠產地為南京市江寧區青龍山動物繁殖場，體重為24～26 g，取其肱骨進行離體拍攝實驗。實驗中所用真實數據由實驗室已有的高分辨型顯微CT拍攝獲得。重建參數，見表1。

表1 重建參數

2.1 仿真結果

使用標準的三維shepp-logan模型進行了平移臺移動仿真實驗，平移臺一共選取了四個位置，分別為P0(0,0,0)、P1(γ,-γ,γ)、P2(-γ,γ,-γ)、P3(γ,-γ,-γ)，其中代表步長。本次仿真實驗分為四組子實驗a、b、c、d，去除環狀偽影后的重建圖像對應圖8a～d，圖8e～h分別為圖8a～d的局部放大圖像。實驗中，分別設置a、b、c、d的步長分別為0、1、5、10。當等于0時等同于平移臺沒有移動，如圖8a、e所示重建圖像帶有明顯環狀偽影；當等于1、5的時候，如圖8b、c、f、g所示重建圖像環狀偽影去除效果明顯，且不影響圖像質量；當等于10時，如圖8d、h所示，雖然環狀偽影也被消除了，但是重建圖像出現了一定程度的條紋狀偽影。實驗說明了位移大小應當盡量控制在10以內，位移過大，會產生影響圖像質量的偽影。

圖8 仿真結果注：a～d為移動步長r分別為0、1、5、10的重建圖像；e～h為a～d的局部放大圖像。

2.2 實際結果

本文對小鼠肱骨進行了高分辨顯微CT成像。使用重復定位精度良好的三維平移臺進行固定選取了四個位置，經過亞像素三維配準算法計算，得到這四個位置分別為(0.00,0.00,0.00)、(3.23,4.23,0.10)、(-4.10,2.30,-0.34)、(4.89,2.01,1.40)。圖9a為平移臺沒有移動的重建結果，可以看到環狀偽影明顯。圖9b為基于三維運動軌跡的重建結果，環狀偽影得到了很好的去除，圖像中沒有出現因為位移而產生其他偽影。實際實驗結果明顯可見環狀偽影減少、去環效果顯著，小鼠肱骨本身結構像也符合其微結構特征，因此論證了此環狀偽影去除方法在實際實驗中的有效性以及可重復性。

圖9 小鼠肱骨圖像注：a. 環狀偽影；b. 基于三維運動軌跡去除環狀偽影圖像。

3 結論

本文提出了一種基于三維運動軌跡的高分辨顯微CT環狀偽影去除方法。基于位移矩陣的FDK重建算法實現過程也可以采用先反投影，再平移的方式。但是，這種方法存在兩個缺點，首先此方法需要額外的內存空間用于保存反投影后的臨時重建數據，特別是對于高分辨顯微CT這種數據量大的系統而言，內存消耗過大。其次平移的過程中多做了插值操作，損失了分辨率。仿真實驗探討了位移大小對于圖像質量的影響，結果證明基于位移矩陣的FDK重建算法在小位移的情況下，對圖像質量基本沒有影響。由仿真實驗以及真實實驗結果可知，本文提出的基于三維運動軌跡的高分辨顯微CT環狀偽影去除方法在不額外增加采圖時間的條件下，對于環狀偽影去除效果優秀、重建速度快，對于具有復雜成像系統的透鏡耦合式高分辨CT非常適用。