go home Home | Main Page | Modules | Namespace List | Class Hierarchy | Alphabetical List | Data Structures | File List | Namespace Members | Data Fields | Globals | Related Pages
itkTransformRigidityPenaltyTerm.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 itkTransformRigidityPenaltyTerm_h
19#define itkTransformRigidityPenaltyTerm_h
20
22
26
28#include "itkNeighborhood.h"
29#include "itkImageRegionIterator.h"
30#include "itkNeighborhoodOperatorImageFilter.h"
31#include "itkNeighborhoodIterator.h"
32
34#include "itkGrayscaleDilateImageFilter.h"
35#include "itkBinaryBallStructuringElement.h"
36#include "itkImageRegionIterator.h"
37
38namespace itk
39{
70template <class TFixedImage, class TScalarType>
71class ITK_TEMPLATE_EXPORT TransformRigidityPenaltyTerm : public TransformPenaltyTerm<TFixedImage, TScalarType>
72{
73public:
75
79 using Pointer = SmartPointer<Self>;
80 using ConstPointer = SmartPointer<const Self>;
81
83 itkNewMacro(Self);
84
87
89 using typename Superclass::CoordinateRepresentationType;
90 using typename Superclass::MovingImageType;
91 using typename Superclass::MovingImagePixelType;
92 using typename Superclass::MovingImagePointer;
93 using typename Superclass::MovingImageConstPointer;
94 using typename Superclass::FixedImageType;
95 using typename Superclass::FixedImagePointer;
96 using typename Superclass::FixedImageConstPointer;
97 using typename Superclass::FixedImageRegionType;
98 using typename Superclass::TransformType;
99 using typename Superclass::TransformPointer;
100 using typename Superclass::InputPointType;
101 using typename Superclass::OutputPointType;
102 using typename Superclass::TransformParametersType;
103 using typename Superclass::TransformJacobianType;
104 using typename Superclass::InterpolatorType;
105 using typename Superclass::InterpolatorPointer;
106 using typename Superclass::RealType;
107 using typename Superclass::GradientPixelType;
108 using typename Superclass::GradientImageType;
109 using typename Superclass::GradientImagePointer;
110 using typename Superclass::GradientImageFilterType;
111 using typename Superclass::GradientImageFilterPointer;
112 using typename Superclass::FixedImageMaskType;
113 using typename Superclass::FixedImageMaskPointer;
114 using typename Superclass::MovingImageMaskType;
115 using typename Superclass::MovingImageMaskPointer;
116 using typename Superclass::MeasureType;
117 using typename Superclass::DerivativeType;
118 using typename Superclass::DerivativeValueType;
119 using typename Superclass::ParametersType;
120 using typename Superclass::FixedImagePixelType;
123 using typename Superclass::ScalarType;
124
133
135 using typename Superclass::SpatialJacobianType;
137 using typename Superclass::SpatialHessianType;
139 using typename Superclass::InternalMatrixType;
140
142 itkStaticConstMacro(FixedImageDimension, unsigned int, FixedImageType::ImageDimension);
143 itkStaticConstMacro(MovingImageDimension, unsigned int, FixedImageType::ImageDimension);
144 itkStaticConstMacro(ImageDimension, unsigned int, FixedImageType::ImageDimension);
145
147 void
148 Initialize() override;
149
155 using CoefficientImagePointer = typename CoefficientImageType::Pointer;
156 using CoefficientImageSpacingType = typename CoefficientImageType::SpacingType;
157
159 using NeighborhoodType = Neighborhood<ScalarType, Self::FixedImageDimension>;
160 using NeighborhoodSizeType = typename NeighborhoodType::SizeType;
161 using CoefficientImageIteratorType = ImageRegionIterator<CoefficientImageType>;
162 using NOIFType = NeighborhoodOperatorImageFilter<CoefficientImageType, CoefficientImageType>;
163 using NeighborhoodIteratorType = NeighborhoodIterator<CoefficientImageType>;
164 using RadiusType = typename NeighborhoodIteratorType::RadiusType;
165
168 using RigidityImagePointer = typename RigidityImageType::Pointer;
169 using RigidityPixelType = typename RigidityImageType::PixelType;
170 using RigidityImageRegionType = typename RigidityImageType::RegionType;
171 using RigidityImageIndexType = typename RigidityImageType::IndexType;
172 using RigidityImagePointType = typename RigidityImageType::PointType;
173 using RigidityImageIteratorType = ImageRegionIterator<RigidityImageType>;
174 using StructuringElementType = BinaryBallStructuringElement<RigidityPixelType, Self::FixedImageDimension>;
175 using SERadiusType = typename StructuringElementType::RadiusType;
176 using DilateFilterType = GrayscaleDilateImageFilter<RigidityImageType, RigidityImageType, StructuringElementType>;
177 using DilateFilterPointer = typename DilateFilterType::Pointer;
178
180 void
182
184 MeasureType
185 GetValue(const ParametersType & parameters) const override;
186
188 void
189 GetDerivative(const ParametersType & parameters, DerivativeType & derivative) const override;
190
192 void
193 BeforeThreadedGetValueAndDerivative(const TransformParametersType & parameters) const override;
194
196 void
197 GetValueAndDerivative(const ParametersType & parameters,
198 MeasureType & value,
199 DerivativeType & derivative) const override;
200
205
207 // itkSetObjectMacro( RigidityCoefficientImage, RigidityImageType );
208
210 itkSetClampMacro(LinearityConditionWeight, ScalarType, 0.0, NumericTraits<ScalarType>::max());
211 itkGetConstMacro(LinearityConditionWeight, ScalarType);
212
214 itkSetClampMacro(OrthonormalityConditionWeight, ScalarType, 0.0, NumericTraits<ScalarType>::max());
215 itkGetConstMacro(OrthonormalityConditionWeight, ScalarType);
216
218 itkSetClampMacro(PropernessConditionWeight, ScalarType, 0.0, NumericTraits<ScalarType>::max());
219 itkGetConstMacro(PropernessConditionWeight, ScalarType);
220
222 itkSetMacro(UseLinearityCondition, bool);
223
225 itkSetMacro(UseOrthonormalityCondition, bool);
226
228 itkSetMacro(UsePropernessCondition, bool);
229
233 itkSetMacro(CalculateLinearityCondition, bool);
234
238 itkSetMacro(CalculateOrthonormalityCondition, bool);
239
243 itkSetMacro(CalculatePropernessCondition, bool);
244
246 itkGetConstReferenceMacro(LinearityConditionValue, MeasureType);
247
249 itkGetConstReferenceMacro(OrthonormalityConditionValue, MeasureType);
250
252 itkGetConstReferenceMacro(PropernessConditionValue, MeasureType);
253
255 itkGetConstReferenceMacro(LinearityConditionGradientMagnitude, MeasureType);
256
258 itkGetConstReferenceMacro(OrthonormalityConditionGradientMagnitude, MeasureType);
259
261 itkGetConstReferenceMacro(PropernessConditionGradientMagnitude, MeasureType);
262
264 // itkGetConstReferenceMacro( RigidityPenaltyTermValue, MeasureType );
265
267 itkSetMacro(DilateRigidityImages, bool);
268
270 itkSetClampMacro(DilationRadiusMultiplier,
271 CoordinateRepresentationType,
272 0.1,
273 NumericTraits<CoordinateRepresentationType>::max());
274
276 itkSetObjectMacro(FixedRigidityImage, RigidityImageType);
277
279 itkSetObjectMacro(MovingRigidityImage, RigidityImageType);
280
282 itkSetMacro(UseFixedRigidityImage, bool);
283
285 itkSetMacro(UseMovingRigidityImage, bool);
286
288 void
289 FillRigidityCoefficientImage(const ParametersType & parameters) const;
290
291protected:
295 ~TransformRigidityPenaltyTerm() override = default;
296
298 void
299 PrintSelf(std::ostream & os, Indent indent) const override;
300
301private:
303 virtual void
305
307 void
309 const std::string & whichF,
310 const unsigned int WhichDimension,
311 const CoefficientImageSpacingType & spacing) const;
312
314 void
315 CreateNDOperator(NeighborhoodType & F, const std::string & whichF, const CoefficientImageSpacingType & spacing) const;
316
319 FilterSeparable(const CoefficientImageType *, const std::vector<NeighborhoodType> & Operators) const;
320
326
327 mutable MeasureType m_RigidityPenaltyTermValue;
328 mutable MeasureType m_LinearityConditionValue;
330 mutable MeasureType m_PropernessConditionValue;
334
341
343 CoordinateRepresentationType m_DilationRadiusMultiplier;
349 std::vector<DilateFilterPointer> m_FixedRigidityImageDilation;
350 std::vector<DilateFilterPointer> m_MovingRigidityImageDilation;
355};
356
357} // end namespace itk
358
359#ifndef ITK_MANUAL_INSTANTIATION
360# include "itkTransformRigidityPenaltyTerm.hxx"
361#endif
362
363#endif // #ifndef itkTransformRigidityPenaltyTerm_h
Deformable transform using a B-spline representation.
This class combines two transforms: an 'initial transform' with a 'current transform'.
typename BSplineOrder1TransformType::Pointer BSplineOrder1TransformPointer
typename ImageSamplerType::OutputVectorContainerPointer ImageSampleContainerPointer
typename ImageSamplerType::OutputVectorContainerType ImageSampleContainerType
typename BSplineOrder3TransformType::Pointer BSplineOrder3TransformPointer
typename BSplineOrder2TransformType::Pointer BSplineOrder2TransformPointer
A cost function that calculates a penalty term on a transformation.
typename TransformType::JacobianOfSpatialJacobianType JacobianOfSpatialJacobianType
typename Superclass::AdvancedTransformType TransformType
typename TransformType::SpatialHessianType SpatialHessianType
typename TransformType::SpatialJacobianType SpatialJacobianType
typename TransformType::InternalMatrixType InternalMatrixType
typename TransformType::JacobianOfSpatialHessianType JacobianOfSpatialHessianType
A cost function that calculates a rigidity penalty term.
itkStaticConstMacro(FixedImageDimension, unsigned int, FixedImageType::ImageDimension)
void FillRigidityCoefficientImage(const ParametersType &parameters) const
std::vector< DilateFilterPointer > m_MovingRigidityImageDilation
CoordinateRepresentationType m_DilationRadiusMultiplier
NeighborhoodOperatorImageFilter< CoefficientImageType, CoefficientImageType > NOIFType
CoefficientImagePointer FilterSeparable(const CoefficientImageType *, const std::vector< NeighborhoodType > &Operators) const
MeasureType GetValue(const ParametersType &parameters) const override
void BeforeThreadedGetValueAndDerivative(const TransformParametersType &parameters) const override
typename BSplineTransformType::SpacingType GridSpacingType
itkStaticConstMacro(ImageDimension, unsigned int, FixedImageType::ImageDimension)
typename CoefficientImageType::SpacingType CoefficientImageSpacingType
typename NeighborhoodType::SizeType NeighborhoodSizeType
itkStaticConstMacro(MovingImageDimension, unsigned int, FixedImageType::ImageDimension)
typename RigidityImageType::RegionType RigidityImageRegionType
GrayscaleDilateImageFilter< RigidityImageType, RigidityImageType, StructuringElementType > DilateFilterType
void Create1DOperator(NeighborhoodType &F, const std::string &whichF, const unsigned int WhichDimension, const CoefficientImageSpacingType &spacing) const
void CreateNDOperator(NeighborhoodType &F, const std::string &whichF, const CoefficientImageSpacingType &spacing) const
Neighborhood< ScalarType, Self::FixedImageDimension > NeighborhoodType
void GetDerivative(const ParametersType &parameters, DerivativeType &derivative) const override
NeighborhoodIterator< CoefficientImageType > NeighborhoodIteratorType
typename DilateFilterType::Pointer DilateFilterPointer
BinaryBallStructuringElement< RigidityPixelType, Self::FixedImageDimension > StructuringElementType
typename BSplineTransformType::ImageType CoefficientImageType
typename BSplineTransformType::Pointer BSplineTransformPointer
void GetValueAndDerivative(const ParametersType &parameters, MeasureType &value, DerivativeType &derivative) const override
ImageRegionIterator< RigidityImageType > RigidityImageIteratorType
typename RigidityImageType::Pointer RigidityImagePointer
typename RigidityImageType::PointType RigidityImagePointType
std::vector< DilateFilterPointer > m_FixedRigidityImageDilation
void PrintSelf(std::ostream &os, Indent indent) const override
typename StructuringElementType::RadiusType SERadiusType
typename CoefficientImageType::Pointer CoefficientImagePointer
typename RigidityImageType::PixelType RigidityPixelType
ITK_DISALLOW_COPY_AND_MOVE(TransformRigidityPenaltyTerm)
typename RigidityImageType::IndexType RigidityImageIndexType
~TransformRigidityPenaltyTerm() override=default
ImageRegionIterator< CoefficientImageType > CoefficientImageIteratorType
typename NeighborhoodIteratorType::RadiusType RadiusType


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