go home Home | Main Page | Modules | Namespace List | Class Hierarchy | Alphabetical List | Data Structures | File List | Namespace Members | Data Fields | Globals | Related Pages
itkParzenWindowHistogramImageToImageMetric.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 itkParzenWindowHistogramImageToImageMetric_h
19#define itkParzenWindowHistogramImageToImageMetric_h
20
23#include <vector>
24
25
26namespace itk
27{
74template <class TFixedImage, class TMovingImage>
76 : public AdvancedImageToImageMetric<TFixedImage, TMovingImage>
77{
78public:
80
84 using Pointer = SmartPointer<Self>;
85 using ConstPointer = SmartPointer<const Self>;
86
89
91 using typename Superclass::CoordinateRepresentationType;
92 using typename Superclass::MovingImageType;
93 using typename Superclass::MovingImagePixelType;
94 using typename Superclass::MovingImageConstPointer;
95 using typename Superclass::FixedImageType;
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;
122 using typename Superclass::ImageSamplerType;
123 using typename Superclass::ImageSamplerPointer;
131 using typename Superclass::ThreaderType;
132 using typename Superclass::ThreadInfoType;
133
135 itkStaticConstMacro(FixedImageDimension, unsigned int, FixedImageType::ImageDimension);
136
138 itkStaticConstMacro(MovingImageDimension, unsigned int, MovingImageType::ImageDimension);
139
146 void
147 Initialize() override;
148
153 void
154 GetDerivative(const ParametersType & parameters, DerivativeType & Derivative) const override;
155
161 void
162 GetValueAndDerivative(const ParametersType & parameters,
163 MeasureType & value,
164 DerivativeType & derivative) const override;
165
172 itkSetClampMacro(NumberOfFixedHistogramBins, unsigned long, 4, NumericTraits<unsigned long>::max());
173 itkGetConstMacro(NumberOfFixedHistogramBins, unsigned long);
181 itkSetClampMacro(NumberOfMovingHistogramBins, unsigned long, 4, NumericTraits<unsigned long>::max());
182 itkGetConstMacro(NumberOfMovingHistogramBins, unsigned long);
185 itkSetClampMacro(FixedKernelBSplineOrder, unsigned int, 0, 3);
186 itkGetConstMacro(FixedKernelBSplineOrder, unsigned int);
187
189 itkSetClampMacro(MovingKernelBSplineOrder, unsigned int, 0, 3);
190 itkGetConstMacro(MovingKernelBSplineOrder, unsigned int);
191
195 itkSetMacro(UseExplicitPDFDerivatives, bool);
196 itkGetConstReferenceMacro(UseExplicitPDFDerivatives, bool);
197 itkBooleanMacro(UseExplicitPDFDerivatives);
198
202 itkSetMacro(UseDerivative, bool);
203 itkGetConstMacro(UseDerivative, bool);
204
208 itkSetMacro(UseFiniteDifferenceDerivative, bool);
209 itkGetConstMacro(UseFiniteDifferenceDerivative, bool);
210
215 itkSetMacro(FiniteDifferencePerturbation, double);
216 itkGetConstMacro(FiniteDifferencePerturbation, double);
217
218protected:
221
224
226 void
227 PrintSelf(std::ostream & os, Indent indent) const override;
228
232 using typename Superclass::FixedImageIndexType;
234 using OffsetValueType = typename FixedImageType::OffsetValueType;
236 using typename Superclass::FixedImagePointType;
243
247 using MarginalPDFType = Array<PDFValueType>;
249 using JointPDFPointer = typename JointPDFType::Pointer;
251 using JointPDFDerivativesPointer = typename JointPDFDerivativesType::Pointer;
253 using IncrementalMarginalPDFPointer = typename IncrementalMarginalPDFType::Pointer;
254 using JointPDFIndexType = JointPDFType::IndexType;
255 using JointPDFRegionType = JointPDFType::RegionType;
256 using JointPDFSizeType = JointPDFType::SizeType;
257 using JointPDFDerivativesIndexType = JointPDFDerivativesType::IndexType;
258 using JointPDFDerivativesRegionType = JointPDFDerivativesType::RegionType;
259 using JointPDFDerivativesSizeType = JointPDFDerivativesType::SizeType;
260 using IncrementalMarginalPDFIndexType = IncrementalMarginalPDFType::IndexType;
261 using IncrementalMarginalPDFRegionType = IncrementalMarginalPDFType::RegionType;
262 using IncrementalMarginalPDFSizeType = IncrementalMarginalPDFType::SizeType;
263 using ParzenValueContainerType = Array<PDFValueType>;
264
268
272 mutable double m_Alpha;
273 mutable DerivativeType m_PerturbedAlphaRight;
274 mutable DerivativeType m_PerturbedAlphaLeft;
275
287 mutable JointPDFRegionType m_JointPDFWindow; // no need for mutable anymore?
294
299
301 void
303
305 void
306 ThreadedComputePDFs(ThreadIdType threadId);
307
309 void
311
313 static ITK_THREAD_RETURN_FUNCTION_CALL_CONVENTION
315
317 void
319
326 void
327 EvaluateParzenValues(double parzenWindowTerm,
328 OffsetValueType parzenWindowIndex,
329 const KernelFunctionType * kernel,
330 ParzenValueContainerType & parzenValues) const;
331
335 virtual void
336 UpdateJointPDFAndDerivatives(const RealType & fixedImageValue,
337 const RealType & movingImageValue,
338 const DerivativeType * imageJacobian,
339 const NonZeroJacobianIndicesType * nzji,
340 JointPDFType * jointPDF) const;
341
352 virtual void
353 UpdateJointPDFAndIncrementalPDFs(RealType fixedImageValue,
354 RealType movingImageValue,
355 RealType movingMaskValue,
356 const DerivativeType & movingImageValuesRight,
357 const DerivativeType & movingImageValuesLeft,
358 const DerivativeType & movingMaskValuesRight,
359 const DerivativeType & movingMaskValuesLeft,
360 const NonZeroJacobianIndicesType & nzji) const;
361
367 void
369 double factor,
370 const DerivativeType & imageJacobian,
371 const NonZeroJacobianIndicesType & nzji) const;
372
374 void
375 NormalizeJointPDF(JointPDFType * pdf, const double factor) const;
376
378 void
379 NormalizeJointPDFDerivatives(JointPDFDerivativesType * pdf, const double factor) const;
380
385 void
386 ComputeMarginalPDF(const JointPDFType * jointPDF, MarginalPDFType & marginalPDF, const unsigned int direction) const;
387
391 virtual void
393 IncrementalMarginalPDFType * fixedIncrementalMarginalPDF,
394 IncrementalMarginalPDFType * movingIncrementalMarginalPDF) const;
395
405 virtual void
406 ComputePDFsAndPDFDerivatives(const ParametersType & parameters) const;
407
431 virtual void
432 ComputePDFsAndIncrementalPDFs(const ParametersType & parameters) const;
433
442 virtual void
443 ComputePDFsSingleThreaded(const ParametersType & parameters) const;
444
445 virtual void
446 ComputePDFs(const ParametersType & parameters) const;
447
449 virtual void
451
452 virtual void
454
459 virtual void
460 GetValueAndAnalyticDerivative(const ParametersType & itkNotUsed(parameters),
461 MeasureType & itkNotUsed(value),
462 DerivativeType & itkNotUsed(derivative)) const
463 {}
464
469 virtual void
470 GetValueAndFiniteDifferenceDerivative(const ParametersType & itkNotUsed(parameters),
471 MeasureType & itkNotUsed(value),
472 DerivativeType & itkNotUsed(derivative)) const
473 {}
474
475private:
477 mutable std::vector<JointPDFPointer> m_ThreaderJointPDFs;
478
482 struct ParzenWindowHistogramMultiThreaderParameterType // can't we use the one from AdvancedImageToImageMetric ?
483 {
485 };
487
489 {
492 };
493 itkPadStruct(ITK_CACHE_LINE_ALIGNMENT,
495 PaddedParzenWindowHistogramGetValueAndDerivativePerThreadStruct);
496 itkAlignedTypedef(ITK_CACHE_LINE_ALIGNMENT,
497 PaddedParzenWindowHistogramGetValueAndDerivativePerThreadStruct,
498 AlignedParzenWindowHistogramGetValueAndDerivativePerThreadStruct);
499 mutable std::vector<AlignedParzenWindowHistogramGetValueAndDerivativePerThreadStruct>
501
511};
512
513} // end namespace itk
514
515#ifndef ITK_MANUAL_INSTANTIATION
516# include "itkParzenWindowHistogramImageToImageMetric.hxx"
517#endif
518
519#endif // end #ifndef itkParzenWindowHistogramImageToImageMetric_h
An extension of the ITK ImageToImageMetric. It is the intended base class for all elastix metrics.
typename TransformType::OutputPointType MovingImagePointType
typename ImageSamplerType::OutputVectorContainerPointer ImageSampleContainerPointer
typename MovingImageType::RegionType MovingImageRegionType
typename AdvancedTransformType::NonZeroJacobianIndicesType NonZeroJacobianIndicesType
typename FixedImageType::PixelType FixedImagePixelType
typename DerivativeType::ValueType DerivativeValueType
GradientImageFilter< MovingImageType, RealType, RealType > CentralDifferenceGradientFilterType
typename FixedImageIndexType::IndexValueType FixedImageIndexValueType
typename MovingImageType::IndexType MovingImageIndexType
typename ImageSamplerType::OutputVectorContainerType ImageSampleContainerType
BSplineInterpolateImageFunction< MovingImageType, CoordinateRepresentationType, double > BSplineInterpolatorType
typename ThreaderType::WorkUnitInfo ThreadInfoType
typename FixedImageType::IndexType FixedImageIndexType
typename BSplineInterpolatorType::CovariantVectorType MovingImageDerivativeType
typename MovingImageLimiterType::OutputType MovingImageLimiterOutputType
typename TransformType::InputPointType FixedImagePointType
typename FixedImageLimiterType::OutputType FixedImageLimiterOutputType
typename ImageSamplerType::Pointer ImageSamplerPointer
typename InterpolatorType::ContinuousIndexType MovingImageContinuousIndexType
This class is a base class for any image sampler.
Kernel used for density estimation and nonparameteric regression.
SmartPointer< Self > Pointer
Base class for all ITK limiter function objects.
A base class for image metrics based on a joint histogram computed using Parzen Windowing.
void InitializeThreadingParameters() const override
ITK_DISALLOW_COPY_AND_MOVE(ParzenWindowHistogramImageToImageMetric)
~ParzenWindowHistogramImageToImageMetric() override=default
virtual void ComputePDFsAndIncrementalPDFs(const ParametersType &parameters) const
typename AdvancedTransformType::NonZeroJacobianIndicesType NonZeroJacobianIndicesType
typename JointPDFDerivativesType::Pointer JointPDFDerivativesPointer
void NormalizeJointPDFDerivatives(JointPDFDerivativesType *pdf, const double factor) const
std::vector< AlignedParzenWindowHistogramGetValueAndDerivativePerThreadStruct > m_ParzenWindowHistogramGetValueAndDerivativePerThreadVariables
void ThreadedComputePDFs(ThreadIdType threadId)
virtual void ComputePDFsSingleThreaded(const ParametersType &parameters) const
itkPadStruct(ITK_CACHE_LINE_ALIGNMENT, ParzenWindowHistogramGetValueAndDerivativePerThreadStruct, PaddedParzenWindowHistogramGetValueAndDerivativePerThreadStruct)
static ITK_THREAD_RETURN_FUNCTION_CALL_CONVENTION ComputePDFsThreaderCallback(void *arg)
virtual void ComputePDFsAndPDFDerivatives(const ParametersType &parameters) const
void ComputeMarginalPDF(const JointPDFType *jointPDF, MarginalPDFType &marginalPDF, const unsigned int direction) const
void GetDerivative(const ParametersType &parameters, DerivativeType &Derivative) const override
itkAlignedTypedef(ITK_CACHE_LINE_ALIGNMENT, PaddedParzenWindowHistogramGetValueAndDerivativePerThreadStruct, AlignedParzenWindowHistogramGetValueAndDerivativePerThreadStruct)
void UpdateJointPDFDerivatives(const JointPDFIndexType &pdfIndex, double factor, const DerivativeType &imageJacobian, const NonZeroJacobianIndicesType &nzji) const
virtual void GetValueAndAnalyticDerivative(const ParametersType &, MeasureType &, DerivativeType &) const
virtual void ComputePDFs(const ParametersType &parameters) const
typename IncrementalMarginalPDFType::Pointer IncrementalMarginalPDFPointer
virtual void UpdateJointPDFAndIncrementalPDFs(RealType fixedImageValue, RealType movingImageValue, RealType movingMaskValue, const DerivativeType &movingImageValuesRight, const DerivativeType &movingImageValuesLeft, const DerivativeType &movingMaskValuesRight, const DerivativeType &movingMaskValuesLeft, const NonZeroJacobianIndicesType &nzji) const
void NormalizeJointPDF(JointPDFType *pdf, const double factor) const
itkStaticConstMacro(FixedImageDimension, unsigned int, FixedImageType::ImageDimension)
void GetValueAndDerivative(const ParametersType &parameters, MeasureType &value, DerivativeType &derivative) const override
void PrintSelf(std::ostream &os, Indent indent) const override
virtual void UpdateJointPDFAndDerivatives(const RealType &fixedImageValue, const RealType &movingImageValue, const DerivativeType *imageJacobian, const NonZeroJacobianIndicesType *nzji, JointPDFType *jointPDF) const
ParzenWindowHistogramMultiThreaderParameterType m_ParzenWindowHistogramThreaderParameters
virtual void ComputeIncrementalMarginalPDFs(const JointPDFDerivativesType *incrementalPDF, IncrementalMarginalPDFType *fixedIncrementalMarginalPDF, IncrementalMarginalPDFType *movingIncrementalMarginalPDF) const
void EvaluateParzenValues(double parzenWindowTerm, OffsetValueType parzenWindowIndex, const KernelFunctionType *kernel, ParzenValueContainerType &parzenValues) const
virtual void GetValueAndFiniteDifferenceDerivative(const ParametersType &, MeasureType &, DerivativeType &) const
itkStaticConstMacro(MovingImageDimension, unsigned int, MovingImageType::ImageDimension)


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