一种大视场CT成像重建方法与流程

本发明涉及辐射成像技术领域，具体涉及一种大视场ct成像重建方法。

背景技术：

锥束ct检测的视场取决于探测器的大小。在实践中，有很多需要ct检测的物体直径大，超出标准ct扫描能够提供的最大视场，这些物体以航空领域为多。中国专利文献(专利申请号分别为：cn200610012217.x，cn200710175588.4，cn200910091282.x，cn104597061a)分别对锥束大视场重建进行了公开。cn200610012217.x公开了一种能达到三倍标准扫描视场的ct成像方法，该方法采用重排技术将锥束数据处理为平行束数据，降低了成像的分辨率，占用内存多，计算速度慢，更为重要的是锥束数据并不能重排成若干共面的倾斜平行束数据。专利cn200710175588.4中的方法，由于采用微分运算，增大了图像的噪声，而且不能在纵向上得到更大的视场。cn200910091282.x公开了一种探测器偏置的大视野锥束x射线倾斜扫描三维数字成像方法，最大能将视场扩大为标准扫描视场的两倍。cn104597061a公开了一种基于虚拟探测器的大视场ct成像方法，该方法是目前唯一在实践中可行的可以在横向和纵向上都可以扩大视场的方法，但它并未能在非螺旋扫描下，在纵向上达到最大视场，而且对探测器数据的处理不够实用，重建未充分考虑探测器和转台的位置对成像质量的影响。

技术实现要素：

针对上述存在的问题，提供了一种探测器的布置方式灵活，对探测数据的处理适应实践的需要，重建方式具有广泛适应性的大视场ct成像重建方法。

本发明采用的技术方案如下：一种大视场ct成像重建方法，包括以下步骤：

s1、通过将物理探测器设置在不同位置，物体放在转台上，物理探测器在每一个位置采集得到相应位于成像平面上的数据，将得到的数据拼接融合为成像平面上的一块矩形区域drrect内的数据；物理探测器在每个位置采集数据时，转台做m周旋转，每个位置上对应的转台的旋转模式相同；螺旋模式下，m大于1，非螺旋模式下，m等于1；

s2、若光源与转台连线到成像面的投影点到所述矩形区域左边界的距离等于到所述矩形区域右边界的距离，则直接对矩形区域drrect内的数据进行fdk重建；

若投影点到所述矩形区域左边界的距离大于到所述矩形区域右边界的距离，将矩形区域向右拓展，再对拓展后的矩形区域drrect内的数据进行fdk重建；

若投影点到所述矩形区域左边界的距离小于到所述矩形区域右边界的距离，将矩形区域向左拓展，再对拓展后的矩形区域drrect内的数据进行fdk重建。

进一步的，步骤1中，探测器在数据采集前进行校正处理，具体包括：探测器在整个采集数据过程中仅进行一次暗场校正和坏像素校正，探测器在每个位置单独都进行亮场校正。

进一步的，所述探测器位置设置于成像平面上得到成像平面上的数据或设置于其他位置通过反透射变换得到相应位于成像平面上的数据。

进一步的，探测器采集的数据在拼接融合前进行矩形化处理，所述矩形化处理具体包括：对每次探测器采集变换到成像面上的数据取矩形区域drrectaux，得到探测器采集的矩形区域drrectaux数据，所述矩形区域面积小于矩形区域drrect，且四条边与矩形区域drrect的四条边平行。

进一步的，对于探测器采集的数据矩形化处理前还需对物理探测器上的数据进行切边处理，所述切边处理为将探测器边缘的像素切除，切掉的像素列数或行数大于等于0。

进一步的，所述步骤1中，在探测器采集的数据进行拼接融合后还需进行硬化校正。

进一步的，所述拼接融合的具体方法为：相邻的矩形数据存在重叠区域和非重叠区域，非重叠区域的数据直接采用原数据，重叠区域数据进行拼接融合处理；两两相邻的重叠部分的数据值为左右相邻重叠区域数据值的加权和，两加权系数和为1，由左至右，左边探测器数据的加权系数从1连续单调变为0，右边探测器数据的加权系数从0连续单调变为1；上下相邻重叠区域的拼接融合和左右相邻处理方式同理。

