RTK  2.6.0
Reconstruction Toolkit
rtkFourDROOSTERConeBeamReconstructionFilter.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright RTK Consortium
4  *
5  * Licensed under the Apache License, Version 2.0 (the "License");
6  * you may not use this file except in compliance with the License.
7  * You may obtain a copy of the License at
8  *
9  * https://www.apache.org/licenses/LICENSE-2.0.txt
10  *
11  * Unless required by applicable law or agreed to in writing, software
12  * distributed under the License is distributed on an "AS IS" BASIS,
13  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14  * See the License for the specific language governing permissions and
15  * limitations under the License.
16  *
17  *=========================================================================*/
18 #ifndef rtkFourDROOSTERConeBeamReconstructionFilter_h
19 #define rtkFourDROOSTERConeBeamReconstructionFilter_h
20 
23 #ifdef RTK_USE_CUDA
26 #else
29 #endif
35 
37 #include <itkSubtractImageFilter.h>
38 #include <itkAddImageFilter.h>
39 
40 #include <itkResampleImageFilter.h>
42 #include <itkIdentityTransform.h>
43 
44 
45 namespace rtk
46 {
206 template <typename VolumeSeriesType, typename ProjectionStackType>
208  : public rtk::IterativeConeBeamReconstructionFilter<VolumeSeriesType, ProjectionStackType>
209 {
210 public:
211  ITK_DISALLOW_COPY_AND_MOVE(FourDROOSTERConeBeamReconstructionFilter);
212 
217  using VolumeType = ProjectionStackType;
219  itk::CovariantVector<typename VolumeSeriesType::ValueType, VolumeSeriesType::ImageDimension - 1>;
222 
224  using CPUVolumeSeriesType =
226 #ifdef RTK_USE_CUDA
227  typedef
228  typename std::conditional<std::is_same<VolumeSeriesType, CPUVolumeSeriesType>::value,
230  itk::CudaImage<CovariantVectorForSpatialGradient, VolumeSeriesType::ImageDimension>>::type
232  typedef typename std::conditional<std::is_same<VolumeSeriesType, CPUVolumeSeriesType>::value,
234  itk::CudaImage<CovariantVectorForTemporalGradient,
235  VolumeSeriesType::ImageDimension>>::type TemporalGradientImageType;
236  typedef typename std::conditional<std::is_same<VolumeSeriesType, CPUVolumeSeriesType>::value,
238  itk::CudaImage<DVFVectorType, VolumeSeriesType::ImageDimension>>::type
240  typedef
241  typename std::conditional<std::is_same<VolumeSeriesType, CPUVolumeSeriesType>::value,
242  itk::Image<DVFVectorType, VolumeSeriesType::ImageDimension - 1>,
243  itk::CudaImage<DVFVectorType, VolumeSeriesType::ImageDimension - 1>>::type DVFImageType;
244  typedef typename std::conditional<std::is_same<VolumeSeriesType, CPUVolumeSeriesType>::value,
246  CudaAverageOutOfROIImageFilter>::type AverageOutOfROIFilterType;
247  typedef typename std::conditional<std::is_same<VolumeSeriesType, CPUVolumeSeriesType>::value,
249  CudaLastDimensionTVDenoisingImageFilter>::type TemporalTVDenoisingFilterType;
250 #else
254  using DVFImageType = itk::Image<DVFVectorType, VolumeSeriesType::ImageDimension - 1>;
258 #endif
259 
261  itkNewMacro(Self);
262 
264 #ifdef itkOverrideGetNameOfClassMacro
265  itkOverrideGetNameOfClassMacro(FourDROOSTERConeBeamReconstructionFilter);
266 #else
268 #endif
269 
270 
272  void
273  SetInputVolumeSeries(const VolumeSeriesType * VolumeSeries);
274  typename VolumeSeriesType::ConstPointer
275  GetInputVolumeSeries();
277 
279  void
280  SetInputProjectionStack(const ProjectionStackType * Projection);
281  typename ProjectionStackType::Pointer
282  GetInputProjectionStack();
284 
286  void
287  SetMotionMask(const VolumeType * mask);
288  typename VolumeType::Pointer
289  GetMotionMask();
291 
293  void
294  SetDisplacementField(const DVFSequenceImageType * DVFs);
295  void
296  SetInverseDisplacementField(const DVFSequenceImageType * DVFs);
298  GetDisplacementField();
300  GetInverseDisplacementField();
302 
303  using FourDCGFilterType =
309  using WarpSequenceFilterType =
314  using AddFilterType = itk::AddImageFilter<VolumeSeriesType, VolumeSeriesType>;
316  using TNVDenoisingFilterType =
318 
319  using ForwardProjectionType = typename Superclass::ForwardProjectionType;
320  using BackProjectionType = typename Superclass::BackProjectionType;
321 
323  virtual void
324  SetWeights(const itk::Array2D<float> _arg);
325 
327  itkSetMacro(DisableDisplacedDetectorFilter, bool);
328  itkGetMacro(DisableDisplacedDetectorFilter, bool);
330 
331  // Regularization steps to perform
332  itkSetMacro(PerformPositivity, bool);
333  itkGetMacro(PerformPositivity, bool);
334  itkSetMacro(PerformMotionMask, bool);
335  itkGetMacro(PerformMotionMask, bool);
336  itkSetMacro(PerformTVSpatialDenoising, bool);
337  itkGetMacro(PerformTVSpatialDenoising, bool);
338  itkSetMacro(PerformWaveletsSpatialDenoising, bool);
339  itkGetMacro(PerformWaveletsSpatialDenoising, bool);
340  itkSetMacro(PerformWarping, bool);
341  itkGetMacro(PerformWarping, bool);
342  itkSetMacro(PerformTVTemporalDenoising, bool);
343  itkGetMacro(PerformTVTemporalDenoising, bool);
344  itkSetMacro(PerformL0TemporalDenoising, bool);
345  itkGetMacro(PerformL0TemporalDenoising, bool);
346  itkSetMacro(PerformTNVDenoising, bool);
347  itkGetMacro(PerformTNVDenoising, bool);
348  itkSetMacro(ComputeInverseWarpingByConjugateGradient, bool);
349  itkGetMacro(ComputeInverseWarpingByConjugateGradient, bool);
350  itkSetMacro(UseNearestNeighborInterpolationInWarping, bool);
351  itkGetMacro(UseNearestNeighborInterpolationInWarping, bool);
352  itkGetMacro(CudaConjugateGradient, bool);
353  itkSetMacro(CudaConjugateGradient, bool);
354 
356  itkSetMacro(UseCudaCyclicDeformation, bool);
357  itkGetMacro(UseCudaCyclicDeformation, bool);
359 
360  // Regularization parameters
361  itkSetMacro(GammaTVSpace, float);
362  itkGetMacro(GammaTVSpace, float);
363  itkSetMacro(GammaTVTime, float);
364  itkGetMacro(GammaTVTime, float);
365  itkSetMacro(GammaTNV, float);
366  itkGetMacro(GammaTNV, float);
367  itkSetMacro(LambdaL0Time, float);
368  itkGetMacro(LambdaL0Time, float);
369  itkSetMacro(SoftThresholdWavelets, float);
370  itkGetMacro(SoftThresholdWavelets, float);
371  itkSetMacro(PhaseShift, float);
372  itkGetMacro(PhaseShift, float);
373 
375  itkGetMacro(NumberOfLevels, unsigned int);
376  itkSetMacro(NumberOfLevels, unsigned int);
378 
380  itkGetMacro(Order, unsigned int);
381  itkSetMacro(Order, unsigned int);
383 
384  // Iterations
385  itkSetMacro(MainLoop_iterations, int);
386  itkGetMacro(MainLoop_iterations, int);
387  itkSetMacro(CG_iterations, int);
388  itkGetMacro(CG_iterations, int);
389  itkSetMacro(TV_iterations, int);
390  itkGetMacro(TV_iterations, int);
391  itkSetMacro(L0_iterations, int);
392  itkGetMacro(L0_iterations, int);
393 
394  // Geometry
395  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry);
396  itkGetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry);
397 
399  virtual void
400  SetSignal(const std::vector<double> signal);
401 
402 protected:
404  ~FourDROOSTERConeBeamReconstructionFilter() override = default;
405 
407  void
408  VerifyPreconditions() const override;
409 
411  void
412  GenerateData() override;
413 
414  void
415  GenerateOutputInformation() override;
416 
417  void
418  GenerateInputRequestedRegion() override;
419 
420  // Inputs are not supposed to occupy the same physical space,
421  // so there is nothing to verify
422  void
423  VerifyInputInformation() const override
424  {}
425 
438  typename AddFilterType::Pointer m_AddFilter;
441 
442  // Booleans :
443  // should warping be performed ?
444  // should conjugate gradient be performed on GPU ?
445  // should wavelets replace TV in spatial denoising ?
455  bool m_UseNearestNeighborInterpolationInWarping; // Default is false, linear interpolation is used instead
459 
460  // Regularization parameters
463  float m_GammaTNV;
467  bool m_DimensionsProcessedForTVSpace[VolumeSeriesType::ImageDimension];
468  bool m_DimensionsProcessedForTVTime[VolumeSeriesType::ImageDimension];
469 
471 
473  unsigned int m_Order;
474  unsigned int m_NumberOfLevels;
475 
476  // Iterations
481 
482  // Geometry
484 
485  // Signal
486  std::vector<double> m_Signal;
487 };
488 } // namespace rtk
489 
490 
491 #ifndef ITK_MANUAL_INSTANTIATION
492 # include "rtkFourDROOSTERConeBeamReconstructionFilter.hxx"
493 #endif
494 
495 #endif
rtk::ThreeDCircularProjectionGeometry::ConstPointer m_Geometry
Denoises along the last dimension, reducing the L0 norm of the gradient.
typename itk::Image< typename VolumeSeriesType::PixelType, VolumeSeriesType::ImageDimension > CPUVolumeSeriesType
Applies an N-D + time Motion Vector Field to an N-D + time sequence of images.
Finds the image sequence that, once warped, equals the input image sequence.
Applies 3D total variation denoising to a 3D + time sequence of images.
Applies a total variation denoising, only alm_SingularValueThresholdFilterong the dimensions specifie...
Implements 4D RecOnstructiOn using Spatial and TEmporal Regularization (short 4D ROOSTER) ...
Projection geometry for a source and a 2-D flat panel.
itk::AddImageFilter< VolumeSeriesType, VolumeSeriesType > AddFilterType
#define itkSetMacro(name, type)
Applies 3D Daubechies wavelets denoising to a 3D + time sequence of images.
SpatialWaveletsDenoisingFilterType::Pointer m_WaveletsDenoisingSpace
Mother class for cone beam reconstruction filters which need runtime selection of their forward and b...
itk::ImageToImageFilter< VolumeSeriesType, VolumeSeriesType >::Pointer m_DownstreamFilter
Implements part of the 4D reconstruction by conjugate gradient.
Averages along the last dimension if the pixel is outside ROI.