RTK  2.0.1
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:
49  return false;
50  }
52  {
53  return !( *this != other );
54  }
55 
56  inline TOutput operator()( const ThreadIdType itkNotUsed(threadId),
57  const double itkNotUsed(stepLengthInVoxel),
58  const TCoordRepType weight,
59  const TInput *p,
60  const int i ) const
61  {
62  return weight*p[i];
63  }
64 };
65 
73 template< class TInput, class TOutput>
75 {
76 public:
78 
79  SumAlongRay()= default;
80  ~SumAlongRay() = default;
81  bool operator!=( const SumAlongRay & ) const
82  {
83  return false;
84  }
85  bool operator==(const SumAlongRay & other) const
86  {
87  return !( *this != other );
88  }
89 
90  inline TOutput operator()(const ThreadIdType itkNotUsed(threadId),
91  const TInput volumeValue,
92  const VectorType &itkNotUsed(stepInMM))
93  {
94  return volumeValue;
95  }
96 };
97 
105 template< class TInput, class TOutput >
107 {
108 public:
110 
111  ProjectedValueAccumulation() = default;
112  ~ProjectedValueAccumulation() = default;
114  {
115  return false;
116  }
117  bool operator==(const ProjectedValueAccumulation & other) const
118  {
119  return !( *this != other );
120  }
121 
122  inline void operator()( const ThreadIdType itkNotUsed(threadId),
123  const TInput &input,
124  TOutput &output,
125  const TOutput &rayCastValue,
126  const VectorType &stepInMM,
127  const VectorType &itkNotUsed(source),
128  const VectorType &itkNotUsed(sourceToPixel),
129  const VectorType &itkNotUsed(nearestPoint),
130  const VectorType &itkNotUsed(farthestPoint)) const
131  {
132  output = input + rayCastValue * stepInMM.GetNorm();
133  }
134 };
135 
136 } // end namespace Functor
137 
138 
154 template <class TInputImage,
155  class TOutputImage,
159  >
161  public ForwardProjectionImageFilter<TInputImage,TOutputImage>
162 {
163 public:
164  ITK_DISALLOW_COPY_AND_ASSIGN(JosephForwardProjectionImageFilter);
165 
171  using InputPixelType = typename TInputImage::PixelType;
172  using OutputPixelType = typename TOutputImage::PixelType;
173  using OutputImageRegionType = typename TOutputImage::RegionType;
174  using CoordRepType = double;
176 
178  itkNewMacro(Self);
179 
182 
184  TInterpolationWeightMultiplication & GetInterpolationWeightMultiplication() { return m_InterpolationWeightMultiplication; }
185  const TInterpolationWeightMultiplication & GetInterpolationWeightMultiplication() const { return m_InterpolationWeightMultiplication; }
186  void SetInterpolationWeightMultiplication(const TInterpolationWeightMultiplication & _arg)
187  {
188  if ( m_InterpolationWeightMultiplication != _arg )
189  {
190  m_InterpolationWeightMultiplication = _arg;
191  this->Modified();
192  }
193  }
195 
198  TProjectedValueAccumulation & GetProjectedValueAccumulation() { return m_ProjectedValueAccumulation; }
199  const TProjectedValueAccumulation & GetProjectedValueAccumulation() const { return m_ProjectedValueAccumulation; }
200  void SetProjectedValueAccumulation(const TProjectedValueAccumulation & _arg)
201  {
202  if ( m_ProjectedValueAccumulation != _arg )
203  {
204  m_ProjectedValueAccumulation = _arg;
205  this->Modified();
206  }
207  }
209 
211  TSumAlongRay & GetSumAlongRay() { return m_SumAlongRay; }
212  const TSumAlongRay & GetSumAlongRay() const { return m_SumAlongRay; }
213  void SetSumAlongRay(const TSumAlongRay & _arg)
214  {
215  if ( m_SumAlongRay != _arg )
216  {
217  m_SumAlongRay = _arg;
218  this->Modified();
219  }
220  }
222 
226  itkGetMacro(InferiorClip, double)
227  itkSetMacro(InferiorClip, double)
228  itkGetMacro(SuperiorClip, double)
229  itkSetMacro(SuperiorClip, double)
230 
231 protected:
233  ~JosephForwardProjectionImageFilter() override = default;
234 
235  void ThreadedGenerateData( const OutputImageRegionType& outputRegionForThread, ThreadIdType threadId ) override;
236 
239 #if ITK_VERSION_MAJOR<5
240  void VerifyInputInformation() override {}
241 #else
242  void VerifyInputInformation() const override {}
243 #endif
244 
245 
246  inline OutputPixelType BilinearInterpolation(const ThreadIdType threadId,
247  const double stepLengthInVoxel,
248  const InputPixelType *pxiyi,
249  const InputPixelType *pxsyi,
250  const InputPixelType *pxiys,
251  const InputPixelType *pxsys,
252  const double x,
253  const double y,
254  const int ox,
255  const int oy);
256 
257  inline OutputPixelType BilinearInterpolationOnBorders(const ThreadIdType threadId,
258  const double stepLengthInVoxel,
259  const InputPixelType *pxiyi,
260  const InputPixelType *pxsyi,
261  const InputPixelType *pxiys,
262  const InputPixelType *pxsys,
263  const double x,
264  const double y,
265  const int ox,
266  const int oy,
267  const double minx,
268  const double miny,
269  const double maxx,
270  const double maxy);
271 
272 private:
273  // Functors
274  TInterpolationWeightMultiplication m_InterpolationWeightMultiplication;
275  TProjectedValueAccumulation m_ProjectedValueAccumulation;
276  TSumAlongRay m_SumAlongRay;
277  double m_InferiorClip{0.};
278  double m_SuperiorClip{1.};
279 };
280 
281 } // end namespace rtk
282 
283 #ifndef ITK_MANUAL_INSTANTIATION
284 #include "rtkJosephForwardProjectionImageFilter.hxx"
285 #endif
286 
287 #endif
void SetInterpolationWeightMultiplication(const TInterpolationWeightMultiplication &_arg)
bool operator==(const InterpolationWeightMultiplication &other) const
Base class for forward projection, i.e. accumulation along x-ray lines.
bool operator!=(const ProjectedValueAccumulation &) const
TOutput operator()(const ThreadIdType, const double, const TCoordRepType weight, const TInput *p, const int i) const
Function to accumulate the ray casting on the projection.
Function to compute the attenuation correction on the projection.
TOutput operator()(const ThreadIdType, const TInput volumeValue, const VectorType &)
bool operator!=(const SumAlongRay &) const
Function to multiply the interpolation weights with the projected volume values.
void SetProjectedValueAccumulation(const TProjectedValueAccumulation &_arg)
TInterpolationWeightMultiplication m_InterpolationWeightMultiplication
unsigned int ThreadIdType
const TInterpolationWeightMultiplication & GetInterpolationWeightMultiplication() const
TInterpolationWeightMultiplication & GetInterpolationWeightMultiplication()
const TProjectedValueAccumulation & GetProjectedValueAccumulation() const
bool operator==(const ProjectedValueAccumulation &other) const
bool operator!=(const InterpolationWeightMultiplication &) const
void operator()(const ThreadIdType, const TInput &input, TOutput &output, const TOutput &rayCastValue, const VectorType &stepInMM, const VectorType &, const VectorType &, const VectorType &, const VectorType &) const
#define itkSetMacro(name, type)
bool operator==(const SumAlongRay &other) const