一种提高锥形束SPECT成像分辨率的方法与流程

本发明涉及一种spect成像方式，特别涉及一种在不改变锥形束spect检测器准直器孔径大小的前提下提高成像分辨率的方法，属于医学图像重建领域。

背景技术：

单光子发射计算机断层成像(singlephotonemissioncomputedtomography,spect)是核医学的一项重要技术，目前在临床上有广泛的应用。在单光子发射断层成像中，放射性示踪剂被注入到病人体内，根据示踪剂在人体内分布，重建其在人体内的密度分布图像，该图像能够直观反映人体器官新陈代谢水平的高低，从而实现功能成像的目的，对于疾病的早期诊断具有重要的意义。

在spect中，放射性示踪剂通过辐射γ光子，穿过人体后被检测器检测，最终重建得到断层图像。在spect中，检测器的准直器孔径是限制图像空间分辨率和灵敏度的关键因素。准直器孔径越小，图像的空间分辨率越高，但接收到的γ光子数量越少，灵敏度越低，成像质量会变差。准直器孔径越大，接收到的γ光子数量越多，灵敏度越高，但图像的空间分辨率越低，成像质量也会变差。

超分辨率重建是从多幅低分辨率图像重建出高分辨率图像。作为一种有效的提高图像分辨率的方法，它已经被广泛应用于汽车牌照、遥感卫星，医学成像等领域。近年来超分辨率重建被应用到磁共振、ct、pet、spect成像领域，以提高图像分辨率，但这些方法主要应用在图像域或平行及扇形检测器投影重建领域。

锥形投影spect系统与传统的平行束及扇形束spect相比，具有高分辨率的特点，有利于对小目标进行检测及重建，锥形束投影重建和超分辨率重建相结合，将会极大提高spect的成像分辨率。

技术实现要素：

为了解决单光子发射计算机断层成像(singlephotonemissioncomputedtomography,spect)分辨率低的问题，本发明公开的一种提高锥形束spect成像分辨率的方法要解决的技术问题是提供一种提高成像分辨率的方法，在不改变检测器中准直器孔径大小的前提下，提高spect图像的分辨率。

本发明是通过下述技术方案实现的。

本发明公开的一种提高锥形束spect成像分辨率的方法，包括如下步骤：在不改变检测器中准直器孔径大小的前提下进行低分辨率投影图采集；由于多幅低分辨率图像都对每个角度同一目标进行成像，但相互之间又包含不同的投影数据信息，通过对采集到的多幅低分辨率投影图像的数据信息进行处理，得到每个角度下的高分辨率投影数据；由此能够从多幅低分辨率投影图像中重建出高分辨率spect图像，即提高重建spect图像的分辨率。

所述的重建出高分辨率spect图像方法可选用解析法或迭代法。

本发明公开的一种提高锥形束spect成像分辨率的方法，包括如下步骤：

步骤一、低分辨率投影图采集。在不改变检测器中准直器孔径大小的前提下，检测器绕物体进行旋转，每次旋转角度，n为检测器绕物体旋转的次数。在每个角度下，分别在其所在位置平面测得一组投影数据，每组投影数据包含m个投影数据。单个旋转角度测量每组m个投影数据结束后，检测器绕物体继续旋转角度，进行下一个角度测量，直至2π范围内的n个角度都测量完毕。

步骤一中单个旋转角度测量每组m个投影数据的具体采集方法通过方法1或方法2实现。

方法1：检测器在其所在的平面内绕检测器中心旋转，共需要旋转m次，每次旋转后采集单个投影数据，并记录每次旋转的角度。

所述的检测器在其所在的平面内绕检测器中心旋转的旋转角度为固定角度均匀旋转或不固定角度非均匀旋转，优选在设定的角度θ范围内均匀旋转m个角度，即每次旋转角度，为避免像素点之间重合，设定的角度范围θ进一步优选为

方法2：检测器中心每次需平移预设的距离，共需要平移m次，每次平移后采集单个投影数据，并记录每次平移的距离。

