RTK  2.6.0
Reconstruction Toolkit
rtkJosephForwardProjectionImageFilter.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  * https://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 rtkJosephForwardProjectionImageFilter_h
20 #define rtkJosephForwardProjectionImageFilter_h
21 
22 #include "rtkConfiguration.h"
24 #include "rtkMacro.h"
25 #include <itkPixelTraits.h>
26 
28 
29 #include <itkVectorImage.h>
30 namespace rtk
31 {
32 namespace Functor
33 {
42 template <class TInput, class TCoordRepType, class TOutput = TInput>
43 class ITK_TEMPLATE_EXPORT InterpolationWeightMultiplication
44 {
45 public:
48  bool
50  {
51  return false;
52  }
53  bool
55  {
56  return !(*this != other);
57  }
58 
59  inline TOutput
60  operator()(const ThreadIdType itkNotUsed(threadId),
61  const double itkNotUsed(stepLengthInVoxel),
62  const TCoordRepType weight,
63  const TInput * p,
64  const int i) const
65  {
66  return weight * p[i];
67  }
68 };
69 
77 template <class TInput, class TOutput>
78 class ITK_TEMPLATE_EXPORT SumAlongRay
79 {
80 public:
82 
83  SumAlongRay() = default;
84  ~SumAlongRay() = default;
85  bool
86  operator!=(const SumAlongRay &) const
87  {
88  return false;
89  }
90  bool
91  operator==(const SumAlongRay & other) const
92  {
93  return !(*this != other);
94  }
95 
96  inline void
97  operator()(const ThreadIdType itkNotUsed(threadId),
98  TOutput & sumValue,
99  const TInput volumeValue,
100  const VectorType & itkNotUsed(stepInMM))
101  {
102  sumValue += static_cast<TOutput>(volumeValue);
103  }
104 };
105 
113 template <class TInput, class TOutput>
114 class ITK_TEMPLATE_EXPORT ProjectedValueAccumulation
115 {
116 public:
118 
119  ProjectedValueAccumulation() = default;
120  ~ProjectedValueAccumulation() = default;
121  bool
123  {
124  return false;
125  }
126  bool
128  {
129  return !(*this != other);
130  }
131 
132  inline void
133  operator()(const ThreadIdType itkNotUsed(threadId),
134  const TInput & input,
135  TOutput & output,
136  const TOutput & rayCastValue,
137  const VectorType & stepInMM,
138  const VectorType & itkNotUsed(source),
139  const VectorType & itkNotUsed(sourceToPixel),
140  const VectorType & itkNotUsed(nearestPoint),
141  const VectorType & itkNotUsed(farthestPoint)) const
142  {
143  output = input + rayCastValue * stepInMM.GetNorm();
144  }
145 };
146 
147 } // end namespace Functor
148 
149 
165 template <class TInputImage,
166  class TOutputImage,
167  class TInterpolationWeightMultiplication = Functor::InterpolationWeightMultiplication<
168  typename TInputImage::PixelType,
170  class TProjectedValueAccumulation =
171  Functor::ProjectedValueAccumulation<typename TInputImage::PixelType, typename TOutputImage::PixelType>,
172  class TSumAlongRay = Functor::SumAlongRay<typename TInputImage::PixelType, typename TOutputImage::PixelType>>
173 class ITK_TEMPLATE_EXPORT JosephForwardProjectionImageFilter
174  : public ForwardProjectionImageFilter<TInputImage, TOutputImage>
175 {
176 public:
177  ITK_DISALLOW_COPY_AND_MOVE(JosephForwardProjectionImageFilter);
178 
184  using InputPixelType = typename TInputImage::PixelType;
185  using OutputPixelType = typename TOutputImage::PixelType;
186  using OutputImageRegionType = typename TOutputImage::RegionType;
187  using CoordRepType = double;
191 
193  itkNewMacro(Self);
194 
196 #ifdef itkOverrideGetNameOfClassMacro
197  itkOverrideGetNameOfClassMacro(JosephForwardProjectionImageFilter);
198 #else
200 #endif
201 
202 
204  TInterpolationWeightMultiplication &
206  {
207  return m_InterpolationWeightMultiplication;
208  }
209  const TInterpolationWeightMultiplication &
211  {
212  return m_InterpolationWeightMultiplication;
213  }
214  void
215  SetInterpolationWeightMultiplication(const TInterpolationWeightMultiplication & _arg)
216  {
217  if (m_InterpolationWeightMultiplication != _arg)
218  {
219  m_InterpolationWeightMultiplication = _arg;
220  this->Modified();
221  }
222  }
224 
227  TProjectedValueAccumulation &
229  {
230  return m_ProjectedValueAccumulation;
231  }
232  const TProjectedValueAccumulation &
234  {
235  return m_ProjectedValueAccumulation;
236  }
237  void
238  SetProjectedValueAccumulation(const TProjectedValueAccumulation & _arg)
239  {
240  if (m_ProjectedValueAccumulation != _arg)
241  {
242  m_ProjectedValueAccumulation = _arg;
243  this->Modified();
244  }
245  }
247 
249  TSumAlongRay &
251  {
252  return m_SumAlongRay;
253  }
254  const TSumAlongRay &
256  {
257  return m_SumAlongRay;
258  }
259  void
260  SetSumAlongRay(const TSumAlongRay & _arg)
261  {
262  if (m_SumAlongRay != _arg)
263  {
264  m_SumAlongRay = _arg;
265  this->Modified();
266  }
267  }
269 
273  void
274  SetInferiorClipImage(const TClipImageType * inferiorClipImage)
275  {
276  // Process object is not const-correct so the const casting is required.
277  this->SetInput("InferiorClipImage", const_cast<TClipImageType *>(inferiorClipImage));
278  }
279  typename TClipImageType::ConstPointer
281  {
282  return static_cast<const TClipImageType *>(this->itk::ProcessObject::GetInput("InferiorClipImage"));
283  }
285 
289  void
290  SetSuperiorClipImage(const TClipImageType * superiorClipImage)
291  {
292  // Process object is not const-correct so the const casting is required.
293  this->SetInput("SuperiorClipImage", const_cast<TClipImageType *>(superiorClipImage));
294  }
295  typename TClipImageType::ConstPointer
297  {
298  return static_cast<const TClipImageType *>(this->itk::ProcessObject::GetInput("SuperiorClipImage"));
299  }
301 
305  itkGetMacro(InferiorClip, double);
306  itkSetMacro(InferiorClip, double);
307  itkGetMacro(SuperiorClip, double);
308  itkSetMacro(SuperiorClip, double);
310 
311 protected:
313  ~JosephForwardProjectionImageFilter() override = default;
314 
316  void
317  GenerateInputRequestedRegion() override;
318 
319  void
320  ThreadedGenerateData(const OutputImageRegionType & outputRegionForThread, ThreadIdType threadId) override;
321 
324  void
325  VerifyInputInformation() const override;
326 
327  inline OutputPixelType
328  BilinearInterpolation(const ThreadIdType threadId,
329  const double stepLengthInVoxel,
330  const InputPixelType * pxiyi,
331  const InputPixelType * pxsyi,
332  const InputPixelType * pxiys,
333  const InputPixelType * pxsys,
334  const double x,
335  const double y,
336  const int ox,
337  const int oy);
338 
339  inline OutputPixelType
340  BilinearInterpolationOnBorders(const ThreadIdType threadId,
341  const double stepLengthInVoxel,
342  const InputPixelType * pxiyi,
343  const InputPixelType * pxsyi,
344  const InputPixelType * pxiys,
345  const InputPixelType * pxsys,
346  const double x,
347  const double y,
348  const int ox,
349  const int oy,
350  const double minx,
351  const double miny,
352  const double maxx,
353  const double maxy);
354 
355 private:
356  // Functors
357  TInterpolationWeightMultiplication m_InterpolationWeightMultiplication;
358  TProjectedValueAccumulation m_ProjectedValueAccumulation;
359  TSumAlongRay m_SumAlongRay;
360  double m_InferiorClip{ 0. };
361  double m_SuperiorClip{ 1. };
362 };
363 
364 } // end namespace rtk
365 
366 #ifndef ITK_MANUAL_INSTANTIATION
367 # include "rtkJosephForwardProjectionImageFilter.hxx"
368 #endif
369 
370 #endif
void SetInterpolationWeightMultiplication(const TInterpolationWeightMultiplication &_arg)
Base class for forward projection, i.e. accumulation along x-ray lines.
void operator()(const ThreadIdType, TOutput &sumValue, const TInput volumeValue, const VectorType &)
Function to accumulate the ray casting on the projection.
const TProjectedValueAccumulation & GetProjectedValueAccumulation() const
Function to compute the attenuation correction on the projection.
bool operator==(const SumAlongRay &other) const
const TInterpolationWeightMultiplication & GetInterpolationWeightMultiplication() const
bool operator!=(const SumAlongRay &) const
typename TPixelType::ValueType ValueType
#define itkSetMacro(name, type)
bool operator!=(const ProjectedValueAccumulation &) const
TOutput operator()(const ThreadIdType, const double, const TCoordRepType weight, const TInput *p, const int i) const
bool operator==(const InterpolationWeightMultiplication &other) const
Function to multiply the interpolation weights with the projected volume values.
typename OutputImageType::RegionType OutputImageRegionType
void SetProjectedValueAccumulation(const TProjectedValueAccumulation &_arg)
TInterpolationWeightMultiplication m_InterpolationWeightMultiplication
void SetSuperiorClipImage(const TClipImageType *superiorClipImage)
bool operator==(const ProjectedValueAccumulation &other) const
DataObject * GetInput(const DataObjectIdentifierType &key)
bool operator!=(const InterpolationWeightMultiplication &) const
TInterpolationWeightMultiplication & GetInterpolationWeightMultiplication()
unsigned int ThreadIdType
void operator()(const ThreadIdType, const TInput &input, TOutput &output, const TOutput &rayCastValue, const VectorType &stepInMM, const VectorType &, const VectorType &, const VectorType &, const VectorType &) const
void SetInferiorClipImage(const TClipImageType *inferiorClipImage)