ENH: default projectors for all iterative filters

SimonRit · Apr 14, 2020 · 2c8ae75 · 2c8ae75
1 parent bc6d102
commit 2c8ae75
Show file tree

Hide file tree

Showing 26 changed files with 119 additions and 377 deletions.
diff --git a/include/rtkADMMTotalVariationConeBeamReconstructionFilter.h b/include/rtkADMMTotalVariationConeBeamReconstructionFilter.h
@@ -174,15 +174,6 @@ class ADMMTotalVariationConeBeamReconstructionFilter
   using DisplacedDetectorFilterType = rtk::DisplacedDetectorImageFilter<TOutputImage>;
   using GatingWeightsFilterType = rtk::MultiplyByVectorImageFilter<TOutputImage>;
 
-  /** Pass the ForwardProjection filter to the conjugate gradient operator */
-  void
-  SetForwardProjectionFilter(ForwardProjectionType _arg) override;
-
-  /** Pass the backprojection filter to the conjugate gradient operator and to the back projection filter generating the
-   * B of AX=B */
-  void
-  SetBackProjectionFilter(BackProjectionType _arg) override;
-
   /** Pass the geometry to all filters needing it */
   itkSetObjectMacro(Geometry, ThreeDCircularProjectionGeometry);
 

diff --git a/include/rtkADMMTotalVariationConeBeamReconstructionFilter.hxx b/include/rtkADMMTotalVariationConeBeamReconstructionFilter.hxx
@@ -95,33 +95,6 @@ ADMMTotalVariationConeBeamReconstructionFilter<TOutputImage,
   m_DisplacedDetectorFilter->ReleaseDataFlagOn();
 }
 
-template <typename TOutputImage, typename TGradientOutputImage>
-void
-ADMMTotalVariationConeBeamReconstructionFilter<TOutputImage, TGradientOutputImage>::SetForwardProjectionFilter(
-  ForwardProjectionType _arg)
-{
-  if (_arg != this->GetForwardProjectionFilter())
-  {
-    Superclass::SetForwardProjectionFilter(_arg);
-    m_ForwardProjectionFilter = this->InstantiateForwardProjectionFilter(_arg);
-    m_CGOperator->SetForwardProjectionFilter(m_ForwardProjectionFilter);
-  }
-}
-
-template <typename TOutputImage, typename TGradientOutputImage>
-void
-ADMMTotalVariationConeBeamReconstructionFilter<TOutputImage, TGradientOutputImage>::SetBackProjectionFilter(
-  BackProjectionType _arg)
-{
-  if (_arg != this->GetBackProjectionFilter())
-  {
-    Superclass::SetBackProjectionFilter(_arg);
-    m_BackProjectionFilter = this->InstantiateBackProjectionFilter(_arg);
-    m_BackProjectionFilterForConjugateGradient = this->InstantiateBackProjectionFilter(_arg);
-    m_CGOperator->SetBackProjectionFilter(m_BackProjectionFilterForConjugateGradient);
-  }
-}
-
 template <typename TOutputImage, typename TGradientOutputImage>
 void
 ADMMTotalVariationConeBeamReconstructionFilter<TOutputImage, TGradientOutputImage>::SetBetaForCurrentIteration(int iter)
@@ -177,6 +150,19 @@ template <typename TOutputImage, typename TGradientOutputImage>
 void
 ADMMTotalVariationConeBeamReconstructionFilter<TOutputImage, TGradientOutputImage>::GenerateOutputInformation()
 {
+  // Set forward projection filter
+  m_ForwardProjectionFilter = this->InstantiateForwardProjectionFilter(this->m_CurrentForwardProjectionConfiguration);
+  // Pass the ForwardProjection filter to the conjugate gradient operator
+  m_CGOperator->SetForwardProjectionFilter(m_ForwardProjectionFilter);
+
+  // Set back projection filter
+  m_BackProjectionFilter = this->InstantiateBackProjectionFilter(this->m_CurrentBackProjectionConfiguration);
+  // Pass the backprojection filter to the conjugate gradient operator and to the back projection filter generating the
+  // B of AX=B
+  m_BackProjectionFilterForConjugateGradient =
+    this->InstantiateBackProjectionFilter(this->m_CurrentBackProjectionConfiguration);
+  m_CGOperator->SetBackProjectionFilter(m_BackProjectionFilterForConjugateGradient);
+
   // Set runtime connections
   m_GradientFilter1->SetInput(this->GetInput(0));
   m_ZeroMultiplyVolumeFilter->SetInput1(this->GetInput(0));

diff --git a/include/rtkADMMWaveletsConeBeamReconstructionFilter.h b/include/rtkADMMWaveletsConeBeamReconstructionFilter.h
@@ -176,15 +176,6 @@ class ADMMWaveletsConeBeamReconstructionFilter
   using ForwardProjectionType = typename Superclass::ForwardProjectionType;
   using BackProjectionType = typename Superclass::BackProjectionType;
 
-  /** Pass the ForwardProjection filter to the conjugate gradient operator */
-  void
-  SetForwardProjectionFilter(ForwardProjectionType _arg) override;
-
-  /** Pass the backprojection filter to the conjugate gradient operator and to the back projection filter generating the
-   * B of AX=B */
-  void
-  SetBackProjectionFilter(BackProjectionType _arg) override;
-
   /** Pass the geometry to all filters needing it */
   itkSetObjectMacro(Geometry, ThreeDCircularProjectionGeometry);
 

diff --git a/include/rtkADMMWaveletsConeBeamReconstructionFilter.hxx b/include/rtkADMMWaveletsConeBeamReconstructionFilter.hxx
@@ -72,31 +72,6 @@ ADMMWaveletsConeBeamReconstructionFilter<TOutputImage>::ADMMWaveletsConeBeamReco
   m_DisplacedDetectorFilter->ReleaseDataFlagOn();
 }
 
