投影旋轉中心（COR）精準定位是確保計算機斷層成像（CT）重建圖像質量的關鍵要素，經典的互相關匹配算法在投影角度為 0～180° 時難以滿足高質量 CT 成像要求，需進行改進創新。本文根據正弦圖上 0° 與 180° 投影數據翻轉后有對應性的特點提出基于這兩行數據平移匹配的 COR 校正算法，該算法利用 OTSU 進行閾值分割以減少背景噪聲影響，通過 L1 范數量化 COR 最小偏移得到準確校正值后進行 CT 重建。分別采用加入隨機梯度噪聲和高斯噪聲的 Sheep-Logan 模型和雄性 SD 大鼠樣本的同質肝臟與異質牙齒圖像驗證新算法的有效性，并將新算法與互相關匹配算法做性能對比。結果表明：新算法運算量少、簡便快速且具有良好的抗噪魯棒性，校正精度高（稀疏采樣投影數據在 10%～50% 時也能很好地校正 COR 值），CT 重建圖像質量有顯著改進，效果優于互相關算法。
引用本文： 趙琦, 趙雨晴, 叢長虹, 冀東江, 秦莉莉, 陳曉冬, 胡春紅. 基于正弦圖平移匹配的投影旋轉中心校正. 生物醫學工程學雜志, 2018, 35(4): 598-605. doi: 10.7507/1001-5515.201707065 復制
引言
X 射線計算機斷層成像（computed tomography，CT）是無損檢測領域的重要分支，并廣泛應用于航天、醫學影像等重要領域。在二維 CT 掃描系統中，投影旋轉中心（center of rotation，COR）的準確性將直接影響重建圖像的質量。但在實驗中，由于受到儀器精度、主觀誤差等因素的影響，COR 偏離探測器的中心，即沒有位于正弦圖（sinogram，sino 圖）的列中心，從而導致重建圖像產生嚴重的偽影[1]。因此，在 CT 重建前進行 COR 校正是必不可少的。
國內外的研究人員對基于 sino 圖的 COR 校正方法進行了研究。1990 年，Azevedo 等[2]提出通過計算投影數據質心的偏移得到 COR 的偏移獲得實際 COR 位置的方法，但該方法不能很好地處理衍射、折射和非線性噪聲干擾且受限于成像視野的大小。2006 年，Liu 等[3]對選擇兩對相近的且投影角度相差 180° 的數據求得插值系數，進而得到 COR 的對角插值法的性能進行研究，該方法對探測器性能與圖像質量要求高且掃描角度多（0～360°）。2014 年，Vo 等[4]提出了利用 sino 圖的傅里葉變換系數對稱性計算 COR 的方法，該方法對帶有隨機梯度噪聲的 sino 圖校正效果不佳。1974 年，Pratt[5]提出了基于互相關的匹配算法。隨著對互相關算法研究的不斷深入，2012 年 Yang 等[6]基于扇束模型將所得 sino 圖中相差 180° 的投影數據均進行互相關匹配來獲得實際 COR 位置，該方法具有一定的抗噪性且能夠準確校正 COR，已被廣泛關注和應用[7-10]。對于平行束成像模型，投影角度范圍通常是 0～180°，即只有一對相差 180° 的投影數據，應用互相關法計算受背景噪聲的影響較大，這將會影響 COR 的校正精度。因此，利用 0° 和 180° 投影數據進行高精度校正 COR 是改進創新的重要方向。
本文選用經典互相關匹配算法作為對比算法，根據 sino 圖中 0° 和 180° 投影數據方向相反數值相同的性質，提出了只對這兩行投影數據進行迭代匹配的校正算法，同時引入 OTSU 分割圖像以減少背景噪聲的干擾，精確 COR 偏移量，并通過L1 范數量化進一步提高 COR 的測量精度，并以此測量值作為 sino 圖列的中心進行校正，然后對校正的 sino 圖進行濾波反投影（filtered back projection，FBP）重建。最后，通過仿真和真實實驗對新算法進行驗證，并討論該算法的實際應用價值。
1 方法
如圖 1 所示，圖像坐標系為（x，y），旋轉坐標系為（t，s），其中t 為射線的位置，s 為射線的投影方向，虛線為部分射線，極坐標下的（r，φ）為圖像中的某一個像素點，θ 為旋轉坐標系的旋轉角度，即投影角度。
圖1
坐標系示意圖
Figure1.
Schematic diagram of the coordinate system
在旋轉坐標系上確定射線的投影地址，即
 |
