ITK  5.0.0
Insight Segmentation and Registration Toolkit
itkGradientImageFilter.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright Insight Software 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 #ifndef itkGradientImageFilter_h
19 #define itkGradientImageFilter_h
20 
21 #include "itkImageToImageFilter.h"
22 #include "itkCovariantVector.h"
23 #include "itkImageRegionIterator.h"
24 
25 namespace itk
26 {
27 
28 
29 template <typename TPixelType, unsigned int VImageDimension > class VectorImage;
30 
31 
54 template< typename TInputImage,
55  typename TOperatorValueType = float,
56  typename TOutputValueType = float,
57  typename TOutputImageType = Image< CovariantVector< TOutputValueType,
58  TInputImage::ImageDimension >,
59  TInputImage::ImageDimension > >
60 class ITK_TEMPLATE_EXPORT GradientImageFilter:
61  public ImageToImageFilter< TInputImage, TOutputImageType >
62 {
63 public:
64  ITK_DISALLOW_COPY_AND_ASSIGN(GradientImageFilter);
65 
67  static constexpr unsigned int InputImageDimension = TInputImage::ImageDimension;
68  static constexpr unsigned int OutputImageDimension = TOutputImageType::ImageDimension;
69 
71  using InputImageType = TInputImage;
72  using InputImagePointer = typename InputImageType::Pointer;
73  using OutputImageType = TOutputImageType;
74  using OutputImagePointer = typename OutputImageType::Pointer;
75 
81 
83  itkNewMacro(Self);
84 
87 
89  using InputPixelType = typename InputImageType::PixelType;
90  using OperatorValueType = TOperatorValueType;
91  using OutputValueType = TOutputValueType;
92  using OutputPixelType = typename OutputImageType::PixelType;
94  OutputValueType, Self::OutputImageDimension >;
96 
103  void GenerateInputRequestedRegion() override;
104 
108  { this->SetUseImageSpacing(true); }
109 
113  { this->SetUseImageSpacing(false); }
114 
117  itkSetMacro(UseImageSpacing, bool);
118  itkGetConstMacro(UseImageSpacing, bool);
119  itkBooleanMacro(UseImageSpacing);
120 
122  void OverrideBoundaryCondition(ImageBoundaryCondition< TInputImage >* boundaryCondition);
123 
124 #ifdef ITK_USE_CONCEPT_CHECKING
125  // Begin concept checking
126  itkConceptMacro( InputConvertibleToOutputCheck,
128  itkConceptMacro( OutputHasNumericTraitsCheck,
130  // End concept checking
131 #endif
132 
143  itkSetMacro(UseImageDirection, bool);
144  itkGetConstMacro(UseImageDirection, bool);
145  itkBooleanMacro(UseImageDirection);
147 
148 protected:
150  ~GradientImageFilter() override;
151  void PrintSelf(std::ostream & os, Indent indent) const override;
152 
163  void DynamicThreadedGenerateData(const OutputImageRegionType & outputRegionForThread) override;
164 
165 
166 private:
167  void GenerateOutputInformation() override;
168 
169  // An overloaded method which may transform the gradient to a
170  // physical vector and converts to the correct output pixel type.
171  template <typename TValue>
173  {
174  if ( this->m_UseImageDirection )
175  {
176  CovariantVectorType physicalGradient;
177  it.GetImage()->TransformLocalVectorToPhysicalVector( gradient, physicalGradient );
178  it.Set( OutputPixelType( physicalGradient.GetDataPointer(), InputImageDimension, false ) );
179  }
180  else
181  {
182  it.Set( OutputPixelType( gradient.GetDataPointer(), InputImageDimension, false ) );
183  }
184  }
185 
186  template <typename T >
188  {
189  // This uses the more efficient set by reference method
190  if ( this->m_UseImageDirection )
191  {
192  it.GetImage()->TransformLocalVectorToPhysicalVector( gradient, it.Value() );
193  }
194  else
195  {
196  it.Value() = gradient;
197  }
198  }
199 
200 
202 
203  // flag to take or not the image direction into account
204  // when computing the derivatives.
206 
207  // allow setting the the m_BoundaryCondition
209 };
210 } // end namespace itk
211 
212 #ifndef ITK_MANUAL_INSTANTIATION
213 #include "itkGradientImageFilter.hxx"
214 #endif
215 
216 #endif
typename InputImageType::Pointer InputImagePointer
ImageBoundaryCondition< TInputImage, TInputImage > * m_BoundaryCondition
Light weight base class for most itk classes.
Templated n-dimensional vector image class.
typename OutputImageType::PixelType OutputPixelType
typename InputImageType::PixelType InputPixelType
void SetOutputPixel(ImageRegionIterator< VectorImage< TValue, OutputImageDimension > > &it, CovariantVectorType &gradient)
typename OutputImageType::RegionType OutputImageRegionType
void SetOutputPixel(ImageRegionIterator< T > &it, CovariantVectorType &gradient)
typename OutputImageType::Pointer OutputImagePointer
A virtual base object that defines an interface to a class of boundary condition objects for use by n...
TOperatorValueType OperatorValueType
Base class for filters that take an image as input and produce an image as output.
Control indentation during Print() invocation.
Definition: itkIndent.h:49
ValueType * GetDataPointer()
#define itkConceptMacro(name, concept)
A templated class holding a n-Dimensional covariant vector.
A multi-dimensional iterator templated over image type that walks a region of pixels.
const ImageType * GetImage() const
Computes the gradient of an image using directional derivatives.