<br>Hi Amardeep,<br><br>Thanks for posting the source code and the image that you are using.<br><br>We reformatted your code and were able to reproduce what you are reporting:<br><br>That the directions of the gradient (of the gradient magnitude) <b>seem</b> to be <br>
pointing in the wrong direction.<br><br>Please find attached the two versions that we reformat of your code.<br><br>After a lot of head scratching we found that the apparent problem, is that<br>your image has DIRECTION cosines that are not those of a unit matrix.<br>
<br>The direction cosines of your image actually are:<br><br> -1 0 0 <br> 0 -1 0 <br> 0 0 1<br><br>Which means that the X and Y axis are reflected with respect to<br>what you would usually expect.<br>
<br>Note that ITK gradient filters *DO* take direction into account <br>when they are computing the gradient vectors. since ITK 3.10,<br>we have set the itkImage to consider Direction by default.<br><br>Therefore, at the end, this is simply a case of using inadecuate <br>
Visualization tools. The viewer that you are using doesn't take<br>the image direction into account.<br><br>For your convenience, we have added the explicit calls for the<br>gradient to take or not to take direction into account. <br>
Please note however, that ignoring the direction will result<br>in inconsistent computations.<br><br><br>Also, you will better understand what is going on, if you <br>use an asymmetric image. The square that you are using<br>
doen't help to disambiguate between the +Y and -Y directions<br>of the axis.<br><br><br><br> Regards,<br><br><br> Luis<br><br>-----------------------------------------------------<br><div class="gmail_quote">On Thu, Feb 12, 2009 at 7:58 AM, Amardeep Singh <span dir="ltr"><<a href="mailto:amar.singh@gmx.de">amar.singh@gmx.de</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="border-left: 1px solid rgb(204, 204, 204); margin: 0pt 0pt 0pt 0.8ex; padding-left: 1ex;">Dear ITK Users<br>
<br>
I am still having problems with the gradient vector flow and so I went a step back and looked at the gradients<br>
of my image. At the moment, I am just dealing with a synthetic image of a square (see attachment).<br>
I calculate the gradient magnitude with the itkGradientMagnitudeImageFilter and then the gradient with the<br>
itkGradientImageFilter. I save the output of the gradient filter as a *.vtk image file and visualize it in Paraview.<br>
I find that the vectors are pointing outwards, whereas I would expect them to point inwards, to the middle of<br>
the edge.<br>
Could anyone tell me whats going wrong, please? Or is everything correct and I misunderstood something?<br>
Thank you for any help.<br>
<br>
Best regards<br>
Amar<br>
<br>
P.S. I have attached my code below:<br>
<br>
/*<br>
* gvfOnSyntheticImage.cc<br>
*<br>
* Created on: 11-Feb-2009<br>
* */<br>
<br>
<br>
<br>
//ITK Headers<br>
#include "itkImage.h"<br>
<br>
#include "itkImageFileReader.h"<br>
#include "itkImageFileWriter.h"<br>
#include "itkCovariantVector.h"<br>
#include "itkGradientVectorFlowImageFilter.h"<br>
//#include "itkLaplacianImageFilter.h"<br>
#include "itkGradientMagnitudeRecursiveGaussianImageFilter.h"<br>
#include "itkGradientMagnitudeImageFilter.h"<br>
#include "itkGradientImageFilter.h"<br>
#include "itkCastImageFilter.h"<br>
<br>
#include "vtkImageWriter.h"<br>
<br>
#include <iostream><br>
<br>
<br>
int main( int argc, char *argv[] )<br>
{<br>
<br>
typedef float InternalPixelType;<br>
typedef float OutputPixelType;<br>
typedef unsigned char BinaryPixelType;<br>
<br>
const unsigned int Dimension = 3;<br>
<br>
typedef itk::Image< InternalPixelType, Dimension > InternalImageType;<br>
typedef itk::Image< OutputPixelType, Dimension > OutputImageType;<br>
typedef itk::Image< BinaryPixelType, Dimension > BinaryImageType;<br>
<br>
<br>
typedef itk::GradientMagnitudeImageFilter< InternalImageType, OutputImageType > GradientMagnitudeFilterType;<br>
GradientMagnitudeFilterType::Pointer gradientMagnFilter = GradientMagnitudeFilterType::New();<br>
// typedef itk::GradientMagnitudeRecursiveGaussianImageFilter<InternalImageType, OutputImageType> GradientGaussianFilterType;<br>
<br>
typedef itk::ImageFileReader< InternalImageType > ReaderType;<br>
//typedef itk::ImageFileReader< BinaryImageType > BinaryReaderType;<br>
typedef itk::ImageFileWriter< OutputImageType > WriterType;<br>
// typedef itk::ImageFileWriter< BinaryImageType > BinaryWriterType;<br>
ReaderType::Pointer imageReader = ReaderType::New();<br>
<br>
argv[1] = "/square.nii";<br>
<br>
imageReader -> SetFileName( argv[1] );<br>
<br>
typedef itk::CovariantVector< float, Dimension > GradientPixelType;<br>
typedef itk::Image< GradientPixelType, Dimension > GradientImageType;<br>
//typedef itk::GradientVectorFlowImageFilter< GradientImageType, GradientImageType > GVFFilterType;<br>
//typedef itk::LaplacianImageFilter<InternalImageType, InternalImageType> LaplacianFilterType;<br>
typedef itk::GradientImageFilter<InternalImageType, float, float> GradientImageFilterType;<br>
typedef itk::GradientMagnitudeImageFilter<InternalImageType, InternalImageType> GradientMagnitudeImageFilterType;<br>
<br>
GradientMagnitudeImageFilterType::Pointer gmFilter = GradientMagnitudeImageFilterType::New();<br>
gmFilter->SetInput(imageReader->GetOutput());<br>
<br>
//GradientGaussianFilterType::Pointer gaussianFilter = GradientGaussianFilterType::New();<br>
//gaussianFilter->SetSigma( 0.5 );<br>
//gaussianFilter->SetInput( imageReader->GetOutput() );<br>
<br>
GradientImageFilterType::Pointer gFilter = GradientImageFilterType::New();<br>
gFilter->SetInput(gmFilter->GetOutput());<br>
<br>
try<br>
{<br>
// gmFilter->Update();<br>
//gaussianFilter->Update();<br>
gFilter->Update();<br>
}<br>
catch( itk::ExceptionObject & excep )<br>
{<br>
std::cerr << "Exception caught !" << std::endl;<br>
std::cerr << excep << std::endl;<br>
}<br>
<br>
<br>
WriterType::Pointer gmWriter = WriterType::New();<br>
gmWriter->SetFileName("/gradientMagnitudeSquare.nii");<br>
gmWriter->SetInput(gmFilter->GetOutput());<br>
<br>
//WriterType::Pointer gaussianWriter = WriterType::New();<br>
//gaussianWriter->SetFileName("/gradientMagnitudeGaussianSquare.nii");<br>
//gaussianWriter->SetInput(gaussianFilter->GetOutput());<br>
<br>
typedef itk::ImageFileWriter< GradientImageType > FieldWriterType;<br>
FieldWriterType::Pointer gradientWriter = FieldWriterType::New();<br>
gradientWriter->SetInput(gFilter->GetOutput());<br>
gradientWriter->SetFileName("/gradientSquare.vtk");<br>
try<br>
{<br>
gmWriter->Update();<br>
gradientWriter->Update();<br>
}<br>
catch( itk::ExceptionObject & excep )<br>
{<br>
std::cerr << "Exception caught !" << std::endl;<br>
std::cerr << excep << std::endl;<br>
}<br>
<br>
<br>
}<br>
<br>
<br>_____________________________________<br>
Powered by <a href="http://www.kitware.com" target="_blank">www.kitware.com</a><br>
<br>
Visit other Kitware open-source projects at<br>
<a href="http://www.kitware.com/opensource/opensource.html" target="_blank">http://www.kitware.com/opensource/opensource.html</a><br>
<br>
Please keep messages on-topic and check the ITK FAQ at: <a href="http://www.itk.org/Wiki/ITK_FAQ" target="_blank">http://www.itk.org/Wiki/ITK_FAQ</a><br>
<br>
Follow this link to subscribe/unsubscribe:<br>
<a href="http://www.itk.org/mailman/listinfo/insight-users" target="_blank">http://www.itk.org/mailman/listinfo/insight-users</a><br>
<br></blockquote></div><br>