Process Image Chunks#

Synopsis#

Process an image in small chunks.

Results#

Input image

Input image#

Output image

Output image#

Code#

C++#

#include "itkImageFileReader.h"
#include "itkImageFileWriter.h"
#include "itkMedianImageFilter.h"

int
main(int argc, char * argv[])
{
  if (argc != 3)
  {
    std::cerr << "Usage: " << std::endl;
    std::cerr << argv[0];
    std::cerr << " <InputFileName> <OutputFileName>";
    std::cerr << std::endl;
    return EXIT_FAILURE;
  }

  const char * inputFileName = argv[1];
  const char * outputFileName = argv[2];

  const int xDiv = 6;
  const int yDiv = 4;
  const int zDiv = 1; // 1 for 2D input

  constexpr unsigned int Dimension = 3;

  using PixelType = unsigned char;
  using ImageType = itk::Image<PixelType, Dimension>;

  const auto input = itk::ReadImage<ImageType>(inputFileName);

  ImageType::RegionType fullRegion = input->GetLargestPossibleRegion();
  ImageType::SizeType   fullSize = fullRegion.GetSize();
  // when reading image from file, start index is always 0

  ImageType::IndexType start;
  ImageType::IndexType end;
  ImageType::SizeType  size;

  using MedianType = itk::MedianImageFilter<ImageType, ImageType>;
  auto median = MedianType::New();
  median->SetInput(input);
  median->SetRadius(2);

  using WriterType = itk::ImageFileWriter<ImageType>;
  auto writer = WriterType::New();
  writer->SetFileName(outputFileName);

  // Use for loops to split the image into chunks.
  // This way of splitting gives us easy control over it.
  // For example, we could make the regions always be of equal size
  // in order to create samples for a deep learning algorithm.
  // ImageRegionSplitterMultidimensional has a similar
  // functionality to what is implemented below.
  for (int z = 0; z < zDiv; ++z)
  {
    start[2] = fullSize[2] * double(z) / zDiv;
    end[2] = fullSize[2] * (z + 1.0) / zDiv;
    size[2] = end[2] - start[2];

    for (int y = 0; y < yDiv; ++y)
    {
      start[1] = fullSize[1] * double(y) / yDiv;
      end[1] = fullSize[1] * (y + 1.0) / yDiv;
      size[1] = end[1] - start[1];

      for (int x = 0; x < xDiv; ++x)
      {
        start[0] = fullSize[0] * double(x) / xDiv;
        end[0] = fullSize[0] * (x + 1.0) / xDiv;
        size[0] = end[0] - start[0];

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

        // now that we have our chunk, request the filter to only process that
        median->GetOutput()->SetRequestedRegion(region);
        median->Update();
        ImageType::Pointer result = median->GetOutput();

        // only needed in case of further manual processing
        result->DisconnectPipeline();

        // possible additional non-ITK pipeline processing

        writer->SetInput(result);

        // convert region into IO region
        itk::ImageIORegion ioRegion(Dimension);
        ioRegion.SetIndex(0, start[0]);
        ioRegion.SetIndex(1, start[1]);
        ioRegion.SetIndex(2, start[2]);
        ioRegion.SetSize(0, region.GetSize()[0]);
        ioRegion.SetSize(1, region.GetSize()[1]);
        ioRegion.SetSize(2, region.GetSize()[2]);
        // tell the writer this is only a chunk of a larger image
        writer->SetIORegion(ioRegion);

        try
        {
          writer->Update();
        }
        catch (const itk::ExceptionObject & error)
        {
          std::cerr << "Exception for chunk: " << x << ' ' << y << ' ' << z << error << std::endl;
          return EXIT_FAILURE;
        }
      }
    }
  }

  return EXIT_SUCCESS;
}
#!/usr/bin/env python

import itk
import argparse

itk.auto_progress(2)

parser = argparse.ArgumentParser(description="Process Image Chunks.")
parser.add_argument("input_image")
parser.add_argument("output_image")
args = parser.parse_args()

xDiv = 6
yDiv = 4
zDiv = 1
# 1 for 2D input

Dimension = 3

PixelType = itk.UC
ImageType = itk.Image[PixelType, Dimension]

