ITK  6.0.0
Insight Toolkit
itkGradientRecursiveGaussianImageFilter.h
Go to the documentation of this file.
1/*=========================================================================
2 *
3 * Copyright NumFOCUS
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#ifndef itkGradientRecursiveGaussianImageFilter_h
19#define itkGradientRecursiveGaussianImageFilter_h
20
23#include "itkImage.h"
24#include "itkCovariantVector.h"
28#include "itkVectorImage.h"
29#include <vector>
30
31namespace itk
32{
54template <
55 typename TInputImage,
56 typename TOutputImage = Image<
57 CovariantVector<typename NumericTraits<typename TInputImage::PixelType>::RealType, TInputImage::ImageDimension>,
58 TInputImage::ImageDimension>>
59class ITK_TEMPLATE_EXPORT GradientRecursiveGaussianImageFilter : public ImageToImageFilter<TInputImage, TOutputImage>
60{
61public:
62 ITK_DISALLOW_COPY_AND_MOVE(GradientRecursiveGaussianImageFilter);
63
69
71 using InputImageType = TInputImage;
72 using PixelType = typename TInputImage::PixelType;
75
81
83 static constexpr unsigned int ImageDimension = TInputImage::ImageDimension;
84
87
92
93
101
104
107
110
113
116
119
121 using OutputImageType = TOutputImage;
122 using OutputPixelType = typename OutputImageType::PixelType;
125
127 itkNewMacro(Self);
128
130 itkOverrideGetNameOfClassMacro(GradientRecursiveGaussianImageFilter);
131
133 void
134 SetSigmaArray(const SigmaArrayType & sigma);
135 void
136 SetSigma(ScalarRealType sigma);
140 GetSigmaArray() const;
141
144 GetSigma() const;
145
149 void
150 SetNormalizeAcrossScale(bool normalize);
151 itkGetConstMacro(NormalizeAcrossScale, bool);
159 void
160 GenerateInputRequestedRegion() override;
161
172 itkSetMacro(UseImageDirection, bool);
173 itkGetConstMacro(UseImageDirection, bool);
174 itkBooleanMacro(UseImageDirection);
177#ifdef ITK_USE_CONCEPT_CHECKING
178 // Begin concept checking
179 // Does not seem to work with wrappings, disabled
180 // itkConceptMacro( InputHasNumericTraitsCheck,
181 // ( Concept::HasNumericTraits< PixelType > ) );
182 // End concept checking
183#endif
184
185protected:
187 ~GradientRecursiveGaussianImageFilter() override = default;
188 void
189 PrintSelf(std::ostream & os, Indent indent) const override;
190
192 void
193 GenerateData() override;
194
195 // Override since the filter produces the entire dataset
196 void
197 EnlargeOutputRequestedRegion(DataObject * output) override;
198
199 void
200 GenerateOutputInformation() override;
201
202private:
203 template <typename TValue>
204 void
206 {
207 // To transform Variable length vector we need to convert to and
208 // fro the CovariantVectorType
209 const CovariantVectorType gradient(it.Get().GetDataPointer());
210 CovariantVectorType physicalGradient;
211 it.GetImage()->TransformLocalVectorToPhysicalVector(gradient, physicalGradient);
212 it.Set(OutputPixelType(physicalGradient.GetDataPointer(), ImageDimension, false));
213 }
214
215 template <typename T>
216 void
217 TransformOutputPixel(ImageRegionIterator<T> & it)
218 {
219 OutputPixelType correctedGradient{};
220 const OutputPixelType & gradient = it.Get();
221
222 const unsigned int nComponents = NumericTraits<OutputPixelType>::GetLength(gradient) / ImageDimension;
223
224 for (unsigned int nc = 0; nc < nComponents; ++nc)
225 {
226 GradientVectorType componentGradient;
227 GradientVectorType correctedComponentGradient;
228 for (unsigned int dim = 0; dim < ImageDimension; ++dim)
229 {
230 componentGradient[dim] =
231 DefaultConvertPixelTraits<OutputPixelType>::GetNthComponent(nc * ImageDimension + dim, gradient);
232 }
233 it.GetImage()->TransformLocalVectorToPhysicalVector(componentGradient, correctedComponentGradient);
234 for (unsigned int dim = 0; dim < ImageDimension; ++dim)
235 {
237 nc * ImageDimension + dim, correctedGradient, correctedComponentGradient[dim]);
238 }
239 }
240 it.Set(correctedGradient);
241 }
242
243 template <template <typename, unsigned int> class P, class T, unsigned int VDimension>
244 void
245 TransformOutputPixel(ImageRegionIterator<Image<P<T, VDimension>, VDimension>> & it)
246 {
247 const OutputPixelType gradient = it.Get();
248 // This uses the more efficient set by reference method
249 it.GetImage()->TransformLocalVectorToPhysicalVector(gradient, it.Value());
250 }
251
253 std::vector<GaussianFilterPointer> m_SmoothingFilters{};
254 DerivativeFilterPointer m_DerivativeFilter{};
255 OutputImageAdaptorPointer m_ImageAdaptor{};
256
258 bool m_NormalizeAcrossScale{};
259
261 bool m_UseImageDirection{ true };
262
264 SigmaArrayType m_Sigma{};
265};
266} // end namespace itk
267
268#ifndef ITK_MANUAL_INSTANTIATION
269# include "itkGradientRecursiveGaussianImageFilter.hxx"
270#endif
271
272#endif
A templated class holding a n-Dimensional covariant vector.
Base class for all data objects in ITK.
static ComponentType GetNthComponent(int c, const PixelType &pixel)
static void SetNthComponent(int c, PixelType &pixel, const ComponentType &v)
ValueType * GetDataPointer()
Computes the gradient of an image by convolution with the first derivative of a Gaussian.
typename NumericTraits< PixelType >::RealType RealType
typename NumericTraits< OutputPixelType >::ValueType OutputComponentType
typename NumericTraits< InternalRealType >::ValueType InternalScalarRealType
typename NumericTraits< RealType >::FloatType InternalRealType
typename OutputImageAdaptorType::Pointer OutputImageAdaptorPointer
typename NumericTraits< PixelType >::ScalarRealType ScalarRealType
typename DerivativeFilterType::Pointer DerivativeFilterPointer
const ImageType * GetImage() const
A multi-dimensional iterator templated over image type that walks a region of pixels.
void Set(const PixelType &value) const
Base class for filters that take an image as input and produce an image as output.
Templated n-dimensional image class.
Definition: itkImage.h:89
Control indentation during Print() invocation.
Definition: itkIndent.h:50
Light weight base class for most itk classes.
Presents an image as being composed of the N-th element of its pixels.
Define additional traits for native types such as int or float.
static unsigned int GetLength()
Base class for computing IIR convolution with an approximation of a Gaussian kernel.
Templated n-dimensional vector image class.
SmartPointer< Self > Pointer
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....