一种X射线超分辨成像方法及其应用与流程

本发明属于x射线成像技术领域，具体涉及一种x射线超分辨成像方法及其应用。

背景技术：

x射线透射成像系统是无损检测物体的一大类方法，随着x射线透射成像系统在各个领域应用的日渐深入，对x射线透射成像系统的空间分辨率自然会提出更高的要求，比如对于细微结构的成像。

就目前技术而言，提高x射线透射成像系统的分辨率可以通过获得更小的探测器像素大小，更薄的闪烁体，更小的x射线源焦点等这些硬件条件来实现，也可以调整系统参数，比如放大倍数。然而，优秀的微焦点x射线源造价过于昂贵、闪烁体加工困难，容易产生响应不一致的像素点，在ct图像域中引入环状伪影、而探测器硬件工艺存在极限，过小的像素大小也会使得每个像素获得的光子数量减少而牺牲探测器其它性能指标，比如信噪比；此外，放大倍数也会受到系统模糊，焦点大小的限制而不能进行无限放大。因而如何进一步提高x射线透射成像系统图像的空间分辨率，提高图像质量成为了本领域的一大技术难点。

正如我们上述提到，在现有条件下得到更高的分辨率，这是本领域的一大技术难点和技术挑战。1984年，通过利用存在亚像素移动的原始分辨率图像序列，tsai和huang的开创性研究开启了综合利用时间、空间信息进行超分辨率研究的大门。图像的超分辨重建在数学上是一个病态逆问题，而且具有高度的计算复杂性，一直是图像处理领域中吸引人们不断深入研究的一项颇具挑战性的课题，若将其合理充分的运用，一定可以被证明是在探测器和光学系统限制下提高成像系统分辨力的有效方法。

就目前x射线成像而言，探测器部分本身光路复杂，对探测器本身进行亚像素偏移非常困难，因此很难获得存在亚像素移动的原始分辨率图像序列。结合目前x射线成像领域的实际发展，如何将图像的超分辨重建合理运用于x射线成像，提高更清晰准确的x成像，已成为当前急需解决和克服的一大技术问题。

技术实现要素：

本发明正是针对现有技术中探测器像素分辨率低的问题，提供了一种x射线超分辨成像方法，利用亚像素成像运动装置的亚像素级移动能力获取亚像素偏移投影成像，多次反复直至获得所需的亚像素偏移原始分辨图像序列为止，利用亚像素偏移分辨率提升算法对原始分辨图像进行重建，获得最终的投影域超分辨图像，同时，该方法也可以运用于ct图像域，不仅大大提高了像素分辨率，还能减少ct图像域环状伪影，进一步提高图像质量。

为了实现上述目的，本发明采用的技术方案是：一种x射线超分辨成像方法，包括以下步骤：

步骤1：待测对象、亚像素成像运动装置、高精度旋转台、及探测器分别固定好，保证位移过程不发生相对移动；

步骤2：进行模型的计算简化，忽略轻微数据损耗，准备多张投影图像；

步骤3：亚像素成像运动装置校准：获得实际操作中亚像素成像运动装置坐标系a3、ct图像域重建数据块中三维坐标系b3以及探测器平面的二维坐标系c2，a3、b3、c2上面的点a、b、c之间的关系是：

a＝tb

c＝kb

t是一个4*4的变换矩阵，将坐标系b3中的点转换到坐标系a3中；

k是一个4*3的投影矩阵，将一个坐标系b3中的点转换到坐标系c2中；

步骤4：确认需求的n*n张原始分辨率图像的偏移网格；

步骤5：确认坐标系c2上的所需齐次坐标，求得这一点在a3上面的坐标；

步骤6：超分辨重建：利用同一角度具有亚像素位移的多张原始分辨图像进行超分辨重建。

作为本发明的一种改进，所述步骤6中x射线原始分辨率图像的亚像素位移方式可以是移动被观测物体本身，也可以是探测器本身。

作为本发明的另一种改进，所述步骤4中偏移网格为：

(0,0)，(0,1/n)，(0,2/n)……(0，(n-1)/n)；

(1/n,0)，(1/n,1/n)，(1/n,2/n)……(1/n，(n-1)/n)；

……

((n-1)/n,0)，((n-1)/n,1/n)，((n-1)/n,2/n)……((n-1)/n，(n-1)/n)。

