itkAdvancedImageToImageMetric.h

Go to the documentation of this file.

/*=========================================================================
 *
 *  Copyright UMC Utrecht and contributors
 *
 *  Licensed under the Apache License, Version 2.0 (the "License");
 *  you may not use this file except in compliance with the License.
 *  You may obtain a copy of the License at
 *
 *        http://www.apache.org/licenses/LICENSE-2.0.txt
 *
 *  Unless required by applicable law or agreed to in writing, software
 *  distributed under the License is distributed on an "AS IS" BASIS,
 *  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 *  See the License for the specific language governing permissions and
 *  limitations under the License.
 *
 *=========================================================================*/
#ifndef __itkAdvancedImageToImageMetric_h
#define __itkAdvancedImageToImageMetric_h
 
#include "itkImageToImageMetric.h"
 
#include "itkImageSamplerBase.h"
#include "itkGradientImageFilter.h"
#include "itkBSplineInterpolateImageFunction.h"
#include "itkReducedDimensionBSplineInterpolateImageFunction.h"
#include "itkAdvancedLinearInterpolateImageFunction.h"
#include "itkLimiterFunctionBase.h"
#include "itkFixedArray.h"
#include "itkAdvancedTransform.h"
#include "vnl/vnl_sparse_matrix.h"
 
#include "itkImageMaskSpatialObject.h"
 
// Needed for checking for B-spline for faster implementation
#include "itkAdvancedBSplineDeformableTransform.h"
#include "itkAdvancedCombinationTransform.h"
 
#include "itkPlatformMultiThreader.h"
 
