From 3675cd1fa48febf8b8ef78d0c6e8d7841a7182d0 Mon Sep 17 00:00:00 2001 From: Niels Dekker Date: Thu, 25 Nov 2021 13:45:24 +0100 Subject: [PATCH] PERF: Make AdvancedNormalizedCorrelationImageToImageMetric faster Added a new member function to `AdvancedImageToImageMetric`, `FastEvaluateMovingImageValueAndDerivative`, which calls the multi-threaded overload of ITK's `itk::BSplineInterpolateImageFunction::EvaluateValueAndDerivativeAtContinuousIndex`, indirectly. (It does so via `EvaluateMovingImageValueAndDerivativeWithOptionalThreadId`, another new member function.) Made `AdvancedNormalizedCorrelationImageToImageMetric::ThreadedGetValueAndDerivative` faster, by calling this new `FastEvaluateMovingImageValueAndDerivative`. A large performance improvement was observed for GoogleTest unit test `itkElastixRegistrationMethod.EulerDiscRotation2D` (which uses the "AdvancedNormalizedCorrelation" metric), from ~1.5 second before this commit down to ~0.9 second after this commit, using Visual Studio 2019, Release configuration. (For a Debug configuration even from ~15 seconds before, down to ~4 seconds after this commit.) --- .../itkAdvancedImageToImageMetric.h | 23 +++++++++++++++++- .../itkAdvancedImageToImageMetric.hxx | 24 ++++++++++++++----- ...ormalizedCorrelationImageToImageMetric.hxx | 3 ++- 3 files changed, 42 insertions(+), 8 deletions(-) diff --git a/Common/CostFunctions/itkAdvancedImageToImageMetric.h b/Common/CostFunctions/itkAdvancedImageToImageMetric.h index 6c8c560b6..de0eb909d 100644 --- a/Common/CostFunctions/itkAdvancedImageToImageMetric.h +++ b/Common/CostFunctions/itkAdvancedImageToImageMetric.h @@ -507,7 +507,21 @@ class ITK_TEMPLATE_EXPORT AdvancedImageToImageMetric : public ImageToImageMetric virtual bool EvaluateMovingImageValueAndDerivative(const MovingImagePointType & mappedPoint, RealType & movingImageValue, - MovingImageDerivativeType * gradient) const; + MovingImageDerivativeType * gradient) const + { + return EvaluateMovingImageValueAndDerivativeWithOptionalThreadId(mappedPoint, movingImageValue, gradient); + } + + /* A faster version of `EvaluateMovingImageValueAndDerivative`: Non-virtual, using multithreading, and doing less + * dynamic memory allocation/decallocation operations, internally. */ + bool + FastEvaluateMovingImageValueAndDerivative(const MovingImagePointType & mappedPoint, + RealType & movingImageValue, + MovingImageDerivativeType * gradient, + const ThreadIdType threadId) const + { + return EvaluateMovingImageValueAndDerivativeWithOptionalThreadId(mappedPoint, movingImageValue, gradient, threadId); + } /** Computes the inner product of transform Jacobian with moving image gradient. * The results are stored in imageJacobian, which is supposed @@ -571,6 +585,13 @@ class ITK_TEMPLATE_EXPORT AdvancedImageToImageMetric : public ImageToImageMetric void operator=(const Self &) = delete; + template + bool + EvaluateMovingImageValueAndDerivativeWithOptionalThreadId(const MovingImagePointType & mappedPoint, + RealType & movingImageValue, + MovingImageDerivativeType * gradient, + const TOptionalThreadId... optionalThreadId) const; + /** Private member variables. */ bool m_UseImageSampler{ false }; bool m_UseFixedImageLimiter{ false }; diff --git a/Common/CostFunctions/itkAdvancedImageToImageMetric.hxx b/Common/CostFunctions/itkAdvancedImageToImageMetric.hxx index 21f93ee0d..7bcb59540 100644 --- a/Common/CostFunctions/itkAdvancedImageToImageMetric.hxx +++ b/Common/CostFunctions/itkAdvancedImageToImageMetric.hxx @@ -112,6 +112,15 @@ AdvancedImageToImageMetric::Initialize(void) if (this->m_UseMultiThread) { this->InitializeThreadingParameters(); + + const auto setNumberOfWorkUnitsIfNotNull = [this](const auto bsplineInterpolator) { + if (!bsplineInterpolator.IsNull()) + { + bsplineInterpolator->SetNumberOfWorkUnits(this->Superclass::GetNumberOfWorkUnits()); + } + }; + setNumberOfWorkUnitsIfNotNull(m_BSplineInterpolator); + setNumberOfWorkUnitsIfNotNull(m_BSplineInterpolatorFloat); } } // end Initialize() @@ -501,15 +510,17 @@ AdvancedImageToImageMetric::CheckForBSplineTransform( /** - * ******************* EvaluateMovingImageValueAndDerivative ****************** + * ******************* EvaluateMovingImageValueAndDerivativeWithOptionalThreadId ****************** */ template +template bool -AdvancedImageToImageMetric::EvaluateMovingImageValueAndDerivative( +AdvancedImageToImageMetric::EvaluateMovingImageValueAndDerivativeWithOptionalThreadId( const MovingImagePointType & mappedPoint, RealType & movingImageValue, - MovingImageDerivativeType * gradient) const + MovingImageDerivativeType * gradient, + const TOptionalThreadId... optionalThreadId) const { /** Check if mapped point inside image buffer. */ MovingImageContinuousIndexType cindex; @@ -523,13 +534,14 @@ AdvancedImageToImageMetric::EvaluateMovingImageValueA if (this->m_InterpolatorIsBSpline && !this->GetComputeGradient()) { /** Compute moving image value and gradient using the B-spline kernel. */ - this->m_BSplineInterpolator->EvaluateValueAndDerivativeAtContinuousIndex(cindex, movingImageValue, *gradient); + this->m_BSplineInterpolator->EvaluateValueAndDerivativeAtContinuousIndex( + cindex, movingImageValue, *gradient, optionalThreadId...); } else if (this->m_InterpolatorIsBSplineFloat && !this->GetComputeGradient()) { /** Compute moving image value and gradient using the B-spline kernel. */ this->m_BSplineInterpolatorFloat->EvaluateValueAndDerivativeAtContinuousIndex( - cindex, movingImageValue, *gradient); + cindex, movingImageValue, *gradient, optionalThreadId...); } else if (this->m_InterpolatorIsReducedBSpline && !this->GetComputeGradient()) { @@ -606,7 +618,7 @@ AdvancedImageToImageMetric::EvaluateMovingImageValueA return sampleOk; -} // end EvaluateMovingImageValueAndDerivative() +} // end EvaluateMovingImageValueAndDerivativeWithOptionalThreadId() /** diff --git a/Components/Metrics/AdvancedNormalizedCorrelation/itkAdvancedNormalizedCorrelationImageToImageMetric.hxx b/Components/Metrics/AdvancedNormalizedCorrelation/itkAdvancedNormalizedCorrelationImageToImageMetric.hxx index 29d7fc840..45c0e8fe8 100644 --- a/Components/Metrics/AdvancedNormalizedCorrelation/itkAdvancedNormalizedCorrelationImageToImageMetric.hxx +++ b/Components/Metrics/AdvancedNormalizedCorrelation/itkAdvancedNormalizedCorrelationImageToImageMetric.hxx @@ -590,7 +590,8 @@ AdvancedNormalizedCorrelationImageToImageMetric::Thre */ if (sampleOk) { - sampleOk = this->EvaluateMovingImageValueAndDerivative(mappedPoint, movingImageValue, &movingImageDerivative); + sampleOk = this->FastEvaluateMovingImageValueAndDerivative( + mappedPoint, movingImageValue, &movingImageDerivative, threadId); } if (sampleOk)