[ VIGRA Homepage | Function Index | Class Index | Namespaces | File List | Main Page ]

details Common Filters VIGRA

Functions

template<... >
void convolveImage (...)
 Convolve an image with the given kernel(s). More...
 
template<... >
void convolveImageWithMask (...)
 Deprecated name of 2-dimensional normalized convolution, i.e. convolution with a mask image. More...
 
template<... >
void gaussianGradient (...)
 Calculate the gradient vector by means of a 1st derivatives of Gaussian filter. More...
 
template<... >
void gaussianGradientMagnitude (...)
 Calculate the gradient magnitude by means of a 1st derivatives of Gaussian filter. More...
 
template<... >
void gaussianSharpening (...)
 Perform sharpening function with gaussian filter. More...
 
template<... >
void gaussianSmoothing (...)
 Perform isotropic Gaussian convolution. More...
 
template<... >
void hessianMatrixOfGaussian (...)
 Filter image with the 2nd derivatives of the Gaussian at the given scale to get the Hessian matrix. More...
 
template<... >
void laplacianOfGaussian (...)
 Filter image with the Laplacian of Gaussian operator at the given scale. More...
 
template<... >
void normalizedConvolveImage (...)
 Performs a 2-dimensional normalized convolution, i.e. convolution with a mask image. More...
 
template<... >
void rieszTransformOfLOG (...)
 Calculate Riesz transforms of the Laplacian of Gaussian. More...
 
template<... >
void simpleSharpening (...)
 Perform simple sharpening function. More...
 
template<... >
void structureTensor (...)
 Calculate the Structure Tensor for each pixel of and image, using Gaussian (derivative) filters. More...
 

Detailed Description

These functions calculate common filters by appropriate sequences of calls to separableConvolveX() and separableConvolveY() or explicit 2-dimensional convolution.

Function Documentation

void vigra::rieszTransformOfLOG (   ...)

Calculate Riesz transforms of the Laplacian of Gaussian.

The Riesz transforms of the Laplacian of Gaussian have the following transfer functions (defined in a polar coordinate representation of the frequency domain):

\[ F_{\sigma}(r, \phi)=(i \cos \phi)^n (i \sin \phi)^m r^2 e^{-r^2 \sigma^2 / 2} \]

where n = xorder and m = yorder determine th e order of the transform, and sigma > 0 is the scale of the Laplacian of Gaussian. This function computes a good spatial domain approximation of these transforms for xorder + yorder <= 2. The filter responses may be used to calculate the monogenic signal or the boundary tensor.

Declarations:

pass 2D array views:

namespace vigra {
template <class T1, class S1,
class T2, class S2>
void
rieszTransformOfLOG(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, T2, S2> dest,
double scale, unsigned int xorder, unsigned int yorder);
}

show deprecated declarations

Usage:

#include <vigra/boundarytensor.hxx>
Namespace: vigra

MultiArrayView<2, double> impulse(17,17), res(17, 17);
impulse(8,8) = 1.0;
// calculate the impulse response of the first order Riesz transform in x-direction
rieszTransformOfLOG(impulse, res, 2.0, 1, 0);
void vigra::convolveImage (   ...)

Convolve an image with the given kernel(s).

If you pass vigra::Kernel2D to this function, it will perform an explicit 2-dimensional convolution. If you pass a single vigra::Kernel1D, it performs a separable convolution, i.e. it concatenates two 1D convolutions (along the x-axis and along the y-axis) with the same kernel via internal calls to separableConvolveX() and separableConvolveY(). If two 1D kernels are specified, separable convolution uses different kernels for the x- and y-axis.

All border treatment modes are supported.

The input pixel type T1 must be a linear space over the kernel's value_type T, i.e. addition of source values, multiplication with kernel values, and NumericTraits must be defined. The kernel's value_type must be an algebraic field, i.e. the arithmetic operations (+, -, *, /) and NumericTraits must be defined. Typically, you will use double for the kernel type.

Declarations:

pass 2D array views:

namespace vigra {
// use the same 1D kernel for all axes
template <class T1, class S1,
class T2, class S2,
class T>
void
convolveImage(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, T2, S2> dest,
Kernel1D<T> const & k);
// use a different kernel for each axis
template <class T1, class S1,
class T2, class S2,
class T>
void
convolveImage(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, T2, S2> dest,
Kernel1D<T> const & kx, Kernel1D<T> const & ky);
// use a non-separable 2D kernel
template <class T1, class S1,
class T2, class S2,
class T3>
void
convolveImage(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, T2, S2> dest,
Kernel2D<T3> const & kernel);
}

show deprecated declarations

Usage:

#include <vigra/convolution.hxx>
Namespace: vigra

MultiArray<2, float> src(w,h), dest1(w,h), dest2(w,h);
...
// create horizontal sobel filter (symmetric difference in x-direction, smoothing in y direction)
Kernel1D<double> kx, ky;
kx.initSymmetricDifference();
ky.initBinomial(1);
// calls separable convolution with the two 1D kernels
convolveImage(src, dest1, kx, ky);
// create a 3x3 Laplacian filter
Kernel2D<double> laplace;
laplace.initExplicitly(Diff2D(-1,-1), Diff2D(1,1)) =
0.375, 0.25, 0.375,
0.25, -2.5, 0.25,
0.375, 0.25, 0.375;
// calls 2D convolution
convolveImage(src, dest2, laplace);

show deprecated examples

Preconditions:

The image must be larger than the kernel radius.

  • For 1D kernels, w > std::max(xkernel.right(), -xkernel.keft()) and h > std::max(ykernel.right(), -ykernel.left()) are required.
  • For 2D kernels, w > std::max(kernel.lowerRight().x, -kernel.upperLeft().x) and h > std::max(kernel.lowerRight().y, -kernel.upperLeft().y) are required.

If BORDER_TREATMENT_CLIP is requested: the sum of kernel elements must be != 0.

Examples:
smooth_convolve.cxx.
void vigra::simpleSharpening (   ...)

Perform simple sharpening function.

This function uses convolveImage() with the following 3x3 filter:

-sharpening_factor/16.0, -sharpening_factor/8.0, -sharpening_factor/16.0,
-sharpening_factor/8.0, 1.0+sharpening_factor*0.75, -sharpening_factor/8.0,
-sharpening_factor/16.0, -sharpening_factor/8.0, -sharpening_factor/16.0;

and uses BORDER_TREATMENT_REFLECT as border treatment mode.

Preconditions:

1. sharpening_factor >= 0
2. scale >= 0

Declarations:

pass 2D array views:

namespace vigra {
template <class T1, class S1,
class T2, class S2>
void
simpleSharpening(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, T2, S2> dest,
double sharpening_factor);
}

show deprecated declarations

Usage:

#include <vigra/convolution.hxx>
Namespace: vigra

MultiArray<2, float> src(w,h), dest(w,h);
...
// sharpening with sharpening_factor = 0.1
vigra::simpleSharpening(src, dest, 0.1);

show deprecated examples

void vigra::gaussianSharpening (   ...)

Perform sharpening function with gaussian filter.

This function uses gaussianSmoothing() at the given scale to create a temporary image 'smooth' and than blends the original and smoothed image according to the formula

dest = (1 + sharpening_factor)*src - sharpening_factor*smooth

Preconditions:

1. sharpening_factor >= 0
2. scale >= 0

Declarations:

pass 2D array views:

namespace vigra {
template <class T1, class S1,
class T2, class S2>
void
gaussianSharpening(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, T2, S2> dest,
double sharpening_factor,
double scale);
}

show deprecated declarations

Usage:

#include <vigra/convolution.hxx>
Namespace: vigra

MultiArray<2, float> src(w,h), dest(w,h);
...
// sharpening with sharpening_factor = 3.0
// smoothing with scale = 0.5
gaussianSharpening(src, dest, 3.0, 0.5);

show deprecated examples

