ITK 6.0.0
Insight Toolkit
 
Loading...
Searching...
No Matches
itkMatrixOffsetTransformBase.h
Go to the documentation of this file.
1/*=========================================================================
2 *
3 * Copyright NumFOCUS
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 * https://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 itkMatrixOffsetTransformBase_h
19#define itkMatrixOffsetTransformBase_h
20
21
22#include "itkMacro.h"
23#include "itkMatrix.h"
24#include "itkTransform.h"
25
26#include <iostream>
27
28namespace itk
29{
30
31/* MatrixOrthogonalityTolerance is a utility to
32 * allow setting the tolerance limits used for
33 * checking if a matrix meet the orthogonality
34 * constraints of being a rigid rotation matrix.
35 * The tolerance needs to be different for
36 * matrices of type float vs. double.
37 */
38template <typename T>
40
41template <>
42class ITK_TEMPLATE_EXPORT MatrixOrthogonalityTolerance<double>
43{
44public:
45 static double
47 {
48 return 1e-10;
49 }
50};
51
52template <>
53class ITK_TEMPLATE_EXPORT MatrixOrthogonalityTolerance<float>
54{
55public:
56 static float
58 {
59 return 1e-5f;
60 }
61};
62
104
105template <typename TParametersValueType = double, unsigned int VInputDimension = 3, unsigned int VOutputDimension = 3>
106class ITK_TEMPLATE_EXPORT MatrixOffsetTransformBase
107 : public Transform<TParametersValueType, VInputDimension, VOutputDimension>
108{
109public:
110 ITK_DISALLOW_COPY_AND_MOVE(MatrixOffsetTransformBase);
111
115
118
120 itkOverrideGetNameOfClassMacro(MatrixOffsetTransformBase);
121
123 itkNewMacro(Self);
124
126 static constexpr unsigned int InputSpaceDimension = VInputDimension;
127 static constexpr unsigned int OutputSpaceDimension = VOutputDimension;
128 static constexpr unsigned int ParametersDimension = VOutputDimension * (VInputDimension + 1);
131 using typename Superclass::FixedParametersType;
145 using typename Superclass::ScalarType;
146
151
155
158
166
168
170 using InputVnlVectorType = vnl_vector_fixed<TParametersValueType, Self::InputSpaceDimension>;
171 using OutputVnlVectorType = vnl_vector_fixed<TParametersValueType, Self::OutputSpaceDimension>;
172
178
182
185
187
190
192
194
198 using InverseTransformBasePointer = typename InverseTransformBaseType::Pointer;
199
203 friend class MatrixOffsetTransformBase<TParametersValueType, VOutputDimension, VInputDimension>;
204
208 virtual void
210
215 GetTransformCategory() const override
216 {
217 return Self::TransformCategoryEnum::Linear;
218 }
219
231 virtual void
232 SetMatrix(const MatrixType & matrix)
233 {
234 m_Matrix = matrix;
235 this->ComputeOffset();
237 m_MatrixMTime.Modified();
238 this->Modified();
239 return;
240 }
241
248
249 virtual const MatrixType &
250 GetMatrix() const
251 {
252 return m_Matrix;
253 }
254
263 void
265 {
266 m_Offset = offset;
267 this->ComputeTranslation();
268 this->Modified();
269 return;
270 }
271
277 const OutputVectorType &
278 GetOffset() const
279 {
280 return m_Offset;
281 }
282
305 void
306 SetCenter(const InputPointType & center)
307 {
308 m_Center = center;
309 this->ComputeOffset();
310 this->Modified();
311 return;
312 }
313
320 const InputPointType &
321 GetCenter() const
322 {
323 return m_Center;
324 }
325
332 void
333 SetTranslation(const OutputVectorType & translation)
334 {
335 m_Translation = translation;
336 this->ComputeOffset();
337 this->Modified();
338 return;
339 }
340
347 const OutputVectorType &
349 {
350 return m_Translation;
351 }
352
357 void
358 SetParameters(const ParametersType & parameters) override;
359
361 const ParametersType &
362 GetParameters() const override;
363
365 void
367
369 const FixedParametersType &
370 GetFixedParameters() const override;
371
383 void
384 Compose(const Self * other, bool pre = false);
385
393
395 TransformPoint(const InputPointType & point) const override;
396
398
400 TransformVector(const InputVectorType & vect) const override;
401
403 TransformVector(const InputVnlVectorType & vect) const override;
404
406 TransformVector(const InputVectorPixelType & vect) const override;
407
409
412
414 TransformCovariantVector(const InputVectorPixelType & vect) const override;
415
417
420
422 TransformDiffusionTensor3D(const InputVectorPixelType & tensor) const override;
423
427
429 TransformSymmetricSecondRankTensor(const InputVectorPixelType & inputTensor) const override;
430
431
441 void
443
444
448 void
451
455 void
457 InverseJacobianPositionType & jac) const override;
459
478 bool
480
483 GetInverseTransform() const override;
484
490 bool
491 IsLinear() const override
492 {
493 return true;
494 }
495
496protected:
499 const InverseMatrixType &
501
502protected:
511#if !defined(ITK_LEGACY_REMOVE)
512 [[deprecated("Removed unused constructor")]] MatrixOffsetTransformBase(const MatrixType & matrix,
513 const OutputVectorType & offset);
514#endif
515 explicit MatrixOffsetTransformBase(unsigned int paramDims = ParametersDimension);
517
519 ~MatrixOffsetTransformBase() override = default;
520
522 void
523 PrintSelf(std::ostream & os, Indent indent) const override;
524
525 const InverseMatrixType &
527 {
528 return m_InverseMatrix;
529 }
530 void
532 {
533 m_InverseMatrix = matrix;
534 m_InverseMatrixMTime.Modified();
535 }
536 bool
538 {
540 {
541 return true;
542 }
543
544 return false;
545 }
546
547 virtual void
549
550 virtual void
552
553 void
554 SetVarMatrix(const MatrixType & matrix)
555 {
556 m_Matrix = matrix;
557 m_MatrixMTime.Modified();
558 }
559
560 virtual void
562
563 void
565 {
566 m_Translation = translation;
567 }
568
569 virtual void
571
572 void
574 {
575 m_Offset = offset;
576 }
577
578 void
580 {
581 m_Center = center;
582 }
583
584 itkGetConstMacro(Singular, bool);
585
586private:
587 MatrixType m_Matrix{ MatrixType::GetIdentity() }; // Matrix of the transformation
588 OutputVectorType m_Offset{}; // Offset of the transformation
590 mutable bool m_Singular{ false }; // Is m_Inverse singular?
591
594
598}; // class MatrixOffsetTransformBase
599} // namespace itk
600
601#ifndef ITK_MANUAL_INSTANTIATION
602# include "itkMatrixOffsetTransformBase.hxx"
603#endif
604
605#endif /* itkMatrixOffsetTransformBase_h */
A templated class holding a n-Dimensional covariant vector.
Control indentation during Print() invocation.
Definition itkIndent.h:50
vnl_vector_fixed< TParametersValueType, Self::InputSpaceDimension > InputVnlVectorType
const OutputVectorType & GetTranslation() const
SymmetricSecondRankTensor< TParametersValueType, VInputDimension > InputSymmetricSecondRankTensorType
void PrintSelf(std::ostream &os, Indent indent) const override
void ComputeJacobianWithRespectToParameters(const InputPointType &p, JacobianType &jacobian) const override
void SetParameters(const ParametersType &parameters) override
OutputVectorPixelType TransformSymmetricSecondRankTensor(const InputVectorPixelType &inputTensor) const override
void SetVarMatrix(const MatrixType &matrix)
const OutputVectorType & GetOffset() const
virtual const MatrixType & GetMatrix() const
void SetVarCenter(const InputPointType &center)
vnl_matrix_fixed< ParametersValueType, VOutputDimension, VInputDimension > JacobianPositionType
InverseTransformBasePointer GetInverseTransform() const override
void SetVarInverseMatrix(const InverseMatrixType &matrix) const
void SetVarOffset(const OutputVectorType &offset)
OutputVectorType TransformVector(const InputVectorType &vect) const override
MatrixOffsetTransformBase< TParametersValueType, VOutputDimension, VInputDimension > InverseTransformType
void SetTranslation(const OutputVectorType &translation)
void ComputeJacobianWithRespectToPosition(const InputPointType &x, JacobianPositionType &jac) const override
vnl_vector_fixed< TParametersValueType, Self::OutputSpaceDimension > OutputVnlVectorType
const InverseMatrixType & GetVarInverseMatrix() const
Matrix< TParametersValueType, Self::InputSpaceDimension, Self::OutputSpaceDimension > InverseMatrixType
const ParametersType & GetParameters() const override
CovariantVector< TParametersValueType, Self::InputSpaceDimension > InputCovariantVectorType
~MatrixOffsetTransformBase() override=default
Transform< TParametersValueType, VInputDimension, VOutputDimension > Superclass
void SetCenter(const InputPointType &center)
const InverseMatrixType & GetInverseMatrix() const
CovariantVector< TParametersValueType, InputDiffusionTensor3DType::Dimension > InputTensorEigenVectorType
OutputVectorPixelType TransformDiffusionTensor3D(const InputVectorPixelType &tensor) const override
CovariantVector< TParametersValueType, Self::OutputSpaceDimension > OutputCovariantVectorType
bool GetInverse(InverseTransformType *inverse) const
const FixedParametersType & GetFixedParameters() const override
OutputPointType TransformPoint(const InputPointType &point) const override
TransformCategoryEnum GetTransformCategory() const override
OutputVnlVectorType TransformVector(const InputVnlVectorType &vect) const override
MatrixOffsetTransformBase(unsigned int paramDims=ParametersDimension)
void ComputeInverseJacobianWithRespectToPosition(const InputPointType &x, InverseJacobianPositionType &jac) const override
OutputVectorPixelType TransformCovariantVector(const InputVectorPixelType &vect) const override
void Compose(const Self *other, bool pre=false)
void SetFixedParameters(const FixedParametersType &) override
OutputVectorPixelType TransformVector(const InputVectorPixelType &vect) const override
void SetOffset(const OutputVectorType &offset)
OutputDiffusionTensor3DType TransformDiffusionTensor3D(const InputDiffusionTensor3DType &tensor) const override
const InputPointType & GetCenter() const
virtual void SetMatrix(const MatrixType &matrix)
OutputCovariantVectorType TransformCovariantVector(const InputCovariantVectorType &vec) const override
OutputSymmetricSecondRankTensorType TransformSymmetricSecondRankTensor(const InputSymmetricSecondRankTensorType &inputTensor) const override
SymmetricSecondRankTensor< TParametersValueType, VOutputDimension > OutputSymmetricSecondRankTensorType
Matrix< TParametersValueType, Self::OutputSpaceDimension, Self::InputSpaceDimension > MatrixType
void SetVarTranslation(const OutputVectorType &translation)
vnl_matrix_fixed< ParametersValueType, VInputDimension, VOutputDimension > InverseJacobianPositionType
A templated class holding a M x N size Matrix.
Definition itkMatrix.h:53
virtual void Modified() const
A templated class holding a geometric point in n-Dimensional space.
Definition itkPoint.h:54
Implements transparent reference counting.
Generate a unique, increasing time value.
OptimizerParameters< ParametersValueType > ParametersType
TParametersValueType ParametersValueType
OptimizerParameters< FixedParametersValueType > FixedParametersType
TransformBaseTemplateEnums::TransformCategory TransformCategoryEnum
VariableLengthVector< TParametersValueType > OutputVectorPixelType
virtual void ComputeInverseJacobianWithRespectToPosition(const InputPointType &pnt, InverseJacobianPositionType &jacobian) const
Transform< TParametersValueType, VOutputDimension, VInputDimension > InverseTransformBaseType
virtual void ComputeJacobianWithRespectToPosition(const InputPointType &x, JacobianPositionType &jacobian) const
OptimizerParameters< FixedParametersValueType > FixedParametersType
VariableLengthVector< TParametersValueType > InputVectorPixelType
virtual OutputVectorType TransformVector(const InputVectorType &) const
TransformBaseTemplateEnums::TransformCategory TransformCategoryEnum
DiffusionTensor3D< TParametersValueType > OutputDiffusionTensor3DType
virtual OutputSymmetricSecondRankTensorType TransformSymmetricSecondRankTensor(const InputSymmetricSecondRankTensorType &inputTensor, const InputPointType &point) const
virtual OutputDiffusionTensor3DType TransformDiffusionTensor3D(const InputDiffusionTensor3DType &tensor) const
virtual OutputCovariantVectorType TransformCovariantVector(const InputCovariantVectorType &) const
A templated class holding a n-Dimensional vector.
Definition itkVector.h:63
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....