itkImageBase.h

Go to the documentation of this file.

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

  Program:   Insight Segmentation & Registration Toolkit
  Module:    $RCSfile: itkImageBase.h,v $
  Language:  C++
  Date:      $Date: 2009-02-18 17:40:55 $
  Version:   $Revision: 1.76 $

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

  Portions of this code are covered under the VTK copyright.
  See VTKCopyright.txt or http://www.kitware.com/VTKCopyright.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.

=========================================================================*/
#ifndef __itkImageBase_h
#define __itkImageBase_h

#include "itkDataObject.h"

#include "itkImageRegion.h"
#include "itkIndex.h"
#include "itkObjectFactory.h"
#include "itkOffset.h"
#include "itkPoint.h"
#include "itkSize.h"
#include "itkFixedArray.h"
#include "itkPoint.h"
#include "itkMatrix.h"
#include "itkContinuousIndex.h"
#include "itkImageHelper.h"
#include <vnl/vnl_matrix_fixed.txx>

#include "itkImageRegion.h"

#ifdef ITK_USE_TEMPLATE_META_PROGRAMMING_LOOP_UNROLLING
#include "itkImageTransformHelper.h"
#endif

namespace itk
{

template <typename TImage>
struct GetImageDimension
{
  itkStaticConstMacro(ImageDimension, unsigned int, TImage::ImageDimension);
};

template<unsigned int VImageDimension=2>
class ITK_EXPORT ImageBase : public DataObject
{
public:
  typedef ImageBase                     Self;
  typedef DataObject                    Superclass;
  typedef SmartPointer<Self>            Pointer;
  typedef SmartPointer<const Self>      ConstPointer;

  itkNewMacro(Self);

  itkTypeMacro(ImageBase, DataObject);

  itkStaticConstMacro(ImageDimension, unsigned int, VImageDimension );

  typedef Index<VImageDimension>                  IndexType;
  typedef typename IndexType::IndexValueType      IndexValueType;

  typedef Offset<VImageDimension>                 OffsetType;
  typedef typename OffsetType::OffsetValueType    OffsetValueType;

  typedef Size<VImageDimension>                   SizeType;
  typedef typename SizeType::SizeValueType        SizeValueType;

  typedef ImageRegion<VImageDimension>            RegionType;

  typedef Vector<double, VImageDimension>         SpacingType;

  typedef Point<double, VImageDimension> PointType;

  typedef Matrix<double, VImageDimension, VImageDimension> DirectionType;

  void Initialize();

  static unsigned int GetImageDimension()
    { return VImageDimension; }

  itkSetMacro(Origin, PointType);
  virtual void SetOrigin( const double origin[VImageDimension] );
  virtual void SetOrigin( const float origin[VImageDimension] );

  virtual void SetDirection( const DirectionType direction );

  itkGetConstReferenceMacro(Direction, DirectionType);

  itkGetConstReferenceMacro(Spacing, SpacingType);

  itkGetConstReferenceMacro(Origin, PointType);

  virtual void SetLargestPossibleRegion(const RegionType &region);

  virtual const RegionType& GetLargestPossibleRegion() const
    { return m_LargestPossibleRegion;};

  virtual void SetBufferedRegion(const RegionType &region);

  virtual const RegionType& GetBufferedRegion() const
  { return m_BufferedRegion;};

  virtual void SetRequestedRegion(const RegionType &region);

  virtual void SetRequestedRegion(DataObject *data);

  virtual const RegionType& GetRequestedRegion() const
  { return m_RequestedRegion;};