作为本发明的一种改进，所述步骤4中偏移网格在二维网格内n*n张图像亚像素偏移位置可任意排列。

作为本发明的又一种改进，所述步骤6中超分辨重建算法采用非均匀插值分辨率提升算法或者正则化分辨率提升算法。

为了实现上述目的，本发明采用的技术方案，一种x射线超分辨成像方法在ct图像域的应用，包括以下步骤：

步骤1：待测对象、亚像素成像运动装置、高精度旋转台、及探测器分别固定好，保证位移过程不发生相对移动；

步骤2：进行模型的计算简化，忽略轻微数据损耗，准备多张投影图像；

步骤3：亚像素成像运动装置校准：获得实际操作中亚像素成像运动装置坐标系a3、ct图像域重建数据块中三维坐标系b3，a3、b3、上面的点a、b之间的关系是：

a＝tb

t是一个4*4的变换矩阵，将坐标系b3中的点转换到坐标系a3中；

步骤4：确认需求的n*n*n张原始分辨率图像的偏移网格；

步骤5：确认坐标系b3上的所需齐次坐标，求得这一点在a3上面的坐标；

步骤6：超分辨重建：利用旋转一周的具有亚像素位移的原始分辨投影图像进行超分辨重建。

作为本发明的一种改进，所述步骤4中偏移网格为

(0,0,0)，……(0,(n-1)/n,0)；

……

((n-1)/n,0，0)，……((n-1)/n,(n-1)/n,0)；

.

.

.

(0,0,(n-1)/n)，……(0,(n-1)/n,(n-1)/n)；

……

((n-1)/n,0，(n-1)/n)，……((n-1)/n,(n-1)/n,(n-1)/n)。

作为本发明的一种改进，所述步骤4中偏移网格在三维网格内n*n*n张图像亚像素偏移位置可任意排列。

作为本发明的又一种改进，所述步骤6中超分辨重建采用投影域分辨率提升算法与ct重建的简单融合法或者ct分辨率提升迭代重建算法。

为了实现上述目的，本发明还采用了一种技术方案，一种x射线超分辨成像系统，包括x射线源、高精度旋转台、亚像素成像运动装置、探测器、整体机架结构及后处理计算机系统，

所述高精度旋转台，用于拍摄投影图，对待测对象进行旋转采图；

所述探测器，用来采集x射线源发出并穿过待测对象的x射线投影数据；

所述亚像素成像运动装置用来对待测对象或探测器进行亚像素级移动，产生亚像素偏移成像序列，

所述后处理计算机系统，用于进行亚像素偏移分辨率提升算法。

作为本发明的又一种改进，所述x射线源是锥束的x射线源。

与现有技术相比，本发明提出了一种x射线的成像方法及其在ct图像域的应用，解决了现有技术中探测器像素分辨率低的问题，基于探测器的亚像素移动硬件实施难，克服了现有技术领域的一大难题；同时由于多次进行拍摄，超分辨x射线投影图像泊松噪声少，信噪比比单张原始分辨率图像更高；在ct图像域的运用中，ct重建图像时，不仅能够提高图像域分辨率，而且还能减少ct图像域环状伪影问题，进一步提高图像质量。

附图说明

图1是x射线超分辨成像系统示意图；

图2是x射线投影到探测器上的投影示意图；

图3是x射线投影到探测器上具有(dx，dy)亚像素偏移的投影示意图；

图4是多幅原始分辨率图像的二维偏移网格示意图；

图5是多幅原始分辨率图像的三维偏移网格示意图；

图6是平移台装置的位移平面跟探测器平面完全平行的投影情况，平移台装置的位移平面不跟探测器平面完全平行的投影情况；

图7是亚像素成像运动装置位移方向、距离校准流程图；

图8是非均匀插值法的算法图解；

图9是独立的x射线超分辨成像方法实例图像。

图10是独立的x射线超分辨成像方法实例图像。

具体实施方式

以下将结合附图和实施例，对本发明进行较为详细的说明。

实施例1

一种x射线超分辨成像方法，我们将具体步骤分为系统搭建，实验操作过程及超分辨成像算法三个过程。

