CN1556501A - 图像配准方法 - Google Patents

Abstract

一种图像配准方法，属于图像处理技术领域。具体如下：步骤一：首先基于图像灰度计算出处理对象与背景之间的阈值；步骤二：根据步骤一确定的处理对象与背景的阈值，使用种子填充法去掉参考图像和待配准图像中背景部分；步骤三：采用快速相关法对去处背景后的参考图像和待配准图像进行粗配准；步骤四：使用最大互信息法对参考图像和由步骤三得出的图像进行精细配准。本发明极大地提高了整个配准方法的速度，使用梯度坡降法寻找最快加速方向，使用powell搜索算法解决梯度坡降法的局部极值问题，使用黄金分割的方法可以在最短时间内得到最优解。本发明方法适用于对多幅图像进行配准为进一步研究处理做准备。

Description

图像配准方法

技术领域

本发明涉及一种图像处理方法，具体是一种图像配准方法。属于图像处理技术领域。

技术背景

图像的配准是把来自不同模态或不同时间的多幅图像进行排列，使这些图像达到空间一致，为图像的进一步处理提供保证。如序列图像的三维重建和压缩，医学图像信息的融合。图像配准的目的就是确定两个图像之间的几何变换关系，通过这个变换可以把一个图像中的任意一点的坐标变换到另一图像中的对应位置的坐标。

经文献检索发现，作者B.Likar，F.Pernus在《Image and Vision Computing》第2001年19期的33-44页上发表文章“A hierarchical approach to elasticregistration based on mutual information”(《图像和视图计算》，“基于互信息的分级弹性配准”)，该文中提到图像配准的方法为：首先对图像进行分块，再分别在各个块图像上用互信息进行本地配准，然后通过插值，从新组成图像，最后进行全局配准。该项技术有以下缺点：(1)如果两幅图像在位置上差别较大，通过分块，对应的小块图像信息差别很大，本地配准结果就很难准确，通过插值，很难保证最后的全局配准结果的正确性。(2)该技术虽然对本地与全局采取不同精度的配准方法来提高配准速度，但如果图像差别较大的情况下，其计算速度依然不能保证。(3)该技术的最优化程序采用二分法查找，虽然能够找到最优变换，但其计算速度不快。

发明内容

本发明的目的在于克服背景技术中存在的不足，提供一种图像配准方法。该方法在使用最大互信息配准前使用最大相关性进行全局粗配准的方法，减少了计算量，提高了配准速度。在最优化程序的选择上，采用梯度坡降法和Powell搜索算法相结合的方法，在很大程度上解决了局部极值的问题，又提高了原Powell搜索算法的计算效率。

本发明通过以下技术方案实现的，本发明方法步骤如下：

步骤一：为了减少配准算法的计算量，使图像免除噪音的干扰，首先基于图像灰度计算出处理对象与背景之间的阈值。

步骤二：根据步骤一确定的处理对象与背景的阈值，使用种子填充法去掉参考图像和待配准图像中背景部分，这样可以避免背景灰度因素在以下步骤中可能产生的干扰。

步骤三：采用快速相关法对去处背景后的参考图像和待配准图像进行粗配准。

步骤四：使用最大互信息法对参考图像和由步骤三得出的图像进行精细配准。

以下对本发明内容作进一步的说明，具体内容如下：

所述的步骤三，采用快速相关法对去处背景后的参考图像和待配准图像进行粗配准具体为：

(1)通过快速相关法找到最优变换。

(2)在查找最优变换时，使用刚性变换(平移变换加旋转变换)。

(3)采用梯度坡降法和Powell搜索算法相结合的方法得出最大相似性。

(4)使用黄金分割查找一维最优变换。

所述的步骤四，使用最大互信息法对参考图像和由步骤三得出的图像进行配准，具体为：

(1)通过最大互信息方法找到最优变换。

(2)在查找最优变换时，使用刚性变换(平移变换加旋转变换)。

(3)采用梯度坡降法和Powell搜索算法相结合的方法得出最大互信息。

(4)使用黄金分割查找一维最优变换。

由于不同的成像设备有不同的成像特点，包含不同的信息，对图象配准技术的研究可以对不同成像设备在同一时期产生的图像信息进行综合，获得更多的可用信息；同时，通过对同一成像设备不同时期的图像配准问题，观察前后不同组成成分的影像学变化情况，得到有用信息。在对骨肉瘤诊断方面，骨肉瘤对化疗药物的敏感程度决定了患者预后，提高对化疗不敏感的骨肉瘤患者的疗效，是临床亟待解决的重要课题。通过对同一医学成像设备不同时期的图像配准问题，观察骨肉瘤化疗前后不同组成成分的影像学变化情况，以肿瘤细胞坏死率测定为金标准，可以很好地解决这一问题。尽管目前已提出了许多配准算法，但都各有其应用的局限性，还没有哪一种算法是通用的，实际用于临床的还不太多，且大多集中于对头部的研究，对骨肉瘤图像的配准应用还很少。本发明方法可以应用在这种场合。

本发明能对图像进行自动配准，主要具有以下特点和作用：(1)本发明在使用互信息进行精细配准前先使用快速相关法对图像进行粗配准。由于快速相关法计算量小，速度快，极大地提高了整个配准方法的速度。(2)采用梯度坡降法和Powell搜索算法相结合的方法得出最大相关性(最大互信息)。使用梯度坡降法寻找最快加速方向。使用powell搜索算法解决梯度坡降法的局部极值问题。(3)本发明使用黄金分割的方法在一维上查找变换T，使得待配准图像B通过T变换，和参考图像A具有最大相似性。由于对给定的一个初始区间，黄金分割是求解一维曲线方程最大值所在位置的最有效方法，因此使用黄金分割的方法可以在最短时间内得到最优解(T变换)。

本发明方法适用于对多幅图像进行配准为进一步研究处理做准备，如对核磁共振设备采集的时间序列图像信息的研究，计算机断层X线摄影图像和磁共振图像的信息融合。

具体实施方式

以下提供实施例来进一步说明本发明方法的应用：

MR图像为某某医院核磁共振影像科室利用马可尼医疗设备公司1.5T核磁共振影像设备(Eclipse 1.5T，Marconi Medical Systems，Inc，)采集的17岁男孩的腿部骨肉瘤的图像。

运行本发明方法的设备为普通个人计算机，其具体配置为P41.8G，256M内存，80G硬盘。发明方法的具体实现平台为Windows2000下的Visual C++计算平台。

具体实施情况如下：

步骤一：为了减少配准算法的计算量，使图像免除噪音的干扰，首先基于图像灰度计算出处理对象与背景之间的阈值。

(1)求出图像中的最大和最小灰度Z₁和Z_k，令阈值初始值为

$T_0=\frac{Z_1+Z_k}{2}$

(2)根据阈值T_k将图像分割成R₁和R₂两部分，分别求出两部分的平均灰度值Z₀和Z_B：

$Z_o=\frac{\underset{z(i,j)<T_k}{\mathrm{\&Sigma;}}z(i,j)\&times;N(i,j)}{\underset{z(i,j)<T_k}{\mathrm{\&Sigma;}}}$

$Z_B=\frac{\underset{z(i,j)>T_k}{\mathrm{\&Sigma;}}z(i,j)\&times;N(i,j)}{\underset{z(i,j)>T_k}{\mathrm{\&Sigma;}}z(i,j)}$

公式中Z(i，j)是图像上(i，j)点的灰度值，N(i，j)是(i，j)点的权重系数，这里取N(i，j)＝1.0。

(3)求出新的阈值

$T_{k+1}=\frac{Z_o+Z_B}{2}$

(4)如果T_k＝T_k+1，则结束，否则k＝k+1，迭代执行上述步骤。

步骤二：根据步骤一确定的处理对象与背景的阈值，使用种子填充法去掉参考图像和待配准图像中背景部分，这样可以避免背景灰度因素在以下步骤中可能产生的干扰。

