CN102068270B - 一种使用静止环状射线阻挡阵列的锥束ct扫描和散射校正方法 - Google Patents

Abstract

本发明公开了一种使用静止环状射线阻挡阵列(Beam Stop Array，BSA)的锥束CT扫描和散射校正方法，使用静止环状BSA测量锥束投影中的散射分量，对因散射分量测量而被遮挡的原发射线分量，使用基于投影相关性的角度内插(Projection Correlation based View Interpolation，PC-VI)方法恢复，同时提出了使PC-VI性能有效发挥的静止环状BSA的设计方法。本发明方法只需单次扫描即可同时获得散射测量值和原发射线分量，形成一种具有可观实用前景的含有散射校正的锥束CT扫描方法，具有剂量低、扫描结果准确稳定、易于实现和操控以及成本低廉的优点。

Description

一种使用静止环状射线阻挡阵列的锥束CT扫描和散射校正方法

技术领域：

本发明属于X射线锥束CT中基于散射测量进行散射校正的扫描方式，具体涉及一种使用静止的环状射线阻挡阵列进行一次扫描同时获得散射分量和原发射线分量的方法。

背景技术：

散射校正是锥束CT中的难题。基于平面射线阻挡阵列(Beam Stop Array，BSA)进行散射测量的校正方法因其能够提供准确稳定的散射测量结果，在千伏级诊断CT和兆伏级放疗CT中获得了广泛研究。由于采用铅质BSA进行射线阻挡，探测器上对应于BSA阴影的信号全部来自散射光子，因此散射测量较为准确，利用散射光通量在面探测器上主要是光滑的低频成分的特点，对BSA测量值进行空间二维空域插值，即可得到当前投影角度的散射光通量。而与此同时，为获取因BSA遮挡而丢失的原发光通量信号，需要额外进行一次扫描，极大的加重了病人的剂量辐射。可见，基于BSA的散射测量方法的关键问题在于：在不增加或少增加剂量的前提下，恢复因BSA遮挡而丢失的原发光通量信号。学术界就这一问题的相关研究如下：2004年Ning R.等学者提出除正常扫描外，只需加装BSA进行一次稀疏角度扫描，利用散射光通量在角度上低频变化的特点，采用角度上的三次样条内插，即可估计出稀疏角度之外其他角度的散射光通量(Ning R，Tang X，and Conover D X-ray scattercorrection algorithm for cone beam CT imaging 2004Med.Phys.311195-202)；为了进一步降低剂量，学者Zhu L.等提出仅使用移动的BSA进行一次扫描的散射校正方法(Zhu L，Strobel N and Fahrig R X-Ray ScatterCorrection for Cone Beam CT Using Moving Blocker Array 2004SPIE 5745251-58)，对于因BSA阻挡而丢失的原发光通量信号，该方法仍然使用空域插值来进行恢复。然而，由于空域插值只适用于低频信号，所以对于存在边缘信息(例如，骨和软组织或空气的边缘)的原发光通量来说，空域插值表现并不理想。正因如此，该方法特别设计了光栅式移动(Raster-moving)的BSA移动模式，使得对原发信号的阻挡位置不断改变，从而降低空域插值误差在特定位置的积累。尽管如此，从该方法的重建图像中，尤其是某些切片图像中，仍能观察到比较明显的高频误差导致的噪声增加和条纹伪影。最近，基于移动BSA的方法在实验型研究中得到了一定的扩展(Jin J，Ren L，Liu Q et al“Combining scatter reduction and correction to improve image quail incone-beam computed tomography(CBCT)”，2010，Med.Phys.375634-5644；Wang J，Mao W，and Solber T“Scatter correction for cone-beam computedtomography using moving blocker strips：A preliminary study”，2010，Med.Phys.375792-5800)，但是从理论角度而言，这些方法仍未在恢复原发信号方面取得本质进展。