namespace itk
{
 
template< class TFixedImage, class TMovingImage >
class AdvancedImageToImageMetric :
  public ImageToImageMetric< TFixedImage, TMovingImage >
{
public:
 
  typedef AdvancedImageToImageMetric                      Self;
  typedef ImageToImageMetric< TFixedImage, TMovingImage > Superclass;
  typedef SmartPointer< Self >                            Pointer;
  typedef SmartPointer< const Self >                      ConstPointer;
 
  itkTypeMacro( AdvancedImageToImageMetric, ImageToImageMetric );
 
  itkStaticConstMacro( MovingImageDimension, unsigned int,
    TMovingImage::ImageDimension );
  itkStaticConstMacro( FixedImageDimension, unsigned int,
    TFixedImage::ImageDimension );
 
  typedef typename Superclass::CoordinateRepresentationType CoordinateRepresentationType;
  typedef typename Superclass::MovingImageType              MovingImageType;
  typedef typename Superclass::MovingImagePixelType         MovingImagePixelType;
  typedef typename MovingImageType::Pointer                 MovingImagePointer;
  typedef typename Superclass::MovingImageConstPointer      MovingImageConstPointer;
  typedef typename Superclass::FixedImageType               FixedImageType;
  typedef typename FixedImageType::Pointer                  FixedImagePointer;
  typedef typename Superclass::FixedImageConstPointer       FixedImageConstPointer;
  typedef typename Superclass::FixedImageRegionType         FixedImageRegionType;
  typedef typename Superclass::TransformType                TransformType;
  typedef typename Superclass::TransformPointer             TransformPointer;
  typedef typename Superclass::InputPointType               InputPointType;
  typedef typename Superclass::OutputPointType              OutputPointType;
  typedef typename Superclass::TransformParametersType      TransformParametersType;
  typedef typename Superclass::TransformJacobianType        TransformJacobianType;
  typedef typename Superclass::InterpolatorType             InterpolatorType;
  typedef typename Superclass::InterpolatorPointer          InterpolatorPointer;
  typedef typename Superclass::RealType                     RealType;
  typedef typename Superclass::GradientPixelType            GradientPixelType;
  typedef typename Superclass::GradientImageType            GradientImageType;
  typedef typename Superclass::GradientImagePointer         GradientImagePointer;
  typedef typename Superclass::GradientImageFilterType      GradientImageFilterType;
  typedef typename Superclass::GradientImageFilterPointer   GradientImageFilterPointer;
  typedef typename Superclass::FixedImageMaskType           FixedImageMaskType;
  typedef typename Superclass::FixedImageMaskPointer        FixedImageMaskPointer;
  typedef typename Superclass::MovingImageMaskType          MovingImageMaskType;
  typedef typename Superclass::MovingImageMaskPointer       MovingImageMaskPointer;
  typedef typename Superclass::MeasureType                  MeasureType;
  typedef typename Superclass::DerivativeType               DerivativeType;
  typedef typename DerivativeType::ValueType                DerivativeValueType;
  typedef typename Superclass::ParametersType               ParametersType;
 
  typedef ImageMaskSpatialObject< itkGetStaticConstMacro( FixedImageDimension ) > FixedImageMaskSpatialObject2Type;
  typedef ImageMaskSpatialObject< itkGetStaticConstMacro( MovingImageDimension ) > MovingImageMaskSpatialObject2Type;
 
  typedef typename FixedImageType::PixelType               FixedImagePixelType;
  typedef typename MovingImageType::RegionType             MovingImageRegionType;
  typedef FixedArray< double, Self::MovingImageDimension > MovingImageDerivativeScalesType;
 
  typedef ImageSamplerBase< FixedImageType >                      ImageSamplerType;
  typedef typename ImageSamplerType::Pointer                      ImageSamplerPointer;
  typedef typename ImageSamplerType::OutputVectorContainerType    ImageSampleContainerType;
  typedef typename ImageSamplerType::OutputVectorContainerPointer ImageSampleContainerPointer;
 
  typedef LimiterFunctionBase< RealType, FixedImageDimension >  FixedImageLimiterType;
  typedef typename FixedImageLimiterType::Pointer               FixedImageLimiterPointer;
  typedef typename FixedImageLimiterType::OutputType            FixedImageLimiterOutputType;
  typedef LimiterFunctionBase< RealType, MovingImageDimension > MovingImageLimiterType;
  typedef typename MovingImageLimiterType::Pointer              MovingImageLimiterPointer;
  typedef typename MovingImageLimiterType::OutputType           MovingImageLimiterOutputType;
 
  typedef typename TransformType::ScalarType ScalarType;
  typedef AdvancedTransform<
    ScalarType, FixedImageDimension, MovingImageDimension >      AdvancedTransformType;
  typedef typename AdvancedTransformType::NumberOfParametersType NumberOfParametersType;
 
  typedef AdvancedCombinationTransform< ScalarType, FixedImageDimension >          CombinationTransformType;
  typedef AdvancedBSplineDeformableTransform< ScalarType, FixedImageDimension, 1 > BSplineOrder1TransformType;
  typedef AdvancedBSplineDeformableTransform< ScalarType, FixedImageDimension, 2 > BSplineOrder2TransformType;
  typedef AdvancedBSplineDeformableTransform< ScalarType, FixedImageDimension, 3 > BSplineOrder3TransformType;
  typedef typename BSplineOrder1TransformType::Pointer                             BSplineOrder1TransformPointer;
  typedef typename BSplineOrder2TransformType::Pointer                             BSplineOrder2TransformPointer;
  typedef typename BSplineOrder3TransformType::Pointer                             BSplineOrder3TransformPointer;
 
  typedef typename DerivativeType::ValueType    HessianValueType;
  typedef vnl_sparse_matrix< HessianValueType > HessianType;
 
  typedef itk::PlatformMultiThreader                      ThreaderType;
  typedef typename ThreaderType::WorkUnitInfo ThreadInfoType;
 
  virtual void SetTransform( AdvancedTransformType * arg )
  {
    this->Superclass::SetTransform( arg );
    if( this->m_AdvancedTransform != arg )
    {
      this->m_AdvancedTransform = arg;
      this->Modified();
    }
  }
 
 
  const AdvancedTransformType * GetTransform( void ) const override
  {
    return this->m_AdvancedTransform.GetPointer();
  }
 
 
  itkSetObjectMacro( ImageSampler, ImageSamplerType );
  virtual ImageSamplerType * GetImageSampler( void ) const
  {
    return this->m_ImageSampler.GetPointer();
  }
 
 
  itkGetConstMacro( UseImageSampler, bool );
 
  itkSetMacro( RequiredRatioOfValidSamples, double );
  itkGetConstMacro( RequiredRatioOfValidSamples, double );
 
  itkSetObjectMacro( MovingImageLimiter, MovingImageLimiterType );
  itkGetConstObjectMacro( MovingImageLimiter, MovingImageLimiterType );
  itkSetObjectMacro( FixedImageLimiter, FixedImageLimiterType );
  itkGetConstObjectMacro( FixedImageLimiter, FixedImageLimiterType );
 
  itkSetMacro( MovingLimitRangeRatio, double );
  itkGetConstMacro( MovingLimitRangeRatio, double );
  itkSetMacro( FixedLimitRangeRatio, double );
  itkGetConstMacro( FixedLimitRangeRatio, double );
 
  itkGetConstMacro( UseFixedImageLimiter, bool );
  itkGetConstMacro( UseMovingImageLimiter, bool );
 
  itkSetMacro( UseMovingImageDerivativeScales, bool );
  itkGetConstMacro( UseMovingImageDerivativeScales, bool );
 
  itkSetMacro( ScaleGradientWithRespectToMovingImageOrientation, bool );
  itkGetConstMacro( ScaleGradientWithRespectToMovingImageOrientation, bool );
 
  itkSetMacro( MovingImageDerivativeScales, MovingImageDerivativeScalesType );
  itkGetConstReferenceMacro( MovingImageDerivativeScales, MovingImageDerivativeScalesType );
 
  void Initialize( void ) override;
 
  virtual void GetSelfHessian( const TransformParametersType & parameters, HessianType & H ) const;
 
  virtual void SetNumberOfWorkUnits( ThreadIdType numberOfThreads );
 
  itkSetMacro( UseMetricSingleThreaded, bool );
  itkGetConstReferenceMacro( UseMetricSingleThreaded, bool );
  itkBooleanMacro( UseMetricSingleThreaded );
 
  // \todo: maybe these can be united, check base class.
  itkSetMacro( UseMultiThread, bool );
  itkGetConstReferenceMacro( UseMultiThread, bool );
  itkBooleanMacro( UseMultiThread );
 
  virtual void BeforeThreadedGetValueAndDerivative(
    const TransformParametersType & parameters ) const;
 
protected:
 
  AdvancedImageToImageMetric();
 
  ~AdvancedImageToImageMetric() override;
 
  void PrintSelf( std::ostream & os, Indent indent ) const override;
 
  typedef typename FixedImageType::IndexType             FixedImageIndexType;
  typedef typename FixedImageIndexType::IndexValueType   FixedImageIndexValueType;
  typedef typename MovingImageType::IndexType            MovingImageIndexType;
  typedef typename TransformType::InputPointType         FixedImagePointType;
  typedef typename TransformType::OutputPointType        MovingImagePointType;
  typedef typename InterpolatorType::ContinuousIndexType MovingImageContinuousIndexType;
 
  typedef BSplineInterpolateImageFunction<
    MovingImageType, CoordinateRepresentationType, double >      BSplineInterpolatorType;
  typedef typename BSplineInterpolatorType::Pointer BSplineInterpolatorPointer;
  typedef BSplineInterpolateImageFunction<
    MovingImageType, CoordinateRepresentationType, float >       BSplineInterpolatorFloatType;
  typedef typename BSplineInterpolatorFloatType::Pointer BSplineInterpolatorFloatPointer;
  typedef ReducedDimensionBSplineInterpolateImageFunction<
    MovingImageType, CoordinateRepresentationType, double >      ReducedBSplineInterpolatorType;
  typedef typename ReducedBSplineInterpolatorType::Pointer ReducedBSplineInterpolatorPointer;
  typedef AdvancedLinearInterpolateImageFunction<
    MovingImageType, CoordinateRepresentationType >              LinearInterpolatorType;
  typedef typename LinearInterpolatorType::Pointer              LinearInterpolatorPointer;
  typedef typename BSplineInterpolatorType::CovariantVectorType MovingImageDerivativeType;
  typedef GradientImageFilter<
    MovingImageType, RealType, RealType >                        CentralDifferenceGradientFilterType;
  typedef typename CentralDifferenceGradientFilterType::Pointer CentralDifferenceGradientFilterPointer;
 
  typedef typename
    AdvancedTransformType::NonZeroJacobianIndicesType NonZeroJacobianIndicesType;
 
  mutable ImageSamplerPointer m_ImageSampler;
 
  bool                                   m_InterpolatorIsLinear;
  bool                                   m_InterpolatorIsBSpline;
  bool                                   m_InterpolatorIsBSplineFloat;
  bool                                   m_InterpolatorIsReducedBSpline;
  LinearInterpolatorPointer              m_LinearInterpolator;
  BSplineInterpolatorPointer             m_BSplineInterpolator;
  BSplineInterpolatorFloatPointer        m_BSplineInterpolatorFloat;
  ReducedBSplineInterpolatorPointer      m_ReducedBSplineInterpolator;
 
  CentralDifferenceGradientFilterPointer m_CentralDifferenceGradientFilter;
 
  bool m_TransformIsAdvanced;
  typename AdvancedTransformType::Pointer m_AdvancedTransform;
  mutable bool m_TransformIsBSpline;
 
  FixedImageLimiterPointer     m_FixedImageLimiter;
  MovingImageLimiterPointer    m_MovingImageLimiter;
  FixedImagePixelType          m_FixedImageTrueMin;
  FixedImagePixelType          m_FixedImageTrueMax;
  MovingImagePixelType         m_MovingImageTrueMin;
  MovingImagePixelType         m_MovingImageTrueMax;
  FixedImageLimiterOutputType  m_FixedImageMinLimit;
  FixedImageLimiterOutputType  m_FixedImageMaxLimit;
  MovingImageLimiterOutputType m_MovingImageMinLimit;
  MovingImageLimiterOutputType m_MovingImageMaxLimit;
 
  virtual inline void ThreadedGetValue( ThreadIdType threadID ){}
 
  virtual inline void AfterThreadedGetValue( MeasureType & value ) const {}
 
  static ITK_THREAD_RETURN_FUNCTION_CALL_CONVENTION GetValueThreaderCallback( void * arg );
 
  void LaunchGetValueThreaderCallback( void ) const;
 
  virtual inline void ThreadedGetValueAndDerivative(
    ThreadIdType threadID ){}
 
  virtual inline void AfterThreadedGetValueAndDerivative(
    MeasureType & value, DerivativeType & derivative ) const {}
 
  static ITK_THREAD_RETURN_FUNCTION_CALL_CONVENTION GetValueAndDerivativeThreaderCallback( void * arg );
 
  void LaunchGetValueAndDerivativeThreaderCallback( void ) const;
 
  static ITK_THREAD_RETURN_FUNCTION_CALL_CONVENTION AccumulateDerivativesThreaderCallback( void * arg );
 
  bool m_UseMetricSingleThreaded;
  bool m_UseMultiThread;
  bool m_UseOpenMP;
 
  struct MultiThreaderParameterType
  {
    // To give the threads access to all members.
    AdvancedImageToImageMetric * st_Metric;
    // Used for accumulating derivatives
    DerivativeValueType * st_DerivativePointer;
    DerivativeValueType   st_NormalizationFactor;
  };
  mutable MultiThreaderParameterType m_ThreaderMetricParameters;
 
  // test per thread struct with padding and alignment
  struct GetValuePerThreadStruct
  {
    SizeValueType st_NumberOfPixelsCounted;
    MeasureType   st_Value;
  };
  itkPadStruct( ITK_CACHE_LINE_ALIGNMENT, GetValuePerThreadStruct,
    PaddedGetValuePerThreadStruct );
  itkAlignedTypedef( ITK_CACHE_LINE_ALIGNMENT, PaddedGetValuePerThreadStruct,
    AlignedGetValuePerThreadStruct );
  mutable AlignedGetValuePerThreadStruct * m_GetValuePerThreadVariables;
  mutable ThreadIdType                     m_GetValuePerThreadVariablesSize;
 
  // test per thread struct with padding and alignment
  struct GetValueAndDerivativePerThreadStruct
  {
    SizeValueType  st_NumberOfPixelsCounted;
    MeasureType    st_Value;
    DerivativeType st_Derivative;
  };
  itkPadStruct( ITK_CACHE_LINE_ALIGNMENT, GetValueAndDerivativePerThreadStruct,
    PaddedGetValueAndDerivativePerThreadStruct );
  itkAlignedTypedef( ITK_CACHE_LINE_ALIGNMENT, PaddedGetValueAndDerivativePerThreadStruct,
    AlignedGetValueAndDerivativePerThreadStruct );
  mutable AlignedGetValueAndDerivativePerThreadStruct * m_GetValueAndDerivativePerThreadVariables;
  mutable ThreadIdType                                  m_GetValueAndDerivativePerThreadVariablesSize;
 
  virtual void InitializeThreadingParameters( void ) const;
 
  virtual void InitializeImageSampler( void );
 
  itkSetMacro( UseImageSampler, bool );
 
  virtual void CheckNumberOfSamples(
    unsigned long wanted, unsigned long found ) const;
 
  virtual void CheckForBSplineInterpolator( void );
 
  virtual bool EvaluateMovingImageValueAndDerivative(
    const MovingImagePointType & mappedPoint,
    RealType & movingImageValue,
    MovingImageDerivativeType * gradient ) const;
 
  virtual void EvaluateTransformJacobianInnerProduct(
    const TransformJacobianType & jacobian,
    const MovingImageDerivativeType & movingImageDerivative,
    DerivativeType & imageJacobian ) const;
 
  virtual void CheckForAdvancedTransform( void );
 
  virtual void CheckForBSplineTransform( void ) const;
 
  virtual bool TransformPoint(
    const FixedImagePointType & fixedImagePoint,
    MovingImagePointType & mappedPoint ) const;
 
  virtual bool EvaluateTransformJacobian(
    const FixedImagePointType & fixedImagePoint,
    TransformJacobianType & jacobian,
    NonZeroJacobianIndicesType & nzji ) const;
 
  virtual bool IsInsideMovingMask( const MovingImagePointType & point ) const;
 
  virtual void InitializeLimiters( void );
 
  itkSetMacro( UseFixedImageLimiter, bool );
  itkSetMacro( UseMovingImageLimiter, bool );
 
  double m_FixedLimitRangeRatio;
  double m_MovingLimitRangeRatio;
 
private:
 
  AdvancedImageToImageMetric( const Self & ); // purposely not implemented
  void operator=( const Self & );             // purposely not implemented
 
  bool   m_UseImageSampler;
  bool   m_UseFixedImageLimiter;
  bool   m_UseMovingImageLimiter;
  double m_RequiredRatioOfValidSamples;
  bool   m_UseMovingImageDerivativeScales;
  bool   m_ScaleGradientWithRespectToMovingImageOrientation;
 
  MovingImageDerivativeScalesType m_MovingImageDerivativeScales;
 
};
 
} // end namespace itk
 
#ifndef ITK_MANUAL_INSTANTIATION
#include "itkAdvancedImageToImageMetric.hxx"
#endif
 
#endif // end #ifndef __itkAdvancedImageToImageMetric_h