CN106651982B - 一种基于阵列x射线源和探测器的ct图像重建方法 - Google Patents
一种基于阵列x射线源和探测器的ct图像重建方法 Download PDFInfo
- Publication number
- CN106651982B CN106651982B CN201611170283.XA CN201611170283A CN106651982B CN 106651982 B CN106651982 B CN 106651982B CN 201611170283 A CN201611170283 A CN 201611170283A CN 106651982 B CN106651982 B CN 106651982B
- Authority
- CN
- China
- Prior art keywords
- mrow
- msub
- msubsup
- detector
- radiographic source
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T11/00—2D [Two Dimensional] image generation
- G06T11/003—Reconstruction from projections, e.g. tomography
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N23/00—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00
- G01N23/02—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by transmitting the radiation through the material
- G01N23/04—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by transmitting the radiation through the material and forming images of the material
- G01N23/046—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by transmitting the radiation through the material and forming images of the material using tomography, e.g. computed tomography [CT]
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N2223/00—Investigating materials by wave or particle radiation
- G01N2223/40—Imaging
- G01N2223/419—Imaging computed tomograph
Landscapes
- Theoretical Computer Science (AREA)
- Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Physics & Mathematics (AREA)
- Health & Medical Sciences (AREA)
- Radiology & Medical Imaging (AREA)
- Pulmonology (AREA)
- Life Sciences & Earth Sciences (AREA)
- Chemical & Material Sciences (AREA)
- Analytical Chemistry (AREA)
- Biochemistry (AREA)
- General Health & Medical Sciences (AREA)
- Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
- Immunology (AREA)
- Pathology (AREA)
- Apparatus For Radiation Diagnosis (AREA)
- Measurement Of Radiation (AREA)
Abstract
本发明公开了一种基于阵列X射线源和探测器的CT图像重建方法，目的在于，解决了传统CT因单锥角原因远离成像对象的缺点，能够更加贴近成像对象，避免不必要的剂量辐射，扫描快速，成本低，小巧灵便，无机械移动，易于实现，使用多个阵列上多个X射线源发射锥型射线，使用探测器接收光子，当分时点亮阵列上的每一个X射线源时，分时发射的每个射线源对应的一个锥角只能覆盖部分成像对象，得到分时的投影图像，使用分时阵列多源CT的重建方法进行重建；当同时点亮阵列上的所有X射线源时，同时发出的射线将成像对象全部覆盖，多个射线源扫描区域存在覆盖，探测器接收的投影会存在混叠部分，基于混叠的阵列多源CT的重建方法进行重建。
Description
技术领域
本发明涉及X射线成像技术领域，具体涉及一种基于阵列X射线源和探测器的CT图像重建方法。
背景技术
随着X射线成像技术的广泛应用，人们对降低X射线成像剂量、提高成像分辨率和成像速度的要求越来越高。传统CT成像技术受限于热阴极X射线管和机械扫描速度的限制，迫使人们寻找新的成像方式替代。值得鼓舞的是，阵列X射线源正开始朝着冷阴极和平板化的方向发展，有望实现X射线源从点光源到可寻址平板光源的变革，同时，由于X射线源能够与被测目标紧密耦合，避免对相关区域的辐射，有望依此降低剂量并缩小成像系统体积，这对于发展新型X射线成像方法，降低X射线检查带来的危害，造福人类的健康有重要的意义。
自2001年日本名古屋工业大学最先报导了以碳纳米管(CNT)作为电子源的X射线管以来，采用CNT等一维纳米材料冷阴极的分立型X射线管成为研究热点，其研究主要集中在微焦X射线源及其在动态成像和分布式X光源CT系统的应用，其中美国北卡罗来纳大学研制的分布式冷阴极X射线管乳腺CT原理型样机已经进入临床验证阶段，韩国VSI公司和日本大阪大学分别报道了非寻址的平板X射线源，美国加州大学洛杉矶分校提出采用热释电晶体阴极制作可寻找的平板X射线源的想法。而到了2015年，中山大学已报导了较大面积的氧化锌纳米冷阴极平板X射线源，并实现了小于25微米的静态成像。这一系列研究成果都为利用X射线源阵列探索新的成像方法与CT图像重建方法铺平了道路。
较于传统CT，阵列X射线源CT每个角度使用多个射线源，相同锥角下，设备体积 更小；阵列X射线源CT还固定射线源位置，扫描时无机械移动，不用考虑高速移动下带来的投影数据矫正问题；此外，无论是分时投影与重建，还是同时投影与重建，相比传统CT，其扫描速度更快，剂量更小，造成的剂量泄露也更少。
发明内容
为了解决现有技术中的问题，本发明提出一种基于阵列X射线源和探测器的CT图像重建方法，使用多个阵列上的多个X射线源进行锥束CT扫描，解决了传统CT因单锥角原因远离成像对象的缺点，能够更加贴近成像对象，避免不必要的剂量辐射，具有扫描快速，成本低，小巧灵便，无机械移动及易于实现的优点。
为了实现以上目的，本发明所采用的技术方案为：包括以下步骤：
1)采用正多边形状的阵列X射线源与探测器对目标进行分时或同时扫描获得光子投影数据，正多边形状的阵列X射线源包括n个阵列，每个阵列有Q个射线源；
2)根据分时扫描光子投影数据建立目标函数：
其中，A表示从图像x到投影数据b间的线性变换，是一个M×N矩阵，且b∈RM为观察数据，x∈RN为重建图像，R(x)表示对图像的罚项约束，β表示正则化因子；
同时扫描光子投影数据建立目标函数：
其中，apq表示第q个射线源的对应第p个探测器的系统矩阵，yp表示第p个探测器测量到的信号与对应的散射分量的和，Ipq表示射线源q对应于探测器p方向的入射光子强度；
3)对分时扫描或同时扫描目标函数进行求解，从而得到CT重建图像。
所述步骤1)中分时扫描包括：首先点亮所有阵列的第1个射线源，探测器收集到光子投影数据后再点亮所有阵列的第2个射线源，直至点亮所有阵列的第Q个射线源，探 测器并收集光子投影数据；或者首先点亮所有阵列的第1、m+1、…、tm+1(tm+1≤Q)个射线源，探测器收集到光子投影数据后再点亮所有阵列的第2、m+2、…、tm+2(tm+2≤Q)个射线源，直至点亮所有阵列的第m、2m、…、tm+m(tm+m≤Q)个射线源，探测器接收光子投影数据；
同时扫描包括：同时使用所有阵列上的所有X线源，各个阵列上的所有射线源发出锥束射线，探测器接收到混叠的光子投影数据。
所述步骤1)中X射线源呈阵列设置，X射线源有相同的张角和发射能谱。
所述步骤3)中分时扫描目标函数采用经典统计迭代法重建框架，设目标函数为f(x)，根据迭代求解公式：xn+1＝xn+λf′(xn)或迭代求解目标函数，得到CT重建图像，其中xn+1表示第n+1次迭代结果，xn表示第n次迭代结果，λ表示迭代步长，f′(xn)表示目标函数的一阶导，f″(xn)表示目标函数的二阶导。
所述步骤3)中同时扫描目标函数采用经典统计迭代法重建框架，设目标函数为：
其中对dp求解对应一阶导数为：
二阶导数为：
并根据代求解公式：xn+1＝xn+λf′(xn)或迭代求解目标函数，得到CT重建图像，其中，apq表示射线源q对应于探测器p的系统矩阵，yp表示探测器p收到的总光子强度，Ipi表示从第i个射线源到第p个探测器的光子强度，api表示从第i个射线 源到第p个探测器对应的系统矩阵，Ipk表示从第k个射线源到第p个探测器的光子强度，apk表示从第k个射线源到第p个探测器对应的系统矩阵，x表示重建图像。
所述同时扫描目标函数采用泰勒展开求解：
Ipqexp(-apqx)经泰勒展开为：
代入公式得到：
dp一阶导数为：
二阶导数为：
将dp一阶导数和二阶导数代入求解公式：xn+1＝xn+λf′(xn)或迭代求解目标函数，得到CT重建图像，其中，x0表示泰勒展开前的重建图像，T表示转置。
与现有技术相比，本发明使用多个阵列上多个X射线源发射锥型射线，再使用探测器进行接收光子投影数据。当分时点亮阵列上的每一个X射线源时，分时发射的每个射线源对应的一个锥角只能覆盖部分成像对象，得到分时的投影图像，使用分时阵列多源CT的重建方法进行重建CT图像；当同时点亮阵列上的所有X射线源时，同时发出的射线将成像对象全部覆盖，但由于多个射线源扫描区域存在覆盖，探测器接收的投影会存在混叠部分，于是提出了基于混叠的阵列多源CT的重建方法进行重建。本发明使用多个阵列上的多个X射线源进行锥束CT扫描，解决了传统CT因单锥角原因远离成像对象的缺点，能够更加贴近成像对象，避免不必要的剂量辐射，具有扫描快速，成本低，小巧灵便，无机械移动及易于实现的优点。
附图说明
图1为阵列X射线源与探测器结构示意图，1为成像目标，2为隔离源与探测器的铅板，3为阵列射线源，4为阵列探测器；
图2为阵列X射线源与探测器剖面示意图；
图3为阵列X射线源与探测器的摆放模型；
图4a为分时扫描下的第一组射线源扫描示意图；图4b为第二组射线源扫描示意图；
图5为同时扫描下的扫描示意图；
图6为阵列X射线源上的多个射线源分时投影流程图；
图7为阵列X射线源上的多个射线源同时投影流程图；
图8为阵列X射线源上的多个射线源同时投影下的重建算法流程图；
图9为示例所用原始模体，外层大圆为人体组织，方形为软骨，两小圆为硬骨；
图10a为使用分时扫描并重建的图像，图10b为使用同时扫描并重建的图像。
具体实施方式
下面结合具体实施例和说明书附图对本发明作进一步的解释说明。
参见图1，使用图1所示的阵列X射线源与探测器CT扫描目标内部结构时，物体置于正中，周围是n个阵列构成围绕。扫描阵列由上下两部分组成，上半部分为Q个成矩形阵列的X射线源，下半部分为P个探测器组。(a)使用分时模式扫描时，所有阵列上的同一对应位置的第i个射线源发出X射线，探测器接受到光子后，所有阵列上的同一对应位置的第i+1个源发出X射线，直至所有X射线源扫描完成；或者首先点亮所有阵列的第1、m+1、…、tm+1(tm+1≤Q)个射线源，探测器收集到光子投影数据后再点亮所有阵列的第2、m+2、…、tm+2(tm+2≤Q)个射线源，直至点亮所有阵列的第m、2m、…、tm+m(tm+m≤Q)个射线源，探测器接收光子投影数据；(b)使用混合模式扫描时，即同时扫描，所有阵列上的多个源同时发出X射线，为保证对目标全扫描，相邻源的扫描 区域会相互覆盖，使对面的阵列上的探测器接收到光子投影，得到混叠的投影数据，最后再使用本发明提出的CT图像重建方法重建出对应图像并输出，选择如图9所示的原始模体，外层大圆为人体组织，方形为软骨，两小圆为硬骨。多个X射线源进行锥束扫描，射线源呈阵列摆放，具有固定的相同的张角和发射能谱，为保证对成像对象进行全扫描，相邻锥束间扫描范围存在重叠覆盖。使用阵列射线源与探测器组合模块时可以根据检测对象形状规格的不同，调节射线源阵列单元与探测器组合模块单元的数量及位置，实现对不同大小规格待检测区域的扫描覆盖。
下面结合图1、图2、图3、图4a、图4b和图6详述分时投影对应的重建方法：
(a)首先，根据扫描物体的具体大小和所需重建图像的内容，布置n个阵列成正n边形分布，每个阵列上半部分为Q个X射线源成面阵或线阵排列，下半部分为P个探测器成面阵排列，n个阵列将扫描物体包围在中央；
(b)使用分时扫描方法一：首先点亮所有阵列的第1个射线源，收集到光子投影数据后再点亮所有阵列的第2个射线源，直至点亮所有阵列的第Q个射线源并收集数据；
分时扫描方法二：在保证光子投影不产生混叠的情况下，首先点亮所有阵列的第1、m+1、…、tm+1(tm+1≤Q)个射线源，收集到光子投影数据后再点亮所有阵列的第2、m+2、…、tm+2(tm+2≤Q)个射线源，直至点亮所有阵列的第m、2m、…、tm+m(tm+m≤Q)个射线源并收集数据；
(c)对于任意一种分时扫描方法，阵列上每个X射线源发出锥束X射线，对于在对面第p个探测器接收到的光子强度：
其中q表示射线源的位置，gpq表示从射线源q至探测器p的线积分，即 apqn表示由射线源q和探测器p形成的X线光束与体素xn之间交集的体积，Ipq表示射线源q对应探测器p方向的入射光子强度，npq表示第q个射线源在探测器p上造 成的噪声强度，根据Bear定理得到含噪投影数据而又根据Radon变换有考虑到各个阵列各个源对应的投影数据有模型b＝Ax+n，再建立目标函数：
其中，A表示从图像x到投影数据b间的线性变换，是一个M×N矩阵，且b∈RM为观察数据，x∈RN为重建图像，R(x)表示对图像的罚项约束，β表示正则化因子；
接下来再对目标函数进行求解，采用经典统计迭代法重建框架，设目标函数为f(x)，根据迭代求解公式：xn+1＝xn+λf′(xn)或迭代求解目标函数，得到如图10a所示的CT重建图像，其中，xn+1表示第n+1次迭代结果，xn表示第n次迭代结果，λ表示迭代步长，f′(xn)表示目标函数的一阶导，f″(xn)表示目标函数的二阶导。
上述使用分时扫描方案获得投影的过程，获得的投影图像由一系列来自不同射线源或不同射线源组的投影图像，每个射线源发出的锥束射线仅穿过被检测物体的一部分，所有射线源发出的锥束在空间上完全覆盖完整的被检测部分，其与传统CT需要被检测物与射线源保持足够远距离才能成像的模式不同，本发明可以明显缩小成像系统体积，降低辐射泄露，扩大X成像的应用场景。同时X射线源阵列分组扫描的方式使得投影图像的获取更加快速，减小了被检测物体可能的移动为成像结果带来的影响。图像重建时根据需要添加不同的正则项也能弥补角度数不足的缺点。
下面结合图1、图2、图3、图5、图7和图8详述同时投影及其对应重建方法：
(a)使用混合模式同时扫描时，每个X射线源阵列上的Q个射线源单元，在不放置被探测物体时(空扫描)发出锥束X射线，其中第q个射线源在探测器单元p上得到光子强度Ipq，在放置被检测物体时，同时发出锥束，探测器单元p上得到光子强度yp，为各射线源发出射线穿过不同衰减路径之后到达p处光子强度之和，在上述扫描模式下，得到的yp与目前单X线点源成像系统不同，除去散射和噪声等成分，yp是一个由多点源产生的多条射线积分的总和，是一个模糊混叠的X线成像，而单点源系统的成像只与一条射线积分相关，是一个有相对清晰结构的图像，使用同时扫描方案可以只一次扫描快速的得到被检测物的X射线投影图像，极大的提升了成像的速度，避免了成像过程中被检测物体移动对成像质量的影响，但是由于图像的混叠模糊，进行CT图像重建需要新的方法。对于在对面第p个探测器接收到的光子强度：
形成最终的模型：
其中apq表示射线源q对应于探测器p的系统矩阵，np表示位置为p的探测器接收的噪声分量；
此时的目标函数为：
(b)使用第一种混合重建的方法直接求取一阶导数和二阶导数：
对每一项对x的一阶导为：
二阶导为：
最后迭代求解公式为：
其中p＝1,2,3,…,P，apq表示射线源q对应于探测器p的系统矩阵，yp表示探测器p收到的总光子强度，Ipi表示从第i个射线源到第p个探测器的光子强度，api表示从第i个射线源到第p个探测器对应的系统矩阵，Ipk表示从第k个射线源到第p个探测器的光子强度，apk表示从第k个射线源到第p个探测器对应的系统矩阵，x表示重建图像，如此循环迭代求得重建图像。
使用第二种混合重建的方法直接求取一阶导数和二阶导数：
Ipqexp(-apqx)在x0处经泰勒展开为：
代入公式(4)有：
其一阶导数为：
其二阶导为：
最后的迭代公式同公式(7)；
添加先验信息，修改目标函数为：
重建方法仍使用牛顿迭代法，参照公式(6)推导其迭代公式为：
以此迭代求解目标函数，x0表示泰勒展开前的重建图像，T表示转置，得到如图10b所示的同时投影模型下的CT重建图像。
使用以上所述的两种投影模型与对应重建方法，作为范例，选取一组简化的特定参数：简化锥束为扇束模型，共n＝20(或n＝10)个阵列，每个阵列上摆放一行共S＝8个源(相距d＝5mm)的X射线源发射扇束，每个扇束扇θ＝4.57°，阵列间相距100mm，中间放置的目标像素大小为256×256，每个像素长宽均为0.25mm，探测器为P＝750×1，每个探测器单元长宽均为0.1mm，以此参数来评估算法的有效性。
定义重建图像的误差为：
其中，m,n为图像长和宽，I为真实目标图像，I′为重建目标图像。
实验结果如下表所示：
由表可见，使用本文方法提出的重建方法在分时投影与同时投影下基本可以达到令人满意的图像重建结果。其中分时模式重建效果更好，而同时扫描模式扫描速度更快。
本发明通过使用一种的无机械移动的阵列X射线源与探测器模型，在这种系统框架下，由于固定存在多个射线源与探测器，因此整个系统在扫描过程中其无需移动，其实现的结构更加简单，成本更加低廉，操作性及维护性更佳，而且探测器与源可移动，即扫描前可以通过改变探测器与源距离目标的远近更个性化地修正扫描参数，探测器与源也可拆卸，即可以更方便地根据需求进行改变。
以上内容是结合具体的实施方式对本发明所作的进一步详细说明，不能认定本发明的具体实施方式仅限于此，对于本发明所属技术领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干简单的推演或替换，都应当视为属于本发明由所提交的权利要求书确定专利保护范围。
Claims (5)
1.一种基于阵列X射线源和探测器的CT图像重建方法，其特征在于，包括以下步骤：
1)采用正多边形状的阵列X射线源与探测器对目标进行分时或同时扫描获得光子投影数据，正多边形状的阵列X射线源包括n个阵列，每个阵列有Q个射线源；
2)根据分时扫描光子投影数据建立目标函数：
<mrow>
<munder>
<mrow>
<mi>m</mi>
<mi>i</mi>
<mi>n</mi>
</mrow>
<mi>x</mi>
</munder>
<mo>|</mo>
<mo>|</mo>
<mi>A</mi>
<mi>x</mi>
<mo>-</mo>
<mi>b</mi>
<mo>|</mo>
<msubsup>
<mo>|</mo>
<mn>2</mn>
<mn>2</mn>
</msubsup>
<mo>+</mo>
<mi>&beta;</mi>
<mi>R</mi>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>)</mo>
</mrow>
</mrow>
其中，A表示从图像x到投影数据b间的线性变换，是一个M×N矩阵，且b∈RM为观察数据，x∈RN为重建图像，R(x)表示对图像的罚项约束，β表示正则化因子；
同时扫描光子投影数据建立目标函数：
<mrow>
<munder>
<mrow>
<mi>m</mi>
<mi>i</mi>
<mi>n</mi>
</mrow>
<mi>x</mi>
</munder>
<mo>|</mo>
<mo>|</mo>
<msubsup>
<mi>&Sigma;</mi>
<mrow>
<mi>p</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>P</mi>
</msubsup>
<mrow>
<mo>(</mo>
<msubsup>
<mi>&Sigma;</mi>
<mrow>
<mi>q</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>Q</mi>
</msubsup>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>exp</mi>
<mo>(</mo>
<mrow>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>x</mi>
</mrow>
<mo>)</mo>
<mo>-</mo>
<msub>
<mi>y</mi>
<mi>p</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>|</mo>
<msubsup>
<mo>|</mo>
<mn>2</mn>
<mn>2</mn>
</msubsup>
<mo>+</mo>
<mi>&beta;</mi>
<mi>R</mi>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>)</mo>
</mrow>
</mrow>
其中，apq表示第q个射线源对应第p个探测器的系统矩阵，yp表示第p个探测器测量到的信号与对应的散射分量的和，Ipq表示射线源q对应于探测器p方向的入射光子强度；
3)对上述分时扫描或同时扫描目标函数进行求解，从而得到CT重建图像；
同时扫描目标函数采用经典统计迭代法重建框架，设目标函数为：
<mrow>
<mi>f</mi>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mo>|</mo>
<mo>|</mo>
<msubsup>
<mi>&Sigma;</mi>
<mrow>
<mi>p</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>P</mi>
</msubsup>
<mrow>
<mo>(</mo>
<msubsup>
<mi>&Sigma;</mi>
<mrow>
<mi>q</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>Q</mi>
</msubsup>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>exp</mi>
<mo>(</mo>
<mrow>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>x</mi>
</mrow>
<mo>)</mo>
<mo>-</mo>
<msub>
<mi>y</mi>
<mi>p</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>|</mo>
<msubsup>
<mo>|</mo>
<mn>2</mn>
<mn>2</mn>
</msubsup>
<mo>+</mo>
<mi>&beta;</mi>
<mi>R</mi>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>p</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>P</mi>
</munderover>
<msub>
<mi>d</mi>
<mi>p</mi>
</msub>
<mo>+</mo>
<mi>&beta;</mi>
<mi>R</mi>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>)</mo>
</mrow>
</mrow>
其中，对dp求解对应一阶导数为：
<mrow>
<mfrac>
<mrow>
<mi>d</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>d</mi>
<mi>p</mi>
</msub>
<mo>)</mo>
</mrow>
</mrow>
<mrow>
<mi>d</mi>
<mi>x</mi>
</mrow>
</mfrac>
<mo>=</mo>
<mn>2</mn>
<msubsup>
<mi>&Sigma;</mi>
<mrow>
<mi>i</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>Q</mi>
</msubsup>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>i</mi>
</mrow>
</msub>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>i</mi>
<mi>r</mi>
</mrow>
</msub>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>i</mi>
</mrow>
</msub>
<mi>x</mi>
<mo>)</mo>
</mrow>
<mo>*</mo>
<mrow>
<mo>(</mo>
<msubsup>
<mi>&Sigma;</mi>
<mrow>
<mi>q</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>Q</mi>
</msubsup>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>exp</mi>
<mo>(</mo>
<mrow>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>x</mi>
</mrow>
<mo>)</mo>
<mo>-</mo>
<msub>
<mi>y</mi>
<mi>p</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>;</mo>
</mrow>
二阶导数为：
<mrow>
<mtable>
<mtr>
<mtd>
<mrow>
<mfrac>
<mrow>
<mi>d</mi>
<msup>
<mrow>
<mo>(</mo>
<msub>
<mi>d</mi>
<mi>p</mi>
</msub>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
</mrow>
<mrow>
<msup>
<mi>dx</mi>
<mn>2</mn>
</msup>
</mrow>
</mfrac>
<mo>=</mo>
<mn>2</mn>
<msubsup>
<mi>&Sigma;</mi>
<mrow>
<mi>i</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>Q</mi>
</msubsup>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>i</mi>
</mrow>
</msub>
<mo>|</mo>
<mo>|</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>i</mi>
</mrow>
</msub>
<mo>|</mo>
<msup>
<mo>|</mo>
<mn>2</mn>
</msup>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>i</mi>
</mrow>
</msub>
<mi>x</mi>
<mo>)</mo>
</mrow>
<mo>*</mo>
<mrow>
<mo>(</mo>
<msubsup>
<mi>&Sigma;</mi>
<mrow>
<mi>q</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>Q</mi>
</msubsup>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>exp</mi>
<mo>(</mo>
<mrow>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>x</mi>
</mrow>
<mo>)</mo>
<mo>-</mo>
<msub>
<mi>y</mi>
<mi>p</mi>
</msub>
<mo>+</mo>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>i</mi>
</mrow>
</msub>
<mi>exp</mi>
<mo>(</mo>
<mrow>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>i</mi>
</mrow>
</msub>
<mi>x</mi>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
<mo>+</mo>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mn>4</mn>
<msubsup>
<mi>&Sigma;</mi>
<mrow>
<mi>i</mi>
<mo>&NotEqual;</mo>
<mi>k</mi>
</mrow>
<mi>Q</mi>
</msubsup>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>i</mi>
</mrow>
</msub>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>k</mi>
</mrow>
</msub>
<mo>|</mo>
<mo>|</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>i</mi>
</mrow>
</msub>
<mo>|</mo>
<mo>|</mo>
<mo>|</mo>
<mo>|</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>k</mi>
</mrow>
</msub>
<mo>|</mo>
<mo>|</mo>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>i</mi>
</mrow>
</msub>
<mi>x</mi>
<mo>)</mo>
</mrow>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>k</mi>
</mrow>
</msub>
<mi>x</mi>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
</mtable>
<mo>;</mo>
</mrow>
并根据代求解公式：xn+1＝xn+λf′(xn)或迭代求解目标函数，得到CT重建图像，其中，apq表示射线源q对应于探测器p的系统矩阵，yp表示探测器p收到的总光子强度，Ipi表示从第i个射线源到第p个探测器的光子强度，api表示从第i个射线源到第p个探测器对应的系统矩阵，Ipk表示从第k个射线源到第p个探测器的光子强度，apk表示从第k个射线源到第p个探测器对应的系统矩阵，x表示重建图像。
2.根据权利要求1所述的一种基于阵列X射线源和探测器的CT图像重建方法，其特征在于，所述步骤1)中分时扫描包括：首先点亮所有阵列的第1个射线源，探测器收集到光子投影数据后再点亮所有阵列的第2个射线源，直至点亮所有阵列的第Q个射线源，探测器并收集光子投影数据；或者首先点亮所有阵列的第1、m+1、…、tm+1(tm+1≤Q)个射线源，探测器收集到光子投影数据后再点亮所有阵列的第2、m+2、…、tm+2(tm+2≤Q)个射线源，直至点亮所有阵列的第m、2m、…、tm+m(tm+m≤Q)个射线源，探测器接收光子投影数据；
同时扫描包括：同时使用所有阵列上的所有X线源，各个阵列上的所有射线源发出锥束射线，探测器接收到混叠的光子投影数据。
3.根据权利要求2所述的一种基于阵列X射线源和探测器的CT图像重建方法，其特征在于，所述步骤1)中X射线源呈阵列设置，X射线源有相同的张角和发射能谱。
4.根据权利要求1所述的一种基于阵列X射线源和探测器的CT图像重建方法，其特征在于，所述步骤3)中分时扫描目标函数采用经典统计迭代法重建框架，设目标函数为f(x)，根据迭代求解公式：xn+1＝xn+λf′(xn)或迭代求解目标函数，得到CT重建图像，其中xn+1表示第n+1次迭代结果，xn表示第n次迭代结果，λ表示迭代步长，f′(xn)表示目标函数的一阶导，f″(xn)表示目标函数的二阶导。
5.根据权利要求1所述的一种基于阵列X射线源和探测器的CT图像重建方法，其特征在于，所述同时扫描目标函数采用泰勒展开求解：
Ipqexp(-apqx)在x0处经泰勒展开为：
<mrow>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>x</mi>
<mo>)</mo>
</mrow>
<mo>&cong;</mo>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
<mo>)</mo>
</mrow>
<mo>+</mo>
<msup>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>-</mo>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
<mo>)</mo>
</mrow>
<mi>T</mi>
</msup>
<mrow>
<mo>(</mo>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>exp</mi>
<mo>(</mo>
<mrow>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
</mrow>
代入公式得到：
<mrow>
<msub>
<mi>d</mi>
<mi>p</mi>
</msub>
<mo>=</mo>
<mo>|</mo>
<mo>|</mo>
<msubsup>
<mi>&Sigma;</mi>
<mrow>
<mi>q</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>Q</mi>
</msubsup>
<msup>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>-</mo>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
<mo>)</mo>
</mrow>
<mi>T</mi>
</msup>
<mrow>
<mo>(</mo>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>exp</mi>
<mo>(</mo>
<mrow>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mrow>
<mo>(</mo>
<msubsup>
<mi>&Sigma;</mi>
<mrow>
<mi>q</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>Q</mi>
</msubsup>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>exp</mi>
<mo>(</mo>
<mrow>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>x</mi>
</mrow>
<mo>)</mo>
<mo>-</mo>
<msub>
<mi>y</mi>
<mi>p</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>|</mo>
<msubsup>
<mo>|</mo>
<mn>2</mn>
<mn>2</mn>
</msubsup>
</mrow>
而dp一阶导数为：
<mrow>
<mtable>
<mtr>
<mtd>
<mrow>
<mfrac>
<mrow>
<mi>d</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>d</mi>
<mi>p</mi>
</msub>
<mo>)</mo>
</mrow>
</mrow>
<mrow>
<mi>d</mi>
<mi>x</mi>
</mrow>
</mfrac>
<mo>=</mo>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mn>2</mn>
<mrow>
<mo>(</mo>
<msubsup>
<mi>&Sigma;</mi>
<mrow>
<mi>q</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>Q</mi>
</msubsup>
<msup>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>-</mo>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
<mo>)</mo>
</mrow>
<mi>T</mi>
</msup>
<mo>(</mo>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>exp</mi>
<mo>(</mo>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mrow>
<mo>(</mo>
<msubsup>
<mi>&Sigma;</mi>
<mrow>
<mi>q</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>Q</mi>
</msubsup>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>exp</mi>
<mo>(</mo>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>x</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msub>
<mi>y</mi>
<mi>p</mi>
</msub>
<mo>)</mo>
<mo>)</mo>
<mrow>
<mo>(</mo>
<msubsup>
<mi>&Sigma;</mi>
<mrow>
<mi>q</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>Q</mi>
</msubsup>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<msub>
<mi>I</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<mi>exp</mi>
<mo>(</mo>
<mo>-</mo>
<msub>
<mi>a</mi>
<mrow>
<mi>p</mi>
<mi>q</mi>
</mrow>
</msub>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
<mo>)</mo>
</mrow>
<mo>)</mo>
</mrow>
</mtd>
</mtr>
</mtable>
<mo>;</mo>
</mrow>
二阶导数为：
将dp一阶导数和二阶导数代入求解公式：xn+1＝xn+λf′(xn)或迭代求解目标函数，得到CT重建图像，其中，x0表示泰勒展开前的重建图像，T表示转置。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201611170283.XA CN106651982B (zh) | 2016-12-16 | 2016-12-16 | 一种基于阵列x射线源和探测器的ct图像重建方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201611170283.XA CN106651982B (zh) | 2016-12-16 | 2016-12-16 | 一种基于阵列x射线源和探测器的ct图像重建方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN106651982A CN106651982A (zh) | 2017-05-10 |
| CN106651982B true CN106651982B (zh) | 2018-04-17 |
Family
ID=58823197
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201611170283.XA Active CN106651982B (zh) | 2016-12-16 | 2016-12-16 | 一种基于阵列x射线源和探测器的ct图像重建方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN106651982B (zh) |
Families Citing this family (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108577876B (zh) * | 2018-02-28 | 2020-10-27 | 西安交通大学 | 一种多边形静止ct及其工作方法 |
| CN109991247A (zh) * | 2018-11-27 | 2019-07-09 | 姚智伟 | 基于平板x射线源阵列的x射线成像系统及扫描成像方法 |
| CN111265231B (zh) * | 2019-04-15 | 2021-08-31 | 清华大学 | 分布式光源ct图像重建方法与系统 |
| CN111991015B (zh) * | 2020-08-13 | 2024-04-26 | 上海联影医疗科技股份有限公司 | 三维图像拼接方法、装置、设备、系统和存储介质 |
| CN113520415A (zh) * | 2020-04-20 | 2021-10-22 | 上海联影医疗科技股份有限公司 | 一种x射线影像获取方法和系统 |
| CN113520416A (zh) * | 2020-04-21 | 2021-10-22 | 上海联影医疗科技股份有限公司 | 一种用于生成对象二维图像的方法和系统 |
| CN111914392B (zh) * | 2020-06-23 | 2024-05-24 | 上海联影医疗科技股份有限公司 | X射线成像设备和x射线图像的建模方法、装置 |
| CN111982939B (zh) * | 2020-07-27 | 2024-03-22 | 清华大学 | 可移动式多段直线光源ct成像系统及方法 |
| CN113393549A (zh) * | 2021-05-30 | 2021-09-14 | 南方医科大学 | 一种基于平板x射线源的ct图像重建方法 |
| CN113962922B (zh) * | 2021-07-14 | 2023-08-22 | 西安交通大学 | 一种基于平面阵列射线源的断层成像方法及系统 |
| CN114199907B (zh) * | 2021-12-06 | 2024-02-09 | 清华大学 | 多分辨率ct成像系统及方法 |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101567090A (zh) * | 2009-04-08 | 2009-10-28 | 西北工业大学 | 一种自适应的锥束ct三维图像快速重建方法 |
| CN103413338A (zh) * | 2013-05-29 | 2013-11-27 | 中国工程物理研究院流体物理研究所 | 一种基于广义变分最小化的少量投影ct图像重建方法 |
| CN103760180A (zh) * | 2014-01-27 | 2014-04-30 | 东南大学 | 基于x射线源阵列的实时在线工业ct检测系统 |
-
2016
- 2016-12-16 CN CN201611170283.XA patent/CN106651982B/zh active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101567090A (zh) * | 2009-04-08 | 2009-10-28 | 西北工业大学 | 一种自适应的锥束ct三维图像快速重建方法 |
| CN103413338A (zh) * | 2013-05-29 | 2013-11-27 | 中国工程物理研究院流体物理研究所 | 一种基于广义变分最小化的少量投影ct图像重建方法 |
| CN103760180A (zh) * | 2014-01-27 | 2014-04-30 | 东南大学 | 基于x射线源阵列的实时在线工业ct检测系统 |
Non-Patent Citations (1)
| Title |
|---|
| X线CT不完备投影数据统计重建研究;许琼;《中国博士学位论文全文数据库信息科技辑》;20120101;正文第1节 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN106651982A (zh) | 2017-05-10 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN106651982B (zh) | 一种基于阵列x射线源和探测器的ct图像重建方法 | |
| CN105361900B (zh) | 静态实时ct成像系统及其成像控制方法 | |
| CN108577876B (zh) | 一种多边形静止ct及其工作方法 | |
| CN101346102B (zh) | X射线ct装置 | |
| CN1106579C (zh) | 医学成象系统中多形态数据的互补重合 | |
| CN102793552B (zh) | Ct图像采集装置及ct扫描成像系统 | |
| CN104284627A (zh) | 锥束计算机断层扫描体积成像系统 | |
| CN102525536B (zh) | 堆叠平板x射线探测器部件及其制作方法 | |
| CN108811488A (zh) | 适应大视野要求的静态实时ct成像系统及其成像方法 | |
| EP2273257A1 (en) | Imaging system using a straight-line trajectory scan and method thereof | |
| CN105640582B (zh) | 一种深部组织x射线激发多光谱断层成像系统及方法 | |
| CN101057788A (zh) | 核素与荧光双模一体的小动物分子成像的数据获取系统 | |
| CN1424925A (zh) | 带有平板成像器的锥面束计算机x线断层扫描 | |
| CN103220978A (zh) | 具有单一探测器的pet-ct系统 | |
| CN1475000A (zh) | 使用圆加多弧轨迹进行圆锥射束体计算机断层摄影的系统和方法 | |
| CN102715914B (zh) | G型臂x光机的三维图像生成方法及装置与g型臂x光机 | |
| CN100376211C (zh) | 运动检查对象的计算机断层造影中设置图像的方法及设备 | |
| CN105899137A (zh) | 辐射检测器和使用辐射检测器的计算机断层扫描设备 | |
| CN103385732A (zh) | 一种静态ct扫描仪 | |
| CN109561865A (zh) | 一种通过时空重叠x射线来重建三维图像的方法和系统 | |
| CN102764137A (zh) | 一种静态ct扫描仪及其散射x光子校正方法 | |
| RU2013137832A (ru) | 4d компьютерная томография с контрастным усилением (кт) | |
| CN106706681B (zh) | 一种基于x射线源阵列成像的投影图像恢复方法 | |
| CN104237278A (zh) | 基于同步辐射的快速三维荧光ct系统及方法 | |
| CN110051378B (zh) | Ct成像方法、ct数据采集方法、系统、设备及存储介质 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| CB03 | Change of inventor or designer information |
Inventor after: Mou Xuanqin Inventor after: Qian Qinrong Inventor after: Cheng Haitao Inventor after: Wang Kai Inventor before: Mou Xuanqin Inventor before: Cheng Haitao Inventor before: Wang Kai |
|
| CB03 | Change of inventor or designer information | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |