Go to the documentation of this file.00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017 #ifndef __itkKernelTransform2_h
00018 #define __itkKernelTransform2_h
00019
00020 #include "itkAdvancedTransform.h"
00021 #include "itkPoint.h"
00022 #include "itkVector.h"
00023 #include "itkMatrix.h"
00024 #include "itkPointSet.h"
00025 #include <deque>
00026 #include <math.h>
00027 #include "vnl/vnl_matrix_fixed.h"
00028 #include "vnl/vnl_matrix.h"
00029 #include "vnl/vnl_vector.h"
00030 #include "vnl/vnl_vector_fixed.h"
00031 #include "vnl/vnl_sample.h"
00032 #include "vnl/algo/vnl_svd.h"
00033 #include "vnl/algo/vnl_qr.h"
00034
00035
00036 namespace itk
00037 {
00038
00078 template <class TScalarType,
00079 unsigned int NDimensions>
00080 class KernelTransform2
00081 : public AdvancedTransform<TScalarType, NDimensions,NDimensions>
00082 {
00083 public:
00084
00086 typedef KernelTransform2 Self;
00087 typedef AdvancedTransform<
00088 TScalarType, NDimensions, NDimensions > Superclass;
00089 typedef SmartPointer<Self> Pointer;
00090 typedef SmartPointer<const Self> ConstPointer;
00091
00093 itkTypeMacro( KernelTransform2, AdvancedTransform );
00094
00096 itkNewMacro( Self );
00097
00099 itkStaticConstMacro( SpaceDimension, unsigned int, NDimensions );
00100
00102 typedef typename Superclass::ScalarType ScalarType;
00103 typedef typename Superclass::ParametersType ParametersType;
00104 typedef typename Superclass::JacobianType JacobianType;
00105 typedef typename Superclass::InputPointType InputPointType;
00106 typedef typename Superclass::OutputPointType OutputPointType;
00107 typedef typename Superclass::InputVectorType InputVectorType;
00108 typedef typename Superclass::OutputVectorType OutputVectorType;
00109
00111 typedef typename Superclass
00112 ::NonZeroJacobianIndicesType NonZeroJacobianIndicesType;
00113 typedef typename Superclass::SpatialJacobianType SpatialJacobianType;
00114 typedef typename Superclass
00115 ::JacobianOfSpatialJacobianType JacobianOfSpatialJacobianType;
00116 typedef typename Superclass::SpatialHessianType SpatialHessianType;
00117 typedef typename Superclass
00118 ::JacobianOfSpatialHessianType JacobianOfSpatialHessianType;
00119 typedef typename Superclass::InternalMatrixType InternalMatrixType;
00120
00124 typedef DefaultStaticMeshTraits< TScalarType,
00125 NDimensions, NDimensions, TScalarType, TScalarType> PointSetTraitsType;
00126 typedef PointSet<InputPointType, NDimensions,
00127 PointSetTraitsType> PointSetType;
00128 typedef typename PointSetType::Pointer PointSetPointer;
00129 typedef typename PointSetType::PointsContainer PointsContainer;
00130 typedef typename PointSetType::PointsContainerIterator PointsIterator;
00131 typedef typename PointSetType::PointsContainerConstIterator PointsConstIterator;
00132
00134 typedef VectorContainer<unsigned long,InputVectorType> VectorSetType;
00135 typedef typename VectorSetType::Pointer VectorSetPointer;
00136
00138 typedef vnl_matrix_fixed<TScalarType, NDimensions, NDimensions> IMatrixType;
00139
00141 virtual unsigned int GetNumberOfParameters(void) const
00142 {
00143 return ( this->m_SourceLandmarks->GetNumberOfPoints() * SpaceDimension );
00144 }
00145
00147 itkGetObjectMacro( SourceLandmarks, PointSetType );
00148
00150 virtual void SetSourceLandmarks( PointSetType * );
00151
00153 itkGetObjectMacro( TargetLandmarks, PointSetType );
00154
00156 virtual void SetTargetLandmarks( PointSetType * );
00157
00161 itkGetObjectMacro( Displacements, VectorSetType );
00162
00164 void ComputeWMatrix( void );
00165
00167 void ComputeLInverse( void );
00168
00170 virtual OutputPointType TransformPoint( const InputPointType & thisPoint ) const;
00171
00173 virtual const JacobianType & GetJacobian( const InputPointType & point ) const;
00174
00176 virtual void GetJacobian(
00177 const InputPointType &,
00178 JacobianType &,
00179 NonZeroJacobianIndicesType & ) const;
00180
00182 virtual void SetIdentity( void );
00183
00189 virtual void SetParameters( const ParametersType & );
00190
00196 virtual void SetFixedParameters( const ParametersType & );
00197
00199 virtual void UpdateParameters( void );
00200
00202 virtual const ParametersType & GetParameters( void ) const;
00203
00205 virtual const ParametersType & GetFixedParameters( void ) const;
00206
00217 virtual void SetStiffness( double stiffness )
00218 {
00219 this->m_Stiffness = stiffness > 0 ? stiffness : 0.0;
00220 this->m_LMatrixComputed = false;
00221 this->m_LInverseComputed = false;
00222 this->m_WMatrixComputed = false;
00223 }
00224 itkGetMacro( Stiffness, double );
00225
00232 virtual void SetAlpha( TScalarType itkNotUsed( Alpha ) ) {};
00233 virtual TScalarType GetAlpha( void ) const { return -1.0; }
00234
00241 itkSetMacro( PoissonRatio, TScalarType );
00242 virtual const TScalarType GetPoissonRatio( void ) const
00243 {
00244 return this->m_PoissonRatio;
00245 };
00246
00248 itkSetMacro( MatrixInversionMethod, std::string );
00249 itkGetConstReferenceMacro( MatrixInversionMethod, std::string );
00250
00251 protected:
00252 KernelTransform2();
00253 virtual ~KernelTransform2();
00254 void PrintSelf( std::ostream& os, Indent indent ) const;
00255
00256 public:
00258 typedef vnl_matrix_fixed<TScalarType, NDimensions, NDimensions> GMatrixType;
00259
00261 typedef vnl_matrix<TScalarType> LMatrixType;
00262
00264 typedef vnl_matrix<TScalarType> KMatrixType;
00265
00267 typedef vnl_matrix<TScalarType> PMatrixType;
00268
00270 typedef vnl_matrix<TScalarType> YMatrixType;
00271
00273 typedef vnl_matrix<TScalarType> WMatrixType;
00274
00276 typedef vnl_matrix<TScalarType> DMatrixType;
00277
00279 typedef vnl_matrix_fixed<TScalarType,NDimensions,NDimensions> AMatrixType;
00280
00282 typedef vnl_vector_fixed<TScalarType,NDimensions> BMatrixType;
00283
00285 typedef vnl_matrix_fixed<TScalarType, 1, NDimensions> RowMatrixType;
00286
00288 typedef vnl_matrix_fixed<TScalarType, NDimensions, 1> ColumnMatrixType;
00289
00291 PointSetPointer m_SourceLandmarks;
00292
00294 PointSetPointer m_TargetLandmarks;
00295
00296 protected:
00297
00305 virtual void ComputeG( const InputVectorType & landmarkVector,
00306 GMatrixType & GMatrix ) const;
00307
00315 virtual void ComputeReflexiveG( PointsIterator, GMatrixType & GMatrix ) const;
00316
00320 virtual void ComputeDeformationContribution(
00321 const InputPointType & inputPoint,
00322 OutputPointType & result ) const;
00323
00325 void ComputeK( void );
00326
00328 void ComputeL( void );
00329
00331 void ComputeP( void );
00332
00334 void ComputeY( void );
00335
00337 void ComputeD( void );
00338
00343 void ReorganizeW( void );
00344
00346 double m_Stiffness;
00347
00351 VectorSetPointer m_Displacements;
00352
00354 LMatrixType m_LMatrix;
00355
00357 LMatrixType m_LMatrixInverse;
00358
00360 KMatrixType m_KMatrix;
00361
00363 PMatrixType m_PMatrix;
00364
00366 YMatrixType m_YMatrix;
00367
00369 WMatrixType m_WMatrix;
00370
00376 DMatrixType m_DMatrix;
00377
00379 AMatrixType m_AMatrix;
00380
00382 BMatrixType m_BVector;
00383
00389
00390
00392 bool m_WMatrixComputed;
00394 bool m_LMatrixComputed;
00396 bool m_LInverseComputed;
00398 bool m_LMatrixDecompositionComputed;
00399
00406 typedef vnl_svd< ScalarType > SVDDecompositionType;
00407 typedef vnl_qr< ScalarType > QRDecompositionType;
00408
00409 SVDDecompositionType * m_LMatrixDecompositionSVD;
00410 QRDecompositionType * m_LMatrixDecompositionQR;
00411
00413 IMatrixType m_I;
00414
00416 NonZeroJacobianIndicesType m_NonZeroJacobianIndices;
00417
00419 mutable NonZeroJacobianIndicesType m_NonZeroJacobianIndicesTemp;
00420
00424 bool m_FastComputationPossible;
00425
00426 private:
00427 KernelTransform2(const Self&);
00428 void operator=(const Self&);
00429
00430 TScalarType m_PoissonRatio;
00431
00433 std::string m_MatrixInversionMethod;
00434
00435 };
00436
00437 }
00438
00439 #ifndef ITK_MANUAL_INSTANTIATION
00440 #include "itkKernelTransform2.txx"
00441 #endif
00442
00443 #endif // __itkKernelTransform2_h