CN1952980A - 三维图像的快速配准方法 - Google Patents

Abstract

本发明公开了一种三维图像的快速配准方法，该方法包括以下步骤：1.分别读入待配准的目标图像和浮动图像；2.分别对三维图像数据进行灰度二值分割，得到图像中待配准物体的向量表示；3.计算物体的质心和协方差矩阵，得到主成分变换表达式和当量子午面；4.利用主成分变换将图像中的物体分别变换到原始图像的中间层，完成粗配准；5.将一幅图像作为目标图像固定，另外一幅作为浮动图像，以变换后的目标图像的当量子午面与目标图像的中间层的互信息量为目标函数，采用Powell算法搜索目标函数的最大值，当互信息量达到最大时，完成精配准。本发明将三维数据的配准减化为二维数据的配准，在保证精度的前提下，减少了配准所需时间。

Description

三维图像的快速配准方法

技术领域

本发明涉及一种在图像平面内的图形图像转换，具体涉及三维图像的快速配准方法，该方法适用于医学、气象卫星、航空航等图形图像片的配准。

背景技术

医学图像配准是指对于一幅医学图像寻求一种(或一系列)空间变换，使它与另一幅图像上的对应点达到空间上的一致。这种一致是指人体上的同一解剖点在两张图像上有相同的空间位置，配准的结果应使两幅图像上所有的解剖点匹配。多模态医学图像配准是目前医学图像处理领域中的一个热点问题，对于临床诊断和治疗有重要意义。医学图像可以提供有病变组织或器官的大小、形状、空间关系等详细信息，比如CT图像可以显示骨骼结构和组织密度分布情况；MR图像和超声(US)图像提供的则是软组织的信息；PET，SPECT能反映人体的功能和代谢信息。单一模态的图像往往不能提供医生所需要的足够的信息，通常要将不同模态的图像融合在一起，得到更丰富的信息以便了解病变组织或器官的综合情况，从而做出准确的诊断或制订合适的治疗方案。

目前已有多种图像配准方法，大体可以分为基于图像特征和基于图像灰度的配准方法。基于图像特征的方法一般通过寻找图像中比较明显的解剖结构特征来计算变换参数，提取的特征包括点、线、边缘和轮廓等。该方法计算效率比较高。特别是基于轮廓的方法，尽管不同模态图像像素的灰度分布特性之间有很大的差异，但物体的一些明显轮廓在两幅图像中均能得较好的保持，这些轮廓特征可以用来作为参照来进行多模图像间的配准。但是基于图像特征的方法的配准精度取决于特征提取的准确性与否。大部分情况下，医学图像特征点的位置比较复杂，很难进行准确地提取。基于图像灰度的配准方法不需要提取特征点，而是利用图像的灰度信息进行配准。由于充分利用了全部的灰度信息，因此一般更容易得到较精确的配准结果。其中互相关法(Cross-Correlation，CC)、最小平方差法(least square method，LSM)是几种常用的方法。互相关法一般实现起来简单，但是计算代价庞大；最小平方差法对图像数据的部分缺失和图像的背景噪声不是非常敏感，具有较强的鲁棒性，但由于需要通过迭代搜寻最优参数，故其执行效率也不是很高，而且其解严重依赖于初始值的选择，不恰当的初始值往往会使解陷入局部极值。

互信息(Mutual Information，MI)是信息论中的一个测度，用来度量两个随机变量之间的相似性。在多模态医学图像配准问题中，虽然两幅图像来源于不同的成像设备，但是它们基于人体共同的解剖信息，所以当两幅图像的空间位置完全一致时，它们的对应象素的灰度互信息达到最大值。作为一种相似性测度，互信息量取得了巨大的成功，特别在医学图像领域。互信息配准方法目前被公认为是配准精度和鲁棒性最好的回溯性配准方法之一。但传统的基于互信息量的方法需要考虑整个三维数据的信息，计算复杂度大，无法满足临床需要。能否综合基于特征方法和基于灰度的优点而设计出一种既快又准的方法一直是国内外学者努力的方向。

发明内容

