计算机断层扫描(computed tomography，CT)是利用计算机技术对被检测物体断层扫描图像进行重建获得三维断层图像的扫描方式，其不受被检测物体材料、形状和表面状况等限制，可实现物体内部结构的无损观测，是一种先进的非接触式检测技术，最早应用于医学影像诊断学的检查，目前已逐步扩展至工业、农业和林业等领域(齐子诚等，2018；张京平等，2014；葛浙东等，2016a)。木材CT断层成像系统可重建木材三切面图像，真实反映裂纹、节子和生长轮等信息，对木材内部结构观测具有重要的科学意义和应用价值(曹正彬等，2017；葛浙东等，2016b)。获取木材断层图像主要分2步：首先，采集X射线穿过不同旋转角度下被检测木材的投影数据；然后，基于反投影图像重建算法生成木材扫描切面断层图像。由于投影数据的准确获取是图像重建的基础，所以为完成木材断层图像的正确重建，必须保持被检测木材在旋转过程中的稳定性，并且系统转台旋转中心(center of rotation, COR)应在射线源焦点与探测器中点的连线上；但在实际操作中，很难准确判断射线源焦点、旋转中心和探测器中点三者是否位于一条直线上，加之机械制造和装配误差的存在以及负载或测试环境的变化，转台在扫描过程中总是存在微小距离偏移，找准偏移位置费力费时。若仍以理论上的旋转中心为基准重建图像，则会引起重建图像畸变(王召巴，2001)，形成向周边扩散的环形伪迹，对后期基于图像的尺寸测量、面积测量、体积测量、密度测量、缺陷识别以及CT图像的精度都有严重影响(庞彦伟等，2001)，因此，在图像重建过程中必须进行旋转中心校正。

扇形束扫描和锥形束扫描是木材CT断层成像系统采用的主流扫描方式(戚玉涵等，2016)，旋转中心确定主要有几何法(Liu et al., 2006)、迭代法(Donath et al., 2006)等算法。几何法需要特定模型完成旋转中心定位，方法简单，但校正精度低；迭代法校正旋转中心偏移量受系统调焦设备机械精度影响，需要进行多次试探调节，很难满足现代各种CT检查设备对自动化程度和检测精度的要求。Azevedo等(1990)提出了一种利用投影正弦图修正投影中心的方法，但该方法对测量数据中的噪声较为敏感，灵敏度相对偏小，修正效果不明显。张俊等(2013)采用一种适用于锥束CT系统旋转平移(RT)轨迹的几何参数标定方法对存在偏差的系统进行仿真，结果发现该方法能有效抑制因几何参数误差造成的几何伪影，提高图像质量。傅健等(2015)基于样品投影质心公式与最小二乘正弦拟合的校正技术，精确求解出了系统旋转中心在每个角度的随机偏移量，可实现旋转中心随机偏移校正，抑制重建伪影。王敬雨等(2017)提出了一种高精度迭代标定方法，一次可标定锥束CT成像系统所有几何参数，并用校正后系统重建Shepp-Logan体模和PCB体模图像，重建图像效果良好，未见几何伪影。以上技术在医学或工业CT成像系统上获得了较好的校正效果，但木材是一种生物质材料，内部结构复杂，密度差异性较大，木材CT断层成像系统在图像重建过程中易受噪声和边界像素定位精度影响，几何参数标定繁琐复杂。

鉴于此，本研究从图像重建参数间的几何关系出发，推导出一种基于正弦图的转台旋转中心校正方法，并将其应用于木材CT断层成像系统，对木段横切面进行图像重建。通过设定不同旋转中心偏移量校正数值，定量计算、分析、修正扫描参数偏差，进而改善重建图像质量，以期为木材CT断层成像系统的稳定性和可靠性提供技术支持。

1 木材CT断层成像原理

木材CT断层成像系统通过X射线发射器发射出具有一定能量和穿透能力的射线，射线沿不同路径穿过物体后，其强度会不同程度衰减，X射线探测器记录沿不同路径衰减后的X射线能量，基于此计算衰减系数在木材断层上的二维分布矩阵，并将其转换为人眼可识别的灰度图像，即木材CT断层图像。

