[Insight-users] Gradient of Vector Image : ImageAdaptor

[Luis Ibanez]

This is a multi-part message in MIME format.
--------------090708050306040108080106
Content-Type: text/plain; charset=us-ascii; format=flowed
Content-Transfer-Encoding: 7bit


Hi Bing,


My mistake,

Instead of defining the input image type of the gradient
filter using the ExternalType of the PixelAccessor, you
should simple USE the OffsetImageAdaptor directly as a
input image type.


Like:



  typedef
     itk::GradientRecursiveGaussianImageFilter<
                                OffsetImageAdaptorType,
                                VectorImageType
                                    > GradientFilterType;



That is, you think of the "OffsetImageAdaptor"
as an "ImageType".


Please find attached a modified version of
the Example ImageAdaptor3.cxx that will do
exactly what you want.

This code is compiling fine on VC++7.0, you will
have to plug in your Offset image in order to be
able to run it.


Regards,


    Luis



--------------------
Bing Jian wrote:
> Hi Luis,
> 
> I am still confusing.
> 
> Below is what I got from gradient->SetInput(adaptor)
> ===================================================================
> cannot convert parameter 1 from 'class itk::SmartPointer<class
> itk::ImageAdaptor<class itk::Image<class itk::Vector<float,3>,3>,class
> VectorPixelAccessor> >' to 'const class itk::Image<float,3> *'
> ===============================================================
> In class PixelAccessor I have
>     typedef itk::Vector<float,3>   InternalType;
>     typedef               float      ExternalType;
> 
> In the cxx file, I also have
> //typedef float ScalarPixelType;
> typedef VectorPixelAccessor::ExternalType ScalarPixelType;
> typedef itk::Image<ScalarPixelType, 3> InputScalarImageType;
> typedef InputScalarImageType::Pointer InputScalarImagePointer;
> 
> typedef itk::Vector<ScalarPixelType, 3> InputVectorPixelType;
> typedef itk::Image<InputVectorPixelType, 3> InputVectorImageType;
> typedef InputVectorImageType::Pointer InputVectorImagePointer;
> 
> typedef itk::CovariantVector<float, 3> GradientVectorPixelType;
> typedef itk::Image<GradientVectorPixelType, 3>
> GradientVectorImageType;
> typedef GradientVectorImageType::Pointer GradientVectorImagePointer;
> 
> typedef itk::Matrix<float, 3, 3> JacobianMatrixPixelType;
> typedef itk::Image<JacobianMatrixPixelType, 3>
> JacobianMatrixImageType;
> typedef JacobianMatrixImageType::Pointer JacobianMatrixImagePointer;
> 
> The filter is defined as follows:
> 	typedef
> itk::GradientRecursiveGaussianImageFilter<InputScalarImageType,GradientVectorImageType>
> GradientFilterType;
> 	typedef GradientFilterType::Pointer GradientFilterPointer;
> 	GradientFilterPointer gradient = GradientFilterType::New();
> 
> The adaptor is defined as follows:
> 	typedef itk::ImageAdaptor<  InputVectorImageType,
> 		VectorPixelAccessor > ImageAdaptorType;
> 	ImageAdaptorType::Pointer adaptor = ImageAdaptorType::New();
> 
> So both the input type of gradient filter and externel type of
> VectorPixelAccessor are InputScalarImageType (float). and both
> the input type of adaptor and internal type of VectorPixelAccessor
> are vector<float,3>. But I still got above error.
> 
> Thanks a lot.
> 


=====================================================================



--------------090708050306040108080106
Content-Type: text/plain;
 name="ImageAdaptor.cxx"
Content-Transfer-Encoding: 7bit
Content-Disposition: inline;
 filename="ImageAdaptor.cxx"

/*=========================================================================

  Program:   Insight Segmentation & Registration Toolkit
  Module:    $RCSfile: ImageAdaptor3.cxx,v $
  Language:  C++
  Date:      $Date: 2003/09/10 14:29:51 $
  Version:   $Revision: 1.10 $

  Copyright (c) Insight Software Consortium. All rights reserved.
  See ITKCopyright.txt or http://www . itk . org/HTML/Copyright . htm for details.

     This software is distributed WITHOUT ANY WARRANTY; without even 
     the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR 
     PURPOSE.  See the above copyright notices for more information.

=========================================================================*/

#if defined(_MSC_VER)
#pragma warning ( disable : 4786 )
#endif