每次平移的距离优选检测器中心每次平移高分辨率图像一个像素大小的距离。

步骤二、对步骤一中检测器绕物体进行旋转，每次旋转角度，每个角度下采集到的m幅低分辨率图像进行处理，由于多幅低分辨率图像都对同一目标进行成像，但相互之间又包含不同的投影数据信息，由此得到每个角度下的高分辨率投影数据。

步骤三、根据步骤二得到的高分辨率投影数据重建出高分辨率spect图像。

由于步骤一中多幅低分辨率投影图像都对同一目标进行成像，但相互之间又包含不同的投影数据信息，由此能够得到高分辨率投影图像，相应的重建图像的分辨率也会相应的提高。由此能够从多幅低分辨率投影图像中重建出高分辨率spect图像，因此能够在不改变检测器中准直器孔径大小的前提下，提高spect图像的分辨率。

步骤三中重建高分辨率spect图像方法选用解析法或迭代法。

步骤一中所述的检测器绕物体旋转的次数n，以及在每个角度下分别在其所在位置平面测得的投影数据个数m与实际低分辨率图像大小和需获得的高分辨率图像大小有关。若低分辨率检测器分辨率为ls×ls个像素，目标获得的高分辨率像素为hs×hs个像素，则m、n、hs、ls可由公式(1)、(2)、(3)给出。

n≥hs(1)

m≥de×de(3)

步骤二中所述的对每个角度下采集到的m幅低分辨率图像进行处理，根据步骤一中检测器绕物体进行旋转，单个旋转角度测量每组m个投影数据的两种采集方法，有针对两种采集方法的不同实现步骤，即分别针对方法1、方法2采集到的m幅低分辨率图像进行处理的方法。

针对方法1中对检测器绕物体旋转后采集到的m幅低分辨率图像进行处理的方法包括如下步骤：

记为检测器绕自身旋转的第i个角度下检测器测量得到的低分辨率投影数据，phi为待求的检测器绕自身旋转的第i个角度下的高分辨率投影数据。

步骤2.1：以检测器所在的平面建立坐标系，坐标系的中心为检测器平面中心。对phi按建立的坐标系进行旋转操作，共需要旋转m次。设旋转第i个角度后的高分辨率投影为ph(i+1)，则：

其中(x,y)为phi中像素点的坐标，(x',y')为ph(i+1)中像素点的坐标，γ为ph(i+1)相对于phi旋转的角度。根据公式(4)算得的ph(i+1)中像素点的坐标(x',y')，由双线性插值公式求得ph(i+1)中每个点的像素值。

步骤2.2：将每个低分辨率区域对应的高分辨率区域块的所有像素点按照公式(5)叠加，即对ph(i+1)进行降采样得到pl(i+1)。

步骤2.3：将经步骤2.2降采样得到的低分辨率图像与直接采集到的低分辨率图像按照公式(6)进行比较，求得真实低分辨率投影与目前低分辨率率投影的差值δpl(i+1)。

步骤2.4：将经步骤2.3得到的对应像素点的差值乘以每次调节的步长作为调节权重，得到对高分辨率图像的调节权重如公式(7)所示。

x'＝x×de+i(i＝0,1,2…(de-1))

y'＝y×de+i(i＝0,1,2…(de-1))

求得ph(i+1)中每个点的调节权重δph(i+1)，其中step为每次的调节步长，调节步长一般在0到1之间取值。

步骤2.5：由步骤2.4得到的调节权重δph(i+1)按照公式(8)更新ph(i+1)。

ph(i+1)＝ph(i+1)+δph(i+1)(8)

步骤2.6；回到步骤2.1进行下一次旋转角度调节或进行下一次循环，直到达到要求的迭代结束条件，得到满足投影分辨率要求的高分辨率投影图像，结束迭代过程，即得到检测器绕物体旋转时单个角度下的高分辨率投影数据。