X射线源发射扇形射线束，射线穿过被检测物体后由等距平板探测器接收。如图 1所示，在扫描断层内建立固定不动的直角坐标系xOy和极坐标系r-θ 。E点是物体扫描断层平面上任意一点，用直角坐标表示其位置为(x, y)，用极坐标表示其位置为(r, θ)，即OE=r。过E点作EP垂直于中心射线束SO于P点。D₁D₂是实际探测器位置，C₀T₀是射线穿过E点后落在实际探测器上的探测点到探测器中点的距离。平行移动D₁D₂至穿过旋转中心位置，即D′₁D′₂为虚拟探测器，OF是射线穿过E点后落在虚拟探测器上的探测点到旋转中心的距离。φ是虚拟探测器不同时刻的旋转角度，即D′₁D′₂与x轴的夹角，φ∈(0，360]。令系统旋转中心O到射线源焦点S的距离SO为d₁，到探测器中点C₀的距离OC₀为d₂，即SO=d₁，OC₀=d₂，SC₀= d₁+ d₂，则：

$ PE = r\cos \left( {\theta - \varphi } \right); $

(1)

$ SP = SO - PO = {d_1} - r\sin \left( {\theta - \varphi } \right)。$

(2)

根据变量的几何关系，计算OF长度，用s′表示，则：

$ \frac{{PE}}{{OF}} = \frac{{SP}}{{SO}}, $

(3)

即

$ s' = \frac{{{d_1}r\cos \left( {\theta - \varphi } \right)}}{{{d_1} - r\sin \left( {\theta - \varphi } \right)}}。$

根据变量的几何关系，计算C₀T₀长度，用s表示，则：

$ \frac{{OF}}{{{C_0}{T_0}}} = \frac{{SO}}{{S{C_0}}}。$

(4)

即

$ s = \frac{{\left( {{d_1} + {d_2}} \right)r\cos \left( {\theta - \varphi } \right)}}{{{d_1} - r\sin \left( {\theta - \varphi } \right)}}。$

由此可知，穿过扫描断层内任意一点的投影地址s与旋转角度φ存在一定函数关系，投影地址轨迹在扫描过程中是一条与旋转角度相关的正弦曲线，等距射线扇形数据重建公式可表示为：

$ f\left( {x,y} \right) = \int_0^{2{\rm{ \mathsf{ π} }}} {P\left( {s,\theta } \right) \times h\left( s \right){\rm{d}}\theta } 。$

(5)

式中：f(x, y)为E点反投影数值；h(s)为滤波函数。

理想情况下，射线源焦点S、旋转中心O和探测器中点C₀位于一条直线上。图 2中实线对应的是xOy直角坐标系，此时系统旋转中心投影为探测器中点C₀，在图像重建时将根据理想旋转中心位置重建图像。实际情况下，由于机械加工精度和环境因素影响，CT系统转台旋转中心与理想位置很可能存在偏差，射线源焦点、旋转中心和探测器中点三者不一定共线，既可以认为旋转中心不在射线源焦点和探测器中点的连线上，也可以理解为射线源焦点不在旋转中心和探测器中点的连线上。为求出实际旋转中心与理想旋转中心的偏移量，将问题转换为求射线源焦点偏离理想焦点的偏移量，再通过几何关系计算系统旋转中心偏移量。

图 2中虚线对应的x′Oy′直角坐标系表示射线源焦点、旋转中心和探测器中点之间的关系。理想情况下，被检测物体上任意一点E经过射线源照射，探测器T₀点收到经过该点后衰减的X射线能量信息。而实际情况下，由于旋转中心偏移，同一位置点E经过射线源照射，探测器T点接收到经过该点后衰减的X射线能量信息。

xOy坐标系与x′Oy′坐标系之间具有如下映射关系：