void vigra::gaussianSmoothing (   ...)

Perform isotropic Gaussian convolution.

This function is a shorthand for the concatenation of a call to separableConvolveX() and separableConvolveY() with a Gaussian kernel of the given scale. If two scales are provided, smoothing in x and y direction will have different strength. The function uses BORDER_TREATMENT_REFLECT.

Function gaussianSmoothMultiArray() performs the same filter operation on arbitrary dimensional arrays.

Declarations:

pass 2D array views:

namespace vigra {
template <class T1, class S1,
class T2, class S2>
void
gaussianSmoothing(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, T2, S2> dest,
double scale_x, double scale_y = scale_x);
}

show deprecated declarations

Usage:

#include <vigra/convolution.hxx>
Namespace: vigra

MultiArray<2, float> src(w,h), dest(w,h);
...
// smooth with scale = 3.0
gaussianSmoothing(src, dest, 3.0);

show deprecated examples

Examples:
smooth.cxx.
void vigra::gaussianGradient (   ...)

Calculate the gradient vector by means of a 1st derivatives of Gaussian filter.

This function is a shorthand for the concatenation of a call to separableConvolveX() and separableConvolveY() with the appropriate kernels at the given scale. Note that this function can either produce two separate result images for the x- and y-components of the gradient, or write into a vector valued image (with at least two components).

Function gaussianGradientMultiArray() performs the same filter operation on arbitrary dimensional arrays.

Declarations:

pass 2D array views:

namespace vigra {
// write x and y component of the gradient into separate images
template <class T1, class S1,
class T2X, class S2X,
class T2Y, class S2Y>
void
gaussianGradient(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, T2X, S2X> destx,
MultiArrayView<2, T2Y, S2Y> desty,
double scale);
// write x and y component of the gradient into a vector-valued image
template <class T1, class S1,
class T2, class S2>
void
gaussianGradient(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, TinyVector<T2, 2>, S2> dest,
double scale);
}

show deprecated declarations

Usage:

#include <vigra/convolution.hxx>
Namespace: vigra

MultiArray<2, float> src(w,h), gradx(w,h), grady(w,h);
...
// calculate gradient vector at scale = 3.0
gaussianGradient(src, gradx, grady, 3.0);
// likewise, but use a vector image to store the gradient
MultiArray<2, TinyVector<float, 2> > dest(w,h);
gaussianGradient(src, dest, 3.0);

show deprecated examples

void vigra::gaussianGradientMagnitude (   ...)

Calculate the gradient magnitude by means of a 1st derivatives of Gaussian filter.

This function calls gaussianGradient() and returns the pixel-wise magnitude of the resulting gradient vectors. If the original image has multiple bands, the squared gradient magnitude is computed for each band separately, and the return value is the square root of the sum of these squared magnitudes.

Anisotropic data should be provided with appropriate vigra::ConvolutionOptions to adjust the filter sizes for the resolution of each axis. Otherwise, the parameter opt is optional unless the parameter sigma is omitted.

If you pass vigra::BlockwiseConvolutionOptions instead, the algorithm will be executed in parallel on data blocks of a certain size. The block size can be customized via BlockwiseConvolutionOptions::blockShape(), but the defaults usually work reasonably. By default, the number of threads equals the capabilities of your hardware, but you can change this via BlockwiseConvolutionOptions::numThreads().

Declarations:

use arbitrary-dimensional arrays:

namespace vigra {
// pass filter scale explicitly
template <unsigned int N, class T1, class S1,
class T2, class S2>
void
gaussianGradientMagnitude(MultiArrayView<N, T1, S1> const & src,
MultiArrayView<N, T2, S2> dest,
double sigma,
ConvolutionOptions<N> opt = ConvolutionOptions<N>());
// likewise, but sum the contributions of each band
template <unsigned int N, class MT, class S1,
class T2, class S2>
void
gaussianGradientMagnitude(MultiArrayView<N+1, Multiband<MT>, S1> const & src,
MultiArrayView<N, T2, S2> dest,
double sigma,
ConvolutionOptions<N> opt = ConvolutionOptions<N>());
// pass filter scale(s) in option object
template <unsigned int N, class T1, class S1,
class T2, class S2>
void
gaussianGradientMagnitude(MultiArrayView<N, T1, S1> const & src,
MultiArrayView<N, T2, S2> dest,
ConvolutionOptions<N> const & opt);
// likewise, but execute algorithm in parallel
template <unsigned int N, class T1, class S1,
class T2, class S2>
void
gaussianGradientMagnitude(MultiArrayView<N, T1, S1> const & src,
MultiArrayView<N, T2, S2> dest,
BlockwiseConvolutionOptions<N> const & opt);
// pass filter scale(s) in option object and sum the contributions of each band
template <unsigned int N, class MT, class S1,
class T2, class S2>
void
gaussianGradientMagnitude(MultiArrayView<N+1, Multiband<MT>, S1> const & src,
MultiArrayView<N, T2, S2> dest,
ConvolutionOptions<N> const & opt);
}

Here, the input element types T1 and MT can be arbitrary scalar types, and T1 may also be TinyVector or RGBValue. The output element type T2 should be the corresponding norm type (see NormTraits). In the Multiband<MT>-version, the input array's right-most dimension is interpreted as a channel axis, therefore it must have one dimension more than the output array.

show deprecated declarations

Usage:

#include <vigra/multi_convolution.hxx> (sequential version)
#include <vigra/multi_blockwise.hxx> (parallel version)
#include <vigra/convolution.hxx> (deprecated API version)
Namespace: vigra

// example 1
{
// use a 3-dimensional float array
MultiArray<3, float> volume(Shape3(w, h, d)), grad(volume.shape());
...
// calculate gradient magnitude at scale = 3.0
gaussianGradientMagnitude(volume, grad, 3.0);
}
// example 2
{
// use a 2-dimensional RGB array
MultiArray<2, RGBValue<float> > rgb(Shape2(w, h));
MultiArray<2, float> grad(rgb.shape());
...
// calculate the color gradient magnitude at scale = 3.0
gaussianGradientMagnitude(rgb, grad, 3.0);
}
// example 3
{
// use a 3-dimensional array whose right-most axis is interpreted as
// a multi-spectral axis with arbitrary many channels
MultiArray<3, Multiband<float> > spectral(Shape3(w, h, channelCount));
MultiArray<2, float> grad(Shape2(w, h));
...
// calculate the multi-channel gradient magnitude at scale = 3.0
// (note that the template parameter N (number of spatial dimensions)
// must be provided explicitly as gaussianGradientMagnitude<2>(...) )
MultiArrayView<3, Multiband<float> > view(spectral);
gaussianGradientMagnitude<2>(view, grad, 3.0);
}

show deprecated examples

Examples:
watershed.cxx.
void vigra::laplacianOfGaussian (   ...)

Filter image with the Laplacian of Gaussian operator at the given scale.

This function calls separableConvolveX() and separableConvolveY() with the appropriate 2nd derivative of Gaussian kernels in x- and y-direction and then sums the results to get the Laplacian.

Function laplacianOfGaussianMultiArray() performs the same filter operation on arbitrary dimensional arrays.

Declarations:

pass 2D array views:

namespace vigra {
template <class T1, class S1,
class T2, class S2>
void
laplacianOfGaussian(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, T2, S2> dest,
double scale);
}

show deprecated declarations

Usage:

#include <vigra/convolution.hxx>
Namespace: vigra

MultiArray<2, float> src(w,h), dest(w,h);
...
// calculate Laplacian of Gaussian at scale = 3.0
laplacianOfGaussian(src, dest, 3.0);

show deprecated examples

void vigra::hessianMatrixOfGaussian (   ...)

Filter image with the 2nd derivatives of the Gaussian at the given scale to get the Hessian matrix.

The Hessian matrix is a symmetric matrix defined as:

