Iterate Region in Image With Access to Index Without Write Access#

Synopsis#

Iterate over a region of an image with efficient access to the current index (without write access).

Results#

../../../../_images/Yinyang4.png

Yinyang.png#

../../../../_images/IterateRegionWithAccessToIndexWithoutWriteAccess.png

Yinyang.png In VTK Window#

Output:

An extensive list of the neighborhood will be printed to the screen.

Code#

C++#

#include "itkImage.h"
#include "itkImageFileReader.h"
#include "itkImageRegionConstIteratorWithIndex.h"

#ifdef ENABLE_QUICKVIEW
#  include "QuickView.h"
#endif

int
main(int argc, char * argv[])
{
  if (argc < 2)
  {
    std::cerr << "Required: filename" << std::endl;
    return EXIT_FAILURE;
  }

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

  ImageType::Pointer image = itk::ReadImage<ImageType>(argv[1]);

  ImageType::SizeType regionSize;
  regionSize[0] = 5;
  regionSize[1] = 4;

  ImageType::IndexType regionIndex;
  regionIndex[0] = 0;
  regionIndex[1] = 0;

  ImageType::RegionType region;
  region.SetSize(regionSize);
  region.SetIndex(regionIndex);

  itk::ImageRegionConstIteratorWithIndex<ImageType> imageIterator(image, region);

  while (!imageIterator.IsAtEnd())
  {
    std::cout << "Index: " << imageIterator.GetIndex() << " value: " << (int)imageIterator.Get() << std::endl;

    ++imageIterator;
  }

#ifdef ENABLE_QUICKVIEW
  // Visualize
  QuickView viewer;
  viewer.AddImage<ImageType>(image);
  viewer.Visualize();
#endif
  return EXIT_SUCCESS;
}

Classes demonstrated#

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

A multi-dimensional iterator templated over image type that walks an image region and is specialized to keep track of its index location.

The “WithIndex” family of iteators was designed for algorithms that use both the values and locations of image pixels in calculations. Unlike ImageRegionIterator, which calculates an index only when requested, ImageRegionIteratorWithIndex maintains its index location as a member variable that is updated during increment and decrement operations. Iteration speed is penalized, but index queries become more efficient.

ImageRegionConstIteratorWithIndex is a multi-dimensional iterator, requiring 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, ImageRegionConstIteratorWithIndex holds a reference to the image over which it is traversing.

ImageRegionConstIteratorWithIndex 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.

operator++ 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.Begin();

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.End();

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

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.
MORE INFORMATION

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

See

ImageRegionRange

See

ImageRegionIndexRange

ITK Sphinx Examples:

Subclassed by itk::FrequencyFFTLayoutImageRegionConstIteratorWithIndex< TImage >, itk::FrequencyHalfHermitianFFTLayoutImageRegionConstIteratorWithIndex< TImage >, itk::FrequencyImageRegionConstIteratorWithIndex< TImage >, itk::FrequencyShiftedFFTLayoutImageRegionConstIteratorWithIndex< TImage >, itk::ImageRegionExclusionConstIteratorWithIndex< TImage >, itk::ImageRegionIteratorWithIndex< TImage >

See itk::ImageRegionConstIteratorWithIndex for additional documentation.