WO2018035905A1 - 一种锥束ct三维重建方法及系统 - Google Patents

一种锥束ct三维重建方法及系统 Download PDF

Info

Publication number
WO2018035905A1
WO2018035905A1 PCT/CN2016/099565 CN2016099565W WO2018035905A1 WO 2018035905 A1 WO2018035905 A1 WO 2018035905A1 CN 2016099565 W CN2016099565 W CN 2016099565W WO 2018035905 A1 WO2018035905 A1 WO 2018035905A1
Authority
WO
WIPO (PCT)
Prior art keywords
dimensional
frame image
voxel
axis
projection
Prior art date
Application number
PCT/CN2016/099565
Other languages
English (en)
French (fr)
Inventor
陈垦
熊璟
王澄
秦文健
谢耀钦
Original Assignee
深圳先进技术研究院
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 深圳先进技术研究院 filed Critical 深圳先进技术研究院
Publication of WO2018035905A1 publication Critical patent/WO2018035905A1/zh

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T11/00—2D [Two Dimensional] image generation
    • 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
    • G06T2210/00—Indexing scheme for image generation or computer graphics
    • G06T2210/41—Medical

Definitions

  • the invention relates to a reconstruction algorithm technology in the field of X-ray CT technology, in particular to a cone beam CT three-dimensional reconstruction method and system.
  • Patent CN104899903A, CN102609978A, CN102279970A provide three-dimensional reconstruction based on FDK algorithm.
  • a filtered back projection accumulation algorithm although the calculation speed is fast and the system resource requirements are low, the reconstructed image quality is poor, there are artifacts, which affect the visualization effect and doctor diagnosis.
  • the ART algorithm needs to calculate the projection matrix.
  • the algorithm proposed by Siddon et al. is usually used.
  • the algorithm is a point-by-point calculation method, that is, the position of the upper and lower points of the ray is determined according to the relative positional relationship between a certain point on the ray and the voxel to be reconstructed. Larger, the computational efficiency is lower, especially in the case where the volume data to be reconstructed is three-dimensional, the defects are more obvious, which limits the clinical application of the ART algorithm.
  • there are some improved algorithms most of them use approximate assumptions. For example, only to judge whether the ray crosses the voxel, the traversal is recorded as 1, otherwise the 01 approximation, which is recorded as 0, affects the accuracy of the projection matrix.
  • the filtered back projection algorithm represented by the FDK algorithm has a fast calculation speed and low requirements on system resources, but the reconstructed image quality is poor, and there are artifacts, which affect the visualization effect and doctor diagnosis.
  • ART The algorithm needs to calculate the projection matrix.
  • the algorithm proposed by Siddon et al. is usually adopted.
  • the algorithm is a point-by-point calculation method, that is, the position of the upper and lower points of the ray is determined according to the relative positional relationship between a certain point on the ray and the voxel to be reconstructed. Large, the computational efficiency is low, especially in the case where the volume data to be reconstructed is three-dimensional, the defects are more obvious, which limits the clinical application of the ART algorithm.
  • there are some improved algorithms most of them use approximate assumptions. For example, only to judge whether the ray crosses the voxel, the traversal is recorded as 1, otherwise the 01 approximation, which is recorded as 0, affects the accuracy of the projection matrix.
  • a cone beam CT three-dimensional reconstruction method includes the following steps:
  • Step S110 collecting projection data
  • 2 + ⁇ R(x), where A is a projection matrix, x is a three-dimensional volume data to be reconstructed, b is an actual acquired projection image, b [ b 0 , b 1 ... b i ... b N ] T , R is a constraint term, ⁇ is an adjustment coefficient, and for each frame image b i , the three-dimensional position coordinates of the ray source corresponding to the frame image are calculated;
  • Step S140 Steps S120 to S130 are repeated until the termination condition is reached.
  • the projection data is stored in an N*H*W three-dimensional array with unsigned short type data, where N is the number of projection data frames, H is the projection data height, and W is the projection data width.
  • step S130 is performed after step S120 is completed, the following steps are further included:
  • Step S121 Calculate a thread for each pixel point of the left half axis of each frame image b i X axis, calculate the coordinates of the pixel point in space, and find the coordinates of the source corresponding to the frame image corresponding to the frame image.
  • Step S122 determining whether the line intersects with each voxel located in the left half axis of the X axis in the three-dimensional volume data to be reconstructed;
  • the three-dimensional position (x, y, z) of the voxel in the three-dimensional volume data is recorded, the length value is set to 0, and the corresponding length of the voxel with respect to the axisymmetric position of the X axis is also set to zero.
  • Step S123 When all the threads are executed, the projection matrix A ij corresponding to all the pixel points b ij of the current frame image is obtained, and the corresponding component is calculated. M is the number of pixels in the image.
  • the present invention also provides a cone beam CT three-dimensional reconstruction system, comprising:
  • a projection data acquisition unit for collecting projection data
  • 2 + ⁇ R(x), wherein A is a projection matrix, x is a three-dimensional volume data to be reconstructed, b is an actual acquired projection image, b [b 0 ,b 1 ...b i ...b N ] T , R is the constraint term, ⁇ is the adjustment coefficient, and for each frame image b i , the three-dimensional position coordinates of the ray source corresponding to the frame image are calculated. ;
  • the iterative unit repeats the first calculation unit and the second calculation unit until a termination condition is reached.
  • the method further includes a thread allocation unit between the first computing unit and the second computing unit, including:
  • a first calculation module assigning a thread to each pixel point of the left half axis of each frame image b i X axis, calculating coordinates of the pixel point in space, and obtaining a coordinate of the source corresponding to the frame image corresponding to the frame image Linear equation
  • Judging module determining whether the straight line intersects each voxel located in the left half axis of the X axis in the three-dimensional volume data to be reconstructed;
  • the three-dimensional position (x, y, z) of the voxel in the three-dimensional volume data is recorded, the length value is set to 0, and the corresponding length of the voxel with respect to the axisymmetric position of the X axis is also set to zero.
  • the second calculation module when all the threads are executed, the projection matrix A ij corresponding to all the pixel points b ij of the current frame image is obtained, and the corresponding component is calculated.
  • M is the number of pixels in the image.
  • the cone beam CT three-dimensional reconstruction method and system provided by the invention adopts the GPU acceleration method, and simultaneously reduces the calculation amount by using geometric symmetry, and establishes a thread for the connection between each pixel and the radiation source, and calculates the connection in three dimensions.
  • the length of each voxel in the volume data, and then the x optimal solution is obtained by an iterative method, thereby achieving the purpose of efficient three-dimensional reconstruction, and the reconstructed image quality is guaranteed while greatly improving the calculation efficiency.
  • FIG. 1 is a flow chart of steps of a three-dimensional reconstruction method of cone beam CT according to an embodiment of the present invention.
  • FIG. 2 is a three-dimensional volume data coordinate system and a projection data coordinate system set according to an embodiment of the present invention.
  • FIG. 3 is a schematic structural diagram of a cone beam CT three-dimensional reconstruction system according to an embodiment of the present invention.
  • FIG. 4 is a schematic structural diagram of a thread allocation unit of a cone beam CT three-dimensional reconstruction system according to an embodiment of the present invention.
  • a three-dimensional reconstruction method of cone beam CT includes the following steps:
  • Step S110 collecting projection data
  • the projection data is stored in an N*H*W three-dimensional array in unsigned short type data.
  • N is the number of projection data frames
  • H is the projection data height
  • W is the projection data width.
  • the three-dimensional volume data coordinate system and the projection data coordinate system are set.
  • the light source rotates around the X axis.
  • 2 + ⁇ R(x), where A is a projection matrix, x is a three-dimensional volume data to be reconstructed, b is an actual acquired projection image, b [ b 0 , b 1 ... b i ... b N ] T , R is a constraint term, ⁇ is an adjustment coefficient, and for each frame image b i , the three-dimensional position coordinates of the ray source corresponding to the frame image are calculated;
  • step S130 the following steps are further performed:
  • Step S121 assigning a thread to each pixel point of the left half axis of each frame image b i X axis, calculating the coordinates of the pixel point in space, and obtaining a linear equation connecting the coordinates with the ray source coordinates corresponding to the frame image. ;
  • Step S122 determining whether the line intersects with each voxel located in the left half axis of the X axis in the three-dimensional volume data to be reconstructed;
  • the three-dimensional position (x, y, z) of the voxel in the three-dimensional volume data is recorded, the length value is set to 0, and the corresponding length of the voxel with respect to the axisymmetric position of the X axis is also set to zero.
  • Step S123 When all the threads are executed, the projection matrix A ij corresponding to all the pixel points b ij of the current frame image is obtained, and the corresponding component is calculated. M is the number of pixels in the image.
  • Step S140 Steps S120 to S130 are repeated until the termination condition is reached.
  • FIG. 3 is a schematic structural diagram of a cone beam CT three-dimensional reconstruction system according to an embodiment of the present invention, including:
  • the projection data collecting unit 110 is configured to collect projection data
  • b [b 0 , b 1 ... b i ... b N ] T
  • R is a constraint term
  • ⁇ is an adjustment coefficient
  • ⁇ is an adjustment coefficient
  • the iteration unit 140 repeats the first calculation unit and the second calculation unit until a termination condition is reached.
  • FIG. 4 is a schematic structural diagram of a thread allocation unit 150 of a cone beam CT three-dimensional reconstruction system according to an embodiment of the present invention.
  • the thread allocation unit 150 is located between the first computing unit 120 and the second computing unit 130, and includes:
  • a first calculating module 151 assigning a thread to each pixel point of the left half axis of each frame image b i X axis, calculating coordinates of the pixel point in space, and obtaining a coordinate of the ray source corresponding to the coordinate of the frame image Linear equation;
  • the determining module 152 determining whether the straight line intersects each voxel located in the left half axis of the X axis in the three-dimensional volume data to be reconstructed;
  • the three-dimensional position (x, y, z) of the voxel in the three-dimensional volume data is recorded, the length value is set to 0, and the length corresponding to the axis-symmetric position of the voxel with respect to the X-axis is also set to zero.
  • the second calculating module 153 when all threads are executed, the projection matrix A ij corresponding to all the pixel points b ij of the current frame image is obtained, and the corresponding component is calculated. M is the number of pixels in the image.
  • the cone beam CT three-dimensional reconstruction method and system provided by the invention adopts the GPU acceleration method, and simultaneously reduces the calculation amount by using geometric symmetry, and establishes a thread for the connection between each pixel and the radiation source, and calculates the connection in three dimensions.
  • the length of each voxel in the volume data, and then the x optimal solution is obtained by an iterative method, thereby achieving the purpose of efficient three-dimensional reconstruction, and the reconstructed image quality is guaranteed while greatly improving the calculation efficiency.
  • the fast calculation method of the cone beam CT projection matrix of the present invention may also have various transformations and modifications, and is not limited to the specific structure of the above embodiment.
  • the scope of the present invention should include such modifications or substitutions and modifications as would be apparent to those skilled in the art.

