[Insight-users] Copy filter output

[Silvano Agliozzo]

Hi Luis, 

thank you for your answer. However I still have a problem.

I've used the duplicator as you suggested, because I have to process 
large 3D images (CT scans). Nevertheless I still get the same results 
for all the different images even though I set the 
RecursiveGaussianImageFilters to calculate different derivatives.
Moreover the values of the filter outputs do not correspond to the 
expected one. 
If you don't mind, I wite hereafter some line of the code:

  ....

  double spacing[ ImageType::ImageDimension ];
  spacing[0] = 1.0; // along X direction
  spacing[1] = 1.0; // along Y direction
  spacing[2] = 1.0; // along Z direction

  double origin[ ImageType::ImageDimension ];
  origin[0] = 0.0; // X coordinate
  origin[1] = 0.0; // Y coordinate
  origin[2] = 0.0; // Z coordinate

  // vol = data to be processedA
  vol->SetRegions(region);
  vol->Allocate();
  vol->FillBuffer( 0.0 );
  vol->SetSpacing( spacing );
  vol->SetOrigin( origin );

  pdxx->SetRegions(region);// second order derivative
  ...
  ...

  the same for pdyy, pdzz, pdxy, pdxz, pdyz  

  typedef itk::RecursiveGaussianImageFilter< ImageType, ImageType > 
FilterType;

  FilterType::Pointer filterX = FilterType::New();
  FilterType::Pointer filterY = FilterType::New();
  FilterType::Pointer filterZ = FilterType::New();

  filterX->SetNormalizeAcrossScale(false);
  filterY->SetNormalizeAcrossScale(false);
  filterZ->SetNormalizeAcrossScale(false);

  filterX->SetDirection(0);
  filterY->SetDirection(1);
  filterZ->SetDirection(2);

  filterX->SetInput(vol);
  filterY->SetInput(filterX->GetOutput());
  filterZ->SetInput(filterY->GetOutput());

  filterX->SetSigma(1.0);
  filterY->SetSigma(1.0);
  filterZ->SetSigma(1.0);

  typedef itk::ImageDuplicator< ImageType > ImageDuplicatorType;
  ImageDuplicatorType::Pointer duplicator = ImageDuplicatorType::New();

  std::cout<<"Getting pdxx ...";
  filterX->SetOrder(FilterType::SecondOrder);
  filterY->SetOrder(FilterType::ZeroOrder);
  filterZ->SetOrder(FilterType::ZeroOrder);
   
  std::cout<<"filter updating ...";
  filterZ->Update();

  duplicator->SetInputImage( filterZ->GetOutput() );
  std::cout<<"duplicator updating ...";
  duplicator->Update();
  
  pdxx = duplicator->GetOutput();
 
  pdyy ...
 
  std::cout<<"Getting pdzz ...";
  filterX->SetOrder(FilterType::ZeroOrder);
  filterY->SetOrder(FilterType::SecondOrder);
  filterZ->SetOrder(FilterType::ZeroOrder);
  
  std::cout<<"filter updating ...";
  filterZ->Update();

  std::cout<<"duplicator updating ...";
  duplicator->Update();
  
  pdzz = duplicator->GetOutput();

  std::cout<<"Getting pdxy ...";
  filterX->SetOrder(FilterType::FirstOrder);
  filterY->SetOrder(FilterType::FirstOrder);
  filterZ->SetOrder(FilterType::ZeroOrder);
   
  std::cout<<"filter updating ...";
  filterZ->Update();

  std::cout<<"duplicator updating ...";
  duplicator->Update();
  
  pdxy = duplicator->GetOutput();

  ....and so on for the rest ...(pdxz, pdyz)


  but pdxx, pdyy, ..., pdyz are the same.
  I also tried to insert the method Modified() before updating the filter 
without get any changes. 


thank you in advance,

a+

Silvano
On Wed, 14 Jan 2004, Luis Ibanez wrote:

> 
> Hi Silvano,
> 
> Welcome to ITK !
> 
>  From your message it looks like you are trying to
> create one RecursiveGaussianImageFilter and
> reuse it in order to compute multiple second
> derivatives.
> 
> Reusing filters is a bad practice in a data
> pipepline system as ITK or VTK.
> 
> Don't be affraid to be generous in creating ITK
> filters (unless your images are really large).
> 
> 
> 
> I assume that you are trying to compute the full
> set of second derivatives in a 3D image, which
> is equivalent to compute the Hessian Matrix for
> every pixel
> 
> 
>             Ixx  Ixy  Ixz
>             Iyx  Iyy  Iyz
>             Izx  Izy  Izz
> 
> 
> given that this is a symmetric matrix, you only
> need to compute the 6 images:
> 
>     {  Ixx   Iyy  Izz   Ixy  Ixz  Iyz  }
> 
> 
> Let's take Izz first
> 
> 
> You need 3 RecursiveGaussianImageFilters in order
> to compute the Izz image.
> 
> FilterType::Pointer ga = FilterType::New();
> FilterType::Pointer gb = FilterType::New();
> FilterType::Pointer gc = FilterType::New();
> 
> ga->SetDirection( 0 );
> gb->SetDirection( 1 );
> gc->SetDirection( 2 );
> 
> ga->SetZeroOrder();
> gb->SetZeroOrder();
> gc->SetSecondOrder();
> 
> ga->SetInput( inputImage );
> gb->SetInput( ga->GetOutput() );
> gc->SetInput( gb->GetOutput() );
> 
> gc->Update();
> 
> At this point you have Izz as the output of the
> filter "gc".
> 
> You can copy this output out of the pipeline
> by using the ImageDuplicator class
> http://www.itk.org/Insight/Doxygen/html/classitk_1_1ImageDuplicator.html
> 
> DuplicatorType::Pointer duplicator = DuplicatorType::New();
> 
> duplicator->SetInput(  gc->GetOutput() );
> duplicator->Update();
> 
> SecondDerivativeImageType::Pointer Izz = duplicator->GetOutput();
> 
> 
> Note the "duplicator" is NOT an ITK filter.
> It doesn't to make part of the pipeline,
> despite the fact that its API look like
> a filter.
> 
> 
> 
> 
> Now, you can readjust the pipeline in order
> to compute Iyy, just do
> 
> gc->SetDirection( 1 );  // gc now works along Y
> gb->SetDirection( 2 );  // gb now works along Z
> 
> gc->Update();
> duplicator->Update();
> 
> SecondDerivativeImageType::Pointer Iyy = duplicator->GetOutput();
> 
> 
> 
> 
> In order to get Ixx you do
> 
> gc->SetDirection( 0 );  // gc now works along X
> ga->SetDirection( 1 );  // ga now works along Y
> 
> gc->Update();
> duplicator->Update();
> 
> SecondDerivativeImageType::Pointer Ixx = duplicator->GetOutput();
> 
> 
> 
> 
> 
> 
> For the cross derivatives you do
> 
> ga->SetDirection( 0 );
> gb->SetDirection( 1 );
> gc->SetDirection( 2 );
> 
> ga->SetZeroOrder();
> gb->SetSecondOrder();
> gc->SetSecondOrder();
> 
> gc->Update();
> duplicator->Update();
> 
> SecondDerivativeImageType::Pointer Iyz = duplicator->GetOutput();
> 
> 
> ga->SetDirection( 1 );
> gb->SetDirection( 0 );
> gc->SetDirection( 2 );
> 
> ga->SetZeroOrder();
> gb->SetSecondOrder();
> gc->SetSecondOrder();
> 
> gc->Update();
> duplicator->Update();
> 
> SecondDerivativeImageType::Pointer Ixz = duplicator->GetOutput();
> 
> 
> ga->SetDirection( 2 );
> gb->SetDirection( 0 );
> gc->SetDirection( 1 );
> 
> ga->SetZeroOrder();
> gb->SetSecondOrder();
> gc->SetSecondOrder();
> 
> gc->Update();
> duplicator->Update();
> 
> SecondDerivativeImageType::Pointer Ixy = duplicator->GetOutput();
> 
> 
> 
> 
> 
> Regards,
> 
> 
>     Luis
> 
> 
> ------------------------
> Silvano Agliozzo wrote:
> 
> > Dear all,
> > 
> > I'am a newcomer to the Itk libraries, so I could ask a trivial and I 
> > apologize in advance.
> > 
> > I need to use the RecursiveGuassianImageFilter in a 3D image in order to 
> > obtain second order derivatives, but I need to store the output of the filter 
> > to different 3D images. I do not want to get the pointer of the filter 
> > output since each time the pipeline of the filter is updated, the voxel 
> > values of all the images previously calculated, are updated too. 
> > 
> > I would like to know how  can I perform this task ? 
> > 
> > thank you very much,
> > 
> > Silvano  
> > 
> > _______________________________________________
> > Insight-users mailing list
> > Insight-users at itk.org
> > http://www.itk.org/mailman/listinfo/insight-users
> > 
> 
> 
>