Andreas Tille · Andreas Tille · Andreas Tille · Andreas Tille · Andreas Tille · Andreas Tille
--- a/.editorconfig
+++ b/.editorconfig
+# EditorConfig: https://EditorConfig.org
+
+# This is the top-most EditorConfig file
+root = true
+
+# Apply to any source file
+[*.*]
+
+# Indentation of two spaces 
+indent_style = space
+indent_size = 2
--- a/.gitignore
+++ b/.gitignore
@@ -37,3 +37,4 @@ cmake*
 Makefile
 CMakeFiles/
 CTestTestfile.cmake
+.idea/
--- a/.travis.yml
+++ b/.travis.yml
-language: cpp
-
-env:
-  global:
-    ELASTIX_SOURCE_DIR=${TRAVIS_BUILD_DIR}
-    ELASTIX_BUILD_DIR=${HOME}/elastix-build
-    ITK_URL=https://github.com/InsightSoftwareConsortium/ITK
-    ITK_SOURCE_DIR=${HOME}/ITK-source
-    ITK_BUILD_DIR=${HOME}/ITK-build
-    OMP_NUM_THREADS=2
-
-cache:
-  directories:
-    - ${ITK_SOURCE_DIR}
-    - ${ITK_BUILD_DIR}
-    
-matrix:
-  include:
-    - os: linux
-      dist: trusty
-      sudo: required
-      compiler: clang
-    - os: linux
-      dist: trusty
-      sudo: required
-      compiler: gcc
-    - os: osx
-      compiler: clang
-    - os: osx
-      compiler: gcc
-
-before_script:
-  - mkdir -p ${ELASTIX_BUILD_DIR}
-
-script:
-  - if [[ ! -e ${ITK_BUILD_DIR}/ITKConfig.cmake ]]; then
-      rm -rf ${ITK_SOURCE_DIR} {ITK_BUILD_DIR} &&
-      git clone ${ITK_URL} ${ITK_SOURCE_DIR} &&
-      cd ${ITK_SOURCE_DIR} &&
-      git checkout v4.11.1 &&
-      mkdir -p ${ITK_BUILD_DIR} &&
-      cd ${ITK_BUILD_DIR} &&
-      cmake -DBUILD_EXAMPLES=OFF -DBUILD_TESTING:BOOL=OFF -DITK_BUILD_DEFAULT_MODULES:BOOL=OFF -DITKGroup_IO:BOOL=ON -DModule_ITKDistanceMap:BOOL=ON -DModule_ITKImageFusion:BOOL=ON -DModule_ITKLabelMap:BOOL:=ON -DModule_ITKMathematicalMorphology:BOOL=ON -DModule_ITKOptimizers:BOOL=ON -DModule_ITKOptimizersv4:BOOL=ON -DModule_ITKRegistrationCommon:BOOL=ON -DModule_ITKVideoIO:BOOL=ON-DCMAKE_BUILD_TYPE=Release ${ITK_SOURCE_DIR} &&
-      make --jobs=4 && 
-      touch ${HOME}/built_cache; fi
-  - if [[ ! -e ${HOME}/built_cache ]]; then
-      cd ${ELASTIX_BUILD_DIR} &&
-      cmake -DELASTIX_USE_OPENMP=OFF -DITK_DIR=${ITK_BUILD_DIR} -DCMAKE_BUILD_TYPE=Release ${ELASTIX_SOURCE_DIR} &&
-      make --jobs=4 && 
-      ctest --jobs=4 --verbose; fi
--- a/CMake/FindSuiteSparse.cmake
+++ b/CMake/FindSuiteSparse.cmake
--- a/CMake/elastixExportTarget.cmake
+++ b/CMake/elastixExportTarget.cmake
@@ -2,7 +2,7 @@ function(elastix_export_target tgt)
  # Remove the build tree's ElastixTargets file if this is the first call:
  get_property(first_time GLOBAL PROPERTY ELASTIX_FIRST_EXPORTED_TARGET)
  if(NOT first_time)
-    file(REMOVE ${CMAKE_BINARY_DIR}/ElastixTargets.cmake)
+    file(REMOVE ${elastix_BINARY_DIR}/ElastixTargets.cmake)
    set_property(GLOBAL PROPERTY ELASTIX_FIRST_EXPORTED_TARGET 1)
  endif()

@@ -17,6 +17,6 @@ function(elastix_export_target tgt)
  endif()

  export(TARGETS ${tgt}
-    APPEND FILE "${CMAKE_BINARY_DIR}/ElastixTargets.cmake"
+    APPEND FILE "${elastix_BINARY_DIR}/ElastixTargets.cmake"
  )
-endfunction()
\ No newline at end of file
+endfunction()
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
 cmake_minimum_required( VERSION 2.8 )

 project( elastix )
+set(CMAKE_CXX_STANDARD 11)
+
+if( BUILD_SHARED_LIBS )
+  message(FATAL_ERROR "Elastix does not support BUILD_SHARED_LIBS")
+endif()

 #---------------------------------------------------------------------
 cmake_policy( SET CMP0012 NEW )
@@ -61,7 +66,7 @@ endif()

 #---------------------------------------------------------------------
 # Find ITK.
-find_package( ITK REQUIRED )
+find_package( ITK 5.0 REQUIRED )
 include( ${ITK_USE_FILE} )

 #---------------------------------------------------------------------
@@ -135,11 +140,32 @@ endif()
 # before declaring any targets.

 if( USE_PreconditionedGradientDescent )
- set( SUITESPARSE_BINARY_DIR "" CACHE PATH "Path to SuiteSparse" )
- set( SUITESPARSE_USE_FILE ${SUITESPARSE_BINARY_DIR}/UseSuiteSparse.cmake )
- if( EXISTS ${SUITESPARSE_USE_FILE} )
-  include( ${SUITESPARSE_USE_FILE} )
- endif()
+# ------------------------------------------------------------------
+## CMake file to locate SuiteSparse and its useful composite projects
+## The first developpement of this file was made fro Windows users who
+## use:
+##  https://github.com/jlblancoc/suitesparse-metis-for-windows
+## Anyway, it chould be work also on linux (tested on fedora 17 when you installed suitesparse from yum)
+# Detect SuiteSparse libraries:
+# If not found automatically, set SuiteSparse_DIR in CMake to the
+# directory where SuiteSparse was built.
+# ------------------------------------------------------------------
+  list( APPEND CMAKE_MODULE_PATH "${elastix_SOURCE_DIR}/.." ) # Add the directory where FindSuiteSparse.cmake module can be found.
+
+  set( SuiteSparse_USE_LAPACK_BLAS ON )
+  find_package( SuiteSparse QUIET NO_MODULE )  # 1st: Try to locate the *config.cmake file.
+  if( NOT SuiteSparse_FOUND )
+    set( SuiteSparse_VERBOSE ON )
+    find_package( SuiteSparse REQUIRED )       # 2nd: Use FindSuiteSparse.cmake module
+    include_directories( ${SuiteSparse_INCLUDE_DIRS} )
+  else()
+    message( STATUS "Find SuiteSparse : include(${USE_SuiteSparse})" )
+    include( ${USE_SuiteSparse} )
+  endif()
+  message( STATUS "SuiteSparse_LIBS: ${SuiteSparse_LIBRARIES}" )
+# ------------------------------------------------------------------
+#   End of SuiteSparse detection
+# ------------------------------------------------------------------
 endif()

 #---------------------------------------------------------------------
@@ -389,6 +415,11 @@ if( ${BUILD_DOCUMENTATION} )
      endif()
    endif()