-template <typename TOutputImage>
-void
-ADMMWaveletsConeBeamReconstructionFilter<TOutputImage>::SetForwardProjectionFilter(ForwardProjectionType _arg)
-{
-  if (_arg != this->GetForwardProjectionFilter())
-  {
-    Superclass::SetForwardProjectionFilter(_arg);
-    m_ForwardProjectionFilterForConjugateGradient = this->InstantiateForwardProjectionFilter(_arg);
-    m_CGOperator->SetForwardProjectionFilter(m_ForwardProjectionFilterForConjugateGradient);
-  }
-}
-
-template <typename TOutputImage>
-void
-ADMMWaveletsConeBeamReconstructionFilter<TOutputImage>::SetBackProjectionFilter(BackProjectionType _arg)
-{
-  if (_arg != this->GetBackProjectionFilter())
-  {
-    Superclass::SetBackProjectionFilter(_arg);
-    m_BackProjectionFilter = this->InstantiateBackProjectionFilter(_arg);
-    m_BackProjectionFilterForConjugateGradient = this->InstantiateBackProjectionFilter(_arg);
-    m_CGOperator->SetBackProjectionFilter(m_BackProjectionFilterForConjugateGradient);
-  }
-}
-
 template <class TOutputImage>
 void
 ADMMWaveletsConeBeamReconstructionFilter<TOutputImage>::VerifyPreconditions() ITKv5_CONST
@@ -132,6 +107,20 @@ template <typename TOutputImage>
 void
 ADMMWaveletsConeBeamReconstructionFilter<TOutputImage>::GenerateOutputInformation()
 {
+  // Set forward projection filter
+  m_ForwardProjectionFilterForConjugateGradient =
+    this->InstantiateForwardProjectionFilter(this->m_CurrentForwardProjectionConfiguration);
+  // Pass the ForwardProjection filter to the conjugate gradient operator
+  m_CGOperator->SetForwardProjectionFilter(m_ForwardProjectionFilterForConjugateGradient);
+
+  // Set back projection filter
+  m_BackProjectionFilter = this->InstantiateBackProjectionFilter(this->m_CurrentBackProjectionConfiguration);
+  // Pass the backprojection filter to the conjugate gradient operator and to the back projection filter generating the
+  // B of AX=B
+  m_BackProjectionFilterForConjugateGradient =
+    this->InstantiateBackProjectionFilter(this->m_CurrentBackProjectionConfiguration);
+  m_CGOperator->SetBackProjectionFilter(m_BackProjectionFilterForConjugateGradient);
+
   // Set runtime connections
   m_ZeroMultiplyFilter->SetInput1(this->GetInput(0));
   m_CGOperator->SetInput(1, this->GetInput(1)); // The projections (the conjugate gradient operator needs them)

diff --git a/include/rtkConjugateGradientConeBeamReconstructionFilter.h b/include/rtkConjugateGradientConeBeamReconstructionFilter.h
@@ -170,15 +170,6 @@ class ConjugateGradientConeBeamReconstructionFilter : public IterativeConeBeamRe
   using ConstantImageSourceType = ConstantImageSource<TOutputImage>;
 #endif
 
-  /** Pass the ForwardProjection filter to the conjugate gradient operator */
-  void
-  SetForwardProjectionFilter(ForwardProjectionType _arg) override;
-
-  /** Pass the backprojection filter to the conjugate gradient operator and to the back projection filter generating the
-   * B of AX=B */
-  void
-  SetBackProjectionFilter(BackProjectionType _arg) override;
-
   /** Set the support mask, if any, for support constraint in reconstruction */
   void
   SetSupportMask(const TSingleComponentImage * SupportMask);

diff --git a/include/rtkConjugateGradientConeBeamReconstructionFilter.hxx b/include/rtkConjugateGradientConeBeamReconstructionFilter.hxx
@@ -123,33 +123,6 @@ ConjugateGradientConeBeamReconstructionFilter<TOutputImage, TSingleComponentImag
   return m_ConjugateGradientFilter->GetResidualCosts();
 }
 
