CN112215953B - 图像重建方法、装置与电子设备 - Google Patents

图像重建方法、装置与电子设备 Download PDF

Info

Publication number
CN112215953B
CN112215953B CN202011244604.2A CN202011244604A CN112215953B CN 112215953 B CN112215953 B CN 112215953B CN 202011244604 A CN202011244604 A CN 202011244604A CN 112215953 B CN112215953 B CN 112215953B
Authority
CN
China
Prior art keywords
voxel
projection point
point
axis
projection
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
Application number
CN202011244604.2A
Other languages
English (en)
Other versions
CN112215953A (zh
Inventor
郑玉爽
许琼
魏存峰
刘双全
王燕芳
史戎坚
魏龙
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Institute of High Energy Physics of CAS
Original Assignee
Institute of High Energy Physics of CAS
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Institute of High Energy Physics of CAS filed Critical Institute of High Energy Physics of CAS
Priority to CN202011244604.2A priority Critical patent/CN112215953B/zh
Publication of CN112215953A publication Critical patent/CN112215953A/zh
Application granted granted Critical
Publication of CN112215953B publication Critical patent/CN112215953B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T17/00—Three dimensional [3D] modelling, e.g. data description of 3D objects
    • 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
    • G06T11/005—Specific pre-processing for tomographic reconstruction, e.g. calibration, source positioning, rebinning, scatter correction, retrospective gating
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/60—Analysis of geometric attributes
    • G06T7/62—Analysis of geometric attributes of area, perimeter, diameter or volume
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10072—Tomographic images
    • G06T2207/10081—Computed x-ray tomography [CT]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/20—Special algorithmic details
    • G06T2207/20068—Projection on vertical or horizontal image axis

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Geometry (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Computer Graphics (AREA)
  • Software Systems (AREA)
  • Apparatus For Radiation Diagnosis (AREA)

Abstract

本公开提供一种图像重建方法、装置与电子设备。方法包括：建立三维坐标系；确定在照射位置n下，被测物体的一个体素的体素值对探测器的一个探元的探测值的第一贡献因子，以及探测器的一个探元的探测值对被测物体的一个体素的体素值的第二贡献因子；根据被测物体的每个体素的假设体素值和每个照射位置对应的第一贡献因子更新每个探元的假设探测值，根据每个探元的假设探测值、实际探测值和每个照射位置对应的第二贡献因子更新每个体素的假设体素值，迭代计算直至达到预设条件时，将全部体素的假设体素值确定为被测物体的重建体素值；根据重建体素值重建被测物体。本公开实施例可以提高重建图像的清晰度。

Description

图像重建方法、装置与电子设备
技术领域
本公开涉及电子计算机断层扫描技术领域，具体而言，涉及一种图像重建方法、装置与电子设备。
背景技术
在传统的X射线CT(Computed Tomography，电子计算机断层扫描)成像技术中，射线源和探测器通常是按照预设的圆轨迹或螺旋轨迹等标准轨迹进行运动，图像重建算法也基于标准轨迹及其采样位置推导而来。然而，在实际操作中，存在一些特殊成像应用(如机械臂CT成像、板状物体CL(Computed Laminography，计算机层析成像技术)成像等)，在这些应用中，扫描轨迹可能不再是标准的圆轨迹或者螺旋轨迹。另外，即使在标准轨迹的CT成像中，受机械系统老化和运动精度的影响，实际的扫描轨迹和采样位置也常常与预设值之间存在一定的偏差。
对于采集轨迹不是标准轨迹的CT应用，相关技术的图像重建算法通常将扫描数据通过插值重排为标准轨迹数据，按照标准轨迹进行图像重建，然而插值操作会导致分辨率下降。对于运动精度导致的实际轨迹与预设轨迹出现微小偏差的情况，相关技术中的图像重建算法常常直接忽略而直接使用预设轨迹和参数进行图像重建，或者更新机械运动系统以提高运动精度，但是忽略偏差会导致成像分辨率不高，更新机械运动系统会造成改造成本过大。
因此，需要一种能够针对非标准轨迹或者运动精度不足的系统进行精准成像的图像重建方法。
需要说明的是，在上述背景技术部分公开的信息仅用于加强对本公开的背景的理解，因此可以包括不构成对本领域普通技术人员已知的现有技术的信息。
发明内容
本公开的目的在于提供一种图像重建方法、装置与电子设备，用于至少在一定程度上克服由于相关技术的限制和缺陷而导致的针对非标准轨迹运动的CT成像系统或者运动精度不足的CT成像系统图像重建体素分辨率不足的问题。
根据本公开实施例的第一方面，提供一种图像重建方法，包括：建立三维坐标系，将所述三维坐标系的原点设置在被测物体的几何中心，其中，所述被测物体的相对两侧分别设置有射线源和探测器；根据所述射线源和所述探测器在所述三维坐标系中的坐标，确定在每个照射位置下，被测物体的一个体素的体素值对所述探测器的一个探元的探测值的第一贡献因子，以及所述探测器的一个探元的探测值对所述被测物体的一个体素的体素值的第二贡献因子；根据所述被测物体的每个体素的假设体素值和每个照射位置对应的所述第一贡献因子更新每个所述照射位置下每个所述探元的假设探测值，根据每个所述照射位置下每个所述探元的假设探测值、实际探测值和所述第二贡献因子更新每个所述体素的假设体素值，迭代计算直至达到预设条件时，将全部所述体素的所述假设体素值确定为所述被测物体的重建体素值；根据所述重建体素值重建所述被测物体。
在本公开的一种示例性实施例中，所述确定在每个照射位置下，被测物体的一个体素的体素值对所述探测器的一个探元的探测值的第一贡献因子，以及所述探测器的一个探元的探测值对所述被测物体的一个体素的体素值的第二贡献因子包括：获取一个照射位置下的射线源坐标以及探测器中心点坐标，所述射线源坐标和所述探测器中心点坐标均包括第一轴坐标、第二轴坐标、第三轴坐标；获取所述射线源坐标中的第一轴坐标与所述探测器中心点坐标中的第一轴坐标的第一差值，以及所述射线源坐标中的第二轴坐标与所述探测器中心点坐标中的第二轴坐标的第二差值；在所述第一差值的绝对值小于等于所述第二差值的绝对值时，根据所述射线源坐标确定所述被测物体的第i体素在第一基准平面的第一投影区域，根据所述射线源坐标、所述探测器的探测器中心点坐标，确定所述探测器的第j探元在所述第一基准平面的第二投影区域，i≥1，j≥1，所述第一基准平面为所述三维坐标系中的第一轴和第三轴所在的平面；在所述第一差值的绝对值大于所述第二差值的绝对值时，根据所述射线源坐标确定所述被测物体的第i体素在第二基准平面的第三投影区域，根据所述射线源坐标、所述探测器的探测器中心点坐标，确定所述探测器的第j探元在所述第二基准平面的第四投影区域，所述第二基准平面为所述三维坐标系中的第二轴和第三轴所在的平面；确定所述第一投影区域和所述第二投影区域的第一重叠面积，或者所述第三投影区域和所述第四投影区域的第二重叠面积；根据所述第一重叠面积与所述第二投影区域的面积的比值、或者所述第二重叠面积与所述第四投影区域的面积的比值确定在所述照射位置下所述第i体素对所述第j探元的所述第一贡献因子，根据所述第一重叠面积与所述第一投影区域的面积的比值、或者所述第二重叠面积与所述第三投影区域的面积的比值确定在所述照射位置下所述第j探元对所述第i体素的所述第二贡献因子。
在本公开的一种示例性实施例中，所述根据所述射线源坐标确定所述被测物体的第i体素在第一基准平面的第一投影区域包括：获取所述被测物体的体素的总行数、总列数、体素边界长度；根据所述第i体素在所述被测物体中的行序号、列序号以及所述总行数、总列数、所述体素边界长度，确定所述第i体素在第一轴方向上的两个边界在第三基准平面上的第一边界投影点和第二边界投影点，以及所述第i体素在第三轴方向上的两个边界在第二基准平面上的第三边界投影点和第四边界投影点，所述第三基准平面为所述第一轴和所述第二轴所在的平面；确定所述射线源坐标在所述第三基准平面的第一源投影点和在所述第二基准平面的第二源投影点；根据所述第一源投影点和所述第一边界投影点的第一连线与所述第一轴的交点，或所述第一连线的延长线与所述第一轴的交点确定第一投影点，根据所述第一源投影点和所述第二边界投影点的第二连线与所述第一轴的交点，或所述第二连线的延长线与所述第一轴的交点确定第二投影点；根据所述第二源投影点和所述第三边界投影点的第三连线与所述第三轴的交点，或所述第三连线的延长线与所述第三轴的交点确定第三投影点，根据所述第二源投影点和所述第四边界投影点的第四连线与所述第三轴的交点，或所述第四连线的延长线与所述第三轴的交点确定第四投影点；根据所述第一投影点、所述第二投影点、所述第三投影点、所述第四投影点确定所述第一投影区域。
在本公开的一种示例性实施例中，所述根据所述射线源坐标、所述探测器的探测器中心点坐标，确定所述探测器的第j探元在所述第一基准平面的第二投影区域包括：获取所述探测器的探元长度和探元数量；根据所述射线源坐标、所述探测器中心点坐标、所述探元长度和所述探元数量确定第j探元在第一轴方向上的两个边界在所述第三基准平面上的第五边界投影点和第六边界投影点，以及所述第j探元在第三轴方向上的两个边界在所述第二基准平面上的第七边界投影点和第八边界投影点；根据所述第一源投影点和所述第五边界投影点的第五连线与所述第一轴的交点确定第五投影点，根据所述第一源投影点和所述第六边界投影点的第六连线与所述第一轴的交点确定第六投影点，根据所述第二源投影点和所述第七边界投影点的第七连线与所述第三轴的交点确定第七投影点，根据所述第二源投影点和所述第八边界投影点的第八连线与所述第三轴的交点确定第八投影点；根据所述第五投影点、所述第六投影点、所述第七投影点、所述第八投影点确定所述第二投影区域。
在本公开的一种示例性实施例中，所述根据所述射线源坐标确定所述被测物体的第i体素在第二基准平面的第三投影区域包括：获取所述被测物体的体素的总行数、总列数、体素边界长度；根据所述第i体素在所述被测物体中的行序号、列序号以及所述总行数、总列数、所述体素边界长度，确定所述第i体素在第二轴方向上的两个边界在第三基准平面上的第九边界投影点和第十边界投影点，以及所述第i体素在第三轴方向上的两个边界在第一基准平面上的第十一边界投影点和第十二边界投影点，所述第三基准平面为所述第一轴和所述第二轴所在的平面；确定所述射线源坐标在所述第三基准平面的第三源投影点和在所述第一基准平面的第四源投影点；根据所述第三源投影点和所述第九边界投影点的第九连线与所述第二轴的交点，或所述第九连线的延长线与所述第二轴的交点确定第九投影点，根据所述第三源投影点和所述第十边界投影点的第十连线与所述第二轴的交点，或所述第十连线的延长线与所述第二轴的交点确定第十投影点；根据所述第四源投影点和所述第十一边界投影点的第十一连线与所述第三轴的交点，或所述第十一连线的延长线与所述第三轴的交点确定第十一投影点，根据所述第四源投影点和所述第十二边界投影点的第十二连线与所述第三轴的交点，或所述第十二连线的延长线与所述第三轴的交点确定第十二投影点；根据所述第九投影点、所述第十投影点、所述第十一投影点、所述第十二投影点确定所述第三投影区域。
在本公开的一种示例性实施例中，所述根据所述射线源坐标、所述探测器的探测器中心点坐标，确定所述探测器的第j探元在所述第二基准平面的第四投影区域包括：获取所述探测器的探元长度和探元数量；根据所述射线源坐标、所述探测器中心点坐标、所述探元长度和所述探元数量确定第j探元在第二轴方向上的两个边界在所述第三基准平面上的第十三边界投影点和第十四边界投影点，以及所述第j探元在第三轴方向上的两个边界在所述第一基准平面上的第十五边界投影点和第十六边界投影点；根据所述第三源投影点和所述第十三边界投影点的第十三连线与所述第二轴的交点确定第十三投影点，根据所述第三源投影点和所述第十四边界投影点的第十四连线与所述第二轴的交点确定第十四投影点，根据所述第四源投影点和所述第十五边界投影点的第十五连线与所述第三轴的交点确定第十五投影点，根据所述第四源投影点和所述第十六边界投影点的第十六连线与所述第三轴的交点确定第十六投影点；根据所述第十三投影点、所述第十四投影点、所述第十五投影点、所述第十六投影点确定所述第四投影区域。
在本公开的一种示例性实施例中，所述根据所述被测物体的每个体素的假设体素值和每个照射位置对应的所述第一贡献因子更新每个所述照射位置下每个所述探元的假设探测值，根据每个所述照射位置下每个所述探元的假设探测值、实际探测值和所述第二贡献因子更新每个所述体素的假设体素值包括：设置所述被测物体中每个体素的假设体素值；将第i体素的所述假设体素值与照射位置n下所述第i体素对第j探元的所述第一贡献因子相乘得到所述照射位置n下所述第i体素对所述第j探元的体素贡献值；将所述照射位置n下所述被测物体中的全部体素对所述第j探元的所述体素贡献值相加，得到所述第j探元在所述照射位置n下的假设探测值；获取所述第j探元在所述照射位置n下的假设探测值与所述第j探元在所述照射位置n下的实际探测值的差值，将所述差值与所述照射位置n下所述第j探元对第i体素的第二贡献因子相乘以得到所述照射位置n下所述第j探元对所述第i体素的探元贡献值；将全部所述照射位置下全部所述探元对所述第i体素的所述探元贡献值相加得到第一值，将全部所述第一贡献因子构成的第一矩阵与全部第二贡献因子构成的第二矩阵相乘得到第二值，将所述第一值与所述第二值相除得到第三值；将所述第i体素的假设体素值与所述第三值的差值更新为所述第i体素的假设体素值。
根据本公开实施例的第二方面，提供一种图像重建装置，包括：坐标系建立模块，设置为建立三维坐标系，将所述三维坐标系的原点设置在被测物体的几何中心，其中，所述被测物体的相对两侧分别设置有射线源和探测器；贡献因子计算模块，设置为根据所述射线源和所述探测器在所述三维坐标系中的坐标，确定在每个照射位置下，被测物体的一个体素的体素值对所述探测器的一个探元的探测值的第一贡献因子，以及所述探测器的一个探元的探测值对所述被测物体的一个体素的体素值的第二贡献因子；重建体素值计算模块，设置为根据所述被测物体的每个体素的假设体素值和每个照射位置对应的所述第一贡献因子更新每个所述照射位置下每个所述探元的假设探测值，根据每个所述照射位置下每个所述探元的假设探测值、实际探测值和所述第二贡献因子更新每个所述体素的假设体素值，迭代计算直至达到预设条件时，将全部所述体素的所述假设体素值确定为所述被测物体的重建体素值；图像重建模块，设置为根据所述重建体素值重建所述被测物体。
根据本公开的第三方面，提供一种电子设备，包括：存储器；以及耦合到所述存储器的处理器，所述处理器被配置为基于存储在所述存储器中的指令，执行如上述任意一项所述的方法。
根据本公开的第四方面，提供一种计算机可读存储介质，其上存储有程序，该程序被处理器执行时实现如上述任意一项所述的图像重建方法。
本公开实施例提供的基于空间坐标的用于CT成像的三维图像重建方法，通过根据X射线源和探测器的精确空间坐标计算参数，确定在每个照射位置下被测物体的每个体素的体素值对每个探元的探测值的第一贡献因子和每个探元的探测值对每个体素的体素值的第二贡献因子，进而通过迭代算法根据被测物体的假设探测值确定探元的假设探测值，根据探元的假设探测值与实际探测值的差值反推被测物体的假设探测值，直至得到一个稳定的被测物体的假设探测值，即可得到被测物体的准确体素值。由于仅利用射线源和探测器的精确空间坐标计算被测物体的体素值，计算过程与射线源和探测器的连续运行轨迹无关，可以避免运行轨迹不规则或者运行轨迹出现误差对重建图像的清晰度造成的影响，极大地提高重建图像的清晰度。
应当理解的是，以上的一般描述和后文的细节描述仅是示例性和解释性的，并不能限制本公开。
附图说明
此处的附图被并入说明书中并构成本说明书的一部分，示出了符合本公开的实施例，并与说明书一起用于解释本公开的原理。显而易见地，下面描述中的附图仅仅是本公开的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1是本公开示例性实施例中图像重建方法的流程图。
图2是本公开实施例提供的三维坐标系的示意图。
图3是本公开一个实施例中步骤S2的子流程图。
图4是图3所示实施例中确定第一投影区域的流程图。
图5是图4所示实施例中确定第一投影点和第二投影点的示意图。
图6是图4所示实施例中确定第三投影点和第四投影点的示意图。
图7是图3所示实施例中确定第二投影区域的流程图。
图8是图7所示实施例中确定第五投影点和第六投影点的示意图。
图9是图7所示实施例中确定第七投影点和第八投影点的示意图。
图10是图3所示实施例中确定第一重叠面积的示意图。
图11是图3所示实施例中确定第三投影区域的流程图。
图12是图11所示实施例中确定第九投影点和第十投影点的示意图。
图13是图11所示实施例中确定第十一投影点和第十二投影点的示意图。
图14是图3所示实施例中确定第四投影区域的流程图。
图15是图14所示实施例中确定第十三投影点和第十四投影点的示意图。
图16是图14所示实施例中确定第十五投影点和第十六投影点的示意图。
图17是图3所示实施例中确定第二重叠面积的示意图。
图18是步骤S3中确定假设探测值的子流程图。
图19A是按照非标准轨迹生成的重建图像。
图19B是使用本公开实施例提供的方法生成的重建图像。
图20是本公开一个示例性实施例中一种图像重建装置的方框图。
图21是本公开一个示例性实施例中一种电子设备的方框图。
具体实施方式
现在将参考附图更全面地描述示例实施方式。然而，示例实施方式能够以多种形式实施，且不应被理解为限于在此阐述的范例；相反，提供这些实施方式使得本公开将更加全面和完整，并将示例实施方式的构思全面地传达给本领域的技术人员。所描述的特征、结构或特性可以以任何合适的方式结合在一个或更多实施方式中。在下面的描述中，提供许多具体细节从而给出对本公开的实施方式的充分理解。然而，本领域技术人员将意识到，可以实践本公开的技术方案而省略所述特定细节中的一个或更多，或者可以采用其它的方法、组元、装置、步骤等。在其它情况下，不详细示出或描述公知技术方案以避免喧宾夺主而使得本公开的各方面变得模糊。
此外，附图仅为本公开的示意性图解，图中相同的附图标记表示相同或类似的部分，因而将省略对它们的重复描述。附图中所示的一些方框图是功能实体，不一定必须与物理或逻辑上独立的实体相对应。可以采用软件形式来实现这些功能实体，或在一个或多个硬件模块或集成电路中实现这些功能实体，或在不同网络和/或处理器装置和/或微控制器装置中实现这些功能实体。
针对上述问题，本发明提出了一种基于空间坐标的高精度三维CT迭代重建算法，通过位置反馈系统(包括外置的光栅尺或机械系统的定位码表等)，获取射线源和探测器在每个照射位置下的精确空间坐标，然后进行迭代重建。在迭代过程中，利用三维距离驱动算法，由空间坐标计算迭代过程中的正、反投影，并对正、反投影过程，进行CUDA加速优化，提升重建速度。
下面结合附图对本公开示例实施方式进行详细说明。
图1是本公开示例性实施例中图像重建方法的流程图。参考图1，图像重建方法100可以包括：
步骤S1，建立三维坐标系，将所述三维坐标系的原点设置在被测物体的几何中心，其中，所述被测物体的相对两侧分别设置有射线源和探测器；
步骤S2，根据所述射线源和所述探测器在所述三维坐标系中的坐标，确定在每个照射位置下，被测物体的一个体素的体素值对所述探测器的一个探元的探测值的第一贡献因子，以及所述探测器的一个探元的探测值对所述被测物体的一个体素的体素值的第二贡献因子；
步骤S3，根据所述被测物体的每个体素的假设体素值和每个照射位置对应的所述第一贡献因子更新每个所述照射位置下每个所述探元的假设探测值，根据每个所述照射位置下每个所述探元的假设探测值、实际探测值和所述第二贡献因子更新每个所述体素的假设体素值，迭代计算直至达到预设条件时，将全部所述体素的所述假设体素值确定为所述被测物体的重建体素值；
步骤S4，根据所述重建体素值重建所述被测物体。
本公开实施例提供的基于空间坐标的CT三维图像重建方法，通过根据X射线源和探测器的精确空间坐标计算参数，确定在每个照射位置下被测物体的每个体素的体素值对每个探元的探测值的第一贡献因子和每个探元的探测值对每个体素的体素值的第二贡献因子，进而通过迭代算法根据被测物体的假设体素值确定探元的假设探测值，根据探元的假设探测值与实际探测值的差值反推被测物体的假设体素值，直至得到一个稳定的被测物体的假设体素值，即可得到被测物体的准确体素值。由于仅利用射线源和探测器的精确空间坐标计算被测物体的体素值，计算过程与射线源和探测器的连续运行轨迹无关，可以避免运行轨迹不规则或者运行轨迹出现误差对重建图像的清晰度造成的影响，极大地提高重建图像的清晰度。
下面，对图像重建方法100的各步骤进行详细说明。
在步骤S1，建立三维坐标系，将所述三维坐标系的原点设置在被测物体的几何中心，其中，所述被测物体的相对两侧分别设置有射线源和探测器。
图2是本公开实施例提供的三维坐标系的示意图。
参考图2，在本公开实施例中，选择被测物体的几何中心作为三维坐标系的原点，然后通过位置反馈系统获取照射位置n下射线源11和探测器12中心的精确坐标，将射线源11在三维坐标系中的坐标记录(xs,ys,zs)，将探测器12的中心坐标记录为(xd,yd,zd)。通过引入射线源和探测器中心的坐标(xs,ys,zs)和(xd,yd,zd)替换标准轨迹的几何参量进行运算，可以避免X射线源和探测器的运行轨迹不标准、有误差对计算结果的影响。在本公开的一个实施例中，将X轴称为第一轴，将Y轴称为第二轴，将Z轴称为第三轴，射线源11的x坐标、y坐标、z坐标分别被称为第一轴坐标、第二轴坐标、第三轴坐标。在本公开的其他实施例中，还可以将Y轴、Z轴称为第一轴，顺次修改其他轴的指代关系。以上坐标轴标注仅为示例，本领域技术人员可以根据实际情况自行修改。
在步骤S2，根据所述射线源和所述探测器在所述三维坐标系中的坐标，确定在每个照射位置下，被测物体的一个体素的体素值对所述探测器的一个探元的探测值的第一贡献因子，以及所述探测器的一个探元的探测值对所述被测物体的一个体素的体素值的第二贡献因子。
图3是本公开一个实施例中步骤S2的子流程图。
参考图3，在一个实施例中，步骤S2可以包括：
步骤S21，获取一个照射位置下的射线源坐标以及探测器中心点坐标，所述射线源坐标和所述探测器中心点坐标均包括第一轴坐标、第二轴坐标、第三轴坐标；
步骤S22，获取所述射线源坐标中的第一轴坐标与所述探测器中心点坐标中的第一轴坐标的第一差值，以及所述射线源坐标中的第二轴坐标与所述探测器中心点坐标中的第二轴坐标的第二差值；
步骤S23，在所述第一差值的绝对值小于等于所述第二差值的绝对值时，根据所述射线源坐标确定所述被测物体的第i体素在第一基准平面的第一投影区域，根据所述射线源坐标、所述探测器的探测器中心点坐标，确定所述探测器的第j探元在所述第一基准平面的第二投影区域，i≥1，j≥1，所述第一基准平面为所述三维坐标系中的第一轴和第三轴所在的平面；
步骤S24，在所述第一差值的绝对值大于所述第二差值的绝对值时，根据所述射线源坐标确定所述被测物体的第i体素在第二基准平面的第三投影区域，根据所述射线源坐标、所述探测器的探测器中心点坐标，确定所述探测器的第j探元在所述第二基准平面的第四投影区域，所述第二基准平面为所述三维坐标系中的第二轴和第三轴所在的平面；
步骤S25，确定所述第一投影区域和所述第二投影区域的第一重叠面积，或者所述第三投影区域和所述第四投影区域的第二重叠面积；
步骤S26，根据所述第一重叠面积与所述第二投影区域的面积的比值、或者所述第二重叠面积与所述第四投影区域的面积的比值确定在所述照射位置下所述第i体素对所述第j探元的所述第一贡献因子，根据所述第一重叠面积与所述第一投影区域的面积的比值、或者所述第二重叠面积与所述第三投影区域的面积的比值确定在所述照射位置下所述第j探元对所述第i体素的所述第二贡献因子。
在图3所示实施例中，根据射线源11在三维坐标系中的坐标特征，选择不同的基准平面进行投影计算。当本公开实施例中第一轴、第二轴、第三轴的设置为X轴、Y轴、Z轴时，如图2所示，射线源坐标中的第一轴坐标为xs，第二轴坐标为ys，第三轴坐标为zs，探测器中心坐标中的第一轴坐标为xd，第二轴坐标为yd，第三轴坐标为zd。当|xs-xd|≤|ys-yd|时，本公开实施例选择在X轴和Z轴所在的XOZ平面(第一基准平面)进行投影计算，当|xs-xd|＞|ys-yd|时，本公开实施例选择在Y轴和Z轴所在的YOZ平面(第二基准平面)进行投影计算。
当本公开其他实施例中的第一轴、第二轴、第三轴设置与上述设置不同时，基准平面的选择也可以与本公开以下附图中的实施例不同，此时各基准平面按照轴设置进行推算即可，本公开不以此为限。
图4是图3所示实施例中确定第一投影区域的流程图。
参考图3，步骤S23可以包括：
步骤S230，获取所述被测物体的体素的总行数、总列数、体素边界长度；
步骤S231，根据所述第i体素在所述被测物体中的行序号、列序号以及所述总行数、总列数、所述体素边界长度，确定所述第i体素在第一轴方向上的两个边界在第三基准平面上的第一边界投影点和第二边界投影点，以及所述第i体素在第三轴方向上的两个边界在第二基准平面上的第三边界投影点和第四边界投影点，所述第三基准平面为所述第一轴和所述第二轴所在的平面；
步骤S232，确定所述射线源坐标在所述第三基准平面的第一源投影点和在所述第二基准平面的第二源投影点；
步骤S233，根据所述第一源投影点和所述第一边界投影点的第一连线与所述第一轴的交点，或所述第一连线的延长线与所述第一轴的交点确定第一投影点，根据所述第一源投影点和所述第二边界投影点的第二连线与所述第一轴的交点，或所述第二连线的延长线与所述第一轴的交点确定第二投影点；
步骤S234，根据所述第二源投影点和所述第三边界投影点的第三连线与所述第三轴的交点，或所述第三连线的延长线与所述第三轴的交点确定第三投影点，根据所述第二源投影点和所述第四边界投影点的第四连线与所述第三轴的交点，或所述第四连线的延长线与所述第三轴的交点确定第四投影点；
步骤S235，根据所述第一投影点、所述第二投影点、所述第三投影点、所述第四投影点确定所述第一投影区域。
图5是图4所示实施例中确定第一投影点和第二投影点的示意图。
在本公开实施例中，可以使用三维距离驱动算法来计算迭代重建中的正投影(被测物体的体素值对探测器的探测值的作用)和反投影(探测器的探测值对被测物体的体素值的作用)。正、反投影的计算原理相同，下面以正投影的计算为例进行推导。
参考图5，假设第一轴为X轴，第二轴为Y轴，第三轴为Z轴。第一基准平面为XOZ平面，则如需计算第i体素在XOZ平面的第一投影区域，则首先需要计算第i体素在X轴上的投影点。由于被测物体13的几何中心为三维坐标系的原点，因此任何照射位置下，被测物体13的几何中心均为原点。
将被测物体13进行离散化，被测物体13的X、Y、Z方向体素数目分别为nx、ny、nz，总体素数量M＝nx*ny*nz，单个体素大小为px*py*pz。
为了计算第i体素的边界在X轴上的投影区域，可以首先将射线源11和被测物体13均映射到XOY平面上(即第三基准平面)，在同一平面进行计算，如图5所示。
记第nx体素在X轴方向上的左边界和右边界在第二基准平面(XOY平面)的投影点分别为Pnx和Pnx+1，当Pnx和Pnx+1映射到X轴上时，体素长度由px变为px′，px与px′比例恒定为f，d为当前行的体素到第一基准平面的垂直距离，由三角形相似可知：
px′＝px*f (2)
根据点Pnx+1与射线源11在第三基准平面上的第一源投影点A(第一源投影坐标为(xs，ys))的第一连线的延长线与X轴的交点，设确定点Pnx+1映射到X轴上的点P′nx+1的坐标为则有：
则第i体素在X轴方向上的左边界的第一投影坐标为：
将公式(4)带入公式(5)，可得：
同理可得第i体素在X轴方向上的右边界的第二投影坐标
根据上述公式，可得该行体素在X轴方向上的边界投影坐标数组为
图6是图4所示实施例中确定第三投影点和第四投影点的示意图。
参考图6，与图5对应，一个体素沿Z轴方向的边界长度pz与映射到Z轴的边界长度pz′的比例仍为f，即：
pz′＝pz*f (7)
Z轴方向体素边界的坐标计算方式与X轴方向计算体素边界的投影坐标的方式相同，射线源11在YOZ平面(第二基准平面)的投影记为第二源投影点B(具有第二源投影坐标(zs,ys))，根据公式(6)可类推第i体素在Z轴方向上的左边界在Z轴的投影点(即第三投影点)的坐标为：
同理可得第i体素在Z轴方向上的右边界在Z轴的投影点的坐标即第四投影点。
该行体素在Z轴上的边界坐标数组记为
如此，即可得到第i体素在XOZ平面的第一投影区域的角坐标为
图7是图3所示实施例中确定第二投影区域的流程图。
参考图7，步骤S23可以包括：
步骤S236，获取所述探测器的探元长度和探元数量；
步骤S237，根据所述射线源坐标、所述探测器中心点坐标、所述探元长度和所述探元数量确定第j探元在第一轴方向上的两个边界在所述第三基准平面上的第五边界投影点和第六边界投影点，以及所述第j探元在第三轴方向上的两个边界在所述第二基准平面上的第七边界投影点和第八边界投影点；
步骤S238，根据所述第一源投影点和所述第五边界投影点的第五连线与所述第一轴的交点确定第五投影点，根据所述第一源投影点和所述第六边界投影点的第六连线与所述第一轴的交点确定第六投影点，根据所述第二源投影点和所述第七边界投影点的第七连线与所述第三轴的交点确定第七投影点，根据所述第二源投影点和所述第八边界投影点的第八连线与所述第三轴的交点确定第八投影点；
步骤S239，根据所述第五投影点、所述第六投影点、所述第七投影点、所述第八投影点确定所述第二投影区域。
图8是图7所示实施例中确定第五投影点和第六投影点的示意图。
参考图8，设探测器12的行、列数分别为dnz和dnx，总探元数N＝dnz*dnx，单个探元的大小为dz*dx。
设探测器12的第j探元在X轴方向上的左边界为DXj，计算第j探元在XOZ平面的投影时，该边界在X轴方向的投影坐标为：
其中，是Dj与射线源11的连线与探测器12的中心与射线源11的连线的夹角，有：
sdd为射线源11在XOY(第三基准平面)上的投影点A(具有第一源投影坐标(xs,ys))到探测器12的中心点在XOY平面的投影点C(坐标为(xd,yd))的距离：
θ是CT的旋转角：
根据公式(9)～(13)，可以计算出探测器12中第j探元的左边界按X轴方向映射到公共平面上的坐标该行探元的边界的坐标数组为由此可得第j探元在X轴方向上的两个边界在X轴上的第五投影点和第六投影点的坐标分别为/>和/>
即当第j探元映射到XOZ平面(第一基准平面)时，X轴方向的探元长度变为：
dx′随j变化。
图9是图7所示实施例中确定第七投影点和第八投影点的示意图。
参考图9，探测器12的中心在Z轴方向的坐标为zd，射线源11在Z轴方向的坐标为zs，当探测器12的中点映射到第一基准平面，即XOZ平面时，其Z轴方向坐标的为zd′，则根据三角形相似原理可得：
第j探元在Z轴方向上的左边界DZj映射到XOZ平面上的点记为DZ′j。
探元沿Z轴方向的长度dz映射到XOZ平面后的大小dz′，
求得探测器中心在XOZ平面上的Z方向坐标zd′以及探元在XOZ平面上的长度dz′后，可求探测器每列体素边界映射到XOZ平面后的坐标：
将(16)、(18)带入(19)式，可得：
由此可得探元边界沿Z轴方向映射到公共平面的坐标数组
即可以得到第j探元的两个边界在Z轴方向上的第七投影点和第八投影点的Z坐标分别为/>和/>
图10是图3所示实施例中确定第一重叠面积的示意图。
参考图10，被测物体的第i体素在X方向的两个边界坐标分别为及/>在Z方向的两个边界坐标分别为/>及/>探测器的第j探元在X方向的两个边界坐标为/>及在Z方向的两个边界坐标分别为/>及/>即，第i体素在第一基准平面(XOZ平面)的第一投影区域的四个角坐标为/>第j探元在第一基准平面(XOZ平面)的第二投影区域的四个角坐标为/>
在第一基准平面上，第一投影区域与第二投影区域在X轴方向的重叠长度Δxij为：
第一投影区域与第二投影区域在Z轴方向重叠长度Δzij为：
/>
重叠区域Q1近似为矩形，则其面积为：
SQ1＝Δxij*Δzij (23)
对于正投影过程，记aij为第i体素对第j探元的第一贡献因子，有：
dx′、dz′分别为第二投影区域在X轴方向上的长度和在Z轴方向上的长度，可根据第二投影区域的端点坐标得到。
对于反投影过程，记bji为第j探元对第i体素的第二贡献因子，有：
px′、pz′分别为第一投影区域在X轴方向上的长度和在Z轴方向上的长度，可根据第一投影区域的端点坐标得到。
图11是本公开另一个实施例中步骤S24确定第三投影区域的子流程图。
参考图11，在一个实施例中，步骤S24中确定第三投影区域的过程可以包括：
步骤S240，获取所述被测物体的体素的总行数、总列数、体素边界长度；
步骤S241，根据所述第i体素在所述被测物体中的行序号、列序号以及所述总行数、总列数、所述体素边界长度，确定所述第i体素在第二轴方向上的两个边界在第三基准平面上的第九边界投影点和第十边界投影点，以及所述第i体素在第三轴方向上的两个边界在第一基准平面上的第十一边界投影点和第十二边界投影点，所述第三基准平面为所述第一轴和所述第二轴所在的平面；
步骤S242，确定所述射线源坐标在所述第三基准平面的第三源投影点和在所述第一基准平面的第四源投影点；
步骤S243，根据所述第三源投影点和所述第九边界投影点的第九连线与所述第二轴的交点，或所述第九连线的延长线与所述第二轴的交点确定第九投影点，根据所述第三源投影点和所述第十边界投影点的第十连线与所述第二轴的交点，或所述第十连线的延长线与所述第二轴的交点确定第十投影点；
步骤S244，根据所述第四源投影点和所述第十一边界投影点的第十一连线与所述第三轴的交点，或所述第十一连线的延长线与所述第三轴的交点确定第十一投影点，根据所述第四源投影点和所述第十二边界投影点的第十二连线与所述第三轴的交点，或所述第十二连线的延长线与所述第三轴的交点确定第十二投影点；
步骤S245，根据所述第九投影点、所述第十投影点、所述第十一投影点、所述第十二投影点确定所述第三投影区域
图12是图11所示实施例中第九投影点和第十投影点的示意图。
参考图12，计算YOZ平面(第二基准平面)上的第三投影区域的方法与计算XOZ平面(第一基准平面)上的第一投影区域的方法相同，差别在于，第九投影点P′ny和第十投影点P′ny+1是像素边界Pny和Pny+1在Y轴上的投影点，同样根据射线源11在XOY平面上的投影点A的坐标得到。
图13是图11所示实施例中第十一投影点和第十二投影点的示意图。
计算第三投影区域的第二步是计算体素边界在第三轴上的第十一投影点和第十二投影点，与计算第三投影点和第四投影点的方式相同，区别在于，图13所示实施例是根据射线源11在XOZ平面的投影点E的坐标计算像素边界Hnz和Hnz+1在Z轴上的第十一投影点H′nz和第十二投影点H′nz+1。
由于计算第九投影点、第十投影点、第十一投影点、第十二投影点的原理计算第一投影点、第二投影点、第三投影点、第四投影点的原理完全相同，可以参见图11所示流程，本公开于此不再赘述。
图14是本公开一个实施例中步骤S24确定第四投影区域的子流程图。
参考图14，步骤S24在确定第四投影区域时包括：
步骤S246，获取所述探测器的探元长度和探元数量；
步骤S247，根据所述射线源坐标、所述探测器中心点坐标、所述探元长度和所述探元数量确定第j探元在第二轴方向上的两个边界在所述第三基准平面上的第十三边界投影点和第十四边界投影点，以及所述第j探元在第三轴方向上的两个边界在所述第一基准平面上的第十五边界投影点和第十六边界投影点；
步骤S248，根据所述第三源投影点和所述第十三边界投影点的第十三连线与所述第二轴的交点确定第十三投影点，根据所述第三源投影点和所述第十四边界投影点的第十四连线与所述第二轴的交点确定第十四投影点，根据所述第四源投影点和所述第十五边界投影点的第十五连线与所述第三轴的交点确定第十五投影点，根据所述第四源投影点和所述第十六边界投影点的第十六连线与所述第三轴的交点确定第十六投影点；
步骤S249，根据所述第十三投影点、所述第十四投影点、所述第十五投影点、所述第十六投影点确定所述第四投影区域。
图15是图14所示实施例中第十三投影点和第十四投影点的示意图。
当投影平面为YOZ平面时，与计算第五投影点和第六投影点的方式相似，首先可以根据射线源11在XOY平面的投影点A与探测器中心点在XOY平面上的投影点C的坐标计算探元边界DYj和DYj+1在Y轴上的第十三投影点DYj′和第十四投影点DY′j+1。
图16是图14所示实施例中第十五投影点和第十六投影点的示意图。
与计算第七投影点和第八投影点的方式相似，可以利用射线源11在XOZ平面上的投影点E的坐标和探测器中心点在XOZ平面上的投影点F、F点与E点的连线与Z轴的交点F′的坐标计算第j探元的左边界FZj与E点的连线与Z轴的交点FZj′的坐标(即第十五投影点的Z轴坐标)，从而推算出第十六投影点FZj+1′的坐标。
由于计算第十三投影点、第十四投影点、第十五投影点、第十六投影点的原理计算第五投影点、第六投影点、第七投影点、第八投影点的原理完全相同，可以参见图14所示流程，本公开于此不再赘述。
图17是图3所示实施例中第二重叠面积的示意图。
与第一重叠区域的处理相似，第二重叠区域Q2的面积根据第三投影区域的四个角点坐标和第四投影区域的四个角点坐标/>在Z轴上的重叠长度Δzij和在Y轴上的重叠长度Δyij得出，本公开于此不再赘述。
在本公开的其他实施例中，确定第一贡献因子和第二贡献因子的方法还可以为其他计算方法，例如射线驱动算法等，本公开不以此为限。
在步骤S3，根据所述被测物体的每个体素的假设体素值和每个照射位置对应的所述第一贡献因子更新每个所述照射位置下每个所述探元的假设探测值，根据每个所述照射位置下每个所述探元的假设探测值、实际探测值和所述第二贡献因子更新每个所述体素的假设体素值，迭代计算直至达到预设条件时，将全部所述体素的所述假设体素值确定为所述被测物体的重建体素值。
图18是步骤S3中确定假设探测值的子流程图。
参考图18，步骤S3可以包括：
步骤S31，设置所述被测物体中每个体素的假设体素值；
步骤S32，将第i体素的所述假设体素值与照射位置n下所述第i体素对第j探元的所述第一贡献因子相乘得到所述照射位置n下所述第i体素对所述第j探元的体素贡献值；
步骤S33，将所述照射位置n下所述被测物体中的全部体素对所述第j探元的所述体素贡献值相加，得到所述第j探元在所述照射位置n下的假设探测值；
步骤S34，获取所述第j探元在所述照射位置n下的假设探测值与所述第j探元在所述照射位置n下的实际探测值的差值，将所述差值与所述照射位置n下所述第j探元对第i体素的第二贡献因子相乘以得到所述照射位置n下所述第j探元对所述第i体素的探元贡献值；
步骤S35，将全部所述照射位置下全部所述探元对所述第i体素的所述探元贡献值相加得到第一值，将全部所述第一贡献因子构成的第一矩阵与全部第二贡献因子构成的第二矩阵相乘得到第二值，将所述第一值与所述第二值相除得到第三值；
步骤S36，将所述第i体素的假设体素值与所述第三值的差值更新为所述第i体素的假设体素值。
在本公开实施例中，可以将每个体素的假设体素值的初始值均设置为非负值，在一些实施例中，可以均设置为零。
图18中，步骤S32～步骤S33可以称为正投影运算，步骤S34～步骤S35可以称为反投影运算，步骤S36为假设体素值更新。
将每个体素的假设体素值与其对第j探元的第一贡献因子相乘可以得到多个体素贡献值，对这些体素贡献值求和，即得到第j探元的假设探测值Pj，如下式所示：
μi为第i体素的假设体素值，Pj为第j探元的假设探测值。
将每个探元的假设探测值与其对第i体素的第二贡献因子相乘可得到多个假设探测值，对多个假设探测值累加，即得第i体素最新的假设体素值，如下式所示：
μi为第i体素的假设体素值，Pj为第j探元的假设探测值。
在图像重建过程中，可选用如下迭代公式：
式中，t为当前迭代次数(t＝1,2,3,…)，μt-1为上一次迭代后的重建图像的假设体素值矩阵，μ0为重建图像的初始假设体素值，可设为全0。下标i表示μt-1中的第i体素。A是由aij组成的矩阵，Aμt-1表示对μt-1使用前述正投影计算得到探测器的假设探测值矩阵Pt-1,[Aμt-1]j表示第j探元的假设探测值。为探测器的实际探测值矩阵，/>代表第j探元的实际探测值。
将[Aμt-1]j与做差，差值记为ΔPj。公式(28)等号右侧的分子部分可改写为：该部分代表反投影计算，与公式(27)式对应，计算得到的值可记为Δμi。分母部分为常数因子，表示对Δμi进行归一化。将第i体素的当前假设体素值/>减去归一化后的Δμi，得到第i体素的最新假设体素值。对全部体素进行计算以进行图像的更新，更新后的图像再继续上述操作，直至达到预设条件，此时的μ即为重建得到的被测物体的体素值。
在本公开实施例中，结束迭代过程的预设条件例如可以为ΔPj的值不再变化(可以设置在相邻两次迭代过程中计算的ΔPj的值之差小于预设值时判断ΔPj的值不再变化)、达到预设迭代次数、相邻两次迭代过程中计算的假设体素值不再变化(即Δμi不再变化，可以设置为判断相邻两次迭代过程中计算的假设体素值之差小于预设阈值)等等。预设条件的设置原则是在计算过程实现收敛时结束计算，因此本领域技术人员可以自行选择收敛判断标准来设置用于结束迭代过程的预设条件。
在本公开的其他实施例中，还可以选用其他一般迭代方法实现迭代过程，例如SIRT(Simultaneous Iterative Reconstruction Technique，联立迭代重建方法)等，在迭代中还可通过使用惩罚函数进一步提高重建质量。由于其他迭代方法和惩罚函数均为一般方法，本公开于此不再赘述。
在步骤S4，根据所述重建体素值重建所述被测物体。
图19A是按照非标准轨迹生成的重建图像。
图19B是使用本公开实施例提供的方法生成的重建图像。在图19B的实验中，使用1027*1027尺寸的探测器，绕物体一周，采集360个角度下的投影。
由图19A和图19B的对比可知，根据本公开实施例提供的方法进行图像重建，可以有效提高图像清晰度。
在图像重建过程中，每次迭代时的正投影、反投影计算耗时最长，为了提升重建速度，本公开实施例的方法可以使用CUDA(Compute Unified Device Architecture，统一计算设备架构)加速优化，使用CUDA进行正投影、反投影计算。CUDA是由NVIDIA(英伟达)开发的一种GPU(Graphics Processing Unit，图形处理器)编程模型，可实现并行计算。CUDA的并行计算依靠线程实现，线程的层次依次为grid、block、thread。在计算过程中，根据探测器尺寸划分线程，使每一个线程thread分别对应探测器中的一个探元，从而使得各探元上的计算依托线程而同时进行，从而提升计算速度。
本公开实施例首先使用位置反馈系统获取射线源和探测器中心的精确空间坐标，使用精确的空间坐标来表征CT系统的扫描轨迹，既能够精准刻画自由轨迹，进行非标准扫描轨迹的描述，也能反映由与机械系统老化及运动精度造成的实际轨迹与预设轨迹之间的偏差。通过使用迭代重建方法进行图像重建，在迭代过程中，利用三维空间坐标代替标准轨迹的几何参数进行正、反投影计算，避免了由于运动精度下降造成的实际轨迹与预设轨迹之间的偏差，从而实现了对运动精度下降导致的重建图像清晰度不足的校正，在不符合标准轨迹的应用场景中也能精确重建图像，无需对CT机械系统进行更换或校正。
对应于上述方法实施例，本公开还提供一种图像重建装置，可以用于执行上述方法实施例。
图20示意性示出本公开一个示例性实施例中一种图像重建装置的方框图。
参考图20，图像重建装置2000可以包括：
坐标系建立模块201，设置为建立三维坐标系，将所述三维坐标系的原点设置在被测物体的几何中心，其中，所述被测物体的相对两侧分别设置有射线源和探测器；
贡献因子计算模块202，设置为根据所述射线源和所述探测器在所述三维坐标系中的坐标，确定在每个照射位置下，被测物体的一个体素的体素值对所述探测器的一个探元的探测值的第一贡献因子，以及所述探测器的一个探元的探测值对所述被测物体的一个体素的体素值的第二贡献因子；
重建体素值计算模块203，设置为根据所述被测物体的每个体素的假设体素值和每个照射位置对应的所述第一贡献因子更新每个所述照射位置下每个所述探元的假设探测值，根据每个所述照射位置下每个所述探元的假设探测值、实际探测值和所述第二贡献因子更新每个所述体素的假设体素值，迭代计算直至达到预设条件时，将全部所述体素的所述假设体素值确定为所述被测物体的重建体素值；
图像重建模块204，设置为根据所述重建体素值重建所述被测物体。
装置2000可以执行如图1～图18所示实施例的方法，各模块的详细功能与图1～图18所示实施例对方法100中各步骤的描述对应，由于装置2000的各功能已在其对应的方法实施例中予以详细说明，本公开于此不再赘述。
应当注意，尽管在上文详细描述中提及了用于动作执行的设备的若干模块或者单元，但是这种划分并非强制性的。实际上，根据本公开的实施方式，上文描述的两个或更多模块或者单元的特征和功能可以在一个模块或者单元中具体化。反之，上文描述的一个模块或者单元的特征和功能可以进一步划分为由多个模块或者单元来具体化。
在本公开的示例性实施例中，还提供了一种能够实现上述方法的电子设备。
所属技术领域的技术人员能够理解，本发明的各个方面可以实现为系统、方法或程序产品。因此，本发明的各个方面可以具体实现为以下形式，即：完全的硬件实施方式、完全的软件实施方式(包括固件、微代码等)，或硬件和软件方面结合的实施方式，这里可以统称为“电路”、“模块”或“系统”。
下面参照图21来描述根据本发明的这种实施方式的电子设备2100。图21显示的电子设备2100仅仅是一个示例，不应对本发明实施例的功能和使用范围带来任何限制。
如图21所示，电子设备2100以通用计算设备的形式表现。电子设备2100的组件可以包括但不限于：上述至少一个处理单元2110、上述至少一个存储单元2120、连接不同系统组件(包括存储单元2120和处理单元2110)的总线2130。
其中，所述存储单元存储有程序代码，所述程序代码可以被所述处理单元2110执行，使得所述处理单元2110执行本说明书上述“示例性方法”部分中描述的根据本发明各种示例性实施方式的步骤。例如，所述处理单元2110可以执行如图1中所示的步骤。
存储单元2120可以包括易失性存储单元形式的可读介质，例如随机存取存储单元(RAM)21201和/或高速缓存存储单元21202，还可以进一步包括只读存储单元(ROM)21203。
存储单元2120还可以包括具有一组(至少一个)程序模块21205的程序/实用工具21204，这样的程序模块21205包括但不限于：操作系统、一个或者多个应用程序、其它程序模块以及程序数据，这些示例中的每一个或某种组合中可能包括网络环境的实现。
总线2130可以为表示几类总线结构中的一种或多种，包括存储单元总线或者存储单元控制器、外围总线、图形加速端口、处理单元或者使用多种总线结构中的任意总线结构的局域总线。
电子设备2100也可以与一个或多个外部设备2200(例如键盘、指向设备、蓝牙设备等)通信，还可与一个或者多个使得用户能与该电子设备2100交互的设备通信，和/或与使得该电子设备2100能与一个或多个其它计算设备进行通信的任何设备(例如路由器、调制解调器等等)通信。这种通信可以通过输入/输出(I/O)接口2150进行。并且，电子设备2100还可以通过网络适配器2160与一个或者多个网络(例如局域网(LAN)，广域网(WAN)和/或公共网络，例如因特网)通信。如图所示，网络适配器2160通过总线2130与电子设备2100的其它模块通信。应当明白，尽管图中未示出，可以结合电子设备2100使用其它硬件和/或软件模块，包括但不限于：微代码、设备驱动器、冗余处理单元、外部磁盘驱动阵列、RAID系统、磁带驱动器以及数据备份存储系统等。
通过以上的实施方式的描述，本领域的技术人员易于理解，这里描述的示例实施方式可以通过软件实现，也可以通过软件结合必要的硬件的方式来实现。因此，根据本公开实施方式的技术方案可以以软件产品的形式体现出来，该软件产品可以存储在一个非易失性存储介质(可以是CD-ROM，U盘，移动硬盘等)中或网络上，包括若干指令以使得一台计算设备(可以是个人计算机、服务器、终端装置、或者网络设备等)执行根据本公开实施方式的方法。
在本公开的示例性实施例中，还提供了一种计算机可读存储介质，其上存储有能够实现本说明书上述方法的程序产品。在一些可能的实施方式中，本发明的各个方面还可以实现为一种程序产品的形式，其包括程序代码，当所述程序产品在终端设备上运行时，所述程序代码用于使所述终端设备执行本说明书上述方法部分中描述的根据本发明各种示例性实施方式的步骤。
根据本发明的实施方式的用于实现上述方法的程序产品可以采用便携式紧凑盘只读存储器(CD-ROM)并包括程序代码，并可以在终端设备，例如个人电脑上运行。然而，本发明的程序产品不限于此，在本文件中，可读存储介质可以是任何包含或存储程序的有形介质，该程序可以被指令执行系统、装置或者器件使用或者与其结合使用。
所述程序产品可以采用一个或多个可读介质的任意组合。可读介质可以是可读信号介质或者可读存储介质。可读存储介质例如可以为但不限于电、磁、光、电磁、红外线、或半导体的系统、装置或器件，或者任意以上的组合。可读存储介质的更具体的例子(非穷举的列表)包括：具有一个或多个导线的电连接、便携式盘、硬盘、随机存取存储器(RAM)、只读存储器(ROM)、可擦式可编程只读存储器(EPROM或闪存)、光纤、便携式紧凑盘只读存储器(CD-ROM)、光存储器件、磁存储器件、或者上述的任意合适的组合。
计算机可读信号介质可以包括在基带中或者作为载波一部分传播的数据信号，其中承载了可读程序代码。这种传播的数据信号可以采用多种形式，包括但不限于电磁信号、光信号或上述的任意合适的组合。可读信号介质还可以是可读存储介质以外的任何可读介质，该可读介质可以发送、传播或者传输用于由指令执行系统、装置或者器件使用或者与其结合使用的程序。
可读介质上包含的程序代码可以用任何适当的介质传输，包括但不限于无线、有线、光缆、RF等等，或者上述的任意合适的组合。
可以以一种或多种程序设计语言的任意组合来编写用于执行本发明操作的程序代码，所述程序设计语言包括面向对象的程序设计语言—诸如Java、C++等，还包括常规的过程式程序设计语言—诸如“C”语言或类似的程序设计语言。程序代码可以完全地在用户计算设备上执行、部分地在用户设备上执行、作为一个独立的软件包执行、部分在用户计算设备上部分在远程计算设备上执行、或者完全在远程计算设备或服务器上执行。在涉及远程计算设备的情形中，远程计算设备可以通过任意种类的网络，包括局域网(LAN)或广域网(WAN)，连接到用户计算设备，或者，可以连接到外部计算设备(例如利用因特网服务提供商来通过因特网连接)。
此外，上述附图仅是根据本发明示例性实施例的方法所包括的处理的示意性说明，而不是限制目的。易于理解，上述附图所示的处理并不表明或限制这些处理的时间顺序。另外，也易于理解，这些处理可以是例如在多个模块中同步或异步执行的。
本领域技术人员在考虑说明书及实践这里公开的发明后，将容易想到本公开的其它实施方案。本申请旨在涵盖本公开的任何变型、用途或者适应性变化，这些变型、用途或者适应性变化遵循本公开的一般性原理并包括本公开未公开的本技术领域中的公知常识或惯用技术手段。说明书和实施例仅被视为示例性的，本公开的真正范围和构思由权利要求指出。

Claims (10)

1.一种图像重建方法，其特征在于，包括：
建立三维坐标系，将所述三维坐标系的原点设置在被测物体的几何中心，其中，所述被测物体的相对两侧分别设置有射线源和探测器；
根据所述射线源和所述探测器在所述三维坐标系中的坐标，确定在每个照射位置下，被测物体的一个体素的体素值对所述探测器的一个探元的探测值的第一贡献因子，以及所述探测器的一个探元的探测值对所述被测物体的一个体素的体素值的第二贡献因子；
根据所述被测物体的每个体素的假设体素值和每个照射位置对应的所述第一贡献因子更新每个所述照射位置下每个所述探元的假设探测值，根据每个所述照射位置下每个所述探元的假设探测值、实际探测值和所述第二贡献因子更新每个所述体素的假设体素值，迭代计算直至达到预设条件时，将全部所述体素的所述假设体素值确定为所述被测物体的重建体素值；
根据所述重建体素值重建所述被测物体。
2.如权利要求1所述的图像重建方法，其特征在于，所述确定在每个照射位置下，被测物体的一个体素的体素值对所述探测器的一个探元的探测值的第一贡献因子，以及所述探测器的一个探元的探测值对所述被测物体的一个体素的体素值的第二贡献因子包括：
获取一个照射位置下的射线源坐标以及探测器中心点坐标，所述射线源坐标和所述探测器中心点坐标均包括第一轴坐标、第二轴坐标、第三轴坐标；
获取所述射线源坐标中的第一轴坐标与所述探测器中心点坐标中的第一轴坐标的第一差值，以及所述射线源坐标中的第二轴坐标与所述探测器中心点坐标中的第二轴坐标的第二差值；
在所述第一差值的绝对值小于等于所述第二差值的绝对值时，根据所述射线源坐标确定所述被测物体的第i体素在第一基准平面的第一投影区域，根据所述射线源坐标、所述探测器的探测器中心点坐标，确定所述探测器的第j探元在所述第一基准平面的第二投影区域，i≥1，j≥1，所述第一基准平面为所述三维坐标系中的第一轴和第三轴所在的平面；
在所述第一差值的绝对值大于所述第二差值的绝对值时，根据所述射线源坐标确定所述被测物体的第i体素在第二基准平面的第三投影区域，根据所述射线源坐标、所述探测器的探测器中心点坐标，确定所述探测器的第j探元在所述第二基准平面的第四投影区域，所述第二基准平面为所述三维坐标系中的第二轴和第三轴所在的平面；
确定所述第一投影区域和所述第二投影区域的第一重叠面积，或者所述第三投影区域和所述第四投影区域的第二重叠面积；
根据所述第一重叠面积与所述第二投影区域的面积的比值或者所述第二重叠面积与所述第四投影区域的面积的比值确定在所述照射位置下所述第i体素对所述第j探元的所述第一贡献因子，根据所述第一重叠面积与所述第一投影区域的面积的比值或者所述第二重叠面积与所述第三投影区域的面积的比值确定在所述照射位置下所述第j探元对所述第i体素的所述第二贡献因子。
3.如权利要求2所述的图像重建方法，其特征在于，所述根据所述射线源坐标确定所述被测物体的第i体素在第一基准平面的第一投影区域包括：
获取所述被测物体的体素的总行数、总列数、体素边界长度；
根据所述第i体素在所述被测物体中的行序号、列序号以及所述总行数、总列数、所述体素边界长度，确定所述第i体素在第一轴方向上的两个边界在第三基准平面上的第一边界投影点和第二边界投影点，以及所述第i体素在第三轴方向上的两个边界在第二基准平面上的第三边界投影点和第四边界投影点，所述第三基准平面为所述第一轴和所述第二轴所在的平面；
确定所述射线源坐标在所述第三基准平面的第一源投影点和在所述第二基准平面的第二源投影点；
根据所述第一源投影点和所述第一边界投影点的第一连线与所述第一轴的交点，或所述第一连线的延长线与所述第一轴的交点确定第一投影点，根据所述第一源投影点和所述第二边界投影点的第二连线与所述第一轴的交点，或所述第二连线的延长线与所述第一轴的交点确定第二投影点；
根据所述第二源投影点和所述第三边界投影点的第三连线与所述第三轴的交点，或所述第三连线的延长线与所述第三轴的交点确定第三投影点，根据所述第二源投影点和所述第四边界投影点的第四连线与所述第三轴的交点，或所述第四连线的延长线与所述第三轴的交点确定第四投影点；
根据所述第一投影点、所述第二投影点、所述第三投影点、所述第四投影点确定所述第一投影区域。
4.如权利要求3所述的图像重建方法，其特征在于，所述根据所述射线源坐标、所述探测器的探测器中心点坐标，确定所述探测器的第j探元在所述第一基准平面的第二投影区域包括：
获取所述探测器的探元长度和探元数量；
根据所述射线源坐标、所述探测器中心点坐标、所述探元长度和所述探元数量确定第j探元在第一轴方向上的两个边界在所述第三基准平面上的第五边界投影点和第六边界投影点，以及所述第j探元在第三轴方向上的两个边界在所述第二基准平面上的第七边界投影点和第八边界投影点；
根据所述第一源投影点和所述第五边界投影点的第五连线与所述第一轴的交点确定第五投影点，根据所述第一源投影点和所述第六边界投影点的第六连线与所述第一轴的交点确定第六投影点，根据所述第二源投影点和所述第七边界投影点的第七连线与所述第三轴的交点确定第七投影点，根据所述第二源投影点和所述第八边界投影点的第八连线与所述第三轴的交点确定第八投影点；
根据所述第五投影点、所述第六投影点、所述第七投影点、所述第八投影点确定所述第二投影区域。
5.如权利要求2所述的图像重建方法，其特征在于，所述根据所述射线源坐标确定所述被测物体的第i体素在第二基准平面的第三投影区域包括：
获取所述被测物体的体素的总行数、总列数、体素边界长度；
根据所述第i体素在所述被测物体中的行序号、列序号以及所述总行数、总列数、所述体素边界长度，确定所述第i体素在第二轴方向上的两个边界在第三基准平面上的第九边界投影点和第十边界投影点，以及所述第i体素在第三轴方向上的两个边界在第一基准平面上的第十一边界投影点和第十二边界投影点，所述第三基准平面为所述第一轴和所述第二轴所在的平面；
确定所述射线源坐标在所述第三基准平面的第三源投影点和在所述第一基准平面的第四源投影点；
根据所述第三源投影点和所述第九边界投影点的第九连线与所述第二轴的交点，或所述第九连线的延长线与所述第二轴的交点确定第九投影点，根据所述第三源投影点和所述第十边界投影点的第十连线与所述第二轴的交点，或所述第十连线的延长线与所述第二轴的交点确定第十投影点；
根据所述第四源投影点和所述第十一边界投影点的第十一连线与所述第三轴的交点，或所述第十一连线的延长线与所述第三轴的交点确定第十一投影点，根据所述第四源投影点和所述第十二边界投影点的第十二连线与所述第三轴的交点，或所述第十二连线的延长线与所述第三轴的交点确定第十二投影点；
根据所述第九投影点、所述第十投影点、所述第十一投影点、所述第十二投影点确定所述第三投影区域。
6.如权利要求5所述的图像重建方法，其特征在于，所述根据所述射线源坐标、所述探测器的探测器中心点坐标，确定所述探测器的第j探元在所述第二基准平面的第四投影区域包括：
获取所述探测器的探元长度和探元数量；
根据所述射线源坐标、所述探测器中心点坐标、所述探元长度和所述探元数量确定第j探元在第二轴方向上的两个边界在所述第三基准平面上的第十三边界投影点和第十四边界投影点，以及所述第j探元在第三轴方向上的两个边界在所述第一基准平面上的第十五边界投影点和第十六边界投影点；
根据所述第三源投影点和所述第十三边界投影点的第十三连线与所述第二轴的交点确定第十三投影点，根据所述第三源投影点和所述第十四边界投影点的第十四连线与所述第二轴的交点确定第十四投影点，根据所述第四源投影点和所述第十五边界投影点的第十五连线与所述第三轴的交点确定第十五投影点，根据所述第四源投影点和所述第十六边界投影点的第十六连线与所述第三轴的交点确定第十六投影点；
根据所述第十三投影点、所述第十四投影点、所述第十五投影点、所述第十六投影点确定所述第四投影区域。
7.如权利要求1所述的图像重建方法，其特征在于，所述根据所述被测物体的每个体素的假设体素值和每个照射位置对应的所述第一贡献因子更新每个所述照射位置下每个所述探元的假设探测值，根据每个所述照射位置下每个所述探元的假设探测值、实际探测值和所述第二贡献因子更新每个所述体素的假设体素值包括：
设置所述被测物体中每个体素的假设体素值；
将第i体素的所述假设体素值与照射位置n下所述第i体素对第j探元的所述第一贡献因子相乘得到所述照射位置n下所述第i体素对所述第j探元的体素贡献值；
将所述照射位置n下所述被测物体中的全部体素对所述第j探元的所述体素贡献值相加，得到所述第j探元在所述照射位置n下的假设探测值；
获取所述第j探元在所述照射位置n下的假设探测值与所述第j探元在所述照射位置n下的实际探测值的差值，将所述差值与所述照射位置n下所述第j探元对第i体素的第二贡献因子相乘以得到所述照射位置n下所述第j探元对所述第i体素的探元贡献值；
将全部所述照射位置下全部所述探元对所述第i体素的所述探元贡献值相加得到第一值，将全部所述第一贡献因子构成的第一矩阵与全部第二贡献因子构成的第二矩阵相乘得到第二值，将所述第一值与所述第二值相除得到第三值；
将所述第i体素的假设体素值与所述第三值的差值更新为所述第i体素的假设体素值。
8.一种图像重建装置，其特征在于，包括：
坐标系建立模块，设置为建立三维坐标系，将所述三维坐标系的原点设置在被测物体的几何中心，其中，所述被测物体的相对两侧分别设置有射线源和探测器；
贡献因子计算模块，设置为根据所述射线源和所述探测器在所述三维坐标系中的坐标，确定在每个照射位置下，被测物体的一个体素的体素值对所述探测器的一个探元的探测值的第一贡献因子，以及所述探测器的一个探元的探测值对所述被测物体的一个体素的体素值的第二贡献因子；
重建体素值计算模块，设置为根据所述被测物体的每个体素的假设体素值和每个照射位置对应的所述第一贡献因子更新每个所述照射位置下每个所述探元的假设探测值，根据每个所述照射位置下每个所述探元的假设探测值、实际探测值和所述第二贡献因子更新每个所述体素的假设体素值，迭代计算直至达到预设条件时，将全部所述体素的所述假设体素值确定为所述被测物体的重建体素值；
图像重建模块，设置为根据所述重建体素值重建所述被测物体。
9.一种电子设备，其特征在于，包括：
存储器；以及
耦合到所述存储器的处理器，所述处理器被配置为基于存储在所述存储器中的指令，执行如权利要求1-7任一项所述的图像重建方法。
10.一种计算机可读存储介质，其上存储有程序，该程序被处理器执行时实现如权利要求1-7任一项所述的图像重建方法。
CN202011244604.2A 2020-11-10 2020-11-10 图像重建方法、装置与电子设备 Active CN112215953B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011244604.2A CN112215953B (zh) 2020-11-10 2020-11-10 图像重建方法、装置与电子设备

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011244604.2A CN112215953B (zh) 2020-11-10 2020-11-10 图像重建方法、装置与电子设备

Publications (2)

Publication Number Publication Date
CN112215953A CN112215953A (zh) 2021-01-12
CN112215953B true CN112215953B (zh) 2023-11-17

Family

ID=74056685

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011244604.2A Active CN112215953B (zh) 2020-11-10 2020-11-10 图像重建方法、装置与电子设备

Country Status (1)

Country Link
CN (1) CN112215953B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114278281B (zh) * 2021-12-24 2023-11-21 北京西华益昌技术开发有限责任公司 测量装置的测量分辨率优化方法、装置、设备及存储介质

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH10192271A (ja) * 1997-01-10 1998-07-28 Toshiba Corp Ｘ線ｃｔ装置及び画像処理装置
CN1584931A (zh) * 2003-08-19 2005-02-23 安捷伦科技有限公司 对被检测对象多个深度层进行并行图像重建的系统和方法
CN1711968A (zh) * 2005-05-26 2005-12-28 西安理工大学 Ct图像的快速渐进式直接体绘制三维重建方法
EP2705332A1 (de) * 2011-05-06 2014-03-12 Werth Messtechnik GmbH BESTIMMUNG VON MAßEN EINES WERKSTÜCKS
CN103654833A (zh) * 2013-11-19 2014-03-26 中国科学院过程工程研究所 Ct探测器偏转角的确定方法和装置
CN106204679A (zh) * 2016-07-18 2016-12-07 上海交通大学 基于可分离足迹函数技术的投影方法、装置及系统
CN107714072A (zh) * 2017-11-20 2018-02-23 中国科学院高能物理研究所 缺失数据的补偿方法、计算机断层扫描成像方法及系统
CN109658465A (zh) * 2018-12-07 2019-04-19 广州华端科技有限公司 图像重建过程中的数据处理、图像重建方法和装置

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US10255696B2 (en) * 2015-12-11 2019-04-09 Shanghai United Imaging Healthcare Co., Ltd. System and method for image reconstruction
US10445905B2 (en) * 2016-03-31 2019-10-15 Shanghai United Imaging Healthcare Co., Ltd. System and method for image reconstruction

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH10192271A (ja) * 1997-01-10 1998-07-28 Toshiba Corp Ｘ線ｃｔ装置及び画像処理装置
CN1584931A (zh) * 2003-08-19 2005-02-23 安捷伦科技有限公司 对被检测对象多个深度层进行并行图像重建的系统和方法
CN1711968A (zh) * 2005-05-26 2005-12-28 西安理工大学 Ct图像的快速渐进式直接体绘制三维重建方法
EP2705332A1 (de) * 2011-05-06 2014-03-12 Werth Messtechnik GmbH BESTIMMUNG VON MAßEN EINES WERKSTÜCKS
CN103654833A (zh) * 2013-11-19 2014-03-26 中国科学院过程工程研究所 Ct探测器偏转角的确定方法和装置
CN106204679A (zh) * 2016-07-18 2016-12-07 上海交通大学 基于可分离足迹函数技术的投影方法、装置及系统
CN107714072A (zh) * 2017-11-20 2018-02-23 中国科学院高能物理研究所 缺失数据的补偿方法、计算机断层扫描成像方法及系统
CN109658465A (zh) * 2018-12-07 2019-04-19 广州华端科技有限公司 图像重建过程中的数据处理、图像重建方法和装置

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
An Accelerated Iterative Cone Beam Computed Tomography Image Reconstruction Approach;Shimaa Abdulsalam Khazal et al.;《Nahrain Journal for Engineering Sciences》;307-314 *
高分辨锥束CT并行重建算法在基于NVDIA GPU 显卡计算平台上的实现;郑海亮 等;《CT理论与应用研究》;805-814 *

Also Published As

Publication number Publication date
CN112215953A (zh) 2021-01-12

Similar Documents

Publication Publication Date Title
CN109521403B (zh) 多线激光雷达的参数标定方法及装置、设备及可读介质
CN110927740B (zh) 一种移动机器人定位方法
US10627520B2 (en) Method and apparatus for constructing reflectance map
WO2021128787A1 (zh) 一种定位的方法及装置
CN109171793B (zh) 一种角度检测和校正方法、装置、设备和介质
CN110780285A (zh) 激光雷达与组合惯导的位姿标定方法、系统及介质
US20070122020A1 (en) Method and device for geometry analysis and calibration of volumetric imaging systems
US8224121B2 (en) System and method for assembling substantially distortion-free images
US9549712B2 (en) Method and apparatus for correcting focus of CT device
CN102028490A (zh) 从在遍历轨迹时拍摄的x射线投影确定图像的方法和设备
US20140093152A1 (en) Methods and devices for locating object in ct imaging
US20090039285A1 (en) Method and device for controlling and monitoring a position of a holding element
WO2023082306A1 (zh) 图像处理方法、装置、电子设备以及计算机可读存储介质
US20150103968A1 (en) Computed Tomography (CT) Image Reconstruction Method
CN112215953B (zh) 图像重建方法、装置与电子设备
CN108520543B (zh) 一种对相对精度地图进行优化的方法、设备及存储介质
CN114049401A (zh) 双目相机标定方法、装置、设备及介质
CN108733211A (zh) 追踪系统、其操作方法、控制器、及电脑可读取记录媒体
CN111220100B (zh) 基于激光束的测量方法、装置、系统、控制设备及介质
CN113034562A (zh) 用于优化深度信息的方法和装置
CN111462321A (zh) 点云地图处理方法、处理装置、电子装置和车辆
Griguletskii et al. TomoSLAM: factor graph optimization for rotation angle refinement in microtomography
CN112589806B (zh) 机器人位姿信息确定方法、装置、设备及存储介质
CN113963058B (zh) 预定轨道ct在线标定方法、装置、电子设备及存储介质
CN113534193B (zh) 确定目标反射点的方法、装置、电子设备及存储介质

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
GR01 Patent grant
GR01 Patent grant