CN106994023B - 锥形束计算机断层成像系统的几何参数确定方法 - Google Patents
锥形束计算机断层成像系统的几何参数确定方法 Download PDFInfo
- Publication number
- CN106994023B CN106994023B CN201710393440.1A CN201710393440A CN106994023B CN 106994023 B CN106994023 B CN 106994023B CN 201710393440 A CN201710393440 A CN 201710393440A CN 106994023 B CN106994023 B CN 106994023B
- Authority
- CN
- China
- Prior art keywords
- projection data
- cbct
- data set
- geometric
- coordinate
- 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.)
- Expired - Fee Related
Links
- 238000007408 cone-beam computed tomography Methods 0.000 title claims abstract description 206
- 238000000034 method Methods 0.000 title claims abstract description 68
- 238000003384 imaging method Methods 0.000 title claims abstract description 37
- 238000005070 sampling Methods 0.000 claims abstract description 68
- 238000012937 correction Methods 0.000 claims description 80
- 239000003550 marker Substances 0.000 claims description 29
- 239000002184 metal Substances 0.000 claims description 11
- 238000004364 calculation method Methods 0.000 claims description 10
- 239000011324 bead Substances 0.000 claims description 9
- 239000000463 material Substances 0.000 claims description 6
- 239000000126 substance Substances 0.000 claims description 4
- 229910000831 Steel Inorganic materials 0.000 claims description 3
- 239000010959 steel Substances 0.000 claims description 3
- 229920005372 Plexiglas® Polymers 0.000 claims 1
- 239000006187 pill Substances 0.000 claims 1
- 239000004926 polymethyl methacrylate Substances 0.000 claims 1
- 239000008188 pellet Substances 0.000 description 12
- 230000005855 radiation Effects 0.000 description 7
- 238000013519 translation Methods 0.000 description 7
- 239000011159 matrix material Substances 0.000 description 6
- 238000012545 processing Methods 0.000 description 6
- 239000013598 vector Substances 0.000 description 5
- 238000013507 mapping Methods 0.000 description 4
- 238000005457 optimization Methods 0.000 description 4
- 238000010586 diagram Methods 0.000 description 3
- 239000011521 glass Substances 0.000 description 3
- 230000000007 visual effect Effects 0.000 description 3
- 238000001514 detection method Methods 0.000 description 2
- 238000003745 diagnosis Methods 0.000 description 2
- 238000009826 distribution Methods 0.000 description 2
- 238000004519 manufacturing process Methods 0.000 description 2
- 238000005259 measurement Methods 0.000 description 2
- 238000009825 accumulation Methods 0.000 description 1
- 230000003044 adaptive effect Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000002591 computed tomography Methods 0.000 description 1
- 238000000354 decomposition reaction Methods 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 230000009191 jumping Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000011218 segmentation Effects 0.000 description 1
Classifications
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B6/00—Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment
- A61B6/02—Arrangements for diagnosis sequentially in different planes; Stereoscopic radiation diagnosis
- A61B6/03—Computed tomography [CT]
- A61B6/032—Transmission computed tomography [CT]
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B6/00—Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment
- A61B6/40—Arrangements for generating radiation specially adapted for radiation diagnosis
- A61B6/4064—Arrangements for generating radiation specially adapted for radiation diagnosis specially adapted for producing a particular type of beam
- A61B6/4085—Cone-beams
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B6/00—Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment
- A61B6/52—Devices using data or image processing specially adapted for radiation diagnosis
- A61B6/5205—Devices using data or image processing specially adapted for radiation diagnosis involving processing of raw data to produce diagnostic data
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B6/00—Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment
- A61B6/52—Devices using data or image processing specially adapted for radiation diagnosis
- A61B6/5258—Devices using data or image processing specially adapted for radiation diagnosis involving detection or reduction of artifacts or noise
Landscapes
- Health & Medical Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- Engineering & Computer Science (AREA)
- Medical Informatics (AREA)
- Radiology & Medical Imaging (AREA)
- Molecular Biology (AREA)
- Biophysics (AREA)
- Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
- Optics & Photonics (AREA)
- Pathology (AREA)
- Physics & Mathematics (AREA)
- Biomedical Technology (AREA)
- Heart & Thoracic Surgery (AREA)
- High Energy & Nuclear Physics (AREA)
- Surgery (AREA)
- Animal Behavior & Ethology (AREA)
- General Health & Medical Sciences (AREA)
- Public Health (AREA)
- Veterinary Medicine (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Pulmonology (AREA)
- Theoretical Computer Science (AREA)
- Apparatus For Radiation Diagnosis (AREA)
Abstract
本发明涉及一种锥形束计算机断层成像系统的几何参数确定方法，包括：分别在CBCT系统探测器未偏置状态以及CBCT系统探测器水平移动至设定偏置状态，获取CBCT系统探测器未偏置状态下的第一投影数据，以及CBCT系统探测器水平移动至设定偏置状态下的第二投影数据；在第一投影数据中获取投影数据集A，计算投影数据集A的中心点坐标，在第二投影数据中获取投影数据集B，计算投影数据集B的中心点坐标；根据投影数据集A的中心点坐标、投影数据集B的中心点坐标、以及投影数据集A以及投影数据集B的之间的对应关系确定估计坐标以及目标函数；计算目标函数取到最小值时的估计坐标，得到优化估计值，根据优化估计值确定CBCT系统在该次采样过程中的几何参数。
Description
技术领域
本发明涉及图像处理技术领域，特别是涉及一种锥形束计算机断层成像系统的几何参数确定方法。
背景技术
锥形束计算机断层成像系统，简称CBCT(Cone-Beam Computerized Tomography)系统，在医学影像诊断领域发挥重要作用。然而，由于CBCT系统的探测器面积有限和物体尺寸较大等因素，被扫描物体在每个视角下的投影信息不能完全包含在探测器范围内，限制了成像范围(如图1所示)。而偏置探测器CBCT系统可以在不增加探测器成本的基础上实现更大的成像范围(如图2所示)，图1为不偏置探测器CBCT系统成像俯视图，图2为偏置探测器CBCT系统成像俯视图。偏置探测器CBCT系统具体原理为，探测器进行水平方向的偏移(假设向左偏移)并固定后进行扫描，每个视角下锥形X射线束只覆盖物体的一半或者一半多一些(物体的左半信息)，利用每个视角下的投影数据进行后续的CBCT图像重建，得到可以反映物体三维结构的CBCT图像。偏置探测器CBCT系统优势在于，第一，实现了更大范围的CBCT成像面积，第二，在大范围CBCT成像面积的前提下又减少了对患者的辐射剂量(减少近一半的剂量)。
精确获取描述CBCT系统的几何参数是重建高质量CBCT图像的前提条件，它影响了医生对患者病情的精确诊断。因此，CBCT系统的几何校正非常必要。
几何校正方法可以分为单角度校正以及全角度校正两大类。所谓单角度校正是指，CBCT系统在各个成像角度下扫描的几何位置均需要进行位置校正，每个角度下获取一套校正参数来描述当前角度下的投照几何关系。全角度校正是指，假定CBCT系统扫描过程是按照严格的圆轨道运行，仅用一套校正参数描述圆轨道扫描的几何位置。
传统的锥束CT圆轨道半覆盖扫描的全角度几何校正方法中，其探测器平移后只覆盖物体的一半信息，容易忽略倾斜角的计算，人为假定值为0，不可以避免带来系统几何校正误差，影响图像重建质量。几何校正实验是在探测器和球管固定、载物台旋转的方式下进行的，这种设计下的机械精度和稳定性在一定程度上保证了可控性。然而，对于病人的扫描，需要探测器和球管的圆周运动，病人保持体位不变，这就导致安装了探测器和球管的机架旋转过程中的机械抖动问题，等价于探测器和球管固定时物体的“跳动”，而传统方案并未考虑这种情况，因而在适用于人体扫描的偏置探测器CBCT系统中，容易使几何参数校正的精度下降，影响CBCT图像重建质量。
发明内容
基于此，有必要针对传统方案容易使几何参数校正的精度下降，影响CBCT图像重建质量的技术问题，提供一种锥形束计算机断层成像系统的几何参数确定方法。
一种锥形束计算机断层成像系统的几何参数确定方法，包括如下步骤：
分别在CBCT系统探测器未偏置状态以及CBCT系统探测器水平移动至设定偏置状态，对CBCT几何参数校正装置中的小球标记物进行扫描，获取CBCT系统探测器未偏置状态下的第一投影数据，以及CBCT系统探测器水平移动至设定偏置状态下的第二投影数据；
在第一投影数据中采样得到投影数据集A，计算投影数据集A的中心点坐标，在第二投影数据中采样得到投影数据集B，计算投影数据集B的中心点坐标；
根据所述投影数据集A的中心点坐标、投影数据集B的中心点坐标、以及所述投影数据集A以及投影数据集B的之间的对应关系确定CBCT几何参数校正装置中各个小球标记物球心投影的估计坐标以及CBCT几何参数校正装置的目标函数；
计算所述目标函数取到最小值时的估计坐标，得到优化估计值，根据所述优化估计值确定CBCT系统在该次采样过程中的几何参数。
上述锥形束计算机断层成像系统的几何参数确定方法，可以分别在CBCT系统探测器未偏置状态以及CBCT系统探测器水平移动至设定偏置状态，对CBCT几何参数校正装置中的小球标记物进行圆周扫描，采样获取CBCT系统探测器未偏置状态下的第一投影数据，以及CBCT系统探测器水平移动至设定偏置状态下的第二投影数据；还可以在第一投影数据中采样得到投影数据集A，在第二投影数据中采样得到的投影数据集B，计算所述投影数据集A的中心点坐标MA(uk,vk)和投影数据集B的中心点坐标MB(um,vm)；以便依据上述投影数据集A的中心点坐标MA(uk,vk)、投影数据集B的中心点坐标MB(um,vm)确定投影数据集A和投影数据集B之间的对应关系，进而确定CBCT几何参数校正装置中各个小球标记物球心投影的估计坐标Est(u'n,v'n)以及CBCT几何参数校正装置的目标函数；计算所述目标函数取到最小值时的估计坐标Est(u'n,v'n)，确定优化估计值，根据所述优化估计值确定CBCT系统第i次采样过程中的几何参数，即上述目标函数取到最小值时的Est(u'n,v'n)所对应的CBCT系统几何参数，使所确定的几何参数具有较高的精度，可以提高依据上述CBCT系统几何参数所重建的CBCT图像的质量。
附图说明
图1为一个实施例的不偏置探测器CBCT系统成像俯视图；
图2为一个实施例的偏置探测器CBCT系统成像俯视图；
图3为一个实施例的锥形束计算机断层成像系统的几何参数确定方法流程图；
图4为一个实施例的基于扫描相机模型的几何关系图；
图5为一个实施例的CBCT几何参数校正装置的外观图；
图6为一个实施例的投影数据集A和投影数据集B的水平平移关系示意图。
具体实施方式
下面结合附图对本发明提供的锥形束计算机断层成像系统的几何参数确定方法的具体实施方式进行详细阐述。
参考图3所示，图3为一个实施例的锥形束计算机断层成像系统的几何参数确定方法流程图，包括如下步骤：
S10，分别在CBCT系统探测器未偏置状态以及CBCT系统探测器水平移动至设定偏置状态，对CBCT几何参数校正装置中的小球标记物进行扫描，获取CBCT系统探测器未偏置状态下的第一投影数据，以及CBCT系统探测器水平移动至设定偏置状态下的第二投影数据；
参考图1所示，上述CBCT系统探测器未偏置状态，CBCT系统可以包括射线源1、探测器2以及全覆盖扫描的成像区域3。参考图2所示，上述CBCT系统探测器未偏置状态，CBCT系统可以包括射线源1、探测器2以及偏移探测器CBCT系统的成像区域4。一般情况下偏移探测器后的成像区域4比不偏移探测器的成像区域3大。
参考图4所示，CBCT系统可以通过11个几何参数θx、θy、θz、Tx、Ty、Tz、u0、v0、f1、f2和dt进行描述，其对应的系统定义可以为：定义右手笛卡尔世界坐标系O-XwYwZw，以被扫描物体中心(图示立方体)为系统坐标原点，其中Yw为被扫描物体的纵轴中轴线，Xw为被扫描物体的横轴中轴线，Zw为分别垂直Xw和Yw的坐标轴，分别为Xw，Yw，Zw轴的单位方向向量；定义图像坐标系I-UiVi，探测器的右下角为坐标原点，u，v分别为探测器上行和列方向，其中和分别是探测器的行和列的方向向量，为探测器的法向量；定义成像坐标系S-XcYcZc，坐标原点为射线源，Xc，Yc分别与图像坐标系u，v平行，Zc垂直于探测器。
S20，在第一投影数据中采样得到投影数据集A，计算投影数据集A的中心点坐标，在第二投影数据中采样得到投影数据集B，计算投影数据集B的中心点坐标；
具体地，上述步骤可以在第一投影数据中获取第i次采样得到的投影数据集A，计算投影数据集A的中心点坐标MA(um,vm)，在第二投影数据中获取第i次采样得到的投影数据集B，计算投影数据集B的中心点坐标MB(um,vm)；其中，m≤N，m为第二投影数据中第i次采样所包括小球标记物像的个数，N为CBCT几何参数校正装置中的小球标记物总数。
上述步骤S20中，可以在计算得到投影数据集A的中心点坐标MA(um,vm)和投影数据集B的中心点坐标MB(um,vm)，进一步确定投影数据集A与投影数据集B之间各个元素之间的对应关系。具体地，可以在获取到第i次采样得到的投影数据集A和投影数据集B后，可以利用投影数据集A或者投影数据集B的水平平移关系，MA(uk,vk)中uk行坐标可以作为标号索引，通过索引可以计算出数投影数据集中小球球心投影坐标MB(um,vm)的标号m，若um落在uk的邻域范围内，则可以确定投影数据集B当前标号m＝k，以此确定投影数据集A与投影数据集B之间各个元素之间的对应关系；上述标号k的作用为找到校正装置第k个小球对应的投影信息。
S30，根据所述投影数据集A的中心点坐标、投影数据集B的中心点坐标、以及所述投影数据集A以及投影数据集B的之间的对应关系确定CBCT几何参数校正装置中各个小球标记物球心投影的估计坐标以及CBCT几何参数校正装置的目标函数；
S40，计算所述目标函数取到最小值时的估计坐标，得到优化估计值，根据所述优化估计值确定CBCT系统在该次采样过程中的几何参数。
本实施例提供的锥形束计算机断层成像系统的几何参数确定方法，可以分别在CBCT系统探测器未偏置状态以及CBCT系统探测器水平移动至设定偏置状态，对CBCT几何参数校正装置中的小球标记物进行圆周扫描，采样获取CBCT系统探测器未偏置状态下的第一投影数据，以及CBCT系统探测器水平移动至设定偏置状态下的第二投影数据；还可以在第一投影数据中采样得到投影数据集A，在第二投影数据中采样得到的投影数据集B，计算所述投影数据集A的中心点坐标MA(uk,vk)和投影数据集B的中心点坐标MB(um,vm)；以便依据上述投影数据集A的中心点坐标MA(uk,vk)、投影数据集B的中心点坐标MB(um,vm)、确定投影数据集A和投影数据集B之间的对应关系，进而确定CBCT几何参数校正装置中各个小球标记物球心投影的估计坐标Est(u'n,v'n)以及CBCT几何参数校正装置的目标函数；计算所述目标函数取到最小值时的估计坐标Est(u'n,v'n)，确定优化估计值，根据所述优化估计值确定CBCT系统第i次采样过程中的几何参数，即上述目标函数取到最小值时的Est(u'n,v'n)所对应的CBCT系统几何参数，使所确定的几何参数具有较高的精度，可以提高依据上述CBCT系统几何参数所重建的CBCT图像的质量。
在一个实施例中，分别在CBCT系统探测器未偏置状态以及CBCT系统探测器水平移动至设定偏置状态，可以对CBCT几何参数校正装置中的小球标记物进行圆周扫描，在扫描过程中，采样获取CBCT系统探测器未偏置状态下的第一投影数据，以及CBCT系统探测器水平移动至设定偏置状态下的第二投影数据；
对获取第一投影数据以及第二投影数据所采用的采样间隔或者采样点等信息完全一致。具体地，CBCT系统探测器未偏置状态，即为CBCT系统中心射线投射至成像探测器中心位置，CBCT几何参数校正装置放置在成像区域范围内，通过激光定位保证CBCT几何参数校正装置中轴线放置于机架旋转轴附近，并加以固定，可以利用75kV/4mA曝光条件下的射线对校正装置进行圆周扫描，得到不同视角下的投影图像，采样间隔可以为0.9度，一共有400个角度下的投影数据A(在第一投影数据中采样得到的投影数据集)，此时每个角度下的投影数据中包含了全部的CBCT几何参数校正装置中N个小球标记物(如17个金属小球)的影像信息，从上至下可以依次标号为第1个小球，第2个小球，…，第N个小球。CBCT系统探测器水平移动至设定偏置状态，即将CBCT系统中的探测器水平移动到设定好的偏置状态，CBCT几何参数校正装置保持不变，仍固定于之前的摆放位置，再次利用同样曝光条件进行圆周扫描，得到不同视角下的投影图像，采样间隔相应的为0.9度，一共有400个角度下的投影数据B，此时每个角度下的投影数据中仅包含了CBCT几何参数校正装置中部分小球的影像信息。上述设定偏置状态可以根据CBCT系统的成像范围(FOV)确定，通常情况下，不同的CBCT系统有不一样的成像要求，需要按实际情况分析。
上述所述在第二投影数据中采样得到投影数据集B，计算投影数据集B的中心点坐标的过程可以包括：
在第二投影数据中获取第i次采样得到的投影数据集B，计算投影数据集B的中心点坐标MB(um,vm)；其中，m≤N，m为第二投影数据中第i次采样所包括小球标记物像的个数，N为CBCT几何参数校正装置中的小球标记物总数；
上述过程可以首先识别投影数据集B中各个像素点的灰度值，依据上述灰度值对投影数据集B进行二值化处理，再对二值化处理后的投影数据集B进行相应计算，确定MB(um,vm)。
CBCT几何参数校正装置中各个小球标记物球心投影的估计坐标为Est(u'n,v'n)，可以根据上述Est(u'n,v'n)以及MB(um,vm)确定CBCT几何参数校正装置的目标函数；
上述小球标记物球心投影的估计坐标Est(u'n,v'n)可以根据相应小球标记物的球心坐标(xn，yn，zn)以及相应的几何参数所确定。通常情况下，上述目标函数可以依据各个小球标记物对应的Est(u'n,v'n)与MB(um,vm)之差的平方和表征。
上述步骤可以根据投影数据集B中小球(小球标记物对应的投影)的中心点坐标和标号，以及N个小球球心位置(xn,yn,zn)，n＝1,2,…，N，再借鉴相应扫描相机模型建立各个角度下校正装置上小球球心空间位置与其投影图像球心位置之间的联系，通过这个关系获取小球球心在投影图像上的估计点(估计坐标)与实际测量点之间的目标函数，并通过CMA-ES(Covariance Matrix Adaptation Evolution Strategy，协方差矩阵自适应进化策略优化算法)等最优化算法求解相应角度下的几何校正参数。
上述计算所述目标函数取到最小值时的估计坐标，得到优化估计值，根据所述优化估计值确定CBCT系统在该次采样过程中的几何参数的过程可以包括：
计算所述目标函数取到最小值时的Est(u'n,v'n)，得到优化估计值，根据所述优化估计值确定CBCT系统探测器第i次采样过程中的几何参数。
上述N个球坐标估计值Est(u'n,v'n)均存在对应的一组CBCT系统几何参数，在得到优化估计值后，将上述得到优化估计值对应的一组CBCT系统几何参数确定为CBCT系统探测器第i次采样过程中的几何参数，便可以有效提高所确定的几何参数的精度。
本实施例提供的锥形束计算机断层成像系统的几何参数确定方法，可以分别在CBCT系统探测器未偏置状态以及CBCT系统探测器水平移动至设定偏置状态，对CBCT几何参数校正装置中的小球标记物进行圆周扫描，采样获取CBCT系统探测器未偏置状态下的第一投影数据，以及CBCT系统探测器水平移动至设定偏置状态下的第二投影数据；还可以在第一投影数据中采样得到投影数据集A，在第二投影数据中获取第i次采样得到的投影数据集B，计算投影数据集A的中心点坐标MA(uk,vk)和投影数据集B的中心点坐标MB(um,vm)，并确定投影数据集A和投影数据集B之间的对应关系，以此确定CBCT几何参数校正装置中各个小球标记物球心投影的估计坐标Est(u'n,v'n)以及CBCT几何参数校正装置的目标函数；计算所述目标函数取到最小值时的Est(u'n,v'n)，确定优化估计值，根据所述优化估计值确定CBCT系统探测器第i次采样过程中的几何参数，即上述目标函数取到最小值时的Est(u'n,v'n)所对应的CBCT系统几何参数，使所确定的几何参数具有较高的精度，可以提高依据上述CBCT系统几何参数所重建的CBCT图像的质量。
在一个实施例中，上述分别在CBCT系统探测器未偏置状态以及CBCT系统探测器水平移动至设定偏置状态，对CBCT几何参数校正装置中的小球标记物进行扫描，获取CBCT系统探测器未偏置状态下的第一投影数据，以及CBCT系统探测器水平移动至设定偏置状态下的第二投影数据的过程之前，还可以包括：
将CBCT几何参数校正装置的圆柱轴线和CBCT系统的旋转轴重合，在所述CBCT系统中设置CBCT几何参数校正装置。
本实施例可以通过激光定位等相关技术将CBCT几何参数校正装置的圆柱轴线和CBCT系统的旋转轴重合，以保证后续所获取的第一投影数据和第二投影数据的质量。
在一个实施例中，上述CBCT几何参数校正装置可以包括空心圆柱以及小球标记物；所述空心圆柱的外表面上设计沿着中心轴方向螺旋轨迹的小孔，小球标记物镶嵌在小孔内，小球数量大于或等于6个。
作为一个实施例，上述空心圆柱的制作材料包括有机玻璃等低密度物质；所述小球标记物的制作材料包括金属等高密度物质。
作为一个实施例，上述空心圆柱的制作材料可以为有机玻璃，所述小球标记物的制作材料可以为钢。
作为一个实施例，上述镶嵌在小孔内的小球标记物按照中心密两边疏的规则进行分布。
具体地，上述CBCT几何参数校正装置的外观图可以参考图5所示，包括：空心圆柱，金属小球(小球标记物)。空心圆柱表面上设计沿着中心轴方向螺旋轨迹的小孔，小球标记物镶嵌在小孔内。中空圆柱的直径可以为50mm(毫米)，小孔和小球数量均为17个，其中16个小孔和小球的直径为1mm，1个小孔和小球的直径为2mm，作为坐标参考点。CBCT几何参数校正装置空心圆柱材料为有机玻璃；校正装置小球标记物材料为钢。
上述CBCT几何参数校正装置上小球分布可以按照中心密两边疏的规则，第7个到第11个范围内的小球，两个小球之间的高度差为2.5mm；第1个到第6个范围内的小球，两个小球之间的高度差为4mm；第12个到第17个范围内的小球，两个小球之间的高度差为4mm。其中第6个到第7个小球之间的高度差为3mm，第11个到第12个小球之间的高度差为3mm。
在一个实施例中，上述在第一投影数据中获取第i次采样得到的投影数据集A，在第二投影数据中获取第i次采样得到的投影数据集B的过程之后，还可以包括：
在第一投影数据中获取第i次采样得到的投影数据集A，识别投影数据集A中各个像素点的灰度值，将投影数据集A中灰度值大于灰度阈值T的像素点灰度置1，将投影数据集A中灰度值小于或者等于灰度阈值T的像素点灰度置0；
识别投影数据集B中各个像素点的灰度值，将投影数据集B中灰度值大于灰度阈值T的像素点灰度置1，将投影数据集B中灰度值小于或者等于灰度阈值T的像素点灰度置0；
根据灰度更新后的投影数据集A计算MA(uk,vk)，根据灰度更新后的投影数据集A计算MB(um,vm)；其中，k≤N，k为第一投影数据中第i次采样所包括小球标记物像的个数。
本实施例对获取第一投影数据以及第二投影数据所采用的采样间隔或者采样点等信息完全一致。具体地，CBCT系统探测器未偏置状态，可以利用75kV/4mA曝光条件下的射线对校正装置进行圆周扫描，得到不同视角下的投影图像，采样间隔可以为0.9度，一共有400个角度下的投影数据A(在第一投影数据中采样得到的投影数据集)，此时每个角度下的投影数据中包含了全部的CBCT几何参数校正装置中N个小球标记物(如17个金属小球)的影像信息，从上至下可以依次标号为第1个小球，第2个小球，…，第17个小球。CBCT系统探测器水平移动至设定偏置状态，即将CBCT系统中的探测器水平移动到设定好的偏置状态，CBCT几何参数校正装置保持不变，仍固定于之前的摆放位置，再次利用同样曝光条件进行圆周扫描，得到不同视角下的投影图像，采样间隔相应的为0.9度，一共有400个角度下的投影数据B，此时每个角度下的投影数据中仅包含了CBCT几何参数校正装置中部分小球的影像信息。
上述灰度阈值T可以依据小球标记物在第一投影数据和第二投影数据上所成像的灰度特征，以及第一投影数据，第二投影数据中除小球标记物像之外的背景特征所确定。若投影数据集A和投影数据集B图像上在坐标(u,v)处的灰度值分别为GA(u,v)和GB(u,v)，对于GA(u,v)≤T的像素，置GA(u,v)＝0，对于GA(u,v)＞T的像素，置GA(u,v)＝1；同样，对于GB(u,v)≤T的像素，置GB(u,v)＝0，对于GB(u,v)>T的像素，置GB(u,v)＝1。这样通过二值化后的图像中小球的投影是黑色的，背景是白色的，使上述投影数据集A和投影数据集B得到灰度更新，便于提取出小球的投影信息。
上述MA(uk,vk)中第一坐标参数uk、第二坐标参数vk，MB(um,vm)中第一坐标参数um、第二坐标参数vm的确定过程可以包括：
式中，符号∑表示累加，u表示投影数据集A或者投影数据集B中像素点的行数，v表示投影数据集A或者投影数据集B中像素点的列数。
作为一个实施例，上述计算所述目标函数取到最小值时的Est(u'n,v'n)，得到优化估计值，根据所述优化估计值确定CBCT系统探测器第i次采样过程中的几何参数的过程之后，还可以包括：
分别获取各次采样所确定的CBCT系统的几何参数；
根据CBCT系统的各组几何参数重建CBCT图像。
本实施例可以利用上述步骤S10至S40所述的锥形束计算机断层成像系统的几何参数确定方法获取各次采样所确定的CBCT系统的几何参数。
本实施例在计算得到各个采样角度下的几何校正参数，可以准确获得重建像素点和对应探测器探测元之间的映射关系，对待成像的物体在探测器偏置CBCT系统下进行扫描，配合CBCT重建算法，重建高精度高分辨率CBCT图像。
作为一个实施例，在根据所述优化估计值确定CBCT系统探测器第i次采样过程中的几何参数的过程之前，还可以确定根据所述MA(uk,vk)和MB(um,vm)确定投影数据集A和投影数据集B之间的对应关系，以便对相应的目标函数进行计算优化处理。
本实施例中，在得到MA(uk,vk)和MB(um,vm)后，可以依据设定偏置状态确定如图6所示的投影数据集A(图示数据集A)和投影数据集B(图示数据集B)的水平平移关系，MA(uk,vk)中uk行坐标可以作为标号索引，通过索引可以计算出数投影数据集中小球球心投影坐标MB(um,vm)的标号m，若um落在uk的邻域范围内，则可以确定投影数据集B当前标号m＝k。标号k的作用为找到校正装置第k个小球对应的投影信息。
在一个实施例中，上述Est(u'n,v'n)根据所述估计坐标计算公式确定；
所述估计坐标计算公式包括：
P3×4＝K[R|t]；
式中，s表示齐次坐标的权重因子，s通常表现为一个常数，可以根据相关坐标深度所确定，xn、yn和zn分别表示第n个金属小球在CBCT几何参数校正装置中的三维坐标参数，u′n和v′n分别表示第n个金属小球对应的估计坐标的二维坐标参数，θx、θy、θz、Tx、Ty、Tz、u0、v0、f1、f2和dt表示第i次采样对应的几何参数，其中，计算第1次采样对应的几何参数的过程中，所采用的相关初始数据可以根据CBCT系统机械装置的配置参数所确定；θx表示世界坐标系围绕成像坐标系Xw轴旋转角度、θy表示表示世界坐标系围绕成像坐标系Yw轴旋转角度、θz表示表示世界坐标系围绕成像坐标系Zw轴旋转角度、Tx表示世界坐标系偏移成像坐标系X方向距离、Ty表示表示世界坐标系偏移成像坐标系Y方向距离、Tz表示表示世界坐标系偏移成像坐标系Z方向距离、u0表示表示射线源垂直入射探测器平面垂足横坐标、v0表示表示射线源垂直入射探测器平面垂足纵坐标、f1表示表示射线源到探测器距离与探测器像素长度比例、f2表示射线源到探测器距离与探测器像素宽度比例，dt表示探测器像素畸变因子；上述K表示内参数矩阵，R表示旋转矩阵，t表示平移向量。其中，
式中，SDD为射线源到探测器的距离，du是探测器像素长度，dv是探测器像素的宽度，通常情况下，探测器像素可以为正方形，dt是探测器像素畸变因子。
本实施例可以第i-1次采样对应的几何参数确定为第i次采样几何参数确定过程中的初始几何参数，使第i次采样几何参数确定过程为依据第i-1次采样所确定的几何参数进行迭代计算的结果，有利于提高第i次采样几何参数确定的效率。若i＝1，则i-1次目标函数的初始几何参数为根据CBCT系统机械装置的配置参数确定。
在一个实施例中，上述目标函数可以为：
式中，f(θx,θy,θz,Tx,Ty,Tz,u0,v0,f1,f2,dt)表示估计坐标Est(u'n,v'n)对应的目标函数；其中目标函数的初始几何参数由i-1次采样对应的几何参数构成，i≥1；若i＝1，则i-1次目标函数的初始几何参数为根据CBCT系统机械装置的配置参数确定。
本实施例中，使上述目标函数f(θx,θy,θz,Tx,Ty,Tz,u0,v0,f1,f2,dt)取到最小值的Est(u'n,v'n)所对应的几何参数，即为相应精度最高的几何参数。
在一个实施例中，以包括17个小球标记物的CBCT几何参数校正装置为例，上述CBCT几何参数校正装置中，中空圆柱的直径为50mm，小孔和小球数量均为17个，其中16个小孔和小球的直径为1mm，1个小孔和小球的直径为2mm，作为坐标参考点，小球分布按照中心密两边疏的规则，第7个到第11个范围内的小球，两个小球之间的高度差为2.5mm；第1个到第6个范围内的小球，两个小球之间的高度差为4mm；第12个到第17个范围内的小球，两个小球之间的高度差为4mm。其中第6个到第7个小球之间的高度差为3mm，第11个到第12个小球之间的高度差为3mm。通过上述锥形束计算机断层成像系统的几何参数确定各个采样角度对应的几何参数，并根据各组几何参数构建高精度高分辨率CBCT图像的过程具体可以包括：
第一步，保持CBCT系统探测器未偏置状态，即中心射线投射至成像探测器中心位置，几何校正装置放置在成像区域范围内，通过激光定位保证校正装置中轴线放置于机架旋转轴附近，并加以固定，利用75kV/4mA曝光条件下的射线对校正装置进行圆周扫描，得到不同视角下的投影图像，采样间隔为0.9度，一共有400个角度下的投影数据集A，此时每个角度下的投影数据中包含了全部的几何校正装置中17个小球的影像信息，从上倒下依次标号为第1个小球，第2个小球，…，第17个小球。
第二步，将CBCT系统中的探测器水平移动到设定偏置状态，几何校正装置不移动仍固定于第一步中的摆放位置，再次利用同样曝光条件进行圆周扫描，得到不同视角下的投影图像，采样间隔为0.9度，一共有400个角度下的投影数据集B，此时每个角度下的投影数据中仅包含了几何校正装置中部分小球的影像信息。
第三步，根据投影图像上灰度差异，选取灰度阈值T，分别对投影数据集A和投影数据集B中不同视角下的图像通过分割算法进行二值化处理。假定投影数据集A和投影数据集B图像上在坐标(u,v)处的灰度值分别为GA(u,v)和GB(u,v)，对于GA(u,v)≤T的像素，置GA(u,v)＝0，对于GA(u,v)＞T的像素，置GA(u,v)＝1；同样，对于GB(u,v)≤T的像素，置GB(u,v)＝0，对于GB(u,v)>T的像素，置GB(u,v)＝1。这样通过二值化后的图像中小球的投影是黑色的，背景是白色的，使上述投影数据集A和投影数据集B得到灰度更新。
第四步，确定更新后(二值化后)的二值化后的中心点坐标分别为MA(uk,vk)、MB(um,vm)。坐标具体计算方法如下公式所示：
第五步，利用投影数据集A或者投影数据集B的水平平移关系，MA(uk,vk)中uk行坐标可以作为标号索引，通过索引可以计算出数投影数据集中小球球心投影坐标MB(um,vm)的标号m，若um落在uk的邻域范围内，则可以确定投影数据集B当前标号m＝k。标号k的作用为找到校正装置第k个小球对应的投影信息。
第六步，借鉴扫描相机的几何模型，根据投影数据集B中小球的中心点坐标MB(um,vm)以及小球球心空间位置(xn，yn，zn)，按照如图4所示的几何模型以及几何参数θx、θy、θz、Tx、Ty、Tz、u0、v0、f1、f2和dt来描述CBCT系统的几何映射关系，通过这个关系获取小球球心投影的估计坐标Est(u'n,v'n)与实际测量坐标MB(um,vm)之间的目标函数f(θx,θy,θz,Tx,Ty,Tz,u0,v0,f1,f2,dt)，具体计算公式如下：
其中Est(u'n,v'n)的确定过程包括：
P3×4＝K[R|t]；
式中，s表示齐次坐标的权重因子，s通常表现为一个常数，xn、yn和zn分别表示第n个金属小球在CBCT几何参数校正装置中的三维坐标参数，u′n和v′n分别表示第n个金属小球对应的估计坐标的二维坐标参数，θx、θy、θz、Tx、Ty、Tz、u0、v0、f1、f2和dt表示第i次采样对应的几何参数，其中，计算第1次采样对应的几何参数的过程中，所采用的相关初始数据可以根据CBCT系统机械装置的配置参数所确定；θx表示世界坐标系围绕成像坐标系Xw轴旋转角度、θy表示表示世界坐标系围绕成像坐标系Yw轴旋转角度、θz表示表示世界坐标系围绕成像坐标系Zw轴旋转角度、Tx表示世界坐标系偏移成像坐标系X方向距离、Ty表示表示世界坐标系偏移成像坐标系Y方向距离、Tz表示表示世界坐标系偏移成像坐标系Z方向距离、u0表示表示射线源垂直入射探测器平面垂足横坐标、v0表示表示射线源垂直入射探测器平面垂足纵坐标、f1表示表示射线源到探测器距离与探测器像素长度比例、f2表示射线源到探测器距离与探测器像素宽度比例，dt表示探测器像素畸变因子；上述K表示内参数矩阵，R表示旋转矩阵，t表示平移向量。其中，
SDD为射线源到探测器的距离，du是探测器像素长度，dv是探测器像素的宽度，通常情况下，探测器像素可以为正方形，dt是探测器像素畸变因子。
第七步，上述CBCT系统探测器尺寸为1024×1024，像素尺寸为0.124mm，射线源到探测器距离约为610mm，射线源到旋转中心距离约为380mm，探测器偏置水平移动距离约为35mm，采用(0,0,0,0,0,-230,90.0240,65.0240,610,610,0)的几何参数作为投影数据集B中第一张图像的初始值，并根据上述所示的目标函数，用CMA-ES算法计算出目标函数最小值所对应的第一张图像的最优几何参数为(-181.800,38.114,2.286,5.886,6.435,-408.341,106.687,67.333,646.843,645.875,2.751)，并以前一张图像的最优几何参数作为下一张图像优化的初始值，如此类推，一直计算到最后一张图像的几何参数。
第八步，利用计算得到各个角度下的几何校正参数，准确获得重建像素点和对应探测器探测元之间的映射关系，对待成像的物体在探测器偏置CBCT系统下进行扫描，配合CBCT重建算法，重建高精度高分辨率CBCT图像。
本实施例中，CBCT几何参数校正装置可以粗糙设计，小球间的相互位置关系可以通过螺旋CT精确获取，以降低装置加工成本；CBCT几何参数校正装置的小球分布按照中心密两边疏的规则，避免了因小球投影面积重叠而无法分割获取小球球心投影点的问题，还可以克服校正装置丢失部分小球信息后其它小球标号的计算问题；利用最优化算法减少解析算法当中矩阵分解计算所带来的误差，可以提高校正算法精度；采用单角度校正方法并借鉴相机模型，精确描述重建体素和探测器像素之间关系。
以上所述实施例的各技术特征可以进行任意的组合，为使描述简洁，未对上述实施例中的各个技术特征所有可能的组合都进行描述，然而，只要这些技术特征的组合不存在矛盾，都应当认为是本说明书记载的范围。
以上所述实施例仅表达了本发明的几种实施方式，其描述较为具体和详细，但并不能因此而理解为对发明专利范围的限制。应当指出的是，对于本领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干变形和改进，这些都属于本发明的保护范围。因此，本发明专利的保护范围应以所附权利要求为准。
Claims (10)
1.一种锥形束计算机断层成像系统的几何参数确定方法，其特征在于，包括如下步骤：
分别在CBCT系统探测器未偏置状态以及CBCT系统探测器水平移动至设定偏置状态，对CBCT几何参数校正装置中的小球标记物进行扫描，获取CBCT系统探测器未偏置状态下的第一投影数据，以及CBCT系统探测器水平移动至设定偏置状态下的第二投影数据；
在第一投影数据中采样得到投影数据集A，计算投影数据集A的中心点坐标，在第二投影数据中采样得到投影数据集B，计算投影数据集B的中心点坐标；
根据所述投影数据集A的中心点坐标、投影数据集B的中心点坐标、以及所述投影数据集A以及投影数据集B的之间的对应关系确定CBCT几何参数校正装置中各个小球标记物球心投影的估计坐标以及CBCT几何参数校正装置的目标函数；
计算所述目标函数取到最小值时的估计坐标，得到优化估计值，根据所述优化估计值确定CBCT系统在采样过程中的几何参数。
2.根据权利要求1所述的锥形束计算机断层成像系统的几何参数确定方法，其特征在于，
对CBCT几何参数校正装置中的小球标记物进行圆周扫描，在扫描过程中，采样获取CBCT系统探测器未偏置状态下的第一投影数据，以及CBCT系统探测器水平移动至设定偏置状态下的第二投影数据；
所述在第二投影数据中采样得到投影数据集B，计算投影数据集B的中心点坐标的过程包括：
在第二投影数据中获取第i次采样得到的投影数据集B，计算投影数据集B的中心点坐标MB(um,vm)；其中，um表示投影数据集B中第一坐标参数，vm表示投影数据集B中第二坐标参数，m≤N，m为第二投影数据中第i次采样所包括小球标记物像的个数，N为CBCT几何参数校正装置中的小球标记物总数；
所述计算所述目标函数取到最小值时的估计坐标，得到优化估计值，根据所述优化估计值确定CBCT系统在采样过程中的几何参数的过程包括：
计算所述目标函数取到最小值时的Est(u'n,v'n)，得到优化估计值，根据所述优化估计值确定CBCT系统第i次采样过程中的几何参数；其中，Est(u′n,v′n)为各个小球标记物球心投影的估计坐标，u′n表示第n个金属小球对应的估计坐标的行坐标，v′n表示第n个金属小球对应的估计坐标的列坐标，n＝1,2，…，N，N为CBCT几何参数校正装置中的小球标记物总数。
3.根据权利要求1所述的锥形束计算机断层成像系统的几何参数确定方法，其特征在于，所述分别在CBCT系统探测器未偏置状态以及CBCT系统探测器水平移动至设定偏置状态，对CBCT几何参数校正装置中的小球标记物进行扫描，获取CBCT系统探测器未偏置状态下的第一投影数据，以及CBCT系统探测器水平移动至设定偏置状态下的第二投影数据的过程之前，还包括：
将CBCT几何参数校正装置的圆柱轴线和CBCT系统的旋转轴重合，在所述CBCT系统中设置CBCT几何参数校正装置。
4.根据权利要求1所述的锥形束计算机断层成像系统的几何参数确定方法，其特征在于，所述CBCT几何参数校正装置包括空心圆柱以及小球标记物；所述空心圆柱的外表面上设计沿着中心轴方向螺旋轨迹的小孔，小球标记物镶嵌在小孔内，小球数量大于或等于6个。
5.根据权利要求4所述的锥形束计算机断层成像系统的几何参数确定方法，其特征在于，所述空心圆柱的制作材料包括低密度物质；所述小球标记物的制作材料包括高密度物质。
6.根据权利要求5所述的锥形束计算机断层成像系统的几何参数确定方法，其特征在于，所述空心圆柱的制作材料为有机玻璃，所述小球标记物的制作材料为钢。
7.根据权利要求1至6任一项所述的锥形束计算机断层成像系统的几何参数确定方法，其特征在于，所述在第一投影数据中采样得到的投影数据集A，在第二投影数据中采样得到的投影数据集B的过程之后，还包括：
在第一投影数据中获取第i次采样得到的投影数据集A，识别投影数据集A中各个像素点的灰度值，将投影数据集A中灰度值大于灰度阈值T的像素点灰度置1，将投影数据集A中灰度值小于或者等于灰度阈值T的像素点灰度置0；
识别投影数据集B中各个像素点的灰度值，将投影数据集B中灰度值大于灰度阈值T的像素点灰度置1，将投影数据集B中灰度值小于或者等于灰度阈值T的像素点灰度置0；
根据灰度更新后的投影数据集A计算MA(uk,vk)，根据灰度更新后的投影数据集B算MB(um,vm)；其中，MA(uk，vk)为投影数据集A中心点坐标，uk表示投影数据集A中第一坐标参数，vk表示投影数据集A中第二坐标参数，k≤N，k为第一投影数据中第i次采样所包括小球标记物像的个数；MB(um，vm)为投影数据集B的中心点坐标，um表示投影数据集B中第一坐标参数，vm表示投影数据集B中第二坐标参数，m≤N，m为第二投影数据中第i次采样所包括小球标记物像的个数，N为CBCT几何参数校正装置中的小球标记物总数。
8.根据权利要求2所述的锥形束计算机断层成像系统的几何参数确定方法，其特征在于，所述计算所述目标函数取到最小值时的Est(u'n,v'n)，得到优化估计值，根据所述优化估计值确定CBCT系统第i次采样过程中的几何参数的过程之后，还包括：
分别获取各次采样所确定的CBCT系统的几何参数；
根据CBCT系统的各组几何参数重建CBCT图像。
9.根据权利要求2所述的锥形束计算机断层成像系统的几何参数确定方法，其特征在于，所述Est(u'n,v'n)根据估计坐标计算公式确定；
所述估计坐标计算公式包括：
<mrow>
<mi>s</mi>
<mfenced open = "[" close = "]">
<mtable>
<mtr>
<mtd>
<msubsup>
<mi>u</mi>
<mi>n</mi>
<mo>&prime;</mo>
</msubsup>
</mtd>
</mtr>
<mtr>
<mtd>
<msubsup>
<mi>v</mi>
<mi>n</mi>
<mo>&prime;</mo>
</msubsup>
</mtd>
</mtr>
<mtr>
<mtd>
<mn>1</mn>
</mtd>
</mtr>
</mtable>
</mfenced>
<mo>=</mo>
<msub>
<mi>P</mi>
<mrow>
<mn>3</mn>
<mo>&times;</mo>
<mn>4</mn>
</mrow>
</msub>
<mfenced open = "[" close = "]">
<mtable>
<mtr>
<mtd>
<msub>
<mi>x</mi>
<mi>n</mi>
</msub>
</mtd>
</mtr>
<mtr>
<mtd>
<msub>
<mi>y</mi>
<mi>n</mi>
</msub>
</mtd>
</mtr>
<mtr>
<mtd>
<msub>
<mi>z</mi>
<mi>n</mi>
</msub>
</mtd>
</mtr>
<mtr>
<mtd>
<mn>1</mn>
</mtd>
</mtr>
</mtable>
</mfenced>
<mo>;</mo>
</mrow>
P3×4＝K[R|t]；
<mrow>
<mi>K</mi>
<mo>=</mo>
<mfenced open = "[" close = "]">
<mtable>
<mtr>
<mtd>
<mrow>
<mi>f</mi>
<mn>1</mn>
</mrow>
</mtd>
<mtd>
<mrow>
<mi>d</mi>
<mi>t</mi>
</mrow>
</mtd>
<mtd>
<msub>
<mi>u</mi>
<mn>0</mn>
</msub>
</mtd>
</mtr>
<mtr>
<mtd>
<mn>0</mn>
</mtd>
<mtd>
<mrow>
<mi>f</mi>
<mn>2</mn>
</mrow>
</mtd>
<mtd>
<msub>
<mi>v</mi>
<mn>0</mn>
</msub>
</mtd>
</mtr>
<mtr>
<mtd>
<mn>0</mn>
</mtd>
<mtd>
<mn>0</mn>
</mtd>
<mtd>
<mn>1</mn>
</mtd>
</mtr>
</mtable>
</mfenced>
<mo>;</mo>
</mrow>
<mrow>
<mi>R</mi>
<mo>=</mo>
<mfenced open = "[" close = "]">
<mtable>
<mtr>
<mtd>
<mrow>
<msub>
<mi>cos&theta;</mi>
<mi>z</mi>
</msub>
</mrow>
</mtd>
<mtd>
<mrow>
<mo>-</mo>
<msub>
<mi>sin&theta;</mi>
<mi>z</mi>
</msub>
</mrow>
</mtd>
<mtd>
<mn>0</mn>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msub>
<mi>sin&theta;</mi>
<mi>z</mi>
</msub>
</mrow>
</mtd>
<mtd>
<mrow>
<msub>
<mi>cos&theta;</mi>
<mi>z</mi>
</msub>
</mrow>
</mtd>
<mtd>
<mn>0</mn>
</mtd>
</mtr>
<mtr>
<mtd>
<mn>0</mn>
</mtd>
<mtd>
<mn>0</mn>
</mtd>
<mtd>
<mn>1</mn>
</mtd>
</mtr>
</mtable>
</mfenced>
<mfenced open = "[" close = "]">
<mtable>
<mtr>
<mtd>
<mrow>
<msub>
<mi>cos&theta;</mi>
<mi>y</mi>
</msub>
</mrow>
</mtd>
<mtd>
<mn>0</mn>
</mtd>
<mtd>
<mrow>
<msub>
<mi>sin&theta;</mi>
<mi>y</mi>
</msub>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mn>0</mn>
</mtd>
<mtd>
<mn>1</mn>
</mtd>
<mtd>
<mn>0</mn>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>-</mo>
<msub>
<mi>sin&theta;</mi>
<mi>y</mi>
</msub>
</mrow>
</mtd>
<mtd>
<mn>0</mn>
</mtd>
<mtd>
<mrow>
<msub>
<mi>cos&theta;</mi>
<mi>y</mi>
</msub>
</mrow>
</mtd>
</mtr>
</mtable>
</mfenced>
<mfenced open = "[" close = "]">
<mtable>
<mtr>
<mtd>
<mn>1</mn>
</mtd>
<mtd>
<mn>0</mn>
</mtd>
<mtd>
<mn>0</mn>
</mtd>
</mtr>
<mtr>
<mtd>
<mn>0</mn>
</mtd>
<mtd>
<mrow>
<msub>
<mi>cos&theta;</mi>
<mi>x</mi>
</msub>
</mrow>
</mtd>
<mtd>
<mrow>
<mo>-</mo>
<msub>
<mi>sin&theta;</mi>
<mi>x</mi>
</msub>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mn>0</mn>
</mtd>
<mtd>
<mrow>
<msub>
<mi>sin&theta;</mi>
<mi>x</mi>
</msub>
</mrow>
</mtd>
<mtd>
<mrow>
<msub>
<mi>cos&theta;</mi>
<mi>x</mi>
</msub>
</mrow>
</mtd>
</mtr>
</mtable>
</mfenced>
<mo>;</mo>
</mrow>
<mrow>
<mi>t</mi>
<mo>=</mo>
<mfenced open = "[" close = "]">
<mtable>
<mtr>
<mtd>
<msub>
<mi>T</mi>
<mi>x</mi>
</msub>
</mtd>
</mtr>
<mtr>
<mtd>
<msub>
<mi>T</mi>
<mi>y</mi>
</msub>
</mtd>
</mtr>
<mtr>
<mtd>
<msub>
<mi>T</mi>
<mi>z</mi>
</msub>
</mtd>
</mtr>
</mtable>
</mfenced>
<mo>;</mo>
</mrow>
式中，s表示齐次坐标的权重因子，xn、yn和zn分别表示第n个金属小球在CBCT几何参数校正装置中的三维坐标参数，u′n和v′n分别表示第n个金属小球对应的估计坐标的二维坐标参数，θx、θy、θz、Tx、Ty、Tz、u0、v0、f1、f2和dt表示第i次采样对应的几何参数。
10.根据权利要求9所述的锥形束计算机断层成像系统的几何参数确定方法，其特征在于，所述目标函数为：
<mrow>
<mtable>
<mtr>
<mtd>
<mrow>
<mi>f</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>&theta;</mi>
<mi>x</mi>
</msub>
<mo>,</mo>
<msub>
<mi>&theta;</mi>
<mi>y</mi>
</msub>
<mo>,</mo>
<msub>
<mi>&theta;</mi>
<mi>z</mi>
</msub>
<mo>,</mo>
<msub>
<mi>T</mi>
<mi>x</mi>
</msub>
<mo>,</mo>
<msub>
<mi>T</mi>
<mi>y</mi>
</msub>
<mo>,</mo>
<msub>
<mi>T</mi>
<mi>z</mi>
</msub>
<mo>,</mo>
<msub>
<mi>u</mi>
<mn>0</mn>
</msub>
<mo>,</mo>
<msub>
<mi>v</mi>
<mn>0</mn>
</msub>
<mo>,</mo>
<mi>f</mi>
<mn>1</mn>
<mo>,</mo>
<mi>f</mi>
<mn>2</mn>
<mo>,</mo>
<mi>d</mi>
<mi>t</mi>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>=</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>n</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>N</mi>
</munderover>
<msub>
<mrow>
<mo>(</mo>
<mi>E</mi>
<mi>s</mi>
<mi>t</mi>
<mo>(</mo>
<mrow>
<msubsup>
<mi>u</mi>
<mi>n</mi>
<mo>&prime;</mo>
</msubsup>
<mo>,</mo>
<msubsup>
<mi>v</mi>
<mi>n</mi>
<mo>&prime;</mo>
</msubsup>
</mrow>
<mo>)</mo>
<mo>-</mo>
<mi>M</mi>
<mi>B</mi>
<mo>(</mo>
<mrow>
<msub>
<mi>u</mi>
<mi>m</mi>
</msub>
<mo>,</mo>
<msub>
<mi>v</mi>
<mi>m</mi>
</msub>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msub>
</mrow>
</mtd>
</mtr>
</mtable>
<mo>;</mo>
</mrow>
式中，f(θx,θy,θz,Tx,Ty,Tz,u0,v0,f1,f2,dt)表示估计坐标Est(u'n,v'n)对应的目标函数。
Priority Applications (2)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710393440.1A CN106994023B (zh) | 2017-05-27 | 2017-05-27 | 锥形束计算机断层成像系统的几何参数确定方法 |
| PCT/CN2017/086851 WO2018218611A1 (zh) | 2017-05-27 | 2017-06-01 | 锥形束计算机断层成像系统的几何参数确定方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710393440.1A CN106994023B (zh) | 2017-05-27 | 2017-05-27 | 锥形束计算机断层成像系统的几何参数确定方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN106994023A CN106994023A (zh) | 2017-08-01 |
| CN106994023B true CN106994023B (zh) | 2018-02-23 |
Family
ID=59436457
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201710393440.1A Expired - Fee Related CN106994023B (zh) | 2017-05-27 | 2017-05-27 | 锥形束计算机断层成像系统的几何参数确定方法 |
Country Status (2)
| Country | Link |
|---|---|
| CN (1) | CN106994023B (zh) |
| WO (1) | WO2018218611A1 (zh) |
Families Citing this family (13)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110850859B (zh) * | 2018-08-01 | 2023-03-07 | 深圳市优必选科技有限公司 | 一种机器人及其避障方法和避障系统 |
| CN109452950A (zh) * | 2018-12-29 | 2019-03-12 | 上海联影医疗科技有限公司 | 一种模体及成像设备的几何标记方法 |
| CN109756717B (zh) * | 2019-03-06 | 2023-02-24 | 成都工业学院 | 一种高分辨率的投影显示装置 |
| CN110084855B (zh) * | 2019-04-19 | 2020-12-15 | 合肥中科离子医学技术装备有限公司 | 一种改进cbct几何参数标定方法 |
| CN110517330B (zh) * | 2019-08-07 | 2021-05-28 | 北京航空航天大学 | 一种偏置扫描模式下的工业锥束ct重建方法 |
| CN110766629B (zh) * | 2019-10-17 | 2022-03-01 | 广州华端科技有限公司 | Cbct系统几何校正方法、装置、计算机设备和存储介质 |
| CN111458739B (zh) * | 2020-03-16 | 2021-11-30 | 北京东软医疗设备有限公司 | 成像方法、装置及系统 |
| CN111513740B (zh) * | 2020-04-13 | 2023-09-12 | 北京东软医疗设备有限公司 | 血管造影机的控制方法、装置、电子设备、存储介质 |
| CN112634353B (zh) * | 2020-12-17 | 2024-03-26 | 广州华端科技有限公司 | Cbct系统几何标定模体自标定方法、装置及介质 |
| CN112529820B (zh) * | 2020-12-29 | 2024-06-14 | 广州开云影像科技有限公司 | Cbct系统半扇束几何校正方法、装置、设备和介质 |
| CN113963058B (zh) * | 2021-09-07 | 2022-11-29 | 于留青 | 预定轨道ct在线标定方法、装置、电子设备及存储介质 |
| CN113749681A (zh) * | 2021-11-10 | 2021-12-07 | 极限人工智能(北京)有限公司 | 评估cbct机图像精度的工装 |
| CN113782174B (zh) * | 2021-11-12 | 2022-02-18 | 极限人工智能有限公司 | Ct机的校验方法、装置、电子设备及可读存储介质 |
Family Cites Families (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US6848827B2 (en) * | 2003-05-13 | 2005-02-01 | General Electric Company | Method and apparatus for calibrating detector spectral response |
| US7613275B2 (en) * | 2005-12-19 | 2009-11-03 | General Electric Company | Method and apparatus for reducing cone beam artifacts using spatially varying weighting functions |
| CN100382763C (zh) * | 2006-03-31 | 2008-04-23 | 北京航空航天大学 | 一种适用于三维ct扫描系统投影坐标原点的标定方法 |
| CN101750430B (zh) * | 2009-06-10 | 2011-10-12 | 中国科学院自动化研究所 | X射线计算机断层成像系统的几何校正方法 |
| CN102598059B (zh) * | 2009-08-06 | 2015-07-01 | 皇家飞利浦电子股份有限公司 | 利用偏移探测器几何结构生成计算断层摄影图像的方法和设备 |
| CN102539460B (zh) * | 2012-01-06 | 2014-07-02 | 公安部第一研究所 | 一种ct系统投影旋转中心定位方法 |
| CN103226113B (zh) * | 2013-03-29 | 2015-10-21 | 中国计量科学研究院 | 锥束3d-ct扫描系统重建体素尺寸的自动标定方法 |
| US9628723B2 (en) * | 2015-02-13 | 2017-04-18 | Toshiba Medical Systems Corporation | Computed tomography scanner calibration with angle correction for scan angle offset |
-
2017
- 2017-05-27 CN CN201710393440.1A patent/CN106994023B/zh not_active Expired - Fee Related
- 2017-06-01 WO PCT/CN2017/086851 patent/WO2018218611A1/zh active Application Filing
Also Published As
| Publication number | Publication date |
|---|---|
| WO2018218611A1 (zh) | 2018-12-06 |
| CN106994023A (zh) | 2017-08-01 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN106994023B (zh) | 锥形束计算机断层成像系统的几何参数确定方法 | |
| CN108122203B (zh) | 一种几何参数的校正方法、装置、设备及系统 | |
| CN105931202B (zh) | 几何校正模体的校正方法和系统 | |
| CN111248934B (zh) | 一种用于cbct系统机械校正的方法及系统 | |
| CN106667512B (zh) | X射线成像设备的几何校正方法、乳腺断层成像设备 | |
| US7844094B2 (en) | Systems and methods for determining geometric parameters of imaging devices | |
| CN105849537B (zh) | 校准设备和计算机断层扫描方法 | |
| CN104107065B (zh) | 3d图像集在不同空间之间的最佳变换 | |
| CN112451093B (zh) | 图像引导机器人微创外科手术的物理空间与图像空间配准方法 | |
| EP3153101B1 (en) | Identification and registration of multi-marker jig | |
| CN109671128B (zh) | 图像重建过程中的数据处理、图像重建方法和装置 | |
| US20200085404A1 (en) | Geometric calibration for cone beam ct using line fiducials | |
| WO2018126335A1 (zh) | 基于小球模体的锥束ct系统几何参数评价及校正方法 | |
| US20210372950A1 (en) | CT Scanner Calibration | |
| CN102488528B (zh) | 一种层析成像几何参数的校准方法 | |
| CN109658465B (zh) | 图像重建过程中的数据处理、图像重建方法和装置 | |
| CN114596222A (zh) | 适用于一般轨迹锥束ct系统几何校正的模体与标定方法 | |
| CN107016655A (zh) | 锥束cl几何全参数迭代校正方法 | |
| Duan et al. | Knowledge-based self-calibration method of calibration phantom by and for accurate robot-based CT imaging systems | |
| US20220130081A1 (en) | Computer-implemented method for determining at least one geometric parameter required for evaluating measurement data | |
| JP2021043046A (ja) | Ｘ線計測装置の校正方法 | |
| Schreiner et al. | Accuracy assessment of a clinical biplane fluoroscope for three-dimensional measurements and targeting | |
| CN113963057B (zh) | 成像几何关系标定方法、装置、电子设备以及存储介质 | |
| CN113729747B (zh) | 一种球形金属标记的锥束ct金属伪影去除系统及去除方法 | |
| Schiffers et al. | Disassemblable fieldwork CT scanner using a 3D-printed calibration phantom |
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 | ||
| CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20180223 |