\[ \mbox{\rm Hessian}(I) = \left( \begin{array}{cc} G_{xx} \ast I & G_{xy} \ast I \\ G_{xy} \ast I & G_{yy} \ast I \end{array} \right) \]

where $G_{xx}, G_{xy}, G_{yy}$ denote 2nd derivatives of Gaussians at the given scale, and $\ast$ is the convolution symbol. This function calls separableConvolveX() and separableConvolveY() with the appropriate 2nd derivative of Gaussian kernels and puts the results in the three destination images. The first destination image will contain the second derivative in x-direction, the second one the mixed derivative, and the third one holds the derivative in y-direction.

Function hessianOfGaussianMultiArray() performs the same filter operation on arbitrary dimensional arrays.

Declarations:

pass 2D array views:

namespace vigra {
template <class T1, class S1,
class T2, class S2>
void
hessianMatrixOfGaussian(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, TinyVector<T2, 3>, S2> dest,
double scale);
}

show deprecated declarations

Usage:

#include <vigra/convolution.hxx>
Namespace: vigra

MultiArray<2, float> src(w,h);
MultiArray<2, TinyVector<float, 3> > hessian(w,h); // will hold the three components of the Hessian
...
// calculate Hessian of Gaussian at scale = 3.0, use a 3-band output image
hessianMatrixOfGaussian(src, hessian, 3.0);

show deprecated examples

void vigra::structureTensor (   ...)

Calculate the Structure Tensor for each pixel of and image, using Gaussian (derivative) filters.

The Structure Tensor is is a smoothed version of the Euclidean product of the gradient vector with itself. I.e. it's a symmetric matrix defined as:

\[ \mbox{\rm StructurTensor}(I) = \left( \begin{array}{cc} G \ast (I_x I_x) & G \ast (I_x I_y) \\ G \ast (I_x I_y) & G \ast (I_y I_y) \end{array} \right) = \left( \begin{array}{cc} A & C \\ C & B \end{array} \right) \]

where $G$ denotes Gaussian smoothing at the outer scale, $I_x, I_y$ are the gradient components taken at the inner scale, $\ast$ is the convolution symbol, and $I_x I_x$ etc. are pixelwise products of the 1st derivative images. This function calls separableConvolveX() and separableConvolveY() with the appropriate Gaussian kernels and puts the results in the three separate destination images (where the first one will contain $G \ast (I_x I_x)$, the second one $G \ast (I_x I_y)$, and the third one holds $G \ast (I_y I_y)$), or into a single 3-band image (where the bands hold the result in the same order as above). The latter form is also applicable when the source image is a multi-band image (e.g. RGB). In this case, tensors are first computed for each band separately, and then summed up to get a single result tensor.

Function structureTensorMultiArray() performs the same filter operation on arbitrary dimensional arrays.

Declarations:

pass 2D array views:

namespace vigra {
// create three separate destination images
template <class T, class S,
class TX, class SX,
class TXY, class SXY,
class TY, class SY>
void
structureTensor(MultiArrayView<2, S, T> const & src,
MultiArrayView<2, TX, SX> destx,
MultiArrayView<2, TXY, SXY> destxy,
MultiArrayView<2, TY, SY> desty,
double inner_scale, double outer_scale);
// create a single 3-band destination image
template <class T1, class S1,
class T2, class S2>
void
structureTensor(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, TinyVector<T2, 3>, S2> dest,
double inner_scale, double outer_scale);
}

show deprecated declarations

Usage:

#include <vigra/convolution.hxx>
Namespace: vigra

MultiArray<2, flost> src(w,h),
stxx(w,h), stxy(w,h), styy(w,h); // use a separate image for each component
...
// calculate Structure Tensor at inner scale = 1.0 and outer scale = 3.0
structureTensor(src, stxx, stxy, styy, 1.0, 3.0);
// likwise with a single 3-band destination image
MultiArray<2, TinyVector<float, 3> > st(w,h);
structureTensor(src, st, 1.0, 3.0);

show deprecated examples

void vigra::normalizedConvolveImage (   ...)

