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 
109 
111  typedef typename Superclass
114  typedef typename Superclass
117  typedef typename Superclass
120 
124  typedef DefaultStaticMeshTraits< TScalarType,
125  NDimensions, NDimensions, TScalarType, TScalarType> PointSetTraitsType;
126  typedef PointSet<InputPointType, NDimensions,
128  typedef typename PointSetType::Pointer PointSetPointer;
129  typedef typename PointSetType::PointsContainer PointsContainer;
130  typedef typename PointSetType::PointsContainerIterator PointsIterator;
131  typedef typename PointSetType::PointsContainerConstIterator PointsConstIterator;
132 
134  typedef VectorContainer<unsigned long,InputVectorType> VectorSetType;
135  typedef typename VectorSetType::Pointer VectorSetPointer;
136 
138  typedef vnl_matrix_fixed<TScalarType, NDimensions, NDimensions> IMatrixType;
139 
141  virtual unsigned int GetNumberOfParameters(void) const
142  {
143  return ( this->m_SourceLandmarks->GetNumberOfPoints() * SpaceDimension );
144  }
145 
147  itkGetObjectMacro( SourceLandmarks, PointSetType );
148 
150  virtual void SetSourceLandmarks( PointSetType * );
151 
153  itkGetObjectMacro( TargetLandmarks, PointSetType );
154 
156  virtual void SetTargetLandmarks( PointSetType * );
157 
161  itkGetObjectMacro( Displacements, VectorSetType );
162 
164  void ComputeWMatrix( void );
165 
167  void ComputeLInverse( void );
168 
170  virtual OutputPointType TransformPoint( const InputPointType & thisPoint ) const;
171 
173  virtual const JacobianType & GetJacobian( const InputPointType & point ) const;
174 
176  virtual void GetJacobian(
177  const InputPointType &,
178  JacobianType &,
179  NonZeroJacobianIndicesType & ) const;
180 
182  virtual void SetIdentity( void );
183 
189  virtual void SetParameters( const ParametersType & );
190 
196  virtual void SetFixedParameters( const ParametersType & );
197 
199  virtual void UpdateParameters( void );
200 
202  virtual const ParametersType & GetParameters( void ) const;
203 
205  virtual const ParametersType & GetFixedParameters( void ) const;
206 
217  virtual void SetStiffness( double stiffness )
218  {
219  this->m_Stiffness = stiffness > 0 ? stiffness : 0.0;
220  this->m_LMatrixComputed = false;
221  this->m_LInverseComputed = false;
222  this->m_WMatrixComputed = false;
223  }
224  itkGetMacro( Stiffness, double );
225 
232  virtual void SetAlpha( TScalarType itkNotUsed( Alpha ) ) {};
233  virtual TScalarType GetAlpha( void ) const { return -1.0; }
234 
241  itkSetMacro( PoissonRatio, TScalarType );
242  virtual const TScalarType GetPoissonRatio( void ) const
243  {
244  return this->m_PoissonRatio;
245  };
246 
248  itkSetMacro( MatrixInversionMethod, std::string );
249  itkGetConstReferenceMacro( MatrixInversionMethod, std::string );
250 
251 protected:
253  virtual ~KernelTransform2();
254  void PrintSelf( std::ostream& os, Indent indent ) const;
255 
256 public:
258  typedef vnl_matrix_fixed<TScalarType, NDimensions, NDimensions> GMatrixType;
259 
261  typedef vnl_matrix<TScalarType> LMatrixType;
262 
264  typedef vnl_matrix<TScalarType> KMatrixType;
265 
267  typedef vnl_matrix<TScalarType> PMatrixType;
268 
270  typedef vnl_matrix<TScalarType> YMatrixType;
271 
273  typedef vnl_matrix<TScalarType> WMatrixType;
274 
276  typedef vnl_matrix<TScalarType> DMatrixType;
277 
279  typedef vnl_matrix_fixed<TScalarType,NDimensions,NDimensions> AMatrixType;
280 
282  typedef vnl_vector_fixed<TScalarType,NDimensions> BMatrixType;
283 
285  typedef vnl_matrix_fixed<TScalarType, 1, NDimensions> RowMatrixType;
286 
288  typedef vnl_matrix_fixed<TScalarType, NDimensions, 1> ColumnMatrixType;
289 
292 
295 
296 protected:
297 
305  virtual void ComputeG( const InputVectorType & landmarkVector,
306  GMatrixType & GMatrix ) const;
307 
315  virtual void ComputeReflexiveG( PointsIterator, GMatrixType & GMatrix ) const;
316 
320  virtual void ComputeDeformationContribution(
321  const InputPointType & inputPoint,
322  OutputPointType & result ) const;
323 
325  void ComputeK( void );
326 
328  void ComputeL( void );
329 
331  void ComputeP( void );
332 
334  void ComputeY( void );
335 
337  void ComputeD( void );
338 
343  void ReorganizeW( void );
344 
346  double m_Stiffness;
347 
352 
355 
358 
361 
364 
367 
370 
377 
380 
383 
389  //GMatrixType m_GMatrix;
390 
399 
406  typedef vnl_svd< ScalarType > SVDDecompositionType;
407  typedef vnl_qr< ScalarType > QRDecompositionType;
408 
411 
414 
417 
420 
425 
426 private:
427  KernelTransform2(const Self&); //purposely not implemented
428  void operator=(const Self&); //purposely not implemented
429 
430  TScalarType m_PoissonRatio;
431 
434 
435 };
436 
437 } // end namespace itk
438 
439 #ifndef ITK_MANUAL_INSTANTIATION
440 #include "itkKernelTransform2.txx"
441 #endif
442 
443 #endif // __itkKernelTransform2_h


Generated on 27-06-2013 for elastix by doxygen 1.8.3.1 elastix logo