CN112233155A - 一种2d-3d影像配准算法 - Google Patents

一种2d-3d影像配准算法 Download PDF

Info

Publication number
CN112233155A
CN112233155A CN202011415198.1A CN202011415198A CN112233155A CN 112233155 A CN112233155 A CN 112233155A CN 202011415198 A CN202011415198 A CN 202011415198A CN 112233155 A CN112233155 A CN 112233155A
Authority
CN
China
Prior art keywords
image
preoperative
matching
gradient
pose
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.)
Granted
Application number
CN202011415198.1A
Other languages
English (en)
Other versions
CN112233155B (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.)
Tuodao Medical Technology Co Ltd
Original Assignee
Nanjing Tuodao Medical Technology Co Ltd
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 Nanjing Tuodao Medical Technology Co Ltd filed Critical Nanjing Tuodao Medical Technology Co Ltd
Priority to CN202011415198.1A priority Critical patent/CN112233155B/zh
Publication of CN112233155A publication Critical patent/CN112233155A/zh
Application granted granted Critical
Publication of CN112233155B publication Critical patent/CN112233155B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/30Determination of transform parameters for the alignment of images, i.e. image registration
    • G06T7/33Determination of transform parameters for the alignment of images, i.e. image registration using feature-based methods
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10072Tomographic images
    • G06T2207/10081Computed x-ray tomography [CT]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10116X-ray image
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/30Subject of image; Context of image processing
    • G06T2207/30004Biomedical image processing

Landscapes

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

Abstract

本发明公开了一种2D‑3D影像配准算法,包括步骤:(1)在术前扫描患者患处得到术前3D影像;在术中分别对患者相对应患处进行正侧位透视扫描得到正侧位2D透视影像;(2)构建术中拍摄正侧位透视影像设备的空间模型,并根据步骤(1)计算得到术前3D影像在空间模型中的初始位姿;(3)采用步骤(2)构建的空间模型生成术前3D影像的正侧位DRR影像,并将其与其对应的术中2D透视影像进行模板匹配,得到术前3D影像的粗匹配变换位姿;(4)使用梯度下降优化算法进行精准优化匹配,得到术前3D影像的最终变换位姿。本发明通过采用模板匹配在进行粗匹配,然后使用梯度下降优化算法得到最优结果,可以在保证配准速度的同时,得到非常精准的配准结果。

Description

一种2D-3D影像配准算法
技术领域
本发明涉及手术导航技术领域,尤其涉及一种3D CT影像和2D透视影像的配准算法。
背景技术
3D CT影像和2D透视影像配准算法计算3D CT影像和2D透视影像之间的位姿变换值,在影像导航置入手术中被用于将术前信息(基于CT的手术计划数据)整合到术中X射线投影影像中,以提高到达目标位置的精度。当前的算法实现通常使用3D CT影像生成2D DRR(数字重建放射)影像,计算DRR影像和透视影像之间的相似度,找到CT的某个位姿使生成的DRR影像和透视影像具有最佳的相似度。由于目标变量有6个自由度(其中3个表示旋转,另外3个表示平移),传统的基于模板匹配的相似度评价的算法或基于随机采样搜索的算法并不能得到精准的结果,且非常耗时,影响算法在实际场景的应用。
发明内容
发明目的:本发明针对上述不足,提出了一种快速、准确的2D-3D影像配准算法。
技术方案:
一种2D-3D影像配准算法,包括步骤:
(1)在术前扫描患者患处得到术前3D影像;在术中分别对患者相对应患处进行正侧位透视扫描得到患者患处的正侧位2D透视影像;
(2)构建术中拍摄正侧位透视影像设备的空间模型,并根据步骤(1)计算得到术前3D影像在空间模型中的初始位姿;
(3)采用步骤(2)构建的空间模型生成术前3D影像的正侧位DRR影像,并将其与其对应的2D透视影像进行模板匹配,得到术前3D影像的粗匹配变换位姿;
(4)使用梯度下降优化算法进行精准优化匹配,得到术前3D影像的最终变换位姿。
所述步骤(3)中,位姿包括三维欧拉角和三维平移向量;具体计算得到术前3D影像的粗匹配变换位姿如下:
通过模板匹配计算生成的DRR影像和其对应的2D透视影像的归一化相关匹配值,得到最佳的匹配位置,据此计算得到术前3D影像粗匹配后的平移向量;
对于初始欧拉角,分别以生成的DRR影像的每个分量角度以预设的偏差取值,并分别计算其与对应的2D透视影像之间的归一化相关匹配值,得到最高匹配值对应的欧拉角;
综上计算得到术前3D影像的粗匹配变换位姿。
所述步骤(4)中:
以术前3D影像的粗匹配变换位姿作为待优化变量;
分别计算正侧位DRR影像的梯度图及经过粗匹配变换位姿变换后的3D影像生成的对应DRR影像的梯度图,并以二者之间的相似度值作为评价函数;
使用梯度下降优化算法计算术前3D影像的最终变换位姿。
使用Scharr算子分别计算所述DRR影像的梯度图和术中2D透视影像的梯度图。
所述DRR影像的梯度图和2D透视影像的梯度图之间的相似度值metric计算如下:
正位DRR影像和正位2D透视影像的单椎骨梯度相似度值
Figure 645788DEST_PATH_IMAGE002
计算方法为:
Figure 850504DEST_PATH_IMAGE004
其中,
Figure 987088DEST_PATH_IMAGE006
表示正位2D透视影像的梯度图,
Figure 593650DEST_PATH_IMAGE008
表示正位DRR影像的梯度图,
Figure 253258DEST_PATH_IMAGE010
是正位透视影像的方差值;侧位相似度
Figure 210850DEST_PATH_IMAGE012
的计算方法类似,最后计算总的相似度值
Figure 885545DEST_PATH_IMAGE014
使用梯度下降优化算法计算术前3D影像的最终变换位姿具体为:
(41)计算当前值T的梯度,梯度由每个分量T[i]的导数组成;
(42)设定偏差Δ,则T[i]的导数为:
Figure 346613DEST_PATH_IMAGE016
(43)沿梯度方向优化当前值T,计算优化后的变量T[i]'为:
Figure 979720DEST_PATH_IMAGE018
其中
Figure 159028DEST_PATH_IMAGE020
为梯度因子,gradientMagnitude为梯度dT的模长,stepLength为优化步长;
(44)计算震荡系数shakeFactor
Figure 371835DEST_PATH_IMAGE022
dPT[i]为上次迭代计算得到的PT每一分量的导数,初始值设为0;若shakeFactor<0,则优化步长stepLength调整为之前的一半,并重新确定偏差,重复步骤(42);否则转至步骤(45);
(45)判断优化步长stepLength是否满足终止条件,若满足则停止迭代,输出最终的T;若不满足,则重新确定偏差Δ,返回步骤(42)。
所述步骤(45)中的终止条件为stepLength<0.03
有益效果:本发明的算法首先采用模板匹配在一定范围内进行粗匹配,然后使用梯度下降优化算法得到最优结果,可以在保证配准速度的同时,得到非常精准的配准结果。
附图说明
图1为本发明的2D-3D影像配准算法框架流程图。
图2为术前CT影像中模拟X光机生成DRR影像的原理图。
图3为在图像x和y方向上进行模板匹配计算的归一化相关匹配值。
具体实施方式
下面结合附图和具体实施例,进一步阐明本发明。
本发明主要应用在人体椎骨的CT透视影像配准中,图1为本发明的流程图,如图1所示,本发明的2D-3D影像配准算法包括步骤:
(1)通过C臂机或CT机在术前扫描患者患处(本实施例为脊椎)得到术前3D影像;对术前3D影像进行自动分割得到椎骨图像,进而得到椎骨的中心点位置(即椎骨中心点在影像坐标系中的坐标)和包围盒(即椎骨的外接长方体);
(2)通过C臂机或CT机在术中分别对患者患处进行正侧位透视扫描,得到患者患处的正侧位术中2D透视影像,并使用透视变换矩阵得到单个椎骨影像,在术中2D透视影像中的椎骨中心点加上标记;
(3)构建术中拍摄正侧位透视影像的C臂机或CT机的空间模型(以X射线发射管中心作为空间模型光心,C臂机成像平面作为空间模型成像平面),根据构建的空间模型中成像平面位姿及步骤(2)标记的术中2D透视影像的椎骨中心点的坐标可以计算出术中2D透视影像中的椎骨中心点在构建的空间模型中的位置坐标;
(4)根据步骤(1)得到的椎骨中心点在影像坐标系中的坐标,可以计算出3D影像在C臂机或CT机空间模型中的初始位姿T 0[r 0 ,t 0],其中r 0是表示旋转的3维欧拉角,t 0是表示平移的3维平移向量;
(5)采用步骤(3)构建的C臂机或CT机空间模型生成术前3D影像正侧位2D DRR影像,生成方法如图2所示:
从X射线源向2D DRR成像平面每个像素的方向采样,由以下公式计算最终的像素值I xy
Figure 952989DEST_PATH_IMAGE024
式中i是术前3D影像的采样值,作为衰减系数。生成算法需要大量计算,并且是可以并行计算的,所以使用CUDA加速计算;
(6)利用模板匹配进行粗配准:
通过模板匹配计算生成的2D DRR影像和其对应的术中2D透视影像的归一化相关匹配值,结果如图3 所示,最亮的像素点对应最佳的匹配位置,这样就可以得到术前3D影像粗匹配后的平移向t 1;对于初始欧拉角r 1,分别以生成的2D DRR影像的每个分量角度以预设的偏差取值,并分别计算其与对应的术中2D透视影像之间的归一化相关匹配值,得到最高匹配值对应的欧拉角r 1,最终得到粗匹配后的术前3D影像的变换位姿T 1[r 1 ,t 1];
(7)粗匹配可以将平移向量的误差控制在0.5毫米以内(透视影像一个像素对应的长度),下面需要使用梯度下降优化算法进行精准优化匹配,具体优化步骤为:
(71)将粗匹配后的变换位姿T 1[r 1 ,t 1]传递给梯度下降优化算法的初始值TT就是优化算法的待优化变量;
(72)由当前的变换位姿T得到变换后的3D影像,并采用步骤(5)生成DRR影像,并使用透视变换矩阵得到单个椎骨影像(脊柱是非刚体,所以算法进行单椎骨配准),使用Scharr算子分别计算该DRR影像的梯度图和透视影像的梯度图,再计算两幅梯度图的相似度值metric
以正位透视为例,正位DRR影像和正位2D透视影像的单椎骨梯度图相似度值
Figure DEST_PATH_IMAGE025
计算方法为:
Figure 429101DEST_PATH_IMAGE026
其中,
Figure DEST_PATH_IMAGE027
表示正位2D透视影像的梯度图,
Figure 33388DEST_PATH_IMAGE008
表示正位DRR影像的梯度图,
Figure 315465DEST_PATH_IMAGE028
是正位透视影像的方差值。该方法计算每个对应像素的相似度
Figure 751126DEST_PATH_IMAGE030
,然后将所有像素的相似度值累加。显然,透视图和DRR图越相似,
Figure 460456DEST_PATH_IMAGE032
越接近0,像素相似度值越接近1,的值越接近梯度的像素数,它们是正相关的;侧位相似度
Figure DEST_PATH_IMAGE033
的计算方法类似,最后计算总的相似度值
Figure 286460DEST_PATH_IMAGE034
本算法框架在计算DRR、梯度图和相似度值都使用CUDA完成,这样不仅可以利用GPU加速计算,还可以减少数据在主机内存和设备内存(GPU显存)间传递所花费的时间,最终只需要将相似度值从显存中拷贝至主机内存。
(73)计算当前6自由度T的梯度,梯度由每个分量T[i]的导数组成;具体为:定义偏差Δ,则T[i]的导数为:
Figure 835210DEST_PATH_IMAGE016
(74)调整优化步长stepLength,具体为计算震荡系数shakeFactor
Figure DEST_PATH_IMAGE035
dPT[i]为上次迭代计算得到的T(定义为PT)每一分量的导数,初始值设为0。如果shakeFactor<0,则stepLength调整为之前的一半,并重新确定偏差,重复步骤(73);否则转至步骤(75);
(75)沿梯度方向优化当前值T,计算优化后的变量T[i]'为:
Figure DEST_PATH_IMAGE036
其中
Figure DEST_PATH_IMAGE037
为梯度因子,gradientMagnitude为梯度dT的模长,stepLength为优化步长;
(76)判断优化步长stepLength是否满足终止条件(即stepLength<0.03),若满足则停止迭代,输出最终的变换位姿T;若不满足,则重新确定偏差Δ,返回步骤(73)。
以上详细描述了本发明的优选实施方式,但是本发明并不限于上述实施方式中的具体细节,在本发明的技术构思范围内,可以对本发明的技术方案进行多种等同变换(如数量、形状、位置等),这些等同变换均属于本发明的保护范围。

Claims (7)

1.一种2D-3D影像配准算法,其特征在于:包括步骤:
(1)在术前扫描患者患处得到术前3D影像;在术中分别对患者相对应患处进行正侧位透视扫描得到患者患处的正侧位2D透视影像;
(2)构建术中拍摄正侧位透视影像设备的空间模型,并根据步骤(1)计算得到术前3D影像在空间模型中的初始位姿;
(3)采用步骤(2)构建的空间模型生成术前3D影像的正侧位DRR影像,并将其与其对应的2D透视影像进行模板匹配,得到术前3D影像的粗匹配变换位姿;
(4)使用梯度下降优化算法进行精准优化匹配,得到术前3D影像的最终变换位姿。
2.根据权利要求1所述的2D-3D影像配准算法,其特征在于:所述步骤(3)中,位姿包括三维欧拉角和三维平移向量;具体计算得到术前3D影像的粗匹配变换位姿如下:
通过模板匹配计算生成的DRR影像和其对应的2D透视影像的归一化相关匹配值,得到最佳的匹配位置,据此计算得到术前3D影像粗匹配后的平移向量;
对于初始欧拉角,分别以生成的DRR影像的每个分量角度以预设的偏差取值,并分别计算其与对应的2D透视影像之间的归一化相关匹配值,得到最高匹配值对应的欧拉角;
综上计算得到术前3D影像的粗匹配变换位姿。
3.根据权利要求1所述的2D-3D影像配准算法,其特征在于:所述步骤(4)中:
以术前3D影像的粗匹配变换位姿作为待优化变量;
分别计算正侧位DRR影像的梯度图及经过粗匹配变换位姿变换后的3D影像生成的对应DRR影像的梯度图,并以二者之间的相似度值作为评价函数;
使用梯度下降优化算法计算术前3D影像的最终变换位姿。
4.根据权利要求3所述的2D-3D影像配准算法,其特征在于:使用Scharr算子分别计算所述DRR影像的梯度图和2D透视影像的梯度图。
5.根据权利要求3所述的2D-3D影像配准算法,其特征在于:所述DRR影像的梯度图和2D透视影像的梯度图之间的相似度值metric计算如下:
正位DRR影像和正位2D透视影像的单椎骨梯度相似度值
Figure 629663DEST_PATH_IMAGE002
计算方法为:
Figure 164550DEST_PATH_IMAGE004
其中,
Figure 543579DEST_PATH_IMAGE006
表示正位2D透视影像的梯度图,
Figure 531126DEST_PATH_IMAGE008
表示正位DRR影像的梯度图,
Figure 817751DEST_PATH_IMAGE010
是正位透视影像的方差值;侧位相似度
Figure 207144DEST_PATH_IMAGE012
的计算方法类似,最后计算总的相似度值
Figure 553812DEST_PATH_IMAGE014
6.根据权利要求3所述的2D-3D影像配准算法,其特征在于:使用梯度下降优化算法计算术前3D影像的最终变换位姿具体为:
(41)计算当前值T的梯度,梯度由每个分量T[i]的导数组成;
(42)设定偏差Δ,则T[i]的导数为:
Figure 28655DEST_PATH_IMAGE016
(43)沿梯度方向优化当前值T,计算优化后的变量T[i]'为:
Figure 322233DEST_PATH_IMAGE018
其中
Figure 831712DEST_PATH_IMAGE020
为梯度因子,gradientMagnitude为梯度dT的模长,stepLength为优化步长;
(44)计算震荡系数shakeFactor
Figure 286964DEST_PATH_IMAGE022
dPT[i]为上次迭代计算得到的PT每一分量的导数,初始值设为0;若shakeFactor<0,则优化步长stepLength调整为之前的一半,并重新确定偏差,重复步骤(42);否则转至步骤(45);
(45)判断优化步长stepLength是否满足终止条件,若满足则停止迭代,输出最终的T;若不满足,则重新确定偏差Δ,返回步骤(42)。
7.根据权利要求6所述的2D-3D影像配准算法,其特征在于:所述步骤(45)中的终止条件为stepLength<0.03
CN202011415198.1A 2020-12-07 2020-12-07 一种2d-3d影像配准算法 Active CN112233155B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011415198.1A CN112233155B (zh) 2020-12-07 2020-12-07 一种2d-3d影像配准算法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011415198.1A CN112233155B (zh) 2020-12-07 2020-12-07 一种2d-3d影像配准算法

Publications (2)

Publication Number Publication Date
CN112233155A true CN112233155A (zh) 2021-01-15
CN112233155B CN112233155B (zh) 2021-02-26

Family

ID=74123566

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011415198.1A Active CN112233155B (zh) 2020-12-07 2020-12-07 一种2d-3d影像配准算法

Country Status (1)

Country Link
CN (1) CN112233155B (zh)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113040908A (zh) * 2021-02-02 2021-06-29 武汉联影智融医疗科技有限公司 手术导航的配准方法、装置、计算机设备和存储介质
CN115082534A (zh) * 2022-07-21 2022-09-20 杭州三坛医疗科技有限公司 双平面图像配准方法、装置及机器人
WO2023006021A1 (zh) * 2021-07-30 2023-02-02 武汉联影智融医疗科技有限公司 一种配准方法和系统

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20130094745A1 (en) * 2011-09-28 2013-04-18 Siemens Corporation Non-rigid 2d/3d registration of coronary artery models with live fluoroscopy images
CN111429491A (zh) * 2020-03-11 2020-07-17 上海嘉奥信息科技发展有限公司 脊柱术前三维影像与术中二维影像配准的方法与系统

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20130094745A1 (en) * 2011-09-28 2013-04-18 Siemens Corporation Non-rigid 2d/3d registration of coronary artery models with live fluoroscopy images
CN111429491A (zh) * 2020-03-11 2020-07-17 上海嘉奥信息科技发展有限公司 脊柱术前三维影像与术中二维影像配准的方法与系统

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
共和国之辉: "《最优化五:梯度法(梯度下降法、最优梯度法、共轭梯度法)》", 《HTTPS://BLOG.CSDN.NET/LITTLEEMPEROR/ARTICLE/DETAILS/105073755》 *
张冉等: "《2D/3D图像配准中的相似性测度和优化算法》", 《激光与红外》 *

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113040908A (zh) * 2021-02-02 2021-06-29 武汉联影智融医疗科技有限公司 手术导航的配准方法、装置、计算机设备和存储介质
WO2023006021A1 (zh) * 2021-07-30 2023-02-02 武汉联影智融医疗科技有限公司 一种配准方法和系统
CN115082534A (zh) * 2022-07-21 2022-09-20 杭州三坛医疗科技有限公司 双平面图像配准方法、装置及机器人