Performs a 2-dimensional normalized convolution, i.e. convolution with a mask image.

This functions computes normalized convolution as defined in Knutsson, H. and Westin, C-F.: Normalized and differential convolution: Methods for Interpolation and Filtering of incomplete and uncertain data. Proc. of the IEEE Conf. on Computer Vision and Pattern Recognition, 1993, 515-523.

The mask image must be binary and encodes which pixels of the original image are valid. It is used as follows: Only pixel under the mask are used in the calculations. Whenever a part of the kernel lies outside the mask, it is ignored, and the kernel is renormalized to its original norm (analogous to the CLIP BorderTreatmentMode). Thus, a useful convolution result is computed whenever at least one valid pixel is within the current window Thus, destination pixels not under the mask still receive a value if they are near the mask. Therefore, this algorithm is useful as an interpolator of sparse input data. If you are only interested in the destination values under the mask, you can perform a subsequent copyImageIf().

The KernelIterator must point to the center of the kernel, and the kernel's size is given by its upper left (x and y of distance <= 0) and lower right (distance >= 0) corners. The image must always be larger than the kernel. At those positions where the kernel does not completely fit into the image, the specified BorderTreatmentMode is applied. Only BORDER_TREATMENT_CLIP and BORDER_TREATMENT_AVOID are currently supported.

The images's pixel type (SrcAccessor::value_type) must be a linear space over the kernel's value_type (KernelAccessor::value_type), i.e. addition of source values, multiplication with kernel values, and NumericTraits must be defined. The kernel's value_type must be an algebraic field, i.e. the arithmetic operations (+, -, *, /) and NumericTraits must be defined.

Declarations:

pass 2D array views:

namespace vigra {
template <class T1, class S1,
class T2, class S2,
class TM, class SM,
class T3>
void
normalizedConvolveImage(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, TM, SM> const & mask,
MultiArrayView<2, T2, S2> dest,
Kernel2D<T3> const & kernel);
}

show deprecated declarations

Usage:

#include <vigra/stdconvolution.hxx>
Namespace: vigra

MultiArray<2, float> src(w,h), dest(w,h);
MultiArray<2, unsigned char> mask(w,h);
...
// define 3x3 binomial filter
vigra::Kernel2D<float> binom;
binom.initExplicitly(Diff2D(-1,-1), Diff2D(1,1)) = // upper left and lower right
0.0625, 0.125, 0.0625,
0.125, 0.25, 0.125,
0.0625, 0.125, 0.0625;
normalizedConvolveImage(src, mask, dest, binom);

show deprecated examples

Preconditions:

  • The image must be longer than the kernel radius: w > std::max(kernel.lowerRight().x, -kernel.upperLeft().x) and h > std::max(kernel.lowerRight().y, -kernel.upperLeft().y).
  • The sum of kernel elements must be != 0.
  • border == BORDER_TREATMENT_CLIP || border == BORDER_TREATMENT_AVOID
void vigra::convolveImageWithMask (   ...)

Deprecated name of 2-dimensional normalized convolution, i.e. convolution with a mask image.

See normalizedConvolveImage() for documentation.

Declarations:

pass 2D array views:

namespace vigra {
template <class SrcIterator, class SrcAccessor,
class MaskIterator, class MaskAccessor,
class DestIterator, class DestAccessor,
class KernelIterator, class KernelAccessor>
void
convolveImageWithMask(SrcIterator src_ul, SrcIterator src_lr, SrcAccessor src_acc,
MaskIterator mul, MaskAccessor am,
DestIterator dest_ul, DestAccessor dest_acc,
KernelIterator ki, KernelAccessor ak,
Diff2D kul, Diff2D klr, BorderTreatmentMode border);
}

show deprecated declarations

© Ullrich Köthe (ullrich.koethe@iwr.uni-heidelberg.de)
Heidelberg Collaboratory for Image Processing, University of Heidelberg, Germany

html generated using doxygen and Python
vigra 1.11.0 (Thu Feb 4 2016)