#include "itkImage.h"
#include "itkImageAdaptor.h"
#include "itkCovariantVector.h"
#include "itkImageRegionIteratorWithIndex.h"
#include "itkGradientRecursiveGaussianImageFilter.h"



//
// ADAPTOR FOR THE OFFSET IMAGES INTO SCALARS
//
class OffsetPixelAccessor  
{
public:
  typedef itk::Offset<3>   InternalType;
  typedef       float      ExternalType;

  void operator=( const OffsetPixelAccessor & vpa )
    {
      m_Index = vpa.m_Index;
    }
  ExternalType Get( const InternalType & input ) const 
    {
      return static_cast<ExternalType>( input[ m_Index ] );
    }
  void SetIndex( unsigned int index )
    {
      m_Index = index;
    }
private:
  unsigned int m_Index;
};





//
// ADAPTOR FOR THE VECTOR IMAGES INTO SCALARS
//
class VectorPixelAccessor  
{
public:
  typedef itk::CovariantVector<float,3>   InternalType;
  typedef                      float      ExternalType;

  void operator=( const VectorPixelAccessor & vpa )
    {
      m_Index = vpa.m_Index;
    }
  ExternalType Get( const InternalType & input ) const 
    {
      return static_cast<ExternalType>( input[ m_Index ] );
    }
  void SetIndex( unsigned int index )
    {
      m_Index = index;
    }
private:
  unsigned int m_Index;
};






int main( int argc, char *argv[] ) 
{

  const   unsigned int   Dimension = 3;

  typedef itk::Offset<Dimension>                     OffsetPixelType;
  typedef itk::Image< OffsetPixelType, Dimension >   OffsetImageType;
  typedef OffsetImageType::Pointer                   OffsetImagePointer;

  OffsetImagePointer  offsetImage = OffsetImageType::New();

  typedef itk::ImageAdaptor<  
                         OffsetImageType, 
                         OffsetPixelAccessor 
                                           >         OffsetImageAdaptorType;

  OffsetImageAdaptorType::Pointer  offsetAdaptor =
                                         OffsetImageAdaptorType::New();

  typedef itk::CovariantVector< float, Dimension  >   VectorPixelType; 
  typedef itk::Image< VectorPixelType, Dimension  >   VectorImageType;
  typedef VectorImageType::Pointer                    VectorImagePointer;

  typedef itk::GradientRecursiveGaussianImageFilter< 
                                        OffsetImageAdaptorType,
                                        VectorImageType> GradientFilterType;

  GradientFilterType::Pointer gradient = GradientFilterType::New();

  typedef itk::ImageAdaptor<  VectorImageType, 
                              VectorPixelAccessor > ImageAdaptorType;

  ImageAdaptorType::Pointer adaptor = ImageAdaptorType::New();

  VectorPixelAccessor  accessor;
  

  // CREATE THE N * N image for the Jacobian of the Offset
  typedef itk::Image< float, Dimension >  ScalarOutputImageType;
  typedef ScalarOutputImageType::Pointer  ScalarOutputImagePointer;

  ScalarOutputImagePointer  outputImages[ 9 ];
 
  for(unsigned int i=0; i<9; i++)
    {
    outputImages[i] = ScalarOutputImageType::New();
    outputImages[i]->SetRegions( offsetImage->GetBufferedRegion() );
    outputImages[i]->Allocate();
    }


  for(unsigned int indexToProcess = 0; indexToProcess< Dimension; indexToProcess++)
    {
    accessor.SetIndex( indexToProcess );
    adaptor->SetPixelAccessor( accessor );
    gradient->SetInput( offsetAdaptor );
    gradient->Update();
    adaptor->SetImage( gradient->GetOutput() );
    VectorImagePointer  vectorImage = gradient->GetOutput();
    for(unsigned int component=0; component < Dimension; component++)
      {
      ScalarOutputImagePointer outputImage = outputImages[ component ];
      typedef itk::ImageRegionIteratorWithIndex< ScalarOutputImageType > ScalarIterator;
      typedef itk::ImageRegionIteratorWithIndex< VectorImageType       > VectorIterator;
      ScalarIterator scalarIter( outputImage, outputImage->GetBufferedRegion() );
      VectorIterator vectorIter( vectorImage, vectorImage->GetBufferedRegion() );
      scalarIter.GoToBegin();
      vectorIter.GoToBegin();
      while( !scalarIter.IsAtEnd() ) 
        {
        scalarIter.Set( vectorIter.Value()[ component ] );
        ++scalarIter;
        ++vectorIter;
        }
      }
    }
 


  return 0;
}




--------------090708050306040108080106--