Landscapes

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

Abstract

一种锥束CT三维重建方法和系统，采用GPU加速的方式，同时利用几何对称性减少计算量，对每一条像素与射线源之间的连线建立一个线程，计算该连线在三维体数据中每一个体素内的长度，进而通过迭代方法得到x最优解，从而达到高效进行三维重建的目的，且使得重建图像质量得到保证的同时大大改善计算效率。

Description

一种锥束CT三维重建方法及系统 技术领域
本发明涉及X射线CT技术领域重建算法技术，尤其是涉及一种锥束CT三维重建方法及系统。
背景技术
相比于传统的CT技术，平板探测器由于有更高的图像分辨率和更大的FOV，因此应用平板探测器技术的锥束CT在成像质量与成像效率上都更有优势。因此如何在锥束CT模型下快速、准确的进行三维重建，成为十分重要的问题。
专利CN104899903A，CN102609978A，CN102279970A提供了基于FDK算法的三维重建。作为滤波反投影累算法，虽然计算速度较快，对系统资源要求低，但重建图像质量差，存在伪影，影响可视化效果和医生诊断。
有文献提出代数迭代算法ART，以解决滤波反投影算法的图像质量问题。ART算法需要计算投影矩阵，目前通常是采用Siddon等提出的算法，该算法是逐点计算方法，即根据射线上某一点和待重建体素间的相对位置关系判断射线上下一点的位置，计算量较大，计算效率较低，尤其是在待重建体数据为三维的情况下缺陷更为明显，限制了ART类算法在临床上的应用。虽然有一些改进算法，但是多数是采用近似假设，如仅判断射线是否穿越体素，穿越记为1，否则记为0的01近似法，会影响投影矩阵的准确度。
以FDK算法为代表的滤波反投影算法，虽然计算速度较快，对系统资源要求低，但重建图像质量差，存在伪影，影响可视化效果和医生诊断。ART 算法需要计算投影矩阵，目前通常是采用Siddon等提出的算法，该算法是逐点计算方法，即根据射线上某一点和待重建体素间的相对位置关系判断射线上下一点的位置，计算量较大，计算效率较低，尤其是在待重建体数据为三维的情况下缺陷更为明显，限制了ART类算法在临床上的应用。虽然有一些改进算法，但是多数是采用近似假设，如仅判断射线是否穿越体素，穿越记为1，否则记为0的01近似法，会影响投影矩阵的准确度。
发明内容
有鉴如此，有必要针对现在技术存在的缺陷，提供一种能够快速进行锥束CT三维重建方法。
为实现上述目的，本发明采用下述技术方案：
一种锥束CT三维重建方法，包括下述步骤：
步骤S110：采集投影数据；
步骤S120：构建第一目标函数，F＝|(Ax-b)|2+λR(x)，其中，A为投影矩阵，x为待重建三维体数据，b为实际采集投影图像，b＝[b0，b1...bi...bN]T，R为约束项，λ为调整系数，对每一帧图像bi，计算对应帧图像的射线源所处三维位置坐标；
步骤S130：构建第二目标函数AT(Axn-b)＝∑Ti，其中，Ti为每一帧图像i对应的分量，并根据
Figure PCTCN2016099565-appb-000001
计算xn+1；
步骤S140：重复步骤S120至S130，直到达到终止条件。
在一些实施例中，步骤S110中，所述投影数据以unsigned short类型数据保存在N*H*W三维数组中，其中N为投影数据帧数，H为投影数据高度，W为投影数据宽度。
在一些实施例中，在完成步骤S120后进行步骤S130前，还包括下述步骤：
步骤S121：计算对每一帧图像biX轴左半轴的每一个像素点分配一个线程，计算该的像素点在空间中的坐标，求得连接该坐标与该帧图像对应射线源坐标的直线方程；
步骤S122：判断该直线与待重建的三维体数据中位于X轴左半轴的每一个体素是否有相交；
若有相交，则计算该直线在该体素内部分的长度值，并记录该体素在三维体数据中的三维位置(x，y，z)，并计算该位置关于X轴的轴对称位置，对应的长度值保持不变；
若无相交，记录该体素在三维体数据中的三维位置(x，y，z)，设长度值为0，且该体素关于X轴的轴对称位置的对应长度也设置为0。
步骤S123：当所有线程执行完成后，即得到对应当前帧图像的所有像素点bij的投影矩阵Aij，并计算对应分量
Figure PCTCN2016099565-appb-000002
M为图像所含像素数。
另外，本发明还提供了一种锥束CT三维重建系统，包括：
投影数据采集单元，用于采集投影数据；
第一计算单元，构建第一目标函数，F＝|(Ax-b)|2+λR(x)，其中，A为投影矩阵，x为待重建三维体数据，b为实际采集投影图像，b＝[b0，b1...bi...bN]T，R为约束项，λ为调整系数，对每一帧图像bi，计算对应帧图像的射线源所处三维位置坐标；
第二计算单元，用于构建第二目标函数AT(Axn-b)＝∑Ti，其中，Ti为每一帧图像i对应的分量，并根据
Figure PCTCN2016099565-appb-000003
计算xn+1；
迭代单元，重复所述第一计算单元及第二计算单元，直到达到终止条件。
在一些实施例中，还包括位于所述第一计算单元和第二计算单元之间的线程分配单元，包括：
第一计算模块：对每一帧图像biX轴左半轴的每一个像素点分配一个线程，计算该像素点在空间中的坐标，求得连接该坐标与该帧图像对应射线源坐标的直线方程；
判断模块：判断该直线与待重建的三维体数据中位于X轴左半轴的每一个体素是否有相交；
若有相交，则计算该直线在该体素内部分的长度值，并记录该体素在三维体数据中的三维位置(x，y，z)，并计算该位置关于X轴的轴对称位置，对应的长度值保持不变；
若无相交，记录该体素在三维体数据中的三维位置(x，y，z)，设长度值为0，且该体素关于X轴的轴对称位置的对应长度也设置为0。
第二计算模块：当所有线程执行完成后，即得到对应当前帧图像的所有像素点bij的投影矩阵Aij，并计算对应分量
Figure PCTCN2016099565-appb-000004
M为图像所含像素数。
本发明采用上述技术方案的优点是：
本发明提供的锥束CT三维重建方法和系统，采用GPU加速的方式，同时利用几何对称性减少计算量，对每一条像素与射线源之间的连线建立一个线程，计算该连线在三维体数据中每一个体素内的长度，进而通过迭代方法得到x最优解，从而达到高效进行三维重建的目的，且使得重建图像质量得到保证的同时大大改善计算效率。
附图说明
图1为本发明实施例提供的锥束CT三维重建方法的步骤流程图。
图2为本发明实施例设定的三维体数据坐标系和投影数据坐标系。
图3为本发明实施例提供的锥束CT三维重建系统的结构示意图。
图4为本发明实施例提供的锥束CT三维重建系统的线程分配单元的结构示意图。
具体实施方式
请参阅图1，本发明实施例提供的一种锥束CT三维重建方法，包括下步骤：
步骤S110：采集投影数据；
优选地，投影数据以unsigned short类型数据保存在N*H*W三维数组中。其中N为投影数据帧数，H为投影数据高度，W为投影数据宽度。
请参阅图2，设定三维体数据坐标系和投影数据坐标系，为后续方便叙述，假设光源绕X轴进行旋转。
步骤S120：构建第一目标函数，F＝|(Ax-b)|2+λR(x)，其中，A为投影矩阵，x为待重建三维体数据，b为实际采集投影图像，b＝[b0，b1...bi...bN]T，R为约束项，λ为调整系数，对每一帧图像bi，计算对应帧图像的射线源所处三维位置坐标；
可以理解，ART类算法的基本方法如下：令目标函数F＝|(Ax-b)|2+λR(x)。其中A为投影矩阵，x为待重建三维体数据，b为投影图像，R为约束项，λ为权重。迭代寻找x，令目标函数F取最小值。
可以理解，基于GPU并行加速技术，对ART类算法进行并行加速，从而达到高效进行三维重建的目的。
步骤S130：构建第二目标函数AT(Axn-b)＝∑Ti，其中，Ti为每一帧图像i对应的分量，并根据
Figure PCTCN2016099565-appb-000005
计算xn+1；
优选地，在进行步骤S130后还进行下述步骤：
步骤S121：对每一帧图像biX轴左半轴的每一个像素点分配一个线程，计算该像素点在空间中的坐标，求得连接该坐标与该帧图像对应射线源坐标的直线方程；
步骤S122：判断该直线与待重建的三维体数据中位于X轴左半轴的每一个体素是否有相交；
若有相交，则计算该直线在该体素内部分的长度值，并记录该体素在三维体数据中的三维位置(x，y，z)，并计算该位置关于X轴的轴对称位置，对应的长度值保持不变；
若无相交，记录该体素在三维体数据中的三维位置(x，y，z)，设长度值为0，且该体素关于X轴的轴对称位置的对应长度也设置为0。
步骤S123：当所有线程执行完成后，即得到对应当前帧图像的所有像素点bij的投影矩阵Aij，并计算对应分量
Figure PCTCN2016099565-appb-000006
M为图像所含像素数。
可以理解，本申请采用最速下降法，即
Figure PCTCN2016099565-appb-000007
其中γ为迭代步长，可以有效降低计算量，提高重建速度。
步骤S140：重复步骤S120至S130，直到达到终止条件。
请参阅图3，为本发明实施例提供的锥束CT三维重建系统的结构示意图，包括：
投影数据采集单元110，用于采集投影数据；
第一计算单元120，构建第一目标函数，F＝|(Ax-b)|2+λR(x)，其中，A为投影矩阵，x为待重建三维体数据，b为实际采集投影图像，b＝[b0，b1...bi...bN]T，R为约束项，λ为调整系数，对每一帧图像bi，计算对应帧图像的射线源所处三维位置坐标；
第二计算单元130，用于构建第二目标函数AT(Axn-b)＝∑Ti，其中，Ti为每一帧图像i对应的分量，并根据
Figure PCTCN2016099565-appb-000008
计算xn+1；
迭代单元140，重复所述第一计算单元及第二计算单元，直到达到终止条件。
请参阅图4，为本发明实施例提供的锥束CT三维重建系统的线程分配单元150的结构示意图。
优选地，线程分配单元150位于所述第一计算单元120和第二计算单元130之间的，包括：
第一计算模块151：对每一帧图像biX轴左半轴的每一个像素点分配一个线程，计算该像素点在空间中的坐标，求得连接该坐标与该帧图像对应射线源坐标的直线方程；
判断模块152：判断该直线与待重建的三维体数据中位于X轴左半轴的每一个体素是否有相交；
若有相交，则计算该直线在该体素内部分的长度值，并记录该体素在三维体数据中的三维位置(x，y，z)，并计算该位置关于X轴的轴对称位置，对应的长度值保持不变；
若无相交，记录该体素在三维体数据中的三维位置(x，y，z)，设长度值为0，且该体素关于X轴的轴对称位置对应的长度也设置为0。
第二计算模块153：当所有线程执行完成后，即得到对应当前帧图像的所有像素点bij的投影矩阵Aij，并计算对应分量
Figure PCTCN2016099565-appb-000009
M为图像所含像素数。
本发明提供的锥束CT三维重建方法和系统，采用GPU加速的方式，同时利用几何对称性减少计算量，对每一条像素与射线源之间的连线建立一个线程，计算该连线在三维体数据中每一个体素内的长度，进而通过迭代方法得到x最优解，从而达到高效进行三维重建的目的，且使得重建图像质量得到保证的同时大大改善计算效率。
当然本发明的锥束CT投影矩阵的快速计算方法还可具有多种变换及改型，并不局限于上述实施方式的具体结构。总之，本发明的保护范围应包括那些对于本领域普通技术人员来说显而易见的变换或替代以及改型。

