【摘要】 針對在FDK算法的反投影過程中，各個體素在探測器上投影分布的特點，本文提出一種新的插值方法。該方法根據(jù)體素投影的特點，采用在重建過程中，根據(jù)其在不同掃描角度下在各個探測器單元上的投影所占面積的加權(quán)和作為反投影值。實際實驗結(jié)果表明，在FDK算法中這個新的插值方法比傳統(tǒng)的插值方法（如：最近鄰插值，雙線性插值）重建出來的圖像邊緣清晰，而且能更好地抑制噪聲。
【關(guān)鍵詞】 體素； 插值； FDK算法
A new interpolation method in the FDK algorithmZHAI Jing, PAN Jinxiao
(National Key Lab for Electronic Measurement and Technology, North University of China,
Taiyuan 030051, China)
Abstract: In the backprojection process of the FDK algorithm， every voxels project distribution in the detector have itself characteristics， This paper presents a new interpolation method. The voxel projection based on the characteristics used in the reconstruction process, According to its different scanning angle detector modules in all the occupied area of the projection of the size as the voxel projection. Actual experimental results show that the FDK algorithm of this new interpolation method give an better result in the reconstructed images verge than the traditional interpolation methods (such as replication interpolation, bilinear interpolation), and this new interpolation method can restrain noise effect.
Key words: voxel; interpolation; FDK algorithm
引言考試大
在FDK算法的反投影過程中，由于數(shù)據(jù)的離散性，會出現(xiàn)象素的投影地址“對不準”現(xiàn)象，一般需要插值運算。插值是指在己知的坐標范圍內(nèi)，一種基于模型的從離散數(shù)據(jù)估計連續(xù)數(shù)據(jù)的方法。經(jīng)典的線性插值技術(shù)包括最鄰近插值(replication)［1］，雙線性插值(bilinear）［2］，雙三次(Bicubic)插值［1，3～6］等。本文考慮了三維重建圖像的每個像素在不同角度時在探測器上的近似面積，（此近似面積大于1）并將在此面積的不同探測器上的投影值的加權(quán)和作為該像素在此掃描角度的投影值。
1 像素投影圖形的取法
通常情狀下，我們認為物體的像素在探測器上的投影是一個點，但在實際中，在某一個掃描角度下，經(jīng)過每個體素的射線在探測器上形成一個幾何圖形。在重建過程中，這個幾何面的形狀不規(guī)則而且其面積很不易求出，因此，我們可近似考慮每個像素的一部分點在探測器上所形成的圖形。具體步驟如下：
如圖1 所示，在錐束圓軌跡掃描結(jié)構(gòu)中，定義錐束的投影角為β，扇角為γ，錐角為κ。設(shè)ABCD—EFGH是要重建圖像的某一個體素(i,j,k)，如圖2所示，A： i－12,j+12,k+12，
B： i－12,j－12,k+12，
C： i+12,j－12,k+12，
D： i+12,j+12,k+12，
E： i－12,j+12,k－12，
F： i－12,j－12,k－12，
G： i+12,j－12,k－12，
H： i+12,j+12,k－12. 當射線源介于x正半軸到y(tǒng)正半軸之間這12π弧度（即0°≤β<90°）時，考慮ACGE這個對角面在探測器上的投影圖形。當射線源介于y正半軸到x負半軸之間這12π弧度（即90°≤β<180°）時，考慮BDHF這個對角面在探測器上的投影圖形。當射線源介于x負半軸到y(tǒng)負半軸之間這12π弧度（即180°≤β<270°）時，考慮ACGE這個對角面在探測器上的投影圖形。當射線源介于y負半軸到x正半軸之間這12π弧度（即270°≤β<360°）時，考慮BDHF這個對角面在探測器上的投影圖形。分別計算在上述條件下，其對角面的投影在探測器上分布的情狀。
2 像素反投影值的計算考試大考試大http:/
一般認為，像素的反投影值是由點的插值取得。有近鄰插值、雙線性插值等。在本文中，我們考慮上訴投影圖形的加權(quán)值。任取某一待重建體素上的頂點(x,y,z)，β為掃描角度，lso為探源到物體中心距離, lso2為探源到探測器距離,那么它在探測器上的落點p的坐標［7］: x′=(lso2/(lso-x×cos(β)+y×sin(β)))
·(x×sin(β)－y×cos(β)); (1)
y′=(lso2/(lso－x×cos(β)+y×sin(β)))×z; (2)因為物體離射線源距離比較遠而離探測器又很近，又由經(jīng)驗可知，當射線源的掃描角度是14π、34π、54π、74π時，投影的幾何面達到，大約就是一個寬是1、長是2的近似矩形，當射線源的掃描角度是12π、π、32π、2π時，投影的幾何面投影的幾何圖形的面達到最小，大約就是一個寬是1、長是1的近似矩形。如圖3所示，其在探測器上的分布共有六種情況，根據(jù)具體情況計算出該投影在每個探測器的探元中的面積，記為Si，設(shè)p′(xi,yi,β)是在掃描角度β時該面積上的投影值，(i,j,k,β)是在β時體素(i,j,k)要取的反投影值，N是投影面占據(jù)探元的個數(shù)，其值是6。我們得到計算像素反投影值的公式：(i,j,k,β)=Ni=1p′(xi,yi,β)*Si.(3)圖3 像素投影在探測器上的近似分布情況示意圖
3 實驗結(jié)果
實驗采用220 kv,10 mA的X射線源。探元的大小為0.127 mm。采用的探測器為 PAXSCAN 2520。工作模式：數(shù)字視頻。數(shù)據(jù)類型：unsigned short。A/D:12bit。射線源—標準件—探測器間距：850 mm～200 mm，旋轉(zhuǎn)一周采樣間隔為1度，某一角度下的投影如圖3示，大小為256*256，分別采用最近鄰插值，雙線性插值以及本文提到的新插值方法重建圖像，圖像大小為256*256。
由實驗結(jié)果表明：新的插值方法比我們通常用的最近鄰和雙線性插值法對圖像有明顯的改良，由圖像的灰度曲線也可表明新方法得到的圖像邊緣要好于前兩種方法。也能看出,由新插值得到圖像的邊緣比較清晰,而且還有抑制噪聲的效果。讓反投影的信息量的值由與周圍差別比較大的某一探元上的獲取改為在周圍的信息量更接近的相鄰幾個探元上綜合獲取，這樣能夠有效地抑制噪聲和孤立點，但是同時圖像的對比度也會有所下降。
4 總結(jié)
通過以上的分析和實驗表明，在基于圓軌跡的錐束CT掃描和重建過程中，如果考慮待重建體素落在探測器上一定的投影面積，并且在反投影重建這個體素點的時候，考慮那些探元上的信息量，將會對重建圖像的質(zhì)量有很大的提高，不僅圖像的邊緣更清晰，而且還能達到抑制噪聲的目的，但是圖像的對比度有所下降。當然，這種新的插值方法不僅適用于圓軌跡的掃描方，而且對所有的錐束掃描方式，并采用平板探測器采集數(shù)據(jù)的各種反投影的錐束重建都是適用的。
【參考文獻】
［1］Parker J A, Kenyon R V, Troxe L D E. Comparison of interpolating methods for image resampling. IEEE Transaction on Medical Imaging. 1983, 2(1)， 31-390.
［2］Jain A K. Fundaments of Digital Image Processing. Englewood Cliffs. NJ: PrenticeHall, 1989.
［3］Chen T C, deFigueiredo R J P. Twodimensional interpolation by generalized spline filters based on partial differential equation image models. IEEE Transaction on Acoustics, Speech,Signal Processing (ASSP). 1985, 33(3)，31-642.
［4］Hou H S, Andrews H C. Cubic splines for image interpolation and digital filtering. IEEE Transaction on Acoustics, Speech, Signal Processing (ASSP). 1978, 26(6)， 508-517.
［5］G Keys R. Cubic convolution interpolation for digital image processing. IEEE Transaction on Acoustics, Speech, Signal Processing (ASSP). 1981，29(6)，1153-11600.
［6］Unser M, Aldroubi A, Eden M. Fast bspline transforms for continuous image representation and1991 interpolation. IEEE Transaction on Pattern Analysis and Machine Intelligence (TPAMI), 13(3)，277-285.［7］孫怡,侯穎,胡家升.體積CT投影數(shù)據(jù)的模擬方法 ［J］. CT理論與應(yīng)用研究,2005,14(1-6).
【關(guān)鍵詞】 體素； 插值； FDK算法
A new interpolation method in the FDK algorithmZHAI Jing, PAN Jinxiao
(National Key Lab for Electronic Measurement and Technology, North University of China,
Taiyuan 030051, China)
Abstract: In the backprojection process of the FDK algorithm， every voxels project distribution in the detector have itself characteristics， This paper presents a new interpolation method. The voxel projection based on the characteristics used in the reconstruction process, According to its different scanning angle detector modules in all the occupied area of the projection of the size as the voxel projection. Actual experimental results show that the FDK algorithm of this new interpolation method give an better result in the reconstructed images verge than the traditional interpolation methods (such as replication interpolation, bilinear interpolation), and this new interpolation method can restrain noise effect.
Key words: voxel; interpolation; FDK algorithm
引言考試大
在FDK算法的反投影過程中，由于數(shù)據(jù)的離散性，會出現(xiàn)象素的投影地址“對不準”現(xiàn)象，一般需要插值運算。插值是指在己知的坐標范圍內(nèi)，一種基于模型的從離散數(shù)據(jù)估計連續(xù)數(shù)據(jù)的方法。經(jīng)典的線性插值技術(shù)包括最鄰近插值(replication)［1］，雙線性插值(bilinear）［2］，雙三次(Bicubic)插值［1，3～6］等。本文考慮了三維重建圖像的每個像素在不同角度時在探測器上的近似面積，（此近似面積大于1）并將在此面積的不同探測器上的投影值的加權(quán)和作為該像素在此掃描角度的投影值。
1 像素投影圖形的取法
通常情狀下，我們認為物體的像素在探測器上的投影是一個點，但在實際中，在某一個掃描角度下，經(jīng)過每個體素的射線在探測器上形成一個幾何圖形。在重建過程中，這個幾何面的形狀不規(guī)則而且其面積很不易求出，因此，我們可近似考慮每個像素的一部分點在探測器上所形成的圖形。具體步驟如下：
如圖1 所示，在錐束圓軌跡掃描結(jié)構(gòu)中，定義錐束的投影角為β，扇角為γ，錐角為κ。設(shè)ABCD—EFGH是要重建圖像的某一個體素(i,j,k)，如圖2所示，A： i－12,j+12,k+12，
B： i－12,j－12,k+12，
C： i+12,j－12,k+12，
D： i+12,j+12,k+12，
E： i－12,j+12,k－12，
F： i－12,j－12,k－12，
G： i+12,j－12,k－12，
H： i+12,j+12,k－12. 當射線源介于x正半軸到y(tǒng)正半軸之間這12π弧度（即0°≤β<90°）時，考慮ACGE這個對角面在探測器上的投影圖形。當射線源介于y正半軸到x負半軸之間這12π弧度（即90°≤β<180°）時，考慮BDHF這個對角面在探測器上的投影圖形。當射線源介于x負半軸到y(tǒng)負半軸之間這12π弧度（即180°≤β<270°）時，考慮ACGE這個對角面在探測器上的投影圖形。當射線源介于y負半軸到x正半軸之間這12π弧度（即270°≤β<360°）時，考慮BDHF這個對角面在探測器上的投影圖形。分別計算在上述條件下，其對角面的投影在探測器上分布的情狀。
2 像素反投影值的計算考試大考試大http:/
一般認為，像素的反投影值是由點的插值取得。有近鄰插值、雙線性插值等。在本文中，我們考慮上訴投影圖形的加權(quán)值。任取某一待重建體素上的頂點(x,y,z)，β為掃描角度，lso為探源到物體中心距離, lso2為探源到探測器距離,那么它在探測器上的落點p的坐標［7］: x′=(lso2/(lso-x×cos(β)+y×sin(β)))
·(x×sin(β)－y×cos(β)); (1)
y′=(lso2/(lso－x×cos(β)+y×sin(β)))×z; (2)因為物體離射線源距離比較遠而離探測器又很近，又由經(jīng)驗可知，當射線源的掃描角度是14π、34π、54π、74π時，投影的幾何面達到，大約就是一個寬是1、長是2的近似矩形，當射線源的掃描角度是12π、π、32π、2π時，投影的幾何面投影的幾何圖形的面達到最小，大約就是一個寬是1、長是1的近似矩形。如圖3所示，其在探測器上的分布共有六種情況，根據(jù)具體情況計算出該投影在每個探測器的探元中的面積，記為Si，設(shè)p′(xi,yi,β)是在掃描角度β時該面積上的投影值，(i,j,k,β)是在β時體素(i,j,k)要取的反投影值，N是投影面占據(jù)探元的個數(shù)，其值是6。我們得到計算像素反投影值的公式：(i,j,k,β)=Ni=1p′(xi,yi,β)*Si.(3)圖3 像素投影在探測器上的近似分布情況示意圖
3 實驗結(jié)果
實驗采用220 kv,10 mA的X射線源。探元的大小為0.127 mm。采用的探測器為 PAXSCAN 2520。工作模式：數(shù)字視頻。數(shù)據(jù)類型：unsigned short。A/D:12bit。射線源—標準件—探測器間距：850 mm～200 mm，旋轉(zhuǎn)一周采樣間隔為1度，某一角度下的投影如圖3示，大小為256*256，分別采用最近鄰插值，雙線性插值以及本文提到的新插值方法重建圖像，圖像大小為256*256。
由實驗結(jié)果表明：新的插值方法比我們通常用的最近鄰和雙線性插值法對圖像有明顯的改良，由圖像的灰度曲線也可表明新方法得到的圖像邊緣要好于前兩種方法。也能看出,由新插值得到圖像的邊緣比較清晰,而且還有抑制噪聲的效果。讓反投影的信息量的值由與周圍差別比較大的某一探元上的獲取改為在周圍的信息量更接近的相鄰幾個探元上綜合獲取，這樣能夠有效地抑制噪聲和孤立點，但是同時圖像的對比度也會有所下降。
4 總結(jié)
通過以上的分析和實驗表明，在基于圓軌跡的錐束CT掃描和重建過程中，如果考慮待重建體素落在探測器上一定的投影面積，并且在反投影重建這個體素點的時候，考慮那些探元上的信息量，將會對重建圖像的質(zhì)量有很大的提高，不僅圖像的邊緣更清晰，而且還能達到抑制噪聲的目的，但是圖像的對比度有所下降。當然，這種新的插值方法不僅適用于圓軌跡的掃描方，而且對所有的錐束掃描方式，并采用平板探測器采集數(shù)據(jù)的各種反投影的錐束重建都是適用的。
【參考文獻】
［1］Parker J A, Kenyon R V, Troxe L D E. Comparison of interpolating methods for image resampling. IEEE Transaction on Medical Imaging. 1983, 2(1)， 31-390.
［2］Jain A K. Fundaments of Digital Image Processing. Englewood Cliffs. NJ: PrenticeHall, 1989.
［3］Chen T C, deFigueiredo R J P. Twodimensional interpolation by generalized spline filters based on partial differential equation image models. IEEE Transaction on Acoustics, Speech,Signal Processing (ASSP). 1985, 33(3)，31-642.
［4］Hou H S, Andrews H C. Cubic splines for image interpolation and digital filtering. IEEE Transaction on Acoustics, Speech, Signal Processing (ASSP). 1978, 26(6)， 508-517.
［5］G Keys R. Cubic convolution interpolation for digital image processing. IEEE Transaction on Acoustics, Speech, Signal Processing (ASSP). 1981，29(6)，1153-11600.
［6］Unser M, Aldroubi A, Eden M. Fast bspline transforms for continuous image representation and1991 interpolation. IEEE Transaction on Pattern Analysis and Machine Intelligence (TPAMI), 13(3)，277-285.［7］孫怡,侯穎,胡家升.體積CT投影數(shù)據(jù)的模擬方法 ［J］. CT理論與應(yīng)用研究,2005,14(1-6).