ReaderType = itk.ImageFileReader[ImageType]
reader = ReaderType.New()
reader.SetFileName(args.input_image)
reader.UpdateOutputInformation()
fullRegion = reader.GetOutput().GetLargestPossibleRegion()
fullSize = fullRegion.GetSize()
# when reading image from file, start index is always 0

# create variables of the proper types
start = itk.Index[Dimension]()
end = itk.Index[Dimension]()
size = itk.Size[Dimension]()

MedianType = itk.MedianImageFilter[ImageType, ImageType]
median = MedianType.New()
median.SetInput(reader.GetOutput())
median.SetRadius(2)

WriterType = itk.ImageFileWriter[ImageType]
writer = WriterType.New()
writer.SetFileName(args.output_image)

# Use for loops to split the image into chunks.
# This way of splitting gives us easy control over it.
# For example, we could make the regions always be of equal size
# in order to create samples for a deep learning algorithm.
# ImageRegionSplitterMultidimensional has a similar
# functionality to what is implemented below.
for z in range(zDiv):
    start[2] = int(fullSize[2] * float(z) / zDiv)
    end[2] = int(fullSize[2] * (z + 1.0) / zDiv)
    size[2] = end[2] - start[2]

    for y in range(yDiv):
        start[1] = int(fullSize[1] * float(y) / yDiv)
        end[1] = int(fullSize[1] * (y + 1.0) / yDiv)
        size[1] = end[1] - start[1]

        for x in range(xDiv):
            start[0] = int(fullSize[0] * float(x) / xDiv)
            end[0] = int(fullSize[0] * (x + 1.0) / xDiv)
            size[0] = end[0] - start[0]

            region = itk.ImageRegion[Dimension]()
            region.SetIndex(start)
            region.SetSize(size)

            # now that we have our chunk, request the filter to only process that
            median.GetOutput().SetRequestedRegion(region)
            median.Update()
            result = median.GetOutput()

            # only needed in case of further manual processing
            result.DisconnectPipeline()

            # possible additional non-ITK pipeline processing

            writer.SetInput(result)

            # convert region into IO region
            ioRegion = itk.ImageIORegion(Dimension)
            ioRegion.SetIndex(0, start[0])
            ioRegion.SetIndex(1, start[1])
            ioRegion.SetIndex(2, start[2])
            ioRegion.SetSize(0, region.GetSize()[0])
            ioRegion.SetSize(1, region.GetSize()[1])
            ioRegion.SetSize(2, region.GetSize()[2])
            # tell the writer this is only a chunk of a larger image
            writer.SetIORegion(ioRegion)

            try:
                writer.Update()
            except Exception as e:
                print("Exception for chunk:", x, y, z, "\n", e)
                sys.exit(1)

Classes demonstrated#

template<typename TInputImage, typename TOutputImage>
class MedianImageFilter : public itk::BoxImageFilter<TInputImage, TOutputImage>

Applies a median filter to an image.

Computes an image where a given pixel is the median value of the the pixels in a neighborhood about the corresponding input pixel.

A median filter is one of the family of nonlinear filters. It is used to smooth an image without being biased by outliers or shot noise.

This filter requires that the input pixel type provides an operator<() (LessThan Comparable).

See

Image

See

Neighborhood

See

NeighborhoodOperator

See

NeighborhoodIterator

ITK Sphinx Examples:

See itk::MedianImageFilter for additional documentation.
template<typename TInputImage>
class ImageFileWriter : public itk::ProcessObject

Writes image data to a single file.

ImageFileWriter writes its input data to a single output file. ImageFileWriter interfaces with an ImageIO class to write out the data. If you wish to write data into a series of files (e.g., a slice per file) use ImageSeriesWriter.

A pluggable factory pattern is used that allows different kinds of writers to be registered (even at run time) without having to modify the code in this class. You can either manually instantiate the ImageIO object and associate it with the ImageFileWriter, or let the class figure it out from the extension. Normally just setting the filename with a suitable suffix (“.png”, “.jpg”, etc) and setting the input to the writer is enough to get the writer to work properly.

See

ImageSeriesReader

See

ImageIOBase

ITK Sphinx Examples:

See itk::ImageFileWriter for additional documentation.