步骤三：采用快速相关法对去处背景后的参考图像和待配准图像进行粗配准。

(1)通过快速相关法找到最优变换。

计算快速相关法测度的方法如下：

对参考图像A和待配准图像B进行逐行扫描(A(i，j)和B(i，j)为图像A和B的灰度值)：

若A(i，j)＝B(i，j)＝背景值，则设该点值为0；

若A(i，j)小于或者大于背景值且B(i，j)小于或者大于背景值，则设该点值为0；

若A(i，j)和B(i，j)中有一个为背景值，另一个非背景值，则设该点的值为1。

计算所有非零点的个数sum，作为参考图像和配准图像的相似性评估，sum越小，则说明两幅图像越相似。

(2)在查找最优变换时，使用刚性变换(平移变换加旋转变换)。

(3)采用梯度坡降法和Powell搜索算法相结合的方法得出最大相似性。

假设有N个相互独立的变量(这里执行的是刚性变换，则有3个变量，x轴变量，y轴变量和旋转变量θ)，则每次Powell迭代需要计算(N+1)次一维查找，每次查找都从起始点开始，前N次查找分别沿着这N个轴进行，第N+1次查找沿着从起始点到前N次查找的结束点的方向(即梯度方向)进行。当迭代找不到更大的相似性时就结束。

使用梯度坡降法来求互信息，一次迭代需要进行N次相似性评估。每次先把一个独立变量x_j增加一个很小的变量ε，求出信息量的增量Δ_j，然后把求得的N维向量{Δ₁，Δ₂，……，Δ_N}化为标准型，用来求得梯度方向。

(4)使用黄金分割查找一维最优变换。

本发明方法使用黄金分割的方法在一维上查找变换T，使得待配准图像B通过T变换，和参考图像A具有最大相似性。为了求得这个变换，假设相似性方程F(此方程是以变换T为变量的)沿着轴方向是一条向上抛物线(也就是有且仅有1个最大值)并且必须找到与这个最大值对应的变换T所在的区间。找到这个区间的方法如下：

a)令x_j-1＝x₀(x₀是起始点)

b)随机找一个数ε，沿着轴方向找到与x_j-1距离为ε的点x_j。

c)在x_j点上求得评估方程F的值，如果F(x_j)＜F(x_j-1)，则停止，否则令x_j-1＝x_j，沿着轴方向找到与点x_j-1距离为ε的点x_j，重复上述步骤，直到方程找到点x_j，使得F(x_j)＜F(x_j-1)。

如果此过程执行3次以上，则只有最后3个点a_j-2，x_j-1，x_j被保存。区间(x_j-2，x_j)即为所求区间。

符合条件的区间找到后，开始在此区间进行黄金分割，求得似的相关性(互信息)为最大值所在的点。步骤如下：

a)令k＝1.0e-4

b)p1＝x_j-2，p4＝x_j，

c)s＝p4-p1，p2＝p1+0.382s，p3＝0.618s

d)如果F(p4)＞F(p1)，则p1＝p2；否则，p4＝p3

e)如果(p4-p1)≥k，则转c)，否则迭代结束，如果F(p4)＞F(p1)，则p4为所求点，否则p1为所求点。

步骤四：使用最大互信息法对参考图像和由步骤三得出的图像进行精细配准。(具体实施方法同步骤三，不同的是这里的相似性测度是互信息)。

(1)在查找最优变换时，使用刚性变换(平移变换加旋转变换)。

(2)采用梯度坡降法和Powell搜索算法相结合的方法得出最大互信息。

(3)使用黄金分割查找一维最优变换。

(4)对待配准图像实施通过最大互信息方法找到的刚性变换，得到最终的准结果。

