RTK  1.4.0
Reconstruction Toolkit
rtkBackProjectionImageFilter.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 rtkBackProjectionImageFilter_h
20 #define rtkBackProjectionImageFilter_h
21 
22 #include "rtkConfiguration.h"
23 
24 #include <itkInPlaceImageFilter.h>
25 #include <itkConceptChecking.h>
26 
28 
29 #include <type_traits>
30 #include <typeinfo>
31 
32 namespace rtk
33 {
34 
49 template <class TInputImage, class TOutputImage>
51  public itk::InPlaceImageFilter<TInputImage,TOutputImage>
52 {
53 public:
59  typedef typename TInputImage::PixelType InputPixelType;
60  typedef typename TInputImage::InternalPixelType InternalInputPixelType;
61  typedef typename TOutputImage::RegionType OutputImageRegionType;
62 
66  typedef itk::Image<InputPixelType, TInputImage::ImageDimension-1> ProjectionImageType;
67  typedef typename ProjectionImageType::Pointer ProjectionImagePointer;
68 
70  itkNewMacro(Self);
71 
74 
76  itkGetConstObjectMacro(Geometry, GeometryType);
77  itkSetConstObjectMacro(Geometry, GeometryType);
79 
81  itkGetMacro(Transpose, bool);
82  itkSetMacro(Transpose, bool);
84 
85 protected:
87  this->SetNumberOfRequiredInputs(2); this->SetInPlace( true );
88  };
89  virtual ~BackProjectionImageFilter() ITK_OVERRIDE {}
90 
92  void GenerateInputRequestedRegion() ITK_OVERRIDE;
93 
94  void BeforeThreadedGenerateData() ITK_OVERRIDE;
95 
96 #if ITK_VERSION_MAJOR<5
97  void ThreadedGenerateData( const OutputImageRegionType& outputRegionForThread, ThreadIdType threadId ) ITK_OVERRIDE;
98 #else
99  void DynamicThreadedGenerateData(const OutputImageRegionType& outputRegionForThread) ITK_OVERRIDE;
100 #endif
101 
103  virtual void CylindricalDetectorCenteredOnSourceBackprojection(const OutputImageRegionType& region,
104  const ProjectionMatrixType& volIndexToProjPP,
106  const ProjectionImagePointer projection);
107 
110  virtual void OptimizedBackprojectionX(const OutputImageRegionType& region, const ProjectionMatrixType& matrix,
111  const ProjectionImagePointer projection);
112 
115  virtual void OptimizedBackprojectionY(const OutputImageRegionType& region, const ProjectionMatrixType& matrix,
116  const ProjectionImagePointer projection);
117 
120 #if ITK_VERSION_MAJOR<5
121  void VerifyInputInformation() ITK_OVERRIDE {}
122 #else
123  void VerifyInputInformation() const ITK_OVERRIDE {}
124 #endif
125 
126 
131  template<class TProjectionImage>
132  typename TProjectionImage::Pointer GetProjection(const unsigned int iProj);
133 
136  ProjectionMatrixType GetIndexToIndexProjectionMatrix(const unsigned int iProj);
137 
138  ProjectionMatrixType GetVolumeIndexToProjectionPhysicalPointMatrix(const unsigned int iProj);
139 
141 
143  GeometryConstPointer m_Geometry;
144 
145 private:
146  BackProjectionImageFilter(const Self&); //purposely not implemented
147  void operator=(const Self&); //purposely not implemented
148 
152 };
153 
154 } // end namespace rtk
155 
156 #ifndef ITK_MANUAL_INSTANTIATION
157 #include "rtkBackProjectionImageFilter.hxx"
158 #endif
159 
160 #endif
itk::Matrix< double, TInputImage::ImageDimension, TInputImage::ImageDimension > GetProjectionPhysicalPointToProjectionIndexMatrix()
void GenerateInputRequestedRegion() override
virtual void CylindricalDetectorCenteredOnSourceBackprojection(const OutputImageRegionType &region, const ProjectionMatrixType &volIndexToProjPP, const itk::Matrix< double, TInputImage::ImageDimension, TInputImage::ImageDimension > &projPPToProjIndex, const ProjectionImagePointer projection)
itk::ImageToImageFilter< TInputImage, TOutputImage > Superclass
void BeforeThreadedGenerateData() override
itk::Image< InputPixelType, TInputImage::ImageDimension-1 > ProjectionImageType
GeometryType::MatrixType ProjectionMatrixType
TOutputImage::RegionType OutputImageRegionType
ProjectionMatrixType GetVolumeIndexToProjectionPhysicalPointMatrix(const unsigned int iProj)
Projection geometry for a source and a 2-D flat panel.
void ThreadedGenerateData(const OutputImageRegionType &outputRegionForThread, ThreadIdType threadId) override
TInputImage::InternalPixelType InternalInputPixelType
virtual void OptimizedBackprojectionY(const OutputImageRegionType &region, const ProjectionMatrixType &matrix, const ProjectionImagePointer projection)
unsigned int ThreadIdType
itk::SmartPointer< const Self > ConstPointer
GeometryType::ConstPointer GeometryConstPointer
ProjectionMatrixType GetIndexToIndexProjectionMatrix(const unsigned int iProj)
void operator=(const Self &)
TProjectionImage::Pointer GetProjection(const unsigned int iProj)
rtk::ThreeDCircularProjectionGeometry GeometryType
#define itkSetMacro(name, type)
ProjectionImageType::Pointer ProjectionImagePointer
virtual void OptimizedBackprojectionX(const OutputImageRegionType &region, const ProjectionMatrixType &matrix, const ProjectionImagePointer projection)