RTK  2.0.1
Reconstruction Toolkit
rtkProjectionsReader.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 
19 #ifndef rtkProjectionsReader_h
20 #define rtkProjectionsReader_h
21 
22 // ITK
23 #include <itkImageSource.h>
24 #include <itkImageIOFactory.h>
26 
27 // RTK
30 
31 // Standard lib
32 #include <vector>
33 #include <string>
34 
35 namespace rtk
36 {
37 
119 template <class TOutputImage>
120 class ITK_EXPORT ProjectionsReader : public itk::ImageSource<TOutputImage>
121 {
122 public:
123  ITK_DISALLOW_COPY_AND_ASSIGN(ProjectionsReader);
124 
129 
131  itkNewMacro(Self);
132 
134  itkTypeMacro(ProjectionsReader, itk::ImageSource);
135 
137  using OutputImageType = TOutputImage;
138  using OutputImagePointer = typename OutputImageType::Pointer;
139  using OutputImageRegionType = typename OutputImageType::RegionType;
140  using OutputImagePixelType = typename OutputImageType::PixelType;
141  using OutputImageDirectionType = typename OutputImageType::DirectionType;
142  using OutputImageSpacingType = typename OutputImageType::SpacingType;
143  using OutputImagePointType = typename OutputImageType::PointType;
144  using OutputImageSizeType = typename OutputImageType::SizeType;
145 
146  using FileNamesContainer = std::vector<std::string>;
149  using WaterPrecorrectionVectorType = std::vector< double >;
150 
155 
157  itkStaticConstMacro(OutputImageDimension, unsigned int,
158  TOutputImage::ImageDimension);
159 
162  void SetFileNames (const FileNamesContainer &name)
163  {
164  if ( m_FileNames != name)
165  {
166  m_FileNames = name;
167  this->Modified();
168  }
169  }
171  {
172  return m_FileNames;
173  }
175 
178  itkGetConstMacro(Origin, OutputImagePointType);
180 
182  itkGetConstMacro(Spacing, OutputImageSpacingType);
183 
185  itkGetConstMacro(Direction, OutputImageDirectionType);
186 
188  itkSetMacro(UpperBoundaryCropSize, OutputImageSizeType);
189  itkGetConstMacro(UpperBoundaryCropSize, OutputImageSizeType);
190  itkSetMacro(LowerBoundaryCropSize, OutputImageSizeType);
191  itkGetConstMacro(LowerBoundaryCropSize, OutputImageSizeType);
193 
195  itkSetMacro(ShrinkFactors, ShrinkFactorsType);
196  itkGetConstReferenceMacro(ShrinkFactors, ShrinkFactorsType);
198 
200  itkSetMacro(MedianRadius, MedianRadiusType);
201  itkGetConstReferenceMacro(MedianRadius, MedianRadiusType);
202  itkGetMacro(ConditionalMedianThresholdMultiplier, double);
203  itkSetMacro(ConditionalMedianThresholdMultiplier, double);
205 
207  itkGetMacro(AirThreshold, double);
208  itkSetMacro(AirThreshold, double);
210 
211  itkGetMacro(ScatterToPrimaryRatio, double);
212  itkSetMacro(ScatterToPrimaryRatio, double);
213 
214  itkGetMacro(NonNegativityConstraintThreshold, double);
215  itkSetMacro(NonNegativityConstraintThreshold, double);
216 
221  itkGetMacro(I0, double);
222  itkSetMacro(I0, double);
224 
228  itkGetMacro(IDark, double);
229  itkSetMacro(IDark, double);
231 
233  itkGetMacro(WaterPrecorrectionCoefficients, WaterPrecorrectionVectorType);
235  {
236  if (this->m_WaterPrecorrectionCoefficients != _arg)
237  {
238  this->m_WaterPrecorrectionCoefficients = _arg;
239  this->Modified();
240  }
241  }
243 
246  itkSetMacro(ComputeLineIntegral, bool);
247  itkGetConstMacro(ComputeLineIntegral, bool);
248  itkBooleanMacro(ComputeLineIntegral);
250 
253  itkGetMacro(VectorComponent, unsigned int)
254  itkSetMacro(VectorComponent, unsigned int)
255 
256 
257  itkGetMacro(ImageIO, itk::ImageIOBase::Pointer);
258 
261  void GenerateOutputInformation(void) override;
262 
263 protected:
265  ~ProjectionsReader() override = default;
266  void PrintSelf(std::ostream& os, itk::Indent indent) const override;
267 
269  void GenerateData() override;
270 
273 
274 private:
278  template<class TInputImage> void PropagateParametersToMiniPipeline();
279  void ConnectElektaRawFilter(itk::ImageBase<OutputImageDimension> **nextInputBase);
280  void PropagateI0(itk::ImageBase<OutputImageDimension> **nextInputBase);
282 
286  itk::ProcessObject::Pointer m_RawDataReader;
287 
289  itk::ProcessObject::Pointer m_VectorComponentSelectionFilter;
290  itk::ProcessObject::Pointer m_ChangeInformationFilter;
291  itk::ProcessObject::Pointer m_ElektaRawFilter;
292  itk::ProcessObject::Pointer m_CropFilter;
293  itk::ProcessObject::Pointer m_ConditionalMedianFilter;
294  itk::ProcessObject::Pointer m_BinningFilter;
295  itk::ProcessObject::Pointer m_ScatterFilter;
296  itk::ProcessObject::Pointer m_I0EstimationFilter;
297 
301 
305 
309 
311  itk::ImageIOBase::Pointer m_ImageIO{nullptr};
312 
323  double m_AirThreshold{32000};
324  double m_ScatterToPrimaryRatio{0.};
325  double m_NonNegativityConstraintThreshold{itk::NumericTraits<double>::NonpositiveMin()};
326  double m_I0{itk::NumericTraits<double>::NonpositiveMin()};
327  double m_IDark{0.};
328  double m_ConditionalMedianThresholdMultiplier{1.};
330  bool m_ComputeLineIntegral{true};
331  unsigned int m_VectorComponent{0};
332 };
334 
335 } //namespace rtk
336 
337 #ifndef ITK_MANUAL_INSTANTIATION
338 #include "rtkProjectionsReader.hxx"
339 #endif
340 
341 #endif // rtkProjectionsReader_h
typename OutputImageType::SpacingType OutputImageSpacingType
typename OutputImageType::RegionType OutputImageRegionType
OutputImageDirectionType m_Direction
Performs the classical water precorrection for beam hardening (Kachelriess, Med. Phys. 2006)
OutputImagePointType m_Origin
typename OutputImageType::Pointer OutputImagePointer
itk::ProcessObject::Pointer m_VectorComponentSelectionFilter
itk::ProcessObject::Pointer m_ScatterFilter
OutputImageSpacingType m_Spacing
OutputImageSizeType m_UpperBoundaryCropSize
WaterPrecorrectionVectorType m_WaterPrecorrectionCoefficients
typename OutputImageType::SizeType OutputImageSizeType
typename OutputImageType::PointType OutputImagePointType
itk::ImageSource< TOutputImage >::Pointer m_RawCastFilter
const FileNamesContainer & GetFileNames() const
FileNamesContainer m_FileNames
itk::ProcessObject::Pointer m_I0EstimationFilter
itk::ProcessObject::Pointer m_BinningFilter
StreamingType::Pointer m_StreamingFilter
MedianRadiusType m_MedianRadius
typename rtk::ConditionalMedianImageFilter< TOutputImage >::MedianRadiusType MedianRadiusType
itk::ProcessObject::Pointer m_ChangeInformationFilter
itk::ImageSource< TOutputImage >::Pointer m_RawToAttenuationFilter
typename OutputImageType::PixelType OutputImagePixelType
OutputImageSizeType m_LowerBoundaryCropSize
typename OutputImageType::DirectionType OutputImageDirectionType
itk::ProcessObject::Pointer m_CropFilter
std::vector< std::string > FileNamesContainer
void SetFileNames(const FileNamesContainer &name)
itk::ProcessObject::Pointer m_ElektaRawFilter
typename itk::ConstNeighborhoodIterator< TInputImage >::RadiusType MedianRadiusType
WaterPrecorrectionType::Pointer m_WaterPrecorrectionFilter
itk::ProcessObject::Pointer m_RawDataReader
std::vector< double > WaterPrecorrectionVectorType
itk::ProcessObject::Pointer m_ConditionalMedianFilter
virtual void SetWaterPrecorrectionCoefficients(const WaterPrecorrectionVectorType _arg)
ShrinkFactorsType m_ShrinkFactors
#define itkSetMacro(name, type)