RTK  2.0.1
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  * http://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 {
195 template< typename VolumeSeriesType, typename ProjectionStackType>
197 {
198 public:
199  ITK_DISALLOW_COPY_AND_ASSIGN(FourDROOSTERConeBeamReconstructionFilter);
200 
205  using VolumeType = ProjectionStackType;
206  using CovariantVectorForSpatialGradient = itk::CovariantVector< typename VolumeSeriesType::ValueType, VolumeSeriesType::ImageDimension - 1>;
209 
211  using CPUVolumeSeriesType = typename itk::Image< typename VolumeSeriesType::PixelType,
212  VolumeSeriesType::ImageDimension>;
213 #ifdef RTK_USE_CUDA
214  typedef typename std::conditional< std::is_same< VolumeSeriesType, CPUVolumeSeriesType >::value,
216  itk::CudaImage<CovariantVectorForSpatialGradient, VolumeSeriesType::ImageDimension> >::type SpatialGradientImageType;
217  typedef typename std::conditional< std::is_same< VolumeSeriesType, CPUVolumeSeriesType >::value,
219  itk::CudaImage<CovariantVectorForTemporalGradient, VolumeSeriesType::ImageDimension> >::type TemporalGradientImageType;
220  typedef typename std::conditional< std::is_same< VolumeSeriesType, CPUVolumeSeriesType >::value,
222  itk::CudaImage<DVFVectorType, VolumeSeriesType::ImageDimension> >::type DVFSequenceImageType;
223  typedef typename std::conditional< std::is_same< VolumeSeriesType, CPUVolumeSeriesType >::value,
224  itk::Image<DVFVectorType, VolumeSeriesType::ImageDimension - 1>,
225  itk::CudaImage<DVFVectorType, VolumeSeriesType::ImageDimension - 1> >::type DVFImageType;
226  typedef typename std::conditional< std::is_same< VolumeSeriesType, CPUVolumeSeriesType >::value,
228  CudaAverageOutOfROIImageFilter >::type AverageOutOfROIFilterType;
229  typedef typename std::conditional< std::is_same< VolumeSeriesType, CPUVolumeSeriesType >::value,
231  CudaLastDimensionTVDenoisingImageFilter >::type TemporalTVDenoisingFilterType;
232 #else
233  using SpatialGradientImageType = itk::Image<CovariantVectorForSpatialGradient, VolumeSeriesType::ImageDimension>;
234  using TemporalGradientImageType = itk::Image<CovariantVectorForTemporalGradient, VolumeSeriesType::ImageDimension>;
235  using DVFSequenceImageType = itk::Image<DVFVectorType, VolumeSeriesType::ImageDimension>;
236  using DVFImageType = itk::Image<DVFVectorType, VolumeSeriesType::ImageDimension - 1>;
237  using AverageOutOfROIFilterType = AverageOutOfROIImageFilter <VolumeSeriesType, VolumeType>;
238  using TemporalTVDenoisingFilterType = TotalVariationDenoisingBPDQImageFilter<VolumeSeriesType, TemporalGradientImageType>;
239 #endif
240 
242  itkNewMacro(Self)
243 
244 
246 
248  void SetInputVolumeSeries(const VolumeSeriesType* VolumeSeries);
249  typename VolumeSeriesType::ConstPointer GetInputVolumeSeries();
251 
253  void SetInputProjectionStack(const ProjectionStackType* Projection);
254  typename ProjectionStackType::Pointer GetInputProjectionStack();
256 
258  void SetMotionMask(const VolumeType* mask);
259  typename VolumeType::Pointer GetMotionMask();
261 
263  void SetDisplacementField(const DVFSequenceImageType* DVFs);
264  void SetInverseDisplacementField(const DVFSequenceImageType* DVFs);
265  typename DVFSequenceImageType::Pointer GetDisplacementField();
266  typename DVFSequenceImageType::Pointer GetInverseDisplacementField();
268 
269  using FourDCGFilterType = rtk::FourDConjugateGradientConeBeamReconstructionFilter<VolumeSeriesType, ProjectionStackType>;
270  using ThresholdFilterType = itk::ThresholdImageFilter<VolumeSeriesType>;
274  using WarpSequenceFilterType = rtk::WarpSequenceImageFilter<VolumeSeriesType, DVFSequenceImageType, VolumeType, DVFImageType>;
275  using UnwarpSequenceFilterType = rtk::UnwarpSequenceImageFilter<VolumeSeriesType, DVFSequenceImageType, VolumeType, DVFImageType>;
276  using SubtractFilterType = itk::SubtractImageFilter<VolumeSeriesType, VolumeSeriesType>;
277  using AddFilterType = itk::AddImageFilter<VolumeSeriesType, VolumeSeriesType>;
279  using TNVDenoisingFilterType = rtk::TotalNuclearVariationDenoisingBPDQImageFilter<VolumeSeriesType, SpatialGradientImageType>;
280 
283 
285  void SetForwardProjectionFilter(ForwardProjectionType fwtype) override;
286 
288  void SetBackProjectionFilter(BackProjectionType bptype) override;
289 
291  virtual void SetWeights(const itk::Array2D<float> _arg);
292 
294  itkSetMacro(DisableDisplacedDetectorFilter, bool)
295  itkGetMacro(DisableDisplacedDetectorFilter, bool)
296 
297  // Regularization steps to perform
298  itkSetMacro(PerformPositivity, bool)
299  itkGetMacro(PerformPositivity, bool)
300  itkSetMacro(PerformMotionMask, bool)
301  itkGetMacro(PerformMotionMask, bool)
302  itkSetMacro(PerformTVSpatialDenoising, bool)
303  itkGetMacro(PerformTVSpatialDenoising, bool)
304  itkSetMacro(PerformWaveletsSpatialDenoising, bool)
305  itkGetMacro(PerformWaveletsSpatialDenoising, bool)
306  itkSetMacro(PerformWarping, bool)
307  itkGetMacro(PerformWarping, bool)
308  itkSetMacro(PerformTVTemporalDenoising, bool)
309  itkGetMacro(PerformTVTemporalDenoising, bool)
310  itkSetMacro(PerformL0TemporalDenoising, bool)
311  itkGetMacro(PerformL0TemporalDenoising, bool)
312  itkSetMacro(PerformTNVDenoising, bool)
313  itkGetMacro(PerformTNVDenoising, bool)
314  itkSetMacro(ComputeInverseWarpingByConjugateGradient, bool)
315  itkGetMacro(ComputeInverseWarpingByConjugateGradient, bool)
316  itkSetMacro(UseNearestNeighborInterpolationInWarping, bool)
317  itkGetMacro(UseNearestNeighborInterpolationInWarping, bool)
318  itkGetMacro(CudaConjugateGradient, bool)
319  itkSetMacro(CudaConjugateGradient, bool)
320 
322  itkSetMacro(UseCudaCyclicDeformation, bool)
323  itkGetMacro(UseCudaCyclicDeformation, bool)
324 
325  // Regularization parameters
326  itkSetMacro(GammaTVSpace, float)
327  itkGetMacro(GammaTVSpace, float)
328  itkSetMacro(GammaTVTime, float)
329  itkGetMacro(GammaTVTime, float)
330  itkSetMacro(GammaTNV, float)
331  itkGetMacro(GammaTNV, float)
332  itkSetMacro(LambdaL0Time, float)
333  itkGetMacro(LambdaL0Time, float)
334  itkSetMacro(SoftThresholdWavelets, float)
335  itkGetMacro(SoftThresholdWavelets, float)
336  itkSetMacro(PhaseShift, float)
337  itkGetMacro(PhaseShift, float)
338 
340  itkGetMacro(NumberOfLevels, unsigned int)
341  itkSetMacro(NumberOfLevels, unsigned int)
342 
344  itkGetMacro(Order, unsigned int)
345  itkSetMacro(Order, unsigned int)
346 
347  // Iterations
348  itkSetMacro(MainLoop_iterations, int)
349  itkGetMacro(MainLoop_iterations, int)
350  itkSetMacro(CG_iterations, int)
351  itkGetMacro(CG_iterations, int)
352  itkSetMacro(TV_iterations, int)
353  itkGetMacro(TV_iterations, int)
354  itkSetMacro(L0_iterations, int)
355  itkGetMacro(L0_iterations, int)
356 
357  // Geometry
358  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
359  itkGetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
360 
362  virtual void SetSignal(const std::vector<double> signal);
363 
364 protected:
365  FourDROOSTERConeBeamReconstructionFilter();
366  ~FourDROOSTERConeBeamReconstructionFilter() override = default;
367 
369  void GenerateData() override;
370 
371  void GenerateOutputInformation() override;
372 
373  void GenerateInputRequestedRegion() override;
374 
375  // Inputs are not supposed to occupy the same physical space,
376  // so there is nothing to verify
377 #if ITK_VERSION_MAJOR<5
378  void VerifyInputInformation() override {}
379 #else
380  void VerifyInputInformation() const override {}
381 #endif
382 
398 
399  // Booleans :
400  // should warping be performed ?
401  // should conjugate gradient be performed on GPU ?
402  // should wavelets replace TV in spatial denoising ?
412  bool m_UseNearestNeighborInterpolationInWarping; //Default is false, linear interpolation is used instead
416 
417  // Regularization parameters
420  float m_GammaTNV;
424  bool m_DimensionsProcessedForTVSpace[VolumeSeriesType::ImageDimension];
425  bool m_DimensionsProcessedForTVTime[VolumeSeriesType::ImageDimension];
426 
428 
430  unsigned int m_Order;
431  unsigned int m_NumberOfLevels;
432 
433  // Iterations
438 
439  // Geometry
441 
442  // Signal
443  std::vector<double> m_Signal;
444 
445 };
446 } //namespace ITK
447 
448 
449 #ifndef ITK_MANUAL_INSTANTIATION
450 #include "rtkFourDROOSTERConeBeamReconstructionFilter.hxx"
451 #endif
452 
453 #endif
rtk::ThreeDCircularProjectionGeometry::ConstPointer m_Geometry
Denoises along the last dimension, reducing the L0 norm of the gradient.
TotalVariationDenoisingBPDQImageFilter< VolumeSeriesType, TemporalGradientImageType > TemporalTVDenoisingFilterType
VolumeSeriesType::ConstPointer GetInputVolumeSeries()
void SetInputVolumeSeries(const VolumeSeriesType *VolumeSeries)
itk::Image< DVFVectorType, VolumeSeriesType::ImageDimension > DVFSequenceImageType
Applies an N-D + time Motion Vector Field to an N-D + time sequence of images.
STL namespace.
Finds the image sequence that, once warped, equals the input image sequence.
AverageOutOfROIImageFilter< VolumeSeriesType, VolumeType > AverageOutOfROIFilterType
Applies 3D total variation denoising to a 3D + time sequence of images.
Applies a total variation denoising, only alm_SingularValueThresholdFilterong the dimensions specifie...
ProjectionStackType::Pointer GetInputProjectionStack()
itk::Image< CovariantVectorForSpatialGradient, VolumeSeriesType::ImageDimension > SpatialGradientImageType
Implements 4D RecOnstructiOn using Spatial and TEmporal Regularization (short 4D ROOSTER) ...
Projection geometry for a source and a 2-D flat panel.
Applies 3D Daubechies wavelets denoising to a 3D + time sequence of images.
SpatialWaveletsDenoisingFilterType::Pointer m_WaveletsDenoisingSpace
itk::CovariantVector< typename VolumeSeriesType::ValueType, VolumeSeriesType::ImageDimension-1 > CovariantVectorForSpatialGradient
itk::Image< CovariantVectorForTemporalGradient, VolumeSeriesType::ImageDimension > TemporalGradientImageType
void SetDisplacementField(const DVFSequenceImageType *DVFs)
typename itk::Image< typename VolumeSeriesType::PixelType, VolumeSeriesType::ImageDimension > CPUVolumeSeriesType
void SetInverseDisplacementField(const DVFSequenceImageType *DVFs)
void SetForwardProjectionFilter(ForwardProjectionType fwtype) override
virtual void SetSignal(const std::vector< double > signal)
Mother class for cone beam reconstruction filters which need runtime selection of their forward and b...
void SetBackProjectionFilter(BackProjectionType bptype) override
itk::ImageToImageFilter< VolumeSeriesType, VolumeSeriesType >::Pointer m_DownstreamFilter
Implements part of the 4D reconstruction by conjugate gradient.
void SetInputProjectionStack(const ProjectionStackType *Projection)
itk::Image< DVFVectorType, VolumeSeriesType::ImageDimension-1 > DVFImageType
DVFSequenceImageType::Pointer GetDisplacementField()
#define itkSetMacro(name, type)
Averages along the last dimension if the pixel is outside ROI.
virtual void SetWeights(const itk::Array2D< float > _arg)
DVFSequenceImageType::Pointer GetInverseDisplacementField()
void SetMotionMask(const VolumeType *mask)