RTK  2.0.1
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:
44  ITK_DISALLOW_COPY_AND_ASSIGN(BackwardDifferenceDivergenceImageFilter);
45 
47  itkStaticConstMacro(InputImageDimension, unsigned int,
48  TInputImage::ImageDimension);
49 
51  using InputImageType = TInputImage;
52 
58 
60  itkNewMacro(Self)
61 
62 
64 
68  { this->SetUseImageSpacing(true); }
69 
73  { this->SetUseImageSpacing(false); }
74 
77  itkSetMacro(UseImageSpacing, bool);
78  itkGetConstMacro(UseImageSpacing, bool);
80 
83  void SetDimensionsProcessed(bool* DimensionsProcessed);
84 
87 
89  using InputPixelType = typename InputImageType::PixelType;
90  using InputImageRegionType = typename InputImageType::RegionType;
91  using InputSizeType = typename InputImageType::SizeType;
93 
94 protected:
97 
98  void GenerateInputRequestedRegion() override;
99 
100  void BeforeThreadedGenerateData() override;
101 
102 #if ITK_VERSION_MAJOR<5
103  void ThreadedGenerateData(const typename InputImageType::RegionType& outputRegionForThread, itk::ThreadIdType itkNotUsed(threadId)) override;
104 #else
105  void DynamicThreadedGenerateData(const typename InputImageType::RegionType& outputRegionForThread) override;
106 #endif
107 
108  void AfterThreadedGenerateData() override;
109 
110 private:
112  typename TInputImage::SpacingType m_InvSpacingCoeffs;
113 
114  // list of the dimensions along which the divergence has
115  // to be computed. The components on other dimensions
116  // are ignored for performance, but the gradient filter
117  // sets them to zero anyway
118  bool m_DimensionsProcessed[TInputImage::ImageDimension];
119 
120  // The default is ConstantBoundaryCondition, but this behavior sometimes needs to be overriden
122  // If so, do not perform boundary processing in AfterThreadedGenerateData
124 };
125 
126 } // end namespace itk
127 
128 #ifndef ITK_MANUAL_INSTANTIATION
129 #include "rtkBackwardDifferenceDivergenceImageFilter.hxx"
130 #endif
131 
132 #endif //__rtkBackwardDifferenceDivergenceImageFilter__
void OverrideBoundaryCondition(itk::ImageBoundaryCondition< TInputImage > *boundaryCondition)
void ThreadedGenerateData(const typename InputImageType::RegionType &outputRegionForThread, itk::ThreadIdType) override
void SetDimensionsProcessed(bool *DimensionsProcessed)
unsigned int ThreadIdType
Computes the backward differences divergence (adjoint of the forward differences gradient) of the inp...
itk::ImageBoundaryCondition< TInputImage, TInputImage > * m_BoundaryCondition
#define itkSetMacro(name, type)