时间戳

首页 > 解决方案 > 即使其中一张图像中的对象略微旋转，如何找到两张图像之间的差异

python - 即使其中一张图像中的对象略微旋转，如何找到两张图像之间的差异

问题描述

给出两个输入图像：i）Master.png（主图像）

ii）query2.png（图像内的旋转对象）

对齐图像后，我得到这个： 与主图像相比，图像有点模糊，并且也形成了黑色边缘

匹配对齐和主图像后的输出是错误的。

如何找到并突出对齐图像和主图像之间的差异。

标签： pythonopencvimage-processing

解决方案

对齐的过程称为注册。您所说的主图像称为Fixed Image。您对齐的图像称为Moving Image。有两种方法可以直观地评估配准：差分法和棋盘格法。在下面的示例中，我将展示使用 ITK 的翻译示例的视觉注册评估。旋转也类似。您可以对 OpenCV 使用类似的概念。您可以忽略注册部分，因为您已经完成了它，但为了完整起见在此处显示。您可以直接转到下面的“两个图像之间的差异”。理想情况下，您应该发布自己的代码 - 您如何注册，以便其他人可以从那里继续。

import itk
from matplotlib import pyplot as plt
import numpy as np
from skimage.morphology import square, closing
from skimage.color import rgb2gray
from skimage.io import imread
from skimage import img_as_ubyte
%matplotlib inline

让我们先看看原始图像

path1 = '../data/3_T1.png'
path2 = '../data/3_T2.png'

PixelType = itk.US    # unsigned short
 
PixelType = itk.F 
Dimension = 2
FixedImageType  = itk.Image[PixelType, Dimension]
MovingImageType = itk.Image[PixelType, Dimension]

FixedImageReaderType  = itk.ImageFileReader[FixedImageType]
MovingImageReaderType = itk.ImageFileReader[MovingImageType]

fixedImageReader  = FixedImageReaderType.New()
movingImageReader = MovingImageReaderType.New()

fixedImageFile = path1
movingImageFile = path2

fixedImageReader.SetFileName(fixedImageFile)
movingImageReader.SetFileName(movingImageFile)
    
fig = plt.figure(figsize=(10, 6))
ax1 = fig.add_subplot(1,2,1)
ax1.imshow(itk.GetArrayViewFromImage(fixedImageReader.GetOutput()), cmap=plt.get_cmap('gray'))
ax1.axis('off')
ax1.set_title('Fixed Image (T1)')
ax2 = fig.add_subplot(1,2,2)
ax2.imshow(itk.GetArrayViewFromImage(movingImageReader.GetOutput()), cmap=plt.get_cmap('gray'))
ax2.set_title('Moving Image (T2)')
ax2.axis('off')
plt.show()

在此处输入图像描述

注册参数和设置

InternalPixelType = itk.F
InternalImageType = itk.Image[InternalPixelType, Dimension]

TransformType    = itk.TranslationTransform[itk.D, Dimension]
OptimizerType    = itk.GradientDescentOptimizer
InterpolatorType = itk.LinearInterpolateImageFunction[InternalImageType, itk.D]
RegistrationType = itk.ImageRegistrationMethod[InternalImageType, InternalImageType]  
MetricType       = itk.MutualInformationImageToImageMetric[InternalImageType, InternalImageType]

transform     = TransformType.New()
optimizer     = OptimizerType.New()
interpolator  = InterpolatorType.New()
registration  = RegistrationType.New()
metric        = MetricType.New()

registration.SetOptimizer(optimizer)
registration.SetTransform(transform)
registration.SetInterpolator(interpolator)
registration.SetMetric(metric)

归一化两个输入图像的统计分布以简化互信息的计算，即重新调整输入图像的强度以产生具有零均值和单位方差的输出图像

FixedNormalizeFilterType  = itk.NormalizeImageFilter[FixedImageType,  InternalImageType]
MovingNormalizeFilterType = itk.NormalizeImageFilter[MovingImageType, InternalImageType]

fixedNormalizer  = FixedNormalizeFilterType.New()
movingNormalizer = MovingNormalizeFilterType.New()

此外，要配准的图像的低通滤波也将增加对噪声的鲁棒性。

GaussianFilterType = itk.DiscreteGaussianImageFilter[ InternalImageType, InternalImageType]
fixedSmoother  = GaussianFilterType.New()
movingSmoother = GaussianFilterType.New()

fixedSmoother.SetVariance( 2.0 )
movingSmoother.SetVariance( 2.0 )

# standard deviation of the Gaussian kernel
metric.SetFixedImageStandardDeviation(0.4)
metric.SetMovingImageStandardDeviation(0.4)

fixedNormalizer.SetInput(  fixedImageReader.GetOutput() )
movingNormalizer.SetInput( movingImageReader.GetOutput() )

fixedSmoother.SetInput( fixedNormalizer.GetOutput() )
movingSmoother.SetInput( movingNormalizer.GetOutput() )

设置用于注册的输入图像

registration.SetFixedImage(    fixedSmoother.GetOutput()    )
registration.SetMovingImage(   movingSmoother.GetOutput()   )
fixedNormalizer.Update()
fixedImageRegion = fixedNormalizer.GetOutput().GetBufferedRegion()
registration.SetFixedImageRegion( fixedImageRegion )

