Stream A Pipeline

Synopsis

Stream a pipeline. The PipelineMonitorImageFilter shows how the image is processed in three sub-regions. Note that the StreamingImageFilter comes at the end of the pipeline and that no calls to Update() occur on any intermediate filters.

Results

The output LargestPossibleRegion is: ImageRegion (0x8dc420)
  Dimension: 2
  Index: [0, 0]
  Size: [3, 3]


Updated RequestedRegions's:
  ImageRegion (0x8dce80)
  Dimension: 2
  Index: [0, 0]
  Size: [3, 1]

  ImageRegion (0x8dcea8)
  Dimension: 2
  Index: [0, 1]
  Size: [3, 1]

  ImageRegion (0x8dced0)
  Dimension: 2
  Index: [0, 2]
  Size: [3, 1]

Code

C++

#include "itkPipelineMonitorImageFilter.h"
#include "itkRandomImageSource.h"
#include "itkStreamingImageFilter.h"

int
main(int argc, char * argv[])
{
  if (argc != 2)
  {
    std::cerr << "Usage: " << argv[0] << " <NumberOfSplits>" << std::endl;
    return EXIT_FAILURE;
  }
  int numberOfSplits = std::stoi(argv[1]);

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

  using SourceType = itk::RandomImageSource<ImageType>;
  SourceType::Pointer source = SourceType::New();
  ImageType::SizeType size;
  size.Fill(numberOfSplits);
  source->SetSize(size);

  using MonitorFilterType = itk::PipelineMonitorImageFilter<ImageType>;
  MonitorFilterType::Pointer monitorFilter = MonitorFilterType::New();
  monitorFilter->SetInput(source->GetOutput());
  // If ITK was built with the Debug CMake configuration, the filter
  // automatically outputs status information to the console
  monitorFilter->DebugOn();

  using StreamingFilterType = itk::StreamingImageFilter<ImageType, ImageType>;
  StreamingFilterType::Pointer streamingFilter = StreamingFilterType::New();
  streamingFilter->SetInput(monitorFilter->GetOutput());
  streamingFilter->SetNumberOfStreamDivisions(numberOfSplits);

  try
  {
    streamingFilter->Update();
  }
  catch (itk::ExceptionObject & error)
  {
    std::cerr << "Error: " << error << std::endl;
    return EXIT_FAILURE;
  }

  std::cout << "The output LargestPossibleRegion is: " << streamingFilter->GetOutput()->GetLargestPossibleRegion()
            << std::endl;
  std::cout << std::endl;

  const MonitorFilterType::RegionVectorType updatedRequestedRegions = monitorFilter->GetUpdatedRequestedRegions();
  std::cout << "Updated RequestedRegions's:" << std::endl;
  for (const auto & updatedRequestedRegion : updatedRequestedRegions)
  {
    std::cout << "  " << updatedRequestedRegion << std::endl;
  }

  return EXIT_SUCCESS;
}

Python

#!/usr/bin/env python
import itk

from distutils.version import StrictVersion as VS
if VS(itk.Version.GetITKVersion()) < VS("4.10.0"):
    print("ITK 4.10.0 is required.")
    sys.exit(1)

if len(sys.argv) != 2:
    print('Usage: ' + sys.argv[0] + ' <NumberOfSplits>')
    sys.exit(1)
numberOfSplits = int(sys.argv[1])

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

source = itk.RandomImageSource[ImageType].New()
size = itk.Size[Dimension]()
size.Fill(numberOfSplits)
source.SetSize(size)

monitorFilter = itk.PipelineMonitorImageFilter[ImageType].New()
monitorFilter.SetInput(source.GetOutput())

streamingFilter = itk.StreamingImageFilter[ImageType, ImageType].New()
streamingFilter.SetInput(monitorFilter.GetOutput())
streamingFilter.SetNumberOfStreamDivisions(numberOfSplits)

streamingFilter.Update()

print('The output LargestPossibleRegion is: ' +
      str(streamingFilter.GetOutput().GetLargestPossibleRegion()))

print('')

updatedRequestedRegions = monitorFilter.GetUpdatedRequestedRegions()
print("Updated RequestedRegion's:")
for ii in range(len(updatedRequestedRegions)):
    print('  ' + str(updatedRequestedRegions[ii]))

Classes demonstrated