值得鼓舞的是，最近，国内学者Yan H.和Mou X.等发现和推导出了锥束投影的相关性(Pr ojection Correlation，PC)，即锥束投影在角度上存在偏微分一致性。该一致性说明，当前投影角度的高频信息可以从其相邻投影角度获取。在此基础上开发的基于锥束投影相关性的角度上的插值算法(ProjectionCorrelation based View Interpolation，PC-VI)，以空域插值结果作为初值，使用投影相关性进行信息的二次恢复，应用于上述使用移动BSA的散射校正方法后，取得了令人印象深刻的优异结果：丢失信号的大部分高频部分得到有效恢复，重建结果基本接近无遮挡的理想重建效果(Yah H，Mou X，Tang S andXu Q and Zankl M，“Projection correlation based view interpolation for conebeam CT：primary fluence restoration in scatter measurement with moving beamstop array”，2010，Phys.Med.Biol.556353-6375)。此外，在PC-VI框架下，移动BSA的移动方式可有更加选择，例如，在Yan H.和Mou X.等的论文中，除既有的光栅式移动模式，还提出了一维横向移动和旋转模式。

移动BSA方法至此获得了在散射测量，原发信号恢复以及低剂量三个方面的完美平衡。然而，即使在具有很强高频信息恢复能力的基于PC-VI的移动BSA方法中，BSA仍然需要移动，其原因也在Yan H和Mou X.等的论文中给出了详细分析：PC-VI算法有效工作的前提是，相邻角度投影被遮挡位置区域不互相重合，即所谓的PC-VI条件。考虑到实际应用中X线球管和探测器本身都在以较高的速度旋转，因此移动BSA将涉及到和球管探测器的运动匹配和校正等很多后续问题，这是目前这类方法的局限点。

发明内容：

本发明的目的在于克服上述现有技术的局限，提供了一种不需要移动的，同时能够满足上述PC-VI条件的环状BSA。在PC-VI算法框架下，该静态环状BSA表现出虚拟的移动BSA的效果，因为不需移动，相比移动BSA，其实现的结构简单，成本低廉，实际可操控性和可维护性更好，综合表现更稳定。总体来看，在基于PC-VI算法恢复被遮挡原发信号的框架下，使用静态环状BSA进行一次扫描同时完成测量散射和获取原发射线分量的方法具有很好的实用前景。

为达到上述目的，本发明采用的技术方案是：

一种使用静止环状BSA的锥束CT扫描和散射校正方法：使用静止环状BSA的扫描机械结构进行锥束CT扫描与散射校正，使用基于投影相关性的角度内插PC-VI方法恢复被静止环状BSA阻挡的原发射线分量。

使用静止环状BSA的扫描机械结构进行锥束CT扫描是指：在通常的X光源点轨迹和平板探测器轨迹间内，加入静止的环状BSA用来测量散射。环状面材料为对X光低衰减的物质，阵列单元材料为铅或其他原子序数较高的对X光强衰减的物质；静止环状BSA的结构以及阵列材料的排列符合给定的要求：使PC-VI方法性能有效发挥的静止环状BSA的设计方法。

使用静止环状BSA进行散射校正是指：

(a)投影测量数据I_A通常包含原发射线分量I_P和散射分量I_S，其中I_P为重建图像所需数据；在使用环状BSA测量I_S时，测量像素分两类，一类为未被BSA遮挡像素，其测量所得信号为I_A ^UB，包含原发射线分量I_P和散射分量I_S；另一类为有遮挡像素，其测量所得信号为I_S ^B，为散射分量I_S在该位置的采样值；

(b)对I_S ^B做空域插值得到任意像素I_S的估计值I_S ^～；对无遮挡像素，从I_A ^UB中减去I_S ^～即可得到原发射线分量计算值I_P ^UB～；对被遮挡像素，其原发射线I_P ^B～完全丢失，利用PC-VI方法恢复被阻挡原发射线部分，记为I_P ^B～。I_P ^{UB ～}与I_P ^B～的集合I_P ^～，即对CT重建所需原发射线投影I_P的计算结果。