$ \left[ {\begin{array}{*{20}{c}} x\\ y \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} {\cos \alpha }&{ - \sin \alpha }\\ {\sin \alpha }&{\cos \alpha } \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {x'}\\ {y'} \end{array}} \right]; $

(6)

同时，存在如下几何关系：

$ {\rm{sin}}\alpha = d/\sqrt {{d^2} + d_1^2} ; $

(7)

$ \cos \alpha = {d_1}/\sqrt {{d^2} + d_1^2} 。$

(8)

当且仅当实际旋转中心位于理想旋转中心位置，即α=0时，C与探测器中点C₀重合，过E点射线的投影地址T与T₀重合，此时重建投影地址正确。当实际旋转中心偏离理想旋转中心，即α≠0时，C与C₀、T与T₀均为不同的投影地址，此时重建图像投影地址错误，无法重建清晰的CT图像，必须对重建地址进行修正。

转台上被检测物体断层上任意一点E在原坐标系中投影地址为T₀，如图 2所示。当旋转中心偏移时，相当于焦点偏移量为d，新的投影地址由T₀变为T，T与T₀之间的距离为t。此时，点E的重建图像正确投影位置为T，为了获取T，有效方法是先计算出t。根据系统几何关系，有：

$ t = d \times \left( {{d_2} + y} \right)/\left( {{d_1} - y} \right)。$

(9)

令T₀与C₀之间的距离为t₀，T与C₀之间的距离为s，则点E正确投影地址T到探测器中点C₀的距离s可表示为：

$ s = {t_0} + t = {t_0} + d \times \left( {{d_2} + y} \right)/\left( {{d_1} - y} \right)。$

(10)

旋转中心校正方法是根据正弦图波峰波谷对称原理，测量E点最左端和最右端的投影地址，进而计算出投影地址对称中心与探测器中点的差值，即旋转中心偏移量。在图像重建过程中，对旋转中心偏移量进行补偿，即可获得任意一点正确的投影地址，成功重建质量良好的CT断层图像。

2 材料与方法 2.1 试验材料

铜丝，直径0.5 mm，长度100 mm。杉木(Cunninghamia lanceolata)木段采自浙江省杭州市临安区，直径205 mm，长度400 mm，气干密度0.37 g·cm^-3，树龄16年。木段表面无皮且光滑，有2处径向裂纹缺陷，从一端延伸至另一端，一端表面有心腐现象，另一端完好。

2.2 试验设备

实验室自行开发木材CT断层成像系统，如图 3所示，转台直径为500 mm，升降位移在-200 mm~+200 mm之间连续可调，旋转中心到X射线源和平板探测器中点的距离均为450 mm。X射线发射器工作电压220 V，最大连续输出功率500 W，射线管工作电压10~160 kV连续可调，稳压误差小于1%，工作电流0~3.12 mA连续可调，稳流误差小于2%。射线源发射扇形X射线束，最大射线张角80°，焦点尺寸0.05 mm，最大数据读出速率4 MB·s^-1。平板探测器有效检测区域512 mm，线性等距排列1 280个探测点，像元尺寸0.4 mm×0.4 mm，用k表示旋转中心偏移的探测点数量，t与k之间存在k=2.5t的关系。

2.3 研究方法

旋转中心校正前，将直径0.5 mm的铜丝竖直固定于转台上偏移旋转中心一定距离的位置。铜丝旋转1周的投影轨迹为正弦曲线，测量波峰、波谷对应的探测点编号，计算出转台旋转中心偏移量。根据CT断层成像原理，在反投影图像重建过程中对偏移量进行校正，重新获得铜丝正确投影轨迹，采集正弦曲线的波峰、波谷横坐标，直至将坐标对称中点坐标校正为640，则校正完毕。对杉木木段断层1进行扫描，设置不同旋转中心校正数值，如 k-2、k-1、k、k+1、k+2，重建断层1的5幅图像。对杉木木段断层2进行扫描，重建旋转中心校正前后断层图像，并对重建结果进行评价与分析。