-template <typename TOutputImage, typename TSingleComponentImage, typename TWeightsImage>
-void
-ConjugateGradientConeBeamReconstructionFilter<TOutputImage, TSingleComponentImage, TWeightsImage>::
-  SetForwardProjectionFilter(ForwardProjectionType _arg)
-{
-  if (_arg != this->GetForwardProjectionFilter())
-  {
-    Superclass::SetForwardProjectionFilter(_arg);
-    m_ForwardProjectionFilter = this->InstantiateForwardProjectionFilter(_arg);
-    m_CGOperator->SetForwardProjectionFilter(m_ForwardProjectionFilter);
-  }
-}
-
-template <typename TOutputImage, typename TSingleComponentImage, typename TWeightsImage>
-void
-ConjugateGradientConeBeamReconstructionFilter<TOutputImage, TSingleComponentImage, TWeightsImage>::
-  SetBackProjectionFilter(BackProjectionType _arg)
-{
-  if (_arg != this->GetBackProjectionFilter())
-  {
-    Superclass::SetBackProjectionFilter(_arg);
-    m_BackProjectionFilter = this->InstantiateBackProjectionFilter(_arg);
-    m_BackProjectionFilterForB = this->InstantiateBackProjectionFilter(_arg);
-    m_CGOperator->SetBackProjectionFilter(m_BackProjectionFilter);
-  }
-}
-
 template <typename TOutputImage, typename TSingleComponentImage, typename TWeightsImage>
 void
 ConjugateGradientConeBeamReconstructionFilter<TOutputImage, TSingleComponentImage, TWeightsImage>::VerifyPreconditions()
@@ -209,6 +182,18 @@ ConjugateGradientConeBeamReconstructionFilter<TOutputImage, TSingleComponentImag
   m_ConjugateGradientFilter->SetA(m_CGOperator.GetPointer());
   m_ConjugateGradientFilter->SetIterationCosts(m_IterationCosts);
 
+  // Set forward projection filter
+  m_ForwardProjectionFilter = this->InstantiateForwardProjectionFilter(this->m_CurrentForwardProjectionConfiguration);
+  // Pass the ForwardProjection filter to the conjugate gradient operator
+  m_CGOperator->SetForwardProjectionFilter(m_ForwardProjectionFilter);
+
+  // Set back projection filter
+  m_BackProjectionFilter = this->InstantiateBackProjectionFilter(this->m_CurrentBackProjectionConfiguration);
+  // Pass the backprojection filter to the conjugate gradient operator and to the back projection filter generating the
+  // B of AX=B
+  m_BackProjectionFilterForB = this->InstantiateBackProjectionFilter(this->m_CurrentBackProjectionConfiguration);
+  m_CGOperator->SetBackProjectionFilter(m_BackProjectionFilter);
+
   // Set runtime connections
   m_ConstantVolumeSource->SetInformationFromImage(this->GetInputVolume());
   m_CGOperator->SetInputProjectionStack(this->GetInputProjectionStack());

diff --git a/include/rtkFourDConjugateGradientConeBeamReconstructionFilter.h b/include/rtkFourDConjugateGradientConeBeamReconstructionFilter.h
@@ -163,14 +163,6 @@ class ITK_EXPORT FourDConjugateGradientConeBeamReconstructionFilter
   typename ProjectionStackType::ConstPointer
   GetInputProjectionStack();
 
-  /** Pass the ForwardProjection filter to the conjugate gradient operator */
-  void
-  SetForwardProjectionFilter(ForwardProjectionType _arg) override;
-
-  /** Pass the backprojection filter to the conjugate gradient operator and to the filter generating the B of AX=B */
-  void
-  SetBackProjectionFilter(BackProjectionType _arg) override;
-
   /** Pass the interpolation weights to subfilters */
   void
   SetWeights(const itk::Array2D<float> _arg);