系统搭建：如图1，高分辨探测器对于精度要求更高，因此使用的亚像素成像运动装置是纳米级平移台。纳米平移台运动自由度为x、y、z，最大行程范围为100μm/轴，分辨率为7nm，重复定位精度为0.05％f.s.，线性度为0.15％f.s.，承载能力为6n，外形尺寸(lxwxh)为70×70×24mm，平移台上侧有负载安装孔(4-m3)，下侧有平台安装固定孔(4-m4)，这是为了防止三维xyz纳米平移台在移动过程中发生相对移动。射线源型号为microfocusx-raysourcel9421-02，焦点大小为5μm，探测器型号为andor_ikon-l_936，实际大小为27.6x27.6mm，像素数目为：2048x2048，像素大小(wxh)为13.5μmx13.5μm，像素有效面积比为100％。射线源、探测器、固定置物台位置固定，射线源到探测器的距离为8cm，固定载物台到探测器距离为0.3cm。

实验操作过程：拍摄几圈模体小球，获得系统校正参数，获得矩阵k，矩阵t，最后得到实验中需要多张原始分辨率图像的三维纳米平移台位置。投影域超分辨重建一共需要拍摄16次，电压为30v，电流为50ma，单次曝光时间为60s，拍摄结束之后进行投影域超分辨重建；ct图像域超分辨重建一共拍摄1000个角度，需要位移4次，每个位置的投影图像各为250张，这4组250间隔排列，此外电压为60v，电流为130ma，单次曝光时间为1500s，拍摄结束之后进行超分辨重建，即超分辨成像算法。

一种x射线超分辨成像方法，包括以下步骤：

步骤1：待测对象、亚像素成像运动装置、高精度旋转台、及探测器分别固定好，保证位移过程不发生相对移动，由于系统需要进行投影域超分辨重建，因此三维xyz纳米平移台需满足：保证射线源到探测器的距离足够大，且三维xyz纳米平移台位置尽量靠近探测器，如图1所示；

步骤2：进行模型的计算简化

1)在放大倍数小、三维平移距离小等这些前提下，不考虑锥束x射线投影过程中造成的投影图像轻微变化问题；

2)进行超分辨重建需要多张具有亚像素平移的原始分辨率投影图像；如图2,图3所示。

步骤3：亚像素成像运动装置校准：获得实际操作中亚像素成像运动装置坐标系a3、ct图像域重建数据块中三维坐标系b3以及探测器平面的二维坐标系c2，a3、b3、c2上面的点a、b、c之间的关系是：

a＝tb

c＝kb

t是一个4*4的变换矩阵，将坐标系b3中的点转换到坐标系a3中；

k是一个4*3的投影矩阵，将一个坐标系b3中的点转换到坐标系c2中；具体做法如下，具体操作流程图见图7：

1)建立重建数据的坐标系b3，如图6所示，是一个成像系统的几何模型，超分辨重建算法依赖于精确的系统空间几何关系，而系统几何参数就表征了这样的空间几何关系，用(x,y,z)轴建立坐标系b3。用z轴表示转台的旋转轴，即射线源和探测器绕z轴作圆周运动。x轴为射线源初始所在的轴，将x射线源的初始位置到z轴的垂足为坐标系原点，y轴则是经过原点并垂直于x轴和z轴的坐标轴。

2)建立探测器平面坐标系c2：在探测器平面中，用u表示探测器的横轴，并用v来表示探测器的纵轴。

3)建立b3坐标系跟c2之间的关系：整个空间几何关系可以用七个参数来表示，即(srd,sdd,θ，φ，η，uc，vc)，其中srd表示的是射线源到旋转轴的距离，sdd表示的是射线源到焦点的距离，θ表示的是探测器的前后偏转角，φ表示的是探测器左右偏转角，η表示的是探测器平面内旋转角，坐标点(uc，vc)表示的是焦点在探测器上的坐标，这7个参数能够精确地表示成像系统的空间关系。

4)利用旋转台拍摄一圈小球模体，用模拟退火算法计算获得3)中的参数值，用这些参数值计算获得4*3的投影矩阵k；

5)ct重建4)中的小球模体数据，获得小球模体在三维坐标系b3的质心坐标(x1，y1，z1)，说明在三维坐标系a3中的点(0，0，0)在三维坐标系b3对应的点坐标是(x1，y1，z1)；将三维纳米平移台移动至(1，0，0)后(在三维坐标系a3)，再次利用旋转台拍摄一圈小球模体并重建小球模体，获得小球质心在三维坐标系b3的坐标(x2，y2，z2)，说明在三维坐标系a3中的x轴向量(1，0，0)在三维坐标系b3对应的向量是(x2-x1，y2-y1，z2-z1)；重复上面的过程，获得三维坐标系a3中的y轴向量(0，1，0)在三维坐标系b3对应的向量是(x3-x1，y3-y1，z3-z1)，获得三维坐标系a3中的z轴向量(0，0，1)在三维坐标系b3对应的向量是(x4-x1，y4-y1，z4-z1)，那么：

