itkMultiBSplineDeformableTransformWithNormal.h

Go to the documentation of this file.

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

  This file is part of the elastix software.

  Copyright (c) University Medical Center Utrecht. All rights reserved.
  See src/CopyrightElastix.txt or http://elastix.isi.uu.nl/legal.php 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 __itkMultiBSplineDeformableTransformWithNormal_h
#define __itkMultiBSplineDeformableTransformWithNormal_h

#include "itkAdvancedBSplineDeformableTransform.h"
#include "itkNearestNeighborInterpolateImageFunction.h"

namespace itk
{
template<
class TScalarType         = double,      // Data type for scalars
unsigned int NDimensions  = 3,           // Number of dimensions
unsigned int VSplineOrder = 3 >
// Spline order
class MultiBSplineDeformableTransformWithNormal :
  public AdvancedTransform< TScalarType, NDimensions, NDimensions >
{
public:

  typedef MultiBSplineDeformableTransformWithNormal Self;
  typedef AdvancedTransform<
    TScalarType, NDimensions, NDimensions >         Superclass;
  typedef SmartPointer< Self >       Pointer;
  typedef SmartPointer< const Self > ConstPointer;

  itkNewMacro( Self );

  itkTypeMacro( MultiBSplineDeformableTransformWithNormal, AdvancedTransform );

  itkStaticConstMacro( SpaceDimension, unsigned int, NDimensions );

  itkStaticConstMacro( SplineOrder, unsigned int, VSplineOrder );

  typedef typename Superclass::ParametersType         ParametersType;
  typedef typename Superclass::NumberOfParametersType NumberOfParametersType;
  typedef typename Superclass::JacobianType           JacobianType;
  typedef typename Superclass::ScalarType             ScalarType;
  typedef typename Superclass::InputPointType         InputPointType;
  typedef typename Superclass::OutputPointType        OutputPointType;
  typedef typename Superclass::InputVectorType        InputVectorType;
  typedef typename Superclass::OutputVectorType       OutputVectorType;
  typedef typename Superclass::InputVnlVectorType     InputVnlVectorType;
  typedef typename Superclass::OutputVnlVectorType    OutputVnlVectorType;
  typedef typename Superclass::InputCovariantVectorType
    InputCovariantVectorType;
  typedef typename Superclass::OutputCovariantVectorType
    OutputCovariantVectorType;

  typedef typename Superclass
    ::NonZeroJacobianIndicesType NonZeroJacobianIndicesType;
  typedef typename Superclass::SpatialJacobianType SpatialJacobianType;
  typedef typename Superclass
    ::JacobianOfSpatialJacobianType JacobianOfSpatialJacobianType;
  typedef typename Superclass::SpatialHessianType SpatialHessianType;
  typedef typename Superclass
    ::JacobianOfSpatialHessianType JacobianOfSpatialHessianType;
  typedef typename Superclass::InternalMatrixType InternalMatrixType;

  typedef BSplineInterpolationWeightFunction2< ScalarType,
    itkGetStaticConstMacro( SpaceDimension ),
    itkGetStaticConstMacro( SplineOrder ) >                  WeightsFunctionType;
  typedef typename WeightsFunctionType::WeightsType WeightsType;

  void SetParameters( const ParametersType & parameters );

  void SetFixedParameters( const ParametersType & parameters );

  void SetParametersByValue( const ParametersType & parameters );

  void SetIdentity( void );

  virtual const ParametersType & GetParameters( void ) const;

  virtual const ParametersType & GetFixedParameters( void ) const;

  typedef typename ParametersType::ValueType PixelType;
  typedef Image< PixelType,
    itkGetStaticConstMacro( SpaceDimension ) >           ImageType;
  typedef typename ImageType::Pointer ImagePointer;

  //virtual ImagePointer * GetCoefficientImage( void )
  //  { return this->m_CoefficientImage; }
  //virtual const ImagePointer * GetCoefficientImage( void ) const
  //  { return this->m_CoefficientImage; }

  // virtual void SetCoefficientImage( ImagePointer images[] );

  typedef ImageRegion< itkGetStaticConstMacro( SpaceDimension ) > RegionType;

  typedef typename RegionType::IndexType    IndexType;
  typedef typename RegionType::SizeType     SizeType;
  typedef typename ImageType::SpacingType   SpacingType;
  typedef typename ImageType::DirectionType DirectionType;
  typedef typename ImageType::PointType     OriginType;
  typedef IndexType                         GridOffsetType;

  virtual void SetGridRegion( const RegionType & region );

  virtual RegionType GetGridRegion( void ) const;

  virtual void SetGridSpacing( const SpacingType & spacing );

  virtual SpacingType GetGridSpacing( void ) const;

  virtual void SetGridDirection( const DirectionType & spacing );

  virtual DirectionType GetGridDirection( void ) const;

  virtual void SetGridOrigin( const OriginType & origin );

  virtual OriginType GetGridOrigin( void ) const;

  typedef Image< unsigned char,
    itkGetStaticConstMacro( SpaceDimension ) >          ImageLabelType;
  typedef typename ImageLabelType::Pointer ImageLabelPointer;

  typedef itk::NearestNeighborInterpolateImageFunction<
    ImageLabelType, TScalarType >                        ImageLabelInterpolator;
  typedef typename ImageLabelInterpolator::Pointer ImageLabelInterpolatorPointer;

  typedef Vector< TScalarType, itkGetStaticConstMacro( SpaceDimension ) > VectorType;
  typedef Vector< VectorType, itkGetStaticConstMacro( SpaceDimension ) >  BaseType;
  typedef Image< VectorType, itkGetStaticConstMacro( SpaceDimension ) >   ImageVectorType;
  typedef typename ImageVectorType::Pointer                               ImageVectorPointer;
  typedef Image< BaseType, itkGetStaticConstMacro( SpaceDimension ) >     ImageBaseType;
  typedef typename ImageBaseType::Pointer                                 ImageBasePointer;

  void SetLabels( ImageLabelType * labels );

  itkGetMacro( Labels, ImageLabelType * );

  itkGetConstMacro( NbLabels, unsigned char );

  void UpdateLocalBases( void );

  itkGetMacro( LocalBases, ImageBaseType * );

  typedef Array< unsigned long > ParameterIndexArrayType;

  virtual OutputVectorType TransformVector( const InputVectorType & ) const
  {
    itkExceptionMacro( << "Method not applicable for deformable transform." );
    return OutputVectorType();
  }


  virtual OutputVnlVectorType TransformVector( const InputVnlVectorType & ) const
  {
    itkExceptionMacro( << "Method not applicable for deformable transform. " );
    return OutputVnlVectorType();
  }


  virtual OutputCovariantVectorType TransformCovariantVector(
    const InputCovariantVectorType & ) const
  {
    itkExceptionMacro( << "Method not applicable for deformable transform. " );
    return OutputCovariantVectorType();
  }


  virtual NumberOfParametersType GetNumberOfParameters( void ) const;

  virtual NumberOfParametersType GetNumberOfParametersPerDimension( void ) const;

  virtual const RegionType & GetValidRegion( void )
  {
    return m_Trans[ 0 ]->GetValidRegion();
  }


  virtual bool IsLinear( void ) const { return false; }

  virtual unsigned long GetNumberOfWeights( void ) const
  {
    return m_Trans[ 0 ]->m_WeightsFunction->GetNumberOfWeights();
  }


  virtual unsigned int GetNumberOfAffectedWeights( void ) const
  {
    return m_Trans[ 0 ]->m_WeightsFunction->GetNumberOfWeights();
  }


  virtual NumberOfParametersType GetNumberOfNonZeroJacobianIndices( void ) const
  {
    return m_Trans[ 0 ]->m_WeightsFunction->GetNumberOfWeights() * SpaceDimension;
  }


  virtual bool GetHasNonZeroSpatialJacobian( void ) const
  {
    return true;
  }


  virtual bool HasNonZeroJacobianOfSpatialJacobian( void ) const
  {
    return true;
  }


  virtual bool GetHasNonZeroSpatialHessian( void ) const
  {
    return true;
  }


  virtual bool HasNonZeroJacobianOfSpatialHessian( void ) const
  {
    return true;
  }


  typedef ContinuousIndex< ScalarType, SpaceDimension > ContinuousIndexType;

  OutputPointType TransformPoint( const InputPointType & point ) const;

  //virtual const JacobianType & GetJacobian( const InputPointType & point ) const;

  virtual void GetJacobian(
    const InputPointType & ipp,
    JacobianType & j,
    NonZeroJacobianIndicesType & ) const;

  virtual void GetSpatialJacobian(
    const InputPointType & ipp,
    SpatialJacobianType & sj ) const;

  virtual void GetJacobianOfSpatialJacobian(
    const InputPointType & ipp,
    JacobianOfSpatialJacobianType & jsj,
    NonZeroJacobianIndicesType & nonZeroJacobianIndices ) const;

  virtual void GetJacobianOfSpatialJacobian(
    const InputPointType &,
    SpatialJacobianType &,
    JacobianOfSpatialJacobianType &,
    NonZeroJacobianIndicesType & ) const;

  virtual void GetSpatialHessian(
    const InputPointType & ipp,
    SpatialHessianType & sh ) const;

  virtual void GetJacobianOfSpatialHessian(
    const InputPointType & ipp,
    JacobianOfSpatialHessianType & jsh,
    NonZeroJacobianIndicesType & nonZeroJacobianIndices ) const
  {
    itkExceptionMacro( << "ERROR: GetJacobianOfSpatialHessian() not yet implemented "
                       << "in the MultiBSplineDeformableTransformWithNormal class." );
  }


  virtual void GetJacobianOfSpatialHessian(
    const InputPointType &,
    SpatialHessianType &,
    JacobianOfSpatialHessianType &,
    NonZeroJacobianIndicesType & ) const;

protected:

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

  MultiBSplineDeformableTransformWithNormal();
  virtual ~MultiBSplineDeformableTransformWithNormal();

  // void WrapAsImages( void );

  /*
  void TransformPointToContinuousGridIndex(
   const InputPointType & point, ContinuousIndexType & index ) const;

  virtual void ComputeNonZeroJacobianIndices(
    NonZeroJacobianIndicesType & nonZeroJacobianIndices,
    const RegionType & supportRegion ) const;
    */

  //virtual bool InsideValidRegion( const ContinuousIndexType& index ) const;

  // BulkTransformPointer  m_BulkTransform;

  // ImagePointer    m_CoefficientImage[ NDimensions ];

  // RegionType          m_GridRegion;
  // SpacingType         m_GridSpacing;
  // DirectionType       m_GridDirection;
  // OriginType          m_GridOrigin;
  // GridOffsetType      m_GridOffsetTable;

  // DirectionType       m_PointToIndexMatrix;
  // SpatialJacobianType m_PointToIndexMatrix2;
  // DirectionType       m_PointToIndexMatrixTransposed;
  // SpatialJacobianType m_PointToIndexMatrixTransposed2;
  // DirectionType       m_IndexToPoint;

  // RegionType      m_ValidRegion;

  // unsigned long   m_Offset;
  // SizeType        m_SupportSize;
  // ContinuousIndexType m_ValidRegionBegin;
  // ContinuousIndexType m_ValidRegionEnd;

  // bool m_SplineOrderOdd;

  const ParametersType * m_InputParametersPointer;

  /*
  typedef typename JacobianType::ValueType      JacobianPixelType;
  typedef Image< JacobianPixelType,
    itkGetStaticConstMacro( SpaceDimension ) >  JacobianImageType;

  typename JacobianImageType::Pointer m_JacobianImage[ NDimensions ];
  */

  // mutable IndexType m_LastJacobianIndex;

  // ImagePointer    m_WrappedImage[ NDimensions ];

  ParametersType m_InternalParametersBuffer;

  typedef AdvancedBSplineDeformableTransform< TScalarType,
    itkGetStaticConstMacro( SpaceDimension ),
    itkGetStaticConstMacro( SplineOrder ) >          TransformType;

  unsigned char                                  m_NbLabels;
  ImageLabelPointer                              m_Labels;
  ImageLabelInterpolatorPointer                  m_LabelsInterpolator;
  ImageVectorPointer                             m_LabelsNormals;
  std::vector< typename TransformType::Pointer > m_Trans;
  std::vector< ParametersType >                  m_Para;
  mutable int                                    m_LastJacobian;
  ImageBasePointer                               m_LocalBases;

private:

  MultiBSplineDeformableTransformWithNormal( const Self & ); // purposely not implemented
  void operator=( const Self & );                            // purposely not implemented

  void DispatchParameters( const ParametersType & parameters );

  void PointToLabel( const InputPointType & p, int & l ) const;

};

} // end namespace itk

#ifndef ITK_MANUAL_INSTANTIATION
#include "itkMultiBSplineDeformableTransformWithNormal.hxx"
#endif

#endif // end __itkMultiBSplineDeformableTransformWithNormal_h