Also Published As

Publication number Publication date
CN112233155B (zh) 2021-02-26

Similar Documents

Publication Publication Date Title
CN112233155B (zh) 一种2d-3d影像配准算法
US11605185B2 (en) System and method for generating partial surface from volumetric data for registration to surface topology image data
JP2022141792A5 (zh)
US6381302B1 (en) Computer assisted 2D adjustment of stereo X-ray images
RU2568635C2 (ru) Регистрация двумерных/трехмерных изображений на основе признаков
US8184886B2 (en) Deformable 2D-3D registration
US20180189966A1 (en) System and method for guidance of laparoscopic surgical procedures through anatomical model augmentation
US20050047544A1 (en) Apparatus and method for registering 2D radiographic images with images reconstructed from 3D scan data
US20150115178A1 (en) Methods and systems for guiding an emission to a target
CN109925054B (zh) 确定靶点路径的辅助方法、装置和系统、可读存储介质
CN104268849A (zh) 图像核对装置以及患者定位装置
CN107680688B (zh) 一种基于3d打印的盆腔仿真微创手术视觉导航验证方法
JP2022534123A (ja) 画像レジストレーション方法及びそれに関係するモデルトレーニング方法、デバイス、装置
CN112289416B (zh) 一种引导针置入精度评价方法
US20170270678A1 (en) Device and method for image registration, and non-transitory recording medium
CN109767458A (zh) 一种半自动分段的顺序优化配准方法
CN109925053B (zh) 手术路径的确定方法、装置和系统、可读存储介质
Maharjan et al. A novel visualization system of using augmented reality in knee replacement surgery: Enhanced bidirectional maximum correntropy algorithm
CN113963057B (zh) 成像几何关系标定方法、装置、电子设备以及存储介质
CN112472293B (zh) 一种术前三维影像与术中透视图像的配准方法
KR102469141B1 (ko) 의료용 화상 처리 장치, 의료용 화상 처리 방법, 및 프로그램
CN114511597A (zh) X光图像与ct图像的配准方法
JP7513980B2 (ja) 医用画像処理装置、治療システム、医用画像処理方法、およびプログラム
EP4144298A1 (en) Object visualisation in x-ray imaging
JP7407831B2 (ja) 介入装置追跡

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
CP03 Change of name, title or address

Address after: 210000 building 3, No. 34, Dazhou Road, Yuhuatai District, Nanjing, Jiangsu Province

Patentee after: Tuodao Medical Technology Co.,Ltd.

Address before: Room 102-86, building 6, 57 Andemen street, Yuhuatai District, Nanjing, Jiangsu 210000

Patentee before: Nanjing Tuodao Medical Technology Co.,Ltd.

CP03 Change of name, title or address