[Insight-users] Using image iterators combined with filters in a loop

[Merrem, Andreas]

-----Ursprüngliche Nachricht-----
Von: Merrem, Andreas 
Gesendet: Mittwoch, 26. Mai 2010 20:08
An: 'Luis Ibanez'
Cc: insight-users at itk.org
Betreff: Using image iterators combined with filters in a loop

 
Hello Luis (or anyone else who can help),

Thanks a lot for your last answer to my registration question. 
Right now I'm stuck with another problem. I am writing a registration program that involves using image iterators and filters within a loop. To test some parts of what I want to do, I have written a test program that should read two images and process them with a loop. I want to multiply all pixel values by a factor and then smooth one of them during each loop iteration.

This way, the program compiles, but it doesn't work the way it should (the first output image is the input image multiplied by the factor only once and then smoothed).
It does work the way it should when I don't use a loop, copy and paste the content of the loop e.g. 3 times, and instantiate a both a new iterator and a new smoothing filter for each "iteration".
However, I can't do that for an arbitrary number of iterations (I would have multiple definitions). 

Is there a way to use one smoothing filter and one iterator for my program? 
My short program is below this message, please ignore the shuffling around with the vector image, it was just a test, and it has nothing to do with the problem. 

Thank you very much for your help,

Andreas
 

//Test Program: smoothing vector images

#include "itkImage.h"
#include "itkImageFileReader.h"
#include "itkImageFileWriter.h"
#include "itkSmoothingRecursiveGaussianImageFilter.h"
#include "itkCovariantVector.h"
#include <iostream>


int main( int argc, char * argv[] )
{


//Read Images

const    unsigned int    Dimension = 2;
typedef  float   PixelType;

typedef itk::Image <PixelType, Dimension> FixedImageType;
typedef itk::ImageFileReader <FixedImageType> FixedImageReaderType;

FixedImageReaderType::Pointer fixedImageReader1  = FixedImageReaderType::New();
FixedImageReaderType::Pointer fixedImageReader2  = FixedImageReaderType::New();


fixedImageReader1-> SetFileName (argv[1]);
fixedImageReader2-> SetFileName (argv[2]);

  try
    {
    fixedImageReader1-> Update();
	fixedImageReader2-> Update();

    }
  catch (itk::ExceptionObject & e)
    {
    std::cerr << "exception in file reader " << std::endl;
    std::cerr << e << std::endl;
    return EXIT_FAILURE;
    }
  FixedImageType::Pointer FixedImage1 = fixedImageReader1-> GetOutput();
  FixedImageType::Pointer FixedImage2 = fixedImageReader2-> GetOutput();

//Create Vector Image

  typedef float VectorComponentType;
  typedef itk::CovariantVector <VectorComponentType, 2> VectorType;
  typedef itk::Image <VectorType, Dimension> VectorImageType;
  VectorImageType::Pointer VectorImage = VectorImageType::New();
  VectorImage-> SetRegions (FixedImage1-> GetLargestPossibleRegion());
  VectorImage-> Allocate();



typedef itk::ImageRegionIterator <VectorImageType> VecIteratorType;
VecIteratorType vecIt (VectorImage, VectorImage-> GetLargestPossibleRegion());

typedef itk::ImageRegionIterator <FixedImageType> FixIteratorType;
FixIteratorType fixIt1 (FixedImage1, FixedImage1-> GetLargestPossibleRegion());
FixIteratorType fixIt2 (FixedImage2, FixedImage2-> GetLargestPossibleRegion());



VectorType pix;

float L=1.3;


  //Smoothing Filter (which doesn't work)

  typedef itk::SmoothingRecursiveGaussianImageFilter <FixedImageType, FixedImageType> SmoothFilterType;
 SmoothFilterType::Pointer Smoother1 = SmoothFilterType::New();
 Smoother1-> SetNormalizeAcrossScale (false);
 Smoother1-> SetSigma (1.0);


  //Output Image from Smoother, and iterator to access it (the one that doesn't work)

 FixedImageType::Pointer SmoothIm = FixedImageType::New();
 SmoothIm=FixedImage1; 
 FixIteratorType smoothIt1 (SmoothIm, SmoothIm->GetLargestPossibleRegion());
// the largest possible region for the smoother output image is the same as for the fixed (input) image

 /////////////////////////////Loop for processing/////////////////////////////

int n=3;

for (int i=0; i<n; i++){

  //Fill Vector Image  
  for ( fixIt1.GoToBegin(), fixIt2.GoToBegin(), vecIt.GoToBegin(); !vecIt.IsAtEnd(); ++fixIt1, ++fixIt2, ++vecIt)
    {
		pix[0] = fixIt1.Get();
		pix[1] = fixIt2.Get();
		vecIt.Set (pix);
    }
 //Vector Image*L 
  for (vecIt.GoToBegin(); !vecIt.IsAtEnd(); ++vecIt)
    {
vecIt.Set (L*vecIt.Get());
    }
// Back to component images

  for ( fixIt1.GoToBegin(), fixIt2.GoToBegin(), vecIt.GoToBegin(); !fixIt1.IsAtEnd(); ++fixIt1, ++fixIt2, ++vecIt)
    {
		fixIt1.Set(vecIt.Get()[0]);
		fixIt2.Set(vecIt.Get()[1]);
    }


  //Copy and smooth first component (this is where the problem is)



 Smoother1-> SetInput (FixedImage1);
 Smoother1-> Update();
 SmoothIm = Smoother1-> GetOutput();

  for ( fixIt1.GoToBegin(), smoothIt1.GoToBegin(); !fixIt1.IsAtEnd(); ++fixIt1, ++smoothIt1)
  {
	  fixIt1.Set(smoothIt1.Get());
  }


}
//////////////////////////////////////////////////////////////////////

   //Write Images to File

  
  
  typedef unsigned char OutputPixelType;
  typedef itk::Image <OutputPixelType,Dimension> OutputImageType;

  typedef itk::ImageFileWriter <OutputImageType> WriterType;
  WriterType::Pointer Writer1 = WriterType::New();
  WriterType::Pointer Writer2 = WriterType::New();
  
  typedef itk::CastImageFilter <FixedImageType,OutputImageType> CastFilterType;
  CastFilterType::Pointer  caster1 =  CastFilterType::New();
  CastFilterType::Pointer  caster2 =  CastFilterType::New();

  caster1-> SetInput (FixedImage1);
  caster2-> SetInput (FixedImage2);
  Writer1-> SetInput (caster1-> GetOutput());
  Writer2-> SetInput (caster2-> GetOutput());
  Writer1-> SetFileName (argv[3]);
  Writer2-> SetFileName (argv[4]);

  Writer1-> Update();
  Writer2-> Update();

  return EXIT_SUCCESS;
  }