基于对称重复的模板校正锥束ct系统几何失真的方法
【技术领域】
[0001] 本发明属于图像处理领域里的计算机层析成像（Computed Tomography，CT)技术 领域，涉及利用具有对称重复的定标模板对锥束CT成像系统进行几何校正的方法。
【背景技术】
[0002] 计算机层析成像（Computed Tomography，CT)系统的几何校正问题是CT研究领域 所要解决的基本问题之一，工业CT系统作为一种新型的无损检测设备，适用于传统检测设 备无法达到的内部尺寸及特性检测。在工业CT系统中，射线源使用锥束射线源，锥束CT成 像系统中要求射线源中心射线穿过旋转轴垂直入射探测器平面的中心，并且旋转轴的投影 落在探测器的中心列上。但是，在工业CT系统的安装及后续使用中，难免会出现几何偏差 问题，由于射线源和旋转台的几何偏差转换成探测器平面的角度偏差和位移偏差，这些偏 差会使重建图像产生条状伪影、阴影伪影、环状伪影等多种伪影，影响图像重建的质量，因 此，在图像重建之前对成像系统进行几何偏差的校正至关重要。
[0003] 由于CT成像系统的几何偏差很难直接测量得到，寻找测量精度较高而且简便易 行的几何校正方法成为大家的目标。CT成像系统的几何校正方法根据是否使用定标模板分 为两大类。不需要使用定标模板的校正方法，直接利用投影数据构造一个代价函数，通过优 化求得这个代价函数的最优解，由于不使用定标模板的方法中，探测器的角度偏差和位移 偏差对代价函数的影响关系不同，求解的参数精度不同，因此还不是目前CT系统几何校正 所采用的主流方法。目前广泛采用的几何校正方法是基于定标模板的校正方法，需要设计 专门的定标模板，依赖模板的先验知识，采集单一角度的投影或在360度范围内采集多个 角度的投影，根据投影轨迹进行求解。目前现有的基于模板的校正方法主要存在以下问题： 假设某些条件是理想的或者忽略某些条件；仅对部分参数进行求解；如果选用简单的校正 模板，则需要采集多角度下的模板投影数据，才能满足求解要求；如果为了减少投影角度， 就必须要增加模板的复杂度，这样就会增加设计成本。
【发明内容】
[0004] 本发明的目的是提供一种基于对称重复的模板校正锥束CT系统几何失真的方 法，该方法利用对称重复的方格纹理作为定标模板，将探测器角度偏差与纹理图像的秩之 间建立联系，并将探测器角度偏差的求解问题等效为一定约束下纹理图像的秩最小化问 题，再对角度偏差进行修正后，利用特定的几何关系对探测器位移偏差进行求解。仅需要采 集单角度的投影数据，模板制作容易，测量方法简单有效，解决了需要采集多角度投影数据 与模板制作复杂的不足。
[0005] 本发明的技术方案是：
[0006] 一、技术方案的原理
[0007] (1)技术方案的基本概念
[0008] 模板：用于实际中测量CT成像系统几何偏差的物体，本发明中指方格纹理。
[0009] 投影：实际场景中物体经过X射线透射在探测器中所成的像。
[0010] 旋转台前进方向：旋转台向射线源移动的方向。
[0011] 源轴距：射线源到旋转轴的距离。
[0012] 源检距：射线源到探测器的距离。
[0013] (2)技术方案的具体原理：
[0014] CT成像系统的几何位置不同，对模板成像得到的投影也不同，本发明根据模板的 投影来确定CT成像系统的几何偏差，包括探测器的角度偏差、位移偏差，源轴距以及源检 距。本发明的基本原理为：首先，调节探测器平面平行于模板平面，当探测器平面与模板平 面不平行时，由于存在投影形变，该纹理图像在原有方向上的对称性和重复性不再保留，从 而不再具备低秩特性。将探测器角度偏差与纹理图像的秩之间建立联系，并将探测器角度 偏差的求解问题等效为一定约束下纹理图像的秩最小化问题，通过求解秩最小化问题，得 到CT成像系统中探测器的角度偏差。由于X射线源发射出的锥形束完全覆盖被测物体， 因此根据光通量相等的原则选取与旋转轴前进方向平行的射线作为中心射线，利用锥束CT 成像系统的特定几何成像关系求解出位移偏差、源检距以及源轴距。下面将对技术方案原 理进行详细介绍，首先介绍锥束CT几何成像系统的定义以及如何利用模板中图案的低秩 特性求解锥束CT成像系统中探测器的角度偏差，最后介绍求解探测器位移偏差、源检距以 及源轴距的求解。
[0015] 步骤一：锥束CT成像几何定义及低秩纹理求解探测器角度偏差
[0016] 如图1所示，射线源装置、探测器装置、转台装置都安装在操作平台上，由于设备 的良好精度与模板的机械制作方式，保证旋转台平面水平，模板垂直放置于旋转台上，并 且转台上的刻度与通孔保证放置时模板平面垂直于旋转台前进方向，模板的纵向中轴与旋 转轴重合。为定义锥束CT几何成像系统，如图2所示，点S为X射线源，建立世界坐标系 Ow-XwYwZw，旋转台上模板平面的中心点为O w，旋转轴为YjiJ，模板平面为X WYW平面，射线源 到XwYw平面的距离为源轴距R，Zdi通过右手定则确定；图像坐标系是以像素为单位的，坐 标原点为0, U轴平行于探测器平面的水平方向，V轴平行于探测器平面的竖直方向，坐标 (u，V)表示该像素在图像中的列数与行数；探测器成像平面坐标系中原点0。取在中心射线 在探测器上的投影位置，X轴和Y轴分别与图像坐标系的u轴和V轴平行，设0。点在图像坐 标系内的坐标为（u。，V。)，探测器平面中心的图像坐标为（u D &)。旋转台通过电机控制进行 Xw轴、Y w轴、Z w轴方向的位移调节，探测器通过螺旋钮实现绕X轴、Y轴、Z ￥轴角度的调节， 并进行三个方向的位移微调。
[0017] 理想系统要求探测器成像平面XY平行于XwYw平面，并且X射线源在世界坐标系的 坐标为（〇, 〇, R)，坐标（u。，V。）与坐标（uD vc)相同，经过人工搭建的锥束CT系统，难以满足 理想系统所需要的精度，存在几何位置偏差。如图3所示，这些偏差归结为探测器的位移偏 差屯火，探测器的角度偏差θ，β，η。其中，如图3(a)所示，探测器成像坐标系的中心0。 与（u D νε)相差的X轴距离为探测器位移偏差dx，探测器成像坐标系的中心0。与（uDv e)相 差的Y轴距离为探测器位移偏差dy;如图3(b)所示，探测器绕Zw轴旋转Tl ;如图3 (c)所 示，探测器绕X轴旋转Θ ;如图3(d)所示，探测器绕Y轴旋转β。
[0018] 空间中的一个点（xw，yw，zw)经过射线源S透射后，在探测器上得到的对应投影点 坐标为（u，v)，用齐次坐标表示两者的映射关系为：
[0020] 其中，s为一个常数因子，N为内参数矩阵，pJP p ￥分别表示探测器u轴方向和v 轴方向的像元尺寸，T是平移矩阵，r是一个3 X 3的旋转矩阵。根据欧拉定理当旋转矩阵用 欧拉角表示时，对于同一个变换，旋转顺序不同，欧拉角也不同，在指定欧拉角时首先约定 旋转顺序。在锥束CT几何成像系统中，欧拉角即为锥束CT成像系统中的探测器角度偏差。 在此，设定旋转顺序为：先绕Z w轴顺时针旋转，再绕Y轴顺时针旋转，最后绕X轴顺时针旋 转。因此，用探测器绕X轴、Y轴、和Zw轴的旋转角θ，β，η表示旋转矩阵 r为：
[0026] 如图4所示，当探测器平面与模板平面之间存在角度偏差时，得到的模板投影 不同于探测器平面平行于模板平面时所得到的模板投影，此时对世界坐标系中同一点 P (Xw，yw，Zw)进行成像得到的投影分别为点P2 (U2, V2)和点P1 (U1, V1)，则有
[0029] 为了之后的推导方便，将图像坐标系的原点建立在探测器成像平面坐标系的中 心，因此内参数矩阵N中的u。= 0，v。= 0,将公式（3)和公式（4)经过代换整理得到公式 (5)：
[0031] 其中 t' = [t' t' 2, t' 3] = [t/R, t2/R, t3/R]为归一化位移矢量。
[0032] 用符号。表示式（5)中（u2, V2)和（Ul，V1)之间运算，则两像素点之间的关系重写 为：
[0033] (u2, ν2) = τ 1 〇 (U1, V1) (6)
[0034] 即倾斜投影I与对正投影1°的关系表示为I 〇 τ = I °，其中转换函数τ为。运 算取决的参数有欧拉角θ，β，η和归一化的位移矢量t'i，t'2，t' 3，需要指出的是，位移矢 量t' i，t' 2, t' 3并不一定是真实的X射线源相对于世界坐标原点的归一化位移。这是因为 当探测器平面平行于模板平面时，只要模板在X射线覆盖范围之内、投影始终不超出探测 器范围，无论X射线源的位置如何变化，获得的模板投影都是低秩的，如图5中所示，X射线 源在不同位置所得到的模板的投影图是相似的，均为低秩图像。因此，式（6)求解τ中 的平移参数t' i，t' 2, t' 3并不唯一。因此，本发明仅通过低秩恢复求取欧拉角，位移求取在 步骤二会详细说明。通常，现实观察到的图像除去空间转换的干扰外，还会受到噪声或者遮 挡的影响，那么加入噪声变量E，则有
[0035] I 〇 τ = (I〇+E) (7)
[0036] 将上述问题转化为单目标优化问题为：
[0038] 其中I |e| |。表示噪声矩阵的范数，即求解出的矩阵中非零元素的个数，用来作为 噪声矩阵稀疏度的表示，常量γ为噪声矩阵稀疏度和平衡纹理矩阵秩两者的权值，在实验 中根据情况