將 r 與 φ 固定，隨著 θ 的變化，t 呈現出正弦曲線。定義投影角度
下射線經過樣品后的投影信息為 Pθ （θ 的取值范圍為 0～180°），那么將 0～180° 投影信息依次排列起來便構成了 sino 圖。將 sino 圖上投影角度 θ 下的投影地址從 1 開始依次向右升序排列，則式（1）等價為
 |
其中，c 為實際的 COR，X 為投影地址。
圖2
0° 與 180° 翻轉的投影地址示意圖
Figure2.
Schematic diagram of the 0° and horizontally flipped 180° projection addresses in sinogram
實際中，COR 偏移如圖 2 所示。其中，探測器的個數為 N，理論 COR 地址為（N+1）/2，
與
分別為投影信息 P0 和 P180 對應的實際 COR 地址，τ1 與 τ2 分別為 P0 和 P180 中同個探測器單元對應的投影地址到 c0 與理論 COR 地址的距離，X0 和 X180 分別為 P0 和 P180 的投影地址，方向相反且與實際 COR 地址的距離均為 τ1。將 P180 以理論 COR 地址為中心進行翻轉變換得到
，設其投影地址為
，則
與理論 COR 地址的距離仍為 τ2。
可表示為
 |
在圖 2 中，offset 表示 P0 和
中同個探測器單元在投影地址上的偏移。從圖 2 可以看出：
 |
 |
由式（4）和（5）可得實際 COR 地址：
 |
由式（6）可以看出，計算實際 COR 地址的關鍵在于得到準確的 offset 值。
由式（2）可得
 |
 |
由式（3）、（7）和（8）可得
 |
由式（6）和（9）可以看出，通過平移
直到與 X0 的投影地址一一對應，則此時對應的平移步長可以用來估計 offset 值，繼而得到實際 COR 地址。本文采用 L1 范數對偏移量 offset 進行量化，其定義如式（10）：
 |
在進行算法實現時，基于 X0 與
在相同投影地址的投影信息差值最小，可得
 |
與
分別為
與
翻轉的投影信息，最終得到 COR 的計算公式：
 |
