ITK  5.0.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:
91  ITK_DISALLOW_COPY_AND_ASSIGN(DisplacementFieldTransform);
92 
98 
100  itkTypeMacro( DisplacementFieldTransform, Transform );
101 
103  itkNewMacro( Self );
104 
106  using InverseTransformBasePointer = typename Superclass::InverseTransformBasePointer;
107 
109  using ScalarType = typename Superclass::ScalarType;
110 
112  using FixedParametersType = typename Superclass::FixedParametersType;
113  using FixedParametersValueType = typename Superclass::FixedParametersValueType;
114  using ParametersType = typename Superclass::ParametersType;
115  using ParametersValueType = typename Superclass::ParametersValueType;
116 
118  using JacobianType = typename Superclass::JacobianType;
119  using JacobianPositionType = typename Superclass::JacobianPositionType;
120  using InverseJacobianPositionType = typename Superclass::InverseJacobianPositionType;
121 
123  using TransformCategoryType = typename Superclass::TransformCategoryType;
124 
126  using NumberOfParametersType = typename Superclass::NumberOfParametersType;
127 
129  using InputPointType = typename Superclass::InputPointType;
130  using OutputPointType = typename Superclass::OutputPointType;
131 
133  using InputVectorType = typename Superclass::InputVectorType;
134  using OutputVectorType = typename Superclass::OutputVectorType;
135 
136  using InputVectorPixelType = typename Superclass::InputVectorPixelType;
137  using OutputVectorPixelType = typename Superclass::OutputVectorPixelType;
138 
140  using InputCovariantVectorType = typename Superclass::InputCovariantVectorType;
141  using OutputCovariantVectorType = typename Superclass::OutputCovariantVectorType;
142 
144  using InputVnlVectorType = typename Superclass::InputVnlVectorType;
145  using OutputVnlVectorType = typename Superclass::OutputVnlVectorType;
146 
148  using InputDiffusionTensor3DType = typename Superclass::InputDiffusionTensor3DType;
149  using OutputDiffusionTensor3DType = typename Superclass::OutputDiffusionTensor3DType;
150 
156 
158  using DerivativeType = typename Superclass::DerivativeType;
159 
161  static constexpr unsigned int Dimension = NDimensions;
162 
167 
170 
179 
182  ScalarType,
185 
190  virtual void SetDisplacementField( DisplacementFieldType* field );
191  itkGetModifiableObjectMacro(DisplacementField, DisplacementFieldType );
193 
196  virtual void SetInverseDisplacementField( DisplacementFieldType * inverseDisplacementField );
197  itkGetModifiableObjectMacro(InverseDisplacementField, DisplacementFieldType );
199 
202  virtual void SetInterpolator( InterpolatorType* interpolator );
203  itkGetModifiableObjectMacro( Interpolator, InterpolatorType );
205 
208  virtual void SetInverseInterpolator( InterpolatorType* interpolator );
209  itkGetModifiableObjectMacro(InverseInterpolator, InterpolatorType );
211 
213  itkGetConstReferenceMacro( DisplacementFieldSetTime, ModifiedTimeType );
214 
217  OutputPointType TransformPoint( const InputPointType& thisPoint ) const override;
218 
220  using Superclass::TransformVector;
222  {
223  itkExceptionMacro( "TransformVector(Vector) unimplemented, use "
224  "TransformVector(Vector,Point)" );
225  }
227 
229  {
230  itkExceptionMacro( "TransformVector(Vector) unimplemented, use "
231  "TransformVector(Vector,Point)" );
232  }
233 
235  {
236  itkExceptionMacro( "TransformVector(Vector) unimplemented, use "
237  "TransformVector(Vector,Point)" );
238  }
239 
241  using Superclass::TransformDiffusionTensor3D;
243  const InputDiffusionTensor3DType & ) const
244  {
245  itkExceptionMacro( "TransformDiffusionTensor(Tensor) unimplemented, use "
246  "TransformDiffusionTensor(Tensor,Point)" );
247  }
249 
251  const
252  {
253  itkExceptionMacro( "TransformDiffusionTensor(Tensor) unimplemented, use "
254  "TransformDiffusionTensor(Tensor,Point)" );
255  }
256 
258  using Superclass::TransformCovariantVector;
260  const InputCovariantVectorType &) const override
261  {
262  itkExceptionMacro( "TransformCovariantVector(CovariantVector) "
263  "unimplemented, use TransformCovariantVector(CovariantVector,Point)" );
264  }
266 
268  const InputVectorPixelType &) const override
269  {
270  itkExceptionMacro( "TransformCovariantVector(CovariantVector) "
271  "unimplemented, use TransformCovariantVector(CovariantVector,Point)" );
272  }
273 
276  void SetParameters(const ParametersType & params) override
277  {
278  if( &(this->m_Parameters) != &params )
279  {
280  if( params.Size() != this->m_Parameters.Size() )
281  {
282  itkExceptionMacro("Input parameters size (" << params.Size()
283  << ") does not match internal size ("
284  << this->m_Parameters.Size() << ").");
285  }
286  // Copy into existing object
287  this->m_Parameters = params;
288  this->Modified();
289  }
290  }
292 
302  void SetFixedParameters( const FixedParametersType & ) override;
303 
326  JacobianType & j) const override
327  {
328  j = this->m_IdentityJacobian;
329  }
331 
339  JacobianType & j) const
340  {
341  j = this->m_IdentityJacobian;
342  }
343 
348  void ComputeJacobianWithRespectToPosition(const InputPointType & x, JacobianPositionType & j ) const override;
349  using Superclass::ComputeJacobianWithRespectToPosition;
350 
355  void ComputeInverseJacobianWithRespectToPosition(const InputPointType & x,
356  InverseJacobianPositionType & j ) const override;
357  using Superclass::ComputeInverseJacobianWithRespectToPosition;
358 
363  virtual void ComputeJacobianWithRespectToPosition(const IndexType & x, JacobianPositionType & j ) const;
364 
376  virtual void GetInverseJacobianOfForwardFieldWithRespectToPosition(const InputPointType & point,
377  JacobianPositionType & jacobian,
378  bool useSVD = false ) const;
379 
391  virtual void GetInverseJacobianOfForwardFieldWithRespectToPosition(const IndexType & index,
392  JacobianPositionType & jacobian,
393  bool useSVD = false ) const;
394 
395  void UpdateTransformParameters( const DerivativeType & update, ScalarType factor = 1.0 ) override;
396 
399  bool GetInverse( Self *inverse ) const;
400 
403  InverseTransformBasePointer GetInverseTransform() const override;
404 
405  virtual void SetIdentity();
406 
409  {
410  return Self::DisplacementField;
411  }
412 
414  {
415  return Dimension;
416  }
417 
425  itkSetMacro(CoordinateTolerance, double);
426  itkGetConstMacro(CoordinateTolerance, double);
428 
436  itkSetMacro(DirectionTolerance, double);
437  itkGetConstMacro(DirectionTolerance, double);
439 
440 protected:
441 
443  ~DisplacementFieldTransform() override = default;
444  void PrintSelf( std::ostream& os, Indent indent ) const override;
445 
449 
453 
456  ModifiedTimeType m_DisplacementFieldSetTime{ 0 };
457 
461 
462 private:
471  virtual void ComputeJacobianWithRespectToPositionInternal(const IndexType & index,
472  JacobianPositionType & jacobian,
473  bool doInverseJacobian) const;
474 
479  virtual void VerifyFixedParametersInformation();
480 
485  virtual void SetFixedParametersFromDisplacementField() const;
486 
489 
490 };
491 
492 } // end namespace itk
493 
494 #ifndef ITK_MANUAL_INSTANTIATION
495 #include "itkDisplacementFieldTransform.hxx"
496 #endif
497 
498 #endif // itkDisplacementFieldTransform_h
void ComputeJacobianWithRespectToParameters(const InputPointType &, JacobianType &j) const override
TPixel PixelType
Definition: itkImage.h:95
typename Superclass::OutputDiffusionTensor3DType OutputDiffusionTensor3DType
Light weight base class for most itk classes.
OutputCovariantVectorType TransformCovariantVector(const InputCovariantVectorType &) const override
typename DisplacementFieldType::RegionType RegionType
typename Superclass::OutputVectorType OutputVectorType
Provides local/dense/high-dimensionaltiy transformation via a a displacement field.
virtual void ComputeJacobianWithRespectToParameters(const IndexType &, JacobianType &j) const
DisplacementFieldType::Pointer m_InverseDisplacementField
typename Superclass::JacobianPositionType JacobianPositionType
DisplacementFieldType::Pointer m_DisplacementField
OutputVnlVectorType TransformVector(const InputVnlVectorType &) const override
typename Superclass::OutputCovariantVectorType OutputCovariantVectorType
typename Superclass::InverseJacobianPositionType InverseJacobianPositionType
typename Superclass::JacobianType JacobianType
typename DisplacementFieldType::SpacingType SpacingType
typename Superclass::InputCovariantVectorType InputCovariantVectorType
typename Superclass::OutputPointType OutputPointType
OutputDiffusionTensor3DType TransformDiffusionTensor(const InputDiffusionTensor3DType &) const
typename Superclass::InverseTransformBasePointer InverseTransformBasePointer
OutputVectorPixelType TransformCovariantVector(const InputVectorPixelType &) const override
TransformCategoryType GetTransformCategory() const override
OutputVectorType TransformVector(const InputVectorType &) const override
typename DisplacementFieldType::SizeType SizeType
typename DisplacementFieldType::PixelType PixelType
Transform points and vectors from an input space to an output space.
Definition: itkTransform.h:83
typename Superclass::InputPointType InputPointType
Class to hold and manage parameters of type Image&lt;Vector&lt;...&gt;,...&gt;, used in Transforms, etc.
Class to hold and manage different parameter types used during optimization.
SizeValueType Size() const
Definition: itkArray.h:120
typename DisplacementFieldType::IndexType IndexType
OutputVectorPixelType TransformVector(const InputVectorPixelType &) const override
typename Superclass::OutputVectorPixelType OutputVectorPixelType
typename Superclass::SpacingType SpacingType
Definition: itkImage.h:143
typename Superclass::OutputVnlVectorType OutputVnlVectorType
unsigned long ModifiedTimeType
Definition: itkIntTypes.h:104
typename Superclass::InputVectorPixelType InputVectorPixelType
typename Superclass::DerivativeType DerivativeType
typename DisplacementFieldType::ConstPointer DisplacementFieldConstPointer
typename DisplacementFieldType::DirectionType DirectionType
NumberOfParametersType GetNumberOfLocalParameters() const override
typename DisplacementFieldType::Pointer DisplacementFieldPointer
typename Superclass::InputDiffusionTensor3DType InputDiffusionTensor3DType
void SetParameters(const ParametersType &params) override
IdentifierType NumberOfParametersType
Control indentation during Print() invocation.
Definition: itkIndent.h:49
Base class for all vector image interpolaters.
typename Superclass::InputVnlVectorType InputVnlVectorType
typename Superclass::ScalarType ScalarType
OutputVectorPixelType TransformDiffusionTensor(const InputVectorPixelType &) const
typename Superclass::InputVectorType InputVectorType
A templated class holding a n-Dimensional covariant vector.
Templated n-dimensional image class.
Definition: itkImage.h:75
typename DisplacementFieldType::PointType PointType
TParametersValueType ParametersValueType