From de18d674281f5fcfc391a15690068e622174c701 Mon Sep 17 00:00:00 2001 From: Mikhail Polkovnikov Date: Fri, 28 Feb 2025 13:06:55 +0300 Subject: [PATCH 1/3] ENH: Use lambda functions instead of functor classes in rtkJosephForwardProjectionFilter Lambdas were added to the classes: JosephForwardProjectionImageFilter, JosephForwardAttenuatedProjectionImageFilter, MaximumIntensityProjectionImageFilter. Python wrappers for classes above have been modified. Similar lambda functions can be added to back projection classes. --- ...phForwardAttenuatedProjectionImageFilter.h | 236 +----------------- ...ForwardAttenuatedProjectionImageFilter.hxx | 128 ++++++---- .../rtkJosephForwardProjectionImageFilter.h | 216 +++++----------- .../rtkJosephForwardProjectionImageFilter.hxx | 154 +++++------- ...rtkMaximumIntensityProjectionImageFilter.h | 121 ++------- ...kMaximumIntensityProjectionImageFilter.hxx | 60 +++++ test/rtkforwardprojectiontest.cxx | 11 + test/rtkmaximumintensityprojectiontest.cxx | 84 ++++++- ...orwardAttenuatedProjectionImageFilter.wrap | 25 -- ...rtkJosephForwardProjectionImageFilter.wrap | 7 - ...MaximumIntensityProjectionImageFilter.wrap | 27 +- 11 files changed, 368 insertions(+), 701 deletions(-) create mode 100644 include/rtkMaximumIntensityProjectionImageFilter.hxx diff --git a/include/rtkJosephForwardAttenuatedProjectionImageFilter.h b/include/rtkJosephForwardAttenuatedProjectionImageFilter.h index 38a4676e4..31b37ebdf 100644 --- a/include/rtkJosephForwardAttenuatedProjectionImageFilter.h +++ b/include/rtkJosephForwardAttenuatedProjectionImageFilter.h @@ -29,213 +29,6 @@ namespace rtk { -namespace Functor -{ -/** \class InterpolationWeightMultiplicationAttenuated - * \brief Function to multiply the interpolation weights with the projected - * volume values and attenuation map. - * - * \author Antoine Robert - * - * \ingroup RTK Functions - */ -template -class ITK_TEMPLATE_EXPORT InterpolationWeightMultiplicationAttenuated -{ -public: - InterpolationWeightMultiplicationAttenuated() - { - for (std::size_t i = 0; i < itk::ITK_MAX_THREADS; i++) - { - m_AttenuationRay[i] = 0; - m_AttenuationPixel[i] = 0; - m_Ex1[i] = 1; - } - } - - ~InterpolationWeightMultiplicationAttenuated() = default; - bool - operator!=(const InterpolationWeightMultiplicationAttenuated &) const - { - return false; - } - - bool - operator==(const InterpolationWeightMultiplicationAttenuated & other) const - { - return !(*this != other); - } - - inline TOutput - operator()(const ThreadIdType threadId, - const double stepLengthInVoxel, - const TCoordinateType weight, - const TInput * p, - const int i) - { - const double w = weight * stepLengthInVoxel; - - m_AttenuationRay[threadId] += w * (p + m_AttenuationMinusEmissionMapsPtrDiff)[i]; - m_AttenuationPixel[threadId] += w * (p + m_AttenuationMinusEmissionMapsPtrDiff)[i]; - return weight * p[i]; - } - - void - SetAttenuationMinusEmissionMapsPtrDiff(std::ptrdiff_t pd) - { - m_AttenuationMinusEmissionMapsPtrDiff = pd; - } - TOutput * - GetAttenuationRay() - { - return m_AttenuationRay; - } - TOutput * - GetAttenuationPixel() - { - return m_AttenuationPixel; - } - TOutput * - GetEx1() - { - return m_Ex1; - } - -private: - std::ptrdiff_t m_AttenuationMinusEmissionMapsPtrDiff; - TInput m_AttenuationRay[itk::ITK_MAX_THREADS]; - TInput m_AttenuationPixel[itk::ITK_MAX_THREADS]; - TInput m_Ex1[itk::ITK_MAX_THREADS]; -}; - -/** \class ComputeAttenuationCorrection - * \brief Function to compute the attenuation correction on the projection. - * - * \author Antoine Robert - * - * \ingroup RTK Functions - */ -template -class ITK_TEMPLATE_EXPORT ComputeAttenuationCorrection -{ -public: - using VectorType = itk::Vector; - - ComputeAttenuationCorrection() = default; - ~ComputeAttenuationCorrection() = default; - bool - operator!=(const ComputeAttenuationCorrection &) const - { - return false; - } - - bool - operator==(const ComputeAttenuationCorrection & other) const - { - return !(*this != other); - } - - inline void - operator()(const ThreadIdType threadId, TOutput & sumValue, const TInput volumeValue, const VectorType & stepInMM) - { - TInput ex2 = exp(-m_AttenuationRay[threadId] * stepInMM.GetNorm()); - TInput wf; - - if (m_AttenuationPixel[threadId] > 0) - { - wf = (m_Ex1[threadId] - ex2) / m_AttenuationPixel[threadId]; - } - else - { - wf = m_Ex1[threadId] * stepInMM.GetNorm(); - } - - m_Ex1[threadId] = ex2; - m_AttenuationPixel[threadId] = 0; - sumValue += wf * volumeValue; - } - - void - SetAttenuationRayVector(TInput * attenuationRayVector) - { - m_AttenuationRay = attenuationRayVector; - } - void - SetAttenuationPixelVector(TInput * attenuationPixelVector) - { - m_AttenuationPixel = attenuationPixelVector; - } - void - SetEx1(TInput * ex1) - { - m_Ex1 = ex1; - } - -private: - TInput * m_AttenuationRay; - TInput * m_AttenuationPixel; - TInput * m_Ex1; -}; - -/** \class ProjectedValueAccumulationAttenuated - * \brief Function to accumulate the ray casting on the projection. - * - * \author Antoine Robert - * - * \ingroup RTK Functions - */ -template -class ITK_TEMPLATE_EXPORT ProjectedValueAccumulationAttenuated -{ -public: - using VectorType = itk::Vector; - - ProjectedValueAccumulationAttenuated() = default; - ~ProjectedValueAccumulationAttenuated() = default; - bool - operator!=(const ProjectedValueAccumulationAttenuated &) const - { - return false; - } - - bool - operator==(const ProjectedValueAccumulationAttenuated & other) const - { - return !(*this != other); - } - - inline void - operator()(const ThreadIdType threadId, - const TInput & input, - TOutput & output, - const TOutput & rayCastValue, - const VectorType & /*stepInMM*/, - const VectorType & itkNotUsed(source), - const VectorType & itkNotUsed(sourceToPixel), - const VectorType & itkNotUsed(nearestPoint), - const VectorType & itkNotUsed(farthestPoint)) - { - output = input + rayCastValue; - m_Attenuation[threadId] = 0; - m_Ex1[threadId] = 1; - } - - void - SetAttenuationVector(TInput * attenuationVector) - { - m_Attenuation = attenuationVector; - } - void - SetEx1(TInput * ex1) - { - m_Ex1 = ex1; - } - -private: - TInput * m_Attenuation; - TInput * m_Ex1; -}; -} // end namespace Functor /** \class JosephForwardAttenuatedProjectionImageFilter * \brief Joseph forward projection. @@ -252,33 +45,16 @@ class ITK_TEMPLATE_EXPORT ProjectedValueAccumulationAttenuated * \ingroup RTK Projector */ -template < - class TInputImage, - class TOutputImage, - class TInterpolationWeightMultiplication = Functor::InterpolationWeightMultiplicationAttenuated< - typename TInputImage::PixelType, - typename itk::PixelTraits::ValueType>, - class TProjectedValueAccumulation = - Functor::ProjectedValueAccumulationAttenuated, - class TSumAlongRay = - Functor::ComputeAttenuationCorrection> +template class ITK_TEMPLATE_EXPORT JosephForwardAttenuatedProjectionImageFilter - : public JosephForwardProjectionImageFilter + : public JosephForwardProjectionImageFilter { public: ITK_DISALLOW_COPY_AND_MOVE(JosephForwardAttenuatedProjectionImageFilter); /** Standard class type alias. */ using Self = JosephForwardAttenuatedProjectionImageFilter; - using Superclass = JosephForwardProjectionImageFilter; + using Superclass = JosephForwardProjectionImageFilter; using Pointer = itk::SmartPointer; using ConstPointer = itk::SmartPointer; using InputPixelType = typename TInputImage::PixelType; @@ -286,6 +62,7 @@ class ITK_TEMPLATE_EXPORT JosephForwardAttenuatedProjectionImageFilter using OutputImageRegionType = typename TOutputImage::RegionType; using CoordinateType = double; using VectorType = itk::Vector; + using WeightCoordinateType = typename itk::PixelTraits::ValueType; /** ImageDimension constants */ static constexpr unsigned int InputImageDimension = TInputImage::ImageDimension; @@ -311,6 +88,11 @@ class ITK_TEMPLATE_EXPORT JosephForwardAttenuatedProjectionImageFilter * to overwrite the method. */ void VerifyInputInformation() const override; + + std::ptrdiff_t m_AttenuationMinusEmissionMapsPtrDiff; + std::array m_AttenuationRay; + std::array m_AttenuationPixel; + std::array m_Ex1; }; } // end namespace rtk diff --git a/include/rtkJosephForwardAttenuatedProjectionImageFilter.hxx b/include/rtkJosephForwardAttenuatedProjectionImageFilter.hxx index 3973c507b..bec8af91a 100644 --- a/include/rtkJosephForwardAttenuatedProjectionImageFilter.hxx +++ b/include/rtkJosephForwardAttenuatedProjectionImageFilter.hxx @@ -31,33 +31,85 @@ namespace rtk { -template -JosephForwardAttenuatedProjectionImageFilter< - TInputImage, - TOutputImage, - TInterpolationWeightMultiplication, - TProjectedValueAccumulation, - TComputeAttenuationCorrection>::JosephForwardAttenuatedProjectionImageFilter() +template +JosephForwardAttenuatedProjectionImageFilter::JosephForwardAttenuatedProjectionImageFilter() + : JosephForwardProjectionImageFilter() { + std::fill(m_AttenuationRay.begin(), m_AttenuationRay.end(), 0.); + std::fill(m_AttenuationPixel.begin(), m_AttenuationPixel.end(), 0.); + std::fill(m_Ex1.begin(), m_Ex1.end(), 1.); + m_AttenuationMinusEmissionMapsPtrDiff = 0; + + /** \brief Function to multiply the interpolation weights with the projected + * volume values and attenuation map. + * + * \author Antoine Robert + */ + auto interpolationWeightMultiplicationAttenuatedFunc = [this](const ThreadIdType threadId, + double stepLengthInVoxel, + const WeightCoordinateType weight, + const InputPixelType * p, + const int i) -> OutputPixelType { + const double w = weight * stepLengthInVoxel; + + this->m_AttenuationRay[threadId] += w * (p + this->m_AttenuationMinusEmissionMapsPtrDiff)[i]; + this->m_AttenuationPixel[threadId] += w * (p + this->m_AttenuationMinusEmissionMapsPtrDiff)[i]; + return weight * p[i]; + }; + this->SetInterpolationWeightMultiplication(interpolationWeightMultiplicationAttenuatedFunc); + + /** \brief Function to compute the attenuation correction on the projection. + * + * \author Antoine Robert + */ + auto computeAttenuationCorrectionFunc = [this](const ThreadIdType threadId, + OutputPixelType & sumValue, + const InputPixelType volumeValue, + const VectorType & stepInMM) { + InputPixelType ex2 = exp(-1. * this->m_AttenuationRay[threadId] * stepInMM.GetNorm()); + InputPixelType wf; + + if (this->m_AttenuationPixel[threadId] > 0) + { + wf = (this->m_Ex1[threadId] - ex2) / this->m_AttenuationPixel[threadId]; + } + else + { + wf = this->m_Ex1[threadId] * stepInMM.GetNorm(); + } + + this->m_Ex1[threadId] = ex2; + this->m_AttenuationPixel[threadId] = 0; + sumValue += wf * volumeValue; + }; + this->SetSumAlongRay(computeAttenuationCorrectionFunc); + + /** \brief Function to accumulate the ray casting on the projection. + * + * \author Antoine Robert + */ + auto projectedValueAccumulationAttenuatedFunc = [this](const ThreadIdType threadId, + const InputPixelType & input, + OutputPixelType & output, + const OutputPixelType & rayCastValue, + const VectorType & itkNotUsed(stepInMM), + const VectorType & itkNotUsed(source), + const VectorType & itkNotUsed(sourceToPixel), + const VectorType & itkNotUsed(nearestPoint), + const VectorType & itkNotUsed(farthestPoint)) { + output = input + rayCastValue; + this->m_AttenuationRay[threadId] = 0; + this->m_Ex1[threadId] = 1; + }; + this->SetProjectedValueAccumulation(projectedValueAccumulationAttenuatedFunc); + this->SetNumberOfRequiredInputs(3); this->DynamicMultiThreadingOff(); } -template +template void -JosephForwardAttenuatedProjectionImageFilter::GenerateInputRequestedRegion() +JosephForwardAttenuatedProjectionImageFilter::GenerateInputRequestedRegion() { Superclass::GenerateInputRequestedRegion(); // Input 2 is the attenuation map relative to the volume @@ -69,17 +121,9 @@ JosephForwardAttenuatedProjectionImageFilterSetRequestedRegion(reqRegion2); } -template +template void -JosephForwardAttenuatedProjectionImageFilter::VerifyInputInformation() const +JosephForwardAttenuatedProjectionImageFilter::VerifyInputInformation() const { Superclass::VerifyInputInformation(); using ImageBaseType = const itk::ImageBase; @@ -160,26 +204,12 @@ JosephForwardAttenuatedProjectionImageFilter +template void -JosephForwardAttenuatedProjectionImageFilter::BeforeThreadedGenerateData() +JosephForwardAttenuatedProjectionImageFilter::BeforeThreadedGenerateData() { - this->GetInterpolationWeightMultiplication().SetAttenuationMinusEmissionMapsPtrDiff( - this->GetInput(2)->GetBufferPointer() - this->GetInput(1)->GetBufferPointer()); - this->GetProjectedValueAccumulation().SetAttenuationVector( - this->GetInterpolationWeightMultiplication().GetAttenuationRay()); - this->GetSumAlongRay().SetAttenuationRayVector(this->GetInterpolationWeightMultiplication().GetAttenuationRay()); - this->GetSumAlongRay().SetAttenuationPixelVector(this->GetInterpolationWeightMultiplication().GetAttenuationPixel()); - this->GetProjectedValueAccumulation().SetEx1(this->GetInterpolationWeightMultiplication().GetEx1()); - this->GetSumAlongRay().SetEx1(this->GetInterpolationWeightMultiplication().GetEx1()); + this->m_AttenuationMinusEmissionMapsPtrDiff = + this->GetInput(2)->GetBufferPointer() - this->GetInput(1)->GetBufferPointer(); } } // end namespace rtk diff --git a/include/rtkJosephForwardProjectionImageFilter.h b/include/rtkJosephForwardProjectionImageFilter.h index 172e2c81f..6447a6bd5 100644 --- a/include/rtkJosephForwardProjectionImageFilter.h +++ b/include/rtkJosephForwardProjectionImageFilter.h @@ -29,123 +29,6 @@ #include namespace rtk { -namespace Functor -{ -/** \class InterpolationWeightMultiplication - * \brief Function to multiply the interpolation weights with the projected - * volume values. - * - * \author Simon Rit - * - * \ingroup RTK Functions - */ -template -class ITK_TEMPLATE_EXPORT InterpolationWeightMultiplication -{ -public: - InterpolationWeightMultiplication() = default; - ~InterpolationWeightMultiplication() = default; - bool - operator!=(const InterpolationWeightMultiplication &) const - { - return false; - } - bool - operator==(const InterpolationWeightMultiplication & other) const - { - return !(*this != other); - } - - inline TOutput - operator()(const ThreadIdType itkNotUsed(threadId), - const double itkNotUsed(stepLengthInVoxel), - const TCoordinateType weight, - const TInput * p, - const int i) const - { - return weight * p[i]; - } -}; - -/** \class SumAlongRay - * \brief Function to compute the attenuation correction on the projection. - * - * \author Antoine Robert - * - * \ingroup RTK Functions - */ -template -class ITK_TEMPLATE_EXPORT SumAlongRay -{ -public: - using VectorType = itk::Vector; - - SumAlongRay() = default; - ~SumAlongRay() = default; - bool - operator!=(const SumAlongRay &) const - { - return false; - } - bool - operator==(const SumAlongRay & other) const - { - return !(*this != other); - } - - inline void - operator()(const ThreadIdType itkNotUsed(threadId), - TOutput & sumValue, - const TInput volumeValue, - const VectorType & itkNotUsed(stepInMM)) - { - sumValue += static_cast(volumeValue); - } -}; - -/** \class ProjectedValueAccumulation - * \brief Function to accumulate the ray casting on the projection. - * - * \author Simon Rit - * - * \ingroup RTK Functions - */ -template -class ITK_TEMPLATE_EXPORT ProjectedValueAccumulation -{ -public: - using VectorType = itk::Vector; - - ProjectedValueAccumulation() = default; - ~ProjectedValueAccumulation() = default; - bool - operator!=(const ProjectedValueAccumulation &) const - { - return false; - } - bool - operator==(const ProjectedValueAccumulation & other) const - { - return !(*this != other); - } - - inline void - operator()(const ThreadIdType itkNotUsed(threadId), - const TInput & input, - TOutput & output, - const TOutput & rayCastValue, - const VectorType & stepInMM, - const VectorType & itkNotUsed(source), - const VectorType & itkNotUsed(sourceToPixel), - const VectorType & itkNotUsed(nearestPoint), - const VectorType & itkNotUsed(farthestPoint)) const - { - output = input + rayCastValue * stepInMM.GetNorm(); - } -}; - -} // end namespace Functor - /** \class JosephForwardProjectionImageFilter * \brief Joseph forward projection. @@ -161,15 +44,7 @@ class ITK_TEMPLATE_EXPORT ProjectedValueAccumulation * * \ingroup RTK Projector */ - -template ::ValueType>, - class TProjectedValueAccumulation = - Functor::ProjectedValueAccumulation, - class TSumAlongRay = Functor::SumAlongRay> +template class ITK_TEMPLATE_EXPORT JosephForwardProjectionImageFilter : public ForwardProjectionImageFilter { @@ -185,78 +60,99 @@ class ITK_TEMPLATE_EXPORT JosephForwardProjectionImageFilter using OutputPixelType = typename TOutputImage::PixelType; using OutputImageRegionType = typename TOutputImage::RegionType; using CoordinateType = double; + using WeightCoordinateType = typename itk::PixelTraits::ValueType; using VectorType = itk::Vector; using TClipImageType = itk::Image; using TClipImagePointer = typename TClipImageType::Pointer; + /** \brief Function to multiply the interpolation weights with the projected + * volume values. + * + * \author Simon Rit + */ + using InterpolationWeightMultiplicationFunc = + std::function; + + /** \brief Function to compute the attenuation correction on the projection. + * + * \author Antoine Robert + */ + using SumAlongRayFunc = + std::function; + + /** \brief Function to accumulate the ray casting on the projection. + * + * \author Simon Rit + */ + using ProjectedValueAccumulationFunc = std::function; + /** Method for creation through the object factory. */ itkNewMacro(Self); /** Run-time type information (and related methods). */ itkOverrideGetNameOfClassMacro(JosephForwardProjectionImageFilter); - /** Get/Set the functor that is used to multiply each interpolation value with a volume value */ - TInterpolationWeightMultiplication & + /** Get/Set the lambda function that is used to multiply each interpolation value with a volume value */ + InterpolationWeightMultiplicationFunc & GetInterpolationWeightMultiplication() { return m_InterpolationWeightMultiplication; } - const TInterpolationWeightMultiplication & + const InterpolationWeightMultiplicationFunc & GetInterpolationWeightMultiplication() const { return m_InterpolationWeightMultiplication; } void - SetInterpolationWeightMultiplication(const TInterpolationWeightMultiplication & _arg) + SetInterpolationWeightMultiplication(const InterpolationWeightMultiplicationFunc & _arg) { - if (m_InterpolationWeightMultiplication != _arg) - { - m_InterpolationWeightMultiplication = _arg; - this->Modified(); - } + m_InterpolationWeightMultiplication = _arg; + this->Modified(); } - /** Get/Set the functor that is used to accumulate values in the projection image after the ray + /** Get/Set the lambda function that is used to accumulate values in the projection image after the ray * casting has been performed. */ - TProjectedValueAccumulation & + ProjectedValueAccumulationFunc & GetProjectedValueAccumulation() { return m_ProjectedValueAccumulation; } - const TProjectedValueAccumulation & + const ProjectedValueAccumulationFunc & GetProjectedValueAccumulation() const { return m_ProjectedValueAccumulation; } void - SetProjectedValueAccumulation(const TProjectedValueAccumulation & _arg) + SetProjectedValueAccumulation(const ProjectedValueAccumulationFunc & _arg) { - if (m_ProjectedValueAccumulation != _arg) - { - m_ProjectedValueAccumulation = _arg; - this->Modified(); - } + m_ProjectedValueAccumulation = _arg; + this->Modified(); } - /** Get/Set the functor that is used to compute the sum along the ray*/ - TSumAlongRay & + /** Get/Set the lambda function that is used to compute the sum along the ray */ + SumAlongRayFunc & GetSumAlongRay() { return m_SumAlongRay; } - const TSumAlongRay & + const SumAlongRayFunc & GetSumAlongRay() const { return m_SumAlongRay; } void - SetSumAlongRay(const TSumAlongRay & _arg) + SetSumAlongRay(const SumAlongRayFunc & _arg) { - if (m_SumAlongRay != _arg) - { - m_SumAlongRay = _arg; - this->Modified(); - } + m_SumAlongRay = _arg; + this->Modified(); } /** Set/Get the inferior clip image. Each pixel of the image @@ -342,12 +238,14 @@ class ITK_TEMPLATE_EXPORT JosephForwardProjectionImageFilter const double maxy); private: - // Functors - TInterpolationWeightMultiplication m_InterpolationWeightMultiplication; - TProjectedValueAccumulation m_ProjectedValueAccumulation; - TSumAlongRay m_SumAlongRay; - double m_InferiorClip{ 0. }; - double m_SuperiorClip{ 1. }; + // lambdas or std::functions + // MUST BE initiated in constructor, or custom functions + // must be set by means of Set methods in application + InterpolationWeightMultiplicationFunc m_InterpolationWeightMultiplication; + SumAlongRayFunc m_SumAlongRay; + ProjectedValueAccumulationFunc m_ProjectedValueAccumulation; + double m_InferiorClip{ 0. }; + double m_SuperiorClip{ 1. }; }; } // end namespace rtk diff --git a/include/rtkJosephForwardProjectionImageFilter.hxx b/include/rtkJosephForwardProjectionImageFilter.hxx index eebc5357d..c6ccddfee 100644 --- a/include/rtkJosephForwardProjectionImageFilter.hxx +++ b/include/rtkJosephForwardProjectionImageFilter.hxx @@ -31,32 +31,37 @@ namespace rtk { -template -JosephForwardProjectionImageFilter::JosephForwardProjectionImageFilter() +template +JosephForwardProjectionImageFilter::JosephForwardProjectionImageFilter() { + this->m_InterpolationWeightMultiplication = [](const ThreadIdType itkNotUsed(threadId), + double itkNotUsed(stepLengthInVoxel), + const WeightCoordinateType weight, + const InputPixelType * p, + int i) -> OutputPixelType { return weight * p[i]; }; + + this->m_SumAlongRay = + [](const ThreadIdType, OutputPixelType & sumValue, const InputPixelType volumeValue, const VectorType &) { + sumValue += static_cast(volumeValue); + }; + + this->m_ProjectedValueAccumulation = [](const ThreadIdType, + const InputPixelType & input, + OutputPixelType & output, + const OutputPixelType & rayCastValue, + const VectorType & stepInMM, + const VectorType &, + const VectorType &, + const VectorType &, + const VectorType &) { output = input + rayCastValue * stepInMM.GetNorm(); }; + this->DynamicMultiThreadingOff(); } -template +template void -JosephForwardProjectionImageFilter::GenerateInputRequestedRegion() +JosephForwardProjectionImageFilter::GenerateInputRequestedRegion() { Superclass::GenerateInputRequestedRegion(); @@ -78,17 +83,9 @@ JosephForwardProjectionImageFilter +template void -JosephForwardProjectionImageFilter::VerifyInputInformation() const +JosephForwardProjectionImageFilter::VerifyInputInformation() const { using ImageBaseType = const itk::ImageBase; @@ -169,19 +166,11 @@ JosephForwardProjectionImageFilter +template void -JosephForwardProjectionImageFilter::ThreadedGenerateData(const OutputImageRegionType & - outputRegionForThread, - ThreadIdType threadId) +JosephForwardProjectionImageFilter::ThreadedGenerateData( + const OutputImageRegionType & outputRegionForThread, + ThreadIdType threadId) { const unsigned int Dimension = TInputImage::ImageDimension; int offsets[3]; @@ -423,30 +412,18 @@ JosephForwardProjectionImageFilter -typename JosephForwardProjectionImageFilter::OutputPixelType -JosephForwardProjectionImageFilter::BilinearInterpolation(const ThreadIdType threadId, - const double stepLengthInVoxel, - const InputPixelType * pxiyi, - const InputPixelType * pxsyi, - const InputPixelType * pxiys, - const InputPixelType * pxsys, - const CoordinateType x, - const CoordinateType y, - const int ox, - const int oy) +template +typename JosephForwardProjectionImageFilter::OutputPixelType +JosephForwardProjectionImageFilter::BilinearInterpolation(const ThreadIdType threadId, + const double stepLengthInVoxel, + const InputPixelType * pxiyi, + const InputPixelType * pxsyi, + const InputPixelType * pxiys, + const InputPixelType * pxsys, + const CoordinateType x, + const CoordinateType y, + const int ox, + const int oy) { int ix = itk::Math::floor(x); int iy = itk::Math::floor(y); @@ -473,34 +450,23 @@ JosephForwardProjectionImageFilter -typename JosephForwardProjectionImageFilter::OutputPixelType -JosephForwardProjectionImageFilter::BilinearInterpolationOnBorders(const ThreadIdType threadId, - const double stepLengthInVoxel, - const InputPixelType * pxiyi, - const InputPixelType * pxsyi, - const InputPixelType * pxiys, - const InputPixelType * pxsys, - const CoordinateType x, - const CoordinateType y, - const int ox, - const int oy, - const CoordinateType minx, - const CoordinateType miny, - const CoordinateType maxx, - const CoordinateType maxy) +template +typename JosephForwardProjectionImageFilter::OutputPixelType +JosephForwardProjectionImageFilter::BilinearInterpolationOnBorders( + const ThreadIdType threadId, + const double stepLengthInVoxel, + const InputPixelType * pxiyi, + const InputPixelType * pxsyi, + const InputPixelType * pxiys, + const InputPixelType * pxsys, + const CoordinateType x, + const CoordinateType y, + const int ox, + const int oy, + const CoordinateType minx, + const CoordinateType miny, + const CoordinateType maxx, + const CoordinateType maxy) { int ix = itk::Math::floor(x); int iy = itk::Math::floor(y); diff --git a/include/rtkMaximumIntensityProjectionImageFilter.h b/include/rtkMaximumIntensityProjectionImageFilter.h index 96d3a6830..ecfdf7b95 100644 --- a/include/rtkMaximumIntensityProjectionImageFilter.h +++ b/include/rtkMaximumIntensityProjectionImageFilter.h @@ -23,95 +23,6 @@ namespace rtk { -namespace Functor -{ - -/** \class MaximumIntensityAlongRay - * \brief Function to compute the maximum intensity (MIP) value along the ray projection. - * - * \author Mikhail Polkovnikov - * - * \ingroup RTK Functions - */ -template -class ITK_TEMPLATE_EXPORT MaximumIntensityAlongRay -{ -public: - using VectorType = itk::Vector; - - MaximumIntensityAlongRay() = default; - ~MaximumIntensityAlongRay() = default; - bool - operator!=(const MaximumIntensityAlongRay &) const - { - return false; - } - bool - operator==(const MaximumIntensityAlongRay & other) const - { - return !(*this != other); - } - - inline void - operator()(const ThreadIdType itkNotUsed(threadId), - TOutput & mipValue, - const TInput volumeValue, - const VectorType & itkNotUsed(stepInMM)) - { - TOutput tmp = static_cast(volumeValue); - if (tmp > mipValue) - { - mipValue = tmp; - } - } -}; - -/** \class MaximumIntensityProjectedValueAccumulation - * \brief Function to calculate maximum intensity step along the ray projection. - * - * \author Mikhail Polkovnikov - * - * \ingroup RTK Functions - */ -template -class ITK_TEMPLATE_EXPORT MaximumIntensityProjectedValueAccumulation -{ -public: - using VectorType = itk::Vector; - - bool - operator!=(const MaximumIntensityProjectedValueAccumulation &) const - { - return false; - } - bool - operator==(const MaximumIntensityProjectedValueAccumulation & other) const - { - return !(*this != other); - } - - inline void - operator()(const ThreadIdType itkNotUsed(threadId), - const TInput & input, - TOutput & output, - const TOutput & rayCastValue, - const VectorType & stepInMM, - const VectorType & itkNotUsed(source), - const VectorType & itkNotUsed(sourceToPixel), - const VectorType & itkNotUsed(nearestPoint), - const VectorType & itkNotUsed(farthestPoint)) const - { - TOutput tmp = static_cast(input); - if (tmp < rayCastValue) - { - tmp = rayCastValue; - } - output = tmp * stepInMM.GetNorm(); - } -}; - -} // end namespace Functor - /** \class MaximumIntensityProjectionImageFilter * \brief MIP filter. @@ -119,27 +30,16 @@ class ITK_TEMPLATE_EXPORT MaximumIntensityProjectedValueAccumulation * Performs a MIP forward projection, i.e. calculation of a maximum intensity * step along the x-ray line. * + * \test rtkmaximumintensityprojectiontest.cxx + * * \author Mikhail Polkovnikov * * \ingroup RTK Projector */ -template ::ValueType>, - class TProjectedValueAccumulation = - Functor::MaximumIntensityProjectedValueAccumulation, - class TSumAlongRay = - Functor::MaximumIntensityAlongRay> +template class ITK_TEMPLATE_EXPORT MaximumIntensityProjectionImageFilter - : public JosephForwardProjectionImageFilter + : public JosephForwardProjectionImageFilter { public: ITK_DISALLOW_COPY_AND_MOVE(MaximumIntensityProjectionImageFilter); @@ -147,6 +47,13 @@ class ITK_TEMPLATE_EXPORT MaximumIntensityProjectionImageFilter /** Standard class type alias. */ using Self = MaximumIntensityProjectionImageFilter; using Pointer = itk::SmartPointer; + using Superclass = JosephForwardProjectionImageFilter; + using ConstPointer = itk::SmartPointer; + using InputPixelType = typename TInputImage::PixelType; + using OutputPixelType = typename TOutputImage::PixelType; + using CoordinateType = double; + using WeightCoordinateType = typename itk::PixelTraits::ValueType; + using VectorType = itk::Vector; /** Method for creation through the object factory. */ itkNewMacro(Self); @@ -155,9 +62,13 @@ class ITK_TEMPLATE_EXPORT MaximumIntensityProjectionImageFilter itkOverrideGetNameOfClassMacro(MaximumIntensityProjectionImageFilter); protected: - MaximumIntensityProjectionImageFilter() = default; + MaximumIntensityProjectionImageFilter(); ~MaximumIntensityProjectionImageFilter() override = default; }; } // end namespace rtk +#ifndef ITK_MANUAL_INSTANTIATION +# include "rtkMaximumIntensityProjectionImageFilter.hxx" +#endif + #endif diff --git a/include/rtkMaximumIntensityProjectionImageFilter.hxx b/include/rtkMaximumIntensityProjectionImageFilter.hxx new file mode 100644 index 000000000..52c7bb8a0 --- /dev/null +++ b/include/rtkMaximumIntensityProjectionImageFilter.hxx @@ -0,0 +1,60 @@ +/*========================================================================= + * + * Copyright RTK Consortium + * + * 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 + * + * https://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 rtkMaximumIntensityProjectionImageFilter_hxx +#define rtkMaximumIntensityProjectionImageFilter_hxx + +namespace rtk +{ + +template +MaximumIntensityProjectionImageFilter::MaximumIntensityProjectionImageFilter() + : JosephForwardProjectionImageFilter() +{ + auto sumAlongRayFunc = + [](const ThreadIdType, OutputPixelType & mipValue, const InputPixelType volumeValue, const VectorType &) { + OutputPixelType tmp = static_cast(volumeValue); + if (tmp > mipValue) + { + mipValue = tmp; + } + }; + this->SetSumAlongRay(sumAlongRayFunc); + + auto projectedValueAccumulationFunc = [](const ThreadIdType, + const InputPixelType & input, + OutputPixelType & output, + const OutputPixelType & rayCastValue, + const VectorType & stepInMM, + const VectorType &, + const VectorType &, + const VectorType &, + const VectorType &) { + OutputPixelType tmp = static_cast(input); + if (tmp < rayCastValue) + { + tmp = rayCastValue; + } + output = tmp * stepInMM.GetNorm(); + }; + this->SetProjectedValueAccumulation(projectedValueAccumulationFunc); +} + +} // end namespace rtk + +#endif diff --git a/test/rtkforwardprojectiontest.cxx b/test/rtkforwardprojectiontest.cxx index 824ce8502..1fb0142cf 100644 --- a/test/rtkforwardprojectiontest.cxx +++ b/test/rtkforwardprojectiontest.cxx @@ -98,6 +98,17 @@ main(int, char **) jfp->SetInput(projInput->GetOutput()); jfp->SetInput(1, volInput->GetOutput()); +#ifndef USE_CUDA + // custom SumAlongRay lambda function (works ONLY without CUDA support) + JFPType::SumAlongRayFunc sumFuncTest = [](const itk::ThreadIdType, + JFPType::OutputPixelType & sum, + const JFPType::InputPixelType in, + const JFPType::VectorType &) -> void { + sum += static_cast(in); + }; + jfp->SetSumAlongRay(sumFuncTest); +#endif + // Ray Box Intersection filter (reference) using RBIType = rtk::RayBoxIntersectionImageFilter; #ifdef USE_CUDA diff --git a/test/rtkmaximumintensityprojectiontest.cxx b/test/rtkmaximumintensityprojectiontest.cxx index 0e8ce027b..c11a5a79b 100644 --- a/test/rtkmaximumintensityprojectiontest.cxx +++ b/test/rtkmaximumintensityprojectiontest.cxx @@ -2,6 +2,7 @@ #include "rtkThreeDCircularProjectionGeometryXMLFile.h" #include "rtkRayBoxIntersectionImageFilter.h" #include "rtkConstantImageSource.h" +#include "rtkJosephForwardProjectionImageFilter.h" #include "rtkMaximumIntensityProjectionImageFilter.h" #include #include @@ -70,21 +71,55 @@ main(int, char **) projInput->SetConstant(0.); projInput->Update(); - // MIP Forward Projection filter - using MIPType = rtk::MaximumIntensityProjectionImageFilter; - MIPType::Pointer mipfp = MIPType::New(); - mipfp->SetInput(projInput->GetOutput()); - mipfp->SetInput(1, volInput->GetOutput()); + // Test-1 + // Joseph Forward Projection filter for MIP projection + using JFPType = rtk::JosephForwardProjectionImageFilter; + JFPType::Pointer jfp = JFPType::New(); + jfp->SetInput(projInput->GetOutput()); + jfp->SetInput(1, volInput->GetOutput()); + + // Function to compute the maximum intensity (MIP) value along the ray projection. + JFPType::SumAlongRayFunc sumAlongFunc = [](const itk::ThreadIdType, + JFPType::OutputPixelType & mipValue, + const JFPType::InputPixelType volumeValue, + const JFPType::VectorType &) -> void { + JFPType::OutputPixelType tmp = static_cast(volumeValue); + if (tmp > mipValue) + { + mipValue = tmp; + } + }; + jfp->SetSumAlongRay(sumAlongFunc); + + // Performs a MIP forward projection, i.e. calculation of a maximum intensity + // step along the x-ray line. + JFPType::ProjectedValueAccumulationFunc projAccumFunc = [](const itk::ThreadIdType itkNotUsed(threadId), + const JFPType::InputPixelType & input, + JFPType::OutputPixelType & output, + const JFPType::OutputPixelType & rayCastValue, + const JFPType::VectorType & stepInMM, + const JFPType::VectorType & itkNotUsed(source), + const JFPType::VectorType & itkNotUsed(sourceToPixel), + const JFPType::VectorType & itkNotUsed(nearestPoint), + const JFPType::VectorType & itkNotUsed(farthestPoint)) { + JFPType::OutputPixelType tmp = static_cast(input); + if (tmp < rayCastValue) + { + tmp = rayCastValue; + } + output = tmp * stepInMM.GetNorm(); + }; + jfp->SetProjectedValueAccumulation(projAccumFunc); using GeometryType = rtk::ThreeDCircularProjectionGeometry; GeometryType::Pointer geometry = GeometryType::New(); geometry->AddProjection(700, 800, 0); - mipfp->SetGeometry(geometry); - mipfp->Update(); + jfp->SetGeometry(geometry); + jfp->Update(); using ConstIteratorType = itk::ImageRegionConstIterator; - ConstIteratorType inputIt(mipfp->GetOutput(), mipfp->GetOutput()->GetRequestedRegion()); + ConstIteratorType inputIt(jfp->GetOutput(), jfp->GetOutput()->GetRequestedRegion()); inputIt.GoToBegin(); @@ -102,6 +137,37 @@ main(int, char **) { return EXIT_FAILURE; } - std::cout << "\n\nTest PASSED! " << std::endl; + std::cout << "\n\nTest-1 PASSED! " << std::endl; + + // Test-2 + // Maximum Intensity Projection (MIP) filter, derived from Joseph Forward Projection filter + using MIPType = rtk::MaximumIntensityProjectionImageFilter; + MIPType::Pointer mipfp = MIPType::New(); + mipfp->SetInput(projInput->GetOutput()); + mipfp->SetInput(1, volInput->GetOutput()); + + mipfp->SetGeometry(geometry); + mipfp->Update(); + + inputIt = ConstIteratorType(mipfp->GetOutput(), mipfp->GetOutput()->GetRequestedRegion()); + + inputIt.GoToBegin(); + + res = false; + while (!inputIt.IsAtEnd()) + { + OutputPixelType pixel = inputIt.Get(); + if (pixel < 4. || pixel > 4.25) + { + res = true; + } + ++inputIt; + } + if (res) + { + return EXIT_FAILURE; + } + + std::cout << "\n\nTest-2 PASSED! " << std::endl; return EXIT_SUCCESS; } diff --git a/wrapping/rtkJosephForwardAttenuatedProjectionImageFilter.wrap b/wrapping/rtkJosephForwardAttenuatedProjectionImageFilter.wrap index 89457c8d3..df8e6dfec 100644 --- a/wrapping/rtkJosephForwardAttenuatedProjectionImageFilter.wrap +++ b/wrapping/rtkJosephForwardAttenuatedProjectionImageFilter.wrap @@ -1,28 +1,3 @@ -set(WRAPPER_AUTO_INCLUDE_HEADERS OFF) -itk_wrap_named_class("rtk::Functor::InterpolationWeightMultiplicationAttenuated" "rtkFunctorInterpolationWeightMultiplicationAttenuatedBackProjectionAttenuated") - foreach(t ${WRAP_ITK_REAL}) - itk_wrap_template("${ITKM_${t}}${ITKM_${t}}${ITKM_${t}}" "${ITKT_${t}}, ${ITKT_${t}}, ${ITKT_${t}}") - endforeach() -itk_end_wrap_class() -itk_wrap_named_class("rtk::Functor::ProjectedValueAccumulationAttenuated" "rtkProjectedValueAccumulationAttenuated") - foreach(t ${WRAP_ITK_REAL}) - itk_wrap_template("${ITKM_${t}}${ITKM_${t}}" "${ITKT_${t}}, ${ITKT_${t}}") - endforeach() -itk_end_wrap_class() -itk_wrap_named_class("rtk::Functor::ComputeAttenuationCorrection" "rtkComputeAttenuationCorrection") - foreach(t ${WRAP_ITK_REAL}) - itk_wrap_template("${ITKM_${t}}${ITKM_${t}}" "${ITKT_${t}}, ${ITKT_${t}}") - endforeach() -itk_end_wrap_class() -set(WRAPPER_AUTO_INCLUDE_HEADERS ON) - -itk_wrap_class("rtk::JosephForwardProjectionImageFilter" POINTER) - foreach(t ${WRAP_ITK_REAL}) - itk_wrap_template("I${ITKM_${t}}3I${ITKM_${t}}3SWM${ITKM_${t}}D${ITKM_${t}}IPC" - "itk::Image<${ITKT_${t}}, 3>, itk::Image< ${ITKT_${t}}, 3>, rtk::Functor::InterpolationWeightMultiplicationAttenuated<${ITKT_${t}}, ${ITKT_${t}}, ${ITKT_${t}}>, rtk::Functor::ProjectedValueAccumulationAttenuated<${ITKT_${t}}, ${ITKT_${t}}>, rtk::Functor::ComputeAttenuationCorrection<${ITKT_${t}}, ${ITKT_${t}}>") - endforeach() -itk_end_wrap_class() - itk_wrap_class("rtk::JosephForwardAttenuatedProjectionImageFilter" POINTER) foreach(t ${WRAP_ITK_REAL}) itk_wrap_template("I${ITKM_${t}}3I${ITKM_${t}}3" "itk::Image<${ITKT_${t}}, 3>, itk::Image<${ITKT_${t}}, 3>") diff --git a/wrapping/rtkJosephForwardProjectionImageFilter.wrap b/wrapping/rtkJosephForwardProjectionImageFilter.wrap index 873afb552..d6820105f 100644 --- a/wrapping/rtkJosephForwardProjectionImageFilter.wrap +++ b/wrapping/rtkJosephForwardProjectionImageFilter.wrap @@ -2,11 +2,4 @@ itk_wrap_class("rtk::JosephForwardProjectionImageFilter" POINTER) foreach(t ${WRAP_ITK_REAL}) itk_wrap_template("I${ITKM_${t}}3I${ITKM_${t}}3" "itk::Image<${ITKT_${t}}, 3>, itk::Image<${ITKT_${t}}, 3>") endforeach() - - foreach(vt ${WRAP_ITK_VECTOR_REAL}) - set(vectorComponents 2 3 4 5) - foreach(nmat ${vectorComponents}) - itk_wrap_template("I${ITKM_${vt}${nmat}}3I${ITKM_${vt}${nmat}}3" "itk::Image<${ITKT_${vt}${nmat}},3>, itk::Image<${ITKT_${vt}${nmat}},3>") - endforeach() - endforeach() itk_end_wrap_class() diff --git a/wrapping/rtkMaximumIntensityProjectionImageFilter.wrap b/wrapping/rtkMaximumIntensityProjectionImageFilter.wrap index b4b1acce6..93b4d4317 100644 --- a/wrapping/rtkMaximumIntensityProjectionImageFilter.wrap +++ b/wrapping/rtkMaximumIntensityProjectionImageFilter.wrap @@ -1,29 +1,4 @@ -set(WRAPPER_AUTO_INCLUDE_HEADERS OFF) -itk_wrap_named_class("rtk::Functor::InterpolationWeightMultiplication" "rtkInterpolationWeightMultiplication") - foreach(t ${WRAP_ITK_REAL}) - itk_wrap_template("${ITKM_${t}}${ITKM_${t}}" "${ITKT_${t}}, ${ITKT_${t}}") - endforeach() -itk_end_wrap_class() -itk_wrap_named_class("rtk::Functor::MaximumIntensityAlongRay" "rtkMaximumIntensityAlongRay") - foreach(t ${WRAP_ITK_REAL}) - itk_wrap_template("${ITKM_${t}}${ITKM_${t}}" "${ITKT_${t}}, ${ITKT_${t}}") - endforeach() -itk_end_wrap_class() -itk_wrap_named_class("rtk::Functor::MaximumIntensityProjectedValueAccumulation" "rtkMaximumIntensityProjectedValueAccumulation") - foreach(t ${WRAP_ITK_REAL}) - itk_wrap_template("${ITKM_${t}}${ITKM_${t}}" "${ITKT_${t}}, ${ITKT_${t}}") - endforeach() -itk_end_wrap_class() -set(WRAPPER_AUTO_INCLUDE_HEADERS ON) - -itk_wrap_class("rtk::JosephForwardProjectionImageFilter" POINTER) - foreach(t ${WRAP_ITK_REAL}) - itk_wrap_template("I${ITKM_${t}}3I${ITKM_${t}}3FWMI${ITKM_${t}}3I${ITKM_${t}}3FVAI${ITKM_${t}}3I${ITKM_${t}}3FIAI${ITKM_${t}}3I${ITKM_${t}}3" - "itk::Image<${ITKT_${t}}, 3>, itk::Image< ${ITKT_${t}}, 3>, rtk::Functor::InterpolationWeightMultiplication<${ITKT_${t}}, ${ITKT_${t}}>, rtk::Functor::MaximumIntensityProjectedValueAccumulation<${ITKT_${t}}, ${ITKT_${t}}>, rtk::Functor::MaximumIntensityAlongRay<${ITKT_${t}}, ${ITKT_${t}}>") - endforeach() -itk_end_wrap_class() - -itk_wrap_class("rtk::MaximumIntensityProjectionImageFilter" POINTER) +itk_wrap_class("rtk::MaximumIntensityProjectionImageFilter" POINTER) foreach(t ${WRAP_ITK_REAL}) itk_wrap_template("I${ITKM_${t}}3I${ITKM_${t}}3" "itk::Image<${ITKT_${t}}, 3>, itk::Image<${ITKT_${t}}, 3>") endforeach() From a7eea4bf62dd724156d547cff133836f026981ee Mon Sep 17 00:00:00 2001 From: Mikhail Polkovnikov Date: Sat, 8 Mar 2025 17:20:33 +0300 Subject: [PATCH 2/3] STYLE: Verbose output in tests --- test/rtkforwardprojectiontest.cxx | 2 +- test/rtkmaximumintensityprojectiontest.cxx | 10 ++++++++-- 2 files changed, 9 insertions(+), 3 deletions(-) diff --git a/test/rtkforwardprojectiontest.cxx b/test/rtkforwardprojectiontest.cxx index 1fb0142cf..53f5a58bb 100644 --- a/test/rtkforwardprojectiontest.cxx +++ b/test/rtkforwardprojectiontest.cxx @@ -99,7 +99,7 @@ main(int, char **) jfp->SetInput(1, volInput->GetOutput()); #ifndef USE_CUDA - // custom SumAlongRay lambda function (works ONLY without CUDA support) + // Adding custom SumAlongRay lambda function doesn't work with CUDA filter JFPType::SumAlongRayFunc sumFuncTest = [](const itk::ThreadIdType, JFPType::OutputPixelType & sum, const JFPType::InputPixelType in, diff --git a/test/rtkmaximumintensityprojectiontest.cxx b/test/rtkmaximumintensityprojectiontest.cxx index c11a5a79b..c8d736cab 100644 --- a/test/rtkmaximumintensityprojectiontest.cxx +++ b/test/rtkmaximumintensityprojectiontest.cxx @@ -73,6 +73,7 @@ main(int, char **) // Test-1 // Joseph Forward Projection filter for MIP projection + std::cout << "\n\n****** Case 1: JosephForwardProjection filter with custom MIP functions ******" << std::endl; using JFPType = rtk::JosephForwardProjectionImageFilter; JFPType::Pointer jfp = JFPType::New(); jfp->SetInput(projInput->GetOutput()); @@ -137,10 +138,11 @@ main(int, char **) { return EXIT_FAILURE; } - std::cout << "\n\nTest-1 PASSED! " << std::endl; + std::cout << "\n\nTest PASSED! " << std::endl; // Test-2 // Maximum Intensity Projection (MIP) filter, derived from Joseph Forward Projection filter + std::cout << "\n\n****** Case 2: MaximumIntensityProjection (MIP) filter ******" << std::endl; using MIPType = rtk::MaximumIntensityProjectionImageFilter; MIPType::Pointer mipfp = MIPType::New(); mipfp->SetInput(projInput->GetOutput()); @@ -167,7 +169,11 @@ main(int, char **) { return EXIT_FAILURE; } + std::cout << "\n\nTest PASSED! " << std::endl; + + std::cout << "\n\n****** Compare two resulted images ******" << std::endl; + CheckImageQuality(jfp->GetOutput(), mipfp->GetOutput(), 0.001, 0.000001, 1.E+19); + std::cout << "\n\nImages are OK! " << std::endl; - std::cout << "\n\nTest-2 PASSED! " << std::endl; return EXIT_SUCCESS; } From e80a2d7fac5dadd89ab4d84721c7849770da26ed Mon Sep 17 00:00:00 2001 From: Mikhail Polkovnikov Date: Thu, 13 Mar 2025 19:07:31 +0300 Subject: [PATCH 3/3] STYLE: Fix review requests --- include/rtkJosephForwardAttenuatedProjectionImageFilter.hxx | 3 --- include/rtkJosephForwardProjectionImageFilter.h | 3 --- include/rtkMaximumIntensityProjectionImageFilter.h | 1 - test/rtkmaximumintensityprojectiontest.cxx | 6 +++--- 4 files changed, 3 insertions(+), 10 deletions(-) diff --git a/include/rtkJosephForwardAttenuatedProjectionImageFilter.hxx b/include/rtkJosephForwardAttenuatedProjectionImageFilter.hxx index bec8af91a..6c3389553 100644 --- a/include/rtkJosephForwardAttenuatedProjectionImageFilter.hxx +++ b/include/rtkJosephForwardAttenuatedProjectionImageFilter.hxx @@ -43,7 +43,6 @@ JosephForwardAttenuatedProjectionImageFilter::JosephF /** \brief Function to multiply the interpolation weights with the projected * volume values and attenuation map. * - * \author Antoine Robert */ auto interpolationWeightMultiplicationAttenuatedFunc = [this](const ThreadIdType threadId, double stepLengthInVoxel, @@ -60,7 +59,6 @@ JosephForwardAttenuatedProjectionImageFilter::JosephF /** \brief Function to compute the attenuation correction on the projection. * - * \author Antoine Robert */ auto computeAttenuationCorrectionFunc = [this](const ThreadIdType threadId, OutputPixelType & sumValue, @@ -86,7 +84,6 @@ JosephForwardAttenuatedProjectionImageFilter::JosephF /** \brief Function to accumulate the ray casting on the projection. * - * \author Antoine Robert */ auto projectedValueAccumulationAttenuatedFunc = [this](const ThreadIdType threadId, const InputPixelType & input, diff --git a/include/rtkJosephForwardProjectionImageFilter.h b/include/rtkJosephForwardProjectionImageFilter.h index 6447a6bd5..c1612cb96 100644 --- a/include/rtkJosephForwardProjectionImageFilter.h +++ b/include/rtkJosephForwardProjectionImageFilter.h @@ -68,21 +68,18 @@ class ITK_TEMPLATE_EXPORT JosephForwardProjectionImageFilter /** \brief Function to multiply the interpolation weights with the projected * volume values. * - * \author Simon Rit */ using InterpolationWeightMultiplicationFunc = std::function; /** \brief Function to compute the attenuation correction on the projection. * - * \author Antoine Robert */ using SumAlongRayFunc = std::function; /** \brief Function to accumulate the ray casting on the projection. * - * \author Simon Rit */ using ProjectedValueAccumulationFunc = std::function::ValueType; using VectorType = itk::Vector; /** Method for creation through the object factory. */ diff --git a/test/rtkmaximumintensityprojectiontest.cxx b/test/rtkmaximumintensityprojectiontest.cxx index c8d736cab..407451c58 100644 --- a/test/rtkmaximumintensityprojectiontest.cxx +++ b/test/rtkmaximumintensityprojectiontest.cxx @@ -117,7 +117,7 @@ main(int, char **) geometry->AddProjection(700, 800, 0); jfp->SetGeometry(geometry); - jfp->Update(); + TRY_AND_EXIT_ON_ITK_EXCEPTION(jfp->Update()); using ConstIteratorType = itk::ImageRegionConstIterator; ConstIteratorType inputIt(jfp->GetOutput(), jfp->GetOutput()->GetRequestedRegion()); @@ -149,7 +149,7 @@ main(int, char **) mipfp->SetInput(1, volInput->GetOutput()); mipfp->SetGeometry(geometry); - mipfp->Update(); + TRY_AND_EXIT_ON_ITK_EXCEPTION(mipfp->Update()); inputIt = ConstIteratorType(mipfp->GetOutput(), mipfp->GetOutput()->GetRequestedRegion()); @@ -173,7 +173,7 @@ main(int, char **) std::cout << "\n\n****** Compare two resulted images ******" << std::endl; CheckImageQuality(jfp->GetOutput(), mipfp->GetOutput(), 0.001, 0.000001, 1.E+19); - std::cout << "\n\nImages are OK! " << std::endl; + std::cout << "\n\nTest PASSED! " << std::endl; return EXIT_SUCCESS; }