ITK  4.13.0
Insight Segmentation and Registration Toolkit
itkDisplacementFieldTransform.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright Insight Software Consortium
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 #ifndef itkDisplacementFieldTransform_h
19 #define itkDisplacementFieldTransform_h
20 
21 #include "itkTransform.h"
22 
23 #include "itkImage.h"
27 
28 namespace itk
29 {
30 
85 template
86 <typename TParametersValueType, unsigned int NDimensions>
87 class ITK_TEMPLATE_EXPORT DisplacementFieldTransform :
88  public Transform<TParametersValueType, NDimensions, NDimensions>
89 {
90 public:
96 
98  itkTypeMacro( DisplacementFieldTransform, Transform );
99 
101  itkNewMacro( Self );
102 
104  typedef typename Superclass::InverseTransformBasePointer InverseTransformBasePointer;
105 
107  typedef typename Superclass::ScalarType ScalarType;
108 
110  typedef typename Superclass::FixedParametersType FixedParametersType;
111  typedef typename Superclass::FixedParametersValueType FixedParametersValueType;
112  typedef typename Superclass::ParametersType ParametersType;
113  typedef typename Superclass::ParametersValueType ParametersValueType;
114 
116  typedef typename Superclass::JacobianType JacobianType;
117 
119  typedef typename Superclass::TransformCategoryType TransformCategoryType;
120 
122  typedef typename Superclass::NumberOfParametersType NumberOfParametersType;
123 
125  typedef typename Superclass::InputPointType InputPointType;
126  typedef typename Superclass::OutputPointType OutputPointType;
127 
129  typedef typename Superclass::InputVectorType InputVectorType;
130  typedef typename Superclass::OutputVectorType OutputVectorType;
131 
132  typedef typename Superclass::InputVectorPixelType InputVectorPixelType;
133  typedef typename Superclass::OutputVectorPixelType OutputVectorPixelType;
134 
136  typedef typename Superclass::InputCovariantVectorType
138  typedef typename Superclass::OutputCovariantVectorType
140 
142  typedef typename Superclass::InputVnlVectorType InputVnlVectorType;
143  typedef typename Superclass::OutputVnlVectorType OutputVnlVectorType;
144 
146  typedef typename Superclass::InputDiffusionTensor3DType
148  typedef typename Superclass::OutputDiffusionTensor3DType
150 
156 
158  typedef typename Superclass::DerivativeType DerivativeType;
159 
161  itkStaticConstMacro( Dimension, unsigned int, NDimensions );
162 
167 
170 
179 
182  ScalarType,
184  Dimension>
186 
191  virtual void SetDisplacementField( DisplacementFieldType* field );
192  itkGetModifiableObjectMacro(DisplacementField, DisplacementFieldType );
194 
197  virtual void SetInverseDisplacementField( DisplacementFieldType * inverseDisplacementField );
198  itkGetModifiableObjectMacro(InverseDisplacementField, DisplacementFieldType );
200 
203  virtual void SetInterpolator( InterpolatorType* interpolator );
204  itkGetModifiableObjectMacro( Interpolator, InterpolatorType );
206 
209  virtual void SetInverseInterpolator( InterpolatorType* interpolator );
210  itkGetModifiableObjectMacro(InverseInterpolator, InterpolatorType );
212 
214  itkGetConstReferenceMacro( DisplacementFieldSetTime, ModifiedTimeType );
215 
218  virtual OutputPointType TransformPoint( const InputPointType& thisPoint ) const ITK_OVERRIDE;
219 
221  using Superclass::TransformVector;
222  virtual OutputVectorType TransformVector(const InputVectorType &) const ITK_OVERRIDE
223  {
224  itkExceptionMacro( "TransformVector(Vector) unimplemented, use "
225  "TransformVector(Vector,Point)" );
226  }
228 
229  virtual OutputVectorPixelType TransformVector(const InputVectorPixelType &) const ITK_OVERRIDE
230  {
231  itkExceptionMacro( "TransformVector(Vector) unimplemented, use "
232  "TransformVector(Vector,Point)" );
233  }
234 
235  virtual OutputVnlVectorType TransformVector(const InputVnlVectorType &) const ITK_OVERRIDE
236  {
237  itkExceptionMacro( "TransformVector(Vector) unimplemented, use "
238  "TransformVector(Vector,Point)" );
239  }
240 
242  using Superclass::TransformDiffusionTensor3D;
244  const InputDiffusionTensor3DType & ) const
245  {
246  itkExceptionMacro( "TransformDiffusionTensor(Tensor) unimplemented, use "
247  "TransformDiffusionTensor(Tensor,Point)" );
248  }
250 
252  const
253  {
254  itkExceptionMacro( "TransformDiffusionTensor(Tensor) unimplemented, use "
255  "TransformDiffusionTensor(Tensor,Point)" );
256  }
257 
259  using Superclass::TransformCovariantVector;
261  const InputCovariantVectorType &) const ITK_OVERRIDE
262  {
263  itkExceptionMacro( "TransformCovariantVector(CovariantVector) "
264  "unimplemented, use TransformCovariantVector(CovariantVector,Point)" );
265  }
267 
269  const InputVectorPixelType &) const ITK_OVERRIDE
270  {
271  itkExceptionMacro( "TransformCovariantVector(CovariantVector) "
272  "unimplemented, use TransformCovariantVector(CovariantVector,Point)" );
273  }
274 
277  virtual void SetParameters(const ParametersType & params) ITK_OVERRIDE
278  {
279  if( &(this->m_Parameters) != &params )
280  {
281  if( params.Size() != this->m_Parameters.Size() )
282  {
283  itkExceptionMacro("Input parameters size (" << params.Size()
284  << ") does not match internal size ("
285  << this->m_Parameters.Size() << ").");
286  }
287  // Copy into existing object
288  this->m_Parameters = params;
289  this->Modified();
290  }
291  }
293 
303  virtual void SetFixedParameters( const FixedParametersType & ) ITK_OVERRIDE;
304 
326  virtual void ComputeJacobianWithRespectToParameters(const InputPointType &,
327  JacobianType & j) const ITK_OVERRIDE
328  {
329  j = this->m_IdentityJacobian;
330  }
332 
340  JacobianType & j) const
341  {
342  j = this->m_IdentityJacobian;
343  }
344 
349  virtual void ComputeJacobianWithRespectToPosition(const InputPointType & x, JacobianType & j ) const ITK_OVERRIDE;
350 
355  virtual void ComputeInverseJacobianWithRespectToPosition(const InputPointType & x, JacobianType & j ) const ITK_OVERRIDE;
356 
361  virtual void ComputeJacobianWithRespectToPosition(const IndexType & x, JacobianType & j ) const;
362 
374  virtual void GetInverseJacobianOfForwardFieldWithRespectToPosition(const InputPointType & point,
375  JacobianType & jacobian,
376  bool useSVD = false ) const;
377 
389  virtual void GetInverseJacobianOfForwardFieldWithRespectToPosition(const IndexType & index, JacobianType & jacobian,
390  bool useSVD = false ) const;
391 
392  virtual void UpdateTransformParameters( const DerivativeType & update, ScalarType factor = 1.0 ) ITK_OVERRIDE;
393 
396  bool GetInverse( Self *inverse ) const;
397 
400  virtual InverseTransformBasePointer GetInverseTransform() const ITK_OVERRIDE;
401 
402  virtual void SetIdentity();
403 
405  virtual TransformCategoryType GetTransformCategory() const ITK_OVERRIDE
406  {
407  return Self::DisplacementField;
408  }
409 
410  virtual NumberOfParametersType GetNumberOfLocalParameters(void) const ITK_OVERRIDE
411  {
412  return Dimension;
413  }
414 
422  itkSetMacro(CoordinateTolerance, double);
423  itkGetConstMacro(CoordinateTolerance, double);
425 
433  itkSetMacro(DirectionTolerance, double);
434  itkGetConstMacro(DirectionTolerance, double);
436 
437 protected:
438 
440  virtual ~DisplacementFieldTransform() ITK_OVERRIDE;
441  void PrintSelf( std::ostream& os, Indent indent ) const ITK_OVERRIDE;
442 
444  typename DisplacementFieldType::Pointer m_DisplacementField;
445  typename DisplacementFieldType::Pointer m_InverseDisplacementField;
446 
448  typename InterpolatorType::Pointer m_Interpolator;
449  typename InterpolatorType::Pointer m_InverseInterpolator;
450 
453  ModifiedTimeType m_DisplacementFieldSetTime;
454 
457  JacobianType m_IdentityJacobian;
458 
459 private:
460  ITK_DISALLOW_COPY_AND_ASSIGN(DisplacementFieldTransform);
461 
470  virtual void ComputeJacobianWithRespectToPositionInternal(const IndexType & index, JacobianType & jacobian,
471  bool doInverseJacobian) const;
472 
477  virtual void VerifyFixedParametersInformation();
478 
483  virtual void SetFixedParametersFromDisplacementField() const;
484 
485  double m_CoordinateTolerance;
486  double m_DirectionTolerance;
487 
488 };
489 
490 } // end namespace itk
491 
492 #ifndef ITK_MANUAL_INSTANTIATION
493 #include "itkDisplacementFieldTransform.hxx"
494 #endif
495 
496 #endif // itkDisplacementFieldTransform_h
CovariantVector< ScalarType, OutputDiffusionTensor3DType::Dimension > OutputTensorEigenVectorType
Superclass::RegionType RegionType
Definition: itkImage.h:137
Superclass::OutputVectorPixelType OutputVectorPixelType
DisplacementFieldType::IndexType IndexType
Transform< TParametersValueType, NDimensions, NDimensions > Superclass
Light weight base class for most itk classes.
Image< OutputVectorType, Dimension > DisplacementFieldType
IdentifierType NumberOfParametersType
virtual OutputVectorType TransformVector(const InputVectorType &) const override
VectorInterpolateImageFunction< DisplacementFieldType, ScalarType > InterpolatorType
DisplacementFieldType::SizeType SizeType
Superclass::NumberOfParametersType NumberOfParametersType
Provides local/dense/high-dimensionaltiy transformation via a a displacement field.
virtual void ComputeJacobianWithRespectToParameters(const IndexType &, JacobianType &j) const
unsigned long ModifiedTimeType
Definition: itkIntTypes.h:164
DisplacementFieldType::PointType PointType
Superclass::ParametersValueType ParametersValueType
Superclass::InputVectorPixelType InputVectorPixelType
OutputDiffusionTensor3DType TransformDiffusionTensor(const InputDiffusionTensor3DType &) const
TPixel PixelType
Definition: itkImage.h:89
Superclass::OutputCovariantVectorType OutputCovariantVectorType
Superclass::InputCovariantVectorType InputCovariantVectorType
ImageVectorOptimizerParametersHelper< ScalarType, OutputVectorType::Dimension, Dimension > OptimizerParametersHelperType
Transform points and vectors from an input space to an output space.
Definition: itkTransform.h:82
virtual OutputVectorPixelType TransformVector(const InputVectorPixelType &) const override
Class to hold and manage parameters of type Image&lt;Vector&lt;...&gt;,...&gt;, used in Transforms, etc.
virtual NumberOfParametersType GetNumberOfLocalParameters(void) const override
Class to hold and manage different parameter types used during optimization.
Superclass::IndexType IndexType
Definition: itkImage.h:119
virtual OutputCovariantVectorType TransformCovariantVector(const InputCovariantVectorType &) const override
DisplacementFieldType::RegionType RegionType
Superclass::SizeType SizeType
Definition: itkImage.h:126
Superclass::FixedParametersType FixedParametersType
virtual OutputVnlVectorType TransformVector(const InputVnlVectorType &) const override
DisplacementFieldType::DirectionType DirectionType
DisplacementFieldType::PixelType PixelType
DisplacementFieldType::SpacingType SpacingType
virtual void SetParameters(const ParametersType &params) override
Superclass::FixedParametersValueType FixedParametersValueType
Superclass::OutputVnlVectorType OutputVnlVectorType
Superclass::InverseTransformBasePointer InverseTransformBasePointer
Superclass::OutputVectorType OutputVectorType
Control indentation during Print() invocation.
Definition: itkIndent.h:49
Base class for all vector image interpolaters.
Superclass::InputVnlVectorType InputVnlVectorType
DisplacementFieldType::ConstPointer DisplacementFieldConstPointer
Superclass::OutputDiffusionTensor3DType OutputDiffusionTensor3DType
Superclass::InputDiffusionTensor3DType InputDiffusionTensor3DType
DisplacementFieldType::Pointer DisplacementFieldPointer
OutputVectorPixelType TransformDiffusionTensor(const InputVectorPixelType &) const
CovariantVector< ScalarType, InputDiffusionTensor3DType::Dimension > InputTensorEigenVectorType
A templated class holding a n-Dimensional covariant vector.
virtual OutputVectorPixelType TransformCovariantVector(const InputVectorPixelType &) const override
Superclass::TransformCategoryType TransformCategoryType
Templated n-dimensional image class.
Definition: itkImage.h:75