diff --git a/include/rtkFourDConjugateGradientConeBeamReconstructionFilter.hxx b/include/rtkFourDConjugateGradientConeBeamReconstructionFilter.hxx
@@ -87,50 +87,6 @@ FourDConjugateGradientConeBeamReconstructionFilter<VolumeSeriesType, ProjectionS
   return static_cast<const ProjectionStackType *>(this->itk::ProcessObject::GetInput(1));
 }
 
-template <class VolumeSeriesType, class ProjectionStackType>
-void
-FourDConjugateGradientConeBeamReconstructionFilter<VolumeSeriesType, ProjectionStackType>::SetForwardProjectionFilter(
-  ForwardProjectionType _arg)
-{
-  if (_arg != this->GetForwardProjectionFilter())
-  {
-    Superclass::SetForwardProjectionFilter(_arg);
-    m_ForwardProjectionFilter = this->InstantiateForwardProjectionFilter(_arg);
-    m_CGOperator->SetForwardProjectionFilter(m_ForwardProjectionFilter);
-  }
-  if (_arg == 2) // The forward projection filter runs on GPU. It is most efficient to also run the interpolation on
-                 // GPU, and to use GPU constant image sources
-  {
-    m_CGOperator->SetUseCudaInterpolation(true);
-    m_CGOperator->SetUseCudaSources(true);
-  }
-}
-
-
-template <class VolumeSeriesType, class ProjectionStackType>
-void
-FourDConjugateGradientConeBeamReconstructionFilter<VolumeSeriesType, ProjectionStackType>::SetBackProjectionFilter(
-  BackProjectionType _arg)
-{
-  if (_arg != this->GetBackProjectionFilter())
-  {
-    Superclass::SetBackProjectionFilter(_arg);
-    m_BackProjectionFilter = this->InstantiateBackProjectionFilter(_arg);
-    m_CGOperator->SetBackProjectionFilter(m_BackProjectionFilter);
-
-    m_BackProjectionFilterForB = this->InstantiateBackProjectionFilter(_arg);
-    m_ProjStackToFourDFilter->SetBackProjectionFilter(m_BackProjectionFilterForB);
-  }
-  if (_arg == 2) // The back projection filter runs on GPU. It is most efficient to also run the splat on GPU, and to
-                 // use GPU constant image sources
-  {
-    m_CGOperator->SetUseCudaSplat(true);
-    m_CGOperator->SetUseCudaSources(true);
-    m_ProjStackToFourDFilter->SetUseCudaSplat(true);
-    m_ProjStackToFourDFilter->SetUseCudaSources(true);
-  }
-}
-
 template <class VolumeSeriesType, class ProjectionStackType>
 void
 FourDConjugateGradientConeBeamReconstructionFilter<VolumeSeriesType, ProjectionStackType>::SetWeights(
@@ -198,6 +154,35 @@ FourDConjugateGradientConeBeamReconstructionFilter<VolumeSeriesType, ProjectionS
   m_DisplacedDetectorFilter->SetDisable(m_DisableDisplacedDetectorFilter);
   m_CGOperator->SetDisableDisplacedDetectorFilter(m_DisableDisplacedDetectorFilter);
 
+  // Set forward projection filter
+  m_ForwardProjectionFilter = this->InstantiateForwardProjectionFilter(this->m_CurrentForwardProjectionConfiguration);
+  // Pass the ForwardProjection filter to the conjugate gradient operator
+  m_CGOperator->SetForwardProjectionFilter(m_ForwardProjectionFilter);
+  if (this->m_CurrentForwardProjectionConfiguration ==
+      ForwardProjectionType::FP_CUDARAYCAST) // The forward projection filter runs on GPU. It is most efficient to also
+                                             // run the interpolation on GPU, and to use GPU constant image sources
+  {
+    m_CGOperator->SetUseCudaInterpolation(true);
+    m_CGOperator->SetUseCudaSources(true);
+  }
+
+  // Set back projection filter
+  m_BackProjectionFilter = this->InstantiateBackProjectionFilter(this->m_CurrentBackProjectionConfiguration);
+  m_CGOperator->SetBackProjectionFilter(m_BackProjectionFilter);
+
+  m_BackProjectionFilterForB = this->InstantiateBackProjectionFilter(this->m_CurrentBackProjectionConfiguration);
+  // Pass the backprojection filter to the conjugate gradient operator and to the filter generating the B of AX=B
+  m_ProjStackToFourDFilter->SetBackProjectionFilter(m_BackProjectionFilterForB);
+  if (this->m_CurrentBackProjectionConfiguration ==
+      BackProjectionType::BP_CUDAVOXELBASED) // The back projection filter runs on GPU. It is most efficient to also run
+                                             // the splat on GPU, and to use GPU constant image sources
+  {
+    m_CGOperator->SetUseCudaSplat(true);
+    m_CGOperator->SetUseCudaSources(true);
+    m_ProjStackToFourDFilter->SetUseCudaSplat(true);
+    m_ProjStackToFourDFilter->SetUseCudaSources(true);
+  }
+
   // Have the last filter calculate its output information
   m_ConjugateGradientFilter->UpdateOutputInformation();
 

diff --git a/include/rtkFourDROOSTERConeBeamReconstructionFilter.h b/include/rtkFourDROOSTERConeBeamReconstructionFilter.h
@@ -310,14 +310,6 @@ class FourDROOSTERConeBeamReconstructionFilter
   using ForwardProjectionType = typename Superclass::ForwardProjectionType;
   using BackProjectionType = typename Superclass::BackProjectionType;
 
-  /** Pass the ForwardProjection filter to SingleProjectionToFourDFilter */
-  void
-  SetForwardProjectionFilter(ForwardProjectionType _arg) override;
-
-  /** Pass the backprojection filter to ProjectionStackToFourD*/
-  void
-  SetBackProjectionFilter(BackProjectionType _arg) override;
-
   /** Pass the interpolation weights to SingleProjectionToFourDFilter */
   virtual void
   SetWeights(const itk::Array2D<float> _arg);

diff --git a/include/rtkFourDROOSTERConeBeamReconstructionFilter.hxx b/include/rtkFourDROOSTERConeBeamReconstructionFilter.hxx
@@ -158,30 +158,6 @@ FourDROOSTERConeBeamReconstructionFilter<VolumeSeriesType, ProjectionStackType>:
   return static_cast<DVFSequenceImageType *>(this->itk::ProcessObject::GetInput("InverseDisplacementField"));
 }
 