  const OffsetValueType *GetOffsetTable() const { return m_OffsetTable; };

#ifdef ITK_USE_TEMPLATE_META_PROGRAMMING_LOOP_UNROLLING
  inline OffsetValueType ComputeOffset(const IndexType &ind) const
    {
    OffsetValueType offset = 0;
    ImageHelper<VImageDimension,VImageDimension>::ComputeOffset(this->GetBufferedRegion().GetIndex(),
                                                                ind,
                                                                m_OffsetTable,
                                                                offset);
    return offset;
    }
#else
  OffsetValueType ComputeOffset(const IndexType &ind) const
    {
    // need to add bounds checking for the region/buffer?
    OffsetValueType offset=0;
    const IndexType &bufferedRegionIndex = this->GetBufferedRegion().GetIndex();

    // data is arranged as [][][][slice][row][col]
    // with Index[0] = col, Index[1] = row, Index[2] = slice
    for (int i=VImageDimension-1; i > 0; i--)
      {
      offset += (ind[i] - bufferedRegionIndex[i])*m_OffsetTable[i];
      }
    offset += (ind[0] - bufferedRegionIndex[0]);

    return offset;
    }
#endif

#ifdef ITK_USE_TEMPLATE_META_PROGRAMMING_LOOP_UNROLLING
  inline IndexType ComputeIndex(OffsetValueType offset) const
    {
    IndexType index;
    const IndexType &bufferedRegionIndex = this->GetBufferedRegion().GetIndex();
    ImageHelper<VImageDimension,VImageDimension>::ComputeIndex(bufferedRegionIndex,
                                                               offset,
                                                               m_OffsetTable,
                                                               index);
    return index;
    }
#else
  IndexType ComputeIndex(OffsetValueType offset) const
    {
    IndexType index;
    const IndexType &bufferedRegionIndex = this->GetBufferedRegion().GetIndex();

    for (int i=VImageDimension-1; i > 0; i--)
      {
      index[i] = static_cast<IndexValueType>(offset / m_OffsetTable[i]);
      offset -= (index[i] * m_OffsetTable[i]);
      index[i] += bufferedRegionIndex[i];
      }
    index[0] = bufferedRegionIndex[0] + static_cast<IndexValueType>(offset);

    return index;
    }
#endif

  virtual void SetSpacing (const SpacingType & spacing);
  virtual void SetSpacing (const double spacing[VImageDimension]);
  virtual void SetSpacing (const float spacing[VImageDimension]);


#ifdef ITK_USE_TEMPLATE_META_PROGRAMMING_LOOP_UNROLLING
  template<class TCoordRep>
  bool TransformPhysicalPointToIndex(
    const Point<TCoordRep, VImageDimension>& point,
    IndexType & index ) const
    {
      ImageTransformHelper<VImageDimension,VImageDimension-1,VImageDimension-1>::TransformPhysicalPointToIndex(
        this->m_PhysicalPointToIndex, this->m_Origin, point, index);

    // Now, check to see if the index is within allowed bounds
    const bool isInside = this->GetLargestPossibleRegion().IsInside( index );
    return isInside;
    }
#else
  template<class TCoordRep>
  bool TransformPhysicalPointToIndex(
            const Point<TCoordRep, VImageDimension>& point,
            IndexType & index                                ) const
    {
    for (unsigned int i = 0; i < VImageDimension; i++)
      {
      TCoordRep sum = NumericTraits<TCoordRep>::Zero;
      for (unsigned int j = 0; j < VImageDimension; j++)
        {
        sum += this->m_PhysicalPointToIndex[i][j] * (point[j] - this->m_Origin[j]);
        }
      index[i] = static_cast< IndexValueType>( sum );
      }

    // Now, check to see if the index is within allowed bounds
    const bool isInside = this->GetLargestPossibleRegion().IsInside( index );

    return isInside;
    }
#endif

  template<class TCoordRep>
  bool TransformPhysicalPointToContinuousIndex(
              const Point<TCoordRep, VImageDimension>& point,
              ContinuousIndex<TCoordRep, VImageDimension>& index   ) const
    {
    Vector<double, VImageDimension> cvector;

    for( unsigned int k = 0; k < VImageDimension; k++ )
      {
      cvector[k] = point[k] - this->m_Origin[k];
      }
    cvector = m_PhysicalPointToIndex * cvector;
    for( unsigned int i = 0; i < VImageDimension; i++ )
      {
      index[i] = static_cast<TCoordRep>(cvector[i]);
      }

    // Now, check to see if the index is within allowed bounds
    const bool isInside = this->GetLargestPossibleRegion().IsInside( index );

    return isInside;
    }