Claims (5)

  1. 一种锥束CT三维重建方法，其特征在于，包括下述步骤：
    步骤S110：采集投影数据；
    步骤S120：构建第一目标函数，F＝|(Ax-b)|2+λR(x)，其中，A为投影矩阵，x为待重建三维体数据，b为实际采集投影图像，b＝[b0，b1...bi...bN]T，R为约束项，λ为调整系数，对每一帧图像bi，计算对应帧图像的射线源所处三维位置坐标；
    步骤S130：构建第二目标函数AT(Axn-b)＝∑Ti，其中，Ti为每一帧图像i对应的分量，并根据
    Figure PCTCN2016099565-appb-100001
    计算xn+1；
    步骤S140：重复步骤S120至S130，直到达到终止条件。
  2. 根据权利要求1所述的锥束CT三维重建方法，其特征在于，步骤S110中，所述投影数据以unsigned short类型数据保存在N*H*W三维数组中，其中N为投影数据帧数，H为投影数据高度，W为投影数据宽度。
  3. 根据权利要求1所述的锥束CT三维重建方法，其特征在于，在完成步骤S120后进行步骤S130前，还包括下述步骤：
    步骤S121：对每一帧图像biX轴左半轴的每一个像素点分配一个线程，计算该像素点在空间中的坐标，求得连接该坐标与该帧图像对应射线源坐标的直线方程；
    步骤S122：判断该直线与待重建的三维体数据中位于X轴左半轴的每一个体素是否有相交；
    若有相交，则计算该直线在该体素内部分的长度值，并记录该体素在三维体数据中的三维位置(x，y，z)，并计算该位置关于X轴的轴对称位置，对应的长度值保持不变；
    若无相交，记录该体素在三维体数据中的三维位置(x，y，z)，设长度值为0，且该体素关于X轴的轴对称位置的对应长度也设置为0；
    步骤S123：当所有线程执行完成后，即得到对应当前帧图像的所有像素点bij的投影矩阵Aij，并计算对应分量
    Figure PCTCN2016099565-appb-100002
    M为图像所含像素数°
  4. 一种锥束CT三维重建系统，其特征在于，包括：
    投影数据采集单元，用于采集投影数据；
    第一计算单元，构建第一目标函数，F＝|(Ax-b)|2+λR(x)，其中，A为投影矩阵，x为待重建三维体数据，b为实际采集投影图像，b＝[b0，b1...bi...bN]T，R为约束项，λ为调整系数，对每一帧图像bi，计算对应帧图像的射线源所处三维位置坐标；
    第二计算单元，用于构建第二目标函数AT(Axn-b)＝∑Ti，其中，Ti为每一帧图像i对应的分量，并根据
    Figure PCTCN2016099565-appb-100003
    计算xn+1；
    迭代单元，重复所述第一计算单元及第二计算单元，直到达到终止条件。
  5. 根据权利要求4所述的锥束CT三维重建系统，其特征在于，还包括位于所述第一计算单元和第二计算单元之间的线程分配单元，包括：
    第一计算模块：对每一帧图像biX轴左半轴的每一个像素点分配一个线程，计算该像素点在空间中的坐标，求得连接该坐标与该帧图像对应射线源坐标的直线方程；
    判断模块：判断该直线与待重建的三维体数据中位于X轴左半轴的每一个体素是否有相交；
    若有相交，则计算该直线在该体素内部分的长度值，并记录该体素在三维体数据中的三维位置(x，y，z)，并计算该位置关于X轴的轴对称位置，对应的长度值保持不变；
    若无相交，记录该体素在三维体数据中的三维位置(x，y，z)，设长度值为0，且该体素关于X轴的轴对称位置的对应长度也设置为0；
    第二计算模块：当所有线程执行完成后，即得到对应当前帧图像的所有像素点bij的投影矩阵Aij，并计算对应分量
    Figure PCTCN2016099565-appb-100004
    M为图像所含像素数。
