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

details Point operators for multi-dimensional arrays. VIGRA

Functions

template<... >
void combineThreeMultiArrays (...)
 Combine three multi-dimensional arrays into one using a ternary function or functor. More...
 
template<... >
void combineTwoMultiArrays (...)
 Combine two multi-dimensional arrays into one using a binary function or functor. More...
 
template<... >
void copyMultiArray (...)
 Copy a multi-dimensional array. More...
 
template<... >
void initMultiArray (...)
 Write a value to every element in a multi-dimensional array. More...
 
template<... >
void initMultiArrayBorder (...)
 Write values to the specified border values in the array. More...
 
template<... >
void inspectMultiArray (...)
 Call an analyzing functor at every element of a multi-dimensional array. More...
 
template<... >
void inspectTwoMultiArrays (...)
 Call an analyzing functor at all corresponding elements of two multi-dimensional arrays. More...
 
template<... >
void tensorDeterminantMultiArray (...)
 Calculate the tensor determinant for every element of a ND tensor array. More...
 
template<... >
void tensorEigenvaluesMultiArray (...)
 Calculate the tensor eigenvalues for every element of a N-D tensor array. More...
 
template<... >
void tensorTraceMultiArray (...)
 Calculate the tensor trace for every element of a N-D tensor array. More...
 
template<... >
void transformMultiArray (...)
 Transform a multi-dimensional array with a unary function or functor. More...
 
template<... >
void vectorToTensorMultiArray (...)
 Calculate the tensor (outer) product of a N-D vector with itself. More...
 

Detailed Description

Copy, transform, and inspect arbitrary dimensional arrays which are represented by iterators compatible to Multi-dimensional Array Iterators. Note that are range is here specified by a pair: an iterator referring to the first point of the array and a shape object specifying the size of the (rectangular) ROI.

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

Function Documentation

void vigra::initMultiArray (   ...)

Write a value to every element in a multi-dimensional array.