本发明的目的是提供一种三维图像的快速配准方法，该方法具有处理图像数据信息量小，配准速度的优点。

本发明实现上述目的的技术解决方案是：

一种三维图像的快速配准方法，该方法由以下步骤组成：

1)分别读入待配准的目标图像和浮动图像，通过线性变换将二者所有像素点的灰度值变换到0～255的范围内；

2)利用FCM算法对步骤1所得到目标图像和浮动图像的灰度值进行二值粗分割，得到目标图像和浮动图像所示物体的向量表示；

3)分别计算目标图像和浮动图像所示物体的质心和协方差矩阵，并对协方差矩阵进行特征值分解，得到特征向量矩阵E和特征值矩阵V，得到主成分变换表达式；再以最大特征值对应的特征向量所指方向为目标图像和浮动图像所示物体的第一主轴方向，次之为第二主轴方向，再次之为第三主轴方向，然后令第一主轴和第二主轴所确定的平面为当量子午面(Equivalent Meridian Plane，EMP)，得到当量子午面；

4)采用双线性插值(Bilinear Interpolation)方法对所到的主成分变换表达式进行主成分变换(Principal Component Transform，PCT)，将目标图像和浮动图像所示物体的质心分别变换到目标图像和浮动图像的中心，并使当量子午面分别与其所在的目标图像和浮动图像的中间层重合，得到变换后的目标图像和浮动图像，完成粗配准；

5)固定变换后的目标图像，对变换后的浮动图像进行微小的调整，其调整方法是，以变换后的目标图像的当量子午面与目标图像的中间层的互信息量为目标函数，采用Powell算法搜索目标函数的最大值，当互信息量达到最大时，完成精配准。

本发明所述的一种三维图像快速配准方法，其中步骤2)的处理过程如下所述：利用FCM算法对步骤1所得到目标图像和浮动图像的灰度值进行二值粗分割，得到阈值T0后，先在阈值T∈[T₀-10，T₀+10]内取不同的T值对目标图像和浮动图像进行分割并计算分割图像与原图像的互信息量，再利用互信息量最大时的阈值T对图像进行二值精分割，得到目标图像和浮动图像所示物体的向量表示。该进一步改进方案使所得到的分割图像边界特征保持完好，虚假信息大大降低，图像边界细腻、连续且定位性能好，可大大提高配准的精度。

本发明所提供的快速配准方法将三维数据的配准减化为二维数据的配准，在保证精度的前提下，减少了配准所需时间。我们知道子午面是地理学中的一个重要概念，它是过地球南北两极所确定的轴线(地轴)的平面。推而广之，对于任意曲线，绕曲线端点所确定的轴线旋转360度都可以形成一个三维体，而过轴线的平面即为该三维体的子午面。对于医学图像，由成像设备产生的断层序列图像也可以重建为三维体，但该三维体一般是非常不规则的，故不能求取一般意义上的子午面。为此本发明提出了当量子午面的概念，即：对于三维物体，分别求出其第一主轴和与之正交的第二主轴，而由第一主轴和第二主轴所确定的平面称之为当量子午面。由于不同模态的图像都有相似的外轮廓，因而当量子午面也基本相同。当量子午面是三维数据中一个非常特别的平面，集合中所有点到该平面的距离的平方和最小，该平面也是唯一确定的。从而两个三维图像的配准可以简化为当量子午面的对准，即当量子午面互信息量最大时就说明图像完全对准。与计算整个三维数据的互信息量相比，这就大大减少了计算量。首先对图像进行二值化，得到待配准图像的向量表示；计算物体的质心和协方差矩阵，对协方差矩阵进行特征值分解，得到特征向量矩阵E和特征值矩阵V，得到当量子午面和主成分变换表达式；利用主成分变换完成粗配准，精细配准时只需要对浮动图像进行微小的调整计算当量子午面的互信息量，这就大大提高了配准速度，减少了陷入局部极值的可能。实验结果表明此方法能准确，快速地处理图像刚性配准问题，特别适用于三维图像的配准。

