ITK  4.1.0
Insight Segmentation and Registration Toolkit
itkMatrixOffsetTransformBase.h
Go to the documentation of this file.
00001 /*=========================================================================
00002  *
00003  *  Copyright Insight Software Consortium
00004  *
00005  *  Licensed under the Apache License, Version 2.0 (the "License");
00006  *  you may not use this file except in compliance with the License.
00007  *  You may obtain a copy of the License at
00008  *
00009  *         http://www.apache.org/licenses/LICENSE-2.0.txt
00010  *
00011  *  Unless required by applicable law or agreed to in writing, software
00012  *  distributed under the License is distributed on an "AS IS" BASIS,
00013  *  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
00014  *  See the License for the specific language governing permissions and
00015  *  limitations under the License.
00016  *
00017  *=========================================================================*/
00018 #ifndef __itkMatrixOffsetTransformBase_h
00019 #define __itkMatrixOffsetTransformBase_h
00020 
00021 #include <iostream>
00022 
00023 #include "itkMatrix.h"
00024 #include "itkTransform.h"
00025 #include "itkMacro.h"
00026 
00027 namespace itk
00028 {
00070 template <
00071   class TScalarType = double,         // Data type for scalars
00072   unsigned int NInputDimensions = 3,  // Number of dimensions in the input space
00073   unsigned int NOutputDimensions = 3>
00074 // Number of dimensions in the output space
00075 class MatrixOffsetTransformBase :
00076   public Transform<TScalarType, NInputDimensions, NOutputDimensions>
00077 {
00078 public:
00080   typedef MatrixOffsetTransformBase Self;
00081   typedef Transform<TScalarType,
00082                     NInputDimensions,
00083                     NOutputDimensions>        Superclass;
00084 
00085   typedef SmartPointer<Self>       Pointer;
00086   typedef SmartPointer<const Self> ConstPointer;
00087 
00089   itkTypeMacro(MatrixOffsetTransformBase, Transform);
00090 
00092   itkNewMacro(Self);
00093 
00095   itkStaticConstMacro(InputSpaceDimension, unsigned int, NInputDimensions);
00096   itkStaticConstMacro(OutputSpaceDimension, unsigned int, NOutputDimensions);
00097   itkStaticConstMacro( ParametersDimension, unsigned int,
00098                        NOutputDimensions * ( NInputDimensions + 1 ) );
00100 
00102   typedef typename Superclass::ParametersType      ParametersType;
00103   typedef typename Superclass::ParametersValueType ParametersValueType;
00104 
00106   typedef typename Superclass::JacobianType JacobianType;
00107 
00109   typedef typename Superclass::ScalarType ScalarType;
00110 
00112   typedef Vector<TScalarType,
00113                  itkGetStaticConstMacro(InputSpaceDimension)>  InputVectorType;
00114   typedef Vector<TScalarType,
00115                  itkGetStaticConstMacro(OutputSpaceDimension)> OutputVectorType;
00116   typedef typename OutputVectorType::ValueType OutputVectorValueType;
00118 
00120   typedef CovariantVector<TScalarType,
00121                           itkGetStaticConstMacro(InputSpaceDimension)>
00122   InputCovariantVectorType;
00123   typedef CovariantVector<TScalarType,
00124                           itkGetStaticConstMacro(OutputSpaceDimension)>
00125   OutputCovariantVectorType;
00127 
00128   typedef typename Superclass::InputVectorPixelType  InputVectorPixelType;
00129   typedef typename Superclass::OutputVectorPixelType OutputVectorPixelType;
00130 
00132   typedef typename Superclass::InputDiffusionTensor3DType
00133   InputDiffusionTensor3DType;
00134   typedef typename Superclass::OutputDiffusionTensor3DType
00135   OutputDiffusionTensor3DType;
00136 
00138   typedef typename Superclass::InputSymmetricSecondRankTensorType
00139   InputSymmetricSecondRankTensorType;
00140   typedef typename Superclass::OutputSymmetricSecondRankTensorType
00141   OutputSymmetricSecondRankTensorType;
00142 
00143   typedef CovariantVector<TScalarType, InputDiffusionTensor3DType::Dimension>
00144   InputTensorEigenVectorType;
00145 
00147   typedef vnl_vector_fixed<TScalarType,
00148                            itkGetStaticConstMacro(InputSpaceDimension)>
00149   InputVnlVectorType;
00150   typedef vnl_vector_fixed<TScalarType,
00151                            itkGetStaticConstMacro(OutputSpaceDimension)>
00152   OutputVnlVectorType;
00154 
00156   typedef Point<TScalarType,
00157                 itkGetStaticConstMacro(InputSpaceDimension)>
00158   InputPointType;
00159   typedef typename InputPointType::ValueType InputPointValueType;
00160   typedef Point<TScalarType,
00161                 itkGetStaticConstMacro(OutputSpaceDimension)>
00162   OutputPointType;
00163   typedef typename OutputPointType::ValueType OutputPointValueType;
00165 
00167   typedef Matrix<TScalarType, itkGetStaticConstMacro(OutputSpaceDimension),
00168                  itkGetStaticConstMacro(InputSpaceDimension)>
00169   MatrixType;
00170   typedef typename MatrixType::ValueType MatrixValueType;
00172 
00174   typedef Matrix<TScalarType, itkGetStaticConstMacro(InputSpaceDimension),
00175                  itkGetStaticConstMacro(OutputSpaceDimension)>
00176   InverseMatrixType;
00177 
00178   typedef InputPointType CenterType;
00179 
00180   typedef OutputVectorType               OffsetType;
00181   typedef typename OffsetType::ValueType OffsetValueType;
00182 
00183   typedef OutputVectorType TranslationType;
00184 
00185   typedef typename TranslationType::ValueType TranslationValueType;
00186 
00189   typedef typename Superclass::InverseTransformBaseType InverseTransformBaseType;
00190   typedef typename InverseTransformBaseType::Pointer    InverseTransformBasePointer;
00191 
00195   virtual void SetIdentity(void);
00196 
00208   virtual void SetMatrix(const MatrixType & matrix)
00209   {
00210     m_Matrix = matrix; this->ComputeOffset();
00211     this->ComputeMatrixParameters();
00212     m_MatrixMTime.Modified(); this->Modified(); return;
00213   }
00215 
00223   virtual const MatrixType & GetMatrix() const
00224   {
00225     return m_Matrix;
00226   }
00227 
00236   void SetOffset(const OutputVectorType & offset)
00237   {
00238     m_Offset = offset; this->ComputeTranslation();
00239     this->Modified(); return;
00240   }
00242 
00248   const OutputVectorType & GetOffset(void) const
00249   {
00250     return m_Offset;
00251   }
00252 
00275   void SetCenter(const InputPointType & center)
00276   {
00277     m_Center = center; this->ComputeOffset();
00278     this->Modified(); return;
00279   }
00281 
00288   const InputPointType & GetCenter() const
00289   {
00290     return m_Center;
00291   }
00292 
00299   void SetTranslation(const OutputVectorType & translation)
00300   {
00301     m_Translation = translation; this->ComputeOffset();
00302     this->Modified(); return;
00303   }
00305 
00312   const OutputVectorType & GetTranslation(void) const
00313   {
00314     return m_Translation;
00315   }
00316 
00321   void SetParameters(const ParametersType & parameters);
00322 
00324   const ParametersType & GetParameters(void) const;
00325 
00327   virtual void SetFixedParameters(const ParametersType &);
00328 
00330   virtual const ParametersType & GetFixedParameters(void) const;
00331 
00343   void Compose(const Self *other, bool pre = 0);
00344 
00353   OutputPointType       TransformPoint(const InputPointType & point) const;
00354 
00355   using Superclass::TransformVector;
00356 
00357   OutputVectorType      TransformVector(const InputVectorType & vector) const;
00358 
00359   OutputVnlVectorType   TransformVector(const InputVnlVectorType & vector) const;
00360 
00361   OutputVectorPixelType TransformVector(const InputVectorPixelType & vector) const;
00362 
00363   using Superclass::TransformCovariantVector;
00364 
00365   OutputCovariantVectorType TransformCovariantVector(const InputCovariantVectorType & vector) const;
00366 
00367   OutputVectorPixelType TransformCovariantVector(const InputVectorPixelType & vector) const;
00368 
00369   using Superclass::TransformDiffusionTensor3D;
00370 
00371   OutputDiffusionTensor3DType TransformDiffusionTensor3D(const InputDiffusionTensor3DType & tensor) const;
00372 
00373   OutputVectorPixelType TransformDiffusionTensor3D(const InputVectorPixelType & tensor ) const;
00374 
00375   using Superclass::TransformSymmetricSecondRankTensor;
00376   OutputSymmetricSecondRankTensorType TransformSymmetricSecondRankTensor( const InputSymmetricSecondRankTensorType & tensor ) const;
00377 
00378   OutputVectorPixelType TransformSymmetricSecondRankTensor( const InputVectorPixelType & tensor ) const;
00379 
00389   virtual void ComputeJacobianWithRespectToParameters(const InputPointType  & x, JacobianType & j) const;
00390 
00394   virtual void ComputeJacobianWithRespectToPosition(const InputPointType  & x, JacobianType & jac) const;
00395 
00399   virtual void ComputeInverseJacobianWithRespectToPosition(const InputPointType  & x, JacobianType & jac) const;
00400 
00419   bool GetInverse(Self *inverse) const;
00421 
00423   virtual InverseTransformBasePointer GetInverseTransform() const;
00424 
00428   const InverseMatrixType & GetInverseMatrix(void) const;
00429 
00435   virtual bool IsLinear() const
00436   {
00437     return true;
00438   }
00439 protected:
00440 
00448   MatrixOffsetTransformBase(const MatrixType & matrix, const OutputVectorType & offset);
00449   MatrixOffsetTransformBase(unsigned int paramDims);
00450   MatrixOffsetTransformBase();
00452 
00454   virtual ~MatrixOffsetTransformBase();
00455 
00457   void PrintSelf(std::ostream & s, Indent indent) const;
00458 
00459   const InverseMatrixType & GetVarInverseMatrix(void) const
00460   {
00461     return m_InverseMatrix;
00462   }
00463   void SetVarInverseMatrix(const InverseMatrixType & matrix) const
00464   {
00465     m_InverseMatrix = matrix; m_InverseMatrixMTime.Modified();
00466   }
00467   bool InverseMatrixIsOld(void) const
00468   {
00469     if( m_MatrixMTime != m_InverseMatrixMTime )
00470       {
00471       return true;
00472       }
00473     else
00474       {
00475       return false;
00476       }
00477   }
00478 
00479   virtual void ComputeMatrixParameters(void);
00480 
00481   virtual void ComputeMatrix(void);
00482 
00483   void SetVarMatrix(const MatrixType & matrix)
00484   {
00485     m_Matrix = matrix; m_MatrixMTime.Modified();
00486   }
00487 
00488   virtual void ComputeTranslation(void);
00489 
00490   void SetVarTranslation(const OutputVectorType & translation)
00491   {
00492     m_Translation = translation;
00493   }
00494 
00495   virtual void ComputeOffset(void);
00496 
00497   void SetVarOffset(const OutputVectorType & offset)
00498   {
00499     m_Offset = offset;
00500   }
00501 
00502   void SetVarCenter(const InputPointType & center)
00503   {
00504     m_Center = center;
00505   }
00506 private:
00507 
00508   MatrixOffsetTransformBase(const Self & other);
00509   const Self & operator=(const Self &);
00510 
00511   MatrixType                m_Matrix;           // Matrix of the transformation
00512   OutputVectorType          m_Offset;           // Offset of the transformation
00513   mutable InverseMatrixType m_InverseMatrix;    // Inverse of the matrix
00514   mutable bool              m_Singular;         // Is m_Inverse singular?
00515 
00516   InputPointType   m_Center;
00517   OutputVectorType m_Translation;
00518 
00520   TimeStamp         m_MatrixMTime;
00521   mutable TimeStamp m_InverseMatrixMTime;
00522 }; // class MatrixOffsetTransformBase
00523 }  // namespace itk
00524 
00525 // Define instantiation macro for this template.
00526 #define ITK_TEMPLATE_MatrixOffsetTransformBase(_, EXPORT, TypeX, TypeY)                         \
00527   namespace itk                                                                                 \
00528   {                                                                                             \
00529   _( 3 ( class EXPORT MatrixOffsetTransformBase<ITK_TEMPLATE_3 TypeX> ) )                     \
00530   namespace Templates                                                                           \
00531   {                                                                                             \
00532   typedef MatrixOffsetTransformBase<ITK_TEMPLATE_3 TypeX> MatrixOffsetTransformBase##TypeY; \
00533   }                                                                                             \
00534   }
00535 
00536 #if ITK_TEMPLATE_EXPLICIT
00537 // template < class TScalarType, unsigned int NInputDimensions, unsigned int
00538 // NOutputDimensions>
00539 //   const unsigned int itk::MatrixOffsetTransformBase< TScalarType,
00540 // NInputDimensions, NOutputDimensions >::InputSpaceDimension;
00541 // template < class TScalarType, unsigned int NInputDimensions, unsigned int
00542 // NOutputDimensions>
00543 //   const unsigned int itk::MatrixOffsetTransformBase< TScalarType,
00544 // NInputDimensions, NOutputDimensions >::OutputSpaceDimension;
00545 // template < class TScalarType, unsigned int NInputDimensions, unsigned int
00546 // NOutputDimensions>
00547 //   const unsigned int itk::MatrixOffsetTransformBase< TScalarType,
00548 // NInputDimensions, NOutputDimensions >::ParametersDimension;
00549 #include "Templates/itkMatrixOffsetTransformBase+-.h"
00550 #endif
00551 
00552 #if ITK_TEMPLATE_TXX
00553 #include "itkMatrixOffsetTransformBase.hxx"
00554 #endif
00555 
00556 #endif /* __itkMatrixOffsetTransformBase_h */
00557