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

vigra/multi_convolution.hxx |

    00001 //-- -*- c++ -*-
    00002 /************************************************************************/
    00003 /*                                                                      */
    00004 /*               Copyright 2003 by Christian-Dennis Rahn                */
    00005 /*                        and Ullrich Koethe                            */
    00006 /*                                                                      */
    00007 /*    This file is part of the VIGRA computer vision library.           */
    00008 /*    The VIGRA Website is                                              */
    00009 /*        http://hci.iwr.uni-heidelberg.de/vigra/                       */
    00010 /*    Please direct questions, bug reports, and contributions to        */
    00011 /*        ullrich.koethe@iwr.uni-heidelberg.de    or                    */
    00012 /*        vigra@informatik.uni-hamburg.de                               */
    00013 /*                                                                      */
    00014 /*    Permission is hereby granted, free of charge, to any person       */
    00015 /*    obtaining a copy of this software and associated documentation    */
    00016 /*    files (the "Software"), to deal in the Software without           */
    00017 /*    restriction, including without limitation the rights to use,      */
    00018 /*    copy, modify, merge, publish, distribute, sublicense, and/or      */
    00019 /*    sell copies of the Software, and to permit persons to whom the    */
    00020 /*    Software is furnished to do so, subject to the following          */
    00021 /*    conditions:                                                       */
    00022 /*                                                                      */
    00023 /*    The above copyright notice and this permission notice shall be    */
    00024 /*    included in all copies or substantial portions of the             */
    00025 /*    Software.                                                         */
    00026 /*                                                                      */
    00027 /*    THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND    */
    00028 /*    EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES   */
    00029 /*    OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND          */
    00030 /*    NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT       */
    00031 /*    HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,      */
    00032 /*    WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING      */
    00033 /*    FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR     */
    00034 /*    OTHER DEALINGS IN THE SOFTWARE.                                   */
    00035 /*                                                                      */
    00036 /************************************************************************/
    00037
    00038 #ifndef VIGRA_MULTI_CONVOLUTION_H
    00039 #define VIGRA_MULTI_CONVOLUTION_H
    00040
    00041 #include "separableconvolution.hxx"
    00042 #include "array_vector.hxx"
    00043 #include "multi_array.hxx"
    00044 #include "accessor.hxx"
    00045 #include "numerictraits.hxx"
    00046 #include "navigator.hxx"
    00047 #include "metaprogramming.hxx"
    00048 #include "multi_pointoperators.hxx"
    00049 #include "functorexpression.hxx"
    00050
    00051 namespace vigra
    00052 {
    00053
    00054
    00055 namespace detail
    00056 {
    00057
    00058 /********************************************************/
    00059 /*                                                      */
    00060 /*        internalSeparableConvolveMultiArray           */
    00061 /*                                                      */
    00062 /********************************************************/
    00063
    00064 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00065           class DestIterator, class DestAccessor, class KernelIterator>
    00066 void
    00067 internalSeparableConvolveMultiArrayTmp(
    00068                       SrcIterator si, SrcShape const & shape, SrcAccessor src,
    00069                       DestIterator di, DestAccessor dest, KernelIterator kit)
    00070 {
    00071     enum { N = 1 + SrcIterator::level };
    00072
    00073     typedef typename NumericTraits<typename DestAccessor::value_type>::RealPromote TmpType;
    00074
    00075     // temporay array to hold the current line to enable in-place operation
    00076     ArrayVector<TmpType> tmp( shape[0] );
    00077
    00078     typedef MultiArrayNavigator<SrcIterator, N> SNavigator;
    00079     typedef MultiArrayNavigator<DestIterator, N> DNavigator;
    00080
    00081     { // only operate on first dimension here
    00082         SNavigator snav( si, shape, 0 );
    00083         DNavigator dnav( di, shape, 0 );
    00084
    00085         for( ; snav.hasMore(); snav++, dnav++ )
    00086         {
    00087              // first copy source to temp for maximum cache efficiency
    00088              copyLine( snav.begin(), snav.end(), src,
    00089                        tmp.begin(), typename AccessorTraits<TmpType>::default_accessor() );
    00090
    00091              convolveLine( srcIterRange(tmp.begin(), tmp.end(),
    00092                                         typename AccessorTraits<TmpType>::default_const_accessor()),
    00093                            destIter( dnav.begin(), dest ),
    00094                            kernel1d( *kit ) );
    00095         }
    00096         ++kit;
    00097     }
    00098
    00099     // operate on further dimensions
    00100     for( int d = 1; d < N; ++d, ++kit )
    00101     {
    00102         DNavigator dnav( di, shape, d );
    00103
    00104         tmp.resize( shape[d] );
    00105
    00106         for( ; dnav.hasMore(); dnav++ )
    00107         {
    00108              // first copy source to temp for maximum cache efficiency
    00109              copyLine( dnav.begin(), dnav.end(), dest,
    00110                        tmp.begin(), typename AccessorTraits<TmpType>::default_accessor() );
    00111
    00112              convolveLine( srcIterRange(tmp.begin(), tmp.end(),
    00113                                         typename AccessorTraits<TmpType>::default_const_accessor()),
    00114                            destIter( dnav.begin(), dest ),
    00115                            kernel1d( *kit ) );
    00116         }
    00117     }
    00118 }
    00119
    00120
    00121 } // namespace detail
    00122
    00123 /** \addtogroup MultiArrayConvolutionFilters Convolution filters for multi-dimensional arrays.
    00124
    00125     These functions realize a separable convolution on an arbitrary dimensional
    00126     array that is specified by iterators (compatible to \ref MultiIteratorPage)
    00127     and shape objects. It can therefore be applied to a wide range of data structures
    00128     (\ref vigra::MultiArrayView, \ref vigra::MultiArray etc.).
    00129 */
    00130 //@{
    00131
    00132 /********************************************************/
    00133 /*                                                      */
    00134 /*             separableConvolveMultiArray              */
    00135 /*                                                      */
    00136 /********************************************************/
    00137
    00138 /** \brief Separated convolution on multi-dimensional arrays.
    00139
    00140     This function computes a separated convolution on all dimensions
    00141     of the given multi-dimensional array. Both source and destination
    00142     arrays are represented by iterators, shape objects and accessors.
    00143     The destination array is required to already have the correct size.
    00144
    00145     There are two variants of this functions: one takes a single kernel
    00146     of type \ref vigra::Kernel1D which is then applied to all dimensions,
    00147     whereas the other requires an iterator referencing a sequence of
    00148     \ref vigra::Kernel1D objects, one for every dimension of the data.
    00149     Then the first kernel in this sequence is applied to the innermost
    00150     dimension (e.g. the x-dimension of an image), while the last is applied to the
    00151     outermost dimension (e.g. the z-dimension in a 3D image).
    00152
    00153     This function may work in-place, which means that <tt>siter == diter</tt> is allowed.
    00154     A full-sized internal array is only allocated if working on the destination
    00155     array directly would cause round-off errors (i.e. if
    00156     <tt>typeid(typename NumericTraits<typename DestAccessor::value_type>::RealPromote)
    00157     != typeid(typename DestAccessor::value_type)</tt>.
    00158
    00159     <b> Declarations:</b>
    00160
    00161     pass arguments explicitly:
    00162     \code
    00163     namespace vigra {
    00164         // apply the same kernel to all dimensions
    00165         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00166                   class DestIterator, class DestAccessor, class T>
    00167         void
    00168         separableConvolveMultiArray(SrcIterator siter, SrcShape const & shape, SrcAccessor src,
    00169                                     DestIterator diter, DestAccessor dest,
    00170                                     Kernel1D<T> const & kernel);
    00171
    00172         // apply each kernel from the sequence 'kernels' in turn
    00173         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00174                   class DestIterator, class DestAccessor, class KernelIterator>
    00175         void
    00176         separableConvolveMultiArray(SrcIterator siter, SrcShape const & shape, SrcAccessor src,
    00177                                     DestIterator diter, DestAccessor dest,
    00178                                     KernelIterator kernels);
    00179     }
    00180     \endcode
    00181
    00182     use argument objects in conjunction with \ref ArgumentObjectFactories :
    00183     \code
    00184     namespace vigra {
    00185         // apply the same kernel to all dimensions
    00186         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00187                   class DestIterator, class DestAccessor, class T>
    00188         void
    00189         separableConvolveMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    00190                                     pair<DestIterator, DestAccessor> const & dest,
    00191                                     Kernel1D<T> const & kernel);
    00192
    00193         // apply each kernel from the sequence 'kernels' in turn
    00194         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00195                   class DestIterator, class DestAccessor, class KernelIterator>
    00196         void
    00197         separableConvolveMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    00198                                     pair<DestIterator, DestAccessor> const & dest,
    00199                                     KernelIterator kernels);
    00200     }
    00201     \endcode
    00202
    00203     <b> Usage:</b>
    00204
    00205     <b>\#include</b> <<a href="multi__convolution_8hxx-source.html">vigra/multi_convolution.hxx</a>>
    00206
    00207     \code
    00208     MultiArray<3, unsigned char>::size_type shape(width, height, depth);
    00209     MultiArray<3, unsigned char> source(shape);
    00210     MultiArray<3, float> dest(shape);
    00211     ...
    00212     Kernel1D<float> gauss;
    00213     gauss.initGaussian(sigma);
    00214     // create 3 Gauss kernels, one for each dimension
    00215     ArrayVector<Kernel1D<float> > kernels(3, gauss);
    00216
    00217     // perform Gaussian smoothing on all dimensions
    00218     separableConvolveMultiArray(srcMultiArrayRange(source), destMultiArray(dest),
    00219                                 kernels.begin());
    00220     \endcode
    00221
    00222     \see vigra::Kernel1D, convolveLine()
    00223 */
    00224 doxygen_overloaded_function(template <...> void separableConvolveMultiArray)
    00225
    00226 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00227           class DestIterator, class DestAccessor, class KernelIterator>
    00228 void
    00229 separableConvolveMultiArray( SrcIterator s, SrcShape const & shape, SrcAccessor src,
    00230                              DestIterator d, DestAccessor dest, KernelIterator kernels )
    00231 {
    00232     typedef typename NumericTraits<typename DestAccessor::value_type>::RealPromote TmpType;
    00233
    00234     if(!IsSameType<TmpType, typename DestAccessor::value_type>::boolResult)
    00235     {
    00236         // need a temporary array to avoid rounding errors
    00237         MultiArray<SrcShape::static_size, TmpType> tmpArray(shape);
    00238         detail::internalSeparableConvolveMultiArrayTmp( s, shape, src,
    00239              tmpArray.traverser_begin(), typename AccessorTraits<TmpType>::default_accessor(), kernels );
    00240         copyMultiArray(srcMultiArrayRange(tmpArray), destIter(d, dest));
    00241     }
    00242     else
    00243     {
    00244         // work directly on the destination array
    00245         detail::internalSeparableConvolveMultiArrayTmp( s, shape, src, d, dest, kernels );
    00246     }
    00247 }
    00248
    00249 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00250           class DestIterator, class DestAccessor, class KernelIterator>
    00251 inline
    00252 void separableConvolveMultiArray(
    00253     triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    00254     pair<DestIterator, DestAccessor> const & dest, KernelIterator kit )
    00255 {
    00256     separableConvolveMultiArray( source.first, source.second, source.third,
    00257                                  dest.first, dest.second, kit );
    00258 }
    00259
    00260 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00261           class DestIterator, class DestAccessor, class T>
    00262 inline void
    00263 separableConvolveMultiArray( SrcIterator s, SrcShape const & shape, SrcAccessor src,
    00264                              DestIterator d, DestAccessor dest,
    00265                              Kernel1D<T> const & kernel )
    00266 {
    00267     ArrayVector<Kernel1D<T> > kernels(shape.size(), kernel);
    00268
    00269     separableConvolveMultiArray( s, shape, src, d, dest, kernels.begin() );
    00270 }
    00271
    00272 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00273           class DestIterator, class DestAccessor, class T>
    00274 inline void
    00275 separableConvolveMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    00276                             pair<DestIterator, DestAccessor> const & dest,
    00277                             Kernel1D<T> const & kernel )
    00278 {
    00279     ArrayVector<Kernel1D<T> > kernels(source.second.size(), kernel);
    00280
    00281     separableConvolveMultiArray( source.first, source.second, source.third,
    00282                                  dest.first, dest.second, kernels.begin() );
    00283 }
    00284
    00285 /********************************************************/
    00286 /*                                                      */
    00287 /*            convolveMultiArrayOneDimension            */
    00288 /*                                                      */
    00289 /********************************************************/
    00290
    00291 /** \brief Convolution along a single dimension of a multi-dimensional arrays.
    00292
    00293     This function computes a convolution along one dimension (specified by
    00294     the parameter <tt>dim</tt> of the given multi-dimensional array with the given
    00295     <tt>kernel</tt>. Both source and destination arrays are represented by
    00296     iterators, shape objects and accessors. The destination array is required to
    00297     already have the correct size.
    00298
    00299     This function may work in-place, which means that <tt>siter == diter</tt> is allowed.
    00300
    00301     <b> Declarations:</b>
    00302
    00303     pass arguments explicitly:
    00304     \code
    00305     namespace vigra {
    00306         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00307                   class DestIterator, class DestAccessor, class T>
    00308         void
    00309         convolveMultiArrayOneDimension(SrcIterator siter, SrcShape const & shape, SrcAccessor src,
    00310                                        DestIterator diter, DestAccessor dest,
    00311                                        unsigned int dim, vigra::Kernel1D<T> const & kernel);
    00312     }
    00313     \endcode
    00314
    00315     use argument objects in conjunction with \ref ArgumentObjectFactories :
    00316     \code
    00317     namespace vigra {
    00318         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00319                   class DestIterator, class DestAccessor, class T>
    00320         void
    00321         convolveMultiArrayOneDimension(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    00322                                        pair<DestIterator, DestAccessor> const & dest,
    00323                                        unsigned int dim, vigra::Kernel1D<T> const & kernel);
    00324     }
    00325     \endcode
    00326
    00327     <b> Usage:</b>
    00328
    00329     <b>\#include</b> <<a href="multi__convolution_8hxx-source.html">vigra/multi_convolution.hxx</a>>
    00330
    00331     \code
    00332     MultiArray<3, unsigned char>::size_type shape(width, height, depth);
    00333     MultiArray<3, unsigned char> source(shape);
    00334     MultiArray<3, float> dest(shape);
    00335     ...
    00336     Kernel1D<float> gauss;
    00337     gauss.initGaussian(sigma);
    00338
    00339     // perform Gaussian smoothing along dimensions 1 (height)
    00340     convolveMultiArrayOneDimension(srcMultiArrayRange(source), destMultiArray(dest), 1, gauss);
    00341     \endcode
    00342
    00343     \see separableConvolveMultiArray()
    00344 */
    00345 doxygen_overloaded_function(template <...> void convolveMultiArrayOneDimension)
    00346
    00347 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00348           class DestIterator, class DestAccessor, class T>
    00349 void
    00350 convolveMultiArrayOneDimension(SrcIterator s, SrcShape const & shape, SrcAccessor src,
    00351                                DestIterator d, DestAccessor dest,
    00352                                unsigned int dim, vigra::Kernel1D<T> const & kernel )
    00353 {
    00354     enum { N = 1 + SrcIterator::level };
    00355     vigra_precondition( dim < N,
    00356                         "convolveMultiArrayOneDimension(): The dimension number to convolve must be smaller "
    00357                         "than the data dimensionality" );
    00358
    00359     typedef typename NumericTraits<typename DestAccessor::value_type>::RealPromote TmpType;
    00360     ArrayVector<TmpType> tmp( shape[dim] );
    00361
    00362     typedef MultiArrayNavigator<SrcIterator, N> SNavigator;
    00363     typedef MultiArrayNavigator<DestIterator, N> DNavigator;
    00364
    00365     SNavigator snav( s, shape, dim );
    00366     DNavigator dnav( d, shape, dim );
    00367
    00368     for( ; snav.hasMore(); snav++, dnav++ )
    00369     {
    00370          // first copy source to temp for maximum cache efficiency
    00371          copyLine( snav.begin(), snav.end(), src,
    00372            tmp.begin(), typename AccessorTraits<TmpType>::default_accessor() );
    00373
    00374          convolveLine( srcIterRange( tmp.begin(), tmp.end(), typename AccessorTraits<TmpType>::default_const_accessor()),
    00375                        destIter( dnav.begin(), dest ),
    00376                        kernel1d( kernel ) );
    00377     }
    00378 }
    00379
    00380 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00381           class DestIterator, class DestAccessor, class T>
    00382 inline void
    00383 convolveMultiArrayOneDimension(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    00384                                pair<DestIterator, DestAccessor> const & dest,
    00385                                unsigned int dim, vigra::Kernel1D<T> const & kernel )
    00386 {
    00387     convolveMultiArrayOneDimension( source.first, source.second, source.third,
    00388                                    dest.first, dest.second, dim, kernel );
    00389 }
    00390
    00391 /********************************************************/
    00392 /*                                                      */
    00393 /*             gaussianSmoothMultiArray                 */
    00394 /*                                                      */
    00395 /********************************************************/
    00396
    00397 /** \brief Isotropic Gaussian smoothing of a multi-dimensional arrays.
    00398
    00399     This function computes an isotropic convolution of the given multi-dimensional
    00400     array with a Gaussian filter at the given standard deviation <tt>sigma</tt>.
    00401     Both source and destination arrays are represented by
    00402     iterators, shape objects and accessors. The destination array is required to
    00403     already have the correct size. This function may work in-place, which means
    00404     that <tt>siter == diter</tt> is allowed. It is implemented by a call to
    00405     \ref separableConvolveMultiArray() with the appropriate kernel.
    00406     If the data are anisotropic (different pixel size along different dimensions)
    00407     you should call \ref separableConvolveMultiArray() directly with the appropriate
    00408     anisotropic Gaussians.
    00409
    00410     <b> Declarations:</b>
    00411
    00412     pass arguments explicitly:
    00413     \code
    00414     namespace vigra {
    00415         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00416                   class DestIterator, class DestAccessor>
    00417         void
    00418         gaussianSmoothMultiArray(SrcIterator siter, SrcShape const & shape, SrcAccessor src,
    00419                                  DestIterator diter, DestAccessor dest,
    00420                                  double sigma);
    00421     }
    00422     \endcode
    00423
    00424     use argument objects in conjunction with \ref ArgumentObjectFactories :
    00425     \code
    00426     namespace vigra {
    00427         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00428                   class DestIterator, class DestAccessor>
    00429         void
    00430         gaussianSmoothMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    00431                                  pair<DestIterator, DestAccessor> const & dest,
    00432                                  double sigma);
    00433     }
    00434     \endcode
    00435
    00436     <b> Usage:</b>
    00437
    00438     <b>\#include</b> <<a href="multi__convolution_8hxx-source.html">vigra/multi_convolution.hxx</a>>
    00439
    00440     \code
    00441     MultiArray<3, unsigned char>::size_type shape(width, height, depth);
    00442     MultiArray<3, unsigned char> source(shape);
    00443     MultiArray<3, float> dest(shape);
    00444     ...
    00445     // perform isotropic Gaussian smoothing at scale 'sigma'
    00446     gaussianSmoothMultiArray(srcMultiArrayRange(source), destMultiArray(dest), sigma);
    00447     \endcode
    00448
    00449     \see separableConvolveMultiArray()
    00450 */
    00451 doxygen_overloaded_function(template <...> void gaussianSmoothMultiArray)
    00452
    00453 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00454           class DestIterator, class DestAccessor>
    00455 void
    00456 gaussianSmoothMultiArray( SrcIterator s, SrcShape const & shape, SrcAccessor src,
    00457                    DestIterator d, DestAccessor dest, double sigma )
    00458 {
    00459     Kernel1D<double> gauss;
    00460     gauss.initGaussian( sigma );
    00461
    00462     separableConvolveMultiArray( s, shape, src, d, dest, gauss);
    00463 }
    00464
    00465 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00466           class DestIterator, class DestAccessor>
    00467 inline void
    00468 gaussianSmoothMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    00469                   pair<DestIterator, DestAccessor> const & dest,
    00470                   double sigma )
    00471 {
    00472     gaussianSmoothMultiArray( source.first, source.second, source.third,
    00473                               dest.first, dest.second, sigma );
    00474 }
    00475
    00476 /********************************************************/
    00477 /*                                                      */
    00478 /*             gaussianGradientMultiArray               */
    00479 /*                                                      */
    00480 /********************************************************/
    00481
    00482 /** \brief Calculate Gaussian gradient of a multi-dimensional arrays.
    00483
    00484     This function computes the Gaussian gradient of the given multi-dimensional
    00485     array with a sequence of first-derivative-of-Gaussian filters at the given
    00486     standard deviation <tt>sigma</tt> (differentiation is applied to each dimension
    00487     in turn, starting with the innermost dimension). Both source and destination arrays
    00488     are represented by iterators, shape objects and accessors. The destination array is
    00489     required to have a vector valued pixel type with as many elements as the number of
    00490     dimensions. This function is implemented by calls to
    00491     \ref separableConvolveMultiArray() with the appropriate kernels.
    00492     If the data are anisotropic (different pixel size along different dimensions)
    00493     you should call \ref separableConvolveMultiArray() directly with the appropriate
    00494     anisotropic Gaussian derivatives.
    00495
    00496     <b> Declarations:</b>
    00497
    00498     pass arguments explicitly:
    00499     \code
    00500     namespace vigra {
    00501         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00502                   class DestIterator, class DestAccessor>
    00503         void
    00504         gaussianGradientMultiArray(SrcIterator siter, SrcShape const & shape, SrcAccessor src,
    00505                                    DestIterator diter, DestAccessor dest,
    00506                                    double sigma);
    00507     }
    00508     \endcode
    00509
    00510     use argument objects in conjunction with \ref ArgumentObjectFactories :
    00511     \code
    00512     namespace vigra {
    00513         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00514                   class DestIterator, class DestAccessor>
    00515         void
    00516         gaussianGradientMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    00517                                    pair<DestIterator, DestAccessor> const & dest,
    00518                                    double sigma);
    00519     }
    00520     \endcode
    00521
    00522     <b> Usage:</b>
    00523
    00524     <b>\#include</b> <<a href="multi__convolution_8hxx-source.html">vigra/multi_convolution.hxx</a>>
    00525
    00526     \code
    00527     MultiArray<3, unsigned char>::size_type shape(width, height, depth);
    00528     MultiArray<3, unsigned char> source(shape);
    00529     MultiArray<3, TinyVector<float, 3> > dest(shape);
    00530     ...
    00531     // compute Gaussian gradient at scale sigma
    00532     gaussianGradientMultiArray(srcMultiArrayRange(source), destMultiArray(dest), sigma);
    00533     \endcode
    00534
    00535     <b> Required Interface:</b>
    00536
    00537     see \ref separableConvolveMultiArray(), in addition:
    00538
    00539     \code
    00540     int dimension = 0;
    00541     VectorElementAccessor<DestAccessor> elementAccessor(0, dest);
    00542     \endcode
    00543
    00544     \see separableConvolveMultiArray()
    00545 */
    00546 doxygen_overloaded_function(template <...> void gaussianGradientMultiArray)
    00547
    00548 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00549           class DestIterator, class DestAccessor>
    00550 void
    00551 gaussianGradientMultiArray(SrcIterator si, SrcShape const & shape, SrcAccessor src,
    00552                            DestIterator di, DestAccessor dest, double sigma )
    00553 {
    00554     typedef typename DestAccessor::value_type DestType;
    00555     typedef typename DestType::value_type     DestValueType;
    00556     typedef typename NumericTraits<DestValueType>::RealPromote KernelType;
    00557
    00558     static const int N = SrcShape::static_size;
    00559
    00560     for(int k=0; k<N; ++k)
    00561         if(shape[k] <=0)
    00562             return;
    00563
    00564     vigra_precondition(N == dest.size(di),
    00565         "gaussianGradientMultiArray(): Wrong number of channels in output array.");
    00566
    00567     vigra_precondition(sigma > 0.0, "gaussianGradientMultiArray(): Scale must be positive.");
    00568
    00569     Kernel1D<KernelType> gauss, derivative;
    00570     gauss.initGaussian(sigma);
    00571
    00572     typedef VectorElementAccessor<DestAccessor> ElementAccessor;
    00573
    00574     // compute gradient components
    00575     for(int d = 0; d < N; ++d )
    00576     {
    00577         ArrayVector<Kernel1D<KernelType> > kernels(N, gauss);
    00578         kernels[d].initGaussianDerivative(sigma, 1);
    00579         separableConvolveMultiArray( si, shape, src, di, ElementAccessor(d, dest), kernels.begin());
    00580     }
    00581 }
    00582
    00583 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00584           class DestIterator, class DestAccessor>
    00585 inline void
    00586 gaussianGradientMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    00587                            pair<DestIterator, DestAccessor> const & dest, double sigma )
    00588 {
    00589     gaussianGradientMultiArray( source.first, source.second, source.third,
    00590                                 dest.first, dest.second, sigma );
    00591 }
    00592
    00593 /********************************************************/
    00594 /*                                                      */
    00595 /*             symmetricGradientMultiArray              */
    00596 /*                                                      */
    00597 /********************************************************/
    00598
    00599 /** \brief Calculate gradient of a multi-dimensional arrays using symmetric difference filters.
    00600
    00601     This function computes the gradient of the given multi-dimensional
    00602     array with a sequence of symmetric difference filters a (differentiation is applied
    00603     to each dimension in turn, starting with the innermost dimension). Both source and
    00604     destination arrays are represented by iterators, shape objects and accessors.
    00605     The destination array is required to have a vector valued pixel type with as many
    00606     elements as the number of dimensions. This function is implemented by calls to
    00607     \ref convolveMultiArrayOneDimension() with the symmetric difference kernel.
    00608
    00609     <b> Declarations:</b>
    00610
    00611     pass arguments explicitly:
    00612     \code
    00613     namespace vigra {
    00614         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00615                   class DestIterator, class DestAccessor>
    00616         void
    00617         symmetricGradientMultiArray(SrcIterator siter, SrcShape const & shape, SrcAccessor src,
    00618                                     DestIterator diter, DestAccessor dest);
    00619     }
    00620     \endcode
    00621
    00622     use argument objects in conjunction with \ref ArgumentObjectFactories :
    00623     \code
    00624     namespace vigra {
    00625         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00626                   class DestIterator, class DestAccessor>
    00627         void
    00628         symmetricGradientMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    00629                                     pair<DestIterator, DestAccessor> const & dest);
    00630     }
    00631     \endcode
    00632
    00633     <b> Usage:</b>
    00634
    00635     <b>\#include</b> <<a href="multi__convolution_8hxx-source.html">vigra/multi_convolution.hxx</a>>
    00636
    00637     \code
    00638     MultiArray<3, unsigned char>::size_type shape(width, height, depth);
    00639     MultiArray<3, unsigned char> source(shape);
    00640     MultiArray<3, TinyVector<float, 3> > dest(shape);
    00641     ...
    00642     // compute gradient
    00643     symmetricGradientMultiArray(srcMultiArrayRange(source), destMultiArray(dest));
    00644     \endcode
    00645
    00646     <b> Required Interface:</b>
    00647
    00648     see \ref convolveMultiArrayOneDimension(), in addition:
    00649
    00650     \code
    00651     int dimension = 0;
    00652     VectorElementAccessor<DestAccessor> elementAccessor(0, dest);
    00653     \endcode
    00654
    00655     \see convolveMultiArrayOneDimension()
    00656 */
    00657 doxygen_overloaded_function(template <...> void symmetricGradientMultiArray)
    00658
    00659 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00660           class DestIterator, class DestAccessor>
    00661 void
    00662 symmetricGradientMultiArray(SrcIterator si, SrcShape const & shape, SrcAccessor src,
    00663                             DestIterator di, DestAccessor dest)
    00664 {
    00665     typedef typename DestAccessor::value_type DestType;
    00666     typedef typename DestType::value_type     DestValueType;
    00667     typedef typename NumericTraits<DestValueType>::RealPromote KernelType;
    00668
    00669     static const int N = SrcShape::static_size;
    00670
    00671     for(int k=0; k<N; ++k)
    00672         if(shape[k] <=0)
    00673             return;
    00674
    00675     vigra_precondition(N == dest.size(di),
    00676         "symmetricGradientMultiArray(): Wrong number of channels in output array.");
    00677
    00678     Kernel1D<KernelType> filter;
    00679     filter.initSymmetricGradient();
    00680
    00681     typedef VectorElementAccessor<DestAccessor> ElementAccessor;
    00682
    00683     // compute gradient components
    00684     for(int d = 0; d < N; ++d )
    00685     {
    00686         convolveMultiArrayOneDimension(si, shape, src,
    00687                                        di, ElementAccessor(d, dest),
    00688                                        d, filter);
    00689     }
    00690 }
    00691
    00692 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00693           class DestIterator, class DestAccessor>
    00694 inline void
    00695 symmetricGradientMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    00696                             pair<DestIterator, DestAccessor> const & dest )
    00697 {
    00698     symmetricGradientMultiArray(source.first, source.second, source.third,
    00699                                 dest.first, dest.second);
    00700 }
    00701
    00702
    00703 /********************************************************/
    00704 /*                                                      */
    00705 /*            laplacianOfGaussianMultiArray             */
    00706 /*                                                      */
    00707 /********************************************************/
    00708
    00709 /** \brief Calculate Laplacian of a N-dimensional arrays using Gaussian derivative filters.
    00710
    00711     This function computes the Laplacian the given N-dimensional
    00712     array with a sequence of second-derivative-of-Gaussian filters at the given
    00713     standard deviation <tt>sigma</tt>. Both source and destination arrays
    00714     are represented by iterators, shape objects and accessors. Both source and destination
    00715     arrays must have scalar value_type. This function is implemented by calls to
    00716     \ref separableConvolveMultiArray() with the appropriate kernels, followed by summation.
    00717
    00718     <b> Declarations:</b>
    00719
    00720     pass arguments explicitly:
    00721     \code
    00722     namespace vigra {
    00723         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00724                   class DestIterator, class DestAccessor>
    00725         void
    00726         laplacianOfGaussianMultiArray(SrcIterator siter, SrcShape const & shape, SrcAccessor src,
    00727                                       DestIterator diter, DestAccessor dest,
    00728                                       double sigma);
    00729     }
    00730     \endcode
    00731
    00732     use argument objects in conjunction with \ref ArgumentObjectFactories :
    00733     \code
    00734     namespace vigra {
    00735         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00736                   class DestIterator, class DestAccessor>
    00737         void
    00738         laplacianOfGaussianMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    00739                                       pair<DestIterator, DestAccessor> const & dest,
    00740                                       double sigma);
    00741     }
    00742     \endcode
    00743
    00744     <b> Usage:</b>
    00745
    00746     <b>\#include</b> <<a href="multi__convolution_8hxx-source.html">vigra/multi_convolution.hxx</a>>
    00747
    00748     \code
    00749     MultiArray<3, float> source(shape);
    00750     MultiArray<3, float> laplacian(shape);
    00751     ...
    00752     // compute Laplacian at scale sigma
    00753     laplacianOfGaussianMultiArray(srcMultiArrayRange(source), destMultiArray(laplacian), sigma);
    00754     \endcode
    00755
    00756     <b> Required Interface:</b>
    00757
    00758     see \ref separableConvolveMultiArray(), in addition:
    00759
    00760     \code
    00761     int dimension = 0;
    00762     VectorElementAccessor<DestAccessor> elementAccessor(0, dest);
    00763     \endcode
    00764
    00765     \see separableConvolveMultiArray()
    00766 */
    00767 doxygen_overloaded_function(template <...> void laplacianOfGaussianMultiArray)
    00768
    00769 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00770           class DestIterator, class DestAccessor>
    00771 void
    00772 laplacianOfGaussianMultiArray(SrcIterator si, SrcShape const & shape, SrcAccessor src,
    00773                               DestIterator di, DestAccessor dest, double sigma )
    00774 {
    00775     using namespace functor;
    00776
    00777     typedef typename DestAccessor::value_type DestType;
    00778     typedef typename NumericTraits<DestType>::RealPromote KernelType;
    00779     typedef typename AccessorTraits<KernelType>::default_accessor DerivativeAccessor;
    00780
    00781     static const int N = SrcShape::static_size;
    00782
    00783     vigra_precondition(sigma > 0.0, "laplacianOfGaussianMultiArray(): Scale must be positive.");
    00784
    00785     Kernel1D<KernelType> gauss;
    00786     gauss.initGaussian(sigma);
    00787
    00788     MultiArray<N, KernelType> derivative(shape);
    00789
    00790     // compute 2nd derivatives and sum them up
    00791     for(int d = 0; d < N; ++d )
    00792     {
    00793         ArrayVector<Kernel1D<KernelType> > kernels(N, gauss);
    00794         kernels[d].initGaussianDerivative(sigma, 2);
    00795         if(d == 0)
    00796         {
    00797             separableConvolveMultiArray( si, shape, src,
    00798                                          di, dest, kernels.begin());
    00799         }
    00800         else
    00801         {
    00802             separableConvolveMultiArray( si, shape, src,
    00803                                          derivative.traverser_begin(), DerivativeAccessor(),
    00804                                          kernels.begin());
    00805             combineTwoMultiArrays(di, shape, dest, derivative.traverser_begin(), DerivativeAccessor(),
    00806                                   di, dest, Arg1() + Arg2() );
    00807         }
    00808     }
    00809 }
    00810
    00811 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00812           class DestIterator, class DestAccessor>
    00813 inline void
    00814 laplacianOfGaussianMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    00815                             pair<DestIterator, DestAccessor> const & dest, double sigma )
    00816 {
    00817     laplacianOfGaussianMultiArray( source.first, source.second, source.third,
    00818                                    dest.first, dest.second, sigma );
    00819 }
    00820
    00821 /********************************************************/
    00822 /*                                                      */
    00823 /*              hessianOfGaussianMultiArray             */
    00824 /*                                                      */
    00825 /********************************************************/
    00826
    00827 /** \brief Calculate Hessian matrix of a N-dimensional arrays using Gaussian derivative filters.
    00828
    00829     This function computes the Hessian matrix the given scalar N-dimensional
    00830     array with a sequence of second-derivative-of-Gaussian filters at the given
    00831     standard deviation <tt>sigma</tt>. Both source and destination arrays
    00832     are represented by iterators, shape objects and accessors. The destination array must
    00833     have a vector valued element type with N*(N+1)/2 elements (it represents the
    00834     upper triangular part of the symmetric Hessian matrix). This function is implemented by calls to
    00835     \ref separableConvolveMultiArray() with the appropriate kernels.
    00836
    00837     <b> Declarations:</b>
    00838
    00839     pass arguments explicitly:
    00840     \code
    00841     namespace vigra {
    00842         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00843                   class DestIterator, class DestAccessor>
    00844         void
    00845         hessianOfGaussianMultiArray(SrcIterator siter, SrcShape const & shape, SrcAccessor src,
    00846                                     DestIterator diter, DestAccessor dest,
    00847                                     double sigma);
    00848     }
    00849     \endcode
    00850
    00851     use argument objects in conjunction with \ref ArgumentObjectFactories :
    00852     \code
    00853     namespace vigra {
    00854         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00855                   class DestIterator, class DestAccessor>
    00856         void
    00857         hessianOfGaussianMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    00858                                     pair<DestIterator, DestAccessor> const & dest,
    00859                                     double sigma);
    00860     }
    00861     \endcode
    00862
    00863     <b> Usage:</b>
    00864
    00865     <b>\#include</b> <<a href="multi__convolution_8hxx-source.html">vigra/multi_convolution.hxx</a>>
    00866
    00867     \code
    00868     MultiArray<3, float> source(shape);
    00869     MultiArray<3, TinyVector<float, 6> > dest(shape);
    00870     ...
    00871     // compute Hessian at scale sigma
    00872     hessianOfGaussianMultiArray(srcMultiArrayRange(source), destMultiArray(dest), sigma);
    00873     \endcode
    00874
    00875     <b> Required Interface:</b>
    00876
    00877     see \ref separableConvolveMultiArray(), in addition:
    00878
    00879     \code
    00880     int dimension = 0;
    00881     VectorElementAccessor<DestAccessor> elementAccessor(0, dest);
    00882     \endcode
    00883
    00884     \see separableConvolveMultiArray(), vectorToTensorMultiArray()
    00885 */
    00886 doxygen_overloaded_function(template <...> void hessianOfGaussianMultiArray)
    00887
    00888 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00889           class DestIterator, class DestAccessor>
    00890 void
    00891 hessianOfGaussianMultiArray(SrcIterator si, SrcShape const & shape, SrcAccessor src,
    00892                             DestIterator di, DestAccessor dest, double sigma )
    00893 {
    00894     typedef typename DestAccessor::value_type DestType;
    00895     typedef typename DestType::value_type     DestValueType;
    00896     typedef typename NumericTraits<DestValueType>::RealPromote KernelType;
    00897
    00898     static const int N = SrcShape::static_size;
    00899     static const int M = N*(N+1)/2;
    00900
    00901     for(int k=0; k<N; ++k)
    00902         if(shape[k] <=0)
    00903             return;
    00904
    00905     vigra_precondition(M == dest.size(di),
    00906         "hessianOfGaussianMultiArray(): Wrong number of channels in output array.");
    00907
    00908     vigra_precondition(sigma > 0.0, "hessianOfGaussianMultiArray(): Scale must be positive.");
    00909
    00910     Kernel1D<KernelType> gauss;
    00911     gauss.initGaussian(sigma);
    00912
    00913     typedef VectorElementAccessor<DestAccessor> ElementAccessor;
    00914
    00915     // compute elements of the Hessian matrix
    00916     for(int b=0, i=0; i<N; ++i)
    00917     {
    00918         for(int j=i; j<N; ++j, ++b)
    00919         {
    00920             ArrayVector<Kernel1D<KernelType> > kernels(N, gauss);
    00921             if(i == j)
    00922             {
    00923                 kernels[i].initGaussianDerivative(sigma, 2);
    00924             }
    00925             else
    00926             {
    00927                 kernels[i].initGaussianDerivative(sigma, 1);
    00928                 kernels[j].initGaussianDerivative(sigma, 1);
    00929             }
    00930             separableConvolveMultiArray(si, shape, src, di, ElementAccessor(b, dest),
    00931                                         kernels.begin());
    00932         }
    00933     }
    00934 }
    00935
    00936 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00937           class DestIterator, class DestAccessor>
    00938 inline void
    00939 hessianOfGaussianMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    00940                             pair<DestIterator, DestAccessor> const & dest, double sigma )
    00941 {
    00942     hessianOfGaussianMultiArray( source.first, source.second, source.third,
    00943                                  dest.first, dest.second, sigma );
    00944 }
    00945
    00946 namespace detail {
    00947
    00948 template<int N, class VectorType>
    00949 struct StructurTensorFunctor
    00950 {
    00951     typedef VectorType result_type;
    00952     typedef typename VectorType::value_type ValueType;
    00953
    00954     template <class T>
    00955     VectorType operator()(T const & in) const
    00956     {
    00957         VectorType res;
    00958         for(int b=0, i=0; i<N; ++i)
    00959         {
    00960             for(int j=i; j<N; ++j, ++b)
    00961             {
    00962                 res[b] = detail::RequiresExplicitCast<ValueType>::cast(in[i]*in[j]);
    00963             }
    00964         }
    00965         return res;
    00966     }
    00967 };
    00968
    00969 } // namespace detail
    00970
    00971 /********************************************************/
    00972 /*                                                      */
    00973 /*               structureTensorMultiArray              */
    00974 /*                                                      */
    00975 /********************************************************/
    00976
    00977 /** \brief Calculate th structure tensor of a multi-dimensional arrays.
    00978
    00979     This function computes the gradient (outer product) tensor for each element
    00980     of the given N-dimensional array with first-derivative-of-Gaussian filters at
    00981     the given <tt>innerScale</tt>, followed by Gaussian smoothing at <tt>outerScale</tt>.
    00982     Both source and destination arrays are represented by iterators, shape objects and
    00983     accessors. The destination array must have a vector valued pixel type with
    00984     N*(N+1)/2 elements (it represents the upper triangular part of the symmetric
    00985     structure tensor matrix). If the source array is also vector valued, the
    00986     resulting structure tensor is the sum of the individual tensors for each channel.
    00987     This function is implemented by calls to
    00988     \ref separableConvolveMultiArray() with the appropriate kernels.
    00989
    00990     <b> Declarations:</b>
    00991
    00992     pass arguments explicitly:
    00993     \code
    00994     namespace vigra {
    00995         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00996                   class DestIterator, class DestAccessor>
    00997         void
    00998         structureTensorMultiArray(SrcIterator siter, SrcShape const & shape, SrcAccessor src,
    00999                                   DestIterator diter, DestAccessor dest,
    01000                                   double innerScale, double outerScale);
    01001     }
    01002     \endcode
    01003
    01004     use argument objects in conjunction with \ref ArgumentObjectFactories :
    01005     \code
    01006     namespace vigra {
    01007         template <class SrcIterator, class SrcShape, class SrcAccessor,
    01008                   class DestIterator, class DestAccessor>
    01009         void
    01010         structureTensorMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    01011                                   pair<DestIterator, DestAccessor> const & dest,
    01012                                   double innerScale, double outerScale);
    01013     }
    01014     \endcode
    01015
    01016     <b> Usage:</b>
    01017
    01018     <b>\#include</b> <<a href="multi__convolution_8hxx-source.html">vigra/multi_convolution.hxx</a>>
    01019
    01020     \code
    01021     MultiArray<3, RGBValue<float> > source(shape);
    01022     MultiArray<3, TinyVector<float, 6> > dest(shape);
    01023     ...
    01024     // compute structure tensor at scales innerScale and outerScale
    01025     structureTensorMultiArray(srcMultiArrayRange(source), destMultiArray(dest), innerScale, outerScale);
    01026     \endcode
    01027
    01028     <b> Required Interface:</b>
    01029
    01030     see \ref separableConvolveMultiArray(), in addition:
    01031
    01032     \code
    01033     int dimension = 0;
    01034     VectorElementAccessor<DestAccessor> elementAccessor(0, dest);
    01035     \endcode
    01036
    01037     \see separableConvolveMultiArray(), vectorToTensorMultiArray()
    01038 */
    01039 doxygen_overloaded_function(template <...> void structureTensorMultiArray)
    01040
    01041 template <class SrcIterator, class SrcShape, class SrcAccessor,
    01042           class DestIterator, class DestAccessor>
    01043 void
    01044 structureTensorMultiArray(SrcIterator si, SrcShape const & shape, SrcAccessor src,
    01045                           DestIterator di, DestAccessor dest,
    01046                           double innerScale, double outerScale)
    01047 {
    01048     static const int N = SrcShape::static_size;
    01049     static const int M = N*(N+1)/2;
    01050
    01051     typedef typename DestAccessor::value_type DestType;
    01052     typedef typename DestType::value_type     DestValueType;
    01053     typedef typename NumericTraits<DestValueType>::RealPromote KernelType;
    01054     typedef TinyVector<KernelType, N> GradientVector;
    01055     typedef typename AccessorTraits<GradientVector>::default_accessor GradientAccessor;
    01056
    01057     for(int k=0; k<N; ++k)
    01058         if(shape[k] <=0)
    01059             return;
    01060
    01061     vigra_precondition(M == dest.size(di),
    01062         "structureTensorMultiArray(): Wrong number of channels in output array.");
    01063
    01064     vigra_precondition(innerScale > 0.0 && outerScale >= 0.0,
    01065          "structureTensorMultiArray(): Scale must be positive.");
    01066
    01067     MultiArray<N, GradientVector> gradient(shape);
    01068     gaussianGradientMultiArray(si, shape, src,
    01069                                gradient.traverser_begin(), GradientAccessor(),
    01070                                innerScale);
    01071
    01072     transformMultiArray(gradient.traverser_begin(), shape, GradientAccessor(),
    01073                         di, dest,
    01074                         detail::StructurTensorFunctor<N, DestType>());
    01075
    01076     gaussianSmoothMultiArray(di, shape, dest, di, dest, outerScale);
    01077 }
    01078
    01079 template <class SrcIterator, class SrcShape, class SrcAccessor,
    01080           class DestIterator, class DestAccessor>
    01081 inline void
    01082 structureTensorMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> const & source,
    01083                           pair<DestIterator, DestAccessor> const & dest,
    01084                           double innerScale, double outerScale)
    01085 {
    01086     structureTensorMultiArray( source.first, source.second, source.third,
    01087                                dest.first, dest.second, innerScale, outerScale );
    01088 }
    01089
    01090 //@}
    01091
    01092 } //-- namespace vigra
    01093
    01094
    01095 #endif        //-- VIGRA_MULTI_CONVOLUTION_H

© 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.7.1 (3 Dec 2010)