CN104537640B - 基于梯度直方图分布匹配确定三维旋转量的方法和系统 - Google Patents

基于梯度直方图分布匹配确定三维旋转量的方法和系统 Download PDF

Info

Publication number
CN104537640B
CN104537640B CN201410708512.3A CN201410708512A CN104537640B CN 104537640 B CN104537640 B CN 104537640B CN 201410708512 A CN201410708512 A CN 201410708512A CN 104537640 B CN104537640 B CN 104537640B
Authority
CN
China
Prior art keywords
dimensional
histogram
gradient
rotation amount
volume data
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
CN201410708512.3A
Other languages
English (en)
Other versions
CN104537640A (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.)
Shenzhen Institute of Advanced Technology of CAS
Original Assignee
Shenzhen Institute of Advanced Technology 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 Shenzhen Institute of Advanced Technology of CAS filed Critical Shenzhen Institute of Advanced Technology of CAS
Priority to CN201410708512.3A priority Critical patent/CN104537640B/zh
Publication of CN104537640A publication Critical patent/CN104537640A/zh
Application granted granted Critical
Publication of CN104537640B publication Critical patent/CN104537640B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • 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/10004Still image; Photographic image
    • G06T2207/10012Stereo images
    • 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

  • Apparatus For Radiation Diagnosis (AREA)
  • Image Processing (AREA)
  • Image Analysis (AREA)

Abstract

本发明提供了一种基于梯度直方图分布匹配确定三维旋转量的方法,其中,所述方法包括:获取待配准的三维体数据;通过三维高斯梯度核函数计算所述待配准的三维体数据对应的三维局部梯度;通过累积所述三维局部梯度的每个体素的幅值和方向值构建所述三维局部梯度对应的幅值加权方向分布直方图;对所述幅值加权方向分布直方图进行直方图归一化匹配确定三维旋转量。上述方法和系统相比传统图像配准技术确定三维旋转量具有简单、快速的优点。

Description

