Dear all,<br>
<br>
I am writing a filter to calculate the gradient structure tensor. This
involves calculating an NxN matrix for each voxel of the image, and
then applying a gaussian filter to each one of the components of the
matrix (i.e., taking one component of the matrix at every voxel to
create a new image, and applying a gaussian filter on this new image).
Based on HessianRecursiveGausianImageFilter, I am using
NthElementImageAdaptor and three instances of
RecursiveGaussianImageFilter to filter in the x, y and z directions. As
I saw in previous submissions to the list, I defined the first of the
gaussian filters to accept the adaptor as an input.<br>
I have been getting errors in the results (basically, the filters
behave as if the pixel spacing were always [1 1 1], though it's not),
so I wrote this very simple piece of code to check that image spacing
has the expected values. I have attached the code below. It works by
reading a tensor image from a file, and applying a gaussian filter to
the first component of the tensor (of course, I would have to define
things as filter direction, sigma, etc., but here I just wanted to
check the spacing). The input image is all zeros, though I don't think
this matters in this case. With an input image with spacing [2 2 2],
the program output is:<br>
<br>
m_ImageAdaptor->GetSpacing(): [2, 2, 2]<br>
m_ComponentFilter->GetOutput()->GetSpacing(): [1, 1, 1]<br>
<br>
which is exactly the same error I got in my structure tensor
calculation program. It seems that the output of the filter has a
different spacing from the input, so when I attach a second filter to
the output of the first one, the result is wrong. Could anybody tell me
what I am doing wrong here? Thanks a lot,<br>
<br>
Vicente<br>
<br>
/**************************** Code follows ***********************/<br>
<br>
#include "itkSymmetricSecondRankTensor.h"<br>
#include "itkImage.h"<br>
#include "itkImageFileReader.h"<br>
#include "itkRecursiveGaussianImageFilter.h"<br>
#include "itkNthElementImageAdaptor.h"<br>
<br>
<br>
<br>
int main( int argc, char *argv[] )<br>
{<br>
const unsigned
int
Dimension = 3;<br>
typedef
double
InputPixelType;<br>
<br>
typedef itk::Image<
itk::SymmetricSecondRankTensor< InputPixelType, Dimension > ,
Dimension > InputImageType;<br>
typedef itk::Image< InputPixelType, Dimension > OutputImageType;<br>
<br>
typedef itk::ImageFileReader< InputImageType > ReaderType;<br>
ReaderType::Pointer reader = ReaderType::New(); <br>
reader->SetFileName( argv[1] );<br>
<br>
typedef itk::NthElementImageAdaptor< InputImageType, InputPixelType > OutputImageAdaptorType;<br>
OutputImageAdaptorType::Pointer m_ImageAdaptor = OutputImageAdaptorType::New(); <br>
m_ImageAdaptor->SetImage( reader->GetOutput() );<br>
m_ImageAdaptor->Allocate();<br>
m_ImageAdaptor->SelectNthElement( 0 );<br>
<br>
typedef itk::RecursiveGaussianImageFilter<
OutputImageAdaptorType, OutputImageType >
ComponentFilterType;<br>
ComponentFilterType::Pointer m_ComponentFilter = ComponentFilterType::New();<br>
m_ComponentFilter->SetInput( m_ImageAdaptor );<br>
<br>
m_ComponentFilter->Update();<br>
<br>
std::cout << "m_ImageAdaptor->GetSpacing():
" << m_ImageAdaptor->GetSpacing() << std::endl;<br>
std::cout << "m_ComponentFilter->GetOutput()->GetSpacing(): " << <br>
m_ComponentFilter->GetOutput()->GetSpacing()
<< std::endl;<br>
<br>
return EXIT_SUCCESS;<br>
}<br>
<br>
<br>
<br>
<br>