RTK  1.4.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  * 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:
203  typedef ProjectionStackType VolumeType;
204  typedef itk::CovariantVector< typename VolumeSeriesType::ValueType, VolumeSeriesType::ImageDimension - 1> CovariantVectorForSpatialGradient;
206  typedef CovariantVectorForSpatialGradient DVFVectorType;
207 
208 #ifdef RTK_USE_CUDA
209  typedef itk::CudaImage<CovariantVectorForSpatialGradient, VolumeSeriesType::ImageDimension> SpatialGradientImageType;
210  typedef itk::CudaImage<CovariantVectorForTemporalGradient, VolumeSeriesType::ImageDimension> TemporalGradientImageType;
211  typedef itk::CudaImage<DVFVectorType, VolumeSeriesType::ImageDimension> DVFSequenceImageType;
212  typedef itk::CudaImage<DVFVectorType, VolumeSeriesType::ImageDimension - 1> DVFImageType;
213 #else
217  typedef itk::Image<DVFVectorType, VolumeSeriesType::ImageDimension - 1> DVFImageType;
218 #endif
219 
221  itkNewMacro(Self)
222 
223 
225 
227  void SetInputVolumeSeries(const VolumeSeriesType* VolumeSeries);
228  typename VolumeSeriesType::ConstPointer GetInputVolumeSeries();
230 
232  void SetInputProjectionStack(const ProjectionStackType* Projection);
233  typename ProjectionStackType::Pointer GetInputProjectionStack();
235 
237  void SetMotionMask(const VolumeType* mask);
238  typename VolumeType::Pointer GetMotionMask();
240 
242  void SetDisplacementField(const DVFSequenceImageType* DVFs);
243  void SetInverseDisplacementField(const DVFSequenceImageType* DVFs);
244  typename DVFSequenceImageType::Pointer GetDisplacementField();
245  typename DVFSequenceImageType::Pointer GetInverseDisplacementField();
247 
248  typedef rtk::FourDConjugateGradientConeBeamReconstructionFilter<VolumeSeriesType, ProjectionStackType> FourDCGFilterType;
249  typedef itk::ThresholdImageFilter<VolumeSeriesType> ThresholdFilterType;
250  typedef itk::ResampleImageFilter<VolumeType, VolumeType> ResampleFilterType;
251  typedef rtk::AverageOutOfROIImageFilter <VolumeSeriesType, VolumeType> AverageOutOfROIFilterType;
254  typedef rtk::WarpSequenceImageFilter<VolumeSeriesType, DVFSequenceImageType, VolumeType, DVFImageType> WarpSequenceFilterType;
255  typedef rtk::TotalVariationDenoisingBPDQImageFilter<VolumeSeriesType, TemporalGradientImageType> TemporalTVDenoisingFilterType;
256  typedef rtk::UnwarpSequenceImageFilter<VolumeSeriesType, DVFSequenceImageType, VolumeType, DVFImageType> UnwarpSequenceFilterType;
257  typedef itk::SubtractImageFilter<VolumeSeriesType, VolumeSeriesType> SubtractFilterType;
258  typedef itk::AddImageFilter<VolumeSeriesType, VolumeSeriesType> AddFilterType;
260  typedef rtk::TotalNuclearVariationDenoisingBPDQImageFilter<VolumeSeriesType, SpatialGradientImageType> TNVDenoisingFilterType;
261 
262  typedef typename Superclass::ForwardProjectionType ForwardProjectionType;
263  typedef typename Superclass::BackProjectionType BackProjectionType;
264 
266  void SetForwardProjectionFilter(ForwardProjectionType fwtype) ITK_OVERRIDE;
267 
269  void SetBackProjectionFilter(BackProjectionType bptype) ITK_OVERRIDE;
270 
272  virtual void SetWeights(const itk::Array2D<float> _arg);
273 
275  itkSetMacro(DisableDisplacedDetectorFilter, bool)
276  itkGetMacro(DisableDisplacedDetectorFilter, bool)
277 
278  // Regularization steps to perform
279  itkSetMacro(PerformPositivity, bool)
280  itkGetMacro(PerformPositivity, bool)
281  itkSetMacro(PerformMotionMask, bool)
282  itkGetMacro(PerformMotionMask, bool)
283  itkSetMacro(PerformTVSpatialDenoising, bool)
284  itkGetMacro(PerformTVSpatialDenoising, bool)
285  itkSetMacro(PerformWaveletsSpatialDenoising, bool)
286  itkGetMacro(PerformWaveletsSpatialDenoising, bool)
287  itkSetMacro(PerformWarping, bool)
288  itkGetMacro(PerformWarping, bool)
289  itkSetMacro(PerformTVTemporalDenoising, bool)
290  itkGetMacro(PerformTVTemporalDenoising, bool)
291  itkSetMacro(PerformL0TemporalDenoising, bool)
292  itkGetMacro(PerformL0TemporalDenoising, bool)
293  itkSetMacro(PerformTNVDenoising, bool)
294  itkGetMacro(PerformTNVDenoising, bool)
295  itkSetMacro(ComputeInverseWarpingByConjugateGradient, bool)
296  itkGetMacro(ComputeInverseWarpingByConjugateGradient, bool)
297  itkSetMacro(UseNearestNeighborInterpolationInWarping, bool)
298  itkGetMacro(UseNearestNeighborInterpolationInWarping, bool)
299  itkGetMacro(CudaConjugateGradient, bool)
300  itkSetMacro(CudaConjugateGradient, bool)
301 
303  itkSetMacro(UseCudaCyclicDeformation, bool)
304  itkGetMacro(UseCudaCyclicDeformation, bool)
305 
306  // Regularization parameters
307  itkSetMacro(GammaTVSpace, float)
308  itkGetMacro(GammaTVSpace, float)
309  itkSetMacro(GammaTVTime, float)
310  itkGetMacro(GammaTVTime, float)
311  itkSetMacro(GammaTNV, float)
312  itkGetMacro(GammaTNV, float)
313  itkSetMacro(LambdaL0Time, float)
314  itkGetMacro(LambdaL0Time, float)
315  itkSetMacro(SoftThresholdWavelets, float)
316  itkGetMacro(SoftThresholdWavelets, float)
317  itkSetMacro(PhaseShift, float)
318  itkGetMacro(PhaseShift, float)
319 
321  itkGetMacro(NumberOfLevels, unsigned int)
322  itkSetMacro(NumberOfLevels, unsigned int)
323 
325  itkGetMacro(Order, unsigned int)
326  itkSetMacro(Order, unsigned int)
327 
328  // Iterations
329  itkSetMacro(MainLoop_iterations, int)
330  itkGetMacro(MainLoop_iterations, int)
331  itkSetMacro(CG_iterations, int)
332  itkGetMacro(CG_iterations, int)
333  itkSetMacro(TV_iterations, int)
334  itkGetMacro(TV_iterations, int)
335  itkSetMacro(L0_iterations, int)
336  itkGetMacro(L0_iterations, int)
337 
338  // Geometry
339  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
340  itkGetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
341 
343  virtual void SetSignal(const std::vector<double> signal);
344 
345 protected:
346  FourDROOSTERConeBeamReconstructionFilter();
347  virtual ~FourDROOSTERConeBeamReconstructionFilter() ITK_OVERRIDE {}
348 
350  void GenerateData() ITK_OVERRIDE;
351 
352  void GenerateOutputInformation() ITK_OVERRIDE;
353 
354  void GenerateInputRequestedRegion() ITK_OVERRIDE;
355 
356  // Inputs are not supposed to occupy the same physical space,
357  // so there is nothing to verify
358  void VerifyInputInformation() ITK_OVERRIDE {}
359 
375 
376  // Booleans :
377  // should warping be performed ?
378  // should conjugate gradient be performed on GPU ?
379  // should wavelets replace TV in spatial denoising ?
389  bool m_UseNearestNeighborInterpolationInWarping; //Default is false, linear interpolation is used instead
393 
394  // Regularization parameters
397  float m_GammaTNV;
401  bool m_DimensionsProcessedForTVSpace[VolumeSeriesType::ImageDimension];
402  bool m_DimensionsProcessedForTVTime[VolumeSeriesType::ImageDimension];
403 
405 
407  unsigned int m_Order;
408  unsigned int m_NumberOfLevels;
409 
410  // Iterations
415 
416  // Geometry
418 
419  // Signal
420  std::vector<double> m_Signal;
421 
422 private:
423  FourDROOSTERConeBeamReconstructionFilter(const Self &); //purposely not implemented
424  void operator=(const Self &); //purposely not implemented
425 
426 };
427 } //namespace ITK
428 
429 
430 #ifndef ITK_MANUAL_INSTANTIATION
431 #include "rtkFourDROOSTERConeBeamReconstructionFilter.hxx"
432 #endif
433 
434 #endif
rtk::ThreeDCircularProjectionGeometry::ConstPointer m_Geometry
Denoises along the last dimension, reducing the L0 norm of the gradient.
VolumeSeriesType::ConstPointer GetInputVolumeSeries()
rtk::IterativeConeBeamReconstructionFilter< VolumeSeriesType, ProjectionStackType > Superclass
void SetInputVolumeSeries(const VolumeSeriesType *VolumeSeries)
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.
Applies 3D total variation denoising to a 3D + time sequence of images.
Applies a total variation denoising, only alm_SingularValueThresholdFilterong the dimensions specifie...
itk::Image< DVFVectorType, VolumeSeriesType::ImageDimension > DVFSequenceImageType
ProjectionStackType::Pointer GetInputProjectionStack()
Implements 4D RecOnstructiOn using Spatial and TEmporal Regularization (short 4D ROOSTER) ...
Projection geometry for a source and a 2-D flat panel.
itk::CovariantVector< typename VolumeSeriesType::ValueType, VolumeSeriesType::ImageDimension-1 > CovariantVectorForSpatialGradient
itk::Image< DVFVectorType, VolumeSeriesType::ImageDimension-1 > DVFImageType
Applies 3D Daubechies wavelets denoising to a 3D + time sequence of images.
SpatialWaveletsDenoisingFilterType::Pointer m_WaveletsDenoisingSpace
void SetDisplacementField(const DVFSequenceImageType *DVFs)
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.
itk::Image< CovariantVectorForSpatialGradient, VolumeSeriesType::ImageDimension > SpatialGradientImageType
itk::Image< CovariantVectorForTemporalGradient, VolumeSeriesType::ImageDimension > TemporalGradientImageType
void SetInputProjectionStack(const ProjectionStackType *Projection)
DVFSequenceImageType::Pointer GetDisplacementField()
#define itkSetMacro(name, type)
Averages along the last dimension if the pixel is outside ROI.
itk::CovariantVector< typename VolumeSeriesType::ValueType, 1 > CovariantVectorForTemporalGradient
virtual void SetWeights(const itk::Array2D< float > _arg)
DVFSequenceImageType::Pointer GetInverseDisplacementField()
void SetMotionMask(const VolumeType *mask)