进一步的，若光源剂量有波动或在采集数据的过程中，探测器的探测性能有变化，需要对相邻的探测器数据进行亮度均衡处理。

进一步的，所述矩形区域drrect向右拓展的具体方法为：定义矩形区域横边长为axisrectlong，竖边长为axisrectshort；拓展后的矩形区域横边长为axisrectlongnew，竖边长不变；光源和转台连线在成像面的投影点到矩形区域drrect左边界距离与到矩形区域drrect右边界距离的差值向下取整为temp；

对于拓展后的矩形区域数据，0-temp区域内的数据直接拷贝矩形区域drrect0-temp区域内的数据，temp-axisrectlong区域内的数据为矩形区域drrect内temp-axisrectlong区域内的数据加权处理之后的数据，axisrectlong-axisrectlongnew区域内的数据值为0；

所述矩形区域drrect向左拓展的方式与向右拓展方式同理。

进一步的，所述矩形区域drrect确定方法为：光源与转台的连线到成像面的投影点在所述矩形区域drrect的内部，投影点到区域drrect边界的距离的最大值对应的边界应能满足，物体成的像和该边界不相交。

与现有技术相比，采用上述技术方案的有益效果为：不仅可以在横向纵向上扩大视场，还可以得到高质量的ct重建图像，重建方法更灵活，将转台偏移重建，探测器偏移重建，标准重建统一到一个框架下处理，充分考虑了成像过程中，光源不稳定，探测器探测效率变化，探测器坏像素，硬化现象，几何参数等因素对成像质量的影响，相比基于重排和bpf的方法，具有重建质量高，重建速度快，占用内存少的优点。

附图说明

图1为本发明的ct扫描流程图。

图2为本发明的ct扫描的结构示意图。

图3为本发明的ct扫描结构的俯视图。

具体实施方式

下面结合附图对本发明做进一步描述。

如图1所示，一种大视场ct成像重建方法，包括以下步骤：

s1、通过将物理探测器设置在不同位置，采集得到相应位于成像平面上的数据，将得到的数据拼接融合为成像平面上的一块矩形区域drrect内的数据；物体放在转台上，物理探测器在每一个位置，转台做m周旋转，每个位置上对应的转台的旋转模式相同；螺旋模式下，m大于1，非螺旋模式下，m等于1；

s2、若光源与转台连线到成像面的投影点到所述矩形区域左边界的距离等于到所述矩形区域右边界的距离，则直接对矩形区域drrect内的数据进行fdk重建；

若投影点到所述矩形区域左边界的距离大于到所述矩形区域右边界的距离，将矩形区域向右拓展，再对拓展后的矩形区域drrect内的数据进行fdk重建；

若投影点到所述矩形区域左边界的距离小于到所述矩形区域右边界的距离，将矩形区域向左拓展，再对拓展后的矩形区域drrect内的数据进行fdk重建。

具体方法如下：

获取矩形区域drrect数据

定义旋转面，法线为转台转轴朝上的方向，光源src在旋转面上；定义成像面，垂直于旋转面；在旋转面上，光源与转台中心rotcenter的连线位于成像面的同一侧；光源与转台中心的连线的延长线srcrotline与成像面相交于o。转台旋转要满足，重建区域的每一个点和光源的连线相交于当前角度下的拓展后的drrect上的点。转台旋转面可以在整个数据采集过程中高度保持不变，也可以上下运动。

矩形区域drrect的四条边，分别平行于旋转轴和旋转面，平行于旋转面的边定义为矩形的横边axisrectlong，中点为dhcenter，另一个边定义为矩形的竖边axisrectshort。

探测器是实际上采集数据的探测器，探测器应尽量放置在成像平面上，否则，物理探测器上数据应该进行反透射变换，得到相应的位于成像平面上的数据，为方便处理，这样得到的成像平面上的数据应取合适的矩形区域drrectaux，这些矩形区域的四条边平行于drrect上的四条边，面积小于矩形区域drrect。以下描述中的物理探测器数据指的是drrectaux上的数据。通过这些矩形数据drdataaux的拼接融合得到矩形区域drrect的数据。

所述矩形区域drrect确定方法为：光源与转台的连线到成像面的投影点在所述矩形区域drrect的内部，投影点到区域drrect边界的距离的最小值越大，成像质量越好，可检物体的尺寸越小，应用中需要选择合适的最小值；距离的最大值对应的边界应能满足，物体成的像和该边界不相交。