可以得到4*4的变换矩阵t^-1为：

6)图像域超分辨需要坐标系a到坐标系b之间的转换关系，即t；投影域超分辨重建需要获得坐标系a到坐标系c之间的转换关系，则联立a＝tb，c＝kb可得：

c＝k*t^-1*a

步骤4：；投影域超分辨重建需要确认坐标系c2中需求的4*4张图像之间的偏移网格：如图4a:

(0,0)，(0,1/4)，(0,2/4)，(0,3/4)；

(1/4,0)，(1/4,1/4)，(1/4,2/4)，(1/4,3/4)；

(2/4,0)，(2/4,1/4)，(2/4,2/4)，(2/4,3/4)；

(3/4,0)，(3/4,1/4)，(3/4,2/4)，(3/4,3/4)。

步骤5：已经确认坐标系c2上所需偏移的齐次坐标，通过：

c＝k*t^-1*a，

可以求得这一偏移在坐标系a3上面的齐次坐标；

这种配准矩阵获取的方法的优点是，它是通过真实坐标直接获得的，准确性要比配准算法精度更高，规避了超分辨重建问题中配准算法不准确的问题。此外纳米平移台重复定位精度高，因此可以一次测量位移矩阵多次使用。

步骤6：超分辨重建：投影域超分辨重建利用同一角度具有亚像素位移的几张原始分辨率图像进行分辨率提升算法重建,本实施了采用非均匀插值分辨率提升算法：

a.原始分辨率图像配准，即运动估计；

b.通过对配准后的原始分辨率图像的非均匀插值生成一张高分辨率图像；

c.对生成的高分辨率图像去模糊。

非均匀插值方法的计算量较小，需要假定原始分辨率图像的噪声和模糊特征都是相同的，这可适用于本系统的图像分辨率提升。重建原理如图8，左边为原始分辨率图像lr(三角形，黑圆点，灰圆点，绿圆点四张低分辨图像)在高分辨网格里面的非均匀采样点分布，可以用直接非均匀插值重建或者迭代重建的方法获得高分辨率hr图像的均匀采样点分布图，最后再对这张图像进行去模糊处理等后处理操作。实验图例如图9所示。

实施例2

本实施例与实施例1的不同之处在于：步骤4中的确定偏移网格，可采用在二维网格内n*n张图像亚像素偏移位置任意排列，这种偏移网格同样也会获得一定的分辨率提升效果，如图4b所示；

步骤6中的超分辨重建同样利用同一角度具有亚像素位移的几张原始分辨率图像进行分辨率提升算法重建，本实施了采用的是正则化分辨率提升算法，一般来说，因为原始分辨率图像数量的不足，模糊矩阵，噪声等不可精确测量的因素的加入会导致它成为一个病态问题，很难直接获得最优解，这时候最好的解决方法就是利用解的先验条件来使得病态问题良态化。通过约束最小化下面这个公式：

这里wk＝dbkmk，α指的是正则化系数，用它来权衡保证项和正则项之间的比例，这里的正则项c是指tv正则项，可以用来保边除噪声，||·||指二范数。

实施例3

本实施例与实施例1不同之处在于：一种x射线超分辨成像方法在ct图像域的应用，包括以下步骤：

步骤1：待测对象、亚像素成像运动装置、高精度旋转台、及探测器分别固定好，保证位移过程不发生相对移动，由于系统需要进行投影域超分辨重建，因此三维xyz纳米平移台需满足：保证射线源到探测器的距离足够大，且三维xyz纳米平移台位置尽量靠近探测器，如图1所示；

步骤2：进行模型的计算简化

1)在放大倍数小、三维平移距离小等这些前提下，不考虑锥束x射线投影过程中造成的投影图像轻微变化问题；

2)进行超分辨重建需要多张具有亚像素平移的原始分辨率投影图像；如图2,图3所示。