PCT/CN2016/099565 2016-08-24 2016-09-21 一种锥束ct三维重建方法及系统 WO2018035905A1 (zh)

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN201610712340.6A CN107784684B (zh) 2016-08-24 2016-08-24 一种锥束ct三维重建方法及系统
CN201610712340.6 2016-08-24

Publications (1)

Publication Number Publication Date
WO2018035905A1 true WO2018035905A1 (zh) 2018-03-01

Family

ID=61245337

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2016/099565 WO2018035905A1 (zh) 2016-08-24 2016-09-21 一种锥束ct三维重建方法及系统

Country Status (2)

Country Link
CN (1) CN107784684B (zh)
WO (1) WO2018035905A1 (zh)

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110969694A (zh) * 2019-11-19 2020-04-07 上海科技大学 无约束式扫描和基于体素的三维实时骨骼成像方法
CN111640054A (zh) * 2020-05-18 2020-09-08 扬州哈工博浩智能科技有限公司 一种基于gpu加速的三维重建方法
CN113963056A (zh) * 2021-09-07 2022-01-21 于留青 Ct图像重建方法、装置、电子设备以及存储介质
CN117726752A (zh) * 2023-12-14 2024-03-19 凝动万生医疗科技(武汉)有限公司 一种快速的三维头颅数据正位校准方法

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116758183B (zh) * 2023-08-18 2023-11-07 有方（合肥）医疗科技有限公司 一种cbct图像重建方法

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20050259780A1 (en) * 2004-05-20 2005-11-24 Eastman Kodak Company Method for reconstructing tomographic images
CN101283913A (zh) * 2008-05-30 2008-10-15 首都师范大学 Ct图像重建的gpu加速方法
CN102279970A (zh) * 2010-06-13 2011-12-14 清华大学 基于gpu的螺旋锥束ct重建方法
CN102779350A (zh) * 2012-06-07 2012-11-14 中国人民解放军信息工程大学 一种锥束ct迭代重建算法投影矩阵构建方法
CN103514629A (zh) * 2012-06-22 2014-01-15 密执安大学评议会 用于迭代重建的方法和设备
CN104574459A (zh) * 2014-12-29 2015-04-29 沈阳东软医疗系统有限公司 一种pet图像的重建方法和装置

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7929747B2 (en) * 2006-04-25 2011-04-19 Wisconsin Alumni Research Foundation System and method for estimating data missing from CT imaging projections
CN102609978B (zh) * 2012-01-13 2014-01-22 中国人民解放军信息工程大学 基于cuda架构的gpu加速锥束ct图像重建的方法
CN104424625A (zh) * 2013-09-04 2015-03-18 中国科学院深圳先进技术研究院 一种gpu加速cbct图像重建方法和装置

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20050259780A1 (en) * 2004-05-20 2005-11-24 Eastman Kodak Company Method for reconstructing tomographic images
CN101283913A (zh) * 2008-05-30 2008-10-15 首都师范大学 Ct图像重建的gpu加速方法
CN102279970A (zh) * 2010-06-13 2011-12-14 清华大学 基于gpu的螺旋锥束ct重建方法
CN102779350A (zh) * 2012-06-07 2012-11-14 中国人民解放军信息工程大学 一种锥束ct迭代重建算法投影矩阵构建方法
CN103514629A (zh) * 2012-06-22 2014-01-15 密执安大学评议会 用于迭代重建的方法和设备
CN104574459A (zh) * 2014-12-29 2015-04-29 沈阳东软医疗系统有限公司 一种pet图像的重建方法和装置

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110969694A (zh) * 2019-11-19 2020-04-07 上海科技大学 无约束式扫描和基于体素的三维实时骨骼成像方法
CN111640054A (zh) * 2020-05-18 2020-09-08 扬州哈工博浩智能科技有限公司 一种基于gpu加速的三维重建方法
CN113963056A (zh) * 2021-09-07 2022-01-21 于留青 Ct图像重建方法、装置、电子设备以及存储介质
CN117726752A (zh) * 2023-12-14 2024-03-19 凝动万生医疗科技(武汉)有限公司 一种快速的三维头颅数据正位校准方法

