[Insight-users] Copy filter output
Luis Ibanez
luis.ibanez at kitware.com
Mon, 19 Jan 2004 13:54:08 -0500
Hi Silvano,
1) The warning message that you get was
put there on purpose with the intention
of keeping in mind that we still have
something to fix in this class.
You can safely disregard the message.
2) The method SetZeroOrder(), SetFirstOrder()
and SetSecondOrder() are entirely equivalent
to SetOrder( ) with arguments Zero, First
and Second. These methods were added in order
to make possible the selection of order from
scripting languages such as Python and Tcl.
The reason is that the enum doesn't have a
counterpart on the scripted languages.
In C++ you can freely use the method of your
preference.
3) About the cross derivatives... you are right !
This is my mistake. I just fixed this in CVS
according to your suggestion.
http://www.itk.org/cgi-bin/cvsweb.cgi/Insight/Examples/Filtering/SecondDerivativeRecursiveGaussianImageFilter.cxx.diff?r1=1.1&r2=1.2&cvsroot=Insight&sortby=date
Thanks for pointing this out.
Please let us know if you find any other problems,
Thanks
Luis
-----------------------------
Silvano Agliozzo wrote:
> Hi Luis,
>
> I update the new version of the ImageDuplicator class.
>
> It works even though it produces the following error each time I
> calculated a derivative.B
>
> WARNING: In /usr/local/itk/Insight/Code/Common/itkImageDuplicator.txx,
> line 52
> Self (0x8099fe0): Input image has not been modified
>
> I suppose this is due to the modification you did.
>
> In the documentation of the RecursiveGaussianImageFilter, the derivative
> order of this filter is set by the method SetOrder(OrderEnumType _arg), but in
> the example
>
> Insight/Examples/Filtering/SecondDerivativeRecursiveGaussianImageFilter.cxx
>
> the order is set differently by the SetZeroOrder(), SetFirstOrder(),
> SetSecondOrder() methods. I checked both way ot setting and I did not
> check any changes, but I would like to know the difference,
>
> Moreover I would like to point out that there could beB a mistake in the
> same example.
>
> For calculating the mixed second order derivative, e.g. Ixy,
> for a 3D image, one should set the filters filterX, filterY, filterZ
> in the following way:
>
> filterX->SetFirstOrder();
> filterY->SetFirstOrder();
> filterZ->SetZeroOrder();
>
> Conversely, in the example this derivative is calculated by setting:
>
> filterX->SetSecondOrder();
> filterY->SetSecondOrder();
> filterZ->SetZeroOrder();
>
> Theoretically, if you set this derivative orders, you get a mixed forth
> order derivative. I did not get insight the implementation of the code,
> and may-be you have an explanation for that.
>
> Concerning the calculated values, they are close to the expected
> ones (I calculated a spherical and cilinder intensity fields), a part the
> second order derivative about z, i.e. Izz, which is always smaller than
> Ixx and Iyy.
>
> a+
>
> Silvano
>
>
> On Thu, 15 Jan 2004, Luis Ibanez wrote:
>
>
>>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
>>>
>>
>>
>>
>>_______________________________________________
>>Insight-users mailing list
>>Insight-users at itk.org
>>http://www.itk.org/mailman/listinfo/insight-users
>>
>
>
>