步骤3：亚像素成像运动装置校准：获得实际操作中亚像素成像运动装置坐标系a3、ct图像域重建数据块中三维坐标系b3，a3、b3、上面的点a、b之间的关系是：

a＝tb

t是一个4*4的变换矩阵，将坐标系b3中的点转换到坐标系a3中；

步骤4：确认需求的n*n*n张原始分辨率图像的偏移网格，如图5所示，这种偏移网格分辨率提升效果最好；

(0,0,0)，……(0,(n-1)/n,0)；

……

((n-1)/n,0，0)，……((n-1)/n,(n-1)/n,0)；

.

.

(0,0,(n-1)/n)，……(0,(n-1)/n,(n-1)/n)；

……

((n-1)/n,0，(n-1)/n)，……((n-1)/n,(n-1)/n,(n-1)/n)；

步骤5：已经确认坐标系b3上所需偏移的齐次坐标，通过：

a＝tb

可以求得这一偏移在坐标系a3上面的齐次坐标。

这种配准矩阵获取的方法的优点是，它是通过真实坐标直接获得的，准确性要比配准算法精度更高，规避了超分辨重建问题中配准算法不准确的问题；

步骤6：超分辨重建：ct图像域超分辨重建利用具有亚像素位移的一圈投影图像进行重建，做法是在ct重建过程中融入分辨率提升算法。

1)简单融合法，第一步、利用fbp或者迭代重建获得不同亚像素位置的重建图像；第二步、对这些不同亚像素位置的重建图像进行超分辨重建，算法同投影域超分辨重建；

2)ct的分辨率提升迭代算法：将分辨率提升重建融入到ct重建中，ct的简单迭代算法的成像模型为：

x＝am0kz

a是正投影矩阵，z是待重建hr图像，x是投影hr图像，m0k是三维位移变换矩阵，这里可通过直接法的前两步获得位移矩阵。考虑探测器成像模型：

yk＝dbkmkx+nkk＝1...p(p为低分辨图像张数)

变换矩阵mk(这里为单位矩阵，可忽略)，系统模糊矩阵bk，降采样矩阵d，加性噪声为nk，将两式联立得到：

yk＝dbkam0kz+nk

将变换矩阵m0k与投影矩阵a结合，得到带变换量的投影矩阵ak，上式可转换为：

yk＝dbkakz+nk

建立一个目标函数：

前面一项为保真项，后面一项为tv正则化项，α为正则化参数。以骨头ct图为例，实例结果见图10。

实施例4

本实施例与实施例3不同之处在于：步骤4中的确定偏移网格，在三维网格内n*n*n张图像亚像素偏移位置可任意排列，这种偏移网格也会获得一定的分辨率提升效果。

实施例5

一种x射线超分辨成像系统，包括x射线源、高精度旋转台、亚像素成像运动装置、探测器、整体机架结构及后处理计算机系统，

所述高精度旋转台，用于拍摄投影图，对待测对象进行旋转采图；

所述探测器，用来采集x射线源发出并穿过待测对象的x射线投影数据；

所述亚像素成像运动装置用来对待测对象或探测器进行亚像素级移动，产生亚像素偏移成像序列，

所述后处理计算机系统，用于进行亚像素偏移分辨率提升算法。

该系统采用亚像素成像运动装置，获取多张具有可控亚像素偏移的x射线原始分辨率投影图像，再利用亚像素偏移分辨率提升算法进行投影域或者ct图像域的超分辨重建，从而实现超分辨的成像，解决了现有技术中探测器像素分辨率低的问题，基于探测器的亚像素移动硬件实施难，克服了现有技术领域的一大难题；亚像素偏移分辨率提升算法的位移矩阵是人为可控的，提高了算法精度，同时由于多次进行拍摄，超分辨x射线投影图像泊松噪声少，信噪比比单张原始分辨率图像更高；在ct图像域的运用中，ct重建图像时，不仅能够提高图像域分辨率，而且还能减少ct图像域环状伪影问题，进一步提高图像质量。

以上显示和描述了本发明的基本原理、主要特征和本发明的优点。本行业的技术人员应该了解，本发明不受上述实例的限制，上述实例和说明书中描述的只是说明本发明的原理，在不脱离本发明精神和范围的前提下本发明还会有各种变化和改进，这些变化和改进都落入要求保护的本发明范围内。本发明要求保护范围由所附的权利要求书及其等同物界定。