go home Home | Main Page | Modules | Namespace List | Class Hierarchy | Alphabetical List | Data Structures | File List | Namespace Members | Data Fields | Globals | Related Pages
itkAdvancedImageToImageMetric.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 __itkAdvancedImageToImageMetric_h
19 #define __itkAdvancedImageToImageMetric_h
20 
21 #include "itkImageToImageMetric.h"
22 
23 #include "itkImageSamplerBase.h"
24 #include "itkGradientImageFilter.h"
25 #include "itkBSplineInterpolateImageFunction.h"
28 #include "itkLimiterFunctionBase.h"
29 #include "itkFixedArray.h"
30 #include "itkAdvancedTransform.h"
31 #include "vnl/vnl_sparse_matrix.h"
32 
33 // Needed for checking for B-spline for faster implementation
36 
37 #include "itkMultiThreader.h"
38 
39 namespace itk
40 {
41 
78 template< class TFixedImage, class TMovingImage >
80  public ImageToImageMetric< TFixedImage, TMovingImage >
81 {
82 public:
83 
86  typedef ImageToImageMetric< TFixedImage, TMovingImage > Superclass;
88  typedef SmartPointer< const Self > ConstPointer;
89 
91  itkTypeMacro( AdvancedImageToImageMetric, ImageToImageMetric );
92 
94  itkStaticConstMacro( MovingImageDimension, unsigned int,
95  TMovingImage::ImageDimension );
96  itkStaticConstMacro( FixedImageDimension, unsigned int,
97  TFixedImage::ImageDimension );
98 
100  typedef typename Superclass::CoordinateRepresentationType CoordinateRepresentationType;
101  typedef typename Superclass::MovingImageType MovingImageType;
102  typedef typename Superclass::MovingImagePixelType MovingImagePixelType;
103  typedef typename MovingImageType::Pointer MovingImagePointer;
104  typedef typename Superclass::MovingImageConstPointer MovingImageConstPointer;
105  typedef typename Superclass::FixedImageType FixedImageType;
106  typedef typename FixedImageType::Pointer FixedImagePointer;
107  typedef typename Superclass::FixedImageConstPointer FixedImageConstPointer;
108  typedef typename Superclass::FixedImageRegionType FixedImageRegionType;
109  typedef typename Superclass::TransformType TransformType;
110  typedef typename Superclass::TransformPointer TransformPointer;
111  typedef typename Superclass::InputPointType InputPointType;
112  typedef typename Superclass::OutputPointType OutputPointType;
113  typedef typename Superclass::TransformParametersType TransformParametersType;
114  typedef typename Superclass::TransformJacobianType TransformJacobianType;
115  typedef typename Superclass::InterpolatorType InterpolatorType;
116  typedef typename Superclass::InterpolatorPointer InterpolatorPointer;
117  typedef typename Superclass::RealType RealType;
118  typedef typename Superclass::GradientPixelType GradientPixelType;
119  typedef typename Superclass::GradientImageType GradientImageType;
120  typedef typename Superclass::GradientImagePointer GradientImagePointer;
121  typedef typename Superclass::GradientImageFilterType GradientImageFilterType;
122  typedef typename Superclass::GradientImageFilterPointer GradientImageFilterPointer;
123  typedef typename Superclass::FixedImageMaskType FixedImageMaskType;
124  typedef typename Superclass::FixedImageMaskPointer FixedImageMaskPointer;
125  typedef typename Superclass::MovingImageMaskType MovingImageMaskType;
126  typedef typename Superclass::MovingImageMaskPointer MovingImageMaskPointer;
127  typedef typename Superclass::MeasureType MeasureType;
128  typedef typename Superclass::DerivativeType DerivativeType;
129  typedef typename DerivativeType::ValueType DerivativeValueType;
130  typedef typename Superclass::ParametersType ParametersType;
131 
133  typedef typename FixedImageType::PixelType FixedImagePixelType;
134  typedef typename MovingImageType::RegionType MovingImageRegionType;
136 
142 
150 
152  typedef typename TransformType::ScalarType ScalarType;
153  typedef AdvancedTransform<
154  ScalarType, FixedImageDimension, MovingImageDimension > AdvancedTransformType;
156 
162 
164  typedef typename DerivativeType::ValueType HessianValueType;
165  typedef vnl_sparse_matrix< HessianValueType > HessianType;
166 
168  typedef itk::MultiThreader ThreaderType;
169  typedef typename ThreaderType::ThreadInfoStruct ThreadInfoType;
170 
174  virtual void SetTransform( AdvancedTransformType * arg )
175  {
176  this->Superclass::SetTransform( arg );
177  if( this->m_AdvancedTransform != arg )
178  {
179  this->m_AdvancedTransform = arg;
180  this->Modified();
181  }
182  }
183 
184 
186  const AdvancedTransformType * GetTransform( void ) const
187  {
188  return this->m_AdvancedTransform.GetPointer();
189  }
190 
191 
193  itkSetObjectMacro( ImageSampler, ImageSamplerType );
194  virtual ImageSamplerType * GetImageSampler( void ) const
195  {
196  return this->m_ImageSampler.GetPointer();
197  }
198 
199 
202  itkGetConstMacro( UseImageSampler, bool );
203 
207  itkSetMacro( RequiredRatioOfValidSamples, double );
208  itkGetConstMacro( RequiredRatioOfValidSamples, double );
209 
212  itkSetObjectMacro( MovingImageLimiter, MovingImageLimiterType );
213  itkGetConstObjectMacro( MovingImageLimiter, MovingImageLimiterType );
214  itkSetObjectMacro( FixedImageLimiter, FixedImageLimiterType );
215  itkGetConstObjectMacro( FixedImageLimiter, FixedImageLimiterType );
216 
223  itkSetMacro( MovingLimitRangeRatio, double );
224  itkGetConstMacro( MovingLimitRangeRatio, double );
225  itkSetMacro( FixedLimitRangeRatio, double );
226  itkGetConstMacro( FixedLimitRangeRatio, double );
227 
230  itkGetConstMacro( UseFixedImageLimiter, bool );
231  itkGetConstMacro( UseMovingImageLimiter, bool );
232 
240  itkSetMacro( UseMovingImageDerivativeScales, bool );
241  itkGetConstMacro( UseMovingImageDerivativeScales, bool );
242 
243  itkSetMacro( ScaleGradientWithRespectToMovingImageOrientation, bool );
244  itkGetConstMacro( ScaleGradientWithRespectToMovingImageOrientation, bool );
245 
246  itkSetMacro( MovingImageDerivativeScales, MovingImageDerivativeScalesType );
247  itkGetConstReferenceMacro( MovingImageDerivativeScales, MovingImageDerivativeScalesType );
248 
257  virtual void Initialize( void ) throw ( ExceptionObject );
258 
262  virtual void GetSelfHessian( const TransformParametersType & parameters, HessianType & H ) const;
263 
265  virtual void SetNumberOfThreads( ThreadIdType numberOfThreads );
266 
268  itkSetMacro( UseMetricSingleThreaded, bool );
269  itkGetConstReferenceMacro( UseMetricSingleThreaded, bool );
270  itkBooleanMacro( UseMetricSingleThreaded );
271 
273  // \todo: maybe these can be united, check base class.
274  itkSetMacro( UseMultiThread, bool );
275  itkGetConstReferenceMacro( UseMultiThread, bool );
276  itkBooleanMacro( UseMultiThread );
277 
284  const TransformParametersType & parameters ) const;
285 
286 protected:
287 
290 
292  virtual ~AdvancedImageToImageMetric();
293 
295  void PrintSelf( std::ostream & os, Indent indent ) const;
296 
300  typedef typename FixedImageType::IndexType FixedImageIndexType;
301  typedef typename FixedImageIndexType::IndexValueType FixedImageIndexValueType;
302  typedef typename MovingImageType::IndexType MovingImageIndexType;
303  typedef typename TransformType::InputPointType FixedImagePointType;
304  typedef typename TransformType::OutputPointType MovingImagePointType;
305  typedef typename InterpolatorType::ContinuousIndexType MovingImageContinuousIndexType;
306 
308  typedef BSplineInterpolateImageFunction<
309  MovingImageType, CoordinateRepresentationType, double > BSplineInterpolatorType;
310  typedef typename BSplineInterpolatorType::Pointer BSplineInterpolatorPointer;
311  typedef BSplineInterpolateImageFunction<
312  MovingImageType, CoordinateRepresentationType, float > BSplineInterpolatorFloatType;
313  typedef typename BSplineInterpolatorFloatType::Pointer BSplineInterpolatorFloatPointer;
315  MovingImageType, CoordinateRepresentationType, double > ReducedBSplineInterpolatorType;
316  typedef typename ReducedBSplineInterpolatorType::Pointer ReducedBSplineInterpolatorPointer;
318  MovingImageType, CoordinateRepresentationType > LinearInterpolatorType;
319  typedef typename LinearInterpolatorType::Pointer LinearInterpolatorPointer;
320  typedef typename BSplineInterpolatorType::CovariantVectorType MovingImageDerivativeType;
321  typedef GradientImageFilter<
322  MovingImageType, RealType, RealType > CentralDifferenceGradientFilterType;
323  typedef typename CentralDifferenceGradientFilterType::Pointer CentralDifferenceGradientFilterPointer;
324 
326  typedef typename
328 
334  mutable ImageSamplerPointer m_ImageSampler;
335 
341  BSplineInterpolatorPointer m_BSplineInterpolator;
342  BSplineInterpolatorFloatPointer m_BSplineInterpolatorFloat;
343  ReducedBSplineInterpolatorPointer m_ReducedBSplineInterpolator;
344  LinearInterpolatorPointer m_LinearInterpolator;
345  CentralDifferenceGradientFilterPointer m_CentralDifferenceGradientFilter;
346 
349  typename AdvancedTransformType::Pointer m_AdvancedTransform;
351 
353  FixedImageLimiterPointer m_FixedImageLimiter;
354  MovingImageLimiterPointer m_MovingImageLimiter;
355  FixedImagePixelType m_FixedImageTrueMin;
356  FixedImagePixelType m_FixedImageTrueMax;
357  MovingImagePixelType m_MovingImageTrueMin;
358  MovingImagePixelType m_MovingImageTrueMax;
359  FixedImageLimiterOutputType m_FixedImageMinLimit;
360  FixedImageLimiterOutputType m_FixedImageMaxLimit;
361  MovingImageLimiterOutputType m_MovingImageMinLimit;
362  MovingImageLimiterOutputType m_MovingImageMaxLimit;
363 
367  virtual inline void ThreadedGetValueAndDerivative(
368  ThreadIdType threadID ){}
369 
372  MeasureType & value, DerivativeType & derivative ) const {}
373 
375  static ITK_THREAD_RETURN_TYPE GetValueAndDerivativeThreaderCallback( void * arg );
376 
379 
381  static ITK_THREAD_RETURN_TYPE AccumulateDerivativesThreaderCallback( void * arg );
382 
387 
392  {
393  // To give the threads access to all members.
395  // Used for accumulating derivatives
396  DerivativeValueType * st_DerivativePointer;
397  DerivativeValueType st_NormalizationFactor;
398  };
399  mutable MultiThreaderParameterType m_ThreaderMetricParameters;
400 
410  // test per thread struct with padding and alignment
412  {
414  MeasureType st_Value;
415  DerivativeType st_Derivative;
416  };
417  itkPadStruct( ITK_CACHE_LINE_ALIGNMENT, GetValueAndDerivativePerThreadStruct,
418  PaddedGetValueAndDerivativePerThreadStruct );
419  itkAlignedTypedef( ITK_CACHE_LINE_ALIGNMENT, PaddedGetValueAndDerivativePerThreadStruct,
420  AlignedGetValueAndDerivativePerThreadStruct );
421  mutable AlignedGetValueAndDerivativePerThreadStruct * m_GetValueAndDerivativePerThreadVariables;
423 
425  virtual void InitializeThreadingParameters( void ) const;
426 
432  virtual void InitializeImageSampler( void ) throw ( ExceptionObject );
433 
436  itkSetMacro( UseImageSampler, bool );
437 
440  virtual void CheckNumberOfSamples(
441  unsigned long wanted, unsigned long found ) const;
442 
447  virtual void CheckForBSplineInterpolator( void );
448 
458  const MovingImagePointType & mappedPoint,
459  RealType & movingImageValue,
460  MovingImageDerivativeType * gradient ) const;
461 
467  const TransformJacobianType & jacobian,
468  const MovingImageDerivativeType & movingImageDerivative,
469  DerivativeType & imageJacobian ) const;
470 
477  virtual void CheckForAdvancedTransform( void );
478 
480  virtual void CheckForBSplineTransform( void );
481 
486  virtual bool TransformPoint(
487  const FixedImagePointType & fixedImagePoint,
488  MovingImagePointType & mappedPoint ) const;
489 
496  virtual bool EvaluateTransformJacobian(
497  const FixedImagePointType & fixedImagePoint,
498  TransformJacobianType & jacobian,
499  NonZeroJacobianIndicesType & nzji ) const;
500 
502  virtual bool IsInsideMovingMask( const MovingImagePointType & point ) const;
503 
509  virtual void ComputeFixedImageExtrema(
510  const FixedImageType * image,
511  const FixedImageRegionType & region );
512 
516  virtual void ComputeMovingImageExtrema(
517  const MovingImageType * image,
518  const MovingImageRegionType & region );
519 
522  virtual void InitializeLimiters( void );
523 
526  itkSetMacro( UseFixedImageLimiter, bool );
527  itkSetMacro( UseMovingImageLimiter, bool );
528 
529 private:
530 
531  AdvancedImageToImageMetric( const Self & ); // purposely not implemented
532  void operator=( const Self & ); // purposely not implemented
533 
543 
544  MovingImageDerivativeScalesType m_MovingImageDerivativeScales;
545 
546 };
547 
548 } // end namespace itk
549 
550 #ifndef ITK_MANUAL_INSTANTIATION
551 #include "itkAdvancedImageToImageMetric.hxx"
552 #endif
553 
554 #endif // end #ifndef __itkAdvancedImageToImageMetric_h
void PrintSelf(std::ostream &os, Indent indent) const
This class combines two transforms: an 'initial transform' with a 'current transform'.
ImageToImageMetric< TFixedImage, TMovingImage > Superclass
virtual bool TransformPoint(const FixedImagePointType &fixedImagePoint, MovingImagePointType &mappedPoint) const
AdvancedTransformType::NumberOfParametersType NumberOfParametersType
const AdvancedTransformType * GetTransform(void) const
MovingImageDerivativeScalesType m_MovingImageDerivativeScales
AdvancedTransform< ScalarType, FixedImageDimension, MovingImageDimension > AdvancedTransformType
Deformable transform using a B-spline representation.
virtual void CheckForAdvancedTransform(void)
virtual bool IsInsideMovingMask(const MovingImagePointType &point) const
LimiterFunctionBase< RealType, FixedImageDimension > FixedImageLimiterType
Superclass::CoordinateRepresentationType CoordinateRepresentationType
Superclass::TransformParametersType TransformParametersType
Superclass::FixedImageMaskPointer FixedImageMaskPointer
An extension of the ITK ImageToImageMetric. It is the intended base class for all elastix metrics...
virtual void CheckForBSplineInterpolator(void)
itkAlignedTypedef(ITK_CACHE_LINE_ALIGNMENT, PaddedGetValueAndDerivativePerThreadStruct, AlignedGetValueAndDerivativePerThreadStruct)
virtual void ComputeFixedImageExtrema(const FixedImageType *image, const FixedImageRegionType &region)
virtual void InitializeLimiters(void)
virtual void ThreadedGetValueAndDerivative(ThreadIdType threadID)
AdvancedCombinationTransform< ScalarType, FixedImageDimension > CombinationTransformType
BSplineInterpolateImageFunction< MovingImageType, CoordinateRepresentationType, float > BSplineInterpolatorFloatType
Superclass::GradientPixelType GradientPixelType
static ITK_THREAD_RETURN_TYPE AccumulateDerivativesThreaderCallback(void *arg)
FixedImageLimiterOutputType m_FixedImageMinLimit
FixedImageLimiterType::OutputType FixedImageLimiterOutputType
Superclass::FixedImageRegionType FixedImageRegionType
AdvancedTransformType::Pointer m_AdvancedTransform
virtual void AfterThreadedGetValueAndDerivative(MeasureType &value, DerivativeType &derivative) const
ImageSamplerType::OutputVectorContainerType ImageSampleContainerType
virtual void BeforeThreadedGetValueAndDerivative(const TransformParametersType &parameters) const
AlignedGetValueAndDerivativePerThreadStruct * m_GetValueAndDerivativePerThreadVariables
FixedArray< double, Self::MovingImageDimension > MovingImageDerivativeScalesType
This class is a base class for any image sampler.
BSplineInterpolateImageFunction< MovingImageType, CoordinateRepresentationType, double > BSplineInterpolatorType
ImageSamplerBase< FixedImageType > ImageSamplerType
FixedImageLimiterOutputType m_FixedImageMaxLimit
Define a front-end to the STL "vector" container that conforms to the IndexedContainerInterface.
MovingImageLimiterOutputType m_MovingImageMaxLimit
Superclass::InterpolatorPointer InterpolatorPointer
BSplineInterpolatorFloatPointer m_BSplineInterpolatorFloat
Transform maps points, vectors and covariant vectors from an input space to an output space...
BSplineInterpolatorPointer m_BSplineInterpolator
Superclass::NumberOfParametersType NumberOfParametersType
AdvancedBSplineDeformableTransform< ScalarType, FixedImageDimension, 3 > BSplineOrder3TransformType
Superclass::MovingImageMaskType MovingImageMaskType
ReducedBSplineInterpolatorPointer m_ReducedBSplineInterpolator
LimiterFunctionBase< RealType, MovingImageDimension > MovingImageLimiterType
Superclass::GradientImageType GradientImageType
Evaluates the B-Spline interpolation of an image. Spline order may be from 0 to 5.
Superclass::TransformPointer TransformPointer
virtual void SetNumberOfThreads(ThreadIdType numberOfThreads)
virtual void EvaluateTransformJacobianInnerProduct(const TransformJacobianType &jacobian, const MovingImageDerivativeType &movingImageDerivative, DerivativeType &imageJacobian) const
Superclass::GradientImageFilterType GradientImageFilterType
MovingImageType::RegionType MovingImageRegionType
MovingImageLimiterOutputType m_MovingImageMinLimit
Superclass::InterpolatorType InterpolatorType
Base class for all ITK limiter function objects.
Superclass::TransformJacobianType TransformJacobianType
virtual void CheckForBSplineTransform(void)
virtual void CheckNumberOfSamples(unsigned long wanted, unsigned long found) const
MovingImageLimiterType::Pointer MovingImageLimiterPointer
AdvancedBSplineDeformableTransform< ScalarType, FixedImageDimension, 1 > BSplineOrder1TransformType
virtual void ComputeMovingImageExtrema(const MovingImageType *image, const MovingImageRegionType &region)
FixedImageLimiterType::Pointer FixedImageLimiterPointer
static ITK_THREAD_RETURN_TYPE GetValueAndDerivativeThreaderCallback(void *arg)
Superclass::MovingImagePixelType MovingImagePixelType
AdvancedBSplineDeformableTransform< ScalarType, FixedImageDimension, 2 > BSplineOrder2TransformType
virtual bool EvaluateMovingImageValueAndDerivative(const MovingImagePointType &mappedPoint, RealType &movingImageValue, MovingImageDerivativeType *gradient) const
virtual void GetSelfHessian(const TransformParametersType &parameters, HessianType &H) const
Superclass::MovingImageMaskPointer MovingImageMaskPointer
MovingImageLimiterType::OutputType MovingImageLimiterOutputType
ThreaderType::ThreadInfoStruct ThreadInfoType
virtual void Initialize(void)
void LaunchGetValueAndDerivativeThreaderCallback(void) const
virtual void InitializeThreadingParameters(void) const
ImageSamplerType::OutputVectorContainerPointer ImageSampleContainerPointer
Superclass::FixedImageMaskType FixedImageMaskType
virtual ImageSamplerType * GetImageSampler(void) const
MultiThreaderParameterType m_ThreaderMetricParameters
CentralDifferenceGradientFilterPointer m_CentralDifferenceGradientFilter
itkStaticConstMacro(MovingImageDimension, unsigned int, TMovingImage::ImageDimension)
Superclass::OutputType OutputType
Superclass::FixedImageConstPointer FixedImageConstPointer
Superclass::MovingImageConstPointer MovingImageConstPointer
virtual void InitializeImageSampler(void)
Linearly interpolate an image at specified positions.
itkPadStruct(ITK_CACHE_LINE_ALIGNMENT, GetValueAndDerivativePerThreadStruct, PaddedGetValueAndDerivativePerThreadStruct)
virtual void SetTransform(AdvancedTransformType *arg)
virtual bool EvaluateTransformJacobian(const FixedImagePointType &fixedImagePoint, TransformJacobianType &jacobian, NonZeroJacobianIndicesType &nzji) const
Superclass::GradientImageFilterPointer GradientImageFilterPointer
vnl_sparse_matrix< HessianValueType > HessianType
Superclass::GradientImagePointer GradientImagePointer


Generated on 04-09-2015 for elastix by doxygen 1.8.9.1 elastix logo