使用PC-VI方法恢复被阻挡原发射线分量是指：

(a)首先，利用未遮挡探测器单元去除散射后的投影I_P ^UB～对遮挡单元作空域内插得到其丢失原发射线信号的初步恢复值I_P ^B～IR，其中上标^IR表示初步恢复，该空域插值可以采用现有的方法，如水平方向的三次样条插值。对所得整个投影(I_P ^UB～和I_P ^B～IR)作对数运算将其转换为线积分并加权，记为g^IR(λ+dλ)，其中λ+dλ表示当前投影角度，乘性加权因子为1/|ξ-η|，上标^IR表示被变换为g之前的投影中含有不够准确的所谓初步恢复部分I_P ^B～IR。根据空域内插的特点，g^UR(λ+dλ)中初步恢复的是丢失信号的低频部分。

投影相关性的核心公式为：

$g_{\mathrm{\λ}}^{*}\≅\frac{j}{r+d}\·(-2\frac{k_1}{k_2}{g}_{k_1}^{*}+k_2{g}_{k_1{k}_2}^{*}-k_1{g}_{k_2{k}_2}^{*}),k_2\≠0---\left(1\right)$

其中g^*为g(λ，u，v)的傅立叶变换g^*(λ，k₁，k₂)的简化表示形式；符号

和

分别表示g^*对变量λ和k₁求偏导；表示g对变量k₁，k₂求偏导；

表示g对变量k₂求两次偏导；涉及的锥束扫描参数罗列如下：X光源为ξ(r，λ，z)，检测器单元为

ξ包含的三个参数中，r代表ξ到旋转中心O的距离，λ为ξ的方位角，z为ξ的纵坐标。η包含的三个参数中，d代表旋转中心到η的距离，u、v分别代表η在面探测器上的横、纵坐标。

定义(1)等号右端部分为G(λ，g^*(λ))，(1)可写成：

$g^{*}(\mathrm{\λ}+\mathrm{d\λ})\≅g^{*}\left(\mathrm{\λ}\right)+\mathrm{d\λ}\·G(\mathrm{\λ},g^{*}\left(\mathrm{\λ}\right)),k_2\≠0---\left(2\right)$

(b)通过如(2)所示的方程，恢复丢失信号频谱中k₂≠0之外的其他频段，尤其是高频部分可以从相邻投影的频谱g^*(λ)计算得到(同理也可用到另外一侧的相邻投影g^*(λ+2dλ))。在k₂＝0处则直接用

的结果。将左右两边相邻角度的结果以权重w相结合(w取值范围为0到1，一般选取为0.5，可根据左右投影角度可提供的有效信息量做调整)，于是得到(λ+dλ)，上标^RR表示细化恢复：