-template <typename VolumeSeriesType, typename ProjectionStackType>
-void
-FourDROOSTERConeBeamReconstructionFilter<VolumeSeriesType, ProjectionStackType>::SetForwardProjectionFilter(
-  ForwardProjectionType _arg)
-{
-  if (_arg != this->GetForwardProjectionFilter())
-  {
-    Superclass::SetForwardProjectionFilter(_arg);
-    m_FourDCGFilter->SetForwardProjectionFilter(_arg);
-  }
-}
-
-template <typename VolumeSeriesType, typename ProjectionStackType>
-void
-FourDROOSTERConeBeamReconstructionFilter<VolumeSeriesType, ProjectionStackType>::SetBackProjectionFilter(
-  BackProjectionType _arg)
-{
-  if (_arg != this->GetBackProjectionFilter())
-  {
-    Superclass::SetBackProjectionFilter(_arg);
-    m_FourDCGFilter->SetBackProjectionFilter(_arg);
-  }
-}
-
 template <typename VolumeSeriesType, typename ProjectionStackType>
 void
 FourDROOSTERConeBeamReconstructionFilter<VolumeSeriesType, ProjectionStackType>::SetWeights(
@@ -266,6 +242,10 @@ FourDROOSTERConeBeamReconstructionFilter<VolumeSeriesType, ProjectionStackType>:
 {
   const int Dimension = VolumeType::ImageDimension;
 
+  // Set projection filters
+  m_FourDCGFilter->SetForwardProjectionFilter(this->m_CurrentForwardProjectionConfiguration);
+  m_FourDCGFilter->SetBackProjectionFilter(this->m_CurrentBackProjectionConfiguration);
+
   // The 4D conjugate gradient filter is the only part that must be in the pipeline
   // whatever was the user wants
   m_FourDCGFilter->SetInputVolumeSeries(this->GetInputVolumeSeries());

diff --git a/include/rtkFourDSARTConeBeamReconstructionFilter.h b/include/rtkFourDSARTConeBeamReconstructionFilter.h
@@ -196,14 +196,6 @@ class ITK_EXPORT FourDSARTConeBeamReconstructionFilter
   itkGetMacro(EnforcePositivity, bool);
   itkSetMacro(EnforcePositivity, bool);
 
-  /** Select the ForwardProjection filter */
-  void
-  SetForwardProjectionFilter(ForwardProjectionType _arg) override;
-
-  /** Select the backprojection filter */
-  void
-  SetBackProjectionFilter(BackProjectionType _arg) override;
-
   /** Pass the interpolation weights to subfilters */
   void
   SetWeights(const itk::Array2D<float> _arg);