背景噪聲會掩蓋 sino 圖的主要部分，使校正的結果存在較大誤差，本文在平移匹配之前，先采用 OTSU[11]分割背景噪聲并將其置零，從而減小噪聲影響，提高校正的準確性。以 SD 雄性大鼠牙齒樣本投影數據為例，圖 3 左列上圖為校正前的 sino 圖，其中紅線代表理論 COR，黃線代表實際 COR，可看出實際 COR 與理論 COR 存在偏移；圖 3 左列下圖代表校正前的 sino 圖中投影信息 P0 和
，黑箭頭所示為背景噪聲。圖 3 右列上圖代表校正后的 sino 圖，其采用插值補列方法使實際 COR 位于 sino 圖列中心以進行 COR 校正；圖 3 左列下圖代表校正后的 sino 圖中投影信息 P0 和
，可以看出，校正后的 P0 和
幾乎完全重合。綜上所述，本文提出新算法的具體步驟及流程圖（見圖 4）如下所示。
（1）提取 sino 圖中投影角度為 0° 和 180° 的投影信息 P0 和 P180；
（2）翻轉 P180 得到投影信息
；
（3）對 P0 和
進行 OTSU 處理，將背景噪聲置零；
（4）對 P0 和
利用 L1 范數量化投影信息最小偏移，即投影地址最佳匹配，然后利用式（12）求得實際 COR 地址，通過插值補列（零填充）的方法使實際 COR 作為新的 sino 圖的列中心以實現 COR 校正。
圖3
校正前后的牙齒樣本 0° 和 180° 投影信息
Figure3.
The 0° and 180° projection data of the dental tissues before and after correction
圖4
新算法流程圖
Figure4.
Flow chart of new algorithm
2 實驗設計
2.1 仿真實驗
仿真實驗采用 512 × 512 的 Sheep-Logan 頭模型，對其進行 radon 變換得到 sino 圖。在 sino 圖的左側或右側，以等差數列的形式（n=4，6，8，
）插入不同列數的零值來模擬 COR 向右或向左偏移，選用互相關算法作為對比算法。實際的 sino 圖含有因光源不穩定和光子散射產生的噪聲，因此，加入隨機梯度噪聲和高斯噪聲（見圖 5）分別進行模擬，并分別測試偏移量的大小和噪聲的強弱對計算 COR 準確性的影響。
圖5
加入隨機梯度噪聲和高斯噪聲的 sino 圖及利用互相關法和新算法校正 COR 后的典型重建圖像
Figure5.
The noisey sinogram with the random gradient and Gaussian noise and typical reconstructed images after correcting COR using cross-correlation method and new algorithm
表1
加入隨機梯度噪聲后不同算法得到的 COR
Table1.
The measurements of COR using different algorithms with the random gradient noises
對加入隨機梯度噪聲的 sino 圖分別使 COR 偏移 4、6、8、10 個像素，測量的 COR 結果如表 1 所示，隨著偏移量增大，新算法比互相關法測量的 COR 更準確，抗噪性更強。在 sino 圖中分別加入方差為 0.2、0.4、0.6、0.8、1.0、1.2 和 1.4 的高斯噪聲并偏移 12 個像素，用上述方法進行 COR 測量，如圖 6 所示。隨著噪聲方差的增大，特別是大于 0.4 時，噪聲干擾嚴重，互相關算法測得的 COR 明顯偏離真實值，而新算法測得的 COR 始終接近真實值，說明本文提出新算法的抗噪性強且校正準確。選擇偏移像素為 10 時和噪聲方差為 1 時兩種算法測得的 COR 分別進行 FBP 重建，并采用同一窗寬成像，如圖 5 所示，互相關法圖像均有明顯拖尾，而新算法無拖尾，圖像邊緣完整，也可直觀地說明新算法的校正效果更好。
圖6
加入不同方差的高斯噪聲用不同方法得到的 COR
Figure6.
The measurements of COR using different algorithms with different Gaussian noise levels
根據上述方法測得的 COR，采用 FBP 算法進行重建，并用相對誤差（relative error，RE）[12]、峰值信噪比（peak signal to noise ratio，PSNR）[12]、通用質量指數（universal quality index，UQI）[13]進行定量分析（見圖 7）。RE 是計算重建圖像與參考圖像誤差的指標，數值越小，則重建圖像越接近參考圖像，其定義為
 |
其中，
和
分別表示重建圖像與參考圖像（標準頭模型）。PSNR 是用于估計重建精度的誤差質量指標，數值越大，抗噪性越好，重建質量越好，其定義如下：
 |
 |
其中，頭模型圖像的灰度范圍是 0～255，Peak 為最大灰度值，即 255。UQI 是全參考圖像質量客觀評價方法之一，用于衡量圖像在相關度、亮度及對比度的失真等方面的質量，數值范圍是［–1，1］，數值越大，重建圖像越接近參考圖像，其定義如下：
 |