$g^{\mathrm{RR}*}(\mathrm{\λ}+\mathrm{d\λ})=\left\{\begin{array}{c}w(g^{*}\left(\mathrm{\λ}\right)+\mathrm{d\λ}\·G(\mathrm{\λ},g^{*}\left(\mathrm{\λ}\right)))\\ +(1-w)(g^{*}(\mathrm{\λ}+2\mathrm{d\λ})-\mathrm{d\λ}\·G(\mathrm{\λ},g^{*}(\mathrm{\λ}+2\mathrm{d\λ}))),k_2\≠0\\ g^{\mathrm{IR}*}(\mathrm{\λ}+\mathrm{d\λ}),k_2=0\end{array}\right.---\left(3\right)$

(c)经过傅立叶逆变换(F^-1)，得到g^RR(λ+dλ)：

$g^{\mathrm{RR}}(\mathrm{\λ}+\mathrm{d\λ})=F^{-1}\left(g^{\mathrm{RR}*}(\mathrm{\λ}+\mathrm{d\λ})\right)---\left(4\right)$

(d)上述(b)(c)的算法迭代进行，以得到更加准确的恢复结果。

(e)记经n次迭代(n为大于1的整数，一般情况下取2或3即可保证精度)，对所得结果去加权(即除以加权因子1/|ξ-η|)得到可直接作为重建的线积分数据，记为g^～(λ+dλ)，也可将其重新变换到投影域，得到被遮挡元发分量的恢复结果I_P ^B～。

使PC-VI方法性能有效发挥的静止环状BSA的设计方法是指：根据该方法设计的静止环状BSA，其机械结构能够保证任意两个相邻投影角度下射线阻挡单元的遮挡区域在投影位置上没有重合或几乎没有重合。

本发明根据PC-VI条件，推导了通用螺旋轨迹CBCT结构下静态环状BSA的机械配置依据，由此依据可以设计出不同参数的静态环状BSA以适应不同配置的CBCT。在此框架下，静态环状BSA和移动BSA在确保PC-VI算法有效工作的意义上是完全等效的，由此得到了基于PC-VI算法的使用静态环状BSA进行散射测量的散射校正方法，从而构成了一种基于静止环状射线阻挡阵列的单次扫描同时获得散射测量值和原发射线分量的锥束CT扫描方法。

附图说明：

图1为本发明提出的锥束扫描方法的示意图。1为X光源点；2为静止环状BSA；3为环状BSA射线阻挡单元；4为环状BSA环面；5为扫描对；6为平板探测器；7为无遮挡像素；8为被遮挡像素。

图2为本发明提出的获取用于锥束重建的每个像素的原发射线投影(I_P ^{UB ～}与I_P ^B～)的流程示意图。

图3为锥束扫描参数图。X光源为ξ(r，λ，z)，检测器单元为

ξ包含的三个参数中，r代表ξ到旋转中心O的距离，λ为ξ的方位角，z为ξ的纵坐标。η包含的三个参数中，d代表旋转中心到η的距离，u、v分别代表η在面探测器上的横、纵坐标。

图4为静止环状BSA机械参数示意图，其中上图为空间立体图，下图为顶视图。X线光源点记为ξ₀，其轨迹为螺旋形，每环绕一周，上升的螺距记为p。X光源点的方位角记为λ。ξ₀到旋转中心o的距离为r，o到平面探测器中心的距离为d。图中o-XYZ形成全局坐标系，而平板探测器所在平面形成二维局部坐标系o′-UV。平板探测器在垂直于Z轴方向上的半宽度为h_U，其扫过区域的最大半径表示为H＝(d²+h_U ²)^0.5。环形BSA整体半径为l，均匀分布于其表面上的铅珠半径为σ，方位角为ψ。

图5为环状BSA具体实例，9，10分别以铅珠和铅条作为射线阻挡单元。

图6为环状BSA在探测器所在平面的投影，其中，竖直白线之间的部分[-ψ_m，ψ_m]表示包含在平板探测器范围内的投影。

图7为对用本文方法设计的一个环状BSA作投影的结果。从上到下依次为：当前角度投影图；相邻角度投影图；上述两投影叠加图；前一图像左上角框内部分放大图。

具体实施方式：

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

参见图1，使用图1所示的静止环状散射测量装置的锥束CT扫描结构时，假如没有BSA存在，投影测量数据I_A包含原发射线分量I_P和散射分量I_S，其中I_P为重建图像所需数据。如图2，在BSA存在的情况下，测量像素分两类，一类为未被BSA遮挡(Un-Blocking，UB)像素，其测量所得信号为I_A ^UB，包含原发射线分量和散射分量；另一类为有遮挡(Blocking，B)像素，其测量所得信号为I_S ^B，为散射分量在该位置的采样值。I_S ^B经空域插值得到任意像素I_S的估计值I_S ^～。对无遮挡像素，从I_A ^UB中减去I_S ^～即可得到原发射线分量计算值I_P ^UB～；对被遮挡像素，其原发射线I_P ^B～完全丢失，利用本发明所提出的PC-VI方法恢复被阻挡原发射线部分，记为I_P ^B～；将I_P ^UB～与I_P ^B～组合起来从而得到每个像素的原发射线投影，所述流程见图2。

下面结合图1～图3详细叙述使用PC-VI方法得到I_P ^B～的过程：

(a)首先，利用未遮挡探测器单元去散射后的投影I_P ^UB～对遮挡单元作空域内插得到其丢失原发射线信号的初步恢复值I_P ^B～IR，其中上标^IR表示初步恢复(Initial Restoration，IR)，该空域插值可以采用现有的方法。对所得整个投影(I_P ^UB～和I_P ^B～IR)作对数运算将其转换为加权的线积分，记为g^IR(λ+dλ)，其中λ+dλ表示当前投影角度，加权权重为1/|ξ-η|，上标^IR表示被变换为g之前的投影中含有不够准确的所谓初步恢复部分I_P ^B～IR。根据空域内插的特点，g^IR(λ+dλ)中初步恢复的是丢失信号的低频部分。

投影相关性的核心公式为：

$g_{\mathrm{\λ}}^{*}\≅\frac{j}{r+d}\·(-2\frac{k_1}{k_2}{g}_{k_1}^{*}+k_2{g}_{k_1{k}_2}^{*}-k_1{g}_{k_2{k}_2}^{*}),k_2\≠0---\left(1\right)$

其中g^*为g(λ，u，v)的傅立叶变换g^*(λ，k₁，k₂)的简化表示形式；符号

和

分别表示g^*对变量λ和k₁求偏导；

表示g对变量k₁，k₂求偏导；

表示g对变量k₂求两次偏导；定义(1)等号右端部分为G(λ，g^*(λ))，(1)可写成：

$g^{*}(\mathrm{\λ}+\mathrm{d\λ})\≅g^{*}\left(\mathrm{\λ}\right)+\mathrm{d\λ}\·G(\mathrm{\λ},g^{*}\left(\mathrm{\λ}\right)),k_2\≠0---\left(2\right)$

(b)通过(2)所示的投影相关性，恢复丢失信号频谱中k₂≠0之外的其他频段，尤其是高频部分可以从相邻投影的频谱g^*(λ)计算得到(同理也可用到另外一侧的相邻投影g^*(λ+2dλ))。在k₂＝0处则直接用

的结果。将左右两边相邻角度的结果以权重w相结合，于是得到上标^RR表示细化恢复(Refined Restoration，RR)：

$g^{\mathrm{RR}*}(\mathrm{\λ}+\mathrm{d\λ})=\left\{\begin{array}{c}w(g^{*}\left(\mathrm{\λ}\right)+\mathrm{d\λ}\·G(\mathrm{\λ},g^{*}\left(\mathrm{\λ}\right)))\\ +(1-w)(g^{*}(\mathrm{\λ}+2\mathrm{d\λ})-\mathrm{d\λ}\·G(\mathrm{\λ},g^{*}(\mathrm{\λ}+2\mathrm{d\λ}))),k_2\≠0\\ g^{\mathrm{IR}*}(\mathrm{\λ}+\mathrm{d\λ}),k_2=0\end{array}\right.---\left(3\right)$

(c)经过傅立叶逆变换(F^-1)，得到g^RR(λ+dλ)：

$g^{\mathrm{RR}}(\mathrm{\λ}+\mathrm{d\λ})=F^{-1}\left(g^{\mathrm{RR}*}(\mathrm{\λ}+\mathrm{d\λ})\right)---\left(4\right)$

(d)对上述(b)(c)的算法迭代进行，以得到更加准确的恢复结果。

(e)记经n次迭代(n一般为3，亦可为其他整数)，对所得结果去加权(除以加权权重1/|ξ-η|)得到可直接重建的线积分数据，记为g^～(λ+dλ)，也可将其重新变换到投影域，对应前述被遮挡信号的恢复结果I_P ^B～。

PC-VI方法非常有效，但其前提是，相邻角度投影被遮挡位置区域不互相重合，即所谓的PC-VI条件。为满足PC-VI条件，BSA必须移动。

申请人注意到目前所有的BSA都是平面式的，此类BSA被安装在X线出线口前，跟随X光源一起旋转。当平面式BSA静止时，其阴影落在探测器上的位置没有变化，因此无法满足PC-VI条件要求的相邻投影被遮挡的区域位置必须彼此不重叠。所以平面BSA在使用时必须运动，才能满足PC-VI条件。

根据这一观察，采用不随X光源旋转的BSA设计，例如，该BSA可以是表面均匀分布有铅珠阵列环状面，如图4所示。只要满足|OH|＜l＜r，即(u²+h_U ²)1/2＜l＜r，该环状BSA即可被固定嵌入X光源轨迹和探测器轨迹之间。需要强调的是，本文方法是用的射线阻挡单元不仅限于铅珠阵列，也可以使用其它阻挡单元，如图5所示。

与平面式跟随光源旋转的BSA不同，环状BSA在相邻角度的投影位置之间会产生位置的偏移，这就成为满足PC-VI条件的基础。在此基础上，我们进一步推导得到了确保满足PC-VI条件的环状BSA的参数设计方法。采用图1所示的参量标记，经推导获得设计方法如下：

$G(r,d,\mathrm{d\λ},h_U,h_V,l,\mathrm{\σ},\mathrm{\ψ},\mathrm{\ζ})$

$=l(\frac{\sin \mathrm{\ψ}}{r-l\cos \mathrm{\ψ}}-\frac{\sin (\mathrm{\ψ}-\mathrm{d\λ})}{r-l\cos (\mathrm{\ψ}-\mathrm{d\λ})})---\left(5\right)$

$-\mathrm{\σ}(\frac{\sqrt{r^2-2\mathrm{rl}\cos \mathrm{\ψ}+l^2+{\mathrm{\ζ}}^2}}{{\left|(r-l\cos \mathrm{\ψ})\right|}^2}+\frac{\sqrt{r^2-2\mathrm{rl}\cos (\mathrm{\ψ}-\mathrm{d\λ})+l^2{(\frac{p}{2\mathrm{\π}}\·\mathrm{d\λ}-\mathrm{\ζ})}^2}}{{|r-l\cos (\mathrm{\ψ}-\mathrm{d\λ})|}^2})>0$

其中，dλ表示相邻角度间隔，h_V表示平板探测器在Z轴方向的半高度值，ψ的取值从-ψ_m到ψ_m。ψ_m计算公式为

${\mathrm{\ψ}}_m=\arccos \frac{{\mathrm{rh}}^2+(r+d)\sqrt{l^2({(r+d)}^2+h^2)-r^2{h}^2}}{l({(r+d)}^2+h^2)}\·---\left(6\right)$

ζ代表铅珠的Z轴方向坐标，其取值从-ζ_m到ζ_m，而且有

${\mathrm{\ζ}}_m=h_V\frac{r-l}{r+d}\·---\left(7\right)$

在如(5)所示方法下，作为范例，设计含有400×8个铅珠的环状BSA，其中400列均匀分布在360度，再选取一个特定参数组合：r＝750，d＝337.5，l＝485，dλ＝2π/880，来进行有效性评估。

图6和图7显示的是模拟的装有环状BSA的锥束CT空扫(无成像物体)的情形。环状BSA在平板探测器所在平面的投影如图6所示(为了显示效果，经取对数运算后显示)，其中白色竖线之间为探测器可探测范围。图7前两图为任意相邻角度的两幅投影，为了测试两者中的被遮挡部分有无重合，特别将两者叠加，见第三图。第四图为第三图左上方四分之一的放大。由图7可见，使用此静止环状BSA，可以达到和移动BSA相同的效果，即任意相邻角度投影中因射线阻挡单元造成的阴影区域互不重叠，从而可以满足PC-VI条件。

在上述扫描模式下，加入扫描对象可视人模体，对其进行胸腔扫描，然后对全部投影进行误差计算，并与Yan H.和Mou X.等原文中的移动BSA情况进行对比。定义每个投影角度的绝对误差如下：

${\left|\mathrm{Error}\left(\mathrm{\λ}\right)\right|}_{\mathrm{PC}-\mathrm{VI}}=\frac{{\mathrm{\Σ}}_{u,v\∈B}|{\mathrm{View}}_{\mathrm{PC}-\mathrm{VI}}(\mathrm{\λ},u,v)-{\mathrm{View}}_0(\mathrm{\λ},u,v)}{{\mathrm{\Σ}}_{u,v\∈B}\left|{\mathrm{View}}_0(\mathrm{\λ},u,v)\right|}---\left(8\right)$

其中，B表示所有被遮挡像素的集合，View₀表示理想投影，View_PC-VI表示被遮挡投影经PC-VI方法恢复的结果。

进行误差对比使用的是平均的绝对误差：

${\left|\mathrm{Error}\right|}_{\mathrm{PC}-\mathrm{VI}}=\underset{\mathrm{\λ}}{\mathrm{mean}}\left\{{\left|\mathrm{Error}\left(\mathrm{\λ}\right)\right|}_{\mathrm{PC}-\mathrm{VI}}\right\}---\left(\text{9}\right)$

其中，mean表示取平均值。

实验结果如表1所示。

由表1可见，使用本文方法设计的环状BSA在PC-VI框架下完全可以达到与移动BSA相似甚至略好的原发投影恢复结果。

以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明，不能认定本发明的具体实施方式仅限于此，对于本发明所属技术领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干简单的推演或替换，都应当视为属于本发明由所提交的权利要求书确定专利保护范围。

表1与采用移动BSA的综合性能比较

Claims (2)

1.一种使用静止环状射线阻挡阵列的锥束CT扫描和散射校正方法，其特征在于：使用静止环状BSA的扫描机械结构进行锥束CT扫描与散射校正，使用基于投影相关性的角度内插PC-VI方法恢复被静止环状BSA阻挡的原发射线分量；

使用静止环状BSA的扫描机械结构进行锥束CT扫描：在通常的X光源点轨迹和平板探测器轨迹间内，加入静止的环状BSA用来测量散射；环状面材料为对X光低衰减的物质，阵列单元材料为铅或其他原子序数较高的对X光强衰减的物质；静止环状BSA的结构以及阵列材料的排列符合给定的要求：使PC-VI方法性能有效发挥的静止环状BSA的设计方法；

使用静止环状BSA进行散射校正：

(a)投影测量数据I_A包含原发射线分量I_P和散射分量I_S，其中I_P为重建图像所需数据；在使用环状BSA测量I_S时，测量像素分两类，一类为未被BSA遮挡像素，其测量所得信号为I_A ^UB，包含原发射线分量I_P和散射分量I_S；另一类为有遮挡像素，其测量所得信号为I_S ^B，为散射分量I_S在该位置的采样值；

(b)对I_S ^B做空域插值得到任意像素I_S的估计值I_S ^～；对无遮挡像素，从I_A ^UB中减去I_S ^～即可得到原发射线分量计算值I_P ^UB～；对被遮挡像素，其原发射线I_P ^B～完全丢失，利用PC-VI方法恢复被阻挡原发射线分量，记为I_P ^B～；I_P ^{UB ～}与I_P ^B～的集合I_P ^～，即对CT重建所需原发射线投影I_P的计算结果；

使用PC-VI方法恢复被静止环状BSA阻挡的原发射线分量：

(a)首先，利用未遮挡探测器单元去除散射后的投影对遮挡单元作空域内插得到其丢失原发射线信号的初步恢复值

其中上标^IR表示初步恢复，该空域插值采用现有的方法，如水平方向的三次样条插值；对所得整个投影

和

作对数运算将其转换为线积分并加权，记为g^IR(λ+dλ)，其中λ+dλ表示当前投影角度，乘性加权因子为1/|ξ-η|，ξ和η为光源和探测器坐标，上标^IR表示被变换为g之前的投影中含有不够准确的所谓初步恢复部分根据空域内插的特点，g^IR(λ+dλ)中初步恢复的是丢失信号的低频部分；

投影相关性的核心公式为：

$g_{\mathrm{\λ}}^{*}\≅\frac{j}{r+d}\·(-2\frac{k_1}{k_2}{g}_{k_1}^{*}+k_2{g}_{k_1{k}_2}^{*}-k_1{g}_{k_2{k}_2}^{*}),k_2\≠0---\left(1\right);$

其中g^*为g(λ，u，v)的傅立叶变换g^*(λ，k₁，k₂)的简化表示形式；符号

和

分别表示g^*对变量λ和k₁求偏导；表示g对变量k₁、k₂求偏导；

表示g对变量k₂求两次偏导；涉及的锥束扫描参数罗列如下：X光源为ξ(r，λ，z)，检测器单元为η(d，u，v)；ξ包含的三个参数中，r代表ξ到旋转中心O的距离，λ为ξ的方位角，z为ξ的纵坐标；η包含的三个参数中，d代表旋转中心到η的距离，u、v分别代表η在面探测器上的横、纵坐标；

定义(1)等号右端部分为G(λ，g^*(λ))，(1)可写成：

$g^{*}(\mathrm{\λ}+\mathrm{d\λ})\≅g^{*}\left(\mathrm{\λ}\right)+\mathrm{d\λ}\·G(\mathrm{\λ},g^{*}\left(\mathrm{\λ}\right)),$ k₂≠0    (2)

(b)通过如(2)所示的方程，恢复丢失信号频谱中k₂≠0之外的其他频段，高频部分可以从相邻投影的频谱g^*(λ)计算得到，同理也可用到另外一侧的相邻投影g^*(λ+2dλ)；在k₂＝0处则直接用g^IR*(λ+dλ)的结果；将左右两边相邻角度的结果以权重w相结合，w的取值范围为大于等于0，小于等于1；于是得到g^RR*(λ+dλ)，上标^RR表示细化恢复：

$g^{\mathrm{RR}*}(\mathrm{\λ}+\mathrm{d\λ})=\left\{\begin{array}{c}w(g^{*}\left(\mathrm{\λ}\right)+\mathrm{d\λ}\·G(\mathrm{\λ},g^{*}\left(\mathrm{\λ}\right)))\\ +(1-w)(g^{*}(\mathrm{\λ}+2\mathrm{d\λ})-\mathrm{d\λ}\·G(\mathrm{\λ},g^{*}(\mathrm{\λ}+2\mathrm{d\λ}))),k_2\≠0\\ g^{\mathrm{IR}*}(\mathrm{\λ}+\mathrm{d\λ}),k_2=0\end{array}\right.---\left(3\right);$

(c)经过傅立叶逆变换(F^-1)，得到g^RR(λ+dλ)：

g^RR(λ+dλ)＝F^-1(g^RR*(λ+dλ))    (4)

(d)对上述(b)(c)的算法迭代进行，以得到更加准确的恢复结果；

(e)记经n次迭代，n的取值范围为大于等于1的正整数；对所得结果除以加权因子1/|ξ-η|得到可直接作为重建的线积分数据，记为g^～(λ+dλ)，或将其重新变换到投影域，得到被遮挡原发射线分量的恢复结果I_P ^B～。

2.如权利要求1所述使用静止环状射线阻挡阵列的锥束CT扫描和散射校正方法；其特征在于：根据该方法设计的静止环状BSA，其机械结构能够保证任意两个相邻投影角度下射线阻挡单元的遮挡区域在投影位置上没有重合或几乎没有重合。

Images