ITK  6.0.0
Insight Toolkit
itkImageRegistrationMethodv4.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 itkImageRegistrationMethodv4_h
19#define itkImageRegistrationMethodv4_h
20
21#include "itkProcessObject.h"
22
33#include "ITKRegistrationMethodsv4Export.h"
34
35#include <vector>
36
37namespace itk
38{
44{
45public:
51 enum class MetricSamplingStrategy : uint8_t
52 {
53 NONE,
54 REGULAR,
55 RANDOM
56 };
57};
58// Define how to print enumeration
59extern ITKRegistrationMethodsv4_EXPORT std::ostream &
61
112template <typename TFixedImage,
113 typename TMovingImage,
115 typename TVirtualImage = TFixedImage,
117class ITK_TEMPLATE_EXPORT ImageRegistrationMethodv4 : public ProcessObject
118{
119public:
120 ITK_DISALLOW_COPY_AND_MOVE(ImageRegistrationMethodv4);
121
127
129 itkNewMacro(Self);
130
132 static constexpr unsigned int ImageDimension = TFixedImage::ImageDimension;
133
135 itkOverrideGetNameOfClassMacro(ImageRegistrationMethodv4);
136
138 using FixedImageType = TFixedImage;
141 using FixedImagesContainerType = std::vector<FixedImageConstPointer>;
142 using MovingImageType = TMovingImage;
145 using MovingImagesContainerType = std::vector<MovingImageConstPointer>;
146
147 using PointSetType = TPointSet;
149 using PointSetsContainerType = std::vector<PointSetConstPointer>;
150
152 using OutputTransformType = TOutputTransform;
154 using RealType = typename OutputTransformType::ScalarType;
155 using DerivativeType = typename OutputTransformType::DerivativeType;
156 using DerivativeValueType = typename DerivativeType::ValueType;
157
160
163
166
168
169 using VirtualImageType = TVirtualImage;
173
177
180 using FixedImageMasksContainerType = std::vector<FixedImageMaskConstPointer>;
183 using MovingImageMasksContainerType = std::vector<MovingImageMaskConstPointer>;
184
193
196
198
201
205 using TransformParametersAdaptorsContainerType = std::vector<TransformParametersAdaptorPointer>;
206
210
213
215#if !defined(ITK_LEGACY_REMOVE)
217 using MetricSamplingStrategyType = MetricSamplingStrategyEnum;
218 static constexpr MetricSamplingStrategyEnum NONE = MetricSamplingStrategyEnum::NONE;
219 static constexpr MetricSamplingStrategyEnum REGULAR = MetricSamplingStrategyEnum::REGULAR;
220 static constexpr MetricSamplingStrategyEnum RANDOM = MetricSamplingStrategyEnum::RANDOM;
221#endif
222
223
225
227 virtual void
229 {
230 this->SetFixedImage(0, image);
231 }
232 virtual const FixedImageType *
234 {
235 return this->GetFixedImage(0);
236 }
237 virtual void
243 virtual void
245 {
246 this->SetMovingImage(0, image);
247 }
248 virtual const MovingImageType *
250 {
251 return this->GetMovingImage(0);
252 }
253 virtual void
259 virtual void
261 {
262 this->SetFixedPointSet(0, pointSet);
263 }
264 virtual const PointSetType *
266 {
267 return this->GetFixedPointSet(0);
268 }
269 virtual void
275 virtual void
277 {
278 this->SetMovingPointSet(0, pointSet);
279 }
280 virtual const PointSetType *
282 {
283 return this->GetMovingPointSet(0);
284 }
285 virtual void
291 itkSetObjectMacro(Optimizer, OptimizerType);
292 itkGetModifiableObjectMacro(Optimizer, OptimizerType);
303 void
305 itkGetConstMacro(OptimizerWeights, OptimizerWeightsType);
309 itkSetObjectMacro(Metric, MetricType);
310 itkGetModifiableObjectMacro(Metric, MetricType);
329 void
331 void
336 void
338
340 virtual void
342 itkGetConstMacro(MetricSamplingPercentagePerLevel, MetricSamplingPercentageArrayType);
346 itkSetGetDecoratedObjectInputMacro(FixedInitialTransform, InitialTransformType);
347
349 itkSetGetDecoratedObjectInputMacro(MovingInitialTransform, InitialTransformType);
350
366 itkSetGetDecoratedObjectInputMacro(InitialTransform, InitialTransformType);
367
369 void
382 void
384 itkGetConstMacro(NumberOfLevels, SizeValueType);
394 void
396 {
397 for (unsigned int level = 0; level < factors.Size(); ++level)
398 {
399 auto shrinkFactors = MakeFilled<ShrinkFactorsPerDimensionContainerType>(factors[level]);
400 this->SetShrinkFactorsPerDimension(level, shrinkFactors);
401 }
402 }
408 ShrinkFactorsPerDimensionContainerType
409 GetShrinkFactorsPerDimension(const unsigned int level) const
410 {
411 if (level >= this->m_ShrinkFactorsPerLevel.size())
412 {
413 itkExceptionMacro("Requesting level greater than the number of levels.");
414 }
415 return this->m_ShrinkFactorsPerLevel[level];
416 }
422 void
424 {
425 if (level >= this->m_ShrinkFactorsPerLevel.size())
426 {
427 this->m_ShrinkFactorsPerLevel.resize(level + 1);
428 }
429 this->m_ShrinkFactorsPerLevel[level] = factors;
430 this->Modified();
431 }
439 itkSetMacro(SmoothingSigmasPerLevel, SmoothingSigmasArrayType);
440 itkGetConstMacro(SmoothingSigmasPerLevel, SmoothingSigmasArrayType);
447 itkSetMacro(SmoothingSigmasAreSpecifiedInPhysicalUnits, bool);
448 itkGetConstMacro(SmoothingSigmasAreSpecifiedInPhysicalUnits, bool);
449 itkBooleanMacro(SmoothingSigmasAreSpecifiedInPhysicalUnits);
454 using Superclass::MakeOutput;
456
460 virtual const DecoratedOutputTransformType *
461 GetOutput() const;
466 {
467 return this->GetOutput();
468 }
469 virtual const DecoratedOutputTransformType *
471 {
472 return this->GetOutput();
473 }
474
475 virtual OutputTransformType *
477 virtual const OutputTransformType *
479
481 itkGetConstMacro(CurrentLevel, SizeValueType);
482
484 itkGetConstReferenceMacro(CurrentIteration, SizeValueType);
485
486 /* Get the current metric value. This is a helper function for reporting observations. */
487 itkGetConstReferenceMacro(CurrentMetricValue, RealType);
488
490 itkGetConstReferenceMacro(CurrentConvergenceValue, RealType);
491
493 itkGetConstReferenceMacro(IsConverged, bool);
494
498 itkSetMacro(InPlace, bool);
499 itkGetConstMacro(InPlace, bool);
500 itkBooleanMacro(InPlace);
508 itkBooleanMacro(InitializeCenterOfLinearOutputTransform);
509 itkSetMacro(InitializeCenterOfLinearOutputTransform, bool);
510 itkGetConstMacro(InitializeCenterOfLinearOutputTransform, bool);
525 void
527
528protected:
530 ~ImageRegistrationMethodv4() override = default;
531 void
532 PrintSelf(std::ostream & os, Indent indent) const override;
533
535 void
536 GenerateData() override;
537
538 virtual void
540
542 virtual void
544
548
550 virtual void
552
553 SizeValueType m_CurrentLevel{};
554 SizeValueType m_NumberOfLevels{ 0 };
555 SizeValueType m_CurrentIteration{};
556 RealType m_CurrentMetricValue{};
557 RealType m_CurrentConvergenceValue{};
558 bool m_IsConverged{};
559
560 FixedImagesContainerType m_FixedSmoothImages{};
561 MovingImagesContainerType m_MovingSmoothImages{};
562 FixedImageMasksContainerType m_FixedImageMasks{};
563 MovingImageMasksContainerType m_MovingImageMasks{};
564 VirtualImagePointer m_VirtualDomainImage{};
565 PointSetsContainerType m_FixedPointSets{};
566 PointSetsContainerType m_MovingPointSets{};
567 SizeValueType m_NumberOfFixedObjects{};
568 SizeValueType m_NumberOfMovingObjects{};
569
570 OptimizerPointer m_Optimizer{};
571 OptimizerWeightsType m_OptimizerWeights{};
572 bool m_OptimizerWeightsAreIdentity{};
573
574 MetricPointer m_Metric{};
575 MetricSamplingStrategyEnum m_MetricSamplingStrategy{};
576 MetricSamplingPercentageArrayType m_MetricSamplingPercentagePerLevel{};
577 SizeValueType m_NumberOfMetrics{};
578 int m_FirstImageMetricIndex{};
579 std::vector<ShrinkFactorsPerDimensionContainerType> m_ShrinkFactorsPerLevel{};
580 SmoothingSigmasArrayType m_SmoothingSigmasPerLevel{};
581 bool m_SmoothingSigmasAreSpecifiedInPhysicalUnits{};
582
583 bool m_ReseedIterator{};
584 int m_RandomSeed{};
585 int m_CurrentRandomSeed{};
586
587
588 TransformParametersAdaptorsContainerType m_TransformParametersAdaptorsPerLevel{};
589
590 CompositeTransformPointer m_CompositeTransform{};
591
592 // TODO: m_OutputTransform should be removed and replaced with a named input parameter for
593 // the pipeline
594 OutputTransformPointer m_OutputTransform{};
595
596
597private:
598 bool m_InPlace{};
599
600 bool m_InitializeCenterOfLinearOutputTransform{};
601
602 // helper function to create the right kind of concrete transform
603 template <typename TTransform>
604 static void
606 {
607 ptr = TTransform::New();
608 }
609
610 static void
612 {
614 }
615};
616} // end namespace itk
617
618#ifndef ITK_MANUAL_INSTANTIATION
619# include "itkImageRegistrationMethodv4.hxx"
620#endif
621
622#endif
enum type for metric sampling strategy
SizeValueType Size() const
Definition: itkArray.h:128
This class contains a list of transforms and concatenates them by composition.
Decorates any subclass of itkObject with a DataObject API.
static Pointer New()
Base class for templated image classes.
Definition: itkImageBase.h:115
Contains all enum classes for ImageRegistrationMethodv4 class.
Interface method for the current registration framework.
virtual void SetFixedImage(SizeValueType, const FixedImageType *)
virtual void InitializeRegistrationAtEachLevel(const SizeValueType)
virtual const DecoratedOutputTransformType * GetTransformOutput() const
typename DecoratedInitialTransformType::Pointer DecoratedInitialTransformPointer
typename MovingImageType::ConstPointer MovingImageConstPointer
typename OutputTransformType::Pointer OutputTransformPointer
virtual void SetMovingPointSet(const PointSetType *pointSet)
std::vector< TransformParametersAdaptorPointer > TransformParametersAdaptorsContainerType
typename DerivativeType::ValueType DerivativeValueType
virtual const DecoratedOutputTransformType * GetOutput() const
virtual const FixedImageType * GetFixedImage() const
virtual void SetMovingImage(SizeValueType, const MovingImageType *)
void PrintSelf(std::ostream &os, Indent indent) const override
typename ImageMetricType::FixedSampledPointSetType MetricSamplePointSetType
typename PointSetType::ConstPointer PointSetConstPointer
typename ShrinkFilterType::ShrinkFactorsType ShrinkFactorsPerDimensionContainerType
std::vector< MovingImageMaskConstPointer > MovingImageMasksContainerType
typename DecoratedOutputTransformType::Pointer DecoratedOutputTransformPointer
typename OutputTransformType::ScalarType RealType
std::vector< FixedImageMaskConstPointer > FixedImageMasksContainerType
virtual VirtualImageBaseConstPointer GetCurrentLevelVirtualDomainImage()
void SetNumberOfLevels(const SizeValueType)
typename VirtualImageType::Pointer VirtualImagePointer
virtual const PointSetType * GetFixedPointSet() const
typename FixedImageType::ConstPointer FixedImageConstPointer
virtual void SetFixedPointSet(const PointSetType *pointSet)
typename FixedImageMaskType::ConstPointer FixedImageMaskConstPointer
static void MakeOutputTransform(SmartPointer< InitialTransformType > &ptr)
typename CompositeTransformType::Pointer CompositeTransformPointer
typename VirtualImageBaseType::ConstPointer VirtualImageBaseConstPointer
virtual DecoratedOutputTransformType * GetOutput()
std::vector< PointSetConstPointer > PointSetsContainerType
void SetOptimizerWeights(OptimizerWeightsType &)
virtual OutputTransformType * GetModifiableTransform()
typename TransformParametersAdaptorType::Pointer TransformParametersAdaptorPointer
~ImageRegistrationMethodv4() override=default
const TransformParametersAdaptorsContainerType & GetTransformParametersAdaptorsPerLevel() const
static void MakeOutputTransform(SmartPointer< TTransform > &ptr)
virtual const MovingImageType * GetMovingImage(SizeValueType) const
typename FixedImageType::Pointer FixedImagePointer
virtual void SetMovingPointSet(SizeValueType, const PointSetType *)
void SetMetricSamplingPercentage(const RealType)
virtual void SetFixedImage(const FixedImageType *image)
virtual const PointSetType * GetFixedPointSet(SizeValueType) const
virtual const FixedImageType * GetFixedImage(SizeValueType) const
std::vector< MovingImageConstPointer > MovingImagesContainerType
void SetShrinkFactorsPerLevel(ShrinkFactorsArrayType factors)
virtual void SetMetricSamplingPercentagePerLevel(const MetricSamplingPercentageArrayType &samplingPercentages)
virtual void SetFixedPointSet(SizeValueType, const PointSetType *)
virtual DecoratedOutputTransformType * GetTransformOutput()
DataObjectPointer MakeOutput(DataObjectPointerArraySizeType) override
Make a DataObject of the correct type to used as the specified output.
virtual const PointSetType * GetMovingPointSet() const
void MetricSamplingReinitializeSeed(int seed)
virtual void SetMetricSamplePoints()
std::vector< FixedImageConstPointer > FixedImagesContainerType
virtual void SetMovingImage(const MovingImageType *image)
typename OptimizerType::Pointer OptimizerPointer
typename MovingImageMaskType::ConstPointer MovingImageMaskConstPointer
typename ImageMetricType::MovingImageMaskType MovingImageMaskType
virtual const MovingImageType * GetMovingImage() const
typename InitialTransformType::Pointer InitialTransformPointer
typename MovingImageType::Pointer MovingImagePointer
typename ImageMetricType::FixedImageMaskType FixedImageMaskType
typename MetricType::Pointer MetricPointer
typename OptimizerType::ScalesType OptimizerWeightsType
virtual const OutputTransformType * GetTransform() const
virtual const PointSetType * GetMovingPointSet(SizeValueType) const
typename OutputTransformType::DerivativeType DerivativeType
void SetShrinkFactorsPerDimension(unsigned int level, ShrinkFactorsPerDimensionContainerType factors)
ShrinkFactorsPerDimensionContainerType GetShrinkFactorsPerDimension(const unsigned int level) const
void SetTransformParametersAdaptorsPerLevel(TransformParametersAdaptorsContainerType &)
Control indentation during Print() invocation.
Definition: itkIndent.h:50
Light weight base class for most itk classes.
Base class for all object-to-object similarity metrics added in ITKv4.
This class takes one ore more ObjectToObject metrics and assigns weights to their derivatives to comp...
Abstract base for object-to-object optimizers.
Generic representation for an optimization method.
Definition: itkOptimizer.h:39
A superclass of the N-dimensional mesh structure; supports point (geometric coordinate and attribute)...
Definition: itkPointSet.h:82
The base class for all process objects (source, filters, mappers) in the Insight data processing pipe...
DataObjectPointerArray::size_type DataObjectPointerArraySizeType
Reduce the size of an image by an integer factor in each dimension.
ObjectType * GetPointer() const noexcept
Implementation of the composite pattern.
Base helper class intended for multi-resolution image registration.
Transform points and vectors from an input space to an output space.
Definition: itkTransform.h:84
A templated class holding a n-Dimensional vector.
Definition: itkVector.h:63
SmartPointer< const Self > ConstPointer
static Pointer New()
SmartPointer< Self > Pointer
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
ITKCommon_EXPORT std::ostream & operator<<(std::ostream &out, typename AnatomicalOrientation::CoordinateEnum value)
class ITK_FORWARD_EXPORT ProcessObject
Definition: itkDataObject.h:41
unsigned long SizeValueType
Definition: itkIntTypes.h:86