Iterate Over Image While Skipping Specific Region#

Synopsis#

Iterator over an image skipping a specified region.

Results#

Output:

0
0
0
0
0
0
Visited 6

Code#

C++#

#include "itkImage.h"
#include "itkImageRegionExclusionConstIteratorWithIndex.h"

namespace
{
using ImageType = itk::Image<int, 2>;
}

static void CreateImage(ImageType::Pointer);

int
main()
{

  ImageType::SizeType exclusionRegionSize;
  exclusionRegionSize.Fill(2);

  ImageType::IndexType exclusionRegionIndex;
  exclusionRegionIndex.Fill(0);

  ImageType::RegionType exclusionRegion(exclusionRegionIndex, exclusionRegionSize);

  auto image = ImageType::New();
  CreateImage(image);

  itk::ImageRegionExclusionConstIteratorWithIndex<ImageType> imageIterator(image, image->GetLargestPossibleRegion());
  imageIterator.SetExclusionRegion(exclusionRegion);

  unsigned int numberVisited = 0;
  while (!imageIterator.IsAtEnd())
  {
    std::cout << imageIterator.Get() << std::endl;
    ++imageIterator;
    ++numberVisited;
  }

  std::cout << "Visited " << numberVisited << std::endl;

  return EXIT_SUCCESS;
}

void
CreateImage(ImageType::Pointer image)
{
  ImageType::SizeType regionSize;
  regionSize.Fill(3);

  ImageType::IndexType regionIndex;
  regionIndex.Fill(0);

  ImageType::RegionType region(regionIndex, regionSize);

  image->SetRegions(region);
  image->Allocate();
  image->FillBuffer(0);
}

Classes demonstrated#

template<typename TImage>
class ImageRegionExclusionConstIteratorWithIndex : public itk::ImageRegionConstIteratorWithIndex<TImage>

A multi-dimensional image iterator that walks an image region, excluding a second region that may partially or completely overlap the first, with read-only access to pixels.

ImageRegionExclusionConstIteratorWithIndex is a class templated over the image type that represents a multi-dimensional iterator. The exclusion region is set after construction. By default the exclusion region is empty and the iterator will behave as the itk::ImageRegionConstIteratorWithIndex with a penalty in performance due to internal bounds checking. There are no restrictions on valid exclusion regions.

As with other ITK image iterators, ImageRegionExclusionConstIteratorWithIndex requires that more information be specified before the iterator can be used than conventional iterators. Whereas the std::vector::iterator from the STL only needs to be passed a pointer to establish the iterator, the multi-dimensional image iterator needs a pointer, the size of the buffer, the size of the region, the start index of the buffer, and the start index of the region. To gain access to this information, ImageRegionExclusionConstIteratorWithIndex holds a reference to the image over which it is traversing.

ImageRegionExclusionConstIteratorWithIndex assumes a particular layout of the image data. The is arranged in a 1D array as if it were [][][][slice][row][col] with Index[0] = col, Index[1] = row, Index[2] = slice, etc.

The operator++ method provides a simple syntax for walking around a region of a multidimensional image. operator++ iterates across a row, constraining the movement to within a region of image. When the iterator reaches the boundary of the region along a row, the iterator automatically wraps to the next row, starting at the first pixel in the row that is part of the region. This allows for simple processing loops of the form:

IteratorType it( image, image->GetRequestedRegion() );

it.SetExclusionRegion( exclusionRegion );
it.GoToBegin();

while( ! it.IsAtEnd() )
{
  it.Set( 100.0 + it.Get() );
  ++it;
}

It also can be used for walking in the reverse direction like

IteratorType it( image, image->GetRequestedRegion() );

it.SetExclusionRegion( exclusionRegion );
it.GoToEnd();

while( !it.IsAtBegin() )
{
  it.Set( 100.0 );
  --it;
}

MORE INFORMATION

For a complete description of the ITK Image Iterators and their API, please see the Iterators chapter in the ITK Software Guide. The ITK Software Guide is available in print and as a free .pdf download from https://www.itk.org.

See

ImageConstIterator

See

ConditionalConstIterator

See

ConstNeighborhoodIterator

See

ConstShapedNeighborhoodIterator

See

ConstSliceIterator

See

CorrespondenceDataStructureIterator

See

FloodFilledFunctionConditionalConstIterator

See

FloodFilledImageFunctionConditionalConstIterator

See

FloodFilledImageFunctionConditionalIterator

See

FloodFilledSpatialFunctionConditionalConstIterator

See

FloodFilledSpatialFunctionConditionalIterator

See

ImageConstIterator

See

ImageConstIteratorWithIndex

See

ImageIterator

See

ImageIteratorWithIndex

See

ImageLinearConstIteratorWithIndex

See

ImageLinearIteratorWithIndex

See

ImageRandomConstIteratorWithIndex

See

ImageRandomIteratorWithIndex

See

ImageRegionConstIterator

See

ImageRegionConstIteratorWithIndex

See

ImageRegionExclusionConstIteratorWithIndex

See

ImageRegionExclusionIteratorWithIndex

See

ImageRegionIterator

See

ImageRegionIteratorWithIndex

See

ImageRegionReverseConstIterator

See

ImageRegionReverseIterator

See

ImageReverseConstIterator

See

ImageReverseIterator

See

ImageSliceConstIteratorWithIndex

See

ImageSliceIteratorWithIndex

See

NeighborhoodIterator

See

PathConstIterator

See

PathIterator

See

ShapedNeighborhoodIterator

See

SliceIterator

See

ImageConstIteratorWithIndex

ITK Sphinx Examples:

Subclassed by itk::ImageRegionExclusionIteratorWithIndex< TImage >

See itk::ImageRegionExclusionConstIteratorWithIndex for additional documentation.