  template<class TCoordRep>
  void TransformContinuousIndexToPhysicalPoint(
            const ContinuousIndex<TCoordRep, VImageDimension>& index,
            Point<TCoordRep, VImageDimension>& point        ) const
    {
    for( unsigned int r=0; r<VImageDimension; r++)
      {
      TCoordRep sum = NumericTraits<TCoordRep>::Zero;
      for( unsigned int c=0; c<VImageDimension; c++ )
        {
        sum += this->m_IndexToPhysicalPoint(r,c) * index[c];
        }
      point[r] = sum + this->m_Origin[r];
      }
    }

#ifdef ITK_USE_TEMPLATE_META_PROGRAMMING_LOOP_UNROLLING
  template<class TCoordRep>
  void TransformIndexToPhysicalPoint(
                      const IndexType & index,
                      Point<TCoordRep, VImageDimension>& point ) const
    {
      ImageTransformHelper<VImageDimension,VImageDimension-1,VImageDimension-1>::TransformIndexToPhysicalPoint(
        this->m_IndexToPhysicalPoint, this->m_Origin, index, point);
    }
#else
  template<class TCoordRep>
  void TransformIndexToPhysicalPoint(
                      const IndexType & index,
                      Point<TCoordRep, VImageDimension>& point ) const
    {
    for (unsigned int i = 0; i < VImageDimension; i++)
      {
      point[i] = this->m_Origin[i];
      for (unsigned int j = 0; j < VImageDimension; j++)
        {
        point[i] += m_IndexToPhysicalPoint[i][j] * index[j];
        }
      }
    }
#endif



  template<class TCoordRep>
  void TransformLocalVectorToPhysicalVector(
    const FixedArray<TCoordRep, VImageDimension> & inputGradient,
          FixedArray<TCoordRep, VImageDimension> & outputGradient ) const
    {
    //
    // This temporary implementation should be replaced with Template MetaProgramming.
    // 
#ifdef ITK_USE_ORIENTED_IMAGE_DIRECTION
    const DirectionType & direction = this->GetDirection();
    for (unsigned int i = 0; i < VImageDimension; i++)
      {
      typedef typename NumericTraits<TCoordRep>::AccumulateType CoordSumType;
      CoordSumType sum = NumericTraits<CoordSumType>::Zero;
      for (unsigned int j = 0; j < VImageDimension; j++)
        {
        sum += direction[i][j] * inputGradient[j];
        }
      outputGradient[i] = static_cast<TCoordRep>( sum );
      }
#else
    for (unsigned int i = 0; i < VImageDimension; i++)
      {
      outputGradient[i] = inputGradient[i];
      }
#endif
    }

  virtual void CopyInformation(const DataObject *data);

  virtual void Graft(const DataObject *data);

  virtual void UpdateOutputInformation();

  virtual void SetRequestedRegionToLargestPossibleRegion();

  virtual bool RequestedRegionIsOutsideOfTheBufferedRegion();

  virtual bool VerifyRequestedRegion();

  virtual unsigned int GetNumberOfComponentsPerPixel() const;
  virtual void SetNumberOfComponentsPerPixel( unsigned int );

protected:
  ImageBase();
  ~ImageBase();
  virtual void PrintSelf(std::ostream& os, Indent indent) const;

  void ComputeOffsetTable();

  virtual void ComputeIndexToPhysicalPointMatrices();

protected:
  SpacingType         m_Spacing;
  PointType           m_Origin;
  DirectionType       m_Direction;

  DirectionType       m_IndexToPhysicalPoint;
  DirectionType       m_PhysicalPointToIndex;

  virtual void InitializeBufferedRegion(void);

private:
  ImageBase(const Self&); //purposely not implemented
  void operator=(const Self&); //purposely not implemented

  OffsetValueType  m_OffsetTable[VImageDimension+1];

  RegionType          m_LargestPossibleRegion;
  RegionType          m_RequestedRegion;
  RegionType          m_BufferedRegion;

};

} // end namespace itk

// Define instantiation macro for this template.
#define ITK_TEMPLATE_ImageBase(_, EXPORT, x, y) namespace itk { \
  _(1(class EXPORT ImageBase< ITK_TEMPLATE_1 x >)) \
  namespace Templates { typedef ImageBase< ITK_TEMPLATE_1 x > ImageBase##y; } \
  }

#if ITK_TEMPLATE_EXPLICIT
# include "Templates/itkImageBase+-.h"
#endif

#if ITK_TEMPLATE_TXX
# include "itkImageBase.txx"
#endif

#endif