[ VIGRA Homepage | Function Index | Class Index | Namespaces | File List | Main Page ]
00001 /************************************************************************/ 00002 /* */ 00003 /* Copyright 1998-2002 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 00037 #ifndef VIGRA_CONVOLUTION_HXX 00038 #define VIGRA_CONVOLUTION_HXX 00039 00040 #include <functional> 00041 #include "stdconvolution.hxx" 00042 #include "separableconvolution.hxx" 00043 #include "recursiveconvolution.hxx" 00044 #include "nonlineardiffusion.hxx" 00045 #include "combineimages.hxx" 00046 00047 /** \page Convolution Functions to Convolve Images and Signals 00048 00049 1D and 2D filters, including separable and recursive convolution, and non-linear diffusion 00050 00051 <b>\#include</b> <<a href="convolution_8hxx-source.html">vigra/convolution.hxx</a>><br> 00052 Namespace: vigra 00053 00054 <UL style="list-style-image:url(documents/bullet.gif)"> 00055 <LI> \ref CommonConvolutionFilters 00056 <BR> <em>Short-hands for the most common 2D convolution filters</em> 00057 <LI> \ref MultiArrayConvolutionFilters 00058 <BR> <em>Convolution filters for arbitrary dimensional arrays (MultiArray etc.)</em> 00059 <LI> \ref ResamplingConvolutionFilters 00060 <BR> <em>Resampling convolution filters</em> 00061 <LI> \ref StandardConvolution 00062 <BR> <em>2D non-separable convolution, with and without ROI mask </em> 00063 <LI> \ref vigra::Kernel2D 00064 <BR> <em>Generic 2-dimensional discrete convolution kernel </em> 00065 <LI> \ref SeparableConvolution 00066 <BR> <em>1D convolution and separable filters in 2 dimensions </em> 00067 <LI> \ref vigra::Kernel1D 00068 <BR> <em>Generic 1-dimensional discrete convolution kernel </em> 00069 <LI> \ref RecursiveConvolution 00070 <BR> <em>Recursive filters (1st and 2nd order)</em> 00071 <LI> \ref NonLinearDiffusion 00072 <BR> <em>Edge-preserving smoothing </em> 00073 <LI> \ref BorderTreatmentMode 00074 <BR> <em>Choose between different border treatment modes </em> 00075 <LI> \ref KernelArgumentObjectFactories 00076 <BR> <em>Factory functions to create argument objects to simplify passing kernels</em> 00077 </UL> 00078 */ 00079 00080 /** \page KernelArgumentObjectFactories Kernel Argument Object Factories 00081 00082 These factory functions allow to create argument objects for 1D 00083 and 2D convolution kernel analogously to 00084 \ref ArgumentObjectFactories for images. 00085 00086 \section Kernel1dFactory kernel1d() 00087 00088 Pass a \ref vigra::Kernel1D to a 1D or separable convolution algorithm. 00089 00090 These factories can be used to create argument objects when we 00091 are given instances or subclasses of \ref vigra::Kernel1D 00092 (analogous to the \ref ArgumentObjectFactories for images). 00093 These factory functions access <TT>kernel.center()</TT>, 00094 <TT>kernel.left()</TT>, <TT>kernel.right()</TT>, <TT>kernel.accessor()</TT>, 00095 and <TT>kernel.borderTreatment()</TT> to obtain the necessary 00096 information. The following factory functions are provided: 00097 00098 <table> 00099 <tr><th bgcolor="#f0e0c0" colspan=2 align=left> 00100 <TT>\ref vigra::Kernel1D "vigra::Kernel1D<SomeType>" kernel;</TT> 00101 </th> 00102 </tr> 00103 <tr><td> 00104 <TT>kernel1d(kernel)</TT> 00105 </td><td> 00106 create argument object from information provided by 00107 kernel 00108 00109 </td></tr> 00110 <tr><td> 00111 <TT>kernel1d(kernel, vigra::BORDER_TREATMENT_CLIP)</TT> 00112 </td><td> 00113 create argument object from information provided by 00114 kernel, but use given border treatment mode 00115 00116 </td></tr> 00117 <tr><td> 00118 <TT>kernel1d(kerneliterator, kernelaccessor,</TT><br> 00119 <TT> kernelleft, kernelright,</TT><br> 00120 <TT> vigra::BORDER_TREATMENT_CLIP)</TT> 00121 </td><td> 00122 create argument object from explicitly given iterator 00123 (pointing to the center of th kernel), accessor, 00124 left and right boundaries, and border treatment mode 00125 00126 </table> 00127 00128 For usage examples see 00129 \ref SeparableConvolution "one-dimensional and separable convolution functions". 00130 00131 \section Kernel2dFactory kernel2d() 00132 00133 Pass a \ref vigra::Kernel2D to a 2D (non-separable) convolution algorithm. 00134 00135 These factories can be used to create argument objects when we 00136 are given instances or subclasses of \ref vigra::Kernel2D 00137 (analogous to the \ref ArgumentObjectFactories for images). 00138 These factory functions access <TT>kernel.center()</TT>, 00139 <TT>kernel.upperLeft()</TT>, <TT>kernel.lowerRight()</TT>, <TT>kernel.accessor()</TT>, 00140 and <TT>kernel.borderTreatment()</TT> to obtain the necessary 00141 information. The following factory functions are provided: 00142 00143 <table> 00144 <tr><th bgcolor="#f0e0c0" colspan=2 align=left> 00145 <TT>\ref vigra::Kernel2D "vigra::Kernel2D<SomeType>" kernel;</TT> 00146 </th> 00147 </tr> 00148 <tr><td> 00149 <TT>kernel2d(kernel)</TT> 00150 </td><td> 00151 create argument object from information provided by 00152 kernel 00153 00154 </td></tr> 00155 <tr><td> 00156 <TT>kernel2d(kernel, vigra::BORDER_TREATMENT_CLIP)</TT> 00157 </td><td> 00158 create argument object from information provided by 00159 kernel, but use given border treatment mode 00160 00161 </td></tr> 00162 <tr><td> 00163 <TT>kernel2d(kerneliterator, kernelaccessor,</TT> 00164 <TT> upperleft, lowerright,</TT> 00165 <TT> vigra::BORDER_TREATMENT_CLIP)</TT> 00166 </td><td> 00167 create argument object from explicitly given iterator 00168 (pointing to the center of th kernel), accessor, 00169 upper left and lower right corners, and border treatment mode 00170 00171 </table> 00172 00173 For usage examples see \ref StandardConvolution "two-dimensional convolution functions". 00174 */ 00175 00176 namespace vigra { 00177 00178 00179 00180 /********************************************************/ 00181 /* */ 00182 /* Common convolution filters */ 00183 /* */ 00184 /********************************************************/ 00185 00186 /** \addtogroup CommonConvolutionFilters Common Filters 00187 00188 These functions calculate common filters by appropriate sequences of calls 00189 to \ref separableConvolveX() and \ref separableConvolveY(). 00190 */ 00191 //@{ 00192 00193 /********************************************************/ 00194 /* */ 00195 /* convolveImage */ 00196 /* */ 00197 /********************************************************/ 00198 00199 /** \brief Apply two separable filters successively, the first in x-direction, 00200 the second in y-direction. 00201 00202 This function is a shorthand for the concatenation of a call to 00203 \ref separableConvolveX() and \ref separableConvolveY() 00204 with the given kernels. 00205 00206 <b> Declarations:</b> 00207 00208 pass arguments explicitly: 00209 \code 00210 namespace vigra { 00211 template <class SrcIterator, class SrcAccessor, 00212 class DestIterator, class DestAccessor, 00213 class T> 00214 void convolveImage(SrcIterator supperleft, 00215 SrcIterator slowerright, SrcAccessor sa, 00216 DestIterator dupperleft, DestAccessor da, 00217 Kernel1D<T> const & kx, Kernel1D<T> const & ky); 00218 } 00219 \endcode 00220 00221 00222 use argument objects in conjunction with \ref ArgumentObjectFactories : 00223 \code 00224 namespace vigra { 00225 template <class SrcIterator, class SrcAccessor, 00226 class DestIterator, class DestAccessor, 00227 class T> 00228 inline void 00229 convolveImage(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00230 pair<DestIterator, DestAccessor> dest, 00231 Kernel1D<T> const & kx, Kernel1D<T> const & ky); 00232 } 00233 \endcode 00234 00235 <b> Usage:</b> 00236 00237 <b>\#include</b> <<a href="convolution_8hxx-source.html">vigra/convolution.hxx</a>> 00238 00239 00240 \code 00241 vigra::FImage src(w,h), dest(w,h); 00242 ... 00243 00244 // implement sobel filter in x-direction 00245 Kernel1D<double> kx, ky; 00246 kx.initSymmetricGradient(); 00247 ky.initBinomial(1); 00248 00249 vigra::convolveImage(srcImageRange(src), destImage(dest), kx, ky); 00250 00251 \endcode 00252 00253 */ 00254 template <class SrcIterator, class SrcAccessor, 00255 class DestIterator, class DestAccessor, 00256 class T> 00257 void convolveImage(SrcIterator supperleft, 00258 SrcIterator slowerright, SrcAccessor sa, 00259 DestIterator dupperleft, DestAccessor da, 00260 Kernel1D<T> const & kx, Kernel1D<T> const & ky) 00261 { 00262 typedef typename 00263 NumericTraits<typename SrcAccessor::value_type>::RealPromote 00264 TmpType; 00265 BasicImage<TmpType> tmp(slowerright - supperleft, SkipInitialization); 00266 00267 separableConvolveX(srcIterRange(supperleft, slowerright, sa), 00268 destImage(tmp), kernel1d(kx)); 00269 separableConvolveY(srcImageRange(tmp), 00270 destIter(dupperleft, da), kernel1d(ky)); 00271 } 00272 00273 template <class SrcIterator, class SrcAccessor, 00274 class DestIterator, class DestAccessor, 00275 class T> 00276 inline void 00277 convolveImage(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00278 pair<DestIterator, DestAccessor> dest, 00279 Kernel1D<T> const & kx, Kernel1D<T> const & ky) 00280 { 00281 convolveImage(src.first, src.second, src.third, 00282 dest.first, dest.second, kx, ky); 00283 } 00284 00285 /********************************************************/ 00286 /* */ 00287 /* simpleSharpening */ 00288 /* */ 00289 /********************************************************/ 00290 00291 /** \brief Perform simple sharpening function. 00292 00293 This function uses \ref convolveImage() with the following filter: 00294 00295 \code 00296 -sharpening_factor/16.0, -sharpening_factor/8.0, -sharpening_factor/16.0, 00297 -sharpening_factor/8.0, 1.0+sharpening_factor*0.75, -sharpening_factor/8.0, 00298 -sharpening_factor/16.0, -sharpening_factor/8.0, -sharpening_factor/16.0; 00299 \endcode 00300 00301 and uses <TT>BORDER_TREATMENT_REFLECT</TT> as border treatment mode. 00302 00303 <b> Preconditions:</b> 00304 \code 00305 1. sharpening_factor >= 0 00306 2. scale >= 0 00307 \endcode 00308 00309 <b> Declarations:</b> 00310 00311 <b> Declarations:</b> 00312 00313 pass arguments explicitly: 00314 \code 00315 namespace vigra { 00316 template <class SrcIterator, class SrcAccessor, 00317 class DestIterator, class DestAccessor> 00318 void simpleSharpening(SrcIterator src_ul, SrcIterator src_lr, SrcAccessor src_acc, 00319 DestIterator dest_ul, DestAccessor dest_acc, double sharpening_factor) 00320 00321 } 00322 \endcode 00323 00324 00325 use argument objects in conjunction with \ref ArgumentObjectFactories : 00326 \code 00327 namespace vigra { 00328 template <class SrcIterator, class SrcAccessor, 00329 class DestIterator, class DestAccessor> 00330 inline 00331 void simpleSharpening(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00332 pair<DestIterator, DestAccessor> dest, double sharpening_factor) 00333 { 00334 simpleSharpening(src.first, src.second, src.third, 00335 dest.first, dest.second, sharpening_factor); 00336 } 00337 00338 } 00339 \endcode 00340 00341 <b> Usage:</b> 00342 00343 <b>\#include</b> <<a href="convolution_8hxx-source.html">vigra/convolution.hxx</a>> 00344 00345 00346 \code 00347 vigra::FImage src(w,h), dest(w,h); 00348 ... 00349 00350 // sharpening with sharpening_factor = 0.1 00351 vigra::simpleSharpening(srcImageRange(src), destImage(dest), 0.1); 00352 00353 \endcode 00354 00355 */ 00356 doxygen_overloaded_function(template <...> void simpleSharpening) 00357 00358 template <class SrcIterator, class SrcAccessor, 00359 class DestIterator, class DestAccessor> 00360 void simpleSharpening(SrcIterator src_ul, SrcIterator src_lr, SrcAccessor src_acc, 00361 DestIterator dest_ul, DestAccessor dest_acc, double sharpening_factor) 00362 { 00363 00364 vigra_precondition(sharpening_factor >= 0.0, 00365 "simpleSharpening(): amount of sharpening must be >= 0."); 00366 00367 Kernel2D<double> kernel; 00368 00369 kernel.initExplicitly(Diff2D(-1,-1), Diff2D(1,1)) = -sharpening_factor/16.0, -sharpening_factor/8.0, -sharpening_factor/16.0, 00370 -sharpening_factor/8.0, 1.0+sharpening_factor*0.75, -sharpening_factor/8.0, 00371 -sharpening_factor/16.0, -sharpening_factor/8.0, -sharpening_factor/16.0; 00372 00373 convolveImage(src_ul, src_lr, src_acc, dest_ul, dest_acc, 00374 kernel.center(), kernel.accessor(), 00375 kernel.upperLeft(), kernel.lowerRight() , BORDER_TREATMENT_REFLECT ); 00376 } 00377 00378 template <class SrcIterator, class SrcAccessor, 00379 class DestIterator, class DestAccessor> 00380 inline 00381 void simpleSharpening(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00382 pair<DestIterator, DestAccessor> dest, double sharpening_factor) 00383 { 00384 simpleSharpening(src.first, src.second, src.third, 00385 dest.first, dest.second, sharpening_factor); 00386 } 00387 00388 00389 /********************************************************/ 00390 /* */ 00391 /* gaussianSharpening */ 00392 /* */ 00393 /********************************************************/ 00394 00395 /** \brief Perform sharpening function with gaussian filter. 00396 00397 00398 This function uses \ref gaussianSmoothing() at the given scale to create a 00399 temporary image 'smooth' and than blends the original and smoothed image 00400 according to the formula 00401 00402 \code 00403 dest = (1 + sharpening_factor)*src - sharpening_factor*smooth 00404 \endcode 00405 00406 <b> Preconditions:</b> 00407 \code 00408 1. sharpening_factor >= 0 00409 2. scale >= 0 00410 \endcode 00411 00412 <b> Declarations:</b> 00413 00414 pass arguments explicitly: 00415 \code 00416 namespace vigra { 00417 template <class SrcIterator, class SrcAccessor, 00418 class DestIterator, class DestAccessor> 00419 void gaussianSharpening(SrcIterator src_ul, SrcIterator src_lr, SrcAccessor src_acc, 00420 DestIterator dest_ul, DestAccessor dest_acc, 00421 double sharpening_factor, double scale) 00422 } 00423 \endcode 00424 00425 00426 use argument objects in conjunction with \ref ArgumentObjectFactories : 00427 \code 00428 namespace vigra { 00429 template <class SrcIterator, class SrcAccessor, 00430 class DestIterator, class DestAccessor> 00431 void gaussianSharpening(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00432 pair<DestIterator, DestAccessor> dest, 00433 double sharpening_factor, double scale) 00434 } 00435 \endcode 00436 00437 <b> Usage:</b> 00438 00439 <b>\#include</b> <<a href="convolution_8hxx-source.html">vigra/convolution.hxx</a>> 00440 00441 00442 \code 00443 vigra::FImage src(w,h), dest(w,h); 00444 ... 00445 00446 // sharpening with sharpening_factor = 3.0 00447 // smoothing with scale = 0.5 00448 vigra::gaussianSmoothing(srcImageRange(src), destImage(dest), 3.0, 0.5); 00449 00450 \endcode 00451 00452 */ 00453 doxygen_overloaded_function(template <...> void gaussianSharpening) 00454 00455 template <class SrcIterator, class SrcAccessor, 00456 class DestIterator, class DestAccessor> 00457 void gaussianSharpening(SrcIterator src_ul, SrcIterator src_lr, SrcAccessor src_acc, 00458 DestIterator dest_ul, DestAccessor dest_acc, double sharpening_factor, 00459 double scale) 00460 { 00461 vigra_precondition(sharpening_factor >= 0.0, 00462 "gaussianSharpening(): amount of sharpening must be >= 0"); 00463 vigra_precondition(scale >= 0.0, 00464 "gaussianSharpening(): scale parameter should be >= 0."); 00465 00466 typedef typename NumericTraits<typename SrcAccessor::value_type>::RealPromote ValueType; 00467 00468 BasicImage<ValueType> tmp(src_lr - src_ul, SkipInitialization); 00469 00470 gaussianSmoothing(src_ul, src_lr, src_acc, tmp.upperLeft(), tmp.accessor(), scale); 00471 00472 SrcIterator i_src = src_ul; 00473 DestIterator i_dest = dest_ul; 00474 typename BasicImage<ValueType>::traverser tmp_ul = tmp.upperLeft(); 00475 typename BasicImage<ValueType>::traverser i_tmp = tmp_ul; 00476 typename BasicImage<ValueType>::Accessor tmp_acc = tmp.accessor(); 00477 00478 for(; i_src.y != src_lr.y ; i_src.y++, i_dest.y++, i_tmp.y++ ) 00479 { 00480 for (;i_src.x != src_lr.x ; i_src.x++, i_dest.x++, i_tmp.x++ ) 00481 { 00482 dest_acc.set((1.0 + sharpening_factor)*src_acc(i_src) - sharpening_factor*tmp_acc(i_tmp), i_dest); 00483 } 00484 i_src.x = src_ul.x; 00485 i_dest.x = dest_ul.x; 00486 i_tmp.x = tmp_ul.x; 00487 } 00488 } 00489 00490 template <class SrcIterator, class SrcAccessor, 00491 class DestIterator, class DestAccessor> 00492 void gaussianSharpening(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00493 pair<DestIterator, DestAccessor> dest, double sharpening_factor, 00494 double scale) 00495 { 00496 gaussianSharpening(src.first, src.second, src.third, 00497 dest.first, dest.second, 00498 sharpening_factor, scale); 00499 } 00500 00501 00502 00503 /********************************************************/ 00504 /* */ 00505 /* gaussianSmoothing */ 00506 /* */ 00507 /********************************************************/ 00508 00509 /** \brief Perform isotropic Gaussian convolution. 00510 00511 This function is a shorthand for the concatenation of a call to 00512 \ref separableConvolveX() and \ref separableConvolveY() with a 00513 Gaussian kernel of the given scale. The function uses 00514 <TT>BORDER_TREATMENT_REFLECT</TT>. 00515 00516 <b> Declarations:</b> 00517 00518 pass arguments explicitly: 00519 \code 00520 namespace vigra { 00521 template <class SrcIterator, class SrcAccessor, 00522 class DestIterator, class DestAccessor> 00523 void gaussianSmoothing(SrcIterator supperleft, 00524 SrcIterator slowerright, SrcAccessor sa, 00525 DestIterator dupperleft, DestAccessor da, 00526 double scale); 00527 } 00528 \endcode 00529 00530 00531 use argument objects in conjunction with \ref ArgumentObjectFactories : 00532 \code 00533 namespace vigra { 00534 template <class SrcIterator, class SrcAccessor, 00535 class DestIterator, class DestAccessor> 00536 inline void 00537 gaussianSmoothing(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00538 pair<DestIterator, DestAccessor> dest, 00539 double scale); 00540 } 00541 \endcode 00542 00543 <b> Usage:</b> 00544 00545 <b>\#include</b> <<a href="convolution_8hxx-source.html">vigra/convolution.hxx</a>> 00546 00547 00548 \code 00549 vigra::FImage src(w,h), dest(w,h); 00550 ... 00551 00552 // smooth with scale = 3.0 00553 vigra::gaussianSmoothing(srcImageRange(src), destImage(dest), 3.0); 00554 00555 \endcode 00556 00557 */ 00558 doxygen_overloaded_function(template <...> void gaussianSmoothing) 00559 00560 template <class SrcIterator, class SrcAccessor, 00561 class DestIterator, class DestAccessor> 00562 void gaussianSmoothing(SrcIterator supperleft, 00563 SrcIterator slowerright, SrcAccessor sa, 00564 DestIterator dupperleft, DestAccessor da, 00565 double scale) 00566 { 00567 typedef typename 00568 NumericTraits<typename SrcAccessor::value_type>::RealPromote 00569 TmpType; 00570 BasicImage<TmpType> tmp(slowerright - supperleft, SkipInitialization); 00571 00572 Kernel1D<double> smooth; 00573 smooth.initGaussian(scale); 00574 smooth.setBorderTreatment(BORDER_TREATMENT_REFLECT); 00575 00576 separableConvolveX(srcIterRange(supperleft, slowerright, sa), 00577 destImage(tmp), kernel1d(smooth)); 00578 separableConvolveY(srcImageRange(tmp), 00579 destIter(dupperleft, da), kernel1d(smooth)); 00580 } 00581 00582 template <class SrcIterator, class SrcAccessor, 00583 class DestIterator, class DestAccessor> 00584 inline void 00585 gaussianSmoothing(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00586 pair<DestIterator, DestAccessor> dest, 00587 double scale) 00588 { 00589 gaussianSmoothing(src.first, src.second, src.third, 00590 dest.first, dest.second, scale); 00591 } 00592 00593 /********************************************************/ 00594 /* */ 00595 /* gaussianGradient */ 00596 /* */ 00597 /********************************************************/ 00598 00599 /** \brief Calculate the gradient vector by means of a 1st derivatives of 00600 Gaussian filter. 00601 00602 This function is a shorthand for the concatenation of a call to 00603 \ref separableConvolveX() and \ref separableConvolveY() with the 00604 appropriate kernels at the given scale. Note that this function can either produce 00605 two separate result images for the x- and y-components of the gradient, or write 00606 into a vector valued image (with at least two components). 00607 00608 <b> Declarations:</b> 00609 00610 pass arguments explicitly: 00611 \code 00612 namespace vigra { 00613 // write x and y component of the gradient into separate images 00614 template <class SrcIterator, class SrcAccessor, 00615 class DestIteratorX, class DestAccessorX, 00616 class DestIteratorY, class DestAccessorY> 00617 void gaussianGradient(SrcIterator supperleft, 00618 SrcIterator slowerright, SrcAccessor sa, 00619 DestIteratorX dupperleftx, DestAccessorX dax, 00620 DestIteratorY dupperlefty, DestAccessorY day, 00621 double scale); 00622 00623 // write x and y component of the gradient into a vector-valued image 00624 template <class SrcIterator, class SrcAccessor, 00625 class DestIterator, class DestAccessor> 00626 void gaussianGradient(SrcIterator supperleft, 00627 SrcIterator slowerright, SrcAccessor src, 00628 DestIterator dupperleft, DestAccessor dest, 00629 double scale); 00630 } 00631 \endcode 00632 00633 00634 use argument objects in conjunction with \ref ArgumentObjectFactories : 00635 \code 00636 namespace vigra { 00637 // write x and y component of the gradient into separate images 00638 template <class SrcIterator, class SrcAccessor, 00639 class DestIteratorX, class DestAccessorX, 00640 class DestIteratorY, class DestAccessorY> 00641 void 00642 gaussianGradient(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00643 pair<DestIteratorX, DestAccessorX> destx, 00644 pair<DestIteratorY, DestAccessorY> desty, 00645 double scale); 00646 00647 // write x and y component of the gradient into a vector-valued image 00648 template <class SrcIterator, class SrcAccessor, 00649 class DestIterator, class DestAccessor> 00650 void 00651 gaussianGradient(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00652 pair<DestIterator, DestAccessor> dest, 00653 double scale); 00654 } 00655 \endcode 00656 00657 <b> Usage:</b> 00658 00659 <b>\#include</b> <<a href="convolution_8hxx-source.html">vigra/convolution.hxx</a>> 00660 00661 00662 \code 00663 vigra::FImage src(w,h), gradx(w,h), grady(w,h); 00664 ... 00665 00666 // calculate gradient vector at scale = 3.0 00667 vigra::gaussianGradient(srcImageRange(src), 00668 destImage(gradx), destImage(grady), 3.0); 00669 00670 \endcode 00671 00672 */ 00673 doxygen_overloaded_function(template <...> void gaussianGradient) 00674 00675 template <class SrcIterator, class SrcAccessor, 00676 class DestIteratorX, class DestAccessorX, 00677 class DestIteratorY, class DestAccessorY> 00678 void gaussianGradient(SrcIterator supperleft, 00679 SrcIterator slowerright, SrcAccessor sa, 00680 DestIteratorX dupperleftx, DestAccessorX dax, 00681 DestIteratorY dupperlefty, DestAccessorY day, 00682 double scale) 00683 { 00684 typedef typename 00685 NumericTraits<typename SrcAccessor::value_type>::RealPromote 00686 TmpType; 00687 BasicImage<TmpType> tmp(slowerright - supperleft, SkipInitialization); 00688 00689 Kernel1D<double> smooth, grad; 00690 smooth.initGaussian(scale); 00691 grad.initGaussianDerivative(scale, 1); 00692 00693 separableConvolveX(srcIterRange(supperleft, slowerright, sa), 00694 destImage(tmp), kernel1d(grad)); 00695 separableConvolveY(srcImageRange(tmp), 00696 destIter(dupperleftx, dax), kernel1d(smooth)); 00697 separableConvolveX(srcIterRange(supperleft, slowerright, sa), 00698 destImage(tmp), kernel1d(smooth)); 00699 separableConvolveY(srcImageRange(tmp), 00700 destIter(dupperlefty, day), kernel1d(grad)); 00701 } 00702 00703 template <class SrcIterator, class SrcAccessor, 00704 class DestIterator, class DestAccessor> 00705 void gaussianGradient(SrcIterator supperleft, 00706 SrcIterator slowerright, SrcAccessor src, 00707 DestIterator dupperleft, DestAccessor dest, 00708 double scale) 00709 { 00710 VectorElementAccessor<DestAccessor> gradx(0, dest), grady(1, dest); 00711 gaussianGradient(supperleft, slowerright, src, 00712 dupperleft, gradx, dupperleft, grady, scale); 00713 } 00714 00715 template <class SrcIterator, class SrcAccessor, 00716 class DestIteratorX, class DestAccessorX, 00717 class DestIteratorY, class DestAccessorY> 00718 inline void 00719 gaussianGradient(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00720 pair<DestIteratorX, DestAccessorX> destx, 00721 pair<DestIteratorY, DestAccessorY> desty, 00722 double scale) 00723 { 00724 gaussianGradient(src.first, src.second, src.third, 00725 destx.first, destx.second, desty.first, desty.second, scale); 00726 } 00727 00728 template <class SrcIterator, class SrcAccessor, 00729 class DestIterator, class DestAccessor> 00730 inline void 00731 gaussianGradient(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00732 pair<DestIterator, DestAccessor> dest, 00733 double scale) 00734 { 00735 gaussianGradient(src.first, src.second, src.third, 00736 dest.first, dest.second, scale); 00737 } 00738 00739 /** \brief Calculate the gradient magnitude by means of a 1st derivatives of 00740 Gaussian filter. 00741 00742 This function calls gaussianGradient() and returns the pixel-wise magnitude of 00743 the resulting gradient vectors. If the original image has multiple bands, 00744 the squared gradient magnitude is computed for each band separately, and the 00745 return value is the square root of the sum of these sqaured magnitudes. 00746 00747 <b> Declarations:</b> 00748 00749 pass arguments explicitly: 00750 \code 00751 namespace vigra { 00752 template <class SrcIterator, class SrcAccessor, 00753 class DestIterator, class DestAccessor> 00754 void gaussianGradientMagnitude(SrcIterator sul, 00755 SrcIterator slr, SrcAccessor src, 00756 DestIterator dupperleft, DestAccessor dest, 00757 double scale); 00758 } 00759 \endcode 00760 00761 00762 use argument objects in conjunction with \ref ArgumentObjectFactories : 00763 \code 00764 namespace vigra { 00765 template <class SrcIterator, class SrcAccessor, 00766 class DestIterator, class DestAccessor> 00767 void 00768 gaussianGradientMagnitude(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00769 pair<DestIterator, DestAccessor> dest, 00770 double scale); 00771 } 00772 \endcode 00773 00774 <b> Usage:</b> 00775 00776 <b>\#include</b> <<a href="convolution_8hxx-source.html">vigra/convolution.hxx</a>> 00777 00778 00779 \code 00780 vigra::FImage src(w,h), grad(w,h); 00781 ... 00782 00783 // calculate gradient magnitude at scale = 3.0 00784 vigra::gaussianGradientMagnitude(srcImageRange(src), destImage(grad), 3.0); 00785 00786 \endcode 00787 00788 */ 00789 doxygen_overloaded_function(template <...> void gaussianGradientMagnitude) 00790 00791 template <class SrcIterator, class SrcAccessor, 00792 class DestIterator, class DestAccessor> 00793 void gaussianGradientMagnitude(SrcIterator sul, 00794 SrcIterator slr, SrcAccessor src, 00795 DestIterator dupperleft, DestAccessor dest, 00796 double scale) 00797 { 00798 typedef typename NumericTraits<typename SrcAccessor::value_type>::RealPromote TmpType; 00799 BasicImage<TmpType> gradx(slr-sul, SkipInitialization), grady(slr-sul, SkipInitialization); 00800 00801 gaussianGradient(srcIterRange(sul, slr, src), 00802 destImage(gradx), destImage(grady), scale); 00803 combineTwoImages(srcImageRange(gradx), srcImage(grady), destIter(dupperleft, dest), 00804 MagnitudeFunctor<TmpType>()); 00805 } 00806 00807 template <class SrcIterator, class SrcAccessor, 00808 class DestIterator, class DestAccessor> 00809 inline void 00810 gaussianGradientMagnitude(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00811 pair<DestIterator, DestAccessor> dest, 00812 double scale) 00813 { 00814 gaussianGradientMagnitude(src.first, src.second, src.third, 00815 dest.first, dest.second, scale); 00816 } 00817 00818 /********************************************************/ 00819 /* */ 00820 /* laplacianOfGaussian */ 00821 /* */ 00822 /********************************************************/ 00823 00824 /** \brief Filter image with the Laplacian of Gaussian operator 00825 at the given scale. 00826 00827 This function calls \ref separableConvolveX() and \ref separableConvolveY() with the appropriate 2nd derivative 00828 of Gaussian kernels in x- and y-direction and then sums the results 00829 to get the Laplacian. 00830 00831 <b> Declarations:</b> 00832 00833 pass arguments explicitly: 00834 \code 00835 namespace vigra { 00836 template <class SrcIterator, class SrcAccessor, 00837 class DestIterator, class DestAccessor> 00838 void laplacianOfGaussian(SrcIterator supperleft, 00839 SrcIterator slowerright, SrcAccessor sa, 00840 DestIterator dupperleft, DestAccessor da, 00841 double scale); 00842 } 00843 \endcode 00844 00845 00846 use argument objects in conjunction with \ref ArgumentObjectFactories : 00847 \code 00848 namespace vigra { 00849 template <class SrcIterator, class SrcAccessor, 00850 class DestIterator, class DestAccessor> 00851 inline void 00852 laplacianOfGaussian(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00853 pair<DestIterator, DestAccessor> dest, 00854 double scale); 00855 } 00856 \endcode 00857 00858 <b> Usage:</b> 00859 00860 <b>\#include</b> <<a href="convolution_8hxx-source.html">vigra/convolution.hxx</a>> 00861 00862 00863 \code 00864 vigra::FImage src(w,h), dest(w,h); 00865 ... 00866 00867 // calculate Laplacian of Gaussian at scale = 3.0 00868 vigra::laplacianOfGaussian(srcImageRange(src), destImage(dest), 3.0); 00869 00870 \endcode 00871 00872 */ 00873 doxygen_overloaded_function(template <...> void laplacianOfGaussian) 00874 00875 template <class SrcIterator, class SrcAccessor, 00876 class DestIterator, class DestAccessor> 00877 void laplacianOfGaussian(SrcIterator supperleft, 00878 SrcIterator slowerright, SrcAccessor sa, 00879 DestIterator dupperleft, DestAccessor da, 00880 double scale) 00881 { 00882 typedef typename 00883 NumericTraits<typename SrcAccessor::value_type>::RealPromote 00884 TmpType; 00885 BasicImage<TmpType> tmp(slowerright - supperleft, SkipInitialization), 00886 tmpx(slowerright - supperleft, SkipInitialization), 00887 tmpy(slowerright - supperleft, SkipInitialization); 00888 00889 Kernel1D<double> smooth, deriv; 00890 smooth.initGaussian(scale); 00891 deriv.initGaussianDerivative(scale, 2); 00892 00893 separableConvolveX(srcIterRange(supperleft, slowerright, sa), 00894 destImage(tmp), kernel1d(deriv)); 00895 separableConvolveY(srcImageRange(tmp), 00896 destImage(tmpx), kernel1d(smooth)); 00897 separableConvolveX(srcIterRange(supperleft, slowerright, sa), 00898 destImage(tmp), kernel1d(smooth)); 00899 separableConvolveY(srcImageRange(tmp), 00900 destImage(tmpy), kernel1d(deriv)); 00901 combineTwoImages(srcImageRange(tmpx), srcImage(tmpy), 00902 destIter(dupperleft, da), std::plus<TmpType>()); 00903 } 00904 00905 template <class SrcIterator, class SrcAccessor, 00906 class DestIterator, class DestAccessor> 00907 inline void 00908 laplacianOfGaussian(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00909 pair<DestIterator, DestAccessor> dest, 00910 double scale) 00911 { 00912 laplacianOfGaussian(src.first, src.second, src.third, 00913 dest.first, dest.second, scale); 00914 } 00915 00916 /********************************************************/ 00917 /* */ 00918 /* hessianMatrixOfGaussian */ 00919 /* */ 00920 /********************************************************/ 00921 00922 /** \brief Filter image with the 2nd derivatives of the Gaussian 00923 at the given scale to get the Hessian matrix. 00924 00925 The Hessian matrix is a symmetric matrix defined as: 00926 00927 \f[ 00928 \mbox{\rm Hessian}(I) = \left( 00929 \begin{array}{cc} 00930 G_{xx} \ast I & G_{xy} \ast I \\ 00931 G_{xy} \ast I & G_{yy} \ast I 00932 \end{array} \right) 00933 \f] 00934 00935 where \f$G_{xx}, G_{xy}, G_{yy}\f$ denote 2nd derivatives of Gaussians 00936 at the given scale, and 00937 \f$\ast\f$ is the convolution symbol. This function calls 00938 \ref separableConvolveX() and \ref separableConvolveY() 00939 with the appropriate 2nd derivative 00940 of Gaussian kernels and puts the results in 00941 the three destination images. The first destination image will 00942 contain the second derivative in x-direction, the second one the mixed 00943 derivative, and the third one holds the derivative in y-direction. 00944 00945 <b> Declarations:</b> 00946 00947 pass arguments explicitly: 00948 \code 00949 namespace vigra { 00950 template <class SrcIterator, class SrcAccessor, 00951 class DestIteratorX, class DestAccessorX, 00952 class DestIteratorXY, class DestAccessorXY, 00953 class DestIteratorY, class DestAccessorY> 00954 void hessianMatrixOfGaussian(SrcIterator supperleft, 00955 SrcIterator slowerright, SrcAccessor sa, 00956 DestIteratorX dupperleftx, DestAccessorX dax, 00957 DestIteratorXY dupperleftxy, DestAccessorXY daxy, 00958 DestIteratorY dupperlefty, DestAccessorY day, 00959 double scale); 00960 } 00961 \endcode 00962 00963 00964 use argument objects in conjunction with \ref ArgumentObjectFactories : 00965 \code 00966 namespace vigra { 00967 template <class SrcIterator, class SrcAccessor, 00968 class DestIteratorX, class DestAccessorX, 00969 class DestIteratorXY, class DestAccessorXY, 00970 class DestIteratorY, class DestAccessorY> 00971 inline void 00972 hessianMatrixOfGaussian(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00973 pair<DestIteratorX, DestAccessorX> destx, 00974 pair<DestIteratorXY, DestAccessorXY> destxy, 00975 pair<DestIteratorY, DestAccessorY> desty, 00976 double scale); 00977 } 00978 \endcode 00979 00980 <b> Usage:</b> 00981 00982 <b>\#include</b> <<a href="convolution_8hxx-source.html">vigra/convolution.hxx</a>> 00983 00984 00985 \code 00986 vigra::FImage src(w,h), hxx(w,h), hxy(w,h), hyy(w,h); 00987 ... 00988 00989 // calculate Hessian of Gaussian at scale = 3.0 00990 vigra::hessianMatrixOfGaussian(srcImageRange(src), 00991 destImage(hxx), destImage(hxy), destImage(hyy), 3.0); 00992 00993 \endcode 00994 00995 */ 00996 doxygen_overloaded_function(template <...> void hessianMatrixOfGaussian) 00997 00998 template <class SrcIterator, class SrcAccessor, 00999 class DestIteratorX, class DestAccessorX, 01000 class DestIteratorXY, class DestAccessorXY, 01001 class DestIteratorY, class DestAccessorY> 01002 void hessianMatrixOfGaussian(SrcIterator supperleft, 01003 SrcIterator slowerright, SrcAccessor sa, 01004 DestIteratorX dupperleftx, DestAccessorX dax, 01005 DestIteratorXY dupperleftxy, DestAccessorXY daxy, 01006 DestIteratorY dupperlefty, DestAccessorY day, 01007 double scale) 01008 { 01009 typedef typename 01010 NumericTraits<typename SrcAccessor::value_type>::RealPromote 01011 TmpType; 01012 BasicImage<TmpType> tmp(slowerright - supperleft, SkipInitialization); 01013 01014 Kernel1D<double> smooth, deriv1, deriv2; 01015 smooth.initGaussian(scale); 01016 deriv1.initGaussianDerivative(scale, 1); 01017 deriv2.initGaussianDerivative(scale, 2); 01018 01019 separableConvolveX(srcIterRange(supperleft, slowerright, sa), 01020 destImage(tmp), kernel1d(deriv2)); 01021 separableConvolveY(srcImageRange(tmp), 01022 destIter(dupperleftx, dax), kernel1d(smooth)); 01023 separableConvolveX(srcIterRange(supperleft, slowerright, sa), 01024 destImage(tmp), kernel1d(smooth)); 01025 separableConvolveY(srcImageRange(tmp), 01026 destIter(dupperlefty, day), kernel1d(deriv2)); 01027 separableConvolveX(srcIterRange(supperleft, slowerright, sa), 01028 destImage(tmp), kernel1d(deriv1)); 01029 separableConvolveY(srcImageRange(tmp), 01030 destIter(dupperleftxy, daxy), kernel1d(deriv1)); 01031 } 01032 01033 template <class SrcIterator, class SrcAccessor, 01034 class DestIteratorX, class DestAccessorX, 01035 class DestIteratorXY, class DestAccessorXY, 01036 class DestIteratorY, class DestAccessorY> 01037 inline void 01038 hessianMatrixOfGaussian(triple<SrcIterator, SrcIterator, SrcAccessor> src, 01039 pair<DestIteratorX, DestAccessorX> destx, 01040 pair<DestIteratorXY, DestAccessorXY> destxy, 01041 pair<DestIteratorY, DestAccessorY> desty, 01042 double scale) 01043 { 01044 hessianMatrixOfGaussian(src.first, src.second, src.third, 01045 destx.first, destx.second, 01046 destxy.first, destxy.second, 01047 desty.first, desty.second, 01048 scale); 01049 } 01050 01051 /********************************************************/ 01052 /* */ 01053 /* structureTensor */ 01054 /* */ 01055 /********************************************************/ 01056 01057 /** \brief Calculate the Structure Tensor for each pixel of 01058 and image, using Gaussian (derivative) filters. 01059 01060 The Structure Tensor is is a smoothed version of the Euclidean product 01061 of the gradient vector with itself. I.e. it's a symmetric matrix defined as: 01062 01063 \f[ 01064 \mbox{\rm StructurTensor}(I) = \left( 01065 \begin{array}{cc} 01066 G \ast (I_x I_x) & G \ast (I_x I_y) \\ 01067 G \ast (I_x I_y) & G \ast (I_y I_y) 01068 \end{array} \right) = \left( 01069 \begin{array}{cc} 01070 A & C \\ 01071 C & B 01072 \end{array} \right) 01073 \f] 01074 01075 where \f$G\f$ denotes Gaussian smoothing at the <i>outer scale</i>, 01076 \f$I_x, I_y\f$ are the gradient components taken at the <i>inner scale</i>, 01077 \f$\ast\f$ is the convolution symbol, and \f$I_x I_x\f$ etc. are pixelwise 01078 products of the 1st derivative images. This function calls 01079 \ref separableConvolveX() and \ref separableConvolveY() with the 01080 appropriate Gaussian kernels and puts the results in 01081 the three separate destination images (where the first one will 01082 contain \f$G \ast (I_x I_x)\f$, the second one \f$G \ast (I_x I_y)\f$, and the 01083 third one holds \f$G \ast (I_y I_y)\f$), or into a single 3-band image (where the bands 01084 hold the result in the same order as above). The latter form is also applicable when 01085 the source image is a multi-band image (e.g. RGB). In this case, tensors are 01086 first computed for each band separately, and then summed up to get a single result tensor. 01087 01088 <b> Declarations:</b> 01089 01090 pass arguments explicitly: 01091 \code 01092 namespace vigra { 01093 // create three separate destination images 01094 template <class SrcIterator, class SrcAccessor, 01095 class DestIteratorX, class DestAccessorX, 01096 class DestIteratorXY, class DestAccessorXY, 01097 class DestIteratorY, class DestAccessorY> 01098 void structureTensor(SrcIterator supperleft, 01099 SrcIterator slowerright, SrcAccessor sa, 01100 DestIteratorX dupperleftx, DestAccessorX dax, 01101 DestIteratorXY dupperleftxy, DestAccessorXY daxy, 01102 DestIteratorY dupperlefty, DestAccessorY day, 01103 double inner_scale, double outer_scale); 01104 01105 // create a single 3-band destination image 01106 template <class SrcIterator, class SrcAccessor, 01107 class DestIterator, class DestAccessor> 01108 void structureTensor(SrcIterator supperleft, 01109 SrcIterator slowerright, SrcAccessor sa, 01110 DestIterator dupperleft, DestAccessor da, 01111 double inner_scale, double outer_scale); 01112 } 01113 \endcode 01114 01115 01116 use argument objects in conjunction with \ref ArgumentObjectFactories : 01117 \code 01118 namespace vigra { 01119 // create three separate destination images 01120 template <class SrcIterator, class SrcAccessor, 01121 class DestIteratorX, class DestAccessorX, 01122 class DestIteratorXY, class DestAccessorXY, 01123 class DestIteratorY, class DestAccessorY> 01124 void 01125 structureTensor(triple<SrcIterator, SrcIterator, SrcAccessor> src, 01126 pair<DestIteratorX, DestAccessorX> destx, 01127 pair<DestIteratorXY, DestAccessorXY> destxy, 01128 pair<DestIteratorY, DestAccessorY> desty, 01129 double nner_scale, double outer_scale); 01130 01131 // create a single 3-band destination image 01132 template <class SrcIterator, class SrcAccessor, 01133 class DestIterator, class DestAccessor> 01134 void 01135 structureTensor(triple<SrcIterator, SrcIterator, SrcAccessor> src, 01136 pair<DestIterator, DestAccessor> dest, 01137 double nner_scale, double outer_scale); 01138 } 01139 \endcode 01140 01141 <b> Usage:</b> 01142 01143 <b>\#include</b> <<a href="convolution_8hxx-source.html">vigra/convolution.hxx</a>> 01144 01145 01146 \code 01147 vigra::FImage src(w,h), stxx(w,h), stxy(w,h), styy(w,h); 01148 vigra::BasicImage<TinyVector<float, 3> > st(w,h); 01149 ... 01150 01151 // calculate Structure Tensor at inner scale = 1.0 and outer scale = 3.0 01152 vigra::structureTensor(srcImageRange(src), 01153 destImage(stxx), destImage(stxy), destImage(styy), 1.0, 3.0); 01154 01155 // dto. with a single 3-band destination image 01156 vigra::structureTensor(srcImageRange(src), destImage(st), 1.0, 3.0); 01157 01158 \endcode 01159 01160 */ 01161 doxygen_overloaded_function(template <...> void structureTensor) 01162 01163 template <class SrcIterator, class SrcAccessor, 01164 class DestIteratorX, class DestAccessorX, 01165 class DestIteratorXY, class DestAccessorXY, 01166 class DestIteratorY, class DestAccessorY> 01167 void structureTensor(SrcIterator supperleft, 01168 SrcIterator slowerright, SrcAccessor sa, 01169 DestIteratorX dupperleftx, DestAccessorX dax, 01170 DestIteratorXY dupperleftxy, DestAccessorXY daxy, 01171 DestIteratorY dupperlefty, DestAccessorY day, 01172 double inner_scale, double outer_scale) 01173 { 01174 typedef typename 01175 NumericTraits<typename SrcAccessor::value_type>::RealPromote 01176 TmpType; 01177 BasicImage<TmpType> tmp(slowerright - supperleft, SkipInitialization), 01178 tmpx(slowerright - supperleft, SkipInitialization), 01179 tmpy(slowerright - supperleft, SkipInitialization); 01180 01181 gaussianGradient(srcIterRange(supperleft, slowerright, sa), 01182 destImage(tmpx), destImage(tmpy), inner_scale); 01183 combineTwoImages(srcImageRange(tmpx), srcImage(tmpx), 01184 destImage(tmp), std::multiplies<TmpType>()); 01185 gaussianSmoothing(srcImageRange(tmp), 01186 destIter(dupperleftx, dax), outer_scale); 01187 combineTwoImages(srcImageRange(tmpy), srcImage(tmpy), 01188 destImage(tmp), std::multiplies<TmpType>()); 01189 gaussianSmoothing(srcImageRange(tmp), 01190 destIter(dupperlefty, day), outer_scale); 01191 combineTwoImages(srcImageRange(tmpx), srcImage(tmpy), 01192 destImage(tmp), std::multiplies<TmpType>()); 01193 gaussianSmoothing(srcImageRange(tmp), 01194 destIter(dupperleftxy, daxy), outer_scale); 01195 } 01196 01197 template <class SrcIterator, class SrcAccessor, 01198 class DestIteratorX, class DestAccessorX, 01199 class DestIteratorXY, class DestAccessorXY, 01200 class DestIteratorY, class DestAccessorY> 01201 inline void 01202 structureTensor(triple<SrcIterator, SrcIterator, SrcAccessor> src, 01203 pair<DestIteratorX, DestAccessorX> destx, 01204 pair<DestIteratorXY, DestAccessorXY> destxy, 01205 pair<DestIteratorY, DestAccessorY> desty, 01206 double inner_scale, double outer_scale) 01207 { 01208 structureTensor(src.first, src.second, src.third, 01209 destx.first, destx.second, 01210 destxy.first, destxy.second, 01211 desty.first, desty.second, 01212 inner_scale, outer_scale); 01213 } 01214 01215 namespace detail { 01216 01217 template <class SrcIterator, class SrcAccessor, 01218 class DestIterator, class DestAccessor> 01219 void structureTensor(SrcIterator supperleft, 01220 SrcIterator slowerright, SrcAccessor src, 01221 DestIterator dupperleft, DestAccessor dest, 01222 double inner_scale, double outer_scale, 01223 VigraTrueType /* isScalar */) 01224 { 01225 typedef VectorElementAccessor<DestAccessor> DA; 01226 structureTensor(supperleft, slowerright, src, 01227 dupperleft, DA(0, dest), 01228 dupperleft, DA(1, dest), 01229 dupperleft, DA(2, dest), 01230 inner_scale, outer_scale); 01231 } 01232 01233 template <class SrcIterator, class SrcAccessor, 01234 class DestIterator, class DestAccessor> 01235 void structureTensor(SrcIterator supperleft, 01236 SrcIterator slowerright, SrcAccessor src, 01237 DestIterator dupperleft, DestAccessor dest, 01238 double inner_scale, double outer_scale, 01239 VigraFalseType /* isScalar */) 01240 { 01241 int bands = src.size(supperleft); 01242 typedef VectorElementAccessor<SrcAccessor> SA; 01243 01244 structureTensor(supperleft, slowerright, SA(0, src), 01245 dupperleft, dest, 01246 inner_scale, outer_scale, 01247 VigraTrueType() /* isScalar */); 01248 01249 BasicImage<typename DestAccessor::value_type> st(slowerright - supperleft, SkipInitialization); 01250 for(int k=1; k < bands; ++k) 01251 { 01252 structureTensor(supperleft, slowerright, SA(k, src), 01253 st.upperLeft(), st.accessor(), 01254 inner_scale, outer_scale, 01255 VigraTrueType() /* isScalar */); 01256 combineTwoImages(srcImageRange(st), srcIter(dupperleft, dest), destIter(dupperleft, dest), 01257 std::plus<typename DestAccessor::value_type>()); 01258 } 01259 } 01260 01261 } // namespace detail 01262 01263 template <class SrcIterator, class SrcAccessor, 01264 class DestIterator, class DestAccessor> 01265 void structureTensor(SrcIterator supperleft, 01266 SrcIterator slowerright, SrcAccessor src, 01267 DestIterator dupperleft, DestAccessor dest, 01268 double inner_scale, double outer_scale) 01269 { 01270 typedef typename 01271 NumericTraits<typename SrcAccessor::value_type>::isScalar isScalar; 01272 detail::structureTensor(supperleft, slowerright, src, 01273 dupperleft, dest, inner_scale, outer_scale, isScalar()); 01274 } 01275 01276 template <class SrcIterator, class SrcAccessor, 01277 class DestIterator, class DestAccessor> 01278 inline void 01279 structureTensor(triple<SrcIterator, SrcIterator, SrcAccessor> src, 01280 pair<DestIterator, DestAccessor> dest, 01281 double inner_scale, double outer_scale) 01282 { 01283 structureTensor(src.first, src.second, src.third, 01284 dest.first, dest.second, 01285 inner_scale, outer_scale); 01286 } 01287 01288 //@} 01289 01290 } // namespace vigra 01291 01292 #endif // VIGRA_CONVOLUTION_HXX
© Ullrich Köthe (ullrich.koethe@iwr.uni-heidelberg.de) |
html generated using doxygen and Python
|