本发明综合了基于特征和基于灰度配准算法快和准的优点，将图像看作一个数据集合，首次提出了当量子午面的概念。并利用主成分分析进行粗配准；精细配准时只需要对浮动图像进行较小的调整计算中间层的互信息量，即将三维数据的配准简化为二维数据的配准，这就大大提高了配准速度，减少了陷入局部极值的可能。本发明在无需人工干预和预处理条件下能自动实现对多模态图像的自动、精确配准，且表现出良好的鲁棒性。

附图说明

图1为本发明所述方法的流程图；

图2为利用主成分变换将当量子午面与XY平面重合的示意图；

图3为人体头部磁共振(Magnetic Resonance，MR)图像和正电子放射断层(PositronEmission Tomography，PET)图像，其中图3(a)、图3(b)、图3(c)分别为MR图像的横断面、矢状面和冠状面图像，图3(d)、图3(e)、图3(f)分别为对应的PET图像。

图4为图3所示图像的分割结果，其中图4(a)、图4(b)、图4(c)分别为图3所示MR图像的横断面、矢状面和冠状面图像的二值精分割结果，图4(d)、图4(e)、图4(f)分别为对应的PET图像的横断面、矢状面和冠状面图像的二值精分割结果。

具体实施方式

下面结合一套头部的MR数据和PET数据(参见图3)，详细阐述本发明的工作步骤。

步骤1，读入MR图像和PET图像，分别采用中值滤波对输入图像进行预处理，以降低噪声对配准的影响。通过线性变换(I-I_min)×255/(I_max-I_min)将所有像素点的灰度值变换到0～255的范围内，其中I为图像灰度值。

步骤2，利用FCM算法对MR图像的灰度值进行二值粗分割得到阈值T0后，先取T∈[T₀-10，T₀+10]，以不同的T值对MR图像进行分割，并计算分割后的MR图像与原MR图像的互信息量，再以互信息量最大时的T值对MR图像进行二值精分割。其中，互信息量计算公式为：

$\mathrm{MI}=\underset{x,y}{\mathrm{\Σ}}p(x,y)\log \frac{p(x,y)}{p\left(x\right)p\left(y\right)},$ 其中分割后的MR图像与原MR图像灰度值出现的联合概率P(x，y)可以用联合直方图估计得到，边际函数P(x)可以分别对P(x，y)的各列求和得到，边际函数P(y)可以分别对P(x，y)的各行求和得到。按照同样的方法可以对PET图像进行二值精分割。MR图像和PET图像的分割结果见图4。从而得到MR图像和PET图像中物体的向量表示X^MR＝{(x_i，y_i，z_i)^T|i＝1，…，n}和X^PET＝{(x_i，y_i，z_i)^T|i＝1，…，n}。

步骤3，计算MR图像中物体的质心 $u^{\mathrm{MR}}=\frac{l}{n}\underset{i=l}{\overset{n}{\mathrm{\Σ}}}{X}_i^{\mathrm{MR}}$ 和协方差矩阵

$C^{\mathrm{MR}}=\frac{l}{n}\underset{i=l}{\overset{n}{\mathrm{\Σ}}}{X}_i^{\mathrm{MR}}{\left(X_i^{\mathrm{MR}}\right)}^T-u^{\mathrm{MR}}{\left(u^{\mathrm{MR}}\right)}^T,$ 对协方差矩阵进行特征值分解，得到特征向量矩阵E^MR和特征值矩阵V^MR。最人特征值对应的特征向量所指方向即为MR图像中物体的第一主轴方向，依次为第二主轴方向和第三主轴方向。第一主轴和第二主轴所确定的平面即为当量子午面，这样就得到了MR图像中物体的当量子午面和主成分变换表达式Y^MR＝E^MR(X^MR-u^MR)；同样计算PET图像中物体的质心 $u^{\mathrm{PET}}=\frac{l}{n}\underset{i=l}{\overset{n}{\mathrm{\Σ}}}{X}_i^{\mathrm{PET}}$ 和协方差矩阵

