[Insight-developers] Image/ImageAdaptor/DataAccessor modification

[Joshua Cates]

Hi all, Luis,

I have run into a problem this week that we have been unable to resolve in
supporting the image adaptors for non-scalar data types.  The problem is a
seemingly unavoidable penalty in working with non-scalar data. The penalty
results because manipulation of non-scalar data cannot be optimized down
to access by reference in many cases.

In trivial cases of scalar Get and Set, the compiler should optimize Get()
and Set() back down to an operator* call on the actual data pointer,
resulting in data access by reference and no loss in performance.  But in
order to modify fields within more complex data types, you have to work
with the copy of the data that Get() returns and then copy back into the
memory location with Set(). (see the examples below)

I wrote a simple test program (shown below) which gave the following
typical results: 

 (1) gcc version 2.95.1 -O2 (Linux)
 Speed of adaptor supporting interator (for reference) 80000
 Speed of interator that does not support adaptors (for reference) 30000 

 a. Modifying scalar image using adaptor api...	110000	 
    compensated = 30000

 b. Modifying scalar image bypassing adaptor api...	50000
    compensated = 20000

 c. Modifying vector image using adaptor api...	1350000	 
    compensated = 1270000

 d. Modifying vector image bypassing adaptor api...	130000
    compensated = 100000


 (2) MIPSpro Compilers: Version 7.3.1.1m -O2 (IRIX)
 Speed of adaptor supporting interator (for reference) 	70000
 Speed of interator that does not support adaptors (for reference) 20000 

 a. Modifying scalar image using adaptor api...	110000	 
 compensated = 40000

 b. Modifying scalar image bypassing adaptor api...	60000
 compensated = 40000

 c. Modifying vector image using adaptor api...	2090000	 
 compensated = 2020000

 d.Modifying vector image bypassing adaptor api...	170000
 compensated = 150000

I have attempted to compensate for the difference in speed of the
different iterators. The correlation between iterator speed and adaptor
support is accidental and does not imply causality.

The compensated values for a and b show no significant difference as we
would expect.  But compensated differences between c and d are significant
in both case 1 and 2.

The problem I've described imposes a serious penalty for certain cases and
merits more testing. Is there a solution within the adaptor framework that
I have not though of? or do we need an API to bypass adaptors in cases
where the trade-off between adaptor benefits and penalties is not
acceptable?

Let's spend some time discussing adaptor issues in the TCON tomorrow.
Particularly I'd like to hear other developers' experiences with
supporting data accessors in their algorithms--other issues they have
encountered or resolved--and hear how folks are using image adaptors.

Josh.


-------------------- itkAdaptorComparisonTest.cxx

#include <iostream>
#include <time.h>
#include "itkImage.h"
#include "itkSimpleImageRegionIterator.h"
#include "itkImageRegionIterator.h"
#include "itkVector.h"

void AdaptorSupportedIteratorSpeed(itk::Image<float, 3> *img)
{
  itk::SimpleImageRegionIterator<itk::Image<float, 3> >
    it (img, img->GetRequestedRegion());

  it.Begin();
  while( ! it.IsAtEnd() )
    {
      ++it;
    }
  
}

void NoAdaptorSupportIteratorSpeed(itk::Image<float, 3> *img)
{
  itk::ImageRegionIterator<float, 3>
    it (img, img->GetRequestedRegion());

  it.Begin();
  while( ! it.IsAtEnd() )
    {
      ++it;
    }
}

void AdaptorSupportedModifyScalars(itk::Image<float, 3> *img)
{
  itk::SimpleImageRegionIterator<itk::Image<float, 3> >
    it (img, img->GetRequestedRegion());

  it.Begin();
  while( ! it.IsAtEnd() )
    {
      // *it += 3.435f;
      it.Set(it.Get() + 3.435f);
      ++it;
    }
}

void NoAdaptorSupportModifyScalars(itk::Image<float, 3> *img)
{
  itk::ImageRegionIterator<float, 3>
    it (img, img->GetRequestedRegion());

  it.Begin();
  while( ! it.IsAtEnd() )
    {
      *it += 3.435f;
      ++it;
    }
}

void AdaptorSupportedModifyVectors(itk::Image<itk::Vector<float, 3>, 3>
*img)
{
  typedef itk::Vector<float, 3> VectorType;
  const unsigned int N = 3;
  unsigned int i;
  VectorType temp_vector;
  
  itk::SimpleImageRegionIterator<itk::Image<VectorType, 3> >
    it (img, img->GetRequestedRegion());
  
  it.Begin();
  while( ! it.IsAtEnd() )
    {
      temp_vector = it.Get();

      for (i = 0; i<N; ++i) temp_vector[i] += 3.435f;

      it.Set(temp_vector);
       ++it;
    }
}

void NoAdaptorSupportModifyVectors(itk::Image<itk::Vector<float, 3>, 3>
*img)
{
  typedef itk::Vector<float, 3> VectorType;
  const unsigned int N = 3;
  unsigned int i;
  
  itk::ImageRegionIterator<VectorType, 3>
    it (img, img->GetRequestedRegion());
  
  it.Begin();
  while( ! it.IsAtEnd() )
    {
      for (i = 0; i<N; ++i)  (*it)[i] += 3.435f;
       ++it;
    }
}

int main()
{
  typedef itk::Image<float, 3> ScalarImageType;
  typedef itk::Image<itk::Vector<float, 3>, 3> VectorImageType;

  clock_t start, stop, no_adaptor_comp, adaptor_comp;
  
  // Set up some images
  itk::ImageRegion<3> region;
  itk::Size<3> size;
   size[0] = 100;
   size[1] = 100;
   size[2] = 100;
  itk::Index<3> index;
   index[0] = 0;
   index[1] = 0;
   index[2] = 0;
  region.SetSize(size);
  region.SetIndex(index);

  ScalarImageType::Pointer scalar_image = ScalarImageType::New();
  VectorImageType::Pointer vector_image = VectorImageType::New();

  scalar_image->SetLargestPossibleRegion(region);
  scalar_image->SetBufferedRegion(region);
  scalar_image->SetRequestedRegion(region);
  scalar_image->Allocate();

  vector_image->SetLargestPossibleRegion(region);
  vector_image->SetBufferedRegion(region);
  vector_image->SetRequestedRegion(region);
  vector_image->Allocate();

  // Time trials
  std::cout << "Speed of adaptor supporting interator (for reference) \t";
  start = clock();
  AdaptorSupportedIteratorSpeed(scalar_image);
  stop = clock();
  adaptor_comp = stop - start;
  std::cout << adaptor_comp << std::endl;

  std::cout << "Speed of interator that does not support adaptors (for
reference) \t";
  start = clock();
  NoAdaptorSupportIteratorSpeed(scalar_image);
  stop = clock();
  no_adaptor_comp = stop - start;
  std::cout << no_adaptor_comp << std::endl;
  
  std::cout << "Modifying scalar image using adaptor api...\t";
  start = clock();
  AdaptorSupportedModifyScalars(scalar_image);
  stop = clock();
  std::cout << (stop - start) << "\t compensated = " << (stop-start) -
    adaptor_comp << std::endl;

  std::cout << "Modifying scalar image bypassing adaptor api...\t";
  start = clock();
  NoAdaptorSupportModifyScalars(scalar_image);
  stop = clock();
  std::cout << (stop - start) << "\t compensated = " << (stop-start) -
    no_adaptor_comp << std::endl;

  std::cout << "Modifying vector image using adaptor api...\t";
  start = clock();
  AdaptorSupportedModifyVectors(vector_image);
  stop = clock();
  std::cout << (stop - start) << "\t compensated = " << (stop-start) -
    adaptor_comp << std::endl;

  std::cout << "Modifying vector image bypassing adaptor api...\t";
  start = clock();
  NoAdaptorSupportModifyVectors(vector_image);
  stop = clock();
  std::cout << (stop - start) << "\t compensated = " << (stop-start) -
    no_adaptor_comp <<std::endl; 

  return 0;
}

______________________________
 Josh Cates			
 School of Computer Science	
 University of Utah
 Email: cates@cs.utah.edu
 Phone: (801) 587-7697
 URL:   www.cs.utk.edu/~cates