所述的迭代结束条件为：计算前后两次迭代过程高分辨率图像的差值的二范数，并与结束迭代阈值条件进行比较，从而确定是否结束迭代。所述的迭代阈值条件根据需要获得的高分辨率投影精度而定。

针对方法2中对检测器绕物体旋转后采集到的m幅低分辨率图像进行处理的方法包括如下步骤：

记为检测器测量得到的低分辨率投影数据，phi为待求的检测器中心第i次平移到预设的距离后高分辨率投影数据。

步骤2.7：以检测器所在的平面建立坐标系，坐标系的中心为检测器平面中心。对phi按建立的坐标系进行平移操作，共需要平移m次。设平移第i次后的高分辨率投影为ph(i+1)，则：

(x',y')＝(x,y)+(δx,δy)(9)

其中(δx,δy)为ph(i+1)中心相对于phi中心的位移。(x',y')和(x,y)分别为ph(i+1)和phi上的点。由公式(10)求得ph(i+1)中每个点的像素值。

ph(i+1)(x',y')＝phi(x0,y0)(10)

步骤2.8：按步骤2.2-步骤2.5对高分辨率图像进行调节。

步骤2.9：回到步骤2.7进行下一次平移调节或进行下一次循环，直到达到要求的迭代结束条件，得到满足投影分辨率要求的高分辨率投影图像，结束迭代过程，即得到检测器绕物体旋转时单个角度下的高分辨率投影数据。

所述的迭代结束条件同步骤2.6里面的迭代结束条件，即所述的迭代结束条件为：计算前后两次迭代过程高分辨率图像的差值的二范数，并与结束迭代阈值条件进行比较，从而确定是否结束迭代。所述的迭代阈值条件根据需要获得的高分辨率投影精度而定。

有益效果：

本发明公开的一种提高锥形束spect成像分辨率的方法，由于多幅低分辨率投影图像都对同一目标进行成像，但相互之间又包含不同的信息，由此能够根据步骤二的处理方法从多幅低分辨率投影图像来获取一幅高分辨率投影图像，相应的重建图像的分辨率也会相应的提高。因此能够在不改变检测器中准直器孔径大小的前提下，提高spect图像的分辨率。

附图说明

图1为锥形spect检测器模型示意图。

图2为锥形检测器绕物体旋转采集示意图。

图3为检测器所在平面坐标系示意图，坐标系的中心为检测器中心。

图4为方法一中检测器平面绕其自身平面中心旋转示意图。

图5为方法二中检测器中心每次移动位置示意图。

图6为理想高分辨率投影图像与实际低分辨率投影图像对应关系示意图。

图7为检测器绕物体旋转时单个角度下由多幅低分辨率投影图像获得高分辨率投影图像流程图。

图8为由低分辨率投影图像获得高分辨率投影图像时单次迭代调整过程示意图。

图9为二维线性插值坐标示意图。

图10为不同层原始图像、真实模拟高分辨率投影图fdk重建图像、由低分辨率投影旋转移动获得高分辨率投影经fdk重建图像、由低分辨率投影平行移动获得高分辨率投影经fdk重建图像以及低分辨率fdk重建图像。其中检测器中心位于物体44层，检测器焦点到检测器中心的距离为380个像素，焦点到物体中心的距离为310个像素。

具体实施方式

下面结合附图和具体实施例对本发明方法的实施方式做详细说明。

实施例1：

检测器绕物体旋转的同时，其自身也需在其所在平面进行移动。本实例以低分辨率检测器分辨率为32×32个像素，目标获得的高分辨率像素为128×128，n＝128，m＝16，de＝4，为例进行说明，但是并不说明该发明仅限于这些条件，本发明中所涉及的范围仅由权利说明书限定。

本实施例公开的一种针对锥形束spect超分辨率成像方法，其具体实施步骤如下：

