RTK  2.0.0
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:
112 
114  typedef ProjectionStackType VolumeType;
115 
117  itkNewMacro(Self)
118 
119 
121 
123  void SetInputVolumeSeries(const VolumeSeriesType* VolumeSeries);
124  typename VolumeSeriesType::ConstPointer GetInputVolumeSeries();
126 
128  void SetInputProjectionStack(const ProjectionStackType* Projections);
129  typename ProjectionStackType::ConstPointer GetInputProjectionStack();
131 
132  typedef rtk::BackProjectionImageFilter< VolumeType, VolumeType > BackProjectionFilterType;
133  typedef itk::ExtractImageFilter< ProjectionStackType, ProjectionStackType > ExtractFilterType;
136  typedef rtk::SplatWithKnownWeightsImageFilter<VolumeSeriesType, VolumeType> SplatFilterType;
137 
139 
141  typedef typename itk::Image< typename VolumeSeriesType::PixelType,
142  VolumeSeriesType::ImageDimension> CPUVolumeSeriesType;
143 #ifdef RTK_USE_CUDA
144  typedef typename std::conditional< std::is_same< VolumeSeriesType, CPUVolumeSeriesType >::value,
146  CudaSplatImageFilter >::type CudaSplatImageFilterType;
147  typedef typename std::conditional< std::is_same< VolumeSeriesType, CPUVolumeSeriesType >::value,
149  CudaConstantVolumeSource >::type CudaConstantVolumeSourceType;
150  typedef typename std::conditional< std::is_same< VolumeSeriesType, CPUVolumeSeriesType >::value,
152  CudaConstantVolumeSeriesSource >::type CudaConstantVolumeSeriesSourceType;
153 #else
154  typedef SplatFilterType CudaSplatImageFilterType;
155  typedef ConstantVolumeSourceType CudaConstantVolumeSourceType;
156  typedef ConstantVolumeSeriesSourceType CudaConstantVolumeSeriesSourceType;
157 #endif
158 
161 
163  itkSetConstObjectMacro(Geometry, GeometryType)
164 
165 
166  itkSetMacro(UseCudaSplat, bool)
167  itkGetMacro(UseCudaSplat, bool)
168  itkSetMacro(UseCudaSources, bool)
169  itkGetMacro(UseCudaSources, bool)
170 
172  itkGetMacro(Weights, itk::Array2D<float>)
173  itkSetMacro(Weights, itk::Array2D<float>)
174 
176  virtual void SetSignal(const std::vector<double> signal);
177 
178 protected:
179  ProjectionStackToFourDImageFilter();
180  virtual ~ProjectionStackToFourDImageFilter() ITK_OVERRIDE {}
181 
183  void GenerateData() ITK_OVERRIDE;
184 
185  void GenerateOutputInformation() ITK_OVERRIDE;
186 
187  void GenerateInputRequestedRegion() ITK_OVERRIDE;
188 
190 
194  typename ExtractFilterType::Pointer m_ExtractFilter;
197 
203  std::vector<double> m_Signal;
204 
205 private:
206  ProjectionStackToFourDImageFilter(const Self &); //purposely not implemented
207  void operator=(const Self &); //purposely not implemented
208 
209 };
210 } //namespace ITK
211 
212 
213 #ifndef ITK_MANUAL_INSTANTIATION
214 #include "rtkProjectionStackToFourDImageFilter.hxx"
215 #endif
216 
217 #endif
rtk::ConstantImageSource< VolumeType > ConstantVolumeSourceType
rtk::ConstantImageSource< VolumeSeriesType > ConstantVolumeSeriesSourceType
Generate an n-dimensional image with constant pixel values.
STL namespace.
Splats (linearly) a 3D volume into a 3D+t sequence of volumes.
itk::ImageToImageFilter< VolumeSeriesType, VolumeSeriesType > Superclass
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)