The initial value can either be a constant of appropriate type (compatible with the destination's value_type), or a functor with compatible result_type. These two cases are automatically distinguished when FunctorTraits<FUNCTOR>::isInitializer yields VigraTrueType. Since the functor is passed by const reference, its operator() must be const, and its internal state may need to be mutable.

Declarations:

pass arbitrary-dimensional array views:

namespace vigra {
template <unsigned int N, class T, class S, class VALUETYPE>
void
initMultiArray(MultiArrayView<N, T, S> s, VALUETYPE const & v);
template <unsigned int N, class T, class S, class FUNCTOR>
void
initMultiArray(MultiArrayView<N, T, S> s, FUNCTOR const & f);
}

show deprecated declarations

Usage:

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

MultiArray<3, unsigned int> array(Shape3(100, 200, 50));
// make an array of all ones
initMultiArray(array, 1);
// equivalent calls:
array = 1;
array.init(1);
// fill the array with random numbers
#include <vigra/random.hxx>

show deprecated examples

void vigra::initMultiArrayBorder (   ...)

Write values to the specified border values in the array.

This functions is similar to initMultiArray(), but it initializes only the array elements whose distance from any array border is at most border_width.

Declarations:

pass arbitrary-dimensional array views:

namespace vigra {
// init equal borders on all array sides
template <unsigned int N, class T, class S,
class VALUETYPE>
void
initMultiArrayBorder( MultiArrayView<N, T, S> array,
MultiArrayIndex border_width, VALUETYPE const & v);
template <unsigned int N, class T, class S,
class FUNCTOR>
void
initMultiArrayBorder( MultiArrayView<N, T, S> array,
MultiArrayIndex border_width, FUNCTOR const & v);
// specify border width individually for all array sides
template <unsigned int N, class T, class S,
class VALUETYPE>
void
initMultiArrayBorder( MultiArrayView<N, T, S> array,
typename MultiArrayShape<N>::type const & lower_border,
typename MultiArrayShape<N>::type const & upper_border,
VALUETYPE const & v);
}

show deprecated declarations

Usage:

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

MultiArray<3, unsigned int> array(Shape3(100, 200, 50));
int border_width = 5;
// init the array interior to 1, the border to 2
initMultiArray(array.subarray(Shape3(border_width), Shape3(-border_width)), 1);
initMultiArrayBorder(array, border_width, 2);
void vigra::copyMultiArray (   ...)

Copy a multi-dimensional array.

This function can be applied in two modes:

Standard Mode:
If the source and destination arrays have the same size, the corresponding array elements are simply copied. If necessary, type conversion takes place.
Expanding Mode:
If the source array has length 1 along some (or even all) dimensions, the source value at index 0 is used for all destination elements in those dimensions. For example, if we have single row of data (column length is 1), we can copy it into a 2D image of the same width: The given row is automatically repeated for every row of the destination image. Again, type conversion os performed if necessary.

The arrays must be represented by iterators compatible with vigra::MultiIterator, and the iteration range is specified by means of shape objects. If only the source shape is given the destination array is assumed to have the same shape, and standard mode is applied. If two shapes are given, the size of corresponding dimensions must be either equal (standard copy), or the source length must be 1 (expanding copy).

Declarations:

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

pass arbitrary-dimensional array views:

namespace vigra {
template <unsigned int N, class T1, class S1,
class T2, class S2>
void
copyMultiArray(MultiArrayView<N, T1, S1> const & source,
MultiArrayView<N, T2, S2> dest);
}

show deprecated declarations

Usage - Standard Mode:

MultiArray<3, int> src(Shape3(100, 200, 50)),
dest(Shape3(100, 200, 50));
...
copyMultiArray(src, dest);
// equivalent to
dest = src;
// copy only the red channel (i.e. channl 0) of an RGB array
MultiArray<3, RGBValue<int> > rgb_src(Shape3(100, 200, 50));
copyMultiArray(rgb_src.bindElementChannel(0), dest);
// equivalent to
dest = rgb_src.bindElementChannel(0);

Usage - Expanding Mode:

The source array is effectively only a 2D image (it has a 3D shape, but 'depth' is a singleton dimension with length 1). Thus, the destination will contain 50 identical copies of this image:

MultiArray<3, int> src(Shape2(100, 200)),
dest(Shape3(100, 200, 50));
...
copyMultiArray(src.insertSingletonDimension(2), dest);
// create an RGB image with three identical color bands
MultiArray<3, RGBValue<int> > rgb_dest(Shape2(100, 200));
copyMultiArray(src.insertSingletonDimension(2), rgb_dest.expandElements(2));

show deprecated examples

void vigra::transformMultiArray (   ...)

Transform a multi-dimensional array with a unary function or functor.

Note: The effect of this function can often be achieved in a simpler and more readable way by means of array expressions.

This function can be applied in three modes:

Standard Mode:
If the source and destination arrays have the same size, the transformation given by the functor is applied to every source element and the result written into the corresponding destination element. Unary functions, unary functors from the STL and the functors specifically defined in Functors to Transform Images can be used in standard mode. Creation of new functors is easiest by using Functor Expressions.
Expanding Mode:
If the source array has length 1 along some (or even all) dimensions, the source value at index 0 is used for all destination elements in those dimensions. In other words, the source index is not incremented along these dimensions, but the transformation functor is applied as usual. So, we can expand a small array (e.g. a single row of data, column length is 1), into a larger one (e.g. a 2D image with the same width): the given values are simply reused as necessary (e.g. for every row of the destination image). The same functors as in standard mode can be applied.
Reducing Mode:
If the destination array has length 1 along some (or even all) dimensions, the source values in these dimensions are reduced to single values by means of a suitable functor (e.g. vigra::ReduceFunctor), which supports two function call operators: one with a single argument to collect the values, and without argument to obtain the final (reduced) result. This behavior is a multi-dimensional generalization of the C++ standard function std::accumulate().

The arrays must be represented by MultiArrayViews. If source and destination shapes match, standard mode is applied. If the shapes differ, the size of corresponding dimensions must either be equal, or the source length must be 1 (expand mode), or the destination length must be 1 (reduce mode). However, reduction and expansion cannot be executed at the same time, so the latter conditions are mutual exclusive, even if they apply to different dimensions.

Declarations:

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

pass arbitrary-dimensional array views:

namespace vigra {
template <unsigned int N, class T1, class S1,
class T2, class S2,
class Functor>
void
transformMultiArray(MultiArrayView<N, T1, S1> const & source,
MultiArrayView<N, T2, S2> dest, Functor const & f);
}

show deprecated declarations

Usage - Standard Mode:

Source and destination array have the same size.

#include <cmath> // for sqrt()
MultiArray<3, float> src(Shape3(100, 200, 50)),
dest(Shape3(100, 200, 50));
...
transformMultiArray(src, dest, &std::sqrt );

Usage - Expand Mode:

The source array is effectively only a 2D image(it has a 3D shape, but depth is a singleton dimension with length 1). Thus, the destination will contain 50 identical copies of the transformed source image.

#include <cmath> // for sqrt()
MultiArray<3, float> src(Shape3(100, 200, 1)),
dest(Shape3(100, 200, 50));
...
transformMultiArray(src, dest, &std::sqrt );

Usage - Reduce Mode:

The destination array is effectively only 1D (it's width and height are singleton dimensions). Thus, it will contain accumulated data for every slice of the source volume (or for every frame, if the source is interpreted as an image sequence). In the example, we use the functor vigra::FindAverage to calculate the average gray value of every slice.

MultiArray<3, float> src(Shape3(100, 200, 50)),
dest(Shape3(1, 1, 50));
...
transformMultiArray(src, dest,
FindAverage<float>() );

Note that the functor must define the appropriate traits described below in order to be recognized as a reduce functor. This is most easily achieved by deriving from UnaryReduceFunctorTag (see vigra::FunctorTraits).

show deprecated examples

Required Interface:

In standard and expand mode, the functor must be a model of UnaryFunction (i.e. support one-argument function call which accepts values of type T1 and a return value that is convertible into T2.

In reduce mode, it must be a model of UnaryAnalyser (i.e. support function call with one argument and no return value functor(arg)) and Initializer (i.e. support function call with no argument, but return value res = functor()). Internally, such functors are recognized by the meta functions FunctorTraits<FUNCTOR>::isUnaryAnalyser and FunctorTraits<FUNCTOR>::isInitializer which must both yield VigraTrueType. Make sure that your functor correctly defines FunctorTraits because otherwise reduce mode will not work. This is most easily achieved by deriving the functor from UnaryReduceFunctorTag (see vigra::FunctorTraits). In addition, the functor must be copy constructible in order to start each reduction with a fresh functor.

See Also
Functors to Transform Images, vigra::multi_math, Functor Expressions
Examples:
graph_agglomerative_clustering.cxx.
void vigra::combineTwoMultiArrays (   ...)

Combine two multi-dimensional arrays into one using a binary function or functor.

Note: The effect of this function can often be achieved in a simpler and more readable way by means of array expressions.

This function can be applied in three modes:

Standard Mode:
If the source and destination arrays have the same size, the transformation given by the functor is applied to every pair of corresponding source elements and the result written into the corresponding destination element. Binary functions, binary functors from the STL and the functors specifically defined in Functors to Combine Images can be used in standard mode. Creation of new functors is easiest by using Functor Expressions.
Expanding Mode:
If the source arrays have length 1 along some (or even all) dimensions, the source values at index 0 are used for all destination elements in those dimensions. In other words, the source index is not incremented along those dimensions, but the transformation functor is applied as usual. So, we can expand small arrays (e.g. a single row of data, column length is 1), into larger ones (e.g. a 2D image with the same width): the given values are simply reused as necessary (e.g. for every row of the destination image). It is not even necessary that the source array shapes are equal. For example, we can combine a small array with one that hase the same size as the destination array. The same functors as in standard mode can be applied.
Reducing Mode:
If the destination array has length 1 along some (or even all) dimensions, the source values in these dimensions are reduced to single values by means of a suitable functor which supports two function call operators: one with two arguments to collect the values, and one without argument to obtain the final (reduced) result. This behavior is a multi-dimensional generalization of the C++ standard function std::accumulate().

The arrays must be represented by MultiArrayViews. If all shapes are identical, standard mode is applied. If the shapes differ, the size of corresponding dimensions must either be equal, or the length of this dimension must be 1 in one or both source arrays (expand mode), or the destination length must be 1 (reduce mode). However, reduction and expansion cannot be executed at the same time, so the latter conditions are mutual exclusive, even if they apply to different dimensions.

Declarations:

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

pass arbitrary-dimensional array views:

namespace vigra {
template <unsigned int N, class T11, class S11,
class T12, class S12,
class T2, class S2,
class Functor>
void
combineTwoMultiArrays(MultiArrayView<N, T11, S11> const & source1,
MultiArrayView<N, T12, S12> const & source2,
MultiArrayView<N, T2, S2> dest,
Functor const & f);
}

show deprecated declarations

Usage - Standard Mode:

Source and destination arrays have the same size.

#include <functional> // for std::plus
MultiArray<3, int> src1(Shape3(100, 200, 50)),
src2(Shape3(100, 200, 50)),
dest(Shape3(100, 200, 50));
...
combineTwoMultiArrays(src1, src2, dest,
std::plus<int>());

Usage - Expand Mode:

One source array is effectively only a 2D image (it has depth 1). This image will be added to every slice of the other source array, and the result is written into the corresponding destination slice.

#include <functional> // for std::plus
MultiArray<3, int> src1(Shape3(100, 200, 1)),
src2(Shape3(100, 200, 50)),
dest(Shape3(100, 200, 50));
...
combineTwoMultiArrays(src1, src2, dest,
std::plus<int>());

Usage - Reduce Mode:

The destination array is only 1D (it's width and height are singleton dimensions). Thus, it will contain accumulated data for every slice of the source volumes (or for every frame, if the sources are interpreted as image sequences). In the example, we use vigra::ReduceFunctor together with a functor expression (see Functor Expressions) to calculate the total absolute difference of the gray values in every pair of source slices.

#include <vigra/functorexpression.hxx>
using namespace vigra::functor;
MultiArray<3, int> src1(Shape3(100, 200, 50)),
src2(Shape3(100, 200, 50)),
dest(Shape3(1, 1, 50));
...
combineTwoMultiArrays(src1, src2, dest,
reduceFunctor(Arg1() + abs(Arg2() - Arg3()), 0) );
// Arg1() is the sum accumulated so far, initialized with 0

Note that the functor must define the appropriate traits described below in order to be recognized as a reduce functor. This is most easily achieved by deriving from BinaryReduceFunctorTag (see vigra::FunctorTraits).

show deprecated examples

Required Interface:

In standard and expand mode, the functor must be a model of BinaryFunction (i.e. support function call with two arguments and a return value which is convertible into T2: T2 res = functor(arg1, arg2)):

In reduce mode, it must be a model of BinaryAnalyser (i.e. support function call with two arguments and no return value functor(arg1, arg2)) and Initializer (i.e. support function call with no argument, but return value res = functor()). Internally, such functors are recognized by the meta functions FunctorTraits<FUNCTOR>::isBinaryAnalyser and FunctorTraits<FUNCTOR>::isInitializer which must both yield VigraTrueType. Make sure that your functor correctly defines FunctorTraits because otherwise reduce mode will not work. This is most easily achieved by deriving the functor from BinaryReduceFunctorTag (see vigra::FunctorTraits). In addition, the functor must be copy constructible in order to start each reduction with a fresh functor.

See Also
Functors to Transform Images, vigra::multi_math, Functor Expressions
void vigra::combineThreeMultiArrays (   ...)

Combine three multi-dimensional arrays into one using a ternary function or functor.

Note: The effect of this function can often be achieved in a simpler and more readable way by means of array expressions.

Except for the fact that it operates on three input arrays, this function is identical to the standard mode of combineTwoMultiArrays() (reduce and expand modes are not supported).

Declarations:

pass arbitrary-dimensional array views:

namespace vigra {
template <unsigned int N, class T11, class S11,
class T12, class S12,
class T13, class S13,
class T2, class S2,
class Functor>
void
combineThreeMultiArrays(MultiArrayView<N, T11, S11> const & source1,
MultiArrayView<N, T12, S12> const & source2,
MultiArrayView<N, T13, S13> const & source3,
MultiArrayView<N, T2, S2> dest,
Functor const & f);
}

show deprecated declarations

Usage:

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

#include <vigra/functorexpression.hxx>
MultiArray<3, int> src1(Shape3(100, 200, 50)),
src2(Shape3(100, 200, 50)),
src3(Shape3(100, 200, 50)),
dest(Shape3(100, 200, 50));
...
using namespace vigra::functor; // activate VIGRA's lambda library
combineThreeMultiArrays(src1, src2, src3, dest,
Arg1()*exp(-abs(Arg2()-Arg3())));

show deprecated examples

See Also
Functors to Transform Images, vigra::multi_math, Functor Expressions
void vigra::inspectMultiArray (   ...)

Call an analyzing functor at every element of a multi-dimensional array.

This function can be used to collect statistics of the array etc. The results must be stored in the functor, which serves as a return value (therefore, it is passed to the function by reference). The array must be represented as a MultiArrayView.

For many common statistics, the use of vigra::acc::extractFeatures() in combination with Feature Accumulators is more convenient.

Declarations:

pass arbitrary-dimensional array views:

namespace vigra {
template <unsigned int N, class T, class S,
class Functor>
void
inspectMultiArray(MultiArrayView<N, T, S> const & s,
Functor & f);
}

show deprecated declarations

Usage:

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

MultiArray<3, int> array(Shape3(100, 200, 50));
... // fill array
// init functor
FindMinMax<int> minmax;
inspectMultiArray(array, minmax);
cout << "Min: " << minmax.min << " Max: " << minmax.max;

The functor must support function call with one argument.

show deprecated examples

void vigra::inspectTwoMultiArrays (   ...)

Call an analyzing functor at all corresponding elements of two multi-dimensional arrays.

This function can be used to collect statistics over tow arrays. For example, one can holde data, and the other region labels or weights. The results must be stored in the functor, which serves as a return value (and is therefore passed by reference). The arrays must be represented by MultiArrayViews.

For many common statistics, the use of vigra::acc::extractFeatures() in combination with Feature Accumulators is more convenient.

Declarations:

pass arbitrary-dimensional array views:

namespace vigra {
template <unsigned int N, class T1, class S1,
class T2, class S2,
class Functor>
void
inspectTwoMultiArrays(MultiArrayView<N, T1, S1> const & s1,
MultiArrayView<N, T2, S2> const & s2,
Functor & f);
}

show deprecated declarations

Usage:

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

MultiArray<3, int> array1(Shape3(100, 200, 50)),
array2(Shape3(100, 200, 50));
// init functor
SomeStatisticsFunctor stats(..);
inspectTwoMultiArrays(array1, array2, stats);

The functor must support function call with two arguments.

show deprecated examples

void vigra::vectorToTensorMultiArray (   ...)

Calculate the tensor (outer) product of a N-D vector with itself.

This function is useful to transform vector arrays into a tensor representation that can be used as input to tensor based processing and analysis functions (e.g. tensor smoothing). When the input array has N dimensions, the input value_type must be a vector of length N, whereas the output value_type mus be vectors of length N*(N-1)/2 which will represent the upper triangular part of the resulting (symmetric) tensor. That is, for 2D arrays the output contains the elements [t11, t12 == t21, t22] in this order, whereas it contains the elements [t11, t12, t13, t22, t23, t33] for 3D arrays.

Currently, N <= 3 is required.

Declarations:

pass arbitrary-dimensional array views:

namespace vigra {
template <unsigned int N, class T1, class S1,
class T2, class S2>
void
vectorToTensorMultiArray(MultiArrayView<N, T1, S1> const & source,
MultiArrayView<N, T2, S2> dest);
}

show deprecated declarations

Usage:

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

MultiArray<3, float> vol(shape);
MultiArray<3, TinyVector<float, 3> > gradient(shape);
MultiArray<3, TinyVector<float, 6> > tensor(shape);
gaussianGradientMultiArray(vol, gradient, 2.0);
vectorToTensorMultiArray(gradient, tensor);
void vigra::tensorTraceMultiArray (   ...)

Calculate the tensor trace for every element of a N-D tensor array.

This function turns a N-D tensor (whose value_type is a vector of length N*(N+1)/2, see vectorToTensorMultiArray()) representing the upper triangular part of a symmetric tensor into a scalar array holding the tensor trace.

Currently, N <= 3 is required.

Declarations:

pass arbitrary-dimensional array views:

namespace vigra {
template <unsigned int N, class T1, class S1,
class T2, class S2>
void
tensorTraceMultiArray(MultiArrayView<N, T1, S1> const & source,
MultiArrayView<N, T2, S2> dest);
}

show deprecated declarations

Usage:

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

MultiArray<3, float> vol(shape);
MultiArray<3, TinyVector<float, 6> > hessian(shape);
MultiArray<3, float> trace(shape);
hessianOfGaussianMultiArray(vol, hessian, 2.0);

Preconditions:

N == 2 or N == 3

void vigra::tensorEigenvaluesMultiArray (   ...)

Calculate the tensor eigenvalues for every element of a N-D tensor array.

This function turns a N-D tensor (whose value_type is a vector of length N*(N+1)/2, see vectorToTensorMultiArray()) representing the upper triangular part of a symmetric tensor into a vector-valued array holding the tensor eigenvalues (thus, the destination value_type must be vectors of length N).

Currently, N <= 3 is required.

Declarations:

pass arbitrary-dimensional array views:

namespace vigra {
template <unsigned int N, class T1, class S1,
class T2, class S2>
void
tensorEigenvaluesMultiArray(MultiArrayView<N, T1, S1> const & source,
MultiArrayView<N, T2, S2> dest);
}

show deprecated declarations

Usage (MultiArrayView API):

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

MultiArray<3, float> vol(shape);
MultiArray<3, TinyVector<float, 6> > hessian(shape);
MultiArray<3, TinyVector<float, 3> > eigenvalues(shape);
hessianOfGaussianMultiArray(vol, hessian, 2.0);
tensorEigenvaluesMultiArray(hessian, eigenvalues);

Preconditions:

N == 2 or N == 3

void vigra::tensorDeterminantMultiArray (   ...)

Calculate the tensor determinant for every element of a ND tensor array.

This function turns a N-D tensor (whose value_type is a vector of length N*(N+1)/2, see vectorToTensorMultiArray()) representing the upper triangular part of a symmetric tensor into the a scalar array holding the tensor determinant.

Currently, N <= 3 is required.

Declarations:

pass arbitrary-dimensional array views:

namespace vigra {
template <unsigned int N, class T1, class S1,
class T2, class S2>
void
tensorDeterminantMultiArray(MultiArrayView<N, T1, S1> const & source,
MultiArrayView<N, T2, S2> dest);
}

show deprecated declarations

Usage (MultiArrayView API):

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

MultiArray<3, float> vol(shape);
MultiArray<3, TinyVector<float, 6> > hessian(shape);
MultiArray<3, float> determinant(shape);
hessianOfGaussianMultiArray(vol, hessian, 2.0);

Preconditions:

N == 2 or N == 3

© 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.1 (Fri May 19 2017)