Create Gaussian Derivative Kernel#

Synopsis#

Create a Gaussian derivative kernel.

Results#

Output:

Size: [3, 3]
Neighborhood:
Radius:[1, 1]
Size:[3, 3]
DataBuffer:NeighborhoodAllocator { this = 0x7ffee86419a8, begin = 0x7fb16bc5cc50, size=9 }
[-1, -1] 0
[0, -1] 0
[1, -1] 0
[-1, 0] 0.208375
[0, 0] 0
[1, 0] -0.208375
[-1, 1] 0
[0, 1] 0
[1, 1] 0

Code#

Python#

#!/usr/bin/env python

import itk

gaussianDerivativeOperator = itk.GaussianDerivativeOperator[itk.F, 2]()
gaussianDerivativeOperator.SetDirection(
    0
)  # Create the operator for the X axis derivative
radius = itk.Size[2]()
radius.Fill(1)
gaussianDerivativeOperator.CreateToRadius(radius)

print("Size: " + str(gaussianDerivativeOperator.GetSize()))

print(gaussianDerivativeOperator)

for i in range(9):
    print(
        str(gaussianDerivativeOperator.GetOffset(i))
        + " "
        + str(gaussianDerivativeOperator.GetElement(i))
    )

C++#

#include <itkGaussianDerivativeOperator.h>

int
main()
{
  using GaussianDerivativeOperatorType = itk::GaussianDerivativeOperator<float, 2>;
  GaussianDerivativeOperatorType gaussianDerivativeOperator;
  gaussianDerivativeOperator.SetDirection(0); // Create the operator for the X axis derivative
  itk::Size<2> radius;
  radius.Fill(1);
  gaussianDerivativeOperator.CreateToRadius(radius);

  std::cout << "Size: " << gaussianDerivativeOperator.GetSize() << std::endl;

  std::cout << gaussianDerivativeOperator << std::endl;

  for (unsigned int i = 0; i < 9; ++i)
  {
    std::cout << gaussianDerivativeOperator.GetOffset(i) << " " << gaussianDerivativeOperator.GetElement(i)
              << std::endl;
  }
  return EXIT_SUCCESS;
}

Classes demonstrated#

template<typename TPixel, unsigned int VDimension = 2, typename TAllocator = NeighborhoodAllocator<TPixel>>
class GaussianDerivativeOperator : public itk::NeighborhoodOperator<TPixel, VDimension, TAllocator>

A NeighborhoodOperator whose coefficients are a one dimensional, discrete derivative Gaussian kernel.

GaussianDerivativeOperator can be used to calculate Gaussian derivatives by taking its inner product with to a Neighborhood (NeighborhooIterator) that is swept across an image region. It is a directional operator. N successive applications oriented along each dimensional direction will calculate separable, efficient, N-D Gaussian derivatives of an image region.

GaussianDerivativeOperator takes three parameters:

(1) The floating-point variance of the desired Gaussian function.

(2) The order of the derivative to be calculated (zero order means it performs only smoothing as a standard itk::GaussianOperator)

(3) The “maximum error” allowed in the discrete Gaussian function. “Maximum errror” is defined as the difference between the area under the discrete Gaussian curve and the area under the continuous Gaussian. Maximum error affects the Gaussian operator size. Care should be taken not to make this value too small relative to the variance lest the operator size become unreasonably large.

References: The Gaussian kernel contained in this operator was described by Tony Lindeberg (Discrete Scale-Space Theory and the Scale-Space Primal Sketch. Dissertation. Royal Institute of Technology, Stockholm, Sweden. May 1991.).

This implementation is derived from the Insight Journal paper:

https://www.insight-journal.org/browse/publication/179
Author

Ivan Macia, Vicomtech, Spain, https://www.vicomtech.org/en

Note

GaussianDerivativeOperator does not have any user-declared “special member function”, following the C++ Rule of Zero: the compiler will generate them if necessary.

See

GaussianOperator

See

NeighborhoodOperator

See

NeighborhoodIterator

See

Neighborhood

ITK Sphinx Examples:

See itk::GaussianDerivativeOperator for additional documentation.