RTK  1.4.0
Reconstruction Toolkit
rtkBackwardDifferenceDivergenceImageFilter.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 rtkBackwardDifferenceDivergenceImageFilter_h
20 #define rtkBackwardDifferenceDivergenceImageFilter_h
21 
22 #include <itkImageToImageFilter.h>
23 #include <itkCastImageFilter.h>
24 
25 namespace rtk
26 {
39 template <typename TInputImage, typename TOutputImage = itk::Image< float, TInputImage::ImageDimension > >
41  public itk::ImageToImageFilter< TInputImage, TOutputImage >
42 {
43 public:
45  itkStaticConstMacro(InputImageDimension, unsigned int,
46  TInputImage::ImageDimension);
47 
49  typedef TInputImage InputImageType;
50 
56 
58  itkNewMacro(Self)
59 
60 
62 
66  { this->SetUseImageSpacing(true); }
67 
71  { this->SetUseImageSpacing(false); }
72 
75  itkSetMacro(UseImageSpacing, bool);
76  itkGetConstMacro(UseImageSpacing, bool);
78 
81  void SetDimensionsProcessed(bool* DimensionsProcessed);
82 
85 
87  typedef typename InputImageType::PixelType InputPixelType;
88  typedef typename InputImageType::RegionType InputImageRegionType;
89  typedef typename InputImageType::SizeType InputSizeType;
91 
92 protected:
94  virtual ~BackwardDifferenceDivergenceImageFilter() ITK_OVERRIDE;
95 
96  void GenerateInputRequestedRegion() ITK_OVERRIDE;
97 
98  void BeforeThreadedGenerateData() ITK_OVERRIDE;
99 
100 #if ITK_VERSION_MAJOR<5
101  void ThreadedGenerateData(const typename InputImageType::RegionType& outputRegionForThread, itk::ThreadIdType itkNotUsed(threadId)) ITK_OVERRIDE;
102 #else
103  void DynamicThreadedGenerateData(const typename InputImageType::RegionType& outputRegionForThread) ITK_OVERRIDE;
104 #endif
105 
106  void AfterThreadedGenerateData() ITK_OVERRIDE;
107 
108 private:
109  BackwardDifferenceDivergenceImageFilter(const Self&); //purposely not implemented
110  void operator=(const Self&); //purposely not implemented
111 
113  typename TInputImage::SpacingType m_InvSpacingCoeffs;
114 
115  // list of the dimensions along which the divergence has
116  // to be computed. The components on other dimensions
117  // are ignored for performance, but the gradient filter
118  // sets them to zero anyway
119  bool m_DimensionsProcessed[TInputImage::ImageDimension];
120 
121  // The default is ConstantBoundaryCondition, but this behavior sometimes needs to be overriden
122  itk::ImageBoundaryCondition< TInputImage, TInputImage >* m_BoundaryCondition;
123  // If so, do not perform boundary processing in AfterThreadedGenerateData
125 };
126 
127 } // end namespace itk
128 
129 #ifndef ITK_MANUAL_INSTANTIATION
130 #include "rtkBackwardDifferenceDivergenceImageFilter.hxx"
131 #endif
132 
133 #endif //__rtkBackwardDifferenceDivergenceImageFilter__
void OverrideBoundaryCondition(itk::ImageBoundaryCondition< TInputImage > *boundaryCondition)
itk::CovariantVector< InputPixelType, InputImageDimension > CovariantVectorType
void ThreadedGenerateData(const typename InputImageType::RegionType &outputRegionForThread, itk::ThreadIdType) override
void SetDimensionsProcessed(bool *DimensionsProcessed)
typename FixedImageType::SpacingType SpacingType
virtual ~BackwardDifferenceDivergenceImageFilter() override
unsigned int ThreadIdType
itk::ImageToImageFilter< InputImageType, TOutputImage > Superclass
Computes the backward differences divergence (adjoint of the forward differences gradient) of the inp...
itk::ImageBoundaryCondition< TInputImage, TInputImage > * m_BoundaryCondition
#define itkSetMacro(name, type)