RTK  2.1.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  * 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 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>
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>
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>
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_EXPORT JosephForwardProjectionImageFilter : public ForwardProjectionImageFilter<TInputImage, TOutputImage>
171 {
172 public:
173 #if ITK_VERSION_MAJOR == 5 && ITK_VERSION_MINOR == 1
174  ITK_DISALLOW_COPY_AND_ASSIGN(JosephForwardProjectionImageFilter);
175 #else
176  ITK_DISALLOW_COPY_AND_MOVE(JosephForwardProjectionImageFilter);
177 #endif
178 
184  using InputPixelType = typename TInputImage::PixelType;
185  using OutputPixelType = typename TOutputImage::PixelType;
186  using OutputImageRegionType = typename TOutputImage::RegionType;
187  using CoordRepType = double;
189 
191  itkNewMacro(Self);
192 
195 
197  TInterpolationWeightMultiplication &
199  {
200  return m_InterpolationWeightMultiplication;
201  }
202  const TInterpolationWeightMultiplication &
204  {
205  return m_InterpolationWeightMultiplication;
206  }
207  void
208  SetInterpolationWeightMultiplication(const TInterpolationWeightMultiplication & _arg)
209  {
210  if (m_InterpolationWeightMultiplication != _arg)
211  {
212  m_InterpolationWeightMultiplication = _arg;
213  this->Modified();
214  }
215  }
217 
220  TProjectedValueAccumulation &
222  {
223  return m_ProjectedValueAccumulation;
224  }
225  const TProjectedValueAccumulation &
227  {
228  return m_ProjectedValueAccumulation;
229  }
230  void
231  SetProjectedValueAccumulation(const TProjectedValueAccumulation & _arg)
232  {
233  if (m_ProjectedValueAccumulation != _arg)
234  {
235  m_ProjectedValueAccumulation = _arg;
236  this->Modified();
237  }
238  }
240 
242  TSumAlongRay &
244  {
245  return m_SumAlongRay;
246  }
247  const TSumAlongRay &
249  {
250  return m_SumAlongRay;
251  }
252  void
253  SetSumAlongRay(const TSumAlongRay & _arg)
254  {
255  if (m_SumAlongRay != _arg)
256  {
257  m_SumAlongRay = _arg;
258  this->Modified();
259  }
260  }
262 
266  itkGetMacro(InferiorClip, double);
267  itkSetMacro(InferiorClip, double);
268  itkGetMacro(SuperiorClip, double);
269  itkSetMacro(SuperiorClip, double);
271 
272 protected:
274  ~JosephForwardProjectionImageFilter() override = default;
275 
276  void
277  ThreadedGenerateData(const OutputImageRegionType & outputRegionForThread, ThreadIdType threadId) override;
278 
281  void
282  VerifyInputInformation() const override
283  {}
284 
285  inline OutputPixelType
286  BilinearInterpolation(const ThreadIdType threadId,
287  const double stepLengthInVoxel,
288  const InputPixelType * pxiyi,
289  const InputPixelType * pxsyi,
290  const InputPixelType * pxiys,
291  const InputPixelType * pxsys,
292  const double x,
293  const double y,
294  const int ox,
295  const int oy);
296 
297  inline OutputPixelType
298  BilinearInterpolationOnBorders(const ThreadIdType threadId,
299  const double stepLengthInVoxel,
300  const InputPixelType * pxiyi,
301  const InputPixelType * pxsyi,
302  const InputPixelType * pxiys,
303  const InputPixelType * pxsys,
304  const double x,
305  const double y,
306  const int ox,
307  const int oy,
308  const double minx,
309  const double miny,
310  const double maxx,
311  const double maxy);
312 
313 private:
314  // Functors
315  TInterpolationWeightMultiplication m_InterpolationWeightMultiplication;
316  TProjectedValueAccumulation m_ProjectedValueAccumulation;
317  TSumAlongRay m_SumAlongRay;
318  double m_InferiorClip{ 0. };
319  double m_SuperiorClip{ 1. };
320 };
321 
322 } // end namespace rtk
323 
324 #ifndef ITK_MANUAL_INSTANTIATION
325 # include "rtkJosephForwardProjectionImageFilter.hxx"
326 #endif
327 
328 #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
bool operator==(const ProjectedValueAccumulation &other) const
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