作为优选，所述探测器为平板探测器，非弧形探测器。

在本实施例中，探测器在采集数据前，需要进行暗场校正，亮场校正，坏像素校正。其中，暗场校正，坏像素校正和探测器的所处位置无关，亮场校正和探测器所处的位置相关；因此，探测器在整个采集数据过程中仅进行一次暗场校正和坏像素校正，而在探测器在每个位置单独进行亮场校正。

对于相邻探测器的数据还需要进行切边处理，即将探测器边缘的像素切除，切掉的像素列数或行数可以为0。

转台在相同的角度下，探测器在不同的n个位置采集n组drdataaux。转台需要较高的旋转定位精度，同一角度下相邻的drdataaux，需要进行拼接融合处理。相邻的drdataaux有重叠区域和非重叠区域，非重叠区域的数据直接采用原数据，重叠区域的数据进行数据融合。

所述融合具体方法为：

两两相邻的重叠部分的数据值为左右相邻重叠区域数据值的加权和，两加权系数和为1，由左至右，左边探测器数据的加权系数从1连续单调变为0，右边探测器数据的加权系数从0连续单调变为1，加权系数的二阶导数不可过大，避免加权系数呈现过大的陡坡结构；上下相邻重叠区域的拼接融合和左右相邻处理方式同理，仅方向改变。

在融合时同时进行探测器探元的像素值计算，用于重建过程中的计算，具体计算过程：

区域的点在相邻的重叠矩形区域drrectaux的平行的两条竖边1(属于右边探测器)，2(属于左边探测器)中间，对应两个探测器上的探元detunit1和detunit2，距离两条边1,2的距离分别为dis1,dis2,dis1+dis2＝l。该点的值pixelvalue＝pixelvalue1*dis2/l+pixelvalue2*dis1/l。pixelvalue1为边2对应的左边探测器上的点detunit1的像素值，pixelvalue2为边1对应的右边探测器上的点detunit2的像素值；上下拼接融合做类似的处理。

在一个优选实施例中，若光源剂量有波动，或在采集数据的过程中，探测器的探测性能有变化，需要对相邻探测器数据进行亮度均衡处理。具体均衡处理如下：

左右相邻的探测器数据，可以逐行进行亮度均衡处理，scale值为一行中重叠部分左右探测器数据的灰度均值的比值，处理方式为b中探测器各行数据乘上各行的scale。

上下相邻的探测器数据，不可以逐列进行亮度均衡处理，需要对探测器数据进行整体的处理，scale值为上下两个探测器的灰度均值的比值，下面探测器像素值都乘上scale。

ct重建

接下来根据所有角度下的drdata数据进行ct重建。其中，所有角度指转台旋转的角度，每个角度执行一次上述获取矩形区域drrect数据的步骤即可得到所有角度下的矩形区域drrect数据。

拼接融合后具体加权拓展drrect的过程如下：

如图2所示，定义光源在成像面的投影为o点，o点与矩形区域drrect左边界垂直距离为disorectleft，o点与矩形区域drrect右边界垂直距离为disorectright；光源到drrect的垂直距离为sdd，sod为光源与转台沿着光源垂直于drrect的直线方向上的距离。下面计算旋转面上drrect一行数据的拓展过程，其它行数据采用同样的方式进行拓展，加权系数，加权范围同旋转面上一行数据的处理方式相同。

以探测器探元大小作为长度单位。ct重建的参数有三个，转台的偏移量rothoffset和探测器的偏移量detectorhoffset，detectorvoffset，这三个值都是相对标准fdk重建的器件布局而言的，向右向上为正方向。detectorvoffset的处理和常规fdk重建的处理方式一样。

其中，min(disorectleft，disorectright)>0，该值越大成像质量越好。

具体拓展过程如下：

矩形区域drrect在横向上拓展，拓展后的矩形区域的横边长为axisrectlongnew。

axisrectlongnew＝axisrectlong+abs(detectorhoffset-rothoffset)*2

一半矩形区域的横边长为hfaxisrectlong＝axisrectlong/2

disorectleft＝hfaxisrectlong-(detectorhoffset-rothoffset*(sdd/sod))

