go home Home | Main Page | Modules | Namespace List | Class Hierarchy | Alphabetical List | Data Structures | File List | Namespace Members | Data Fields | Globals | Related Pages
itkGradientDifferenceImageToImageMetric2.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/*=========================================================================
19
20 Program: Insight Segmentation & Registration Toolkit
21 Module: $RCSfile: itkGradientDifferenceImageToImageMetric2.h,v $
22 Date: $Date: 2011-29-04 14:33 $
23 Version: $Revision: 2.0 $
24
25 Copyright (c) Insight Software Consortium. All rights reserved.
26 See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.
27
28 This software is distributed WITHOUT ANY WARRANTY; without even
29 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
30 PURPOSE. See the above copyright notices for more information.
31
32=========================================================================*/
33#ifndef itkGradientDifferenceImageToImageMetric2_h
34#define itkGradientDifferenceImageToImageMetric2_h
35
37
38#include "itkSobelOperator.h"
39#include "itkNeighborhoodOperatorImageFilter.h"
40#include "itkPoint.h"
41#include "itkCastImageFilter.h"
42#include "itkResampleImageFilter.h"
43#include "itkOptimizer.h"
46
47namespace itk
48{
73template <class TFixedImage, class TMovingImage>
74class ITK_TEMPLATE_EXPORT GradientDifferenceImageToImageMetric
75 : public AdvancedImageToImageMetric<TFixedImage, TMovingImage>
76{
77public:
79
83
84 using Pointer = SmartPointer<Self>;
85 using ConstPointer = SmartPointer<const Self>;
86
88 itkNewMacro(Self);
89
91 itkTypeMacro(GradientDifferenceImageToImageMetric, ImageToImageMetric);
92
95#if defined(_MSC_VER) && (_MSC_VER == 1300)
96 using RealType = double;
97#else
98 using typename Superclass::RealType;
99#endif
100
101 using typename Superclass::TransformType;
102 using ScalarType = typename TransformType::ScalarType;
103 using typename Superclass::TransformPointer;
104 using typename Superclass::TransformParametersType;
105 using typename Superclass::TransformJacobianType;
106 using typename Superclass::InterpolatorType;
107 using InterpolatorPointer = typename InterpolatorType::Pointer;
108 using typename Superclass::MeasureType;
109 using typename Superclass::DerivativeType;
110 using typename Superclass::FixedImageType;
111 using typename Superclass::MovingImageType;
112 using typename Superclass::FixedImageConstPointer;
113 using typename Superclass::MovingImageConstPointer;
114 using FixedImagePixelType = typename TFixedImage::PixelType;
115 using MovedImagePixelType = typename TMovingImage::PixelType;
116 using MovingImageRegionType = typename MovingImageType::RegionType;
117 using ScalesType = typename Optimizer::ScalesType;
118
119 itkStaticConstMacro(FixedImageDimension, unsigned int, FixedImageType::ImageDimension);
120 itkStaticConstMacro(MovedImageDimension, unsigned int, MovingImageType::ImageDimension);
121
123 using CombinationTransformPointer = typename CombinationTransformType::Pointer;
124 using TransformedMovingImageType = itk::Image<FixedImagePixelType, Self::FixedImageDimension>;
125 using TransformMovingImageFilterType = itk::ResampleImageFilter<MovingImageType, TransformedMovingImageType>;
127 using RayCastInterpolatorPointer = typename RayCastInterpolatorType::Pointer;
128 using FixedGradientImageType = itk::Image<RealType, Self::FixedImageDimension>;
129 using CastFixedImageFilterType = itk::CastImageFilter<FixedImageType, FixedGradientImageType>;
130 using CastFixedImageFilterPointer = typename CastFixedImageFilterType::Pointer;
131 using FixedGradientPixelType = typename FixedGradientImageType::PixelType;
132 using MovedGradientImageType = itk::Image<RealType, Self::MovedImageDimension>;
133 using CastMovedImageFilterType = itk::CastImageFilter<TransformedMovingImageType, MovedGradientImageType>;
134 using CastMovedImageFilterPointer = typename CastMovedImageFilterType::Pointer;
135 using MovedGradientPixelType = typename MovedGradientImageType::PixelType;
136
138 void
139 GetDerivative(const TransformParametersType & parameters, DerivativeType & derivative) const override;
140
142 MeasureType
143 GetValue(const TransformParametersType & parameters) const override;
144
146 void
147 GetValueAndDerivative(const TransformParametersType & parameters,
148 MeasureType & Value,
149 DerivativeType & derivative) const override;
150
151 void
152 Initialize() override;
153
155 void
157
159 itkSetMacro(Scales, ScalesType);
160 itkGetConstReferenceMacro(Scales, ScalesType);
161
164 itkSetMacro(DerivativeDelta, double);
165 itkGetConstReferenceMacro(DerivativeDelta, double);
166
167protected:
170 void
171 PrintSelf(std::ostream & os, Indent indent) const override;
172
174 void
176
178 void
180
182 MeasureType
183 ComputeMeasure(const TransformParametersType & parameters, const double * subtractionFactor) const;
184
185 using FixedSobelFilter = NeighborhoodOperatorImageFilter<FixedGradientImageType, FixedGradientImageType>;
186
187 using MovedSobelFilter = NeighborhoodOperatorImageFilter<MovedGradientImageType, MovedGradientImageType>;
188
189private:
191 mutable MovedGradientPixelType m_Variance[FixedImageDimension]{};
192
194 mutable MovedGradientPixelType m_MinMovedGradient[MovedImageDimension]{};
195 mutable MovedGradientPixelType m_MaxMovedGradient[MovedImageDimension]{};
196
198 mutable FixedGradientPixelType m_MinFixedGradient[FixedImageDimension]{};
199 mutable FixedGradientPixelType m_MaxFixedGradient[FixedImageDimension]{};
200
202 typename TransformMovingImageFilterType::Pointer m_TransformMovingImageFilter{
203 TransformMovingImageFilterType::New()
204 };
205
207 CastFixedImageFilterPointer m_CastFixedImageFilter{ CastFixedImageFilterType::New() };
208
209 SobelOperator<FixedGradientPixelType, Self::FixedImageDimension> m_FixedSobelOperators[FixedImageDimension]{};
210
211 typename FixedSobelFilter::Pointer m_FixedSobelFilters[Self::FixedImageDimension]{};
212
213 ZeroFluxNeumannBoundaryCondition<MovedGradientImageType> m_MovedBoundCond{};
214 ZeroFluxNeumannBoundaryCondition<FixedGradientImageType> m_FixedBoundCond{};
215
217 CastMovedImageFilterPointer m_CastMovedImageFilter{ CastMovedImageFilterType::New() };
218
219 SobelOperator<MovedGradientPixelType, Self::MovedImageDimension> m_MovedSobelOperators[MovedImageDimension]{};
220
221 typename MovedSobelFilter::Pointer m_MovedSobelFilters[Self::MovedImageDimension]{};
222
223 ScalesType m_Scales{};
224 double m_DerivativeDelta{ 0.001 };
225 double m_Rescalingfactor{ 1.0 };
226 CombinationTransformPointer m_CombinationTransform{ CombinationTransformType::New() };
227};
228
229} // end namespace itk
230
231#ifndef ITK_MANUAL_INSTANTIATION
232# include "itkGradientDifferenceImageToImageMetric2.hxx"
233#endif
234
235#endif
This class combines two transforms: an 'initial transform' with a 'current transform'.
An extension of the ITK ImageToImageMetric. It is the intended base class for all elastix metrics.
Projective interpolation of an image at specified positions.
Computes similarity between two objects to be registered.
void GetValueAndDerivative(const TransformParametersType &parameters, MeasureType &Value, DerivativeType &derivative) const override
itkStaticConstMacro(FixedImageDimension, unsigned int, FixedImageType::ImageDimension)
~GradientDifferenceImageToImageMetric() override=default
typename CastFixedImageFilterType::Pointer CastFixedImageFilterPointer
typename CombinationTransformType::Pointer CombinationTransformPointer
itkStaticConstMacro(MovedImageDimension, unsigned int, MovingImageType::ImageDimension)
itk::CastImageFilter< TransformedMovingImageType, MovedGradientImageType > CastMovedImageFilterType
NeighborhoodOperatorImageFilter< FixedGradientImageType, FixedGradientImageType > FixedSobelFilter
void PrintSelf(std::ostream &os, Indent indent) const override
typename itk::AdvancedCombinationTransform< ScalarType, FixedImageDimension > CombinationTransformType
typename RayCastInterpolatorType::Pointer RayCastInterpolatorPointer
typename MovedGradientImageType::PixelType MovedGradientPixelType
void GetDerivative(const TransformParametersType &parameters, DerivativeType &derivative) const override
itk::Image< RealType, Self::MovedImageDimension > MovedGradientImageType
MeasureType GetValue(const TransformParametersType &parameters) const override
NeighborhoodOperatorImageFilter< MovedGradientImageType, MovedGradientImageType > MovedSobelFilter
MeasureType ComputeMeasure(const TransformParametersType &parameters, const double *subtractionFactor) const
typename CastMovedImageFilterType::Pointer CastMovedImageFilterPointer
typename FixedGradientImageType::PixelType FixedGradientPixelType
typename itk::AdvancedRayCastInterpolateImageFunction< MovingImageType, ScalarType > RayCastInterpolatorType
itk::CastImageFilter< FixedImageType, FixedGradientImageType > CastFixedImageFilterType
itk::Image< RealType, Self::FixedImageDimension > FixedGradientImageType
itk::ResampleImageFilter< MovingImageType, TransformedMovingImageType > TransformMovingImageFilterType
ITK_DISALLOW_COPY_AND_MOVE(GradientDifferenceImageToImageMetric)
itk::Image< FixedImagePixelType, Self::FixedImageDimension > TransformedMovingImageType


Generated on 2024-07-17 for elastix by doxygen 1.11.0 (9b424b03c9833626cd435af22a444888fbbb192d) elastix logo