其中，Cov（x，y）表示
和
的協方差，
和
分別表示
和
的方差，ux 和 uy 分別表示
和
的平均值。
互相關法和新算法的定量分析結果如圖 7 所示，當噪聲方差大于 0.4 時，噪聲干擾嚴重。總體而言，新算法測得 COR 的重建圖像 RE 更小，PSNR 更高，UQI 更高，質量更好。
圖7
含有不同噪聲程度的重建圖像的質量指標
Figure7.
The quality metrics of the reconstructed images with different noise levels
此外，本研究對 COR 偏移 12 像素的 sino 圖進行稀疏采樣并測量 COR，將 COR 與未采樣時的位置相對應，從而探究本文所提方法在不同探測器采樣率下是否具有魯棒性。如圖 8 所示，探測器采樣率分別為 10%、13%、17%、20%、25%、33%、50% 和 100%，可以看出，新算法測得的 COR 基本接近真實值，具有很好的魯棒性。
圖8
不同探測器采樣率下測得的 COR
Figure8.
The measurements of COR for different sampling rates
2.2 真實實驗
針對基于同質性和異質性生物材料投影信息的不同，本文分別選用肝纖維化雄性 SD 大鼠肝臟（同質性）和健康雄性 SD 大鼠的牙齒（異質性）作為實驗樣本，樣品分別由首都醫科大學附屬北京友誼醫院動物研究中心和北京維通利華實驗動物技術有限公司提供，實驗在上海同步輻射裝置（Shanghai Synchrotron Radiation Facility，SSRF）的 BL13W1 線站進行，采用類同軸相襯成像方法，并利用基于相位恢復的 CT 重建方法獲得樣品斷層圖像。探測器分辨率為 9 μm，樣品與探測器之間的距離為 0.8 m；CT 掃描時，成像能量分別為 24 keV 和 33 keV，樣品在轉臺上分別以 0.35 °/s 和 0.15 °/s 的速度旋轉 180°，總投影數分別為 731 張和 959 張，投影圖像大小分別為 2 928 × 581 和 3 992 × 513，重建圖像大小分別為 1 025 × 1 675 和 720 × 630。
本文選取第 213 層相位恢復后的牙齒數據并得到 sino 圖，用互相關算法與新算法進行COR校正并進行 FBP 重建，采用同一窗寬進行成像。圖 9 右圖中，紅箭頭為牙釉質，白箭頭為牙髓腔，黃箭頭為牙槽骨、藍箭頭指示的黑色區域分布著牙周韌帶，均清晰可見。從綠框放大的圖像中看出，互相關法存在條狀偽影且圖像內容不清晰，空隙（紫箭頭所示）存在拖尾。本文算法校正后所得重建圖像清晰，無條狀偽影，且可以清晰地看到牙槽骨表面的空隙（紫箭頭所示）和軟組織（橙箭頭所示），可應用于正畸期間牙槽骨重塑的研究。
圖9
不同 COR 重建牙組織 CT 及感興趣區域放大圖像
Figure9.
The CT images of the dental tissues with different COR values and the corresponding magnified regions of interest
同樣，選取第 280 層相位恢復后的肝臟數據并得到 sino 圖，用上述兩種方法進行COR校正并重建，采用同一窗寬進行成像。圖 10 右圖中，可觀察到輪廓清晰、無拖尾的肝靜脈（紅箭頭所示），從綠框放大的圖像中也可以看出，互相關法存在拖尾且管腔不清晰，影響圖像質量，而新算法可清晰地觀察到小脈管，重建質量有很大提高，可應用于肝纖維化微脈管的研究。
圖10
不同 COR 重建肝臟組織 CT 及感興趣區域放大圖像
Figure10.
The CT images of the liver tissues with different COR values and the corresponding magnified regions of interest
對上述兩個樣品的 sino 圖進行稀疏采樣并測量 COR，探測器采樣率分別為 10%、13%、17%、20%、25%、33%、50%、100%，如圖 11 所示。其中，黑線和藍線分別為牙齒數據和肝臟數據測得的 COR，可以看出，探測器采樣率不同時，牙齒數據和肝臟數據測得的 COR 值變化幅度都在三個像素以內，說明本文方法具有很好的魯棒性，在投影數據不足的情況下（采樣率為 10%～50%）也能很好地校正。
圖11
不同探測器采樣率下測得的 COR
Figure11.
The measurements of COR wih the different sampling levels
3 結論
本文利用 sino 圖的 0° 和 180° 投影數據翻轉后具有對應性的特點，只對這兩行數據進行平移匹配，采用L1 范數量化 COR 最小偏移，并在該計算方法的基礎上引入 OTSU 分割背景噪聲，從而簡單快速且更準確地進行校正。仿真實驗結果表明，分別加入隨機梯度噪聲和不同噪聲強度的高斯噪聲，新算法均能準確校正 COR，優于互相關算法，具有較好的抗噪性。選擇不同采樣率模擬稀疏采樣計算 COR，新算法在 10%～50% 稀疏采樣投影數據不足時也能很好地校正 COR 值，優于互相關算法，具有較好的魯棒性。實際數據實驗結果表明，經新算法COR校正后，重建圖像能清晰地呈現出牙齒結構和肝臟管腔，無互相關法校正重建后的拖尾、條狀偽影，有效地提高了圖像質量，有助于臨床上對于牙槽骨重塑和肝纖維化脈管的研究。下一步工作是將該算法推廣到扇束和錐束投影模型，增大該算法的普適性。
引言
X 射線計算機斷層成像（computed tomography，CT）是無損檢測領域的重要分支，并廣泛應用于航天、醫學影像等重要領域。在二維 CT 掃描系統中，投影旋轉中心（center of rotation，COR）的準確性將直接影響重建圖像的質量。但在實驗中，由于受到儀器精度、主觀誤差等因素的影響，COR 偏離探測器的中心，即沒有位于正弦圖（sinogram，sino 圖）的列中心，從而導致重建圖像產生嚴重的偽影[1]。因此，在 CT 重建前進行 COR 校正是必不可少的。
國內外的研究人員對基于 sino 圖的 COR 校正方法進行了研究。1990 年，Azevedo 等[2]提出通過計算投影數據質心的偏移得到 COR 的偏移獲得實際 COR 位置的方法，但該方法不能很好地處理衍射、折射和非線性噪聲干擾且受限于成像視野的大小。2006 年，Liu 等[3]對選擇兩對相近的且投影角度相差 180° 的數據求得插值系數，進而得到 COR 的對角插值法的性能進行研究，該方法對探測器性能與圖像質量要求高且掃描角度多（0～360°）。2014 年，Vo 等[4]提出了利用 sino 圖的傅里葉變換系數對稱性計算 COR 的方法，該方法對帶有隨機梯度噪聲的 sino 圖校正效果不佳。1974 年，Pratt[5]提出了基于互相關的匹配算法。隨著對互相關算法研究的不斷深入，2012 年 Yang 等[6]基于扇束模型將所得 sino 圖中相差 180° 的投影數據均進行互相關匹配來獲得實際 COR 位置，該方法具有一定的抗噪性且能夠準確校正 COR，已被廣泛關注和應用[7-10]。對于平行束成像模型，投影角度范圍通常是 0～180°，即只有一對相差 180° 的投影數據，應用互相關法計算受背景噪聲的影響較大，這將會影響 COR 的校正精度。因此，利用 0° 和 180° 投影數據進行高精度校正 COR 是改進創新的重要方向。
本文選用經典互相關匹配算法作為對比算法，根據 sino 圖中 0° 和 180° 投影數據方向相反數值相同的性質，提出了只對這兩行投影數據進行迭代匹配的校正算法，同時引入 OTSU 分割圖像以減少背景噪聲的干擾，精確 COR 偏移量，并通過L1 范數量化進一步提高 COR 的測量精度，并以此測量值作為 sino 圖列的中心進行校正，然后對校正的 sino 圖進行濾波反投影（filtered back projection，FBP）重建。最后，通過仿真和真實實驗對新算法進行驗證，并討論該算法的實際應用價值。
1 方法
如圖 1 所示，圖像坐標系為（x，y），旋轉坐標系為（t，s），其中t 為射線的位置，s 為射線的投影方向，虛線為部分射線，極坐標下的（r，φ）為圖像中的某一個像素點，θ 為旋轉坐標系的旋轉角度，即投影角度。
圖1
坐標系示意圖
Figure1.
Schematic diagram of the coordinate system
在旋轉坐標系上確定射線的投影地址，即
|
將 r 與 φ 固定，隨著 θ 的變化，t 呈現出正弦曲線。定義投影角度
下射線經過樣品后的投影信息為 Pθ （θ 的取值范圍為 0～180°），那么將 0～180° 投影信息依次排列起來便構成了 sino 圖。將 sino 圖上投影角度 θ 下的投影地址從 1 開始依次向右升序排列，則式（1）等價為
|
其中，c 為實際的 COR，X 為投影地址。
圖2
0° 與 180° 翻轉的投影地址示意圖
Figure2.
Schematic diagram of the 0° and horizontally flipped 180° projection addresses in sinogram
實際中，COR 偏移如圖 2 所示。其中，探測器的個數為 N，理論 COR 地址為（N+1）/2，
與
分別為投影信息 P0 和 P180 對應的實際 COR 地址，τ1 與 τ2 分別為 P0 和 P180 中同個探測器單元對應的投影地址到 c0 與理論 COR 地址的距離，X0 和 X180 分別為 P0 和 P180 的投影地址，方向相反且與實際 COR 地址的距離均為 τ1。將 P180 以理論 COR 地址為中心進行翻轉變換得到
，設其投影地址為
，則
與理論 COR 地址的距離仍為 τ2。
可表示為
|
在圖 2 中，offset 表示 P0 和
中同個探測器單元在投影地址上的偏移。從圖 2 可以看出：
|
|
由式（4）和（5）可得實際 COR 地址：
|
由式（6）可以看出，計算實際 COR 地址的關鍵在于得到準確的 offset 值。
由式（2）可得
|
|
由式（3）、（7）和（8）可得
|
由式（6）和（9）可以看出，通過平移
直到與 X0 的投影地址一一對應，則此時對應的平移步長可以用來估計 offset 值，繼而得到實際 COR 地址。本文采用 L1 范數對偏移量 offset 進行量化，其定義如式（10）：
|
在進行算法實現時，基于 X0 與
在相同投影地址的投影信息差值最小，可得
|
與
分別為
與
翻轉的投影信息，最終得到 COR 的計算公式：
|
背景噪聲會掩蓋 sino 圖的主要部分，使校正的結果存在較大誤差，本文在平移匹配之前，先采用 OTSU[11]分割背景噪聲并將其置零，從而減小噪聲影響，提高校正的準確性。以 SD 雄性大鼠牙齒樣本投影數據為例，圖 3 左列上圖為校正前的 sino 圖，其中紅線代表理論 COR，黃線代表實際 COR，可看出實際 COR 與理論 COR 存在偏移；圖 3 左列下圖代表校正前的 sino 圖中投影信息 P0 和
，黑箭頭所示為背景噪聲。圖 3 右列上圖代表校正后的 sino 圖，其采用插值補列方法使實際 COR 位于 sino 圖列中心以進行 COR 校正；圖 3 左列下圖代表校正后的 sino 圖中投影信息 P0 和
，可以看出，校正后的 P0 和
幾乎完全重合。綜上所述，本文提出新算法的具體步驟及流程圖（見圖 4）如下所示。
（1）提取 sino 圖中投影角度為 0° 和 180° 的投影信息 P0 和 P180；
（2）翻轉 P180 得到投影信息
；
（3）對 P0 和
進行 OTSU 處理，將背景噪聲置零；
（4）對 P0 和
利用 L1 范數量化投影信息最小偏移，即投影地址最佳匹配，然后利用式（12）求得實際 COR 地址，通過插值補列（零填充）的方法使實際 COR 作為新的 sino 圖的列中心以實現 COR 校正。
圖3
校正前后的牙齒樣本 0° 和 180° 投影信息
Figure3.
The 0° and 180° projection data of the dental tissues before and after correction
圖4
新算法流程圖
Figure4.
Flow chart of new algorithm
2 實驗設計
2.1 仿真實驗
仿真實驗采用 512 × 512 的 Sheep-Logan 頭模型，對其進行 radon 變換得到 sino 圖。在 sino 圖的左側或右側，以等差數列的形式（n=4，6，8，
）插入不同列數的零值來模擬 COR 向右或向左偏移，選用互相關算法作為對比算法。實際的 sino 圖含有因光源不穩定和光子散射產生的噪聲，因此，加入隨機梯度噪聲和高斯噪聲（見圖 5）分別進行模擬，并分別測試偏移量的大小和噪聲的強弱對計算 COR 準確性的影響。
圖5
加入隨機梯度噪聲和高斯噪聲的 sino 圖及利用互相關法和新算法校正 COR 后的典型重建圖像
Figure5.
The noisey sinogram with the random gradient and Gaussian noise and typical reconstructed images after correcting COR using cross-correlation method and new algorithm
表1
加入隨機梯度噪聲后不同算法得到的 COR
Table1.
The measurements of COR using different algorithms with the random gradient noises
對加入隨機梯度噪聲的 sino 圖分別使 COR 偏移 4、6、8、10 個像素，測量的 COR 結果如表 1 所示，隨著偏移量增大，新算法比互相關法測量的 COR 更準確，抗噪性更強。在 sino 圖中分別加入方差為 0.2、0.4、0.6、0.8、1.0、1.2 和 1.4 的高斯噪聲并偏移 12 個像素，用上述方法進行 COR 測量，如圖 6 所示。隨著噪聲方差的增大，特別是大于 0.4 時，噪聲干擾嚴重，互相關算法測得的 COR 明顯偏離真實值，而新算法測得的 COR 始終接近真實值，說明本文提出新算法的抗噪性強且校正準確。選擇偏移像素為 10 時和噪聲方差為 1 時兩種算法測得的 COR 分別進行 FBP 重建，并采用同一窗寬成像，如圖 5 所示，互相關法圖像均有明顯拖尾，而新算法無拖尾，圖像邊緣完整，也可直觀地說明新算法的校正效果更好。
圖6
加入不同方差的高斯噪聲用不同方法得到的 COR
Figure6.
The measurements of COR using different algorithms with different Gaussian noise levels
根據上述方法測得的 COR，采用 FBP 算法進行重建，并用相對誤差（relative error，RE）[12]、峰值信噪比（peak signal to noise ratio，PSNR）[12]、通用質量指數（universal quality index，UQI）[13]進行定量分析（見圖 7）。RE 是計算重建圖像與參考圖像誤差的指標，數值越小，則重建圖像越接近參考圖像，其定義為
|
其中，
和
分別表示重建圖像與參考圖像（標準頭模型）。PSNR 是用于估計重建精度的誤差質量指標，數值越大，抗噪性越好，重建質量越好，其定義如下：
|
|
其中，頭模型圖像的灰度范圍是 0～255，Peak 為最大灰度值，即 255。UQI 是全參考圖像質量客觀評價方法之一，用于衡量圖像在相關度、亮度及對比度的失真等方面的質量，數值范圍是［–1，1］，數值越大，重建圖像越接近參考圖像，其定義如下：
|
其中，Cov（x，y）表示
和
的協方差，
和
分別表示
和
的方差，ux 和 uy 分別表示
和
的平均值。
互相關法和新算法的定量分析結果如圖 7 所示，當噪聲方差大于 0.4 時，噪聲干擾嚴重。總體而言，新算法測得 COR 的重建圖像 RE 更小，PSNR 更高，UQI 更高，質量更好。
圖7
含有不同噪聲程度的重建圖像的質量指標
Figure7.
The quality metrics of the reconstructed images with different noise levels
此外，本研究對 COR 偏移 12 像素的 sino 圖進行稀疏采樣并測量 COR，將 COR 與未采樣時的位置相對應，從而探究本文所提方法在不同探測器采樣率下是否具有魯棒性。如圖 8 所示，探測器采樣率分別為 10%、13%、17%、20%、25%、33%、50% 和 100%，可以看出，新算法測得的 COR 基本接近真實值，具有很好的魯棒性。
圖8
不同探測器采樣率下測得的 COR
Figure8.
The measurements of COR for different sampling rates
2.2 真實實驗
針對基于同質性和異質性生物材料投影信息的不同，本文分別選用肝纖維化雄性 SD 大鼠肝臟（同質性）和健康雄性 SD 大鼠的牙齒（異質性）作為實驗樣本，樣品分別由首都醫科大學附屬北京友誼醫院動物研究中心和北京維通利華實驗動物技術有限公司提供，實驗在上海同步輻射裝置（Shanghai Synchrotron Radiation Facility，SSRF）的 BL13W1 線站進行，采用類同軸相襯成像方法，并利用基于相位恢復的 CT 重建方法獲得樣品斷層圖像。探測器分辨率為 9 μm，樣品與探測器之間的距離為 0.8 m；CT 掃描時，成像能量分別為 24 keV 和 33 keV，樣品在轉臺上分別以 0.35 °/s 和 0.15 °/s 的速度旋轉 180°，總投影數分別為 731 張和 959 張，投影圖像大小分別為 2 928 × 581 和 3 992 × 513，重建圖像大小分別為 1 025 × 1 675 和 720 × 630。
本文選取第 213 層相位恢復后的牙齒數據并得到 sino 圖，用互相關算法與新算法進行COR校正并進行 FBP 重建，采用同一窗寬進行成像。圖 9 右圖中，紅箭頭為牙釉質，白箭頭為牙髓腔，黃箭頭為牙槽骨、藍箭頭指示的黑色區域分布著牙周韌帶，均清晰可見。從綠框放大的圖像中看出，互相關法存在條狀偽影且圖像內容不清晰，空隙（紫箭頭所示）存在拖尾。本文算法校正后所得重建圖像清晰，無條狀偽影，且可以清晰地看到牙槽骨表面的空隙（紫箭頭所示）和軟組織（橙箭頭所示），可應用于正畸期間牙槽骨重塑的研究。
圖9
不同 COR 重建牙組織 CT 及感興趣區域放大圖像
Figure9.
The CT images of the dental tissues with different COR values and the corresponding magnified regions of interest
同樣，選取第 280 層相位恢復后的肝臟數據并得到 sino 圖，用上述兩種方法進行COR校正并重建，采用同一窗寬進行成像。圖 10 右圖中，可觀察到輪廓清晰、無拖尾的肝靜脈（紅箭頭所示），從綠框放大的圖像中也可以看出，互相關法存在拖尾且管腔不清晰，影響圖像質量，而新算法可清晰地觀察到小脈管，重建質量有很大提高，可應用于肝纖維化微脈管的研究。
圖10
不同 COR 重建肝臟組織 CT 及感興趣區域放大圖像
Figure10.
The CT images of the liver tissues with different COR values and the corresponding magnified regions of interest
對上述兩個樣品的 sino 圖進行稀疏采樣并測量 COR，探測器采樣率分別為 10%、13%、17%、20%、25%、33%、50%、100%，如圖 11 所示。其中，黑線和藍線分別為牙齒數據和肝臟數據測得的 COR，可以看出，探測器采樣率不同時，牙齒數據和肝臟數據測得的 COR 值變化幅度都在三個像素以內，說明本文方法具有很好的魯棒性，在投影數據不足的情況下（采樣率為 10%～50%）也能很好地校正。
圖11
不同探測器采樣率下測得的 COR
Figure11.
The measurements of COR wih the different sampling levels
3 結論
本文利用 sino 圖的 0° 和 180° 投影數據翻轉后具有對應性的特點，只對這兩行數據進行平移匹配，采用L1 范數量化 COR 最小偏移，并在該計算方法的基礎上引入 OTSU 分割背景噪聲，從而簡單快速且更準確地進行校正。仿真實驗結果表明，分別加入隨機梯度噪聲和不同噪聲強度的高斯噪聲，新算法均能準確校正 COR，優于互相關算法，具有較好的抗噪性。選擇不同采樣率模擬稀疏采樣計算 COR，新算法在 10%～50% 稀疏采樣投影數據不足時也能很好地校正 COR 值，優于互相關算法，具有較好的魯棒性。實際數據實驗結果表明，經新算法COR校正后，重建圖像能清晰地呈現出牙齒結構和肝臟管腔，無互相關法校正重建后的拖尾、條狀偽影，有效地提高了圖像質量，有助于臨床上對于牙槽骨重塑和肝纖維化脈管的研究。下一步工作是將該算法推廣到扇束和錐束投影模型，增大該算法的普適性。
