RTK  2.5.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 TOutput
97  operator()(const ThreadIdType itkNotUsed(threadId), const TInput volumeValue, const VectorType & itkNotUsed(stepInMM))
98  {
99  return volumeValue;
100  }
101 };
102 
110 template <class TInput, class TOutput>
111 class ITK_TEMPLATE_EXPORT ProjectedValueAccumulation
112 {
113 public:
115 
116  ProjectedValueAccumulation() = default;
117  ~ProjectedValueAccumulation() = default;
118  bool
120  {
121  return false;
122  }
123  bool
125  {
126  return !(*this != other);
127  }
128 
129  inline void
130  operator()(const ThreadIdType itkNotUsed(threadId),
131  const TInput & input,
132  TOutput & output,
133  const TOutput & rayCastValue,
134  const VectorType & stepInMM,
135  const VectorType & itkNotUsed(source),
136  const VectorType & itkNotUsed(sourceToPixel),
137  const VectorType & itkNotUsed(nearestPoint),
138  const VectorType & itkNotUsed(farthestPoint)) const
139  {
140  output = input + rayCastValue * stepInMM.GetNorm();
141  }
142 };
143 
144 } // end namespace Functor
145 
146 
162 template <class TInputImage,
163  class TOutputImage,
164  class TInterpolationWeightMultiplication = Functor::InterpolationWeightMultiplication<
165  typename TInputImage::PixelType,
167  class TProjectedValueAccumulation =
168  Functor::ProjectedValueAccumulation<typename TInputImage::PixelType, typename TOutputImage::PixelType>,
169  class TSumAlongRay = Functor::SumAlongRay<typename TInputImage::PixelType, typename TOutputImage::PixelType>>
170 class ITK_TEMPLATE_EXPORT JosephForwardProjectionImageFilter
171  : public ForwardProjectionImageFilter<TInputImage, TOutputImage>
172 {
173 public:
174  ITK_DISALLOW_COPY_AND_MOVE(JosephForwardProjectionImageFilter);
175 
181  using InputPixelType = typename TInputImage::PixelType;
182  using OutputPixelType = typename TOutputImage::PixelType;
183  using OutputImageRegionType = typename TOutputImage::RegionType;
184  using CoordRepType = double;
188 
190  itkNewMacro(Self);
191 
193 #ifdef itkOverrideGetNameOfClassMacro
195 #else
197 #endif
198 
199 
201  TInterpolationWeightMultiplication &
203  {
204  return m_InterpolationWeightMultiplication;
205  }
206  const TInterpolationWeightMultiplication &
208  {
209  return m_InterpolationWeightMultiplication;
210  }
211  void
212  SetInterpolationWeightMultiplication(const TInterpolationWeightMultiplication & _arg)
213  {
214  if (m_InterpolationWeightMultiplication != _arg)
215  {
216  m_InterpolationWeightMultiplication = _arg;
217  this->Modified();
218  }
219  }
221 
224  TProjectedValueAccumulation &
226  {
227  return m_ProjectedValueAccumulation;
228  }
229  const TProjectedValueAccumulation &
231  {
232  return m_ProjectedValueAccumulation;
233  }
234  void
235  SetProjectedValueAccumulation(const TProjectedValueAccumulation & _arg)
236  {
237  if (m_ProjectedValueAccumulation != _arg)
238  {
239  m_ProjectedValueAccumulation = _arg;
240  this->Modified();
241  }
242  }
244 
246  TSumAlongRay &
248  {
249  return m_SumAlongRay;
250  }
251  const TSumAlongRay &
253  {
254  return m_SumAlongRay;
255  }
256  void
257  SetSumAlongRay(const TSumAlongRay & _arg)
258  {
259  if (m_SumAlongRay != _arg)
260  {
261  m_SumAlongRay = _arg;
262  this->Modified();
263  }
264  }
266 
270  void
271  SetInferiorClipImage(const TClipImageType * inferiorClipImage)
272  {
273  // Process object is not const-correct so the const casting is required.
274  this->SetInput("InferiorClipImage", const_cast<TClipImageType *>(inferiorClipImage));
275  }
276  typename TClipImageType::ConstPointer
278  {
279  return static_cast<const TClipImageType *>(this->itk::ProcessObject::GetInput("InferiorClipImage"));
280  }
282 
286  void
287  SetSuperiorClipImage(const TClipImageType * superiorClipImage)
288  {
289  // Process object is not const-correct so the const casting is required.
290  this->SetInput("SuperiorClipImage", const_cast<TClipImageType *>(superiorClipImage));
291  }
292  typename TClipImageType::ConstPointer
294  {
295  return static_cast<const TClipImageType *>(this->itk::ProcessObject::GetInput("SuperiorClipImage"));
296  }
298 
302  itkGetMacro(InferiorClip, double);
303  itkSetMacro(InferiorClip, double);
304  itkGetMacro(SuperiorClip, double);
305  itkSetMacro(SuperiorClip, double);
307 
308 protected:
310  ~JosephForwardProjectionImageFilter() override = default;
311 
313  void
314  GenerateInputRequestedRegion() override;
315 
316  void
317  ThreadedGenerateData(const OutputImageRegionType & outputRegionForThread, ThreadIdType threadId) override;
318 
321  void
322  VerifyInputInformation() const override;
323 
324  inline OutputPixelType
325  BilinearInterpolation(const ThreadIdType threadId,
326  const double stepLengthInVoxel,
327  const InputPixelType * pxiyi,
328  const InputPixelType * pxsyi,
329  const InputPixelType * pxiys,
330  const InputPixelType * pxsys,
331  const double x,
332  const double y,
333  const int ox,
334  const int oy);
335 
336  inline OutputPixelType
337  BilinearInterpolationOnBorders(const ThreadIdType threadId,
338  const double stepLengthInVoxel,
339  const InputPixelType * pxiyi,
340  const InputPixelType * pxsyi,
341  const InputPixelType * pxiys,
342  const InputPixelType * pxsys,
343  const double x,
344  const double y,
345  const int ox,
346  const int oy,
347  const double minx,
348  const double miny,
349  const double maxx,
350  const double maxy);
351 
352 private:
353  // Functors
354  TInterpolationWeightMultiplication m_InterpolationWeightMultiplication;
355  TProjectedValueAccumulation m_ProjectedValueAccumulation;
356  TSumAlongRay m_SumAlongRay;
357  double m_InferiorClip{ 0. };
358  double m_SuperiorClip{ 1. };
359 };
360 
361 } // end namespace rtk
362 
363 #ifndef ITK_MANUAL_INSTANTIATION
364 # include "rtkJosephForwardProjectionImageFilter.hxx"
365 #endif
366 
367 #endif
void SetInterpolationWeightMultiplication(const TInterpolationWeightMultiplication &_arg)
Base class for forward projection, i.e. accumulation along x-ray lines.
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
TOutput operator()(const ThreadIdType, const TInput volumeValue, const VectorType &)
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
itkOverrideGetNameOfClassMacro(AddImageFilter)
void SetInferiorClipImage(const TClipImageType *inferiorClipImage)