go home Home | Main Page | Modules | Namespace List | Class Hierarchy | Alphabetical List | Data Structures | File List | Namespace Members | Data Fields | Globals | Related Pages
itkKernelTransform2.h
Go to the documentation of this file.
1 /*=========================================================================
2 
3 Program: Insight Segmentation & Registration Toolkit
4 Module: $RCSfile: itkKernelTransform2.h,v $
5 Language: C++
6 Date: $Date: 2006-11-28 14:22:18 $
7 Version: $Revision: 1.1 $
8 
9 Copyright (c) Insight Software Consortium. All rights reserved.
10 See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.
11 
12 This software is distributed WITHOUT ANY WARRANTY; without even
13 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
14 PURPOSE. See the above copyright notices for more information.
15 
16 =========================================================================*/
17 #ifndef __itkKernelTransform2_h
18 #define __itkKernelTransform2_h
19 
20 #include "itkAdvancedTransform.h"
21 #include "itkPoint.h"
22 #include "itkVector.h"
23 #include "itkMatrix.h"
24 #include "itkPointSet.h"
25 #include <deque>
26 #include <math.h>
27 #include "vnl/vnl_matrix_fixed.h"
28 #include "vnl/vnl_matrix.h"
29 #include "vnl/vnl_vector.h"
30 #include "vnl/vnl_vector_fixed.h"
31 #include "vnl/vnl_sample.h"
32 #include "vnl/algo/vnl_svd.h"
33 #include "vnl/algo/vnl_qr.h"
34 
35 
36 namespace itk
37 {
38 
78 template <class TScalarType, // probably only float and double make sense here
79  unsigned int NDimensions> // Number of dimensions
81  : public AdvancedTransform<TScalarType, NDimensions,NDimensions>
82 {
83 public:
84 
87  typedef AdvancedTransform<
88  TScalarType, NDimensions, NDimensions > Superclass;
89  typedef SmartPointer<Self> Pointer;
90  typedef SmartPointer<const Self> ConstPointer;
91 
93  itkTypeMacro( KernelTransform2, AdvancedTransform );
94 
96  itkNewMacro( Self );
97 
99  itkStaticConstMacro( SpaceDimension, unsigned int, NDimensions );
100 
114 
116  typedef typename Superclass
119  typedef typename Superclass
122  typedef typename Superclass
125 
129  typedef DefaultStaticMeshTraits< TScalarType,
130  NDimensions, NDimensions, TScalarType, TScalarType> PointSetTraitsType;
131  typedef PointSet<InputPointType, NDimensions,
133  typedef typename PointSetType::Pointer PointSetPointer;
134  typedef typename PointSetType::PointsContainer PointsContainer;
135  typedef typename PointSetType::PointsContainerIterator PointsIterator;
136  typedef typename PointSetType::PointsContainerConstIterator PointsConstIterator;
137 
139  typedef VectorContainer<unsigned long,InputVectorType> VectorSetType;
140  typedef typename VectorSetType::Pointer VectorSetPointer;
141 
143  typedef vnl_matrix_fixed<TScalarType, NDimensions, NDimensions> IMatrixType;
144 
147  {
148  return ( this->m_SourceLandmarks->GetNumberOfPoints() * SpaceDimension );
149  }
150 
152  itkGetObjectMacro( SourceLandmarks, PointSetType );
153 
155  virtual void SetSourceLandmarks( PointSetType * );
156 
158  itkGetObjectMacro( TargetLandmarks, PointSetType );
159 
161  virtual void SetTargetLandmarks( PointSetType * );
162 
166  itkGetObjectMacro( Displacements, VectorSetType );
167 
169  void ComputeWMatrix( void );
170 
172  void ComputeLInverse( void );
173 
175  virtual OutputPointType TransformPoint( const InputPointType & thisPoint ) const;
176 
179  {
180  itkExceptionMacro(
181  << "TransformVector(const InputVectorType &) is not implemented "
182  << "for KernelTransform" );
183  }
185  {
186  itkExceptionMacro(
187  << "TransformVector(const InputVnlVectorType &) is not implemented "
188  << "for KernelTransform" );
189  }
191  {
192  itkExceptionMacro(
193  << "TransformCovariantVector(const InputCovariantVectorType &) is not implemented "
194  << "for KernelTransform" );
195  }
196 
198  virtual void GetJacobian(
199  const InputPointType &,
200  JacobianType &,
201  NonZeroJacobianIndicesType & ) const;
202 
204  virtual void SetIdentity( void );
205 
211  virtual void SetParameters( const ParametersType & );
212 
218  virtual void SetFixedParameters( const ParametersType & );
219 
221  virtual void UpdateParameters( void );
222 
224  virtual const ParametersType & GetParameters( void ) const;
225 
227  virtual const ParametersType & GetFixedParameters( void ) const;
228 
239  virtual void SetStiffness( double stiffness )
240  {
241  this->m_Stiffness = stiffness > 0 ? stiffness : 0.0;
242  this->m_LMatrixComputed = false;
243  this->m_LInverseComputed = false;
244  this->m_WMatrixComputed = false;
245  }
246  itkGetMacro( Stiffness, double );
247 
254  virtual void SetAlpha( TScalarType itkNotUsed( Alpha ) ) {};
255  virtual TScalarType GetAlpha( void ) const { return -1.0; }
256 
263  itkSetMacro( PoissonRatio, TScalarType );
264  virtual const TScalarType GetPoissonRatio( void ) const
265  {
266  return this->m_PoissonRatio;
267  };
268 
270  itkSetMacro( MatrixInversionMethod, std::string );
271  itkGetConstReferenceMacro( MatrixInversionMethod, std::string );
272 
274  virtual void GetSpatialJacobian(
275  const InputPointType & ipp, SpatialJacobianType & sj ) const
276  {
277  itkExceptionMacro( << "Not implemented for KernelTransform2" );
278  }
279  virtual void GetSpatialHessian(
280  const InputPointType & ipp, SpatialHessianType & sh ) const
281  {
282  itkExceptionMacro( << "Not implemented for KernelTransform2" );
283  }
286  NonZeroJacobianIndicesType & nonZeroJacobianIndices ) const
287  {
288  itkExceptionMacro( << "Not implemented for KernelTransform2" );
289  }
291  const InputPointType & ipp, SpatialJacobianType & sj,
293  NonZeroJacobianIndicesType & nonZeroJacobianIndices ) const
294  {
295  itkExceptionMacro( << "Not implemented for KernelTransform2" );
296  }
298  const InputPointType & ipp, JacobianOfSpatialHessianType & jsh,
299  NonZeroJacobianIndicesType & nonZeroJacobianIndices ) const
300  {
301  itkExceptionMacro( << "Not implemented for KernelTransform2" );
302  }
304  const InputPointType & ipp, SpatialHessianType & sh,
306  NonZeroJacobianIndicesType & nonZeroJacobianIndices ) const
307  {
308  itkExceptionMacro( << "Not implemented for KernelTransform2" );
309  }
310 
311 protected:
313  virtual ~KernelTransform2();
314  void PrintSelf( std::ostream& os, Indent indent ) const;
315 
316 public:
318  typedef vnl_matrix_fixed<TScalarType, NDimensions, NDimensions> GMatrixType;
319 
321  typedef vnl_matrix<TScalarType> LMatrixType;
322 
324  typedef vnl_matrix<TScalarType> KMatrixType;
325 
327  typedef vnl_matrix<TScalarType> PMatrixType;
328 
330  typedef vnl_matrix<TScalarType> YMatrixType;
331 
333  typedef vnl_matrix<TScalarType> WMatrixType;
334 
336  typedef vnl_matrix<TScalarType> DMatrixType;
337 
339  typedef vnl_matrix_fixed<TScalarType,NDimensions,NDimensions> AMatrixType;
340 
342  typedef vnl_vector_fixed<TScalarType,NDimensions> BMatrixType;
343 
345  typedef vnl_matrix_fixed<TScalarType, 1, NDimensions> RowMatrixType;
346 
348  typedef vnl_matrix_fixed<TScalarType, NDimensions, 1> ColumnMatrixType;
349 
352 
355 
356 protected:
357 
365  virtual void ComputeG( const InputVectorType & landmarkVector,
366  GMatrixType & GMatrix ) const;
367 
375  virtual void ComputeReflexiveG( PointsIterator, GMatrixType & GMatrix ) const;
376 
380  virtual void ComputeDeformationContribution(
381  const InputPointType & inputPoint,
382  OutputPointType & result ) const;
383 
385  void ComputeK( void );
386 
388  void ComputeL( void );
389 
391  void ComputeP( void );
392 
394  void ComputeY( void );
395 
397  void ComputeD( void );
398 
403  void ReorganizeW( void );
404 
406  double m_Stiffness;
407 
412 
415 
418 
421 
424 
427 
430 
437 
440 
443 
449  //GMatrixType m_GMatrix;
450 
459 
466  typedef vnl_svd< ScalarType > SVDDecompositionType;
467  typedef vnl_qr< ScalarType > QRDecompositionType;
468 
471 
474 
477 
480 
485 
486 private:
487  KernelTransform2(const Self&); //purposely not implemented
488  void operator=(const Self&); //purposely not implemented
489 
490  TScalarType m_PoissonRatio;
491 
494 
495 };
496 
497 } // end namespace itk
498 
499 #ifndef ITK_MANUAL_INSTANTIATION
500 #include "itkKernelTransform2.txx"
501 #endif
502 
503 #endif // __itkKernelTransform2_h
Superclass::OutputPointType OutputPointType
AdvancedTransform< TScalarType, NDimensions, NDimensions > Superclass
Matrix< ScalarType, OutputSpaceDimension, InputSpaceDimension > SpatialJacobianType
virtual ~KernelTransform2()
NonZeroJacobianIndicesType m_NonZeroJacobianIndices
virtual void ComputeReflexiveG(PointsIterator, GMatrixType &GMatrix) const
Superclass::OutputCovariantVectorType OutputCovariantVectorType
vnl_matrix_fixed< TScalarType, 1, NDimensions > RowMatrixType
virtual void SetAlpha(TScalarType)
virtual void GetJacobianOfSpatialJacobian(const InputPointType &ipp, SpatialJacobianType &sj, JacobianOfSpatialJacobianType &jsj, NonZeroJacobianIndicesType &nonZeroJacobianIndices) const
vnl_matrix_fixed< TScalarType, NDimensions, NDimensions > AMatrixType
itkStaticConstMacro(SpaceDimension, unsigned int, NDimensions)
vnl_matrix< TScalarType > LMatrixType
PointSetPointer m_TargetLandmarks
vnl_matrix< TScalarType > KMatrixType
virtual OutputVnlVectorType TransformVector(const InputVnlVectorType &) const
vnl_matrix< TScalarType > YMatrixType
vnl_matrix< TScalarType > WMatrixType
void ComputeWMatrix(void)
NonZeroJacobianIndicesType m_NonZeroJacobianIndicesTemp
virtual void GetJacobianOfSpatialHessian(const InputPointType &ipp, JacobianOfSpatialHessianType &jsh, NonZeroJacobianIndicesType &nonZeroJacobianIndices) const
virtual void GetJacobian(const InputPointType &, JacobianType &, NonZeroJacobianIndicesType &) const
virtual const ParametersType & GetParameters(void) const
Superclass::SpatialHessianType SpatialHessianType
vnl_matrix_fixed< TScalarType, NDimensions, NDimensions > GMatrixType
PointSetPointer m_SourceLandmarks
vnl_svd< ScalarType > SVDDecompositionType
Superclass::InputVectorType InputVectorType
PointSetType::PointsContainerConstIterator PointsConstIterator
vnl_vector_fixed< TScalarType, NDimensions > BMatrixType
virtual void ComputeG(const InputVectorType &landmarkVector, GMatrixType &GMatrix) const
PointSet< InputPointType, NDimensions, PointSetTraitsType > PointSetType
vnl_matrix< TScalarType > DMatrixType
virtual void UpdateParameters(void)
virtual void GetSpatialJacobian(const InputPointType &ipp, SpatialJacobianType &sj) const
Superclass::JacobianOfSpatialJacobianType JacobianOfSpatialJacobianType
Superclass::NonZeroJacobianIndicesType NonZeroJacobianIndicesType
virtual void SetStiffness(double stiffness)
virtual void GetJacobianOfSpatialJacobian(const InputPointType &ipp, JacobianOfSpatialJacobianType &jsj, NonZeroJacobianIndicesType &nonZeroJacobianIndices) const
void PrintSelf(std::ostream &os, Indent indent) const
virtual void SetFixedParameters(const ParametersType &)
vnl_matrix_fixed< TScalarType, NDimensions, NDimensions > IMatrixType
Transform maps points, vectors and covariant vectors from an input space to an output space...
virtual void SetTargetLandmarks(PointSetType *)
Superclass::OutputVectorType OutputVectorType
virtual void SetSourceLandmarks(PointSetType *)
DefaultStaticMeshTraits< TScalarType, NDimensions, NDimensions, TScalarType, TScalarType > PointSetTraitsType
virtual void SetParameters(const ParametersType &)
Superclass::NumberOfParametersType NumberOfParametersType
Superclass::JacobianOfSpatialHessianType JacobianOfSpatialHessianType
Superclass::SpatialJacobianType SpatialJacobianType
SmartPointer< Self > Pointer
Superclass::InputPointType InputPointType
virtual const TScalarType GetPoissonRatio(void) const
vnl_matrix< TScalarType > PMatrixType
PointSetType::PointsContainer PointsContainer
Superclass::InternalMatrixType InternalMatrixType
Superclass::OutputVnlVectorType OutputVnlVectorType
VectorSetType::Pointer VectorSetPointer
virtual NumberOfParametersType GetNumberOfParameters(void) const
Superclass::ParametersType ParametersType
PointSetType::PointsContainerIterator PointsIterator
virtual OutputCovariantVectorType TransformCovariantVector(const InputCovariantVectorType &) const
Superclass::ScalarType ScalarType
virtual void ComputeDeformationContribution(const InputPointType &inputPoint, OutputPointType &result) const
QRDecompositionType * m_LMatrixDecompositionQR
virtual void GetJacobianOfSpatialHessian(const InputPointType &ipp, SpatialHessianType &sh, JacobianOfSpatialHessianType &jsh, NonZeroJacobianIndicesType &nonZeroJacobianIndices) const
Superclass::InputVnlVectorType InputVnlVectorType
Superclass::InputCovariantVectorType InputCovariantVectorType
VectorContainer< unsigned long, InputVectorType > VectorSetType
virtual void GetSpatialHessian(const InputPointType &ipp, SpatialHessianType &sh) const
Superclass::JacobianType JacobianType
virtual TScalarType GetAlpha(void) const
void ComputeLInverse(void)
VectorSetPointer m_Displacements
SVDDecompositionType * m_LMatrixDecompositionSVD
virtual OutputVectorType TransformVector(const InputVectorType &) const
SmartPointer< const Self > ConstPointer
PointSetType::Pointer PointSetPointer
void operator=(const Self &)
vnl_qr< ScalarType > QRDecompositionType
vnl_matrix_fixed< TScalarType, NDimensions, 1 > ColumnMatrixType
FixedArray< Matrix< ScalarType, InputSpaceDimension, InputSpaceDimension >, OutputSpaceDimension > SpatialHessianType
virtual OutputPointType TransformPoint(const InputPointType &thisPoint) const
virtual void SetIdentity(void)
virtual const ParametersType & GetFixedParameters(void) const


Generated on 04-01-2014 for elastix by doxygen 1.8.5 elastix logo