本发明方法取得了较好的实施效果：(1)对核磁共振图像进行配准，得到较的效果(见表1)；(2)对于骨肉瘤化疗前后不同组成成分的图像进行配准，得出变化情况，为进一步的诊断提供了依据。

                   表1：仿真实验的配准参数

	旋转变量θ	x轴变量	y轴变量
				快速相关法配准结果	14.8060	-2.2539	2.8416
最大互信息配准结果	0.1335	-0.0570	0.5219
				最终结果	14.99395	-2.3109	3.3635
真实结果	15	-2.5	3

Claims (5)

1、一种图像配准方法，其特征在于，方法步骤如下：

步骤一：首先基于图像灰度计算出处理对象与背景之间的阈值；

步骤二：根据步骤一确定的处理对象与背景的阈值，使用种子填充法去掉参考图像和待配准图像中背景部分；

步骤三：采用快速相关法对去处背景后的参考图像和待配准图像进行粗配准；

步骤四：使用最大互信息法对参考图像和由步骤三得出的图像进行精细配准。

2、根据权利要求1所述的图像配准方法，其特征是，所述的步骤三，采用快速相关法对去处背景后的参考图像和待配准图像进行粗配准，具体为：

(1)通过快速相关法找到最优变换；

(2)使用刚性变换即平移变换加旋转变换查找最优变换；

(3)采用梯度坡降法和Powell搜索算法相结合的方法得出最大相似性；

(4)使用黄金分割查找一维最优变换。

3、根据权利要求1或2所述的图像配准方法，其特征是，所述的快速相关法，计算其测度的方法如下：

对参考图像A和待配准图像B进行逐行扫描，A(i，j)和B(i，j)为图像A和B的灰度值：

若A(i，j)＝B(i，j)＝背景值，则设该点值为0；

若A(i，j)小于或者大于背景值且B(i，j)小于或者大于背景值，则设该点值为0；

若A(i，j)和B(i，j)中有一个为背景值，另一个非背景值，则设该点的值为1；

计算所有非零点的个数sum，作为参考图像和配准图像的相似性评估，sum越小，则说明两幅图像越相似。

4、根据权利要求1所述的图像配准方法，其特征是，所述的步骤四，使用最大互信息法对参考图像和由步骤三得出的图像进行精细配准，具体为：

(1)通过最大互信息方法找到最优变换；

(2)使用刚性变换即平移变换加旋转变换查找最优变换；

(3)采用梯度坡降法和Powell搜索算法相结合的方法得出最大互信息；

(4)使用黄金分割查找一维最优变换。

5、根据权利要求2或4所述的图像配准方法，其特征是，所述的使用黄金分割查找一维最优变换，具体如下：

使用黄金分割的方法在一维上查找变换T，使待配准图像B通过T变换，和参考图像A具有最大相似性，假设相似性方程F沿着轴方向是一条向上抛物线，此方是以变换T为变量，为了求得一维最优变换，必须找到与这个最大值对应的变换T所在的区间，找到这个区间的方法如下：

d)令x_j-1＝x₀，x₀是起始点，

e)随机找一个数ε，沿着轴方向找到与x_j-1距离为ε的点x_j，

f)在x_j点上求得评估方程F的值，如果F(x_j)＜F(x_j-1)，则停止，否则令x_j-1＝x_j，沿着轴方向找到与点x_j-1距离为ε的点x_j，重复上述步骤，直到方程找到点x_j，使得F(x_j)＜F(x_j-1)，

如果此过程执行3次以上，则只有最后3个点x_j-2，x_j-1，x_j被保存，区间(x_j-2，x_j)即为所求区间；

  然后在此区间进行黄金分割，求得相似的相关性或互信息为最大值所在的点，步骤如下：

  f)令k＝1.0e-4，

  g)p1＝x_j-2，p4＝x_j，

  h)s＝p4-p1，p2＝p1+0.382s，p3＝0.618s，

  i)如果F(p4)＞F(p1)，则p1＝p2，否则，p4＝p3，

  如果(p4-p1)≥k，则转c)，否则迭代结束，如果F(p4)＞F(p1)，则p4为所求点，否则p1为所求点。