ITK  6.0.0
Insight Toolkit
itkDiscreteGaussianImageFilter.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 itkDiscreteGaussianImageFilter_h
19#define itkDiscreteGaussianImageFilter_h
20
21#include "itkGaussianOperator.h"
23#include "itkImage.h"
25
26namespace itk
27{
63template <typename TInputImage, typename TOutputImage = TInputImage>
64class ITK_TEMPLATE_EXPORT DiscreteGaussianImageFilter : public ImageToImageFilter<TInputImage, TOutputImage>
65{
66public:
67 ITK_DISALLOW_COPY_AND_MOVE(DiscreteGaussianImageFilter);
68
74
76 itkNewMacro(Self);
77
79 itkOverrideGetNameOfClassMacro(DiscreteGaussianImageFilter);
80
82 using InputImageType = TInputImage;
83 using OutputImageType = TOutputImage;
84
87 using OutputPixelType = typename TOutputImage::PixelType;
88 using OutputInternalPixelType = typename TOutputImage::InternalPixelType;
89 using InputPixelType = typename TInputImage::PixelType;
90 using InputInternalPixelType = typename TInputImage::InternalPixelType;
91
95
98 static constexpr unsigned int ImageDimension = TOutputImage::ImageDimension;
99
104
107#ifndef ITK_FUTURE_LEGACY_REMOVE
108 using InputBoundaryConditionPointerType [[deprecated("Please just use `BoundaryConditionType *` instead!")]] =
110#endif
114
118 using ScalarRealType = double;
119
123
130 itkSetMacro(Variance, ArrayType);
131 itkGetConstMacro(Variance, const ArrayType);
137 itkSetMacro(MaximumError, ArrayType);
138 itkGetConstMacro(MaximumError, const ArrayType);
143 itkGetConstMacro(MaximumKernelWidth, unsigned int);
144 itkSetMacro(MaximumKernelWidth, unsigned int);
152 itkGetConstMacro(FilterDimensionality, unsigned int);
153 itkSetMacro(FilterDimensionality, unsigned int);
157 itkSetMacro(InputBoundaryCondition, BoundaryConditionType *);
158 itkGetConstMacro(InputBoundaryCondition, BoundaryConditionType *);
159 itkSetMacro(RealBoundaryCondition, RealBoundaryConditionPointerType);
160 itkGetConstMacro(RealBoundaryCondition, RealBoundaryConditionPointerType);
165 void
167 {
168 m_Variance.Fill(v);
169 this->Modified();
170 }
173 void
175 {
176 m_MaximumError.Fill(v);
177 this->Modified();
178 }
179
180 void
181 SetVariance(const double * v)
182 {
183 ArrayType dv;
184
185 for (unsigned int i = 0; i < ImageDimension; ++i)
186 {
187 dv[i] = v[i];
188 }
189 this->SetVariance(dv);
190 }
191
192 void
193 SetVariance(const float * v)
194 {
195 ArrayType dv;
196
197 for (unsigned int i = 0; i < ImageDimension; ++i)
198 {
199 dv[i] = v[i];
200 }
201 this->SetVariance(dv);
202 }
203
206 void
207 SetSigma(const ArrayType & sigma)
208 {
209 ArrayType variance;
210 for (unsigned int i = 0; i < ImageDimension; ++i)
211 {
212 variance[i] = sigma[i] * sigma[i];
213 }
214 this->SetVariance(variance);
215 }
216
219 void
220 SetSigmaArray(const ArrayType & sigmas)
221 {
222 this->SetSigma(sigmas);
223 }
224 void
225 SetSigma(double sigma)
226 {
227 this->SetVariance(sigma * sigma);
228 }
232 ArrayType
234 {
235 ArrayType sigmas;
236 for (unsigned int i = 0; i < ImageDimension; ++i)
237 {
238 sigmas[i] = std::sqrt(m_Variance[i]);
239 }
240 return sigmas;
241 }
246 double
247 GetSigma() const
248 {
249 return std::sqrt(m_Variance[0]);
250 }
251
252 void
253 SetMaximumError(const double * v)
254 {
255 ArrayType dv;
256
257 for (unsigned int i = 0; i < ImageDimension; ++i)
258 {
259 dv[i] = v[i];
260 }
261 this->SetMaximumError(dv);
262 }
263
264 void
265 SetMaximumError(const float * v)
266 {
267 ArrayType dv;
268
269 for (unsigned int i = 0; i < ImageDimension; ++i)
270 {
271 dv[i] = v[i];
272 }
273 this->SetMaximumError(dv);
274 }
275
277 unsigned int
278 GetKernelRadius(const unsigned int dimension) const;
279
283
288
294 itkSetMacro(UseImageSpacing, bool);
295 itkGetConstMacro(UseImageSpacing, bool);
296 itkBooleanMacro(UseImageSpacing);
299#if !defined(ITK_FUTURE_LEGACY_REMOVE)
304 void
305 SetUseImageSpacingOn()
306 {
307 this->SetUseImageSpacing(true);
308 }
309
313 void
314 SetUseImageSpacingOff()
315 {
316 this->SetUseImageSpacing(false);
317 }
318#endif
319
329 itkLegacyMacro(unsigned int GetInternalNumberOfStreamDivisions() const;)
330 itkLegacyMacro(void SetInternalNumberOfStreamDivisions(unsigned int);)
331
332#ifdef ITK_USE_CONCEPT_CHECKING
333 // Begin concept checking
334
336
337 // End concept checking
338#endif
339
340protected:
342 {
343 m_Variance.Fill(0.0);
344 m_MaximumError.Fill(0.01);
345 m_MaximumKernelWidth = 32;
346 m_UseImageSpacing = true;
347 m_FilterDimensionality = ImageDimension;
348 m_InputBoundaryCondition = &m_InputDefaultBoundaryCondition;
349 m_RealBoundaryCondition = &m_RealDefaultBoundaryCondition;
350 }
351
352 ~DiscreteGaussianImageFilter() override = default;
353 void
354 PrintSelf(std::ostream & os, Indent indent) const override;
355
362 void
364
370 void
371 GenerateData() override;
372
374 void
375 GenerateKernel(const unsigned int dimension, KernelType & oper) const;
376
380
381private:
384 ArrayType m_Variance{};
385
389 ArrayType m_MaximumError{};
390
393 unsigned int m_MaximumKernelWidth{};
394
396 unsigned int m_FilterDimensionality{};
397
399 bool m_UseImageSpacing{};
400
403 BoundaryConditionType * m_InputBoundaryCondition{};
404
406 InputDefaultBoundaryConditionType m_InputDefaultBoundaryCondition{};
407
409 RealBoundaryConditionPointerType m_RealBoundaryCondition{};
410
412 RealDefaultBoundaryConditionType m_RealDefaultBoundaryCondition{};
413};
414} // end namespace itk
415
416#ifndef ITK_MANUAL_INSTANTIATION
417# include "itkDiscreteGaussianImageFilter.hxx"
418#endif
419
420#endif
Blurs an image by separable convolution with discrete gaussian kernels. This filter performs Gaussian...
void GenerateKernel(const unsigned int dimension, KernelType &oper) const
void GenerateInputRequestedRegion() override
typename NumericTraits< InputPixelType >::ValueType InputPixelValueType
typename TOutputImage::InternalPixelType OutputInternalPixelType
unsigned int GetKernelRadius(const unsigned int dimension) const
void SetVariance(const typename ArrayType::ValueType v)
typename TInputImage::InternalPixelType InputInternalPixelType
typename KernelType::RadiusType RadiusType
typename NumericTraits< RealOutputPixelType >::ValueType RealOutputPixelValueType
typename TInputImage::PixelType InputPixelType
void SetMaximumError(const typename ArrayType::ValueType v)
ArrayType GetKernelRadius() const
ArrayType GetKernelVarianceArray() const
~DiscreteGaussianImageFilter() override=default
typename NumericTraits< OutputPixelType >::ValueType OutputPixelValueType
typename TOutputImage::PixelType OutputPixelType
void PrintSelf(std::ostream &os, Indent indent) const override
itkLegacyMacro(unsigned int GetInternalNumberOfStreamDivisions() const ;) itkLegacyMacro(void SetInternalNumberOfStreamDivisions(unsigned int)
Set/Get number of pieces to divide the input for the internal composite pipeline. The upstream pipeli...
typename NumericTraits< OutputPixelType >::RealType RealOutputPixelType
A NeighborhoodOperator whose coefficients are a one dimensional, discrete Gaussian kernel.
A virtual base object that defines an interface to a class of boundary condition objects for use by n...
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.
Define additional traits for native types such as int or float.
A function object that determines a neighborhood of values at an image boundary according to a Neuman...
#define itkConceptMacro(name, concept)
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....