步骤一、低分辨率投影图采集。在不改变检测器中准直器孔径大小的前提下，检测器绕物体进行旋转，每次旋转角度。图1为锥形束spect检测器模型示意图，锥形束spect检测器为一平板，且其所有准直器孔对应一个相同的焦点。被检测物体位于检测器及焦点中间，从物体发出的伽马射线可以被检测器接收。正常情况下，spect在进行采集时，检测器平板需绕物体均匀旋转一周。在旋转过程中，检测器会采集到不同角度下从物体内部发出的伽马射线，根据各个角度下采集到的伽马射线，我们就可以重建出物体内部的射线源强度分布。图2为检测器绕物体旋转进行光子采集示意图，实线为检测器初始位置，虚线为检测器绕物体旋转一定角度到达的位置，检测器焦点会绕着图中的虚线圆旋转一周。

在每个角度下，分别在检测器所在位置平面测得一组投影数据，每组投影数据包含16个投影数据，单个旋转角度测量每组16个投影数据结束后，检测器绕物体继续旋转角度，进行下一个角度测量，直至2π范围内的128个角度都测量完毕。单个旋转角度测量每组16个投影数据具体采集方法通过方法1或方法2实现。

方法1：检测器在其所在的平面内绕检测器中心旋转，在范围内均匀旋转16个角度，即每次旋转角度，共需要旋转16次，每次旋转后采集单个投影数据，并记录每次旋转的角度。如图3所示为检测器平面及在其上面建立的坐标系示意图，坐标系原点为检测器中心位置。图4为检测器绕其自身中心在角度范围内均匀旋转示意图。图中实线为检测器初始位置，虚线为检测器每次旋转一定角度后到达的位置。

方法2：检测器中心每次平移高分辨率图像一个像素大小的距离，共需要平移16次，每次平移后采集单个投影数据，并记录每次平移的距离。检测器中心移动位置如图5所示。图5中的小黑点为检测器每次移动后的中心位置。

步骤二、对步骤一中检测器绕物体进行旋转，每次旋转角度，每个角度下采集到的16幅低分辨率投影图像进行处理，得到每个角度下的高分辨率投影数据。

步骤二中所述的每个角度下采集到的16幅低分辨率投影图像进行处理，根据步骤一中检测器绕物体进行旋转，单个旋转角度测量每组16个投影数据的两种采集方法，有针对两种采集方法的不同实现步骤，即分别针对方法1、方法2采集到的16幅低分辨率图像进行处理的方法。

针对方法1中对检测器绕物体旋转后采集到的16幅低分辨率投影图像进行处理的方法包括如下步骤：

记为检测器绕自身旋转的第i个角度下检测器测量得到的低分辨率投影数据，phi为待求的检测器绕自身旋转的第i个角度下高分辨率投影数据。

图6为理想高分辨率投影图像与实际低分辨率投影图像之间关系示意图。在实际采集过程中，通过检测器平面绕其自身旋转，可以获得图中右边每个角度下的低分辨率图像。图中标粗黑色方框代表高分辨率像素点与低分辨率像素点之间的对应关系。

图7为由16幅32×32像素的低分辨率投影图像获得1幅128×128像素高分辨率投影图像流程图。

流程图中初始高分辨率投影图像获得说明：

初始高分辨率图像获得选用下述两种方法之一：(1)，对高分辨率图像随机赋值，此种方法在最终迭代后会有大量的随机噪声。(2)，由采集到的低分辨率图像插值成高分辨率图像大小，并将此图像作为原始高分辨率图像，此种方法在最终迭代后获得高分辨率图像较好，且收敛速度较快。

单次迭代说明：根据低分辨率图像采集时每次旋转的角度大小，高分辨率图像绕其中心每次旋转相同的角度，然后进行降采样。此处降采样的方法为：将每个低分辨率区域对应的高分辨率区域块的所有像素叠加(如图8所示：将高分辨率图像每个加粗黑色方框内的所有像素值相加，作为对应低分辨率图像中的一个像素值)。将经降采样得到的低分辨率图像与采集到的低分辨率图像进行比较，将对应像素点的差值乘以每次调节的步长作为调节权重，对高分辨率图像进行调节。

