go home Home | Main Page | Modules | Namespace List | Class Hierarchy | Alphabetical List | Data Structures | File List | Namespace Members | Data Fields | Globals | Related Pages
elxTransformIO.h
Go to the documentation of this file.
1/*=========================================================================
2 *
3 * Copyright UMC Utrecht and contributors
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 elxTransformIO_h
19#define elxTransformIO_h
20
22
23#include <itkCompositeTransform.h>
24#include <itkTransform.h>
25#include <itkTransformBase.h>
26
27#include <cassert>
28#include <string>
29
30namespace elastix
31{
32class BaseComponent;
33class Configuration;
34
36{
37public:
38 static itk::OptimizerParameters<double>
39 GetParameters(const bool fixed, const itk::TransformBase & transform)
40 {
41 return fixed ? transform.GetFixedParameters() : transform.GetParameters();
42 }
43
44 static void
45 SetParameters(const bool fixed, itk::TransformBase & transform, const itk::OptimizerParameters<double> & parameters)
46 {
47 fixed ? transform.SetFixedParameters(parameters) : transform.SetParameters(parameters);
48 }
49
50
54 static std::string
55 ConvertITKNameOfClassToElastixClassName(const std::string & itkNameOfClass);
56
57
60 template <unsigned NDimension>
61 static itk::SmartPointer<itk::CompositeTransform<double, NDimension>>
63 const itk::AdvancedCombinationTransform<double, NDimension> & advancedCombinationTransform)
64 {
65 const auto numberOfTransforms = advancedCombinationTransform.GetNumberOfTransforms();
66
67 if ((numberOfTransforms > 1) && (!advancedCombinationTransform.GetUseComposition()))
68 {
69 // A combination of multiple transforms can only be converted to CompositeTransform when the original combination
70 // uses composition.
71 return nullptr;
72 }
73
74 const auto compositeTransform = itk::CompositeTransform<double, NDimension>::New();
75
76 for (itk::SizeValueType n{}; n < numberOfTransforms; ++n)
77 {
78 const auto nthTransform = advancedCombinationTransform.GetNthTransform(n);
79 const auto singleItkTransform = ConvertToSingleItkTransform(*nthTransform);
80 compositeTransform->AddTransform((singleItkTransform == nullptr) ? nthTransform : singleItkTransform);
81 }
82 return compositeTransform;
83 }
84
85
89 template <unsigned NDimension>
90 static itk::SmartPointer<itk::CompositeTransform<double, NDimension>>
92 const itk::AdvancedCombinationTransform<double, NDimension> & advancedCombinationTransform)
93 {
94 const auto numberOfTransforms = advancedCombinationTransform.GetNumberOfTransforms();
95
96 if ((numberOfTransforms > 1) && (!advancedCombinationTransform.GetUseComposition()))
97 {
98 // A combination of multiple transforms can only be converted to CompositeTransform when the original combination
99 // uses composition.
100 return nullptr;
101 }
102
103 const auto compositeTransform = itk::CompositeTransform<double, NDimension>::New();
104
105 for (itk::SizeValueType n{}; n < numberOfTransforms; ++n)
106 {
107 const auto nthTransform = advancedCombinationTransform.GetNthTransform(n);
108 const auto singleItkTransform = ConvertToSingleItkTransform(*nthTransform);
109
110 if (singleItkTransform == nullptr)
111 {
112 return nullptr;
113 }
114 compositeTransform->AddTransform(singleItkTransform);
115 }
116 return compositeTransform;
117 }
118
119
122 template <unsigned NDimension>
123 static itk::SmartPointer<itk::Transform<double, NDimension, NDimension>>
124 ConvertToSingleItkTransform(const itk::Transform<double, NDimension, NDimension> & elxTransform)
125 {
126 // Do not use this function for elastix combination transforms!
127 using CombinationTransformType = itk::AdvancedCombinationTransform<double, NDimension>;
128 assert(dynamic_cast<const CombinationTransformType *>(&elxTransform) == nullptr);
129
130 return dynamic_cast<itk::Transform<double, NDimension, NDimension> *>(
131 ConvertItkTransformBaseToSingleItkTransform(elxTransform).GetPointer());
132 }
133
134
135 static void
136 Write(const itk::TransformBase & itkTransform, const std::string & fileName);
137
138 static itk::TransformBase::Pointer
139 Read(const std::string & fileName);
140
142 template <typename TElastixTransform>
143 static std::string
144 MakeDeformationFieldFileName(const TElastixTransform & elxTransform)
145 {
146 return MakeDeformationFieldFileName(*(elxTransform.GetConfiguration()),
147 elxTransform.GetElastix()->GetCurrentTransformParameterFileName());
148 }
149
150private:
151 static itk::TransformBase::Pointer
152 ConvertItkTransformBaseToSingleItkTransform(const itk::TransformBase & elxTransform);
153
154 static std::string
155 MakeDeformationFieldFileName(Configuration & configuration, const std::string & transformParameterFileName);
156};
157} // namespace elastix
158
159#endif
A class that deals with user given parameters and command line arguments.
static itk::SmartPointer< itk::Transform< double, NDimension, NDimension > > ConvertToSingleItkTransform(const itk::Transform< double, NDimension, NDimension > &elxTransform)
static std::string MakeDeformationFieldFileName(const TElastixTransform &elxTransform)
Makes the deformation field file name, as used by BSplineTransformWithDiffusion and DeformationFieldT...
static itk::SmartPointer< itk::CompositeTransform< double, NDimension > > ConvertToCompositionOfItkTransforms(const itk::AdvancedCombinationTransform< double, NDimension > &advancedCombinationTransform)
static itk::TransformBase::Pointer Read(const std::string &fileName)
static itk::SmartPointer< itk::CompositeTransform< double, NDimension > > ConvertToItkCompositeTransform(const itk::AdvancedCombinationTransform< double, NDimension > &advancedCombinationTransform)
static void Write(const itk::TransformBase &itkTransform, const std::string &fileName)
static void SetParameters(const bool fixed, itk::TransformBase &transform, const itk::OptimizerParameters< double > &parameters)
static std::string MakeDeformationFieldFileName(Configuration &configuration, const std::string &transformParameterFileName)
static itk::OptimizerParameters< double > GetParameters(const bool fixed, const itk::TransformBase &transform)
static std::string ConvertITKNameOfClassToElastixClassName(const std::string &itkNameOfClass)
static itk::TransformBase::Pointer ConvertItkTransformBaseToSingleItkTransform(const itk::TransformBase &elxTransform)
This class combines two transforms: an 'initial transform' with a 'current transform'.
SizeValueType GetNumberOfTransforms() const
const TransformTypePointer GetNthTransform(SizeValueType n) const
virtual bool GetUseComposition() const


Generated on 2023-01-13 for elastix by doxygen 1.9.6 elastix logo