3 结果与分析 3.1 正弦图旋转中心校正分析

投影正弦图轨迹表示木材扫描断层在360°旋转角度下的线积分投影数据，水平轴横坐标表示平板探测器接收穿过断层某点射线束的探测点编号n，垂直轴纵坐标表示转台旋转角度φ。正弦图可视为各角度下线积分投影一列列叠放起来的数据集，重建图像中任一像素点在不同旋转角度下的投影数据均对应正弦图中的某一条正弦曲线。如图 4所示，铜丝竖立于偏离转台旋转中心的某一位置，转台旋转中心校正前铜丝投影轨迹如图 4a所示，正弦图波峰、波谷坐标分别为(320, 10)、(940, 190)，对称中心横坐标地址为630。由于探测器中心地址为640，可知反映在探测器上的旋转中心偏移量为10个探测点，即k=10，转台旋转中心校正后的铜丝投影轨迹如图 4b所示。

3.2 断层图像质量分析

分别重建杉木木段断层1、2横切面图像，并对获得图像进行分析。图 5a为杉木木段整体图像，在断层1位置扫描，设置不同旋转中心校正数值，并重建断层1横切面图像。图 5b和c分别为系统转台旋转中心校正数值为k-2和k-1时的重建结果，即校正不完全的横切面断层图像。图 5d为系统转台旋转中心校正数值与偏移量相等时的重建结果。图 5e和f分别为系统转台旋转中心校正数值为k+1和k+2时的重建结果，即校正过量的横切面断层图像。可以看出，图 5b—f差异不明显，肉眼难以区分，木段上部均有1处大裂纹，周边分布多处细碎小裂纹且非常清晰，由此可知，5幅断层图像均可以反映木材内部结构特征，但需定量指标区分、评价重建图像质量。

为定量分析重建图像质量，引入均方误差(mean-square error，MSE)和峰值信噪比(peak signal to noise ratio，PSNR)来比较图像的质量和差异性(刘明娜等，2011)。MSE是反映重建图像和原始图像之间差异的一种度量，MSE越低，误差越小；而PSNR表示峰值误差的量度，PSNR越高，重建图像质量越好。计算PSNR之前，需先计算MSE(傅健等，2003)：

$ {\rm{MSE}} = \frac{{\sum\limits_{x = 1}^X {\sum\limits_{y = 1}^Y {{{\left[ {f\left( {x,y} \right) - \bar f\left( {x,y} \right)} \right]}^2}} } }}{{X \times Y}}。$

(11)

式中：f(x, y)表示原始图像像素矩阵；$\tilde f\left({x, y} \right) $表示重建图像像素矩阵；X和Y表示图像尺寸，本研究中重建的木材断层图像尺寸为1 024×1 024。

$ {\rm{PSNR}}\left( {x,y} \right) = 10 \times \lg \left( {\frac{{{\rm{Pea}}{{\rm{k}}^2}}}{{{\rm{MSE}}}}} \right)。$

(12)

式中：Peak为重建图像最大灰度值，本研究中重建的木材断层图像为灰度图像，Peak=255。

为进一步定量分析重建图像质量，选择结构相似度(structural similarity index measurement，SSIM)和相对误差(relative error，RE)2个指标(赵雨晴等，2017；李俊峰等，2014)对图像质量做更详尽的分析：

$ {\rm{SSIM}}\left( {x,y} \right) = \frac{{\left( {2{u_x}{u_y} + {c_1}} \right)\left( {2{\sigma _{xy}} + {c_2}} \right)}}{{\left( {u_x^2 + u_y^2 + {c_1}} \right)\left( {\sigma _x^2 + \sigma _y^2 + {c_2}} \right)}}; $

(13)

$ {\rm{RE}}\left( {x,y} \right) = \frac{{{{\left\| {x - y} \right\|}_2}}}{{{{\left\| x \right\|}_2}}} \times 100\% 。$

(14)