初始化

initialParameters = transform.GetParameters()
initialParameters[0] = 0.0 # Initial offset in mm along X
initialParameters[0] = 0.0 # Initial offset in mm along Y

registration.SetInitialTransformParameters( initialParameters )

我们现在应该定义度量计算中要考虑的空间样本的数量。我们推迟了这个设置，直到我们完成了图像的预处理，因为样本的数量通常被定义为固定图像中像素总数的一小部分。

空间样本的数量通常可以低至固定图像中像素总数的 $1%$。增加样本数量可以提高度量从一次迭代到另一次迭代的平滑度，因此在将此度量与依赖于度量值连续性的优化器结合使用时会有所帮助。当然，权衡是更多的样本会导致每次度量评估的计算时间更长。

numberOfPixels = fixedImageRegion.GetNumberOfPixels()
numberOfSamples = numberOfPixels * 0.01
numberOfSamples = int(numberOfSamples)
metric.SetNumberOfSpatialSamples( numberOfSamples )
metric.ReinitializeSeed( 121212 ) # For consistent results when regression testing

优化器设置 - 由于互信息的较大值比较小的值表示更好的匹配，我们需要最大化成本函数

optimizer.SetNumberOfIterations( 200 )
optimizer.MaximizeOn()
optimizer.SetLearningRate( 15.0 )

执行注册

registration.Update()

提取结果

finalParameters = registration.GetLastTransformParameters()
TranslationAlongX = finalParameters[0]
TranslationAlongY = finalParameters[1]
numberOfIterations = optimizer.GetCurrentIteration()
bestValue = optimizer.GetValue()

print("Result = ")
print(" Translation X = " + str(TranslationAlongX))
print(" Translation Y = " + str(TranslationAlongY))
print(" Iterations    = " + str(numberOfIterations))
print(" Metric value  = " + str(bestValue))
print(" Numb. Samples = " + str(numberOfSamples))

输出：

Result = 
 Translation X = 10.029903702181082
 Translation Y = 9.982473996881163
 Iterations    = 200
 Metric value  = 0.5775022237597298
 Numb. Samples = 530

转换：

finalTransform = TransformType.New()
finalTransform.SetParameters( finalParameters )
finalTransform.SetFixedParameters( transform.GetFixedParameters() )

重采样

ResampleFilterType = itk.ResampleImageFilter[MovingImageType, FixedImageType]
resample = ResampleFilterType.New()
resample.SetTransform( finalTransform )
resample.SetInput( movingImageReader.GetOutput() )

fixedImage = fixedImageReader.GetOutput()

resample.SetSize(    fixedImage.GetLargestPossibleRegion().GetSize() )
resample.SetOutputOrigin(  fixedImage.GetOrigin() )
resample.SetOutputSpacing( fixedImage.GetSpacing() )
resample.SetOutputDirection( fixedImage.GetDirection() )
resample.SetDefaultPixelValue( 100 )

OutputPixelType = itk.ctype('unsigned char')
OutputImageType = itk.Image[OutputPixelType, Dimension]

resample.Update()

两个图像之间的差异

difference = itk.SubtractImageFilter.New(Input1=fixedImage, Input2=resample)
resample.SetDefaultPixelValue(1)
fig = plt.figure(figsize=(10, 12))

ax1 = fig.add_subplot(1,3,1)
ax1.imshow(itk.GetArrayViewFromImage(resample), cmap=plt.get_cmap('gray'))
ax1.axis('off')
ax1.set_title('Registered Image')

ax3 = fig.add_subplot(1,3,2)
ax3.imshow(itk.GetArrayViewFromImage(difference), cmap=plt.get_cmap('gray'))
ax3.set_title('Difference After Registration')
ax3.axis('off')

identityTransform = TransformType.New()
identityTransform.SetIdentity()
resample.SetTransform(identityTransform)
difference.Update()

ax4 = fig.add_subplot(1,3,3)
ax4.imshow(itk.GetArrayViewFromImage(difference), cmap=plt.get_cmap('gray'))
ax4.set_title('Difference Before Registration')
ax4.axis('off')

plt.show()

在此处输入图像描述

固定/移动图像棋盘

CheckerBoardFilterType = itk.CheckerBoardImageFilter[FixedImageType]
checker = CheckerBoardFilterType.New()
checker.SetInput1( fixedImage )
checker.SetInput2( resample.GetOutput() )

fig = plt.figure(figsize=(10, 6))
ax1 = fig.add_subplot(1,2,1)
ax1.imshow(itk.GetArrayViewFromImage(checker), cmap=plt.get_cmap('gray'))
ax1.axis('off')
ax1.set_title('Before Registration')

resample.SetTransform( finalTransform )
checker.SetInput2( resample.GetOutput() )

ax2 = fig.add_subplot(1,2,2)
ax2.imshow(itk.GetArrayViewFromImage(checker), cmap=plt.get_cmap('gray'))
ax2.set_title('After Registration')
ax2.axis('off')
plt.show()

在此处输入图像描述

推荐阅读