[ VIGRA Homepage | Function Index | Class Index | Namespaces | File List | Main Page ]
00001 /************************************************************************/ 00002 /* */ 00003 /* Copyright 1998-2004 by Ullrich Koethe */ 00004 /* Cognitive Systems Group, University of Hamburg, Germany */ 00005 /* */ 00006 /* This file is part of the VIGRA computer vision library. */ 00007 /* The VIGRA Website is */ 00008 /* http://kogs-www.informatik.uni-hamburg.de/~koethe/vigra/ */ 00009 /* Please direct questions, bug reports, and contributions to */ 00010 /* ullrich.koethe@iwr.uni-heidelberg.de or */ 00011 /* vigra@informatik.uni-hamburg.de */ 00012 /* */ 00013 /* Permission is hereby granted, free of charge, to any person */ 00014 /* obtaining a copy of this software and associated documentation */ 00015 /* files (the "Software"), to deal in the Software without */ 00016 /* restriction, including without limitation the rights to use, */ 00017 /* copy, modify, merge, publish, distribute, sublicense, and/or */ 00018 /* sell copies of the Software, and to permit persons to whom the */ 00019 /* Software is furnished to do so, subject to the following */ 00020 /* conditions: */ 00021 /* */ 00022 /* The above copyright notice and this permission notice shall be */ 00023 /* included in all copies or substantial portions of the */ 00024 /* Software. */ 00025 /* */ 00026 /* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND */ 00027 /* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES */ 00028 /* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND */ 00029 /* NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT */ 00030 /* HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, */ 00031 /* WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING */ 00032 /* FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR */ 00033 /* OTHER DEALINGS IN THE SOFTWARE. */ 00034 /* */ 00035 /************************************************************************/ 00036 00037 00038 #ifndef VIGRA_BOUNDARYTENSOR_HXX 00039 #define VIGRA_BOUNDARYTENSOR_HXX 00040 00041 #include <cmath> 00042 #include <functional> 00043 #include "utilities.hxx" 00044 #include "array_vector.hxx" 00045 #include "basicimage.hxx" 00046 #include "combineimages.hxx" 00047 #include "numerictraits.hxx" 00048 #include "convolution.hxx" 00049 00050 namespace vigra { 00051 00052 namespace detail { 00053 00054 /***********************************************************************/ 00055 00056 typedef ArrayVector<Kernel1D<double> > KernelArray; 00057 00058 void 00059 initGaussianPolarFilters1(double std_dev, KernelArray & k) 00060 { 00061 typedef KernelArray::value_type Kernel; 00062 typedef Kernel::iterator iterator; 00063 00064 vigra_precondition(std_dev >= 0.0, 00065 "initGaussianPolarFilter1(): " 00066 "Standard deviation must be >= 0."); 00067 00068 k.resize(4); 00069 00070 int radius = (int)(4.0*std_dev + 0.5); 00071 std_dev *= 1.08179074376; 00072 double f = 1.0 / VIGRA_CSTD::sqrt(2.0 * M_PI) / std_dev; // norm 00073 double a = 0.558868151788 / VIGRA_CSTD::pow(std_dev, 5); 00074 double b = -2.04251639729 / VIGRA_CSTD::pow(std_dev, 3); 00075 double sigma22 = -0.5 / std_dev / std_dev; 00076 00077 00078 for(unsigned int i=0; i<k.size(); ++i) 00079 { 00080 k[i].initExplicitly(-radius, radius); 00081 k[i].setBorderTreatment(BORDER_TREATMENT_REFLECT); 00082 } 00083 00084 int ix; 00085 iterator c = k[0].center(); 00086 for(ix=-radius; ix<=radius; ++ix) 00087 { 00088 double x = (double)ix; 00089 c[ix] = f * VIGRA_CSTD::exp(sigma22 * x * x); 00090 } 00091 00092 c = k[1].center(); 00093 for(ix=-radius; ix<=radius; ++ix) 00094 { 00095 double x = (double)ix; 00096 c[ix] = f * x * VIGRA_CSTD::exp(sigma22 * x * x); 00097 } 00098 00099 c = k[2].center(); 00100 double b2 = b / 3.0; 00101 for(ix=-radius; ix<=radius; ++ix) 00102 { 00103 double x = (double)ix; 00104 c[ix] = f * (b2 + a * x * x) * VIGRA_CSTD::exp(sigma22 * x * x); 00105 } 00106 00107 c = k[3].center(); 00108 for(ix=-radius; ix<=radius; ++ix) 00109 { 00110 double x = (double)ix; 00111 c[ix] = f * x * (b + a * x * x) * VIGRA_CSTD::exp(sigma22 * x * x); 00112 } 00113 } 00114 00115 void 00116 initGaussianPolarFilters2(double std_dev, KernelArray & k) 00117 { 00118 typedef Kernel1D<double>::iterator iterator; 00119 00120 vigra_precondition(std_dev >= 0.0, 00121 "initGaussianPolarFilter2(): " 00122 "Standard deviation must be >= 0."); 00123 00124 k.resize(3); 00125 00126 int radius = (int)(4.0*std_dev + 0.5); 00127 double f = 1.0 / VIGRA_CSTD::sqrt(2.0 * M_PI) / std_dev; // norm 00128 double sigma2 = std_dev*std_dev; 00129 double sigma22 = -0.5 / sigma2; 00130 00131 for(unsigned int i=0; i<k.size(); ++i) 00132 { 00133 k[i].initExplicitly(-radius, radius); 00134 k[i].setBorderTreatment(BORDER_TREATMENT_REFLECT); 00135 } 00136 00137 int ix; 00138 iterator c = k[0].center(); 00139 for(ix=-radius; ix<=radius; ++ix) 00140 { 00141 double x = (double)ix; 00142 c[ix] = f * VIGRA_CSTD::exp(sigma22 * x * x); 00143 } 00144 00145 c = k[1].center(); 00146 double f1 = f / sigma2; 00147 for(ix=-radius; ix<=radius; ++ix) 00148 { 00149 double x = (double)ix; 00150 c[ix] = f1 * x * VIGRA_CSTD::exp(sigma22 * x * x); 00151 } 00152 00153 c = k[2].center(); 00154 double f2 = f / (sigma2 * sigma2); 00155 for(ix=-radius; ix<=radius; ++ix) 00156 { 00157 double x = (double)ix; 00158 c[ix] = f2 * (x * x - sigma2) * VIGRA_CSTD::exp(sigma22 * x * x); 00159 } 00160 } 00161 00162 void 00163 initGaussianPolarFilters3(double std_dev, KernelArray & k) 00164 { 00165 typedef Kernel1D<double>::iterator iterator; 00166 00167 vigra_precondition(std_dev >= 0.0, 00168 "initGaussianPolarFilter3(): " 00169 "Standard deviation must be >= 0."); 00170 00171 k.resize(4); 00172 00173 int radius = (int)(4.0*std_dev + 0.5); 00174 std_dev *= 1.15470053838; 00175 double sigma22 = -0.5 / std_dev / std_dev; 00176 double f = 1.0 / VIGRA_CSTD::sqrt(2.0 * M_PI) / std_dev; // norm 00177 double a = 0.883887052922 / VIGRA_CSTD::pow(std_dev, 5); 00178 00179 for(unsigned int i=0; i<k.size(); ++i) 00180 { 00181 k[i].initExplicitly(-radius, radius); 00182 k[i].setBorderTreatment(BORDER_TREATMENT_REFLECT); 00183 } 00184 00185 //double b = -1.3786348292 / VIGRA_CSTD::pow(std_dev, 3); 00186 00187 int ix; 00188 iterator c = k[0].center(); 00189 for(ix=-radius; ix<=radius; ++ix) 00190 { 00191 double x = (double)ix; 00192 c[ix] = f * VIGRA_CSTD::exp(sigma22 * x * x); 00193 } 00194 00195 c = k[1].center(); 00196 for(ix=-radius; ix<=radius; ++ix) 00197 { 00198 double x = (double)ix; 00199 c[ix] = f * x * VIGRA_CSTD::exp(sigma22 * x * x); 00200 } 00201 00202 c = k[2].center(); 00203 double a2 = 3.0 * a; 00204 for(ix=-radius; ix<=radius; ++ix) 00205 { 00206 double x = (double)ix; 00207 c[ix] = f * a2 * x * x * VIGRA_CSTD::exp(sigma22 * x * x); 00208 } 00209 00210 c = k[3].center(); 00211 for(ix=-radius; ix<=radius; ++ix) 00212 { 00213 double x = (double)ix; 00214 c[ix] = f * a * x * x * x * VIGRA_CSTD::exp(sigma22 * x * x); 00215 } 00216 } 00217 00218 template <class SrcIterator, class SrcAccessor, 00219 class DestIterator, class DestAccessor> 00220 void 00221 evenPolarFilters(SrcIterator supperleft, SrcIterator slowerright, SrcAccessor src, 00222 DestIterator dupperleft, DestAccessor dest, 00223 double scale, bool noLaplacian) 00224 { 00225 vigra_precondition(dest.size(dupperleft) == 3, 00226 "evenPolarFilters(): image for even output must have 3 bands."); 00227 00228 int w = slowerright.x - supperleft.x; 00229 int h = slowerright.y - supperleft.y; 00230 00231 typedef typename 00232 NumericTraits<typename SrcAccessor::value_type>::RealPromote TmpType; 00233 typedef BasicImage<TinyVector<TmpType, 3> > TmpImage; 00234 typedef typename TmpImage::traverser TmpTraverser; 00235 TmpImage t(w, h); 00236 00237 KernelArray k2; 00238 initGaussianPolarFilters2(scale, k2); 00239 00240 // calculate filter responses for even filters 00241 VectorElementAccessor<typename TmpImage::Accessor> tmpBand(0, t.accessor()); 00242 convolveImage(srcIterRange(supperleft, slowerright, src), 00243 destImage(t, tmpBand), k2[2], k2[0]); 00244 tmpBand.setIndex(1); 00245 convolveImage(srcIterRange(supperleft, slowerright, src), 00246 destImage(t, tmpBand), k2[1], k2[1]); 00247 tmpBand.setIndex(2); 00248 convolveImage(srcIterRange(supperleft, slowerright, src), 00249 destImage(t, tmpBand), k2[0], k2[2]); 00250 00251 // create even tensor from filter responses 00252 TmpTraverser tul(t.upperLeft()); 00253 TmpTraverser tlr(t.lowerRight()); 00254 for(; tul.y != tlr.y; ++tul.y, ++dupperleft.y) 00255 { 00256 typename TmpTraverser::row_iterator tr = tul.rowIterator(); 00257 typename TmpTraverser::row_iterator trend = tr + w; 00258 typename DestIterator::row_iterator d = dupperleft.rowIterator(); 00259 if(noLaplacian) 00260 { 00261 for(; tr != trend; ++tr, ++d) 00262 { 00263 TmpType v = 0.5*sq((*tr)[0]-(*tr)[2]) + 2.0*sq((*tr)[1]); 00264 dest.setComponent(v, d, 0); 00265 dest.setComponent(0, d, 1); 00266 dest.setComponent(v, d, 2); 00267 } 00268 } 00269 else 00270 { 00271 for(; tr != trend; ++tr, ++d) 00272 { 00273 dest.setComponent(sq((*tr)[0]) + sq((*tr)[1]), d, 0); 00274 dest.setComponent(-(*tr)[1] * ((*tr)[0] + (*tr)[2]), d, 1); 00275 dest.setComponent(sq((*tr)[1]) + sq((*tr)[2]), d, 2); 00276 } 00277 } 00278 } 00279 } 00280 00281 template <class SrcIterator, class SrcAccessor, 00282 class DestIterator, class DestAccessor> 00283 void 00284 oddPolarFilters(SrcIterator supperleft, SrcIterator slowerright, SrcAccessor src, 00285 DestIterator dupperleft, DestAccessor dest, 00286 double scale, bool addResult) 00287 { 00288 vigra_precondition(dest.size(dupperleft) == 3, 00289 "oddPolarFilters(): image for odd output must have 3 bands."); 00290 00291 int w = slowerright.x - supperleft.x; 00292 int h = slowerright.y - supperleft.y; 00293 00294 typedef typename 00295 NumericTraits<typename SrcAccessor::value_type>::RealPromote TmpType; 00296 typedef BasicImage<TinyVector<TmpType, 4> > TmpImage; 00297 typedef typename TmpImage::traverser TmpTraverser; 00298 TmpImage t(w, h); 00299 00300 detail::KernelArray k1; 00301 detail::initGaussianPolarFilters1(scale, k1); 00302 00303 // calculate filter responses for odd filters 00304 VectorElementAccessor<typename TmpImage::Accessor> tmpBand(0, t.accessor()); 00305 convolveImage(srcIterRange(supperleft, slowerright, src), 00306 destImage(t, tmpBand), k1[3], k1[0]); 00307 tmpBand.setIndex(1); 00308 convolveImage(srcIterRange(supperleft, slowerright, src), 00309 destImage(t, tmpBand), k1[2], k1[1]); 00310 tmpBand.setIndex(2); 00311 convolveImage(srcIterRange(supperleft, slowerright, src), 00312 destImage(t, tmpBand), k1[1], k1[2]); 00313 tmpBand.setIndex(3); 00314 convolveImage(srcIterRange(supperleft, slowerright, src), 00315 destImage(t, tmpBand), k1[0], k1[3]); 00316 00317 // create odd tensor from filter responses 00318 TmpTraverser tul(t.upperLeft()); 00319 TmpTraverser tlr(t.lowerRight()); 00320 for(; tul.y != tlr.y; ++tul.y, ++dupperleft.y) 00321 { 00322 typename TmpTraverser::row_iterator tr = tul.rowIterator(); 00323 typename TmpTraverser::row_iterator trend = tr + w; 00324 typename DestIterator::row_iterator d = dupperleft.rowIterator(); 00325 if(addResult) 00326 { 00327 for(; tr != trend; ++tr, ++d) 00328 { 00329 TmpType d0 = (*tr)[0] + (*tr)[2]; 00330 TmpType d1 = -(*tr)[1] - (*tr)[3]; 00331 00332 dest.setComponent(dest.getComponent(d, 0) + sq(d0), d, 0); 00333 dest.setComponent(dest.getComponent(d, 1) + d0 * d1, d, 1); 00334 dest.setComponent(dest.getComponent(d, 2) + sq(d1), d, 2); 00335 } 00336 } 00337 else 00338 { 00339 for(; tr != trend; ++tr, ++d) 00340 { 00341 TmpType d0 = (*tr)[0] + (*tr)[2]; 00342 TmpType d1 = -(*tr)[1] - (*tr)[3]; 00343 00344 dest.setComponent(sq(d0), d, 0); 00345 dest.setComponent(d0 * d1, d, 1); 00346 dest.setComponent(sq(d1), d, 2); 00347 } 00348 } 00349 } 00350 } 00351 00352 } // namespace detail 00353 00354 /** \addtogroup CommonConvolutionFilters Common Filters 00355 */ 00356 //@{ 00357 00358 /********************************************************/ 00359 /* */ 00360 /* rieszTransformOfLOG */ 00361 /* */ 00362 /********************************************************/ 00363 00364 /** \brief Calculate Riesz transforms of the Laplacian of Gaussian. 00365 00366 The Riesz transforms of the Laplacian of Gaussian have the following transfer 00367 functions (defined in a polar coordinate representation of the frequency domain): 00368 00369 \f[ 00370 F_{\sigma}(r, \phi)=(i \cos \phi)^n (i \sin \phi)^m r^2 e^{-r^2 \sigma^2 / 2} 00371 \f] 00372 00373 where <i>n</i> = <tt>xorder</tt> and <i>m</i> = <tt>yorder</tt> determine th e 00374 order of the transform, and <tt>sigma > 0</tt> is the scale of the Laplacian 00375 of Gaussian. This function computes a good spatial domain approximation of 00376 these transforms for <tt>xorder + yorder <= 2</tt>. The filter responses may be used 00377 to calculate the monogenic signal or the boundary tensor. 00378 00379 <b> Declarations:</b> 00380 00381 pass arguments explicitly: 00382 \code 00383 namespace vigra { 00384 template <class SrcIterator, class SrcAccessor, 00385 class DestIterator, class DestAccessor> 00386 void rieszTransformOfLOG(SrcIterator supperleft, SrcIterator slowerright, SrcAccessor src, 00387 DestIterator dupperleft, DestAccessor dest, 00388 double scale, unsigned int xorder, unsigned int yorder); 00389 } 00390 \endcode 00391 00392 00393 use argument objects in conjunction with \ref ArgumentObjectFactories : 00394 \code 00395 namespace vigra { 00396 template <class SrcIterator, class SrcAccessor, 00397 class DestIterator, class DestAccessor> 00398 void rieszTransformOfLOG(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00399 pair<DestIterator, DestAccessor> dest, 00400 double scale, unsigned int xorder, unsigned int yorder); 00401 } 00402 \endcode 00403 00404 <b> Usage:</b> 00405 00406 <b>\#include</b> <<a href="boundarytensor_8hxx-source.html">vigra/boundarytensor.hxx</a>> 00407 00408 \code 00409 FImage impulse(17,17), res(17, 17); 00410 impulse(8,8) = 1.0; 00411 00412 // calculate the impulse response of the first order Riesz transform in x-direction 00413 rieszTransformOfLOG(srcImageRange(impulse), destImage(res), 2.0, 1, 0); 00414 \endcode 00415 00416 */ 00417 doxygen_overloaded_function(template <...> void rieszTransformOfLOG) 00418 00419 template <class SrcIterator, class SrcAccessor, 00420 class DestIterator, class DestAccessor> 00421 void rieszTransformOfLOG(SrcIterator supperleft, SrcIterator slowerright, SrcAccessor src, 00422 DestIterator dupperleft, DestAccessor dest, 00423 double scale, unsigned int xorder, unsigned int yorder) 00424 { 00425 unsigned int order = xorder + yorder; 00426 00427 vigra_precondition(order <= 2, 00428 "rieszTransformOfLOG(): can only compute Riesz transforms up to order 2."); 00429 vigra_precondition(scale > 0.0, 00430 "rieszTransformOfLOG(): scale must be positive."); 00431 00432 int w = slowerright.x - supperleft.x; 00433 int h = slowerright.y - supperleft.y; 00434 00435 typedef typename NumericTraits<typename SrcAccessor::value_type>::RealPromote TmpType; 00436 typedef BasicImage<TmpType> TmpImage; 00437 00438 switch(order) 00439 { 00440 case 0: 00441 { 00442 detail::KernelArray k2; 00443 detail::initGaussianPolarFilters2(scale, k2); 00444 00445 TmpImage t1(w, h), t2(w, h); 00446 00447 convolveImage(srcIterRange(supperleft, slowerright, src), 00448 destImage(t1), k2[2], k2[0]); 00449 convolveImage(srcIterRange(supperleft, slowerright, src), 00450 destImage(t2), k2[0], k2[2]); 00451 combineTwoImages(srcImageRange(t1), srcImage(t2), 00452 destIter(dupperleft, dest), std::plus<TmpType>()); 00453 break; 00454 } 00455 case 1: 00456 { 00457 detail::KernelArray k1; 00458 detail::initGaussianPolarFilters1(scale, k1); 00459 00460 TmpImage t1(w, h), t2(w, h); 00461 00462 if(xorder == 1) 00463 { 00464 convolveImage(srcIterRange(supperleft, slowerright, src), 00465 destImage(t1), k1[3], k1[0]); 00466 convolveImage(srcIterRange(supperleft, slowerright, src), 00467 destImage(t2), k1[1], k1[2]); 00468 } 00469 else 00470 { 00471 convolveImage(srcIterRange(supperleft, slowerright, src), 00472 destImage(t1), k1[0], k1[3]); 00473 convolveImage(srcIterRange(supperleft, slowerright, src), 00474 destImage(t2), k1[2], k1[1]); 00475 } 00476 combineTwoImages(srcImageRange(t1), srcImage(t2), 00477 destIter(dupperleft, dest), std::plus<TmpType>()); 00478 break; 00479 } 00480 case 2: 00481 { 00482 detail::KernelArray k2; 00483 detail::initGaussianPolarFilters2(scale, k2); 00484 00485 convolveImage(srcIterRange(supperleft, slowerright, src), 00486 destIter(dupperleft, dest), k2[xorder], k2[yorder]); 00487 break; 00488 } 00489 /* for test purposes only: compute 3rd order polar filters */ 00490 case 3: 00491 { 00492 detail::KernelArray k3; 00493 detail::initGaussianPolarFilters3(scale, k3); 00494 00495 TmpImage t1(w, h), t2(w, h); 00496 00497 if(xorder == 3) 00498 { 00499 convolveImage(srcIterRange(supperleft, slowerright, src), 00500 destImage(t1), k3[3], k3[0]); 00501 convolveImage(srcIterRange(supperleft, slowerright, src), 00502 destImage(t2), k3[1], k3[2]); 00503 } 00504 else 00505 { 00506 convolveImage(srcIterRange(supperleft, slowerright, src), 00507 destImage(t1), k3[0], k3[3]); 00508 convolveImage(srcIterRange(supperleft, slowerright, src), 00509 destImage(t2), k3[2], k3[1]); 00510 } 00511 combineTwoImages(srcImageRange(t1), srcImage(t2), 00512 destIter(dupperleft, dest), std::minus<TmpType>()); 00513 break; 00514 } 00515 } 00516 } 00517 00518 template <class SrcIterator, class SrcAccessor, 00519 class DestIterator, class DestAccessor> 00520 inline 00521 void rieszTransformOfLOG(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00522 pair<DestIterator, DestAccessor> dest, 00523 double scale, unsigned int xorder, unsigned int yorder) 00524 { 00525 rieszTransformOfLOG(src.first, src.second, src.third, dest.first, dest.second, 00526 scale, xorder, yorder); 00527 } 00528 //@} 00529 00530 /** \addtogroup TensorImaging Tensor Image Processing 00531 */ 00532 //@{ 00533 00534 /********************************************************/ 00535 /* */ 00536 /* boundaryTensor */ 00537 /* */ 00538 /********************************************************/ 00539 00540 /** \brief Calculate the boundary tensor for a scalar valued image. 00541 00542 These functions calculate a spatial domain approximation of 00543 the boundary tensor as described in 00544 00545 U. Köthe: <a href="http://kogs-www.informatik.uni-hamburg.de/~koethe/papers/abstracts/polarfilters.html"> 00546 <i>"Integrated Edge and Junction Detection with the Boundary Tensor"</i></a>, 00547 in: ICCV 03, Proc. of 9th Intl. Conf. on Computer Vision, Nice 2003, vol. 1, 00548 pp. 424-431, Los Alamitos: IEEE Computer Society, 2003 00549 00550 with the Laplacian of Gaussian as the underlying bandpass filter (see 00551 \ref rieszTransformOfLOG()). The output image must have 3 bands which will hold the 00552 tensor components in the order t11, t12 (== t21), t22. The function 00553 \ref boundaryTensor1() with the same interface implements a variant of the 00554 boundary tensor where the 0th-order Riesz transform has been dropped, so that the 00555 tensor is no longer sensitive to blobs. 00556 00557 <b> Declarations:</b> 00558 00559 pass arguments explicitly: 00560 \code 00561 namespace vigra { 00562 template <class SrcIterator, class SrcAccessor, 00563 class DestIterator, class DestAccessor> 00564 void boundaryTensor(SrcIterator supperleft, SrcIterator slowerright, SrcAccessor src, 00565 DestIterator dupperleft, DestAccessor dest, 00566 double scale); 00567 } 00568 \endcode 00569 00570 use argument objects in conjunction with \ref ArgumentObjectFactories : 00571 \code 00572 namespace vigra { 00573 template <class SrcIterator, class SrcAccessor, 00574 class DestIterator, class DestAccessor> 00575 void boundaryTensor(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00576 pair<DestIterator, DestAccessor> dest, 00577 double scale); 00578 } 00579 \endcode 00580 00581 <b> Usage:</b> 00582 00583 <b>\#include</b> <<a href="boundarytensor_8hxx-source.html">vigra/boundarytensor.hxx</a>> 00584 00585 \code 00586 FImage img(w,h); 00587 FVector3Image bt(w,h); 00588 ... 00589 boundaryTensor(srcImageRange(img), destImage(bt), 2.0); 00590 \endcode 00591 00592 */ 00593 doxygen_overloaded_function(template <...> void boundaryTensor) 00594 00595 template <class SrcIterator, class SrcAccessor, 00596 class DestIterator, class DestAccessor> 00597 void boundaryTensor(SrcIterator supperleft, SrcIterator slowerright, SrcAccessor src, 00598 DestIterator dupperleft, DestAccessor dest, 00599 double scale) 00600 { 00601 vigra_precondition(dest.size(dupperleft) == 3, 00602 "boundaryTensor(): image for even output must have 3 bands."); 00603 vigra_precondition(scale > 0.0, 00604 "boundaryTensor(): scale must be positive."); 00605 00606 detail::evenPolarFilters(supperleft, slowerright, src, 00607 dupperleft, dest, scale, false); 00608 detail::oddPolarFilters(supperleft, slowerright, src, 00609 dupperleft, dest, scale, true); 00610 } 00611 00612 template <class SrcIterator, class SrcAccessor, 00613 class DestIterator, class DestAccessor> 00614 inline 00615 void boundaryTensor(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00616 pair<DestIterator, DestAccessor> dest, 00617 double scale) 00618 { 00619 boundaryTensor(src.first, src.second, src.third, 00620 dest.first, dest.second, scale); 00621 } 00622 00623 /** \brief Boundary tensor variant. 00624 00625 This function implements a variant of the boundary tensor where the 00626 0th-order Riesz transform has been dropped, so that the tensor is no 00627 longer sensitive to blobs. See \ref boundaryTensor() for more detailed 00628 documentation. 00629 00630 <b> Declarations:</b> 00631 00632 <b>\#include</b> <<a href="boundarytensor_8hxx-source.html">vigra/boundarytensor.hxx</a>> 00633 00634 pass arguments explicitly: 00635 \code 00636 namespace vigra { 00637 template <class SrcIterator, class SrcAccessor, 00638 class DestIterator, class DestAccessor> 00639 void boundaryTensor1(SrcIterator supperleft, SrcIterator slowerright, SrcAccessor src, 00640 DestIterator dupperleft, DestAccessor dest, 00641 double scale); 00642 } 00643 \endcode 00644 00645 use argument objects in conjunction with \ref ArgumentObjectFactories : 00646 \code 00647 namespace vigra { 00648 template <class SrcIterator, class SrcAccessor, 00649 class DestIterator, class DestAccessor> 00650 void boundaryTensor1(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00651 pair<DestIterator, DestAccessor> dest, 00652 double scale); 00653 } 00654 \endcode 00655 */ 00656 doxygen_overloaded_function(template <...> void boundaryTensor1) 00657 00658 template <class SrcIterator, class SrcAccessor, 00659 class DestIterator, class DestAccessor> 00660 void boundaryTensor1(SrcIterator supperleft, SrcIterator slowerright, SrcAccessor src, 00661 DestIterator dupperleft, DestAccessor dest, 00662 double scale) 00663 { 00664 vigra_precondition(dest.size(dupperleft) == 3, 00665 "boundaryTensor1(): image for even output must have 3 bands."); 00666 vigra_precondition(scale > 0.0, 00667 "boundaryTensor1(): scale must be positive."); 00668 00669 detail::evenPolarFilters(supperleft, slowerright, src, 00670 dupperleft, dest, scale, true); 00671 detail::oddPolarFilters(supperleft, slowerright, src, 00672 dupperleft, dest, scale, true); 00673 } 00674 00675 template <class SrcIterator, class SrcAccessor, 00676 class DestIterator, class DestAccessor> 00677 inline 00678 void boundaryTensor1(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00679 pair<DestIterator, DestAccessor> dest, 00680 double scale) 00681 { 00682 boundaryTensor1(src.first, src.second, src.third, 00683 dest.first, dest.second, scale); 00684 } 00685 00686 /********************************************************/ 00687 /* */ 00688 /* boundaryTensor3 */ 00689 /* */ 00690 /********************************************************/ 00691 00692 /* Add 3rd order Riesz transform to boundary tensor 00693 ??? Does not work -- bug or too coarse approximation for 3rd order ??? 00694 */ 00695 template <class SrcIterator, class SrcAccessor, 00696 class DestIteratorEven, class DestAccessorEven, 00697 class DestIteratorOdd, class DestAccessorOdd> 00698 void boundaryTensor3(SrcIterator supperleft, SrcIterator slowerright, SrcAccessor sa, 00699 DestIteratorEven dupperleft_even, DestAccessorEven even, 00700 DestIteratorOdd dupperleft_odd, DestAccessorOdd odd, 00701 double scale) 00702 { 00703 vigra_precondition(even.size(dupperleft_even) == 3, 00704 "boundaryTensor3(): image for even output must have 3 bands."); 00705 vigra_precondition(odd.size(dupperleft_odd) == 3, 00706 "boundaryTensor3(): image for odd output must have 3 bands."); 00707 00708 detail::evenPolarFilters(supperleft, slowerright, sa, 00709 dupperleft_even, even, scale, false); 00710 00711 int w = slowerright.x - supperleft.x; 00712 int h = slowerright.y - supperleft.y; 00713 00714 typedef typename 00715 NumericTraits<typename SrcAccessor::value_type>::RealPromote TmpType; 00716 typedef BasicImage<TinyVector<TmpType, 4> > TmpImage; 00717 TmpImage t1(w, h), t2(w, h); 00718 00719 detail::KernelArray k1, k3; 00720 detail::initGaussianPolarFilters1(scale, k1); 00721 detail::initGaussianPolarFilters3(scale, k3); 00722 00723 // calculate filter responses for odd filters 00724 VectorElementAccessor<typename TmpImage::Accessor> tmpBand(0, t1.accessor()); 00725 convolveImage(srcIterRange(supperleft, slowerright, sa), 00726 destImage(t1, tmpBand), k1[3], k1[0]); 00727 tmpBand.setIndex(1); 00728 convolveImage(srcIterRange(supperleft, slowerright, sa), 00729 destImage(t1, tmpBand), k1[1], k1[2]); 00730 tmpBand.setIndex(2); 00731 convolveImage(srcIterRange(supperleft, slowerright, sa), 00732 destImage(t1, tmpBand), k3[3], k3[0]); 00733 tmpBand.setIndex(3); 00734 convolveImage(srcIterRange(supperleft, slowerright, sa), 00735 destImage(t1, tmpBand), k3[1], k3[2]); 00736 00737 tmpBand.setIndex(0); 00738 convolveImage(srcIterRange(supperleft, slowerright, sa), 00739 destImage(t2, tmpBand), k1[0], k1[3]); 00740 tmpBand.setIndex(1); 00741 convolveImage(srcIterRange(supperleft, slowerright, sa), 00742 destImage(t2, tmpBand), k1[2], k1[1]); 00743 tmpBand.setIndex(2); 00744 convolveImage(srcIterRange(supperleft, slowerright, sa), 00745 destImage(t2, tmpBand), k3[0], k3[3]); 00746 tmpBand.setIndex(3); 00747 convolveImage(srcIterRange(supperleft, slowerright, sa), 00748 destImage(t2, tmpBand), k3[2], k3[1]); 00749 00750 // create odd tensor from filter responses 00751 typedef typename TmpImage::traverser TmpTraverser; 00752 TmpTraverser tul1(t1.upperLeft()); 00753 TmpTraverser tlr1(t1.lowerRight()); 00754 TmpTraverser tul2(t2.upperLeft()); 00755 for(; tul1.y != tlr1.y; ++tul1.y, ++tul2.y, ++dupperleft_odd.y) 00756 { 00757 typename TmpTraverser::row_iterator tr1 = tul1.rowIterator(); 00758 typename TmpTraverser::row_iterator trend1 = tr1 + w; 00759 typename TmpTraverser::row_iterator tr2 = tul2.rowIterator(); 00760 typename DestIteratorOdd::row_iterator o = dupperleft_odd.rowIterator(); 00761 for(; tr1 != trend1; ++tr1, ++tr2, ++o) 00762 { 00763 TmpType d11 = (*tr1)[0] + (*tr1)[2]; 00764 TmpType d12 = -(*tr1)[1] - (*tr1)[3]; 00765 TmpType d31 = (*tr2)[0] - (*tr2)[2]; 00766 TmpType d32 = (*tr2)[1] - (*tr2)[3]; 00767 TmpType d111 = 0.75 * d11 + 0.25 * d31; 00768 TmpType d112 = 0.25 * (d12 + d32); 00769 TmpType d122 = 0.25 * (d11 - d31); 00770 TmpType d222 = 0.75 * d12 - 0.25 * d32; 00771 TmpType d2 = sq(d112); 00772 TmpType d3 = sq(d122); 00773 00774 odd.setComponent(0.25 * (sq(d111) + 2.0*d2 + d3), o, 0); 00775 odd.setComponent(0.25 * (d111*d112 + 2.0*d112*d122 + d122*d222), o, 1); 00776 odd.setComponent(0.25 * (d2 + 2.0*d3 + sq(d222)), o, 2); 00777 } 00778 } 00779 } 00780 00781 template <class SrcIterator, class SrcAccessor, 00782 class DestIteratorEven, class DestAccessorEven, 00783 class DestIteratorOdd, class DestAccessorOdd> 00784 inline 00785 void boundaryTensor3(triple<SrcIterator, SrcIterator, SrcAccessor> src, 00786 pair<DestIteratorEven, DestAccessorEven> even, 00787 pair<DestIteratorOdd, DestAccessorOdd> odd, 00788 double scale) 00789 { 00790 boundaryTensor3(src.first, src.second, src.third, 00791 even.first, even.second, odd.first, odd.second, scale); 00792 } 00793 00794 //@} 00795 00796 } // namespace vigra 00797 00798 #endif // VIGRA_BOUNDARYTENSOR_HXX
© Ullrich Köthe (ullrich.koethe@iwr.uni-heidelberg.de) |
html generated using doxygen and Python
|