RTK  1.4.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:
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 
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 
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:
169  typedef typename TInputImage::PixelType InputPixelType;
170  typedef typename TOutputImage::PixelType OutputPixelType;
171  typedef typename TOutputImage::RegionType OutputImageRegionType;
172  typedef double CoordRepType;
174 
176  itkNewMacro(Self);
177 
180 
182  TInterpolationWeightMultiplication & GetInterpolationWeightMultiplication() { return m_InterpolationWeightMultiplication; }
183  const TInterpolationWeightMultiplication & GetInterpolationWeightMultiplication() const { return m_InterpolationWeightMultiplication; }
184  void SetInterpolationWeightMultiplication(const TInterpolationWeightMultiplication & _arg)
185  {
186  if ( m_InterpolationWeightMultiplication != _arg )
187  {
188  m_InterpolationWeightMultiplication = _arg;
189  this->Modified();
190  }
191  }
193 
196  TProjectedValueAccumulation & GetProjectedValueAccumulation() { return m_ProjectedValueAccumulation; }
197  const TProjectedValueAccumulation & GetProjectedValueAccumulation() const { return m_ProjectedValueAccumulation; }
198  void SetProjectedValueAccumulation(const TProjectedValueAccumulation & _arg)
199  {
200  if ( m_ProjectedValueAccumulation != _arg )
201  {
202  m_ProjectedValueAccumulation = _arg;
203  this->Modified();
204  }
205  }
207 
209  TSumAlongRay & GetSumAlongRay() { return m_SumAlongRay; }
210  const TSumAlongRay & GetSumAlongRay() const { return m_SumAlongRay; }
211  void SetSumAlongRay(const TSumAlongRay & _arg)
212  {
213  if ( m_SumAlongRay != _arg )
214  {
215  m_SumAlongRay = _arg;
216  this->Modified();
217  }
218  }
220 
224  itkGetMacro(InferiorClip, double)
225  itkSetMacro(InferiorClip, double)
226  itkGetMacro(SuperiorClip, double)
227  itkSetMacro(SuperiorClip, double)
228 
229 protected:
231  virtual ~JosephForwardProjectionImageFilter() ITK_OVERRIDE {}
232 
233  void ThreadedGenerateData( const OutputImageRegionType& outputRegionForThread, ThreadIdType threadId ) ITK_OVERRIDE;
234 
237  void VerifyInputInformation() ITK_OVERRIDE {}
238 
239  inline OutputPixelType BilinearInterpolation(const ThreadIdType threadId,
240  const double stepLengthInVoxel,
241  const InputPixelType *pxiyi,
242  const InputPixelType *pxsyi,
243  const InputPixelType *pxiys,
244  const InputPixelType *pxsys,
245  const double x,
246  const double y,
247  const int ox,
248  const int oy);
249 
250  inline OutputPixelType BilinearInterpolationOnBorders(const ThreadIdType threadId,
251  const double stepLengthInVoxel,
252  const InputPixelType *pxiyi,
253  const InputPixelType *pxsyi,
254  const InputPixelType *pxiys,
255  const InputPixelType *pxsys,
256  const double x,
257  const double y,
258  const int ox,
259  const int oy,
260  const double minx,
261  const double miny,
262  const double maxx,
263  const double maxy);
264 
265 private:
266  JosephForwardProjectionImageFilter(const Self&); //purposely not implemented
267  void operator=(const Self&); //purposely not implemented
268 
269  // Functors
270  TInterpolationWeightMultiplication m_InterpolationWeightMultiplication;
271  TProjectedValueAccumulation m_ProjectedValueAccumulation;
272  TSumAlongRay m_SumAlongRay;
275 };
276 
277 } // end namespace rtk
278 
279 #ifndef ITK_MANUAL_INSTANTIATION
280 #include "rtkJosephForwardProjectionImageFilter.hxx"
281 #endif
282 
283 #endif
void SetInterpolationWeightMultiplication(const TInterpolationWeightMultiplication &_arg)
bool operator==(const InterpolationWeightMultiplication &other) const
Base class for forward projection, i.e. accumulation along x-ray lines.
ForwardProjectionImageFilter< TInputImage, TOutputImage > Superclass
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
itk::Vector< CoordRepType, TInputImage::ImageDimension > VectorType
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