disorectright＝hfaxisrectlong+(detectorhoffset-rothoffset*(sdd/sod))

若disorectleft＝disorectright，直接对矩形区域drrect内的数据进行重建；

若disorectleft>disorectright，将探测器向右拓展，再对矩形区域drrect内的数据进行重建；

若disorectleft<disorectright，将探测器向左拓展，再对矩形区域drrect内的数据进行重建。

以disorectleft>disorectright情况为例，说明拓展处理过程，即矩形区域drrect数据向右拓展，向左拓展的方式与向右拓展方式同理。光源与drrect左右两边连线与光源到成像面垂线的夹角为：

maxangle＝atan(((hfaxisrectlong+detectorhoffset))/sdd)

minangle＝atan(((-hfaxisrectlong+detectorhoffset))/sdd)

temp＝floor(disorectleft-disorectright)

theta＝maxangle-atan(rothoffset/sod)

对于拓展后区域矩形区域，0到temp数据直接拷贝原来axisrectlong长的探测器上的数据；从temp到axisrectlong的数据需要加权处理；从axisrectlong到axisrectlongnew数据值为0。

其中，加权处理的加权系数为weightfactor；

alpha＝atan(rothoffset/sod)-atan((xindex-hfaxisrectlong+detectorhoffset)/sdd)；

weightfactor＝0.5*(sin(pi*alpha/(2.0*theta))+1)；

其中xindex指数据当前横坐标，即以矩形区域的左端为坐标原点的坐标值，xindex范围为从temp到axisrectlong。

在另一个实施例中，weightfactor可以采用其它的形式，需要保证以o为对称点，两侧对称的点的加权系数和为1，weightfactor的值由1逐步过渡到0，加权系数的二阶导数不可过大，避免加权系数呈现明显的陡坡结构。

若在重建过程中，temp到axisrectlong的范围过于小时，重建断层数据的中心会出现伪影，需要先利用扫描的对称性计算axisrectlong以外的指定范围的数据，并在原始数据和计算出的数据间进行平滑过渡处理。

本发明根据实际探测情况还提供了探测器采集数据的具体设置方式，具体过程如下：

锥束成像的情况，假定旋转平面置于中心位置。在保证成像质量的情况下，非螺旋锥束ct成像，被检物体的最大长度为betamax。物体垂直于旋转平面的方向定义为物体的长度方向。光源与探测器的距离为sdd，光源与物体的距离为sod。

首先分析，纵向上成像的情况。

物理探测器的尺寸，宽为dw，高为dh。物体长度为objlen，物体在探测器上成的像长度为objlenimaging＝objlen*sdd/sod。ct成像有三种情况：

1)：objlenimaging小于dh

探测器横向上重叠采集数据，重叠的程度取决于相邻探测器所采集到的数据的灰度均值的差异。差异大，要求重叠区域大，差异小，要求重叠区域小。转台旋转面固定不变

2)：objlenimaging大于dh，小于betamax

探测器需要像1)中描述的方式采集数据，还需要在纵向上进行采集、拼接融合处理。直到纵向上长度能满足重建的需要。

3)：objlenimaging大于betamax

探测器横向上采用1)的方式采集，处理数据。同时转台相对物体垂直于旋转平面移动。旋转的角度整数倍于2pi。探测器也可以横向上采用2)的方式采集，但是这样的锥角大，容易造成成像质量的下降。

横向上成像的情况。rothoffset接近0，detectorhoffset足够大，成像变为经典的探测器偏置重建成像模式。detectorhoffset接近0，rothoffset足够大，成像变为经典的转台偏置重建模式。

本发明并不局限于前述的具体实施方式。本发明扩展到任何在本说明书中披露的新特征或任何新的组合，以及披露的任一新的方法或过程的步骤或任何新的组合。如果本领域技术人员，在不脱离本发明的精神所做的非实质性改变或改进，都应该属于本发明权利要求保护的范围。

本说明书中公开的所有特征，或公开的所有方法或过程中的步骤，除了互相排斥的特征和/或步骤以外，均可以以任何方式组合。

本说明书中公开的任一特征，除非特别叙述，均可被其他等效或具有类似目的的替代特征加以替换。即，除非特别叙述，每个特征只是一系列等效或类似特征中的一个例子而已。