基于梯度直方图分布匹配确定三维旋转量的方法和系统
技术领域
本发明涉及图像处理技术领域,特别是涉及一种基于梯度直方图分布匹配确定三维旋转量的方法和系统。
背景技术
三维图像配准是将两幅或多幅包含相同场景或目标的图像进行几何对准的过程,其中一幅是参考图像,另外一幅是检测图像。图像配准技术主要应用于医学计算机断层扫描、核磁共振成像、肿瘤监控、治疗验证以及病人的解剖图前后对比分析等。图像配准通常可被分为刚性配准和非刚性配准,其中刚性配准涉及的变换通常只包括平移和旋转,实现的是将两幅图像整体对齐,计算复杂度低而且具有较高的鲁棒性。
在临床医学配准中,一般将刚性配准作为预配准,用于消除影像间较大的全局形变差异,然后采用非刚性配准作为精配准实现局部形变的配准。因此,快速、高效的刚性配准在医学临床中仍然是急需的。特别是在三维体数据中如何实现快速准确确定三维旋转量。
目前实现确定待配准的三维体数据的三维旋转量,一般是预先在病人体内植入标记或基于点、线的特征提取和匹配方式,以及基于相似性度量的迭代优化方式。从而依据标记点实现三维旋转量的度量。在于体内植入标记点需要手术麻醉完成,增加了病患的负担。同时,在医学体数据影像中提取标记点和匹配,往往比较繁琐和费时。
发明内容
基于此,有必要针对上述技术问题,提供一种简单快速确定三维图像配准中三维旋转量的基于梯度直方图分布匹配确定三维旋转量的方法和系统。
一种基于梯度直方图分布匹配确定三维旋转量的方法,所述方法包括:
获取待配准的三维体数据;
通过三维高斯梯度核函数计算所述待配准的三维体数据对应的三维局部梯度;
通过累积所述三维局部梯度的每个体素的幅值和方向值构建所述三维局部梯度对应的幅值加权方向分布直方图;
对所述幅值加权方向分布直方图进行直方图归一化匹配确定三维旋转量。
在其中一个实施例中,所述对所述幅值加权方向分布直方图进行直方图归一化匹配确定三维旋转量的步骤之后,所述方法还包括:
通过将所述三维旋转量进行直接的矩阵相乘构建三维旋转矩阵。
在其中一个实施例中,所述通过三维高斯梯度核函数计算所述待配准的三维体数据对应的三维局部梯度的步骤,包括:
将所述三维体数据与三维高斯梯度核函数进行卷积运算得到三维局部梯度。
在其中一个实施例中,所述对所幅值加权方向分布直方图进行直方图归一化匹配确定三维旋转量的步骤,包括:
将幅值加权方向分布直方图代入预设的最小化代价函数进行计算得到最优平移量,所述最优平移量即为三维旋转量。
在其中一个实施例中,所述三维高斯梯度核函数是三维高斯函数对三维坐标系中三个坐标轴方向分别求一阶导数得到的函数表达式。
一种基于梯度直方图分布匹配确定三维旋转量的系统,所述系统包括:
数据获取模块,用于获取待配准的三维体数据;
梯度计算模块,用于通过三维高斯梯度核函数计算所述待配准的三维体数据对应的三维局部梯度;
直方图构建模块,用于通过累积所述三维局部梯度的每个体素的幅值和方向值构建所述三维局部梯度对应的幅值加权方向分布直方图;
旋转量确定模块,用于对所述幅值加权方向分布直方图进行直方图归一化匹配确定三维旋转量。
在其中一个实施例中,所述系统还包括:
矩阵构建模块,用于通过将所述三维旋转量进行直接的矩阵相乘构建三维旋转矩阵。
在其中一个实施例中,梯度计算模块还用于将所述三维体数据与三维高斯梯度核函数进行卷积运算得到三维局部梯度。
在其中一个实施例中,所述旋转量确定模块包括:
将幅值加权方向分布直方图代入预设的最小化代价函数进行计算得到最优平移量,所述最优平移量即为三维旋转量。
在其中一个实施例中,所述三维高斯梯度核函数是三维高斯函数对三维坐标系中三个坐标轴方向分别求一阶导数得到的函数表达式。
上述基于梯度直方图分布匹配确定三维旋转量的方法和系统,由于在确定待配准的三维体数据的三维旋转量过程中不需要进行任何的迭代优化过程,也不需要在人体内植入任何标记点或在三维体数据中进行三维特征的提取和匹配,相比传统图像配准技术确定三维旋转量具有简单、快速的优点。
附图说明
图1为一个实施例中基于梯度直方图分布匹配确定三维旋转量的方法流程示意图;
图2为一个实施例中三维旋转坐标系统的结构示意图;
图3为一个实施例中用图形来表示三维高斯梯度核函数的矩阵图;
图4为另一个实施例中基于梯度直方图分布匹配确定三维旋转量的方法流程示意图;
图5为一个实施例中将本发明提供的方法应用于医学三维临床数据处理的效果图;
图6为一个实施例中基于梯度直方图分布匹配确定三维旋转量的系统结构示意图;
图7为另一个实施例中基于梯度直方图分布匹配确定三维旋转量的系统结构示意图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白,以下结合附图及实施例,对本发明进行进一步详细说明。应当理解,此处所描述的具体实施例仅仅用以解释本发明,并不用于限定本发明。
如图1所示,在一个实施例中,提供了一种基于梯度直方图分布匹配确定三维旋转量的方法,该方法包括:
步骤101,获取待配准的三维体数据。
待配准三维体数据包括:参考三维体数据和检测三维体数据。本实施例中,三维图像配准是指刚性配准,即参考三维体数据和检测三维体数据之间仅存在旋转和平移的变换。待配准的三维体数据可以通过医学计算机断层扫描成像(CT)、核磁共振成像(MR)等高端医学影像设备获得。
步骤102,通过三维高斯梯度核函数计算待配准的三维体数据对应的三维局部梯度。
三维高斯梯度核函数是三维高斯函数对三维坐标系中三个坐标轴方向分别求一阶导数得到的函数表达式。将三维体数据与三维高斯梯度核函数进行卷积运算得到三维局部梯度。三维局部梯度的计算直接影响对三维旋转量的估计,
本实施例中采用三维高斯梯度核函数分别计算待配准三维体数据在如图2所示的三维坐标系中X,Y和Z三个坐标轴方向的梯度。高斯函数和高斯梯度函数的表达式如下:
其中σ是高斯函数的标准偏差。对X,Y和Z方向的三维各向同性高斯梯度函数定义如下:
其中是标准偏差为σ的三维高斯函数对X,Y,Z轴的一阶导数。归一化生成高斯梯度核函数。
具体的,归一化过程为:将函数离散化,即采用离散化的矩阵替代原始的连续函数,便于计算机处理。并将离散化的矩阵进行归一化处理使得矩阵中的所有元素之和为1生成高斯梯度核函数,高斯梯度核函数其实就是一个三维矩阵。
如图3所示,为一个实施例中,通过图形来表示对进行归一化处理生成的7Χ7Χ7的三维高斯梯度核函数,32表明矩阵中值较大(且>0),而31则表明值较小(且<0)。
将待配准的三维体数据分别与三维高斯梯度核函数进行卷积运算,获得三维局部梯度[Gx,Gy,Gz],即:Gx是由待配准的三维体数据与卷积运算获得,Gy是由待配准的三维体数据与卷积运算获得,Gz是由待配准的三维体数据与卷积运算获得。
步骤103,通过累积三维局部梯度的每个体素的幅值和方向值构建三维局部梯度对应的幅值加权方向分布直方图。
本实施例中,直方图的角度分辨率可以预先设定,为了提高三维旋转量的估计精度,角度分辨率设定较小。此外,为了降低三维体数据间尺度变换的影响,采用了直方图归一化处理。分别创建三个分别针对旋转角度θXYYZ和θZX的幅值加权方向分布直方图HXY,HYZ,和HZX
步骤104,对幅值加权方向分布直方图进行直方图归一化匹配确定三维旋转量。
具体的,将幅值加权方向分布直方图代入预设的最小化代价函数进行计算得到最优平移量,该最优平移量即为三维旋转量。
本实施例中,直方图匹配主要是确定待配准的三维体数据中参考三维体数据与检测三维体数据之间的三维旋转量,若HA和HB分别对应参考三维体数据的幅值加权方向分布直方图和检测三维体数据的幅值加权方向分布直方图,直方图匹配主要是确定HA和HB之间的最优平移量,即公式(5)中最小化代价函数TD
求J使
其中,K是直方图的角度分块数(即:1/角度分辨率后取整的值)。
方向直方图HA和HB之间相关的最佳平移量可以通过公式(7)逆傅里叶变换确定:
其中表示一维傅里叶变换,表示一维逆傅里叶变换,*表示复数共轭,·表示点乘运算。最优平移量可以由公式(7)计算出来的参数C的最大值位置确定,即为公式(6)中TD(J)的最小值位置J。
上述基于梯度直方图分布匹配确定三维旋转量的方法,由于在确定待配准的三维体数据的三维旋转量过程中不需要进行任何的迭代优化过程,也不需要在人体内植入任何标记点或在三维体数据中进行三维特征的提取和匹配,相比传统图像配准技术确定三维旋转量具有简单、快速的优点。
本实施例中,创造性的将三维旋转量的计算转化为了对归一化直方图匹配的最优平移量测度。通过将复杂的三维旋转量计算转为采用归一化直方图匹配确定直方图之间的平移量。减少了三维旋转量的计算量。
如图4所示,在一个实施例中,基于三维图像配准确定三维旋转量的方法还包括:
步骤105,通过将三维旋转量进行直接的矩阵相乘构建三维旋转矩阵。
本实施例中,三维旋转量分别对应为θSI,θAP和θRL,根据以下公式(8)构造三维体数据间的三维旋转矩阵:
三维旋转矩阵的构造,可以通过直接矩阵相乘实现体数据间旋转差异的消除。计算出来的R即为三维旋转矩阵。计算三维旋转矩阵的目的是用于医学临床肿瘤监控(肿瘤的生长或消亡的趋势监控)、治疗验证(放疗或化疗前后的影像变化,用以评估验证治疗效果)以及病人的解剖图前后对比(手术前和手术后的效果对比)等。因此,在待配准的三维体数据的三维旋转量确定后,可以通过包含该测量值的三维旋转矩阵,实现待配准的三维影像的旋转差异的消除,从而实现上述医学临床的监控、治疗对比等目的。
如图5所示为将本发明提供的方法应用于医学三维临床数据处理的效果图。其中51为待配准的三维体数据;52为对待配准的三维体数据构建的幅值加权分布直方图;53为对直方图归一化匹配结果获得的三维旋转量θSI,θAP,θRL分别为9°,-6°和7°。采用本发明提供的方法在计算机的CPU 3.3GHz,8GRAM的Matlab平台下,检测所耗费的时间约为4-6秒。如果采用其它的基于优化迭代的方法(标准ITK提供的函数,C++程序运行),一般耗费时间为30-40秒,而如果采用基于标记点植入或特征匹配对应的方法,一般耗费时间为180-300秒。因此,本发明提供的方法在效率上具有优势。
如图6所示,在一个实施例中,提供了一种基于梯度直方图分布匹配确定三维旋转量的系统,该系统包括:
数据获取模块61,用于获取待配准的三维体数据。
梯度计算模块62,用于通过三维高斯梯度核函数计算待配准的三维体数据对应的三维局部梯度。本实施例中,三维高斯梯度核函数是三维高斯函数对三维坐标系中三个坐标轴方向分别求一阶导数得到的函数表达式。
直方图构建模块63,用于通过累积三维局部梯度的每个体素的幅值和方向值构建三维局部梯度对应的幅值加权方向分布直方图。
旋转量确定模块64,用于对幅值加权方向分布直方图进行直方图归一化匹配确定三维旋转量。
如图7所示,在一个实施例中,基于三维图像配准确定三维旋转量的系统还包括:
矩阵构建模块65,用于通过将三维旋转量进行直接的矩阵相乘构建三维旋转矩阵。
在一个实施例中,梯度计算模块62还用于将三维体数据与三维高斯梯度核函数进行卷积运算得到三维局部梯度。
在一个实施例中,旋转量确定模块64还用于将幅值加权方向分布直方图代入预设的最小化代价函数进行计算得到最优平移量,最优平移量即为三维旋转量。
以上实施例仅表达了本发明的几种实施方式,其描述较为具体和详细,但并不能因此而理解为对本发明专利范围的限制。应当指出的是,对于本领域的普通技术人员来说,在不脱离本发明构思的前提下,还可以做出若干变形和改进,这些都属于本发明的保护范围。因此,本发明专利的保护范围应以所附权利要求为准。

