[Insight-users] Copy filter output
Luis Ibanez
luis.ibanez at kitware.com
Thu, 15 Jan 2004 13:26:20 -0500
Hi Silvano,
It seems that you found an interesting Bug.
For some reason the pipeline is not updating
the modified time of the output image in the
RecursiveGaussianImageFilter.
This has been entered as Bug#517.
http://www.itk.org/Bug/bug.php?op=show&bugid=517&pos=11
In the meantime, the following actions have
been taken to make your life easier:
1) The check for the modified time has been
temporarily removed from the ImageDuplicator
class. So it will always duplicate the image
regardless of the time stamp. This file is in
Insight/Code/Common/itkImageDuplicator.txx
2) The full code example of the Second derivative
computation has been commited under
Insight/Examples/Filtering
SecondDerivativeRecursiveGaussianImageFilter.cxx
You may simply update these two files in your CVS
checkout and give them a try.
Please let us know if you find any further problems,
Thanks
Luis
--------------------------
Silvano Agliozzo wrote:
> 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
>>>
>>
>>
>>
>
>
>
> _______________________________________________
> Insight-users mailing list
> Insight-users at itk.org
> http://www.itk.org/mailman/listinfo/insight-users
>