Squared Difference of Two Images#

Synopsis#

Compute the squared difference of corresponding pixels in two images.

Results#

Code#

C++#

#include "itkImage.h"
#include "itkSquaredDifferenceImageFilter.h"
#include "itkImageRegionIterator.h"

using UnsignedCharImageType = itk::Image<unsigned char, 2>;
using FloatImageType = itk::Image<float, 2>;

static void
CreateImage1(UnsignedCharImageType::Pointer image);
static void
CreateImage2(UnsignedCharImageType::Pointer image);

int
main()
{
  auto image1 = UnsignedCharImageType::New();
  CreateImage1(image1);

  auto image2 = UnsignedCharImageType::New();
  CreateImage2(image2);

  using SquaredDifferenceImageFilterType =
    itk::SquaredDifferenceImageFilter<UnsignedCharImageType, UnsignedCharImageType, FloatImageType>;

  auto squaredDifferenceFilter = SquaredDifferenceImageFilterType::New();
  squaredDifferenceFilter->SetInput1(image1);
  squaredDifferenceFilter->SetInput2(image2);
  squaredDifferenceFilter->Update();

  return EXIT_SUCCESS;
}

void
CreateImage1(UnsignedCharImageType::Pointer image)
{
  UnsignedCharImageType::IndexType start;
  start.Fill(0);

  UnsignedCharImageType::SizeType size;
  size.Fill(10);

  UnsignedCharImageType::RegionType region;
  region.SetSize(size);
  region.SetIndex(start);

  image->SetRegions(region);
  image->Allocate();

  itk::ImageRegionIterator<UnsignedCharImageType> imageIterator(image, region);

  while (!imageIterator.IsAtEnd())
  {
    imageIterator.Set(255);

    ++imageIterator;
  }
}


void
CreateImage2(UnsignedCharImageType::Pointer image)
{
  // Create an image with 2 connected components
  UnsignedCharImageType::IndexType start;
  start.Fill(0);

  UnsignedCharImageType::SizeType size;
  size.Fill(10);

  UnsignedCharImageType::RegionType region;
  region.SetSize(size);
  region.SetIndex(start);

  image->SetRegions(region);
  image->Allocate();

  itk::ImageRegionIterator<UnsignedCharImageType> imageIterator(image, region);

  while (!imageIterator.IsAtEnd())
  {
    imageIterator.Set(100);

    ++imageIterator;
  }
}

Classes demonstrated#

template<typename TInputImage1, typename TInputImage2, typename TOutputImage>
class SquaredDifferenceImageFilter : public itk::BinaryGeneratorImageFilter<TInputImage1, TInputImage2, TOutputImage>

Implements pixel-wise the computation of squared difference.

This filter is parameterized over the types of the two input images and the type of the output image.

Numeric conversions (castings) are done by the C++ defaults.

The filter will walk over all the pixels in the two input images, and for each one of them it will do the following:

  • cast the input 1 pixel value to double

  • cast the input 2 pixel value to double

  • compute the difference of the two pixel values

  • compute the square of the difference

  • cast the double value resulting from sqr() to the pixel type of the output image

  • store the casted value into the output image.

The filter expect all images to have the same dimension (e.g. all 2D, or all 3D, or all ND)

ITK Sphinx Examples:

See itk::SquaredDifferenceImageFilter for additional documentation.