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

vigra/multi_tensorutilities.hxx |

    00001 /************************************************************************/
    00002 /*                                                                      */
    00003 /*               Copyright 2009-2010 by Ullrich Koethe                  */
    00004 /*                                                                      */
    00005 /*    This file is part of the VIGRA computer vision library.           */
    00006 /*    The VIGRA Website is                                              */
    00007 /*        http://hci.iwr.uni-heidelberg.de/vigra/                       */
    00008 /*    Please direct questions, bug reports, and contributions to        */
    00009 /*        ullrich.koethe@iwr.uni-heidelberg.de    or                    */
    00010 /*        vigra@informatik.uni-hamburg.de                               */
    00011 /*                                                                      */
    00012 /*    Permission is hereby granted, free of charge, to any person       */
    00013 /*    obtaining a copy of this software and associated documentation    */
    00014 /*    files (the "Software"), to deal in the Software without           */
    00015 /*    restriction, including without limitation the rights to use,      */
    00016 /*    copy, modify, merge, publish, distribute, sublicense, and/or      */
    00017 /*    sell copies of the Software, and to permit persons to whom the    */
    00018 /*    Software is furnished to do so, subject to the following          */
    00019 /*    conditions:                                                       */
    00020 /*                                                                      */
    00021 /*    The above copyright notice and this permission notice shall be    */
    00022 /*    included in all copies or substantial portions of the             */
    00023 /*    Software.                                                         */
    00024 /*                                                                      */
    00025 /*    THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND    */
    00026 /*    EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES   */
    00027 /*    OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND          */
    00028 /*    NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT       */
    00029 /*    HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,      */
    00030 /*    WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING      */
    00031 /*    FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR     */
    00032 /*    OTHER DEALINGS IN THE SOFTWARE.                                   */
    00033 /*                                                                      */
    00034 /************************************************************************/
    00035
    00036 #ifndef VIGRA_MULTI_TENSORUTILITIES_HXX
    00037 #define VIGRA_MULTI_TENSORUTILITIES_HXX
    00038
    00039 #include <cmath>
    00040 #include "utilities.hxx"
    00041 #include "mathutil.hxx"
    00042 #include "metaprogramming.hxx"
    00043 #include "multi_pointoperators.hxx"
    00044
    00045 namespace vigra {
    00046
    00047 namespace detail {
    00048
    00049 template <int N, class ArgumentVector, class ResultVector>
    00050 class OuterProductFunctor
    00051 {
    00052 public:
    00053     typedef ArgumentVector argument_type;
    00054     typedef ResultVector result_type;
    00055     typedef typename ArgumentVector::value_type ValueType;
    00056
    00057     result_type operator()(argument_type const & in) const
    00058     {
    00059         result_type res;
    00060         for(int b=0, i=0; i<N; ++i)
    00061         {
    00062             for(int j=i; j<N; ++j, ++b)
    00063             {
    00064                 res[b] = detail::RequiresExplicitCast<ValueType>::cast(in[i]*in[j]);
    00065             }
    00066         }
    00067         return res;
    00068     }
    00069 };
    00070
    00071 template <int N, class ArgumentVector>
    00072 class TensorTraceFunctor
    00073 {
    00074 public:
    00075
    00076     typedef ArgumentVector argument_type;
    00077     typedef typename ArgumentVector::value_type result_type;
    00078
    00079     result_type exec(argument_type const & v, MetaInt<1>) const
    00080     {
    00081         return v[0];
    00082     }
    00083
    00084     result_type exec(argument_type const & v, MetaInt<2>) const
    00085     {
    00086         return v[0] + v[2];
    00087     }
    00088
    00089     result_type exec(argument_type const & v, MetaInt<3>) const
    00090     {
    00091         return v[0] + v[3] + v[5];
    00092     }
    00093
    00094     template <int N2>
    00095     void exec(argument_type const & v, result_type & r, MetaInt<N2>) const
    00096     {
    00097         vigra_fail("tensorTraceMultiArray(): Sorry, can only handle dimensions up to 3.");
    00098     }
    00099
    00100     result_type operator()( const argument_type & a ) const
    00101     {
    00102         return exec(a, MetaInt<N>());
    00103     }
    00104 };
    00105
    00106 template <int N, class ArgumentVector, class ResultVector>
    00107 class EigenvaluesFunctor
    00108 {
    00109 public:
    00110
    00111     typedef ArgumentVector argument_type;
    00112     typedef ResultVector result_type;
    00113
    00114     void exec(argument_type const & v, result_type & r, MetaInt<1>) const
    00115     {
    00116         symmetric2x2Eigenvalues(v[0], &r[0]);
    00117     }
    00118
    00119     void exec(argument_type const & v, result_type & r, MetaInt<2>) const
    00120     {
    00121         symmetric2x2Eigenvalues(v[0], v[1], v[2], &r[0], &r[1]);
    00122     }
    00123
    00124     void exec(argument_type const & v, result_type & r, MetaInt<3>) const
    00125     {
    00126         symmetric3x3Eigenvalues(v[0], v[1], v[2], v[3], v[4], v[5], &r[0], &r[1], &r[2]);
    00127     }
    00128
    00129     template <int N2>
    00130     void exec(argument_type const & v, result_type & r, MetaInt<N2>) const
    00131     {
    00132         vigra_fail("tensorEigenvaluesMultiArray(): Sorry, can only handle dimensions up to 3.");
    00133     }
    00134
    00135     result_type operator()( const argument_type & a ) const
    00136     {
    00137         result_type res;
    00138         exec(a, res, MetaInt<N>());
    00139         return res;
    00140     }
    00141 };
    00142
    00143
    00144 template <int N, class ArgumentVector>
    00145 class DeterminantFunctor
    00146 {
    00147 public:
    00148
    00149     typedef ArgumentVector argument_type;
    00150     typedef typename ArgumentVector::value_type result_type;
    00151
    00152     result_type exec(argument_type const & v, MetaInt<1>) const
    00153     {
    00154         return v[0];
    00155     }
    00156
    00157     result_type exec(argument_type const & v, MetaInt<2>) const
    00158     {
    00159         return v[0]*v[2] - sq(v[1]);
    00160     }
    00161
    00162     result_type exec(argument_type const & v, MetaInt<3>) const
    00163     {
    00164         result_type r0, r1, r2;
    00165         symmetric3x3Eigenvalues(v[0], v[1], v[2], v[3], v[4], v[5], &r0, &r1, &r2);
    00166         return r0*r1*r2;
    00167     }
    00168
    00169     template <int N2>
    00170     void exec(argument_type const & v, result_type & r, MetaInt<N2>) const
    00171     {
    00172         vigra_fail("tensorDeterminantMultiArray(): Sorry, can only handle dimensions up to 3.");
    00173     }
    00174
    00175     result_type operator()( const argument_type & a ) const
    00176     {
    00177         return exec(a, MetaInt<N>());
    00178     }
    00179 };
    00180
    00181 } // namespace detail
    00182
    00183
    00184 /** \addtogroup MultiPointoperators
    00185 */
    00186 //@{
    00187
    00188 /********************************************************/
    00189 /*                                                      */
    00190 /*                vectorToTensorMultiArray              */
    00191 /*                                                      */
    00192 /********************************************************/
    00193
    00194 /** \brief Calculate the tensor (outer) product of a N-D vector with itself.
    00195
    00196     This function is useful to transform vector arrays into a tensor representation
    00197     that can be used as input to tensor based processing and analysis functions
    00198     (e.g. tensor smoothing). When the input array has N dimensions, the input value_type
    00199     must be a vector of length N, whereas the output value_type mus be vectors of length
    00200     N*(N-1)/2 which will represent the upper triangular part of the resulting (symmetric)
    00201     tensor. That is, for 2D arrays the output contains the elements
    00202     <tt>[t11, t12 == t21, t22]</tt> in this order, whereas it contains the elements
    00203     <tt>[t11, t12, t13, t22, t23, t33]</tt> for 3D arrays.
    00204
    00205     Currently, <tt>N <= 3</tt> is required.
    00206
    00207     <b> Declarations:</b>
    00208
    00209     pass arguments explicitly:
    00210     \code
    00211     namespace vigra {
    00212         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00213                   class DestIterator, class DestAccessor>
    00214         void
    00215         vectorToTensorMultiArray(SrcIterator  si, SrcShape const & shape, SrcAccessor src,
    00216                                  DestIterator di, DestAccessor dest);
    00217     }
    00218     \endcode
    00219
    00220     use argument objects in conjunction with \ref ArgumentObjectFactories :
    00221     \code
    00222     namespace vigra {
    00223         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00224                   class DestIterator, class DestAccessor>
    00225         void
    00226         vectorToTensorMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> s,
    00227                                  pair<DestIterator, DestAccessor> d);
    00228     }
    00229     \endcode
    00230
    00231     <b> Usage:</b>
    00232
    00233     <b>\#include</b> <vigra/multi_tensorutilities.hxx>
    00234
    00235     \code
    00236     MultiArray<3, float> vol(shape);
    00237     MultiArray<3, TinyVector<float, 3> > gradient(shape);
    00238     MultiArray<3, TinyVector<float, 6> > tensor(shape);
    00239
    00240     gaussianGradientMultiArray(srcMultiArrayRange(vol), destMultiArray(gradient), 2.0);
    00241     vectorToTensorMultiArray(srcMultiArrayRange(gradient), destMultiArray(tensor));
    00242     \endcode
    00243
    00244 */
    00245 doxygen_overloaded_function(template <...> void vectorToTensorMultiArray)
    00246
    00247 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00248           class DestIterator, class DestAccessor>
    00249 void
    00250 vectorToTensorMultiArray(SrcIterator  si, SrcShape const & shape, SrcAccessor src,
    00251                          DestIterator di, DestAccessor dest)
    00252 {
    00253     static const int N = SrcShape::static_size;
    00254     static const int M = N*(N+1)/2;
    00255
    00256     typedef typename SrcAccessor::value_type  SrcType;
    00257     typedef typename DestAccessor::value_type DestType;
    00258
    00259     for(int k=0; k<N; ++k)
    00260         if(shape[k] <=0)
    00261             return;
    00262
    00263     vigra_precondition(N == (int)src.size(si),
    00264         "vectorToTensorMultiArray(): Wrong number of channels in input array.");
    00265     vigra_precondition(M == (int)dest.size(di),
    00266         "vectorToTensorMultiArray(): Wrong number of channels in output array.");
    00267
    00268     transformMultiArray(si, shape, src, di, dest,
    00269                         detail::OuterProductFunctor<N, SrcType, DestType>());
    00270 }
    00271
    00272 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00273           class DestIterator, class DestAccessor>
    00274 inline void
    00275 vectorToTensorMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> s,
    00276                          pair<DestIterator, DestAccessor> d)
    00277 {
    00278     vectorToTensorMultiArray(s.first, s.second, s.third, d.first, d.second);
    00279 }
    00280
    00281 /********************************************************/
    00282 /*                                                      */
    00283 /*                tensorTraceMultiArray                 */
    00284 /*                                                      */
    00285 /********************************************************/
    00286
    00287 /** \brief Calculate the tensor trace for every element of a N-D tensor array.
    00288
    00289     This function turns a N-D tensor (whose value_type is a vector of length N*(N+1)/2,
    00290     see \ref vectorToTensorMultiArray()) representing the upper triangular part of a
    00291     symmetric tensor into a scalar array holding the tensor trace.
    00292
    00293     Currently, <tt>N <= 3</tt> is required.
    00294
    00295     <b> Declarations:</b>
    00296
    00297     pass arguments explicitly:
    00298     \code
    00299     namespace vigra {
    00300         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00301                   class DestIterator, class DestAccessor>
    00302         void
    00303         tensorTraceMultiArray(SrcIterator si,  SrcShape const & shape, SrcAccessor src,
    00304                               DestIterator di, DestAccessor dest);
    00305     }
    00306     \endcode
    00307
    00308
    00309     use argument objects in conjunction with \ref ArgumentObjectFactories :
    00310     \code
    00311     namespace vigra {
    00312         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00313                   class DestIterator, class DestAccessor>
    00314         void
    00315         tensorTraceMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> s,
    00316                               pair<DestIterator, DestAccessor> d);
    00317     }
    00318     \endcode
    00319
    00320     <b> Usage:</b>
    00321
    00322     <b>\#include</b> <vigra/multi_tensorutilities.hxx>
    00323
    00324     \code
    00325     MultiArray<3, float> vol(shape);
    00326     MultiArray<3, TinyVector<float, 6> > hessian(shape);
    00327     MultiArray<3, float> trace(shape);
    00328
    00329     hessianOfGaussianMultiArray(srcMultiArrayRange(vol), destMultiArray(hessian), 2.0);
    00330     tensorTraceMultiArray(srcMultiArrayRange(hessian), destMultiArray(trace));
    00331     \endcode
    00332
    00333 */
    00334 doxygen_overloaded_function(template <...> void tensorTraceMultiArray)
    00335
    00336 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00337           class DestIterator, class DestAccessor>
    00338 void
    00339 tensorTraceMultiArray(SrcIterator si,  SrcShape const & shape, SrcAccessor src,
    00340                       DestIterator di, DestAccessor dest)
    00341 {
    00342     static const int N = SrcShape::static_size;
    00343     typedef typename SrcAccessor::value_type  SrcType;
    00344
    00345     transformMultiArray(si, shape, src, di, dest,
    00346                         detail::TensorTraceFunctor<N, SrcType>());
    00347 }
    00348
    00349 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00350           class DestIterator, class DestAccessor>
    00351 inline void
    00352 tensorTraceMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> s,
    00353                       pair<DestIterator, DestAccessor> d)
    00354 {
    00355     tensorTraceMultiArray(s.first, s.second, s.third, d.first, d.second);
    00356 }
    00357
    00358
    00359 /********************************************************/
    00360 /*                                                      */
    00361 /*             tensorEigenvaluesMultiArray              */
    00362 /*                                                      */
    00363 /********************************************************/
    00364
    00365 /** \brief Calculate the tensor eigenvalues for every element of a N-D tensor array.
    00366
    00367     This function turns a N-D tensor (whose value_type is a vector of length N*(N+1)/2,
    00368     see \ref vectorToTensorMultiArray()) representing the upper triangular part of a
    00369     symmetric tensor into a vector-valued array holding the tensor eigenvalues (thus,
    00370     the destination value_type must be vectors of length N).
    00371
    00372     Currently, <tt>N <= 3</tt> is required.
    00373
    00374     <b> Declarations:</b>
    00375
    00376     pass arguments explicitly:
    00377     \code
    00378     namespace vigra {
    00379         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00380                   class DestIterator, class DestAccessor>
    00381         void
    00382         tensorEigenvaluesMultiArray(SrcIterator si,  SrcShape const & shape, SrcAccessor src,
    00383                                     DestIterator di, DestAccessor dest);
    00384     }
    00385     \endcode
    00386
    00387
    00388     use argument objects in conjunction with \ref ArgumentObjectFactories :
    00389     \code
    00390     namespace vigra {
    00391         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00392                   class DestIterator, class DestAccessor>
    00393         void
    00394         tensorEigenvaluesMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> s,
    00395                                     pair<DestIterator, DestAccessor> d);
    00396     }
    00397     \endcode
    00398
    00399     <b> Usage:</b>
    00400
    00401     <b>\#include</b> <vigra/multi_tensorutilities.hxx>
    00402
    00403     \code
    00404     MultiArray<3, float> vol(shape);
    00405     MultiArray<3, TinyVector<float, 6> > hessian(shape);
    00406     MultiArray<3, TinyVector<float, 3> > eigenvalues(shape);
    00407
    00408     hessianOfGaussianMultiArray(srcMultiArrayRange(vol), destMultiArray(hessian), 2.0);
    00409     tensorEigenvaluesMultiArray(srcMultiArrayRange(hessian), destMultiArray(eigenvalues));
    00410     \endcode
    00411
    00412 */
    00413 doxygen_overloaded_function(template <...> void tensorEigenvaluesMultiArray)
    00414
    00415 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00416           class DestIterator, class DestAccessor>
    00417 void
    00418 tensorEigenvaluesMultiArray(SrcIterator si,  SrcShape const & shape, SrcAccessor src,
    00419                             DestIterator di, DestAccessor dest)
    00420 {
    00421     static const int N = SrcShape::static_size;
    00422     static const int M = N*(N+1)/2;
    00423
    00424     typedef typename SrcAccessor::value_type  SrcType;
    00425     typedef typename DestAccessor::value_type DestType;
    00426
    00427     for(int k=0; k<N; ++k)
    00428         if(shape[k] <=0)
    00429             return;
    00430
    00431     vigra_precondition(M == (int)src.size(si),
    00432         "tensorEigenvaluesMultiArray(): Wrong number of channels in input array.");
    00433     vigra_precondition(N == (int)dest.size(di),
    00434         "tensorEigenvaluesMultiArray(): Wrong number of channels in output array.");
    00435
    00436     transformMultiArray(si, shape, src, di, dest,
    00437                         detail::EigenvaluesFunctor<N, SrcType, DestType>());
    00438 }
    00439
    00440 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00441           class DestIterator, class DestAccessor>
    00442 inline void
    00443 tensorEigenvaluesMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> s,
    00444                             pair<DestIterator, DestAccessor> d)
    00445 {
    00446     tensorEigenvaluesMultiArray(s.first, s.second, s.third, d.first, d.second);
    00447 }
    00448
    00449 /********************************************************/
    00450 /*                                                      */
    00451 /*             tensorDeterminantMultiArray              */
    00452 /*                                                      */
    00453 /********************************************************/
    00454
    00455 /** \brief Calculate the tensor determinant for every element of a ND tensor array.
    00456
    00457     This function turns a N-D tensor (whose value_type is a vector of length N*(N+1)/2,
    00458     see \ref vectorToTensorMultiArray()) representing the upper triangular part of a
    00459     symmetric tensor into the a scalar array holding the tensor determinant.
    00460
    00461     Currently, <tt>N <= 3</tt> is required.
    00462
    00463     <b> Declarations:</b>
    00464
    00465     pass arguments explicitly:
    00466     \code
    00467     namespace vigra {
    00468         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00469                   class DestIterator, class DestAccessor>
    00470         void
    00471         tensorDeterminantMultiArray(SrcIterator si,  SrcShape const & shape, SrcAccessor src,
    00472                                     DestIterator di, DestAccessor dest);
    00473     }
    00474     \endcode
    00475
    00476
    00477     use argument objects in conjunction with \ref ArgumentObjectFactories :
    00478     \code
    00479     namespace vigra {
    00480         template <class SrcIterator, class SrcShape, class SrcAccessor,
    00481                   class DestIterator, class DestAccessor>
    00482         void
    00483         tensorDeterminantMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> s,
    00484                                     pair<DestIterator, DestAccessor> d);
    00485     }
    00486     \endcode
    00487
    00488     <b> Usage:</b>
    00489
    00490     <b>\#include</b> <vigra/multi_tensorutilities.hxx>
    00491
    00492     \code
    00493     MultiArray<3, float> vol(shape);
    00494     MultiArray<3, TinyVector<float, 6> > hessian(shape);
    00495     MultiArray<3, float> determinant(shape);
    00496
    00497     hessianOfGaussianMultiArray(srcMultiArrayRange(vol), destMultiArray(hessian), 2.0);
    00498     tensorDeterminantMultiArray(srcMultiArrayRange(hessian), destMultiArray(determinant));
    00499     \endcode
    00500
    00501 */
    00502 doxygen_overloaded_function(template <...> void tensorDeterminantMultiArray)
    00503
    00504 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00505           class DestIterator, class DestAccessor>
    00506 void
    00507 tensorDeterminantMultiArray(SrcIterator si,  SrcShape const & shape, SrcAccessor src,
    00508                             DestIterator di, DestAccessor dest)
    00509 {
    00510     typedef typename SrcAccessor::value_type  SrcType;
    00511
    00512     static const int N = SrcShape::static_size;
    00513     static const int M = N*(N+1)/2;
    00514
    00515     for(int k=0; k<N; ++k)
    00516         if(shape[k] <=0)
    00517             return;
    00518
    00519     vigra_precondition(M == (int)src.size(si),
    00520         "tensorDeterminantMultiArray(): Wrong number of channels in output array.");
    00521
    00522     transformMultiArray(si, shape, src, di, dest,
    00523                         detail::DeterminantFunctor<N, SrcType>());
    00524 }
    00525
    00526 template <class SrcIterator, class SrcShape, class SrcAccessor,
    00527           class DestIterator, class DestAccessor>
    00528 inline void
    00529 tensorDeterminantMultiArray(triple<SrcIterator, SrcShape, SrcAccessor> s,
    00530                             pair<DestIterator, DestAccessor> d)
    00531 {
    00532     tensorDeterminantMultiArray(s.first, s.second, s.third, d.first, d.second);
    00533 }
    00534
    00535 //@}
    00536
    00537 } // namespace vigra
    00538
    00539 #endif /* VIGRA_MULTI_TENSORUTILITIES_HXX */

© 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.9.0 (Tue Nov 6 2012)