式中：x为原始图像；y为重建图像；u_x为x的平均值；u_y为y的平均值；σ_x²为x的方差；σ_y²为y的方差；σ_xy为x和y的协方差差；c₁=(k₁L)²、c₂=(k₂L)²为用来维持稳定的常数；L为像素灰度值的动态范围，本研究中L=255、k₁=0.01、k₂=0.03。

选择图 5中5幅木段扫描断层1的CT图像作为旋转中心校正前后的评价依据，由于图 5d为旋转中心校正数值与偏移量相等时的重建图像，所以定义其为原始图像，另外用图 5b、c、e和f 4幅CT图像与之比较，计算结果如表 1所示。可以看出，校正数值逐渐接近旋转中心偏移量，PSNR均增大且图像质量提高，SSIM增加且逐渐接近于1，RE均有所减小，图像愈加近似于原始图像。这表明，基于正弦图的系统旋转中心校正方法可有效改善图像质量。

为进一步分析旋转中心校正的重要性和必要性，重建杉木木段断层2横切面图像，同时对旋转中心校正前后的断层图像进行对比和分析。图 6a为系统旋转中心校正前对杉木木段断层2位置扫描后的断层图像重建结果，图像边缘有向周边扩散现象，环形伪影干扰严重，裂纹开口处黏连，心边材区分不明显，生长轮界限不清晰，多处小裂纹伴有假象，断层图像无法反映杉木内部真实结构信息，更无法用于测量与分析。图 6b为系统旋转中心校正后对杉木木段断层2位置扫描后的断层图像重建结果，图像边缘平滑、清晰，无向周边扩散现象，更无环形伪影干扰；心边材界限明显，心材生长轮清晰可见，边材生长轮重叠严重且难以区分；2处大裂纹位于木材上半部，左侧较窄且长裂纹由外向内延伸至木材髓心处，右侧较宽且短裂纹延伸至心材处。系统旋转中心校正后重建的断层图像质量良好，可以准确反映木材内部结构特征。这表明，基于正弦图的旋转中心校正方法在木材CT断层成像系统中具有良好的应用效果。

4 讨论

根据前文分析可知，穿过扫描断层内任意一点的投影地址与旋转角度存在一定函数关系，投影地址轨迹在扫描过程中是一条与旋转角度相关的正弦曲线。系统旋转中心校正是基于平面几何学中正弦图波峰、波谷对称的相关计算方法，将其应用于木材CT断层成像系统，在木材断层图像边界和细节保持方面获得了良好效果。该方法不足之处在于借助极细的铜丝作为旋转中心校正工具，使得点扩散效应下不止1个探测点接收到钢丝投影值，计算旋转中心偏移量的精度通常达到像素级，很难达到纳米CT成像系统的分辨率要求。但尽管如此，基于正弦图的系统旋转中心校正方法可有效解决木材CT成像系统旋转中心随机偏移问题，满足木材断层扫描成像需求，为木材内部缺陷检测和生长结构研究提供了关键技术支撑，工程应用性较强。

本研究更多的是对木材CT成像系统参数调教和几何运算方法的探讨，选择细铜丝的正弦图投影轨迹作为研究对象，有待利用更多投影数据信息进行多元算法分析，这也是下一步的研究重点。另外，基于正弦图的旋转中心校正方法对探测器短时间内性能不稳定所导致的较短且微弱的环状伪影还有一定不足，有待进一步改进和优化，以期构建更加精细、高分辨率和完整智能化的木材CT成像系统旋转中心校正方法。

5 结论

本研究提出一种基于正弦图的系统旋转中心校正方法，根据铜丝原始投影轨迹的平面几何关系，自动确定转台旋转中心位置，实现方法简单、便捷，具有一定的抗噪性，能适应旋转中心随机偏移情况，工程应用性强。通过试验结果可知，基于正弦图的系统旋转中心校正方法可有效解决旋转中心随机偏移问题，抑制环形伪影缺陷，明显提升木材CT断层图像的重建效果，具有可行性和有效性，对木材CT断层成像系统的调试具有较高应用价值。