Also Published As

Publication number Publication date
CN107784684A (zh) 2018-03-09
CN107784684B (zh) 2021-05-25

Similar Documents

Publication Publication Date Title
WO2018035905A1 (zh) 一种锥束ct三维重建方法及系统
US20220292646A1 (en) System and method for image reconstruction
Jia et al. GPU-based iterative cone-beam CT reconstruction using tight frame regularization
US8600132B2 (en) Method and apparatus for motion correcting medical images
US9524567B1 (en) Method and system for iterative computed tomography reconstruction
JP5543448B2 (ja) 高効率コンピュータ断層撮影
US9256967B2 (en) Systems and methods for partial volume correction in PET penalized-likelihood image reconstruction
CN106846465B (zh) 一种ct三维重建方法及系统
KR101493375B1 (ko) 화상처리장치, 화상처리방법, 및 컴퓨터 판독가능한 기억매체
US20170178366A1 (en) System and method for image reconstruction
JP6026214B2 (ja) 連続マルチスケール再構成において詳細画像を補うｘ線コンピュータ断層撮像装置（ｘ線ｃｔ装置）、医用画像処理装置及び医用画像処理方法
JP2018064627A (ja) 放射線撮影装置、放射線撮影システム、放射線撮影方法、及びプログラム
CN103198497A (zh) 确定运动场和利用运动场进行运动补偿重建的方法和系统
CN103027705A (zh) 产生运动补偿的ct 图像数据组的方法和系统
CN103479379B (zh) 一种倾斜螺旋扫描的图像重建方法及装置
WO2015144069A1 (zh) 图像滤波方法及ct系统
US8488850B2 (en) Isotropic resolution image reconstruction
US20120288178A1 (en) Medical image processing apparatus, a medical imaging apparatus, and a method of processing medical images
JP2015231528A (ja) Ｘ線コンピュータ断層撮像装置及び医用画像処理装置
WO2010021387A1 (ja) コンプトンカメラ
US9495770B2 (en) Practical model based CT construction
CN103310471A (zh) Ct 图像生成装置及方法、ct 图像生成系统
JP7382496B2 (ja) Ｐｅｔイメージングシステムのセットアップおよび品質管理のための短軸ファントムを用いた連続寝台移動収集
US9336598B2 (en) Shape-extraction method and shape-extraction system
CN110853113A (zh) 基于bpf的tof-pet图像重建算法及重建系统

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 16913926

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

32PN Ep: public notification in the ep bulletin as address of the adressee cannot be established

Free format text: NOTING OF LOSS OF RIGHTS PURSUANT TO RULE 112(1) EPC (EPO FORM 1205A DATED 17/06/2019)

122 Ep: pct application non-entry in european phase

Ref document number: 16913926

Country of ref document: EP

Kind code of ref document: A1