CN112102142A - 一种基于gpu并行计算的不同分辨率医学图像配准系统 - Google Patents
一种基于gpu并行计算的不同分辨率医学图像配准系统 Download PDFInfo
- Publication number
- CN112102142A CN112102142A CN202010951147.4A CN202010951147A CN112102142A CN 112102142 A CN112102142 A CN 112102142A CN 202010951147 A CN202010951147 A CN 202010951147A CN 112102142 A CN112102142 A CN 112102142A
- Authority
- CN
- China
- Prior art keywords
- image
- registration
- gpu
- affine
- deformation
- 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.)
- Pending
Links
- PXFBZOLANLWPMH-UHFFFAOYSA-N 16-Epiaffinine Natural products C1C(C2=CC=CC=C2N2)=C2C(=O)CC2C(=CC)CN(C)C1C2CO PXFBZOLANLWPMH-UHFFFAOYSA-N 0.000 claims description 42
- 238000000034 method Methods 0.000 claims description 22
- 230000009466 transformation Effects 0.000 claims description 14
- 238000009795 derivation Methods 0.000 claims description 6
- 239000011159 matrix material Substances 0.000 claims description 4
- 238000012897 Levenberg–Marquardt algorithm Methods 0.000 claims description 3
- 238000000354 decomposition reaction Methods 0.000 claims description 3
- 238000012804 iterative process Methods 0.000 claims description 3
- 238000012545 processing Methods 0.000 description 7
- 238000003384 imaging method Methods 0.000 description 6
- 238000005457 optimization Methods 0.000 description 3
- 238000005070 sampling Methods 0.000 description 3
- 210000004204 blood vessel Anatomy 0.000 description 2
- 230000005855 radiation Effects 0.000 description 2
- 238000001356 surgical procedure Methods 0.000 description 2
- 238000002679 ablation Methods 0.000 description 1
- 210000003484 anatomy Anatomy 0.000 description 1
- 238000004422 calculation algorithm Methods 0.000 description 1
- 238000004364 calculation method Methods 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 230000004927 fusion Effects 0.000 description 1
- 238000003709 image segmentation Methods 0.000 description 1
- 238000002347 injection Methods 0.000 description 1
- 239000007924 injection Substances 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000008092 positive effect Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000002945 steepest descent method Methods 0.000 description 1
- 238000012285 ultrasound imaging Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T1/00—General purpose image data processing
- G06T1/20—Processor architectures; Processor configuration, e.g. pipelining
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T3/00—Geometric image transformations in the plane of the image
- G06T3/02—Affine transformations
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T3/00—Geometric image transformations in the plane of the image
- G06T3/40—Scaling of whole images or parts thereof, e.g. expanding or contracting
- G06T3/4007—Scaling of whole images or parts thereof, e.g. expanding or contracting based on interpolation, e.g. bilinear interpolation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/30—Determination of transform parameters for the alignment of images, i.e. image registration
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Image Processing (AREA)
Abstract
本发明涉及一种基于GPU并行计算的不同分辨率医学图像配准系统，第一图像获取模块，用于获取参考图像；第二图像获取模块，用于获取待配准图像；控制点生成模块，用于在空间中生成多个初始变形参数为零的控制节点；配准模块，用于进行不同尺度下全局配准，并通过B样条曲线拟合的方式拟合所述控制节点上的变形参数。本发明能够在不提取特征点的前提下在GPU上实现不同分辨率医学图像的配准。
Description
技术领域
本发明涉及图像处理技术领域，特别是涉及一种基于GPU并行计算的不同分辨率医学图像配准系统。
背景技术
医学图像配准是图像处理研究领域中的一项重要且复杂的任务，由于医学图像会来自不同的采集设备、不同的扫描时间、不同的拍摄角度、不同的研究对象等。因此，往往在进行图像融合、图像分割等后续处理前需要先进行图像配准已达到针对同一解剖结构在不同的采集图像中对应于空间中的同一位置。
传统的图像配准方法主要使用CPU，并使用包括基于特征的方法，该方法往往需要利用边缘、表面、像素等信息进行特征提取，配准的成功与否很大程度上依赖于提取的特征。这就导致其存在注入适用性差、处理速度慢、配准评估不统一等问题。
现有专利文献CN111091567A公开了医学图像配准方法，该方法包括获取目标部位的参考医学图像序列以及待配准医学图像序列；确定待配准医学图像序列的各个图像层面中的血管边界，以得到各个图像层面的边界子图；分别对参考医学图像序列的至少一个图像层面以及待配准医学图像序列的至少一个图像层面进行血管边界的采样，得到第一边界点集合以及第二边界点集合；基于第一边界点集合以及第二边界点集合中对应边界点的偏移量，确定最优变换；根据最优变换对各个边界子图进行调整，以对所述待配准医学图像序列进行配准。由此可见，该方法需要提取血管边界作为特征点，当信息量较少时，且配准精度将受到较大的限制，另外，该方法中待配准图像和参考图像的大小一致，因此其不适用于复杂工作场景，如穿刺导航和不同成像设备间配准。
发明内容
本发明所要解决的技术问题是提供一种基于GPU并行计算的不同分辨率医学图像配准系统，能够在不提取特征点的前提下实现医学图像的配准。
本发明解决其技术问题所采用的技术方案是：提供一种基于GPU并行计算的不同分辨率医学图像配准系统，包括：第一图像获取模块，用于获取参考图像；第二图像获取模块，用于获取待配准图像；控制点生成模块，用于在空间中生成多个初始变形参数为零的控制节点；配准模块，用于进行不同尺度下全局配准并通过曲线B样条拟合的方式拟合所述控制节点上的变形参数。
所述带配准图像的分辨率和所述参考图像的分辨率不同。
所述配准模块包括仿射配准和自由变形配准，所述仿射配准是指通过优化所述参考图像的放射变形参数，最小化所述参考图像和待配准图像差异以消除图像全局仿射变形，所述自由变形配准是指通过B样条曲线拟合所述控制节点上的变形参数。
所述仿射配准具体为：给定初始仿射变换参数x，将所述仿射变换参数x由CPU复制到GPU，对所述待配准图像使用仿射变换参数x进行仿射变换，记为Imoving，在GPU上计算惩罚函数f，f为参考图像Iref和仿射变换后的待配准图像Imoving的差异，表示为：f＝0.5∑(Iref-Imoving)2；在GPU上通过链式求导法则计算惩罚函数f对仿射配准参数x的偏导数，记为J，将其传递到CPU上，再用Levenberg-Marquardt算法求迭代过程仿射配准参数x的增量dh，具体为：(JT*J+mu*I)*dh＝JT*f，其中，mu为迭代中根据惩罚函数f和J进行调整的参数，I为单位矩阵，dh通过Cholesky分解求得，以x＝x+dh更新仿射配准参数x，当满足迭代终止条件，迭代终止。
根据仿射配准计算求得仿射变换参数时，还对待配准图像Imoving进行缩放并通过插值计算，使待配准图像Imoving的分辨率和参考图像Iref一致。
所述自由变形配准具体为：给定常数α0为初始步长；将控制节点上的变形参数β初始化为0，并由CPU复制到GPU，在GPU中对图像进行B样条自由变形，得到变形后的图像，记为Imoving'，在GPU上计算惩罚函数f'，f'为参考图像Iref和变形后图像Imoving'的差异，表示为：f＝0.5∑(Iref-Imoving')2；在GPU上通过链式求导法则计算惩罚函数f'对变形参数β偏导数，记为J'，通过线搜索方法找到最佳步长，记为α，以x'＝x'-αJ更新步长x'，当满足迭代终止条件，迭代终止。
有益效果
由于采用了上述的技术方案，本发明与现有技术相比，具有以下的优点和积极效果：本发明通过控制节点的变形参数直接对原始二维或三维图像配准，不用提取特征点，无信息损失，配准精度高。本发明使用高分辨率术前图像做为基准参考图像，并通过插值配准与术中低分辨率图像进行配准，实现不同分辨率医学图像的配准，可以有效降低了术中成像扫描的CT辐射量和扫描时间。本发明还使用GPU进行并行计算，使得配准效率高，适用于更广泛医学场景。
附图说明
图1是本发明的结构方框图。
具体实施方式
下面结合具体实施例，进一步阐述本发明。应理解，这些实施例仅用于说明本发明而不用于限制本发明的范围。此外应理解，在阅读了本发明讲授的内容之后，本领域技术人员可以对本发明作各种改动或修改，这些等价形式同样落于本申请所附权利要求书所限定的范围。
本发明的实施方式涉及一种基于GPU并行计算的不同分辨率医学图像配准系统，如图1所示，包括：第一图像获取模块、第二图像获取模块、控制点生成模块、配准模块和输出模块。
本实施方式中，第一图像获取模块，用于获取参考图像。该第一图像获取模块通常在对病人术前进行全面扫描时使用，例如可以使用CT，MRI或超声波成像等方式使用较高的分辨率和较小的成像层间距对病人进行全面三维术前扫描，获得二维或三维的参考图像，记为Iref。
第二图像获取模块，用于获取待配准图像。该第二图像获取模块通常在手术中使用，例如在消融/穿刺手术中对病人进行快速扫描或使用较大的成像层间距进行扫描，获得低分辨率的二维或三维图像作为待配准图像，记为Imoving。
所述控制点生成模块前还包括输入模块，所述输入模块用于将所述参考图像和待配准图像输入至GPU中，通过GPU进行并行计算。
控制点生成模块，用于在空间中生成多个初始变形参数为零的控制节点，通过改变这些控制节点的变形参数，能够使得带配准图像和参考图像之间的偏差缩小，从而实现配准的目的。
配准模块，用于进行不同尺度下全局配准并通过曲线拟合的方式拟合所述控制节点上的变形参数。由于本实施方式中，参考图像的分辨率高于待配准图像，因此在配准模块中还需要对参考图像进行插值处理，具体处理方式为：根据仿射变换参数，对待配准图像Imoving进行缩放和插值，使得所述待配准图像和参考图像具有相同的分辨率和矩阵大小。当两者处于相同的分辨率时，使用B样条曲线拟合控制节点上的变形参数，并使用最速下降法和线搜索方法对变形参数进行非线性优化以最小化待配准图像Imoving和参考图像Iref之间的偏差，从而实现配准的目的。
本实施方式中包括仿射配准和自由变形配准，所述仿射配准是指通过优化所述参考图像的放射变形参数，最小化所述参考图像和待配准图像差异以消除图像全局仿射变形，所述自由变形配准是指通过B样条曲线拟合所述控制节点上的变形参数。
所述仿射配准具体为：给定仿射配准参数x，将所述仿射配准参数x由CPU复制到GPU，对所述参考图像和待配准图像采用不同降采样系数进行降采样；对所述待配准图像使用仿射配准参数x进行仿射变换，记为Imoving，在GPU上计算惩罚函数f，f为参考图像Iref和仿射变换后的待配准图像Imoving的差异，表示为：f＝0.5∑(Iref-Imoving)2；在GPU上通过链式求导法则计算惩罚函数f对仿射配准参数x的偏导数，记为J，将其传递到CPU上，再用Levenberg-Marquardt算法求迭代过程仿射配准参数x的增量dh，具体为：(JT*J+mu*I)*dh＝JT*f，其中，mu为迭代中根据惩罚函数f和J进行调整的参数，I为单位矩阵，dh通过Cholesky分解求得，更新仿射配准参数x＝x+dh，当满足迭代终止条件，迭代终止。
所述自由变形配准具体为：给定常数α0为步长x'；将控制节点上的变形参数β初始化为0，并由CPU复制到GPU，在GPU中对图像进行B样条自由变形，得到变形后的图像，记为Imoving'，在GPU上计算惩罚函数f'，f'为参考图像Iref和变形后图像Imoving'的差异，表示为：f＝0.5∑(Iref-Imoving')2；在GPU上通过链式求导法则计算惩罚函数f'对变形参数β偏导数，记为J'，通过黄线搜索方法(包括但不限于Wolfe-Powell方法，黄金分割搜索方法)最佳步长，记为α，更新步长x'＝x'-αJ，当满足迭代终止条件，迭代终止。
值得一提的是，本实施方式还可以使用其他的非线性优化方法实现对变形参数进行非线性优化，比如LM算法和高斯-牛顿法等。
输出模块，用于将控制节点上的变形参数的值和配准后的图像从GPU复制到CPU并进行输出。
不难发现，本发明通过控制节点的变形参数直接对原始二维或三维图像配准，不用提取特征点，无信息损失，配准精度高。本发明使用高分辨率术前图像做为基准参考图像，并通过插值配准与术中低分辨率图像进行配准，实现不同分辨率医学图像的配准，并且降低了术中成像扫描CT辐射量和扫描时间。本发明还使用GPU进行并行计算，使得配准效率高，适用于更广泛医学场景。
Claims (6)
1.一种基于GPU并行计算的不同分辨率医学图像配准系统，其特征在于，包括：第一图像获取模块，用于获取参考图像；第二图像获取模块，用于获取待配准图像；控制点生成模块，用于在空间中生成多个初始变形参数为零的控制节点；配准模块，用于进行不同尺度下全局配准并通过曲线B样条拟合的方式拟合所述控制节点上的变形参数。
2.根据权利要求1所述的基于GPU并行计算的不同分辨率医学图像配准系统，其特征在于，所述带配准图像的分辨率和所述参考图像的分辨率不同。
3.根据权利要求1所述的基于GPU并行计算的不同分辨率医学图像配准系统，其特征在于，所述配准模块包括仿射配准和自由变形配准，所述仿射配准是指通过优化所述参考图像的仿射变形参数，最小化所述参考图像和待配准图像差异以消除图像全局仿射变形影响，所述自由变形配准是指通过B样条曲线拟合所述控制节点上的变形参数。
4.根据权利要求3所述的基于GPU并行计算的不同分辨率医学图像配准系统，其特征在于，所述仿射配准具体为：给定初始仿射变换参数x，将所述仿射变换参数x由CPU复制到GPU，对所述待配准图像使用仿射变换参数x进行仿射变换，记为Imoving，在GPU上计算惩罚函数f，f为参考图像Iref和仿射变换后的待配准图像Imoving的差异，表示为：f＝0.5∑(Iref-Imoving)2；在GPU上通过链式求导法则计算惩罚函数f对仿射配准参数x的偏导数，记为J，将其传递到CPU上，再用Levenberg-Marquardt算法求迭代过程仿射配准参数x的增量dh，具体为：(JT*J+mu*I)*dh＝JT*f，其中，mu为迭代中根据惩罚函数f和J进行调整的参数，I为单位矩阵，dh通过Cholesky分解求得，以x＝x+dh更新仿射配准参数x，当满足迭代终止条件，迭代终止。
5.根据权利要求3所述的基于GPU并行计算的不同分辨率医学图像配准系统，其特征在于，根据仿射配准计算求得仿射变换参数时，还对待配准图像Imoving进行缩放并通过插值计算，使待配准图像Imoving的分辨率和参考图像Iref一致。
6.根据权利要求3所述的基于GPU并行计算的不同分辨率医学图像配准系统，其特征在于，所述自由变形配准具体为：给定常数α0为初始步长；将控制节点上的变形参数β初始化为0，并由CPU复制到GPU，在GPU中对图像进行B样条自由变形，得到变形后的图像，记为Imoving'，在GPU上计算惩罚函数f'，f'为参考图像Iref和变形后图像Imoving'的差异，表示为：f＝0.5∑(Iref-Imoving')2；在GPU上通过链式求导法则计算惩罚函数f'对变形参数β偏导数，记为J'，通过线搜索方法找到最佳步长，记为α，以x'＝x'-αJ更新步长x'，当满足迭代终止条件，迭代终止。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010951147.4A CN112102142A (zh) | 2020-09-11 | 2020-09-11 | 一种基于gpu并行计算的不同分辨率医学图像配准系统 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010951147.4A CN112102142A (zh) | 2020-09-11 | 2020-09-11 | 一种基于gpu并行计算的不同分辨率医学图像配准系统 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN112102142A true CN112102142A (zh) | 2020-12-18 |
Family
ID=73752082
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202010951147.4A Pending CN112102142A (zh) | 2020-09-11 | 2020-09-11 | 一种基于gpu并行计算的不同分辨率医学图像配准系统 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN112102142A (zh) |
Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20080050043A1 (en) * | 2006-08-22 | 2008-02-28 | Siemens Medical Solutions Usa, Inc. | Methods and Systems for Registration of Images |
| CN104574278A (zh) * | 2015-02-03 | 2015-04-29 | 中国人民解放军国防科学技术大学 | 基于局部反馈机制与稀疏全局调整的航拍图像拼接方法 |
| CN106611411A (zh) * | 2015-10-19 | 2017-05-03 | 上海联影医疗科技有限公司 | 一种医学图像中肋骨分割的方法及医学图像处理装置 |
| US20170278244A1 (en) * | 2016-03-24 | 2017-09-28 | The Chinese University Of Hong Kong | Method and a system for non-rigid image registration |
| CN108765474A (zh) * | 2018-04-17 | 2018-11-06 | 天津工业大学 | 一种针对ct与光学扫描牙齿模型的高效配准方法 |
-
2020
- 2020-09-11 CN CN202010951147.4A patent/CN112102142A/zh active Pending
Patent Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20080050043A1 (en) * | 2006-08-22 | 2008-02-28 | Siemens Medical Solutions Usa, Inc. | Methods and Systems for Registration of Images |
| CN104574278A (zh) * | 2015-02-03 | 2015-04-29 | 中国人民解放军国防科学技术大学 | 基于局部反馈机制与稀疏全局调整的航拍图像拼接方法 |
| CN106611411A (zh) * | 2015-10-19 | 2017-05-03 | 上海联影医疗科技有限公司 | 一种医学图像中肋骨分割的方法及医学图像处理装置 |
| US20170278244A1 (en) * | 2016-03-24 | 2017-09-28 | The Chinese University Of Hong Kong | Method and a system for non-rigid image registration |
| CN108765474A (zh) * | 2018-04-17 | 2018-11-06 | 天津工业大学 | 一种针对ct与光学扫描牙齿模型的高效配准方法 |
Non-Patent Citations (7)
| Title |
|---|
| SERDAR K. BALCI ET AL.: "Non-rigid Groupwise Registration using B-Spline Deformation Model", 《RESEARCHGATE》 * |
| 何永富等: "《数值优化计算方法与应用》", 31 October 1994, 成都科技大学出版社 * |
| 信息化理论学术研讨会论文集编委会: "《信息化理论与综合信息系统 中国电子学会电子系统工程分会第十三届信息化理论学术研讨会论文集》", 国防工业出版社 * |
| 徐进: "基于图像配准的心脏冠脉 CTA 模型建立方法研究", 《中国优秀博硕士学位论文全文数据库(硕士) 医药卫生科技辑》 * |
| 杜晓刚: "图像引导放疗中的医学图像配准关键技术研究", 《中国优秀博士学位论文全文数据库》 * |
| 罗四维等: "《计算机视觉检测逆问题导论》", 31 March 2017, 北京交通大学出版社 * |
| 谭兆: "遥感影像与地形图自动配准方法的研究", 《科技创新导报》 * |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| JP7221421B2 (ja) | Ｃｔ画像の椎骨位置決め方法、装置、機器及び媒体 | |
| US7945117B2 (en) | Methods and systems for registration of images | |
| US7565019B2 (en) | Method of volume-panorama imaging processing | |
| US8861891B2 (en) | Hierarchical atlas-based segmentation | |
| CN107610162A (zh) | 一种基于互信息与图像分割的三维多模态医学图像自动配准方法 | |
| CN111091567B (zh) | 医学图像配准方法、医疗设备及存储介质 | |
| JP2009520558A (ja) | ポイント・ベースの適応的弾性画像登録 | |
| US7440628B2 (en) | Method and system for motion correction in a sequence of images | |
| US7062078B2 (en) | Method and device for the registration of images | |
| JP2006127241A (ja) | 超解像処理の高速化方法 | |
| US6597816B1 (en) | Correcting distortion in an imaging system using parametric motion estimation | |
| JP7083015B2 (ja) | 超解像処理装置、方法及びプログラム | |
| JP4274400B2 (ja) | 画像の位置合わせ方法および装置 | |
| CN110930444B (zh) | 一种基于双边优化的点云匹配方法、介质、终端和装置 | |
| CN109961435B (zh) | 脑图像获取方法、装置、设备及存储介质 | |
| CN116363181A (zh) | 一种基于特征的ct图像和超声图像肝脏配准方法 | |
| WO2018133027A1 (zh) | 基于灰度约束的三维数字散斑的整像素搜索方法及装置 | |
| Ashtari et al. | Double line image rotation | |
| CN108076341A (zh) | 一种视频卫星成像在轨实时稳像方法及系统 | |
| CN112102142A (zh) | 一种基于gpu并行计算的不同分辨率医学图像配准系统 | |
| CN112700439A (zh) | 一种基于神经网络的oct人眼图像采集配准方法及系统 | |
| US20230099906A1 (en) | Image registration method, computer device, and storage medium | |
| US20230206477A1 (en) | Image processing method, image processing device, program, and trained model | |
| KR102348863B1 (ko) | 비지도 학습 기법을 이용한 영상데이터의 정합 방법 및 장치 | |
| CN117036428A (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 | ||
| RJ01 | Rejection of invention patent application after publication |
Application publication date: 20201218 |
|
| RJ01 | Rejection of invention patent application after publication |