$C^{\mathrm{PET}}=\frac{l}{n}\underset{i=l}{\overset{n}{\mathrm{\Σ}}}{X}_i^{\mathrm{PET}}{\left(X_i^{\mathrm{PET}}\right)}^T-u^{\mathrm{PET}}{\left(u^{\mathrm{PET}}\right)}^T,$ 对协方差矩阵进行特征值分解，得到特征向量矩阵E^PET和特征值矩阵V^PET，最大特征值对应的特征向量所指方向即为PET图像中物体的第一主轴方向，依次为第二主轴方向和第三主轴方向。第一主轴和第二主轴所确定的平面即为当量子午面，这样就得到了PET图像中物体的当量子午面和主成分变换表达式Y^PET＝E^PET(X^PET-u^PET)。

步骤4，利用主成分变换表达式Y^MR＝E^MR(X^MR-u^MR)对MR图像进行主成分变换，采用双线性插值方法进行插值，将目标图像所示物体的质心变换到目标图像的中心，并使当量子午面与MR图像的中间层重合，得到变换后的图像MR’；再利用主成分变换表达式Y^PET＝E^PET(X^PET-u^PET)对PET图像进行主成分变换，采用双线性插值方法进行插值，将浮动图像所示物体的质心变换到浮动图像的的中心，并使当量子午面与PET图像的中间层重合，得到变换后的图像PET’。这样，MR图像中的物体的当量子午面与PET图像中的物体的当量子午曲重合，完成粗配准。

步骤5，将图像MR’作为目标图像固定不动，图像PET’作为浮动图像而进行微小的调整。具体调整方法是：以变换后的目标图像MR’的当量子午面与目标图像MR的中间层的互信息量为目标函数，采用Powell算法搜索目标函数的最大值。互信息量计算公式为：

$\mathrm{MI}=\underset{x,y}{\mathrm{\Σ}}p(x,y)\log \frac{p(x,y)}{p\left(x\right)p\left(y\right)},$ 其中变换后的目标图像MR’与变换后的浮动图像PET’灰度值出现的联合概率P(x，y)可以用联合直方图估计得到，边际函数P(x)可以分别对P(x，y)的各列求和得到，边际函数P(y)可以分别对P(x，y)的各行求和得到。当目标函数达到最大时，两图像完全配准。

Claims (2)

1、一种三维图像快速配准方法，其特征是：

1)分别读入待配准的目标图像和浮动图像，通过线性变换将二者所有像素点的灰度值变换到0～255的范围内；

2)利用FCM算法对步骤1所得到目标图像和浮动图像的灰度值进行二值粗分割，得到目标图像和浮动图像所示物体的向量表示；

3)分别计算待目标图像和浮动图像所示物体的质心和协方差矩阵，并对协方差矩阵进行特征值分解，得到特征向量矩阵E和特征值矩阵V，得到主成分变换表达式；再以最大特征值对应的特征向量所指方向为目标图像和浮动图像所示物体的第一主轴方向，次之为第二主轴方向，再次之为第三主轴方向，然后令第一主轴和第二主轴所确定的平面为当量子午面，得到当量子午面；

4)采用双线性插值方法对所到的主成分变换表达式进行主成分变换，将目标图像和浮动图像所示物体的质心分别变换到目标图像和浮动图像的中心，并使当量子午面分别与其所在的目标图像和浮动图像的中间层重合，得到变换后的目标图像和浮动图像，完成粗配准；

5)固定变换后的目标图像，对变换后的浮动图像进行微小的调整，其调整方法是，以变换后的目标图像的当量子午面与目标图像的中间层的互信息量为目标函数，采用Powell算法搜索目标函数的最大值，当互信息量达到最大时，完成精配准。

2、根据权利要求1所述的一种三维图像快速配准方法，其特征是步骤2)的处理过程如下所述：利用FCM算法对步骤1所得到目标图像和浮动图像的灰度值进行二值粗分割，得到阈值T0后，先在阈值T∈[T₀-10，T₀+10]内取不同的T值对目标图像和浮动图像进行分割并计算分割图像与原图像的互信息量，再利用互信息量最大时的阈值T对图像进行二值精分割，得到目标图像和浮动图像所示物体的向量表示。

Images