Claims (8)

1.一种基于梯度直方图分布匹配确定三维旋转量的方法,所述方法包括:
获取待配准的三维体数据;
通过三维高斯梯度核函数计算所述待配准的三维体数据对应的三维局部梯度;
通过累积所述三维局部梯度的每个体素的幅值和方向值构建所述三维局部梯度对应的幅值加权方向分布直方图;
对所述幅值加权方向分布直方图进行直方图归一化匹配确定三维旋转量;
所述三维高斯梯度核函数是三维高斯函数对三维坐标系中三个坐标轴方向分别求一阶导数,并对所述三个坐标轴方向的一阶导数进行归一化处理得到的函数表达式;
所述三维高斯函数对三维坐标系中三个坐标轴方向分别求得的一阶导数为:
<mrow> <msubsup> <mi>G</mi> <mi>x</mi> <mo>&amp;prime;</mo> </msubsup> <mrow> <mo>(</mo> <mi>&amp;sigma;</mi> <mo>,</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>=</mo> <mo>-</mo> <mfrac> <mi>x</mi> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mfrac> <msup> <mrow> <mo>(</mo> <mfrac> <mn>1</mn> <mrow> <msqrt> <mrow> <mn>2</mn> <mi>&amp;pi;</mi> </mrow> </msqrt> <mi>&amp;sigma;</mi> </mrow> </mfrac> <mo>)</mo> </mrow> <mn>3</mn> </msup> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>x</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>*</mo> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>y</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>*</mo> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>z</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>;</mo> </mrow>
<mrow> <msubsup> <mi>G</mi> <mi>y</mi> <mo>&amp;prime;</mo> </msubsup> <mrow> <mo>(</mo> <mi>&amp;sigma;</mi> <mo>,</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>=</mo> <mo>-</mo> <mfrac> <mi>y</mi> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mfrac> <msup> <mrow> <mo>(</mo> <mfrac> <mn>1</mn> <mrow> <msqrt> <mrow> <mn>2</mn> <mi>&amp;pi;</mi> </mrow> </msqrt> <mi>&amp;sigma;</mi> </mrow> </mfrac> <mo>)</mo> </mrow> <mn>3</mn> </msup> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>x</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>*</mo> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>y</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>*</mo> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>z</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>;</mo> </mrow>
<mrow> <msubsup> <mi>G</mi> <mi>z</mi> <mo>&amp;prime;</mo> </msubsup> <mrow> <mo>(</mo> <mi>&amp;sigma;</mi> <mo>,</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>=</mo> <mo>-</mo> <mfrac> <mi>z</mi> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mfrac> <msup> <mrow> <mo>(</mo> <mfrac> <mn>1</mn> <mrow> <msqrt> <mrow> <mn>2</mn> <mi>&amp;pi;</mi> </mrow> </msqrt> <mi>&amp;sigma;</mi> </mrow> </mfrac> <mo>)</mo> </mrow> <mn>3</mn> </msup> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>x</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>*</mo> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>y</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>*</mo> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>z</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>;</mo> </mrow>
其中,G′x(σ,x,y,z),G′y(σ,x,y,z),和G′z(σ,x,y,z)分别表示标准偏差为σ的三维高斯函数对X,Y,Z轴的一阶导数;
所述对所述三个坐标轴方向的一阶导数进行归一化处理的步骤包括:将所述一阶导数离散化,获得离散化的矩阵,将所述离散化的矩阵进行归一化处理,使得矩阵中的所有元素之和为1,生成三维高斯梯度核函数;
所述对所述幅值加权方向分布直方图进行直方图归一化匹配确定三维旋转量的步骤包括:
对参考三维体数据的幅值加权方向分布直方图和检测三维体数据的幅值加权方向分布直方图进行逆傅里叶变换,得到最佳平移量;所述逆傅里叶变换的计算公式为:
其中,表示一维傅里叶变换,表示一维逆傅里叶变换,*表示复数共轭,·表示点乘运算,HA表示参考三维体数据的幅值加权方向分布直方图,HB表示检测三维体数据的幅值加权方向分布直方图,C表示参数,最佳平移量由参数C的最大值位置确定。
2.根据权利要求1所述的方法,其特征在于,所述对所述幅值加权方向分布直方图进行直方图归一化匹配确定三维旋转量的步骤之后,所述方法还包括:
通过将所述三维旋转量进行直接的矩阵相乘构建三维旋转矩阵。
3.根据权利要求1所述的方法,其特征在于,所述通过三维高斯梯度核函数计算所述待配准的三维体数据对应的三维局部梯度的步骤,包括:
将所述三维体数据与三维高斯梯度核函数进行卷积运算得到三维局部梯度。
4.根据权利要求1所述的方法,其特征在于,所述对幅值加权方向分布直方图进行直方图归一化匹配确定三维旋转量的步骤,包括:
将幅值加权方向分布直方图代入预设的最小化代价函数进行计算得到最优平移量,所述最优平移量即为三维旋转量。
5.一种基于梯度直方图分布匹配确定三维旋转量的系统,其特征在于,所述系统包括:
数据获取模块,用于获取待配准的三维体数据;
梯度计算模块,用于通过三维高斯梯度核函数计算所述待配准的三维体数据对应的三维局部梯度;
直方图构建模块,用于通过累积所述三维局部梯度的每个体素的幅值和方向值构建所述三维局部梯度对应的幅值加权方向分布直方图;
旋转量确定模块,用于对所述幅值加权方向分布直方图进行直方图归一化匹配确定三维旋转量;
所述基于梯度直方图分布匹配确定三维旋转量的系统,还用于通过三维高斯函数对三维坐标系中三个坐标轴方向分别求一阶导数,并对所述三个坐标轴方向的一阶导数进行归一化处理得到三维高斯梯度核函数;
所述通过三维高斯函数对三维坐标系中三个坐标轴方向分别求得的一阶导数为:
<mrow> <msubsup> <mi>G</mi> <mi>x</mi> <mo>&amp;prime;</mo> </msubsup> <mrow> <mo>(</mo> <mi>&amp;sigma;</mi> <mo>,</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>=</mo> <mo>-</mo> <mfrac> <mi>x</mi> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mfrac> <msup> <mrow> <mo>(</mo> <mfrac> <mn>1</mn> <mrow> <msqrt> <mrow> <mn>2</mn> <mi>&amp;pi;</mi> </mrow> </msqrt> <mi>&amp;sigma;</mi> </mrow> </mfrac> <mo>)</mo> </mrow> <mn>3</mn> </msup> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>x</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>*</mo> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>y</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>*</mo> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>z</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>;</mo> </mrow>
<mrow> <msubsup> <mi>G</mi> <mi>y</mi> <mo>&amp;prime;</mo> </msubsup> <mrow> <mo>(</mo> <mi>&amp;sigma;</mi> <mo>,</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>=</mo> <mo>-</mo> <mfrac> <mi>y</mi> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mfrac> <msup> <mrow> <mo>(</mo> <mfrac> <mn>1</mn> <mrow> <msqrt> <mrow> <mn>2</mn> <mi>&amp;pi;</mi> </mrow> </msqrt> <mi>&amp;sigma;</mi> </mrow> </mfrac> <mo>)</mo> </mrow> <mn>3</mn> </msup> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>x</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>*</mo> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>y</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>*</mo> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>z</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>;</mo> </mrow>
<mrow> <msubsup> <mi>G</mi> <mi>z</mi> <mo>&amp;prime;</mo> </msubsup> <mrow> <mo>(</mo> <mi>&amp;sigma;</mi> <mo>,</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>=</mo> <mo>-</mo> <mfrac> <mi>z</mi> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mfrac> <msup> <mrow> <mo>(</mo> <mfrac> <mn>1</mn> <mrow> <msqrt> <mrow> <mn>2</mn> <mi>&amp;pi;</mi> </mrow> </msqrt> <mi>&amp;sigma;</mi> </mrow> </mfrac> <mo>)</mo> </mrow> <mn>3</mn> </msup> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>x</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>*</mo> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>y</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>*</mo> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msup> <mi>z</mi> <mn>2</mn> </msup> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> </mrow> </msup> <mo>;</mo> </mrow>
其中,G′x(σ,x,y,z),G′y(σ,x,y,z),和G′z(σ,x,y,z)分别表示标准偏差为σ的三维高斯函数对X,Y,Z轴的一阶导数;
所述基于梯度直方图分布匹配确定三维旋转量的系统,还用于将所述一阶导数离散化,获得离散化的矩阵,将所述离散化的矩阵进行归一化处理,使得矩阵中的所有元素之和为1,生成三维高斯梯度核函数;
所述基于梯度直方图分布匹配确定三维旋转量的系统,还用于对参考三维体数据的幅值加权方向分布直方图和检测三维体数据的幅值加权方向分布直方图进行逆傅里叶变换,得到最佳平移量;所述逆傅里叶变换的计算公式为:
其中,表示一维傅里叶变换,表示一维逆傅里叶变换,*表示复数共轭,·表示点乘运算,HA表示参考三维体数据的幅值加权方向分布直方图,HB表示检测三维体数据的幅值加权方向分布直方图,C表示参数,最佳平移量由参数C的最大值位置确定。
6.根据权利要求5所述的系统,其特征在于,所述系统还包括:
矩阵构建模块,用于通过将所述三维旋转量进行直接的矩阵相乘构建三维旋转矩阵。
7.根据权利要求5所述的系统,其特征在于,梯度计算模块还用于将所述三维体数据与三维高斯梯度核函数进行卷积运算得到三维局部梯度。
8.根据权利要求5所述的系统,其特征在于,所述旋转量确定模块包括:
将幅值加权方向分布直方图代入预设的最小化代价函数进行计算得到最优平移量,所述最优平移量即为三维旋转量。
CN201410708512.3A 2014-11-27 2014-11-27 基于梯度直方图分布匹配确定三维旋转量的方法和系统 Active CN104537640B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201410708512.3A CN104537640B (zh) 2014-11-27 2014-11-27 基于梯度直方图分布匹配确定三维旋转量的方法和系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201410708512.3A CN104537640B (zh) 2014-11-27 2014-11-27 基于梯度直方图分布匹配确定三维旋转量的方法和系统

Publications (2)

Publication Number Publication Date
CN104537640A CN104537640A (zh) 2015-04-22
CN104537640B true CN104537640B (zh) 2018-01-09

Family

ID=52853159

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201410708512.3A Active CN104537640B (zh) 2014-11-27 2014-11-27 基于梯度直方图分布匹配确定三维旋转量的方法和系统

Country Status (1)

Country Link
CN (1) CN104537640B (zh)

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103679713A (zh) * 2013-12-04 2014-03-26 华南理工大学 一种针对部分匹配图像的二维图像配准方法
CN104021547A (zh) * 2014-05-17 2014-09-03 清华大学深圳研究生院 肺部 ct 的三维配准方法

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103679713A (zh) * 2013-12-04 2014-03-26 华南理工大学 一种针对部分匹配图像的二维图像配准方法
CN104021547A (zh) * 2014-05-17 2014-09-03 清华大学深圳研究生院 肺部 ct 的三维配准方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
A Novel Technique for Prealignment in Multimodality Medical Image Registration;Wu Zhou et al;《Hindawi Publishing Corporation BioMed Research International》;20140522;第2部分,图1 *
Real-time estimation of prostate tumor rotation and translation with a kV imaging system based on an iterative closest point algorithm;Joubin Nasehi Tehrani et al;《PHYSICS IN MEDICINE AND BIOLOGY》;20131115;第2.4节 *
Volume Registration using the 3-D pseudo-polar Fourier transform;Yosi Keller;《IEEE Trans Signal Proces》;20061231;第54卷(第11期);第1部分 *

Also Published As

Publication number Publication date
CN104537640A (zh) 2015-04-22

Similar Documents

Publication Publication Date Title
Stemkens et al. Image-driven, model-based 3D abdominal motion estimation for MR-guided radiotherapy
Ceritoglu et al. Multi-contrast large deformation diffeomorphic metric mapping for diffusion tensor imaging
US8908944B2 (en) Information processing apparatus, information processing method, and program
CN107123137B (zh) 医学图像处理方法及设备
US10083506B2 (en) Whole body image registration method and method for analyzing images thereof
EP2422318B1 (en) Quantification of medical image data
US8417005B1 (en) Method for automatic three-dimensional segmentation of magnetic resonance images
JP5977214B2 (ja) 画像処理方法および装置並びにプログラム
JP6541334B2 (ja) 画像処理装置、画像処理方法、およびプログラム
CN102622759A (zh) 一种结合灰度与几何信息的医学图像配准方法
Yap et al. TIMER: Tensor image morphing for elastic registration
WO2014066001A2 (en) System and method for diagnosis of focal cortical dysplasia
WO2014176154A1 (en) System and method for image intensity bias estimation and tissue segmentation
Liu et al. Automatic localization of the anterior commissure, posterior commissure, and midsagittal plane in MRI scans using regression forests
Alam et al. Medical image registration: Classification, applications and issues
Richards et al. Evaluating methods for constructing average high-density electrode positions
US20160061923A1 (en) Sheet tractography using diffusion tensor mri
Garlapati et al. Towards measuring neuroimage misalignment
Ferrari et al. Detection of the midsagittal plane in MR images using a sheetness measure from eigenanalysis of local 3D phase congruency responses
CN104537640B (zh) 基于梯度直方图分布匹配确定三维旋转量的方法和系统
Ruppertshofen et al. Multi-level approach for the discriminative generalized hough transform
Prakash et al. Morphologic relationship among the corpus callosum, fornix, anterior commissure, and posterior commissure: MRI-based variability study
Schwing et al. Reliable extraction of the mid-sagittal plane in 3D brain MRI via hierarchical landmark detection
Galdames et al. Registration of renal SPECT and 2.5 D US images
Pinto et al. Initialization of deformable models in 3D magnetic resonance images guided by automatically detected phase congruency point landmarks

Legal Events

Date Code Title Description
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant