RTK  2.6.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  * https://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>>
40 class ITK_TEMPLATE_EXPORT BackwardDifferenceDivergenceImageFilter
41  : public itk::ImageToImageFilter<TInputImage, TOutputImage>
42 {
43 public:
44  ITK_DISALLOW_COPY_AND_MOVE(BackwardDifferenceDivergenceImageFilter);
45 
47  static constexpr unsigned int InputImageDimension = TInputImage::ImageDimension;
48 
50  using InputImageType = TInputImage;
51 
57 
59  itkNewMacro(Self);
60 
62 #ifdef itkOverrideGetNameOfClassMacro
63  itkOverrideGetNameOfClassMacro(BackwardDifferenceDivergenceImageFilter);
64 #else
66 #endif
67 
68 
71  void
73  {
74  this->SetUseImageSpacing(true);
75  }
76 
79  void
81  {
82  this->SetUseImageSpacing(false);
83  }
84 
87  itkSetMacro(UseImageSpacing, bool);
88  itkGetConstMacro(UseImageSpacing, bool);
90 
93  void
94  SetDimensionsProcessed(bool * DimensionsProcessed);
95 
97  void
98  OverrideBoundaryCondition(itk::ImageBoundaryCondition<TInputImage> * boundaryCondition);
99 
101  using InputPixelType = typename InputImageType::PixelType;
102  using InputImageRegionType = typename InputImageType::RegionType;
103  using InputSizeType = typename InputImageType::SizeType;
105 
106 protected:
109 
110  void
111  GenerateInputRequestedRegion() override;
112 
113  void
114  BeforeThreadedGenerateData() override;
115 
116  void
117  DynamicThreadedGenerateData(const typename InputImageType::RegionType & outputRegionForThread) override;
118 
119  void
120  AfterThreadedGenerateData() override;
121 
122 private:
124  typename TInputImage::SpacingType m_InvSpacingCoeffs;
125 
126  // list of the dimensions along which the divergence has
127  // to be computed. The components on other dimensions
128  // are ignored for performance, but the gradient filter
129  // sets them to zero anyway
130  bool m_DimensionsProcessed[TInputImage::ImageDimension];
131 
132  // The default is ConstantBoundaryCondition, but this behavior sometimes needs to be overriden
134  // If so, do not perform boundary processing in AfterThreadedGenerateData
136 };
137 
138 } // namespace rtk
139 
140 #ifndef ITK_MANUAL_INSTANTIATION
141 # include "rtkBackwardDifferenceDivergenceImageFilter.hxx"
142 #endif
143 
144 #endif //__rtkBackwardDifferenceDivergenceImageFilter__
TInputImage InputImageType
#define itkSetMacro(name, type)
itk::ImageBoundaryCondition< TInputImage, TInputImage > * m_BoundaryCondition
typename InputImageType::RegionType InputImageRegionType
Computes the backward differences divergence (adjoint of the forward differences gradient) of the inp...