itkKernelTransform.h

Go to the documentation of this file.

/*=========================================================================
 *
 *  Copyright Insight Software Consortium
 *
 *  Licensed under the Apache License, Version 2.0 (the "License");
 *  you may not use this file except in compliance with the License.
 *  You may obtain a copy of the License at
 *
 *         http://www.apache.org/licenses/LICENSE-2.0.txt
 *
 *  Unless required by applicable law or agreed to in writing, software
 *  distributed under the License is distributed on an "AS IS" BASIS,
 *  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 *  See the License for the specific language governing permissions and
 *  limitations under the License.
 *
 *=========================================================================*/
#ifndef __itkKernelTransform_h
#define __itkKernelTransform_h

#include "itkTransform.h"
#include "itkMatrix.h"
#include "itkPointSet.h"
#include <deque>
#include <math.h>
#include "vnl/vnl_matrix_fixed.h"
#include "vnl/vnl_matrix.h"
#include "vnl/vnl_vector.h"
#include "vnl/vnl_vector_fixed.h"
#include "vnl/algo/vnl_svd.h"
#include "vnl/vnl_sample.h"

namespace itk
{
template <class TScalarType, // probably only float and double make sense here
          unsigned int NDimensions>
// Number of dimensions
class ITK_EXPORT KernelTransform :
  public Transform<TScalarType, NDimensions, NDimensions>
{
public:
  typedef KernelTransform                                  Self;
  typedef Transform<TScalarType, NDimensions, NDimensions> Superclass;
  typedef SmartPointer<Self>                               Pointer;
  typedef SmartPointer<const Self>                         ConstPointer;

  itkTypeMacro(KernelTransform, Transform);

  itkNewMacro(Self);

  itkStaticConstMacro(SpaceDimension, unsigned int, NDimensions);

  typedef typename Superclass::ScalarType ScalarType;

  typedef typename Superclass::ParametersType ParametersType;

  typedef typename Superclass::JacobianType JacobianType;

  typedef typename Superclass::InputPointType  InputPointType;
  typedef typename Superclass::OutputPointType OutputPointType;

  typedef typename Superclass::InputVectorType  InputVectorType;
  typedef typename Superclass::OutputVectorType OutputVectorType;

  typedef typename Superclass::InputCovariantVectorType  InputCovariantVectorType;
  typedef typename Superclass::OutputCovariantVectorType OutputCovariantVectorType;

  typedef typename Superclass::InputVnlVectorType  InputVnlVectorType;
  typedef typename Superclass::OutputVnlVectorType OutputVnlVectorType;

  typedef DefaultStaticMeshTraits<TScalarType, NDimensions, NDimensions, TScalarType, TScalarType> PointSetTraitsType;
  typedef PointSet<InputPointType, NDimensions, PointSetTraitsType>                                PointSetType;

  typedef typename PointSetType::Pointer                      PointSetPointer;
  typedef typename PointSetType::PointsContainer              PointsContainer;
  typedef typename PointSetType::PointsContainerIterator      PointsIterator;
  typedef typename PointSetType::PointsContainerConstIterator PointsConstIterator;
  typedef typename PointSetType::PointIdentifier              PointIdentifier;

  typedef itk::VectorContainer<SizeValueType, InputVectorType> VectorSetType;
  typedef typename VectorSetType::Pointer                      VectorSetPointer;

  itkGetObjectMacro(SourceLandmarks, PointSetType);

  virtual void SetSourceLandmarks(PointSetType *);

  itkGetObjectMacro(TargetLandmarks, PointSetType);

  virtual void SetTargetLandmarks(PointSetType *);

  itkGetObjectMacro(Displacements, VectorSetType);

  void ComputeWMatrix(void);

  virtual OutputPointType TransformPoint(const InputPointType & thisPoint) const;

  using Superclass::TransformVector;
  virtual OutputVectorType TransformVector(const InputVectorType &) const                       \
  {                                                                                             \
    itkExceptionMacro(                                                                          \
      << "TransformVector(const InputVectorType &) is not implemented for KernelTransform");    \
  }

  virtual OutputVnlVectorType TransformVector(const InputVnlVectorType &) const                 \
  {                                                                                             \
    itkExceptionMacro(                                                                          \
      << "TransformVector(const InputVnlVectorType &) is not implemented for KernelTransform"); \
  }

  using Superclass::TransformCovariantVector;
  virtual OutputCovariantVectorType TransformCovariantVector(const InputCovariantVectorType &) const           \
  {                                                                                                            \
    itkExceptionMacro(                                                                                         \
      << "TransformCovariantVector(const InputCovariantVectorType &) is not implemented for KernelTransform"); \
  }

  typedef vnl_matrix_fixed<TScalarType, NDimensions, NDimensions> IMatrixType;

  virtual void ComputeJacobianWithRespectToParameters( const InputPointType  & p, JacobianType & jacobian) const;

  virtual void ComputeJacobianWithRespectToPosition(const InputPointType &,
                                                    JacobianType &) const           \
  {                                                                                 \
    itkExceptionMacro( "ComputeJacobianWithRespectToPosition not yet implemented "  \
                       "for " << this->GetNameOfClass() );                          \
  }

  virtual void SetParameters(const ParametersType &);

  virtual void SetFixedParameters(const ParametersType &);

  virtual void UpdateParameters(void) const;

  virtual const ParametersType & GetParameters(void) const;

  virtual const ParametersType & GetFixedParameters(void) const;

  virtual bool IsLinear() const
  {
    return false;
  }

  itkSetClampMacro( Stiffness, double, 0.0, NumericTraits<double>::max() );
  itkGetConstMacro(Stiffness, double);
protected:
  KernelTransform();
  virtual ~KernelTransform();
  void PrintSelf(std::ostream & os, Indent indent) const;

public:
  typedef vnl_matrix_fixed<TScalarType, NDimensions, NDimensions> GMatrixType;

  typedef vnl_matrix<TScalarType> LMatrixType;

  typedef vnl_matrix<TScalarType> KMatrixType;

  typedef vnl_matrix<TScalarType> PMatrixType;

  typedef vnl_matrix<TScalarType> YMatrixType;

  typedef vnl_matrix<TScalarType> WMatrixType;

  typedef vnl_matrix<TScalarType> DMatrixType;

  typedef vnl_matrix_fixed<TScalarType, NDimensions, NDimensions> AMatrixType;

  typedef vnl_vector_fixed<TScalarType, NDimensions> BMatrixType;

  typedef vnl_matrix_fixed<TScalarType, 1, NDimensions> RowMatrixType;

  typedef vnl_matrix_fixed<TScalarType, NDimensions, 1> ColumnMatrixType;
protected:

  virtual void ComputeG(const InputVectorType & landmarkVector, GMatrixType & gmatrix) const;

  virtual const GMatrixType & ComputeReflexiveG(PointsIterator) const;

  virtual void ComputeDeformationContribution(const InputPointType & inputPoint, OutputPointType & result) const;

  void ComputeK();

  void ComputeL();

  void ComputeP();

  void ComputeY();

  void ComputeD();

  void ReorganizeW(void);

  double m_Stiffness;

  VectorSetPointer m_Displacements;

  LMatrixType m_LMatrix;

  KMatrixType m_KMatrix;

  PMatrixType m_PMatrix;

  YMatrixType m_YMatrix;

  WMatrixType m_WMatrix;

  DMatrixType m_DMatrix;

  AMatrixType m_AMatrix;

  BMatrixType m_BVector;

  mutable GMatrixType m_GMatrix;

  bool m_WMatrixComputed;

  IMatrixType m_I;

  PointSetPointer m_SourceLandmarks;

  PointSetPointer m_TargetLandmarks;
private:

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

};
} // end namespace itk

// Define instantiation macro for this template.
#define ITK_TEMPLATE_KernelTransform(_, EXPORT, TypeX, TypeY)     \
  namespace itk                                                   \
  {                                                               \
  _( 2 ( class EXPORT KernelTransform<ITK_TEMPLATE_2 TypeX> ) ) \
  namespace Templates                                             \
  {                                                               \
  typedef KernelTransform<ITK_TEMPLATE_2 TypeX>                 \
  KernelTransform##TypeY;                                       \
  }                                                               \
  }

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

#if ITK_TEMPLATE_TXX
#include "itkKernelTransform.hxx"
#endif

#endif // __itkKernelTransform_h