RTK  2.0.1
Reconstruction Toolkit
rtkProjectionStackToFourDImageFilter.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 rtkProjectionStackToFourDImageFilter_h
19 #define rtkProjectionStackToFourDImageFilter_h
20 
21 #include <itkExtractImageFilter.h>
22 #include <itkArray2D.h>
23 
26 #include "rtkConstantImageSource.h"
28 
29 #ifdef RTK_USE_CUDA
30  #include "rtkCudaSplatImageFilter.h"
33 #endif
34 
35 namespace rtk
36 {
104 template< typename VolumeSeriesType, typename ProjectionStackType, typename TFFTPrecision=double>
105 class ProjectionStackToFourDImageFilter : public itk::ImageToImageFilter< VolumeSeriesType, VolumeSeriesType >
106 {
107 public:
108  ITK_DISALLOW_COPY_AND_ASSIGN(ProjectionStackToFourDImageFilter);
109 
114 
116  using VolumeType = ProjectionStackType;
117 
119  itkNewMacro(Self)
120 
121 
123 
125  void SetInputVolumeSeries(const VolumeSeriesType* VolumeSeries);
126  typename VolumeSeriesType::ConstPointer GetInputVolumeSeries();
128 
130  void SetInputProjectionStack(const ProjectionStackType* Projections);
131  typename ProjectionStackType::ConstPointer GetInputProjectionStack();
133 
135  using ExtractFilterType = itk::ExtractImageFilter< ProjectionStackType, ProjectionStackType >;
138  using SplatFilterType = rtk::SplatWithKnownWeightsImageFilter<VolumeSeriesType, VolumeType>;
139 
141 
143  using CPUVolumeSeriesType = typename itk::Image< typename VolumeSeriesType::PixelType,
144  VolumeSeriesType::ImageDimension>;
145 #ifdef RTK_USE_CUDA
146  typedef typename std::conditional< std::is_same< VolumeSeriesType, CPUVolumeSeriesType >::value,
148  CudaSplatImageFilter >::type CudaSplatImageFilterType;
149  typedef typename std::conditional< std::is_same< VolumeSeriesType, CPUVolumeSeriesType >::value,
151  CudaConstantVolumeSource >::type CudaConstantVolumeSourceType;
152  typedef typename std::conditional< std::is_same< VolumeSeriesType, CPUVolumeSeriesType >::value,
154  CudaConstantVolumeSeriesSource >::type CudaConstantVolumeSeriesSourceType;
155 #else
156  using CudaSplatImageFilterType = SplatFilterType;
157  using CudaConstantVolumeSourceType = ConstantVolumeSourceType;
158  using CudaConstantVolumeSeriesSourceType = ConstantVolumeSeriesSourceType;
159 #endif
160 
163 
165  itkSetConstObjectMacro(Geometry, GeometryType)
166 
167 
168  itkSetMacro(UseCudaSplat, bool)
169  itkGetMacro(UseCudaSplat, bool)
170  itkSetMacro(UseCudaSources, bool)
171  itkGetMacro(UseCudaSources, bool)
172 
174  itkGetMacro(Weights, itk::Array2D<float>)
175  itkSetMacro(Weights, itk::Array2D<float>)
176 
178  virtual void SetSignal(const std::vector<double> signal);
179 
180 protected:
181  ProjectionStackToFourDImageFilter();
182  ~ProjectionStackToFourDImageFilter() override = default;
183 
185  void GenerateData() override;
186 
187  void GenerateOutputInformation() override;
188 
189  void GenerateInputRequestedRegion() override;
190 
192 
194  typename SplatFilterType::Pointer m_SplatFilter;
197  typename ConstantVolumeSourceType::Pointer m_ConstantVolumeSource;
198  typename ConstantVolumeSeriesSourceType::Pointer m_ConstantVolumeSeriesSource;
199 
202  GeometryType::ConstPointer m_Geometry;
205  std::vector<double> m_Signal;
206 
207 };
208 } //namespace ITK
209 
210 
211 #ifndef ITK_MANUAL_INSTANTIATION
212 #include "rtkProjectionStackToFourDImageFilter.hxx"
213 #endif
214 
215 #endif
Generate an n-dimensional image with constant pixel values.
rtk::ConstantImageSource< VolumeSeriesType > ConstantVolumeSeriesSourceType
typename itk::Image< typename VolumeSeriesType::PixelType, VolumeSeriesType::ImageDimension > CPUVolumeSeriesType
STL namespace.
Splats (linearly) a 3D volume into a 3D+t sequence of volumes.
rtk::SplatWithKnownWeightsImageFilter< VolumeSeriesType, VolumeType > SplatFilterType
Projection geometry for a source and a 2-D flat panel.
void SetInputProjectionStack(const ProjectionStackType *Projections)
VolumeSeriesType::ConstPointer GetInputVolumeSeries()
Implements part of the 4D reconstruction by conjugate gradient.
ProjectionStackType::ConstPointer GetInputProjectionStack()
virtual void SetSignal(const std::vector< double > signal)
void SetInputVolumeSeries(const VolumeSeriesType *VolumeSeries)
ConstantVolumeSeriesSourceType::Pointer m_ConstantVolumeSeriesSource
void SetBackProjectionFilter(const typename BackProjectionFilterType::Pointer _arg)
#define itkSetMacro(name, type)
rtk::ConstantImageSource< VolumeType > ConstantVolumeSourceType