对应过程如图8所示：(a)为高分辨率图像，(b)为由高分辨率图像经降采样得到的图像，(c)为实际采集到的低分辨率图像，(d)为由(b)和(c)相减的差值得到的图像，将得到的(d)作为调节依据，对原始高分辨率图像(a)进行调节，最终得到调节后的高分辨率图像(e)。然后获得的(e)进入下一次调节过程。

迭代结束条件：迭代结束条件的表示有两种方式，这两种方式实质上是一样的，(1)，通过计算图8中(a)与(e)的差值的二范数并与结束迭代阈值条件进行比较，从而确定是否结束迭代。(2)，计算图8中(d)矩阵的二范数并与结束迭代阈值条件进行比较，从而确定是否结束迭代。

步骤2.1：对phi按图(3)所示建立坐标系进行旋转操作，共需要旋转16次。设旋转第i个角度后的高分辨率投影为ph(i+1)，则根据公式(4)算得的ph(i+1)中像素点的坐标(x',y')，由线性插值公式(11)求得ph(i+1)中每个点的像素值。图9为二维线性插值坐标示意图。

步骤2.2：将每个低分辨率区域对应的高分辨率区域块的所有像素按照公式(5)叠加，即对ph(i+1)进行降采样得到pl(i+1)。

步骤2.3：将经步骤2.2降采样得到的低分辨率图像与直接采集到的低分辨率图像按照公式(6)进行比较，求得真实低分辨率投影与目前低分辨率率投影的差值δpl(i+1)。

步骤2.4：将经步骤2.3得到的对应像素点的差值乘以每次调节的步长作为调节权重，得到对高分辨率图像的调节权重如公式(7)所示。求得ph(i+1)中每个点的调节权重δph(i+1)，其中step为每次的调节步长，step＝1。

步骤2.5：由步骤2.4得到的调节权重δph(i+1)按照公式(8)更新ph(i+1)。

步骤2.6：回到步骤2.1进行下一次旋转角度调节或进行下一次循环，直到达到要求的迭代结束条件，得到满足投影分辨率要求的高分辨率投影图像，结束迭代过程，即得到检测器绕物体旋转时单个角度下的高分辨率投影数据。

所述的迭代结束条件为：计算前后两次迭代过程高分辨率图像的差值的二范数，并与结束迭代阈值条件进行比较，从而确定是否结束迭代。所述的迭代阈值条件根据需要获得的高分辨率投影精度而定。

针对方法2中对检测器绕物体旋转后采集到的16幅低分辨率投影图像进行处理的方法包括如下步骤：

记为检测器测量得到的低分辨率投影数据，phi为待求的检测器中心第i次平移到预设的距离后高分辨率投影数据。

步骤2.7：对phi按图(4)所示的坐标系进行平移操作，共需要平移16次。设平移第i次后的高分辨率投影为ph(i+1)，则根据公式(9)，(10)求得ph(i+1)中每个点的像素值。

步骤2.8：按步骤2.2-步骤2.5对高分辨率图像进行调节。

步骤2.9：回到步骤2.7进行下一次平移调节或进行下一次循环，直到达到要求的迭代结束条件，得到满足投影分辨率要求的高分辨率投影图像，结束迭代过程，即得到检测器绕物体旋转时单个角度下的高分辨率投影数据。

所述的迭代结束条件同步骤2.6里面的迭代结束条件。

步骤三、根据步骤二得到的高分辨率投影数据根据fdk重建算法重建出高分辨率spect图像。图10为利用fdk重建算法对投影图像重建所得结果。

以上结合具体实施例对本发明的技术方案和具体实施方式作了说明，但这些说明不能被理解为限制了本发明的范围，这些仅是举例说明，可以对这些实施方式做出多种变更或修改，而不背离本发明的原理和实质。本发明的保护范围由随附的权利要求书限定，任何在本发明权利要求基础上的改动都是本发明的保护范围。