+    # Configure the doxygen configuration
+    configure_file(
+      ${ELASTIX_DOXYGEN_DIR}/doxyfile.in
+      ${ELASTIX_HELP_DIR}/doxyfile.out @ONLY )
+
    # Configure the footer of the help website.
    configure_file(
      ${ELASTIX_DOXYGEN_DIR}/DoxygenFooter.html.in
@@ -408,6 +439,14 @@ if( BUILD_TESTING )
  add_subdirectory( Testing )
 endif()

+mark_as_advanced( ELASTIX_USE_GTEST )
+option( ELASTIX_USE_GTEST "Use GoogleTest to test Elastix implementation" OFF )
+
+if( ELASTIX_USE_GTEST )
+  enable_testing()
+  add_subdirectory( Common/GTesting )
+endif()
+
 #---------------------------------------------------------------------
 # Packaging

@@ -475,7 +514,7 @@ if( ELASTIX_ENABLE_PACKAGER )
  #set( CPACK_PACKAGE_DESCRIPTION_FILE
  #  "${CMAKE_CURRENT_SOURCE_DIR}/ReadMe.txt" )
  set( CPACK_RESOURCE_FILE_LICENSE
-	  "${CMAKE_CURRENT_SOURCE_DIR}/NOTICE" )
+    "${CMAKE_CURRENT_SOURCE_DIR}/NOTICE" )

  # The default install directories: .../elastix_v4.3
  set( CPACK_PACKAGE_INSTALL_DIRECTORY
@@ -563,6 +602,6 @@ foreach( LIB IN LISTS AllComponentLibs )
  elastix_export_target( ${LIB} )
 endforeach()

-configure_file( ${CMAKE_SOURCE_DIR}/ElastixConfig.cmake.in ${CMAKE_BINARY_DIR}/ElastixConfig.cmake @ONLY )
-configure_file( ${CMAKE_SOURCE_DIR}/ElastixConfigVersion.cmake.in ${CMAKE_BINARY_DIR}/ElastixConfigVersion.cmake @ONLY )
-configure_file( ${CMAKE_SOURCE_DIR}/UseElastix.cmake.in ${CMAKE_BINARY_DIR}/UseElastix.cmake @ONLY )
+configure_file( ${elastix_SOURCE_DIR}/ElastixConfig.cmake.in ${elastix_BINARY_DIR}/ElastixConfig.cmake @ONLY )
+configure_file( ${elastix_SOURCE_DIR}/ElastixConfigVersion.cmake.in ${elastix_BINARY_DIR}/ElastixConfigVersion.cmake @ONLY )
+configure_file( ${elastix_SOURCE_DIR}/UseElastix.cmake.in ${elastix_BINARY_DIR}/UseElastix.cmake @ONLY )
--- a/CONTRIBUTING.md
+++ b/CONTRIBUTING.md
@@ -5,13 +5,13 @@

 ### Do you have questions about the source code? ###

-* Ask any question about how to use `elastix` on the _mailing list_. You can subscribe to the mailing list [here](http://lists.bigr.nl/mailman/listinfo/elastix), after which you can email to elastix@bigr.nl.
+* Ask any question about how to use `elastix` on the [mailing list](https://groups.google.com/forum/#!forum/elastix-imageregistration).

 * Do not open an issue on GitHub for general questions, registration questions, or issues you may have while running `elastix`. _GitHub issues are primarily intended for bug reports and fixes._

-* General information about `elastix` can be found on our [website](http://elastix.isi.uu.nl/).
+* General information about `elastix` can be found on our [wiki](https://github.com/SuperElastix/elastix/wiki) or at the [website](http://elastix.isi.uu.nl/).

-* [The manual](http://elastix.isi.uu.nl/download/elastix_manual_v4.8.pdf) explains much of the inner workings of image registration.
+* [The manual](https://github.com/SuperElastix/elastix/releases/download/4.9.0/elastix-4.9.0-manual.pdf) explains much of the inner workings of image registration.

 ### Did you find a bug? ###

@@ -31,9 +31,9 @@

 * Make sure the PR is done with respect to the [develop branch](https://github.com/SuperElastix/elastix/tree/develop).

-* Ensure the PR description (log message) _clearly describes the problem and solution_. Include the relevant issue number if applicable. One-line messages are fine for small changes, but bigger changes should look like this:  
+* Ensure the PR description (log message) _clearly describes the problem and solution_. Include the relevant issue number if applicable. One-line messages are fine for small changes, but bigger changes should look like this:
    $ git commit -m "ENH: A brief summary of the commit
-    > 
+    >
    > A paragraph describing what changed and its impact."

 * We use the following tags for commit messages:
@@ -60,7 +60,7 @@

 Start reading our code and you'll get the hang of it. We optimize for readability:

-* We indent using two spaces (soft tabs). This ensures a consistent layout of the code on different text editors. 
+* We indent using two spaces (soft tabs). This ensures a consistent layout of the code on different text editors.

 * We ALWAYS put spaces after list items and function arguments (`[1, 2, 3]`, not `[1,2,3]`), around operators (`x += 1`, not `x+=1`), etc.

@@ -68,9 +68,9 @@ Start reading our code and you'll get the hang of it. We optimize for readabilit

 * Member variables of classes should start their name with `m_`. For example: `m_NumberOfIterations`.

-* Type definitions should start with a capital. ImageType for example, instead of imageType. 
+* Type definitions should start with a capital. ImageType for example, instead of imageType.

 * This is open source software. Consider the people who will read your code, and make it look nice for them. It's sort of like driving a car: Perhaps you love doing donuts when you're alone, but with passengers the goal is to make the ride as smooth as possible.

-Thanks,  
+Thanks,
 The `elastix` team
--- a/Common/CMakeLists.txt
+++ b/Common/CMakeLists.txt
@@ -20,20 +20,20 @@ set( CommonFiles
  itkAdvancedLinearInterpolateImageFunction.hxx
  itkAdvancedRayCastInterpolateImageFunction.h
  itkAdvancedRayCastInterpolateImageFunction.hxx
+  itkComputeImageExtremaFilter.h
+  itkComputeImageExtremaFilter.hxx
  itkComputeDisplacementDistribution.h
  itkComputeDisplacementDistribution.hxx
  itkComputeJacobianTerms.h
  itkComputeJacobianTerms.hxx
+  itkComputePreconditionerUsingDisplacementDistribution.h
+  itkComputePreconditionerUsingDisplacementDistribution.hxx
  itkErodeMaskImageFilter.h
  itkErodeMaskImageFilter.hxx
  itkGenericMultiResolutionPyramidImageFilter.h
  itkGenericMultiResolutionPyramidImageFilter.hxx
  itkImageFileCastWriter.h
  itkImageFileCastWriter.hxx
-  itkImageMaskSpatialObject2.h
-  itkImageMaskSpatialObject2.hxx
-  itkImageSpatialObject2.h
-  itkImageSpatialObject2.hxx
  itkMeshFileReaderBase.h
  itkMeshFileReaderBase.hxx
  itkMultiOrderBSplineDecompositionImageFilter.h
@@ -94,6 +94,8 @@ set( TransformFiles
  Transforms/itkAdvancedEuler3DTransform.h
  Transforms/itkAdvancedEuler3DTransform.hxx
  Transforms/itkAdvancedIdentityTransform.h
+  Transforms/itkAdvancedImageMomentsCalculator.h
+  Transforms/itkAdvancedImageMomentsCalculator.hxx
  Transforms/itkAdvancedMatrixOffsetTransformBase.h
  Transforms/itkAdvancedMatrixOffsetTransformBase.hxx
  Transforms/itkAdvancedRigid2DTransform.h

--- a/Common/CostFunctions/itkAdvancedImageToImageMetric.h
+++ b/Common/CostFunctions/itkAdvancedImageToImageMetric.h
@@ -30,11 +30,13 @@
 #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 "itkMultiThreader.h"
+#include "itkPlatformMultiThreader.h"

 namespace itk
 {
@@ -129,6 +131,9 @@ public:
  typedef typename DerivativeType::ValueType                DerivativeValueType;
  typedef typename Superclass::ParametersType               ParametersType;

+  typedef ImageMaskSpatialObject< itkGetStaticConstMacro( FixedImageDimension ) > FixedImageMaskSpatialObject2Type;
+  typedef ImageMaskSpatialObject< itkGetStaticConstMacro( MovingImageDimension ) > MovingImageMaskSpatialObject2Type;
+
  /** Some useful extra typedefs. */
  typedef typename FixedImageType::PixelType               FixedImagePixelType;
  typedef typename MovingImageType::RegionType             MovingImageRegionType;
@@ -168,8 +173,8 @@ public:
  typedef vnl_sparse_matrix< HessianValueType > HessianType;

  /** Typedefs for multi-threading. */
-  typedef itk::MultiThreader                      ThreaderType;
-  typedef typename ThreaderType::ThreadInfoStruct ThreadInfoType;
+  typedef itk::PlatformMultiThreader                      ThreaderType;
+  typedef typename ThreaderType::WorkUnitInfo ThreadInfoType;

  /** Public methods ********************/

@@ -186,7 +191,7 @@ public:


  /** Get the advanced transform. */
-  const AdvancedTransformType * GetTransform( void ) const ITK_OVERRIDE
+  const AdvancedTransformType * GetTransform( void ) const override
  {
    return this->m_AdvancedTransform.GetPointer();
  }
@@ -257,7 +262,7 @@ public:
   * \li Check if a B-spline interpolator has been set
   * \li Check if an AdvancedTransform has been set
   */
-  virtual void Initialize( void ) throw ( ExceptionObject ) ITK_OVERRIDE;
+  void Initialize( void ) override;

  /** Experimental feature: compute SelfHessian.
   * This base class just returns an identity matrix of the right size.
@@ -265,7 +270,7 @@ public:
  virtual void GetSelfHessian( const TransformParametersType & parameters, HessianType & H ) const;

  /** Set number of threads to use for computations. */
-  virtual void SetNumberOfThreads( ThreadIdType numberOfThreads );
+  virtual void SetNumberOfWorkUnits( ThreadIdType numberOfThreads );

  /** Switch the function BeforeThreadedGetValueAndDerivative on or off. */
  itkSetMacro( UseMetricSingleThreaded, bool );
@@ -292,10 +297,10 @@ protected:
  AdvancedImageToImageMetric();

  /** Destructor. */
-  virtual ~AdvancedImageToImageMetric();
+  ~AdvancedImageToImageMetric() override;

  /** PrintSelf. */
-  void PrintSelf( std::ostream & os, Indent indent ) const;
+  void PrintSelf( std::ostream & os, Indent indent ) const override;

  /** Protected Typedefs ******************/

@@ -458,7 +463,7 @@ protected:
  /** Methods for image sampler support **********/

  /** Initialize variables related to the image sampler; called by Initialize. */
-  virtual void InitializeImageSampler( void ) throw ( ExceptionObject );
+  virtual void InitializeImageSampler( void );

  /** Inheriting classes can specify whether they use the image sampler functionality
   * Make sure to set it before calling Initialize; default: false. */
@@ -530,22 +535,6 @@ protected:
  /** Convenience method: check if point is inside the moving mask. *****************/
  virtual bool IsInsideMovingMask( const MovingImagePointType & point ) const;

-  /** Methods for the support of gray value limiters. ***************/
-
-  /** Compute the extrema of fixed image over a region
-   * Initializes the m_Fixed[True]{Max,Min}[Limit]
-   * This method is called by InitializeLimiters() and uses the FixedLimitRangeRatio */
-  virtual void ComputeFixedImageExtrema(
-    const FixedImageType * image,
-    const FixedImageRegionType & region );
-
-  /** Compute the extrema of the moving image over a region
-   * Initializes the m_Moving[True]{Max,Min}[Limit]
-   * This method is called by InitializeLimiters() and uses the MovingLimitRangeRatio; */
-  virtual void ComputeMovingImageExtrema(
-    const MovingImageType * image,
-    const MovingImageRegionType & region );
-
  /** Initialize the {Fixed,Moving}[True]{Max,Min}[Limit] and the {Fixed,Moving}ImageLimiter
   * Only does something when Use{Fixed,Moving}Limiter is set to true; */
  virtual void InitializeLimiters( void );
@@ -555,6 +544,9 @@ protected:
  itkSetMacro( UseFixedImageLimiter, bool );
  itkSetMacro( UseMovingImageLimiter, bool );

+  double m_FixedLimitRangeRatio;
+  double m_MovingLimitRangeRatio;
+
 private:

  AdvancedImageToImageMetric( const Self & ); // purposely not implemented
@@ -562,8 +554,6 @@ private:

  /** Private member variables. */
  bool   m_UseImageSampler;
-  double m_FixedLimitRangeRatio;
-  double m_MovingLimitRangeRatio;
  bool   m_UseFixedImageLimiter;
  bool   m_UseMovingImageLimiter;
  double m_RequiredRatioOfValidSamples;

--- a/Common/CostFunctions/itkAdvancedImageToImageMetric.hxx
+++ b/Common/CostFunctions/itkAdvancedImageToImageMetric.hxx
@@ -20,14 +20,15 @@

 #include "itkAdvancedImageToImageMetric.h"

-#include "itkImageRegionConstIterator.h"          // used for extrema computation
-#include "itkImageRegionConstIteratorWithIndex.h" // used for extrema computation
 #include "itkAdvancedRayCastInterpolateImageFunction.h"
+#include "itkComputeImageExtremaFilter.h"

 #ifdef ELASTIX_USE_OPENMP
 #include <omp.h>
 #endif

+#include "itkTimeProbe.h"
+
 namespace itk
 {

@@ -86,14 +87,12 @@ AdvancedImageToImageMetric< TFixedImage, TMovingImage >
  /** Threading related variables. */
  this->m_UseMetricSingleThreaded = true;
  this->m_UseMultiThread = false;
-  this->m_Threader->SetUseThreadPool( false ); // setting to true makes elastix hang
-                                               // at a WaitForSingleMethodThread()

  /** OpenMP related. Switch to on when available */
 #ifdef ELASTIX_USE_OPENMP
  this->m_UseOpenMP = true;

-  const int nthreads = static_cast< int >( this->m_NumberOfThreads );
+  const int nthreads = static_cast< int >( Self::GetNumberOfWorkUnits() );
  omp_set_num_threads( nthreads );
 #else
  this->m_UseOpenMP = false;
@@ -103,9 +102,9 @@ AdvancedImageToImageMetric< TFixedImage, TMovingImage >
  this->m_ThreaderMetricParameters.st_Metric = this;

  // Multi-threading structs
-  this->m_GetValuePerThreadVariables                  = NULL;
+  this->m_GetValuePerThreadVariables                  = nullptr;
  this->m_GetValuePerThreadVariablesSize              = 0;
-  this->m_GetValueAndDerivativePerThreadVariables     = NULL;
+  this->m_GetValueAndDerivativePerThreadVariables     = nullptr;
  this->m_GetValueAndDerivativePerThreadVariablesSize = 0;

 } // end Constructor
@@ -125,21 +124,23 @@ AdvancedImageToImageMetric< TFixedImage, TMovingImage >


 /**
- * ********************* SetNumberOfThreads ****************************
+ * ********************* SetNumberOfWorkUnits ****************************
 */

 template< class TFixedImage, class TMovingImage >
 void
 AdvancedImageToImageMetric< TFixedImage, TMovingImage >
-::SetNumberOfThreads( ThreadIdType numberOfThreads )
+::SetNumberOfWorkUnits( ThreadIdType numberOfThreads )
 {
-  Superclass::SetNumberOfThreads( numberOfThreads );
+  // Note: This is a workaround for ITK5, which renamed NumberOfThreads
+  // to NumberOfWorkUnits
+  Superclass::SetNumberOfWorkUnits( numberOfThreads );

 #ifdef ELASTIX_USE_OPENMP
-  const int nthreads = static_cast< int >( this->m_NumberOfThreads );
+  const int nthreads = static_cast< int >( Self::GetNumberOfWorkUnits() );
  omp_set_num_threads( nthreads );
 #endif
-} // end SetNumberOfThreads()
+} // end SetNumberOfWorkUnits()


 /**
@@ -149,7 +150,7 @@ AdvancedImageToImageMetric< TFixedImage, TMovingImage >
 template< class TFixedImage, class TMovingImage >
 void
 AdvancedImageToImageMetric< TFixedImage, TMovingImage >
-::Initialize( void ) throw ( ExceptionObject )
+::Initialize( void )
 {
  /** Initialize transform, interpolator, etc. */
  Superclass::Initialize();
@@ -187,6 +188,8 @@ void
 AdvancedImageToImageMetric< TFixedImage, TMovingImage >
 ::InitializeThreadingParameters( void ) const
 {
+  const ThreadIdType numberOfThreads = Self::GetNumberOfWorkUnits();
+
  /** Resize and initialize the threading related parameters.
   * The SetSize() functions do not resize the data when this is not
   * needed, which saves valuable re-allocation time.
@@ -198,23 +201,23 @@ AdvancedImageToImageMetric< TFixedImage, TMovingImage >
   */

  /** Only resize the array of structs when needed. */
-  if( this->m_GetValuePerThreadVariablesSize != this->m_NumberOfThreads )
+  if( this->m_GetValuePerThreadVariablesSize != numberOfThreads )
  {
    delete[] this->m_GetValuePerThreadVariables;
-    this->m_GetValuePerThreadVariables     = new AlignedGetValuePerThreadStruct[ this->m_NumberOfThreads ];
-    this->m_GetValuePerThreadVariablesSize = this->m_NumberOfThreads;
+    this->m_GetValuePerThreadVariables     = new AlignedGetValuePerThreadStruct[ numberOfThreads ];
+    this->m_GetValuePerThreadVariablesSize = numberOfThreads;
  }

  /** Only resize the array of structs when needed. */
-  if( this->m_GetValueAndDerivativePerThreadVariablesSize != this->m_NumberOfThreads )
+  if( this->m_GetValueAndDerivativePerThreadVariablesSize != numberOfThreads )
  {
    delete[] this->m_GetValueAndDerivativePerThreadVariables;
-    this->m_GetValueAndDerivativePerThreadVariables     = new AlignedGetValueAndDerivativePerThreadStruct[ this->m_NumberOfThreads ];
-    this->m_GetValueAndDerivativePerThreadVariablesSize = this->m_NumberOfThreads;
+    this->m_GetValueAndDerivativePerThreadVariables     = new AlignedGetValueAndDerivativePerThreadStruct[ numberOfThreads ];
+    this->m_GetValueAndDerivativePerThreadVariablesSize = numberOfThreads;
  }

  /** Some initialization. */
-  for( ThreadIdType i = 0; i < this->m_NumberOfThreads; ++i )
+  for( ThreadIdType i = 0; i < numberOfThreads; ++i )
  {
    this->m_GetValuePerThreadVariables[ i ].st_NumberOfPixelsCounted = NumericTraits< SizeValueType >::Zero;
    this->m_GetValuePerThreadVariables[ i ].st_Value                 = NumericTraits< MeasureType >::Zero;
@@ -229,152 +232,58 @@ AdvancedImageToImageMetric< TFixedImage, TMovingImage >


 /**
- * ****************** ComputeFixedImageExtrema ***************************
+ * ****************** InitializeLimiters *****************************
 */

 template< class TFixedImage, class TMovingImage >
 void
 AdvancedImageToImageMetric< TFixedImage, TMovingImage >
-::ComputeFixedImageExtrema(
-  const FixedImageType * image,
-  const FixedImageRegionType & region )
+::InitializeLimiters( void )
 {
-  /** NB: We can't use StatisticsImageFilterWithMask to do this because
-   * the filter computes the min/max for the largest possible region.
-   * This filter is multi-threaded though.
-   */
-  FixedImagePixelType trueMinTemp = NumericTraits< FixedImagePixelType >::max();
-  FixedImagePixelType trueMaxTemp = NumericTraits< FixedImagePixelType >::NonpositiveMin();
-
-  /** If no mask. */
-  if( this->m_FixedImageMask.IsNull() )
-  {
-    typedef ImageRegionConstIterator< FixedImageType > IteratorType;
-    IteratorType it( image, region );
-    for( it.GoToBegin(); !it.IsAtEnd(); ++it )
-    {
-      const FixedImagePixelType sample = it.Get();
-      trueMinTemp = vnl_math_min( trueMinTemp, sample );
-      trueMaxTemp = vnl_math_max( trueMaxTemp, sample );
-    }
-  }
-  /** Excluded extrema outside the mask.
-   * Because we have to call TransformIndexToPhysicalPoint() and
-   * check IsInside() this way is much (!) slower.
-   */
-  else
+  /** Set up fixed limiter. */
+  if( this->GetUseFixedImageLimiter() )
  {
-    typedef ImageRegionConstIteratorWithIndex< FixedImageType > IteratorType;
-    IteratorType it( image, region );
-
-    for( it.GoToBegin(); !it.IsAtEnd(); ++it )
+    if( this->GetFixedImageLimiter() == 0 )
    {
-      OutputPointType point;
-      image->TransformIndexToPhysicalPoint( it.GetIndex(), point );
-      if( this->m_FixedImageMask->IsInside( point ) )
-      {
-        const FixedImagePixelType sample = it.Get();
-        trueMinTemp = vnl_math_min( trueMinTemp, sample );
-        trueMaxTemp = vnl_math_max( trueMaxTemp, sample );
-      }
+      itkExceptionMacro( << "No fixed image limiter has been set!" );
    }
-  }
-
-  /** Update member variables. */
-  this->m_FixedImageTrueMin = trueMinTemp;
-  this->m_FixedImageTrueMax = trueMaxTemp;
-
-  this->m_FixedImageMinLimit = static_cast< FixedImageLimiterOutputType >(
-    trueMinTemp - this->m_FixedLimitRangeRatio * ( trueMaxTemp - trueMinTemp ) );
-  this->m_FixedImageMaxLimit = static_cast< FixedImageLimiterOutputType >(
-    trueMaxTemp + this->m_FixedLimitRangeRatio * ( trueMaxTemp - trueMinTemp ) );
-
-} // end ComputeFixedImageExtrema()
-

-/**
- * ****************** ComputeMovingImageExtrema ***************************
- */
-
-template< class TFixedImage, class TMovingImage >
-void
-AdvancedImageToImageMetric< TFixedImage, TMovingImage >
-::ComputeMovingImageExtrema(
-  const MovingImageType * image,
-  const MovingImageRegionType & region )
-{
-  /** NB: We can't use StatisticsImageFilter to do this because
-   * the filter computes the min/max for the largest possible region.
-   */
-  MovingImagePixelType trueMinTemp = NumericTraits< MovingImagePixelType >::max();
-  MovingImagePixelType trueMaxTemp = NumericTraits< MovingImagePixelType >::NonpositiveMin();
+    itk::TimeProbe timer;
+    timer.Start();

-  /** If no mask. */
-  if( this->m_MovingImageMask.IsNull() )
-  {
-    typedef ImageRegionConstIterator< MovingImageType > IteratorType;
-    IteratorType iterator( image, region );
-    for( iterator.GoToBegin(); !iterator.IsAtEnd(); ++iterator )
+    typedef typename itk::ComputeImageExtremaFilter<FixedImageType> ComputeFixedImageExtremaFilterType;
+    typename ComputeFixedImageExtremaFilterType::Pointer computeFixedImageExtrema
+      = ComputeFixedImageExtremaFilterType::New();
+    computeFixedImageExtrema->SetInput( this->GetFixedImage() );
+    computeFixedImageExtrema->SetImageRegion( this->GetFixedImageRegion() );
+    if( this->m_FixedImageMask.IsNotNull() )
    {
-      const MovingImagePixelType sample = iterator.Get();
-      trueMinTemp = vnl_math_min( trueMinTemp, sample );
-      trueMaxTemp = vnl_math_max( trueMaxTemp, sample );
-    }
-  }
-  /** Excluded extrema outside the mask.
-   * Because we have to call TransformIndexToPhysicalPoint() and
-   * check IsInside() this way is much (!) slower.
-   */
-  else
-  {
-    typedef ImageRegionConstIteratorWithIndex< MovingImageType > IteratorType;
-    IteratorType it( image, region );
+      computeFixedImageExtrema->SetUseMask( true );

-    for( it.GoToBegin(); !it.IsAtEnd(); ++it )
-    {
-      OutputPointType point;
-      image->TransformIndexToPhysicalPoint( it.GetIndex(), point );
-      if( this->m_MovingImageMask->IsInside( point ) )
+      const FixedImageMaskSpatialObject2Type * fMask
+        = dynamic_cast< const FixedImageMaskSpatialObject2Type * >( this->m_FixedImageMask.GetPointer() );
+      if( fMask )
+      {
+        computeFixedImageExtrema->SetImageSpatialMask( fMask );
+      }
+      else
      {
-        const MovingImagePixelType sample = it.Get();
-        trueMinTemp = vnl_math_min( trueMinTemp, sample );
-        trueMaxTemp = vnl_math_max( trueMaxTemp, sample );
+        computeFixedImageExtrema->SetImageMask( this->GetFixedImageMask() );
      }
    }
-  }
-
-  /** Update member variables. */
-  this->m_MovingImageTrueMin = trueMinTemp;
-  this->m_MovingImageTrueMax = trueMaxTemp;

-  this->m_MovingImageMinLimit = static_cast< MovingImageLimiterOutputType >(
-    trueMinTemp - this->m_MovingLimitRangeRatio * ( trueMaxTemp - trueMinTemp ) );
-  this->m_MovingImageMaxLimit = static_cast< MovingImageLimiterOutputType >(
-    trueMaxTemp + this->m_MovingLimitRangeRatio * ( trueMaxTemp - trueMinTemp ) );
+    computeFixedImageExtrema->Update();
+    timer.Stop();
+    elxout << "  Computing the fixed image extrema took "
+      << static_cast< long >( timer.GetMean() * 1000 ) << " ms." << std::endl;

-} // end ComputeMovingImageExtrema()
+    this->m_FixedImageTrueMax = computeFixedImageExtrema->GetMaximum();
+    this->m_FixedImageTrueMin = computeFixedImageExtrema->GetMinimum();

-
-/**
- * ****************** InitializeLimiter *****************************
- */
-
-template< class TFixedImage, class TMovingImage >
-void
-AdvancedImageToImageMetric< TFixedImage, TMovingImage >
-::InitializeLimiters( void )
-{
-  /** Set up fixed limiter. */
-  if( this->GetUseFixedImageLimiter() )
-  {
-    if( this->GetFixedImageLimiter() == 0 )
-    {
-      itkExceptionMacro( << "No fixed image limiter has been set!" );
-    }
-
-    this->ComputeFixedImageExtrema(
-      this->GetFixedImage(),
-      this->GetFixedImageRegion() );
+    this->m_FixedImageMinLimit = static_cast< FixedImageLimiterOutputType >(
+      this->m_FixedImageTrueMin - this->m_FixedLimitRangeRatio * ( this->m_FixedImageTrueMax - this->m_FixedImageTrueMin ) );
+    this->m_FixedImageMaxLimit = static_cast< FixedImageLimiterOutputType >(
+      this->m_FixedImageTrueMax + this->m_FixedLimitRangeRatio * ( this->m_FixedImageTrueMax - this->m_FixedImageTrueMin ) );

    this->m_FixedImageLimiter->SetLowerThreshold(
      static_cast< RealType >( this->m_FixedImageTrueMin ) );
@@ -394,9 +303,41 @@ AdvancedImageToImageMetric< TFixedImage, TMovingImage >
      itkExceptionMacro( << "No moving image limiter has been set!" );
    }

-    this->ComputeMovingImageExtrema(
-      this->GetMovingImage(),
-      this->GetMovingImage()->GetBufferedRegion() );
+    itk::TimeProbe timer;
+    timer.Start();
+
+    typedef typename itk::ComputeImageExtremaFilter<MovingImageType> ComputeMovingImageExtremaFilterType;
+    typename ComputeMovingImageExtremaFilterType::Pointer computeMovingImageExtrema
+      = ComputeMovingImageExtremaFilterType::New();
+    computeMovingImageExtrema->SetInput( this->GetMovingImage() );
+    computeMovingImageExtrema->SetImageRegion( this->GetMovingImage()->GetBufferedRegion() );
+    if( this->m_MovingImageMask.IsNotNull() )
+    {
+      computeMovingImageExtrema->SetUseMask( true );
+      const MovingImageMaskSpatialObject2Type * mMask
+        = dynamic_cast< const MovingImageMaskSpatialObject2Type * >( this->m_MovingImageMask.GetPointer() );
+      if( mMask )
+      {
+        computeMovingImageExtrema->SetImageSpatialMask( mMask );
+      }
+      else
+      {
+        computeMovingImageExtrema->SetImageMask( this->GetMovingImageMask() );
+      }
+    }
+    computeMovingImageExtrema->Update();
+
+    timer.Stop();
+    elxout << "  Computing the moving image extrema took "
+      << static_cast< long >( timer.GetMean() * 1000 ) << " ms." << std::endl;
+
+    this->m_MovingImageTrueMax = computeMovingImageExtrema->GetMaximum();
+    this->m_MovingImageTrueMin = computeMovingImageExtrema->GetMinimum();
+
+    this->m_MovingImageMinLimit = static_cast< MovingImageLimiterOutputType >(
+      this->m_MovingImageTrueMin - this->m_MovingLimitRangeRatio * ( this->m_MovingImageTrueMax - this->m_MovingImageTrueMin ) );
+    this->m_MovingImageMaxLimit = static_cast< MovingImageLimiterOutputType >(
+      this->m_MovingImageTrueMax + this->m_MovingLimitRangeRatio * ( this->m_MovingImageTrueMax - this->m_MovingImageTrueMin ) );

    this->m_MovingImageLimiter->SetLowerThreshold(
      static_cast< RealType >( this->m_MovingImageTrueMin ) );
@@ -408,7 +349,7 @@ AdvancedImageToImageMetric< TFixedImage, TMovingImage >
    this->m_MovingImageLimiter->Initialize();
  }

-} // end InitializeLimiter()
+} // end InitializeLimiters()


 /**
@@ -418,7 +359,7 @@ AdvancedImageToImageMetric< TFixedImage, TMovingImage >
 template< class TFixedImage, class TMovingImage >
 void
 AdvancedImageToImageMetric< TFixedImage, TMovingImage >
-::InitializeImageSampler( void ) throw ( ExceptionObject )
+::InitializeImageSampler( void )
 {
  if( this->GetUseImageSampler() )
  {
@@ -615,11 +556,11 @@ AdvancedImageToImageMetric< TFixedImage, TMovingImage >
  else if( testPtr_combo )
  {
    /** Check if the current transform is a B-spline transform. */
-    BSplineOrder1TransformType * testPtr_1b = dynamic_cast< BSplineOrder1TransformType * >(
+    const BSplineOrder1TransformType * testPtr_1b = dynamic_cast< const BSplineOrder1TransformType * >(
      testPtr_combo->GetCurrentTransform() );
-    BSplineOrder2TransformType * testPtr_2b = dynamic_cast< BSplineOrder2TransformType * >(
+    const BSplineOrder2TransformType * testPtr_2b = dynamic_cast< const BSplineOrder2TransformType * >(
      testPtr_combo->GetCurrentTransform() );
-    BSplineOrder3TransformType * testPtr_3b = dynamic_cast< BSplineOrder3TransformType * >(
+    const BSplineOrder3TransformType * testPtr_3b = dynamic_cast< const BSplineOrder3TransformType * >(
      testPtr_combo->GetCurrentTransform() );
    if( testPtr_1b || testPtr_2b || testPtr_3b )
    {
@@ -863,7 +804,7 @@ AdvancedImageToImageMetric< TFixedImage, TMovingImage >
  /** If a mask has been set: */
  if( this->m_MovingImageMask.IsNotNull() )
  {
-    return this->m_MovingImageMask->IsInside( point );
+    return this->m_MovingImageMask->IsInsideInWorldSpace( point );
  }

  /** If no mask has been set, just return true. */
@@ -930,14 +871,14 @@ AdvancedImageToImageMetric< TFixedImage, TMovingImage >
 ::GetValueThreaderCallback( void * arg )
 {
  ThreadInfoType * infoStruct = static_cast< ThreadInfoType * >( arg );
-  ThreadIdType     threadID   = infoStruct->ThreadID;
+  ThreadIdType     threadID   = infoStruct->WorkUnitID;

  MultiThreaderParameterType * temp
    = static_cast< MultiThreaderParameterType * >( infoStruct->UserData );

  temp->st_Metric->ThreadedGetValue( threadID );

-  return ITK_THREAD_RETURN_VALUE;
+  return itk::ITK_THREAD_RETURN_DEFAULT_VALUE;

 } // end GetValueThreaderCallback()

@@ -971,14 +912,14 @@ AdvancedImageToImageMetric< TFixedImage, TMovingImage >
 ::GetValueAndDerivativeThreaderCallback( void * arg )
 {
  ThreadInfoType * infoStruct = static_cast< ThreadInfoType * >( arg );
-  ThreadIdType     threadID   = infoStruct->ThreadID;
+  ThreadIdType     threadID   = infoStruct->WorkUnitID;

  MultiThreaderParameterType * temp
    = static_cast< MultiThreaderParameterType * >( infoStruct->UserData );

  temp->st_Metric->ThreadedGetValueAndDerivative( threadID );

-  return ITK_THREAD_RETURN_VALUE;
+  return itk::ITK_THREAD_RETURN_DEFAULT_VALUE;

 } // end GetValueAndDerivativeThreaderCallback()

@@ -1012,15 +953,15 @@ AdvancedImageToImageMetric< TFixedImage, TMovingImage >
 ::AccumulateDerivativesThreaderCallback( void * arg )
 {
  ThreadInfoType * infoStruct  = static_cast< ThreadInfoType * >( arg );
-  ThreadIdType     threadID    = infoStruct->ThreadID;
-  ThreadIdType     nrOfThreads = infoStruct->NumberOfThreads;
+  ThreadIdType     threadID    = infoStruct->WorkUnitID;
+  ThreadIdType     nrOfThreads = infoStruct->NumberOfWorkUnits;

  MultiThreaderParameterType * temp
    = static_cast< MultiThreaderParameterType * >( infoStruct->UserData );

  const unsigned int numPar  = temp->st_Metric->GetNumberOfParameters();
  const unsigned int subSize = static_cast< unsigned int >(
-    vcl_ceil( static_cast< double >( numPar )
+    std::ceil( static_cast< double >( numPar )
    / static_cast< double >( nrOfThreads ) ) );
  const unsigned int jmin = threadID * subSize;
  unsigned int       jmax = ( threadID + 1 ) * subSize;
@@ -1044,7 +985,7 @@ AdvancedImageToImageMetric< TFixedImage, TMovingImage >
    temp->st_DerivativePointer[ j ] = tmp * normalization;
  }

-  return ITK_THREAD_RETURN_VALUE;
+  return itk::ITK_THREAD_RETURN_DEFAULT_VALUE;

 } // end AccumulateDerivativesThreaderCallback()


--- a/Common/CostFunctions/itkExponentialLimiterFunction.h
+++ b/Common/CostFunctions/itkExponentialLimiterFunction.h
@@ -67,18 +67,18 @@ public:
  typedef typename Superclass::DerivativeType      DerivativeType;

  /** Limit the input value */
-  virtual OutputType Evaluate( const InputType & input ) const;
+  OutputType Evaluate( const InputType & input ) const override;

  /** Limit the input value and change the input function derivative accordingly */
-  virtual OutputType Evaluate( const InputType & input, DerivativeType & derivative ) const;
+  OutputType Evaluate( const InputType & input, DerivativeType & derivative ) const override;

  /** Initialize the limiter; calls the ComputeLimiterSettings() function */
-  virtual void Initialize( void ) throw ( ExceptionObject );
+  void Initialize( void ) override;

 protected:

  ExponentialLimiterFunction();
-  ~ExponentialLimiterFunction(){}
+  ~ExponentialLimiterFunction() override{}

  virtual void ComputeLimiterSettings( void );


--- a/Common/CostFunctions/itkExponentialLimiterFunction.hxx
+++ b/Common/CostFunctions/itkExponentialLimiterFunction.hxx
@@ -33,7 +33,7 @@ ExponentialLimiterFunction< TInput, NDimension >
 ::ExponentialLimiterFunction()
 {
  this->ComputeLimiterSettings();
-}   // end Constructor
+} // end Constructor


 /**
@@ -43,10 +43,10 @@ ExponentialLimiterFunction< TInput, NDimension >
 template< class TInput, unsigned int NDimension >
 void
 ExponentialLimiterFunction< TInput, NDimension >
-::Initialize( void ) throw ( ExceptionObject )
+::Initialize( void )
 {
  this->ComputeLimiterSettings();
-}   // end Initialize()
+} // end Initialize()


 /**
@@ -63,7 +63,7 @@ ExponentialLimiterFunction< TInput, NDimension >
  if( diffU > 1e-10 )
  {
    return static_cast< OutputType >(
-      this->m_UTminUB * vcl_exp( this->m_UTminUBinv * diffU ) + this->m_UpperBound );
+      this->m_UTminUB * std::exp( this->m_UTminUBinv * diffU ) + this->m_UpperBound );
  }

  /** Apply a soft limit if the input is smaller than the LowerThreshold */
@@ -71,12 +71,12 @@ ExponentialLimiterFunction< TInput, NDimension >
  if( diffL < -1e-10 )
  {
    return static_cast< OutputType >(
-      this->m_LTminLB * vcl_exp( this->m_LTminLBinv * diffL ) + this->m_LowerBound );
+      this->m_LTminLB * std::exp( this->m_LTminLBinv * diffL ) + this->m_LowerBound );
  }

  /** Leave the value as it is */
  return static_cast< OutputType >( input );
-}   // end Evaluate()
+} // end Evaluate()


 /**
@@ -92,7 +92,7 @@ ExponentialLimiterFunction< TInput, NDimension >
  const double diffU = static_cast< double >( input - this->m_UpperThreshold );
  if( diffU > 1e-10 )
  {
-    const double temp           = this->m_UTminUB * vcl_exp( this->m_UTminUBinv * diffU );
+    const double temp           = this->m_UTminUB * std::exp( this->m_UTminUBinv * diffU );
    const double gradientfactor = this->m_UTminUBinv * temp;
    for( unsigned int i = 0; i < Dimension; ++i )
    {
@@ -105,7 +105,7 @@ ExponentialLimiterFunction< TInput, NDimension >
  const double diffL = static_cast< double >( input - this->m_LowerThreshold );
  if( diffL < -1e-10 )
  {
-    const double temp           = this->m_LTminLB * vcl_exp( this->m_LTminLBinv * diffL );
+    const double temp           = this->m_LTminLB * std::exp( this->m_LTminLBinv * diffL );
    const double gradientfactor = this->m_LTminLBinv * temp;
    for( unsigned int i = 0; i < Dimension; ++i )
    {
@@ -116,7 +116,7 @@ ExponentialLimiterFunction< TInput, NDimension >

  /** Leave the value and derivative as they are */
  return static_cast< OutputType >( input );
-}   // end Evaluate()
+} // end Evaluate()


 /**
@@ -150,7 +150,7 @@ ExponentialLimiterFunction< TInput, NDimension >::ComputeLimiterSettings( void )
    this->m_LTminLB    = 0.0;
    this->m_LTminLBinv = 0.0;
  }
-}   // end ComputeLimiterSettings()
+} // end ComputeLimiterSettings()


 } // end namespace itk

--- a/Common/CostFunctions/itkHardLimiterFunction.h
+++ b/Common/CostFunctions/itkHardLimiterFunction.h
@@ -63,15 +63,15 @@ public:
  typedef typename Superclass::DerivativeType      DerivativeType;

  /** Limit the input value */
-  virtual OutputType Evaluate( const InputType & input ) const;
+  OutputType Evaluate( const InputType & input ) const override;

  /** Limit the input value and change the input function derivative accordingly */
-  virtual OutputType Evaluate( const InputType & input, DerivativeType & derivative ) const;
+  OutputType Evaluate( const InputType & input, DerivativeType & derivative ) const override;

 protected:

  HardLimiterFunction(){}
-  ~HardLimiterFunction(){}
+  ~HardLimiterFunction() override{}

 private:


--- a/Common/CostFunctions/itkHardLimiterFunction.hxx
+++ b/Common/CostFunctions/itkHardLimiterFunction.hxx
@@ -29,9 +29,9 @@ typename HardLimiterFunction< TInput, NDimension >::OutputType
 HardLimiterFunction< TInput, NDimension >
 ::Evaluate( const InputType & input ) const
 {
-  OutputType output = vnl_math_min( static_cast< OutputType >( input ), this->m_UpperBound );
-  return ( vnl_math_max( output, this->m_LowerBound ) );
-}   // end Evaluate()
+  OutputType output = std::min( static_cast< OutputType >( input ), this->m_UpperBound );
+  return ( std::max( output, this->m_LowerBound ) );
+} // end Evaluate()


 template< class TInput, unsigned int NDimension >
@@ -50,7 +50,7 @@ HardLimiterFunction< TInput, NDimension >
    return ( this->m_LowerBound );
  }
  return ( static_cast< OutputType >( input ) );
-}   // end Evaluate()
+} // end Evaluate()


 } // end namespace itk

--- a/Common/CostFunctions/itkImageToImageMetricWithFeatures.h
+++ b/Common/CostFunctions/itkImageToImageMetricWithFeatures.h
@@ -202,7 +202,7 @@ public:


  /** Initialize the metric. */
-  virtual void Initialize( void ) throw ( ExceptionObject );
+  virtual void Initialize( void );

 protected:


--- a/Common/CostFunctions/itkImageToImageMetricWithFeatures.hxx
+++ b/Common/CostFunctions/itkImageToImageMetricWithFeatures.hxx
@@ -34,7 +34,7 @@ ImageToImageMetricWithFeatures< TFixedImage, TMovingImage, TFixedFeatureImage, T
  this->m_NumberOfFixedFeatureImages  = 0;
  this->m_NumberOfMovingFeatureImages = 0;

-}   // end Constructor
+} // end Constructor


 /**
@@ -44,7 +44,7 @@ ImageToImageMetricWithFeatures< TFixedImage, TMovingImage, TFixedFeatureImage, T
 template< class TFixedImage, class TMovingImage, class TFixedFeatureImage, class TMovingFeatureImage >
 void
 ImageToImageMetricWithFeatures< TFixedImage, TMovingImage, TFixedFeatureImage, TMovingFeatureImage >
-::Initialize( void ) throw ( ExceptionObject )
+::Initialize( void )
 {
  /** Call the superclass. */
  this->Superclass::Initialize();
@@ -86,7 +86,7 @@ ImageToImageMetricWithFeatures< TFixedImage, TMovingImage, TFixedFeatureImage, T
  /** Check if the moving feature image interpolators are B-spline interpolators. */
  this->CheckForBSplineFeatureInterpolators();

-}   // end Initialize()
+} // end Initialize()


 /**
@@ -106,7 +106,7 @@ ImageToImageMetricWithFeatures< TFixedImage, TMovingImage, TFixedFeatureImage, T
    this->Modified();
  }

-}   // end SetNumberOfFixedFeatureImages()
+} // end SetNumberOfFixedFeatureImages()


 /**
@@ -134,7 +134,7 @@ ImageToImageMetricWithFeatures< TFixedImage, TMovingImage, TFixedFeatureImage, T
    }
  }

-}   // end SetFixedFeatureImage()
+} // end SetFixedFeatureImage()


 /**
@@ -148,7 +148,7 @@ TFixedFeatureImage, TMovingFeatureImage >::FixedFeatureImageType
 ::GetFixedFeatureImage( unsigned int i ) const
 {
  return this->m_FixedFeatureImages[ i ].GetPointer();
-}   // end GetFixedFeatureImage()
+} // end GetFixedFeatureImage()

 /**
 * ********************* SetFixedFeatureInterpolator ****************************
@@ -175,7 +175,7 @@ ImageToImageMetricWithFeatures< TFixedImage, TMovingImage, TFixedFeatureImage, T
    }
  }

-}   // end SetFixedFeatureInterpolator()
+} // end SetFixedFeatureInterpolator()


 /**
@@ -189,7 +189,7 @@ TFixedFeatureImage, TMovingFeatureImage >::FixedFeatureInterpolatorType
 ::GetFixedFeatureInterpolator( unsigned int i ) const
 {
  return this->m_FixedFeatureInterpolators[ i ].GetPointer();
-}   // end GetFixedFeatureInterpolator()
+} // end GetFixedFeatureInterpolator()

 /**
 * ********************* SetNumberOfMovingFeatureImages ****************************
@@ -208,7 +208,7 @@ ImageToImageMetricWithFeatures< TFixedImage, TMovingImage, TFixedFeatureImage, T
    this->Modified();
  }

-}   // end SetNumberOfMovingFeatureImages()
+} // end SetNumberOfMovingFeatureImages()


 /**
@@ -236,7 +236,7 @@ ImageToImageMetricWithFeatures< TFixedImage, TMovingImage, TFixedFeatureImage, T
    }
  }

-}   // end SetMovingFeatureImage()
+} // end SetMovingFeatureImage()


 /**
@@ -250,7 +250,7 @@ TFixedFeatureImage, TMovingFeatureImage >::MovingFeatureImageType
 ::GetMovingFeatureImage( unsigned int i ) const
 {
  return this->m_MovingFeatureImages[ i ].GetPointer();
-}   // end GetMovingFeatureImage()
+} // end GetMovingFeatureImage()

 /**
 * ********************* SetMovingFeatureInterpolator ****************************
@@ -277,7 +277,7 @@ ImageToImageMetricWithFeatures< TFixedImage, TMovingImage, TFixedFeatureImage, T
    }
  }

-}   // end SetMovingFeatureInterpolator()
+} // end SetMovingFeatureInterpolator()


 /**
@@ -290,7 +290,7 @@ const typename ImageToImageMetricWithFeatures< TFixedImage, TMovingImage, TFixed
 ::GetMovingFeatureInterpolator( unsigned int i ) const
 {
  return this->m_MovingFeatureInterpolators[ i ].GetPointer();
-}   // end GetMovingFeatureInterpolator()
+} // end GetMovingFeatureInterpolator()

 /**
 * ****************** CheckForBSplineFeatureInterpolators **********************
@@ -327,9 +327,9 @@ ImageToImageMetricWithFeatures< TFixedImage, TMovingImage, TFixedFeatureImage, T
      itkDebugMacro( << "Interpolator " << i << " is NOT B-spline." );
      itkExceptionMacro( << "Interpolator " << i << " is NOT B-spline." );
    }
-  }   // end for-loop
+  } // end for-loop

-}   // end CheckForBSplineFeatureInterpolators()
+} // end CheckForBSplineFeatureInterpolators()


 /**
@@ -370,7 +370,7 @@ ImageToImageMetricWithFeatures< TFixedImage, TMovingImage, TFixedFeatureImage, T
       << this->m_MovingFeatureInterpolators[ i ].GetPointer() << std::endl;
  }

-}   // end PrintSelf()
+} // end PrintSelf()


 } // end namespace itk

--- a/Common/CostFunctions/itkLimiterFunctionBase.h
+++ b/Common/CostFunctions/itkLimiterFunctionBase.h
@@ -19,7 +19,7 @@
 #define __itkLimiterFunctionBase_h

 #include "itkFunctionBase.h"
-#include "itkExceptionObject.h"
+#include "itkMacro.h"

 namespace itk
 {
@@ -73,7 +73,7 @@ public:
    itkGetStaticConstMacro( Dimension ) >            DerivativeType;

  /** Limit the input value. */
-  virtual OutputType Evaluate( const InputType & input ) const = 0;
+  OutputType Evaluate( const InputType & input ) const override = 0;

  /** Limit the input value and change the input function derivative accordingly */
  virtual OutputType Evaluate( const InputType & input, DerivativeType & derivative ) const = 0;
@@ -99,7 +99,7 @@ public:
  itkGetConstMacro( LowerThreshold, InputType );

  /** Initialize the limiter */
-  virtual void Initialize( void ) throw ( ExceptionObject ) {}
+  virtual void Initialize( void ) {}

 protected:

@@ -114,7 +114,7 @@ protected:
  }


-  ~LimiterFunctionBase(){}
+  ~LimiterFunctionBase() override{}

  OutputType m_UpperBound;
  OutputType m_LowerBound;

--- a/Common/CostFunctions/itkMultiInputImageToImageMetricBase.h
+++ b/Common/CostFunctions/itkMultiInputImageToImageMetricBase.h
@@ -115,7 +115,7 @@ public:
  virtual void SetFixedImage( const FixedImageType * _arg, unsigned int pos );

  /** Set the first fixed image. */
-  virtual void SetFixedImage( const FixedImageType * _arg )
+  void SetFixedImage( const FixedImageType * _arg ) override
  {
    this->SetFixedImage( _arg, 0 );
  }
@@ -125,7 +125,7 @@ public:
  virtual const FixedImageType * GetFixedImage( unsigned int pos ) const;

  /** Get the first fixed image. */
-  virtual const FixedImageType * GetFixedImage( void ) const
+  const FixedImageType * GetFixedImage( void ) const override
  {
    return this->GetFixedImage( 0 );
  }
@@ -143,7 +143,7 @@ public:
  virtual void SetFixedImageMask( FixedImageMaskType * _arg, unsigned int pos );

  /** Set the first fixed image mask. */
-  virtual void SetFixedImageMask( FixedImageMaskType * _arg )
+  void SetFixedImageMask( FixedImageMaskType * _arg ) override
  {
    this->SetFixedImageMask( _arg, 0 );
  }
@@ -153,7 +153,7 @@ public:
  virtual FixedImageMaskType * GetFixedImageMask( unsigned int pos ) const;

  /** Get the first fixed image mask. */
-  virtual FixedImageMaskType * GetFixedImageMask( void ) const
+  FixedImageMaskType * GetFixedImageMask( void ) const override
  {
    return this->GetFixedImageMask( 0 );
  }
@@ -171,7 +171,7 @@ public:
  virtual void SetFixedImageRegion( const FixedImageRegionType _arg, unsigned int pos );

  /** Set the first fixed image region. */
-  virtual void SetFixedImageRegion( const FixedImageRegionType _arg )
+  void SetFixedImageRegion( const FixedImageRegionType _arg ) override
  {
    this->SetFixedImageRegion( _arg, 0 );
  }
@@ -181,7 +181,7 @@ public:
  virtual const FixedImageRegionType & GetFixedImageRegion( unsigned int pos ) const;

  /** Get the first fixed image region. */
-  virtual const FixedImageRegionType & GetFixedImageRegion( void ) const
+  const FixedImageRegionType & GetFixedImageRegion( void ) const override
  {
    return this->GetFixedImageRegion( 0 );
  }
@@ -199,7 +199,7 @@ public:
  virtual void SetMovingImage( const MovingImageType * _arg, unsigned int pos );

  /** Set the first moving image. */
-  virtual void SetMovingImage( const MovingImageType * _arg )
+  void SetMovingImage( const MovingImageType * _arg ) override
  {
    this->SetMovingImage( _arg, 0 );
  }
@@ -209,7 +209,7 @@ public:
  virtual const MovingImageType * GetMovingImage( unsigned int pos ) const;

  /** Get the first moving image. */
-  virtual const MovingImageType * GetMovingImage( void ) const
+  const MovingImageType * GetMovingImage( void ) const override
  {
    return this->GetMovingImage( 0 );
  }
@@ -227,7 +227,7 @@ public:
  virtual void SetMovingImageMask( MovingImageMaskType * _arg, unsigned int pos );

  /** Set the first moving image mask. */
-  virtual void SetMovingImageMask( MovingImageMaskType * _arg )
+  void SetMovingImageMask( MovingImageMaskType * _arg ) override
  {
    this->SetMovingImageMask( _arg, 0 );
  }
@@ -237,7 +237,7 @@ public:
  virtual MovingImageMaskType * GetMovingImageMask( unsigned int pos ) const;

  /** Get the first moving image mask. */
-  virtual MovingImageMaskType * GetMovingImageMask( void ) const
+  MovingImageMaskType * GetMovingImageMask( void ) const override
  {
    return this->GetMovingImageMask( 0 );
  }
@@ -257,7 +257,7 @@ public:
  virtual void SetInterpolator( InterpolatorType * _arg, unsigned int pos );

  /** Set the first interpolator. */
-  virtual void SetInterpolator( InterpolatorType * _arg )
+  void SetInterpolator( InterpolatorType * _arg ) override
  {
    return this->SetInterpolator( _arg, 0 );
  }
@@ -267,7 +267,7 @@ public:
  virtual InterpolatorType * GetInterpolator( unsigned int pos ) const;

  /** Get the first interpolator. */
-  virtual InterpolatorType * GetInterpolator( void ) const
+  InterpolatorType * GetInterpolator( void ) const override
  {
    return this->GetInterpolator( 0 );
  }
@@ -315,7 +315,7 @@ public:
  /** ******************** Other public functions ******************** */

  /** Initialisation. */
-  virtual void Initialize( void ) throw ( ExceptionObject );
+  void Initialize( void ) override;

 protected:

@@ -323,7 +323,7 @@ protected:
  MultiInputImageToImageMetricBase();

  /** Destructor. */
-  virtual ~MultiInputImageToImageMetricBase() {}
+  ~MultiInputImageToImageMetricBase() override {}

  /** Typedef's from the Superclass. */
  typedef typename Superclass::MovingImagePointType           MovingImagePointType;
@@ -337,7 +337,7 @@ protected:
  typedef std::vector< BSplineInterpolatorPointer >    BSplineInterpolatorVectorType;

  /** Initialize variables related to the image sampler; called by Initialize. */
-  virtual void InitializeImageSampler( void ) throw ( ExceptionObject );
+  void InitializeImageSampler( void ) override;

  /** Check if all interpolators (for the moving image) are of type
   * BSplineInterpolateImageFunction.
@@ -347,14 +347,14 @@ protected:
  /** Check if mappedPoint is inside all moving images.
   * If so, the moving image value and possibly derivative are computed.
   */
-  virtual bool EvaluateMovingImageValueAndDerivative(
+  bool EvaluateMovingImageValueAndDerivative(
    const MovingImagePointType & mappedPoint,
    RealType & movingImageValue,
-    MovingImageDerivativeType * gradient ) const;
+    MovingImageDerivativeType * gradient ) const override;

  /** IsInsideMovingMask: Returns the AND of all moving image masks. */
-  virtual bool IsInsideMovingMask(
-    const MovingImagePointType & mappedPoint ) const;
+  bool IsInsideMovingMask(
+    const MovingImagePointType & mappedPoint ) const override;

  /** Protected member variables. */
  FixedImageVectorType             m_FixedImageVector;

--- a/Common/CostFunctions/itkMultiInputImageToImageMetricBase.hxx
+++ b/Common/CostFunctions/itkMultiInputImageToImageMetricBase.hxx
@@ -118,7 +118,7 @@ MultiInputImageToImageMetricBase< TFixedImage, TMovingImage >

  this->m_InterpolatorsAreBSpline = false;

-}   // end Constructor()
+} // end Constructor()


 /** Set components. */
@@ -137,6 +137,7 @@ itkImplementationGetObjectMacro( MovingImageMask, MovingImageMaskType );
 itkImplementationGetObjectMacro( Interpolator, InterpolatorType );
 itkImplementationGetObjectMacro( FixedImageInterpolator, FixedImageInterpolatorType );

+
 /**
 * ************************ SetFixedImageRegion *************************
 */
@@ -161,7 +162,7 @@ MultiInputImageToImageMetricBase< TFixedImage, TMovingImage >
    this->Modified();
  }

-}   // end SetFixedImageRegion()
+} // end SetFixedImageRegion()


 /**
@@ -181,7 +182,8 @@ const typename MultiInputImageToImageMetricBase< TFixedImage, TMovingImage >

  return this->m_FixedImageRegionVector[ pos ];

-}   // end GetFixedImageRegion()
+} // end GetFixedImageRegion()
+

 /**
 * ****************** CheckForBSplineInterpolators **********************
@@ -217,9 +219,9 @@ MultiInputImageToImageMetricBase< TFixedImage, TMovingImage >
      itkDebugMacro( << "Interpolator " << i << " is NOT B-spline." );
      itkExceptionMacro( << "Interpolator " << i << " is NOT B-spline." );
    }
-  }   // end for-loop
+  } // end for-loop

-}   // end CheckForBSplineInterpolators()
+} // end CheckForBSplineInterpolators()


 /**
@@ -229,7 +231,7 @@ MultiInputImageToImageMetricBase< TFixedImage, TMovingImage >
 template< class TFixedImage, class TMovingImage >
 void
 MultiInputImageToImageMetricBase< TFixedImage, TMovingImage >
-::Initialize( void ) throw ( ExceptionObject )
+::Initialize( void )
 {
  /** Connect the interpolators. */
  for( unsigned int i = 0; i < this->GetNumberOfInterpolators(); ++i )
@@ -249,7 +251,7 @@ MultiInputImageToImageMetricBase< TFixedImage, TMovingImage >
  /** Call the superclass' implementation. */
  this->Superclass::Initialize();

-}   // end Initialize()
+} // end Initialize()


 /**
@@ -259,7 +261,7 @@ MultiInputImageToImageMetricBase< TFixedImage, TMovingImage >
 template< class TFixedImage, class TMovingImage >
 void
 MultiInputImageToImageMetricBase< TFixedImage, TMovingImage >
-::InitializeImageSampler( void ) throw ( ExceptionObject )
+::InitializeImageSampler( void )
 {
  if( this->GetUseImageSampler() )
  {
@@ -288,7 +290,7 @@ MultiInputImageToImageMetricBase< TFixedImage, TMovingImage >
    }
  }

-}   // end InitializeImageSampler()
+} // end InitializeImageSampler()


 /**
@@ -317,7 +319,7 @@ MultiInputImageToImageMetricBase< TFixedImage, TMovingImage >
  return this->Superclass::EvaluateMovingImageValueAndDerivative(
    mappedPoint, movingImageValue, gradient );

-}   // end EvaluateMovingImageValueAndDerivative()
+} // end EvaluateMovingImageValueAndDerivative()


 /**
@@ -340,7 +342,7 @@ MultiInputImageToImageMetricBase< TFixedImage, TMovingImage >
    MovingImageMaskPointer movingImageMask = this->GetMovingImageMask( i );
    if( movingImageMask.IsNotNull() )
    {
-      inside &= movingImageMask->IsInside( mappedPoint );
+      inside &= movingImageMask->IsInsideInWorldSpace( mappedPoint );
    }

    /** If the point falls outside one mask, we can skip the rest. */
@@ -351,7 +353,7 @@ MultiInputImageToImageMetricBase< TFixedImage, TMovingImage >
  }
  return inside;

-}   // end IsInsideMovingMask()
+} // end IsInsideMovingMask()


 } // end namespace itk

--- a/Common/CostFunctions/itkParzenWindowHistogramImageToImageMetric.h
+++ b/Common/CostFunctions/itkParzenWindowHistogramImageToImageMetric.h
@@ -33,7 +33,7 @@ namespace itk
 * where the probability density distribution are estimated using
 * Parzen histograms.
 *
- * One the PDF's have been constructed, the metric value and derivative
+ * Once the PDF's have been constructed, the metric value and derivative
 * can be computed. Inheriting classes should make sure to call
 * the function ComputePDFs(AndPDFDerivatives) before using m_JointPDF and m_Alpha
 * (and m_JointPDFDerivatives).
@@ -143,7 +143,7 @@ public:
   * (3) InitializeKernels()
   * (4) Resize AlphaDerivatives
   */
-  void Initialize( void ) throw ( ExceptionObject );
+  void Initialize( void ) override;

  /** Get the derivatives of the match measure. This method simply calls the
   * the GetValueAndDerivative, since this will be mostly almost as fast
@@ -151,7 +151,7 @@ public:
   */
  void GetDerivative(
    const ParametersType & parameters,
-    DerivativeType & Derivative ) const;
+    DerivativeType & Derivative ) const override;

  /**  Get the value and derivatives for single valued optimizers.
   * This method calls this->GetValueAndAnalyticDerivative or
@@ -159,7 +159,7 @@ public:
   * m_UseFiniteDifferenceDerivative.
   */
  void GetValueAndDerivative( const ParametersType & parameters,
-    MeasureType & value, DerivativeType & derivative ) const;
+    MeasureType & value, DerivativeType & derivative ) const override;

  /** Number of bins to use for the fixed image in the histogram.
   * Typical value is 32.  The minimum value is 4 due to the padding
@@ -221,10 +221,10 @@ protected:
  ParzenWindowHistogramImageToImageMetric();

  /** The destructor. */
-  virtual ~ParzenWindowHistogramImageToImageMetric();
+  ~ParzenWindowHistogramImageToImageMetric() override;

  /** Print Self. */
-  void PrintSelf( std::ostream & os, Indent indent ) const;
+  void PrintSelf( std::ostream & os, Indent indent ) const override;

  /** Protected Typedefs ******************/

@@ -322,7 +322,7 @@ protected:
  mutable ThreadIdType                                                       m_ParzenWindowHistogramGetValueAndDerivativePerThreadVariablesSize;

  /** Initialize threading related parameters. */
-  virtual void InitializeThreadingParameters( void ) const;
+  void InitializeThreadingParameters( void ) const override;

  /** Multi-threaded versions of the ComputePDF function. */
  inline void ThreadedComputePDFs( ThreadIdType threadId );