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  * Copyright UMC Utrecht and contributors
4  *
5  * Licensed under the Apache License, Version 2.0 (the "License");
6  * you may not use this file except in compliance with the License.
7  * You may obtain a copy of the License at
8  *
9  * http://www.apache.org/licenses/LICENSE-2.0.txt
10  *
11  * Unless required by applicable law or agreed to in writing, software
12  * distributed under the License is distributed on an "AS IS" BASIS,
13  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14  * See the License for the specific language governing permissions and
15  * limitations under the License.
16  *
17  *=========================================================================*/
18 /*=========================================================================
19 
20 Program: Insight Segmentation & Registration Toolkit
21 Module: $RCSfile: itkKernelTransform2.h,v $
22 Language: C++
23 Date: $Date: 2006-11-28 14:22:18 $
24 Version: $Revision: 1.1 $
25 
26 Copyright (c) Insight Software Consortium. All rights reserved.
27 See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.
28 
29 This software is distributed WITHOUT ANY WARRANTY; without even
30 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
31 PURPOSE. See the above copyright notices for more information.
32 
33 =========================================================================*/
34 #ifndef __itkKernelTransform2_h
35 #define __itkKernelTransform2_h
36 
37 #include "itkAdvancedTransform.h"
38 #include "itkPoint.h"
39 #include "itkVector.h"
40 #include "itkMatrix.h"
41 #include "itkPointSet.h"
42 #include <deque>
43 #include <math.h>
44 #include "vnl/vnl_matrix_fixed.h"
45 #include "vnl/vnl_matrix.h"
46 #include "vnl/vnl_vector.h"
47 #include "vnl/vnl_vector_fixed.h"
48 #include "vnl/vnl_sample.h"
49 #include "vnl/algo/vnl_svd.h"
50 #include "vnl/algo/vnl_qr.h"
51 
52 namespace itk
53 {
54 
94 template< class TScalarType, // probably only float and double make sense here
95 unsigned int NDimensions >
96 // Number of dimensions
98  public AdvancedTransform< TScalarType, NDimensions, NDimensions >
99 {
100 public:
101 
104  typedef AdvancedTransform<
105  TScalarType, NDimensions, NDimensions > Superclass;
107  typedef SmartPointer< const Self > ConstPointer;
108 
110  itkTypeMacro( KernelTransform2, AdvancedTransform );
111 
113  itkNewMacro( Self );
114 
116  itkStaticConstMacro( SpaceDimension, unsigned int, NDimensions );
117 
131 
133  typedef typename Superclass
136  typedef typename Superclass
139  typedef typename Superclass
142 
146  typedef DefaultStaticMeshTraits< TScalarType,
147  NDimensions, NDimensions, TScalarType, TScalarType > PointSetTraitsType;
148  typedef PointSet< InputPointType, NDimensions,
149  PointSetTraitsType > PointSetType;
150  typedef typename PointSetType::Pointer PointSetPointer;
151  typedef typename PointSetType::PointsContainer PointsContainer;
152  typedef typename PointSetType::PointsContainerIterator PointsIterator;
153  typedef typename PointSetType::PointsContainerConstIterator PointsConstIterator;
154 
156  typedef VectorContainer< unsigned long, InputVectorType > VectorSetType;
157  typedef typename VectorSetType::Pointer VectorSetPointer;
158 
160  typedef vnl_matrix_fixed< TScalarType, NDimensions, NDimensions > IMatrixType;
161 
163  virtual NumberOfParametersType GetNumberOfParameters( void ) const
164  {
165  return ( this->m_SourceLandmarks->GetNumberOfPoints() * SpaceDimension );
166  }
167 
168 
170  itkGetObjectMacro( SourceLandmarks, PointSetType );
171 
173  virtual void SetSourceLandmarks( PointSetType * );
174 
176  itkGetObjectMacro( TargetLandmarks, PointSetType );
177 
179  virtual void SetTargetLandmarks( PointSetType * );
180 
184  itkGetObjectMacro( Displacements, VectorSetType );
185 
187  void ComputeWMatrix( void );
188 
190  void ComputeLInverse( void );
191 
193  virtual OutputPointType TransformPoint( const InputPointType & thisPoint ) const;
194 
196  virtual OutputVectorType TransformVector( const InputVectorType & ) const
197  {
198  itkExceptionMacro(
199  << "TransformVector(const InputVectorType &) is not implemented "
200  << "for KernelTransform" );
201  }
202 
203 
204  virtual OutputVnlVectorType TransformVector( const InputVnlVectorType & ) const
205  {
206  itkExceptionMacro(
207  << "TransformVector(const InputVnlVectorType &) is not implemented "
208  << "for KernelTransform" );
209  }
210 
211 
212  virtual OutputCovariantVectorType TransformCovariantVector( const InputCovariantVectorType & ) const
213  {
214  itkExceptionMacro(
215  << "TransformCovariantVector(const InputCovariantVectorType &) is not implemented "
216  << "for KernelTransform" );
217  }
218 
219 
221  virtual void GetJacobian(
222  const InputPointType &,
223  JacobianType &,
224  NonZeroJacobianIndicesType & ) const;
225 
227  virtual void SetIdentity( void );
228 
234  virtual void SetParameters( const ParametersType & );
235 
241  virtual void SetFixedParameters( const ParametersType & );
242 
244  virtual void UpdateParameters( void );
245 
247  virtual const ParametersType & GetParameters( void ) const;
248 
250  virtual const ParametersType & GetFixedParameters( void ) const;
251 
262  virtual void SetStiffness( double stiffness )
263  {
264  this->m_Stiffness = stiffness > 0 ? stiffness : 0.0;
265  this->m_LMatrixComputed = false;
266  this->m_LInverseComputed = false;
267  this->m_WMatrixComputed = false;
268  }
269 
270 
271  itkGetMacro( Stiffness, double );
272 
279  virtual void SetAlpha( TScalarType itkNotUsed( Alpha ) ) {}
280  virtual TScalarType GetAlpha( void ) const { return -1.0; }
281 
288  itkSetMacro( PoissonRatio, TScalarType );
289  virtual const TScalarType GetPoissonRatio( void ) const
290  {
291  return this->m_PoissonRatio;
292  }
293 
294 
296  itkSetMacro( MatrixInversionMethod, std::string );
297  itkGetConstReferenceMacro( MatrixInversionMethod, std::string );
298 
300  virtual void GetSpatialJacobian(
301  const InputPointType & ipp, SpatialJacobianType & sj ) const
302  {
303  itkExceptionMacro( << "Not implemented for KernelTransform2" );
304  }
305 
306 
307  virtual void GetSpatialHessian(
308  const InputPointType & ipp, SpatialHessianType & sh ) const
309  {
310  itkExceptionMacro( << "Not implemented for KernelTransform2" );
311  }
312 
313 
315  const InputPointType & ipp, JacobianOfSpatialJacobianType & jsj,
316  NonZeroJacobianIndicesType & nonZeroJacobianIndices ) const
317  {
318  itkExceptionMacro( << "Not implemented for KernelTransform2" );
319  }
320 
321 
323  const InputPointType & ipp, SpatialJacobianType & sj,
324  JacobianOfSpatialJacobianType & jsj,
325  NonZeroJacobianIndicesType & nonZeroJacobianIndices ) const
326  {
327  itkExceptionMacro( << "Not implemented for KernelTransform2" );
328  }
329 
330 
332  const InputPointType & ipp, JacobianOfSpatialHessianType & jsh,
333  NonZeroJacobianIndicesType & nonZeroJacobianIndices ) const
334  {
335  itkExceptionMacro( << "Not implemented for KernelTransform2" );
336  }
337 
338 
340  const InputPointType & ipp, SpatialHessianType & sh,
341  JacobianOfSpatialHessianType & jsh,
342  NonZeroJacobianIndicesType & nonZeroJacobianIndices ) const
343  {
344  itkExceptionMacro( << "Not implemented for KernelTransform2" );
345  }
346 
347 
348 protected:
349 
351  virtual ~KernelTransform2();
352  void PrintSelf( std::ostream & os, Indent indent ) const;
353 
354 public:
355 
357  typedef vnl_matrix_fixed< TScalarType, NDimensions, NDimensions > GMatrixType;
358 
360  typedef vnl_matrix< TScalarType > LMatrixType;
361 
363  typedef vnl_matrix< TScalarType > KMatrixType;
364 
366  typedef vnl_matrix< TScalarType > PMatrixType;
367 
369  typedef vnl_matrix< TScalarType > YMatrixType;
370 
372  typedef vnl_matrix< TScalarType > WMatrixType;
373 
375  typedef vnl_matrix< TScalarType > DMatrixType;
376 
378  typedef vnl_matrix_fixed< TScalarType, NDimensions, NDimensions > AMatrixType;
379 
381  typedef vnl_vector_fixed< TScalarType, NDimensions > BMatrixType;
382 
384  typedef vnl_matrix_fixed< TScalarType, 1, NDimensions > RowMatrixType;
385 
387  typedef vnl_matrix_fixed< TScalarType, NDimensions, 1 > ColumnMatrixType;
388 
390  PointSetPointer m_SourceLandmarks;
391 
393  PointSetPointer m_TargetLandmarks;
394 
395 protected:
396 
404  virtual void ComputeG( const InputVectorType & landmarkVector,
405  GMatrixType & GMatrix ) const;
406 
414  virtual void ComputeReflexiveG( PointsIterator, GMatrixType & GMatrix ) const;
415 
419  virtual void ComputeDeformationContribution(
420  const InputPointType & inputPoint,
421  OutputPointType & result ) const;
422 
424  void ComputeK( void );
425 
427  void ComputeL( void );
428 
430  void ComputeP( void );
431 
433  void ComputeY( void );
434 
436  void ComputeD( void );
437 
442  void ReorganizeW( void );
443 
445  double m_Stiffness;
446 
450  VectorSetPointer m_Displacements;
451 
453  LMatrixType m_LMatrix;
454 
456  LMatrixType m_LMatrixInverse;
457 
459  KMatrixType m_KMatrix;
460 
462  PMatrixType m_PMatrix;
463 
465  YMatrixType m_YMatrix;
466 
468  WMatrixType m_WMatrix;
469 
475  DMatrixType m_DMatrix;
476 
478  AMatrixType m_AMatrix;
479 
481  BMatrixType m_BVector;
482 
488  //GMatrixType m_GMatrix;
489 
498 
505  typedef vnl_svd< ScalarType > SVDDecompositionType;
506  typedef vnl_qr< ScalarType > QRDecompositionType;
507 
508  SVDDecompositionType * m_LMatrixDecompositionSVD;
509  QRDecompositionType * m_LMatrixDecompositionQR;
510 
512  IMatrixType m_I;
513 
515  NonZeroJacobianIndicesType m_NonZeroJacobianIndices;
516 
518  mutable NonZeroJacobianIndicesType m_NonZeroJacobianIndicesTemp;
519 
524 
525 private:
526 
527  KernelTransform2( const Self & ); // purposely not implemented
528  void operator=( const Self & ); // purposely not implemented
529 
530  TScalarType m_PoissonRatio;
531 
534 
535 };
536 
537 } // end namespace itk
538 
539 #ifndef ITK_MANUAL_INSTANTIATION
540 #include "itkKernelTransform2.hxx"
541 #endif
542 
543 #endif // __itkKernelTransform2_h
Superclass::OutputPointType OutputPointType
AdvancedTransform< TScalarType, NDimensions, NDimensions > Superclass
Matrix< ScalarType, OutputSpaceDimension, InputSpaceDimension > SpatialJacobianType
vnl_matrix_fixed< TScalarType, 1, NDimensions > RowMatrixType
virtual ~KernelTransform2()
NonZeroJacobianIndicesType m_NonZeroJacobianIndices
virtual void ComputeReflexiveG(PointsIterator, GMatrixType &GMatrix) const
Superclass::OutputCovariantVectorType OutputCovariantVectorType
virtual void SetAlpha(TScalarType)
virtual void GetJacobianOfSpatialJacobian(const InputPointType &ipp, SpatialJacobianType &sj, JacobianOfSpatialJacobianType &jsj, NonZeroJacobianIndicesType &nonZeroJacobianIndices) const
itkStaticConstMacro(SpaceDimension, unsigned int, NDimensions)
PointSetPointer m_TargetLandmarks
vnl_matrix< TScalarType > LMatrixType
virtual OutputVnlVectorType TransformVector(const InputVnlVectorType &) const
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
PointSetPointer m_SourceLandmarks
vnl_svd< ScalarType > SVDDecompositionType
Superclass::InputVectorType InputVectorType
PointSetType::PointsContainerConstIterator PointsConstIterator
virtual void ComputeG(const InputVectorType &landmarkVector, GMatrixType &GMatrix) const
vnl_matrix< TScalarType > WMatrixType
vnl_matrix_fixed< TScalarType, NDimensions, NDimensions > GMatrixType
virtual void UpdateParameters(void)
virtual void GetSpatialJacobian(const InputPointType &ipp, SpatialJacobianType &sj) const
Superclass::JacobianOfSpatialJacobianType JacobianOfSpatialJacobianType
Superclass::NonZeroJacobianIndicesType NonZeroJacobianIndicesType
DefaultStaticMeshTraits< TScalarType, NDimensions, NDimensions, TScalarType, TScalarType > PointSetTraitsType
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 &)
Transform maps points, vectors and covariant vectors from an input space to an output space...
virtual void SetTargetLandmarks(PointSetType *)
vnl_matrix< TScalarType > KMatrixType
Superclass::OutputVectorType OutputVectorType
vnl_matrix_fixed< TScalarType, NDimensions, 1 > ColumnMatrixType
virtual void SetSourceLandmarks(PointSetType *)
virtual void SetParameters(const ParametersType &)
Superclass::NumberOfParametersType NumberOfParametersType
Superclass::JacobianOfSpatialHessianType JacobianOfSpatialHessianType
Superclass::SpatialJacobianType SpatialJacobianType
vnl_matrix< TScalarType > PMatrixType
vnl_matrix< TScalarType > DMatrixType
Superclass::InputPointType InputPointType
virtual const TScalarType GetPoissonRatio(void) const
PointSetType::PointsContainer PointsContainer
vnl_matrix_fixed< TScalarType, NDimensions, NDimensions > AMatrixType
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
SmartPointer< Self > Pointer
Superclass::ScalarType ScalarType
vnl_matrix< TScalarType > YMatrixType
virtual void ComputeDeformationContribution(const InputPointType &inputPoint, OutputPointType &result) const
SmartPointer< const Self > ConstPointer
QRDecompositionType * m_LMatrixDecompositionQR
virtual void GetJacobianOfSpatialHessian(const InputPointType &ipp, SpatialHessianType &sh, JacobianOfSpatialHessianType &jsh, NonZeroJacobianIndicesType &nonZeroJacobianIndices) const
Superclass::InputVnlVectorType InputVnlVectorType
Superclass::InputCovariantVectorType InputCovariantVectorType
virtual void GetSpatialHessian(const InputPointType &ipp, SpatialHessianType &sh) const
Superclass::JacobianType JacobianType
virtual TScalarType GetAlpha(void) const
vnl_matrix_fixed< TScalarType, NDimensions, NDimensions > IMatrixType
void ComputeLInverse(void)
VectorSetPointer m_Displacements
SVDDecompositionType * m_LMatrixDecompositionSVD
virtual OutputVectorType TransformVector(const InputVectorType &) const
VectorContainer< unsigned long, InputVectorType > VectorSetType
PointSetType::Pointer PointSetPointer
void operator=(const Self &)
vnl_qr< ScalarType > QRDecompositionType
PointSet< InputPointType, NDimensions, PointSetTraitsType > PointSetType
FixedArray< Matrix< ScalarType, InputSpaceDimension, InputSpaceDimension >, OutputSpaceDimension > SpatialHessianType
vnl_vector_fixed< TScalarType, NDimensions > BMatrixType
virtual OutputPointType TransformPoint(const InputPointType &thisPoint) const
virtual void SetIdentity(void)
virtual const ParametersType & GetFixedParameters(void) const


Generated on 04-09-2015 for elastix by doxygen 1.8.9.1 elastix logo