[ VIGRA Homepage | Function Index | Class Index | Namespaces | File List | Main Page ]
00001 /************************************************************************/ 00002 /* */ 00003 /* Copyright 1998-2002 by Ullrich Koethe */ 00004 /* Cognitive Systems Group, University of Hamburg, Germany */ 00005 /* */ 00006 /* This file is part of the VIGRA computer vision library. */ 00007 /* ( Version 1.6.0, Aug 13 2008 ) */ 00008 /* The VIGRA Website is */ 00009 /* http://kogs-www.informatik.uni-hamburg.de/~koethe/vigra/ */ 00010 /* Please direct questions, bug reports, and contributions to */ 00011 /* ullrich.koethe@iwr.uni-heidelberg.de or */ 00012 /* vigra@informatik.uni-hamburg.de */ 00013 /* */ 00014 /* Permission is hereby granted, free of charge, to any person */ 00015 /* obtaining a copy of this software and associated documentation */ 00016 /* files (the "Software"), to deal in the Software without */ 00017 /* restriction, including without limitation the rights to use, */ 00018 /* copy, modify, merge, publish, distribute, sublicense, and/or */ 00019 /* sell copies of the Software, and to permit persons to whom the */ 00020 /* Software is furnished to do so, subject to the following */ 00021 /* conditions: */ 00022 /* */ 00023 /* The above copyright notice and this permission notice shall be */ 00024 /* included in all copies or substantial portions of the */ 00025 /* Software. */ 00026 /* */ 00027 /* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND */ 00028 /* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES */ 00029 /* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND */ 00030 /* NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT */ 00031 /* HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, */ 00032 /* WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING */ 00033 /* FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR */ 00034 /* OTHER DEALINGS IN THE SOFTWARE. */ 00035 /* */ 00036 /************************************************************************/ 00037 00038 00039 #ifndef VIGRA_EDGEDETECTION_HXX 00040 #define VIGRA_EDGEDETECTION_HXX 00041 00042 #include <vector> 00043 #include <queue> 00044 #include <cmath> // sqrt(), abs() 00045 #include "utilities.hxx" 00046 #include "numerictraits.hxx" 00047 #include "stdimage.hxx" 00048 #include "stdimagefunctions.hxx" 00049 #include "recursiveconvolution.hxx" 00050 #include "separableconvolution.hxx" 00051 #include "labelimage.hxx" 00052 #include "mathutil.hxx" 00053 #include "pixelneighborhood.hxx" 00054 #include "linear_solve.hxx" 00055 00056 00057 namespace vigra { 00058 00059 /** \addtogroup EdgeDetection Edge Detection 00060 Edge detectors based on first and second derivatives, 00061 and related post-processing. 00062 */ 00063 //@{ 00064 00065 /********************************************************/ 00066 /* */ 00067 /* differenceOfExponentialEdgeImage */ 00068 /* */ 00069 /********************************************************/ 00070 00071 /** \brief Detect and mark edges in an edge image using the Shen/Castan zero-crossing detector. 00072 00073 This operator applies an exponential filter to the source image 00074 at the given <TT>scale</TT> and subtracts the result from the original image. 00075 Zero crossings are detected in the resulting difference image. Whenever the 00076 gradient at a zero crossing is greater than the given <TT>gradient_threshold</TT>, 00077 an edge point is marked (using <TT>edge_marker</TT>) in the destination image on 00078 the darker side of the zero crossing (note that zero crossings occur 00079 <i>between</i> pixels). For example: 00080 00081 \code 00082 sign of difference image resulting edge points (*) 00083 00084 + - - * * . 00085 + + - => . * * 00086 + + + . . . 00087 \endcode 00088 00089 Non-edge pixels (<TT>.</TT>) remain untouched in the destination image. 00090 The result can be improved by the post-processing operation \ref removeShortEdges(). 00091 A more accurate edge placement can be achieved with the function 00092 \ref differenceOfExponentialCrackEdgeImage(). 00093 00094 The source value type 00095 (<TT>SrcAccessor::value_type</TT>) must be a linear algebra, i.e. addition, 00096 subtraction and multiplication of the type with itself, and multiplication 00097 with double and 00098 \ref NumericTraits "NumericTraits" must 00099 be defined. In addition, this type must be less-comparable. 00100 00101 <b> Declarations:</b> 00102 00103 pass arguments explicitly: 00104 \code 00105 namespace vigra { 00106 template <class SrcIterator, class SrcAccessor, 00107 class DestIterator, class DestAccessor, 00108 class GradValue, 00109 class DestValue = DestAccessor::value_type> 00110 void differenceOfExponentialEdgeImage( 00111 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 00112 DestIterator dul, DestAccessor da, 00113 double scale, GradValue gradient_threshold, 00114 DestValue edge_marker = NumericTraits<DestValue>::one()) 00115 } 00116 \endcode 00117 00118 use argument objects in conjunction with \ref ArgumentObjectFactories : 00119 \code 00120 namespace vigra { 00121 template <class SrcIterator, class SrcAccessor, 00122 class DestIterator, class DestAccessor, 00123 class GradValue, 00124 class DestValue = DestAccessor::value_type> 00125 void differenceOfExponentialEdgeImage( 00126 triple<SrcIterator, SrcIterator, SrcAccessor> src, 00127 pair<DestIterator, DestAccessor> dest, 00128 double scale, GradValue gradient_threshold, 00129 DestValue edge_marker = NumericTraits<DestValue>::one()) 00130 } 00131 \endcode 00132 00133 <b> Usage:</b> 00134 00135 <b>\#include</b> <<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>><br> 00136 Namespace: vigra 00137 00138 \code 00139 vigra::BImage src(w,h), edges(w,h); 00140 00141 // empty edge image 00142 edges = 0; 00143 ... 00144 00145 // find edges at scale 0.8 with gradient larger than 4.0, mark with 1 00146 vigra::differenceOfExponentialEdgeImage(srcImageRange(src), destImage(edges), 00147 0.8, 4.0, 1); 00148 \endcode 00149 00150 <b> Required Interface:</b> 00151 00152 \code 00153 SrcImageIterator src_upperleft, src_lowerright; 00154 DestImageIterator dest_upperleft; 00155 00156 SrcAccessor src_accessor; 00157 DestAccessor dest_accessor; 00158 00159 SrcAccessor::value_type u = src_accessor(src_upperleft); 00160 double d; 00161 GradValue gradient_threshold; 00162 00163 u = u + u 00164 u = u - u 00165 u = u * u 00166 u = d * u 00167 u < gradient_threshold 00168 00169 DestValue edge_marker; 00170 dest_accessor.set(edge_marker, dest_upperleft); 00171 \endcode 00172 00173 <b> Preconditions:</b> 00174 00175 \code 00176 scale > 0 00177 gradient_threshold > 0 00178 \endcode 00179 */ 00180 doxygen_overloaded_function(template <...> void differenceOfExponentialEdgeImage) 00181 00182 template <class SrcIterator, class SrcAccessor, 00183 class DestIterator, class DestAccessor, 00184 class GradValue, class DestValue> 00185 void differenceOfExponentialEdgeImage( 00186 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 00187 DestIterator dul, DestAccessor da, 00188 double scale, GradValue gradient_threshold, DestValue edge_marker) 00189 { 00190 vigra_precondition(scale > 0, 00191 "differenceOfExponentialEdgeImage(): scale > 0 required."); 00192 00193 vigra_precondition(gradient_threshold > 0, 00194 "differenceOfExponentialEdgeImage(): " 00195 "gradient_threshold > 0 required."); 00196 00197 int w = slr.x - sul.x; 00198 int h = slr.y - sul.y; 00199 int x,y; 00200 00201 typedef typename 00202 NumericTraits<typename SrcAccessor::value_type>::RealPromote 00203 TMPTYPE; 00204 typedef BasicImage<TMPTYPE> TMPIMG; 00205 00206 TMPIMG tmp(w,h); 00207 TMPIMG smooth(w,h); 00208 00209 recursiveSmoothX(srcIterRange(sul, slr, sa), destImage(tmp), scale / 2.0); 00210 recursiveSmoothY(srcImageRange(tmp), destImage(tmp), scale / 2.0); 00211 00212 recursiveSmoothX(srcImageRange(tmp), destImage(smooth), scale); 00213 recursiveSmoothY(srcImageRange(smooth), destImage(smooth), scale); 00214 00215 typename TMPIMG::Iterator iy = smooth.upperLeft(); 00216 typename TMPIMG::Iterator ty = tmp.upperLeft(); 00217 DestIterator dy = dul; 00218 00219 static const Diff2D right(1, 0); 00220 static const Diff2D bottom(0, 1); 00221 00222 00223 TMPTYPE thresh = (gradient_threshold * gradient_threshold) * 00224 NumericTraits<TMPTYPE>::one(); 00225 TMPTYPE zero = NumericTraits<TMPTYPE>::zero(); 00226 00227 for(y=0; y<h-1; ++y, ++iy.y, ++ty.y, ++dy.y) 00228 { 00229 typename TMPIMG::Iterator ix = iy; 00230 typename TMPIMG::Iterator tx = ty; 00231 DestIterator dx = dy; 00232 00233 for(x=0; x<w-1; ++x, ++ix.x, ++tx.x, ++dx.x) 00234 { 00235 TMPTYPE diff = *tx - *ix; 00236 TMPTYPE gx = tx[right] - *tx; 00237 TMPTYPE gy = tx[bottom] - *tx; 00238 00239 if((gx * gx > thresh) && 00240 (diff * (tx[right] - ix[right]) < zero)) 00241 { 00242 if(gx < zero) 00243 { 00244 da.set(edge_marker, dx, right); 00245 } 00246 else 00247 { 00248 da.set(edge_marker, dx); 00249 } 00250 } 00251 if(((gy * gy > thresh) && 00252 (diff * (tx[bottom] - ix[bottom]) < zero))) 00253 { 00254 if(gy < zero) 00255 { 00256 da.set(edge_marker, dx, bottom); 00257 } 00258 else 00259 { 00260 da.set(edge_marker, dx); 00261 } 00262 } 00263 } 00264 } 00265 00266 typename TMPIMG::Iterator ix = iy; 00267 typename TMPIMG::Iterator tx = ty; 00268 DestIterator dx = dy; 00269 00270 for(x=0; x<w-1; ++x, ++ix.x, ++tx.x, ++dx.x) 00271 { 00272 TMPTYPE diff = *tx - *ix; 00273 TMPTYPE gx = tx[right] - *tx; 00274 00275 if((gx * gx > thresh) && 00276 (diff * (tx[right] - ix[right]) < zero)) 00277 { 00278 if(gx < zero) 00279 { 00280 da.set(edge_marker, dx, right); 00281 } 00282 else 00283 { 00284 da.set(edge_marker, dx); 00285 } 00286 } 00287 } 00288 } 00289 00290 template <class SrcIterator, class SrcAccessor, 00291 class DestIterator, class DestAccessor, 00292 class GradValue> 00293 inline 00294 void differenceOfExponentialEdgeImage( 00295 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 00296 DestIterator dul, DestAccessor da, 00297 double scale, GradValue gradient_threshold) 00298 { 00299 differenceOfExponentialEdgeImage(sul, slr, sa, dul, da, 00300 scale, gradient_threshold, 1); 00301 } 00302 00303 template <class SrcIterator, class SrcAccessor, 00304 class DestIterator, class DestAccessor, 00305 class GradValue, class DestValue> 00306 inline 00307 void differenceOfExponentialEdgeImage( 00308 triple<SrcIterator, SrcIterator, SrcAccessor> src, 00309 pair<DestIterator, DestAccessor> dest, 00310 double scale, GradValue gradient_threshold, 00311 DestValue edge_marker) 00312 { 00313 differenceOfExponentialEdgeImage(src.first, src.second, src.third, 00314 dest.first, dest.second, 00315 scale, gradient_threshold, 00316 edge_marker); 00317 } 00318 00319 template <class SrcIterator, class SrcAccessor, 00320 class DestIterator, class DestAccessor, 00321 class GradValue> 00322 inline 00323 void differenceOfExponentialEdgeImage( 00324 triple<SrcIterator, SrcIterator, SrcAccessor> src, 00325 pair<DestIterator, DestAccessor> dest, 00326 double scale, GradValue gradient_threshold) 00327 { 00328 differenceOfExponentialEdgeImage(src.first, src.second, src.third, 00329 dest.first, dest.second, 00330 scale, gradient_threshold, 1); 00331 } 00332 00333 /********************************************************/ 00334 /* */ 00335 /* differenceOfExponentialCrackEdgeImage */ 00336 /* */ 00337 /********************************************************/ 00338 00339 /** \brief Detect and mark edges in a crack edge image using the Shen/Castan zero-crossing detector. 00340 00341 This operator applies an exponential filter to the source image 00342 at the given <TT>scale</TT> and subtracts the result from the original image. 00343 Zero crossings are detected in the resulting difference image. Whenever the 00344 gradient at a zero crossing is greater than the given <TT>gradient_threshold</TT>, 00345 an edge point is marked (using <TT>edge_marker</TT>) in the destination image 00346 <i>between</i> the corresponding original pixels. Topologically, this means we 00347 must insert additional pixels between the original ones to represent the 00348 boundaries between the pixels (the so called zero- and one-cells, with the original 00349 pixels being two-cells). Within VIGRA, such an image is called \ref CrackEdgeImage. 00350 To allow insertion of the zero- and one-cells, the destination image must have twice the 00351 size of the original (precisely, <TT>(2*w-1)</TT> by <TT>(2*h-1)</TT> pixels). Then the algorithm 00352 proceeds as follows: 00353 00354 \code 00355 sign of difference image insert zero- and one-cells resulting edge points (*) 00356 00357 + . - . - . * . . . 00358 + - - . . . . . . * * * . 00359 + + - => + . + . - => . . . * . 00360 + + + . . . . . . . . * * 00361 + . + . + . . . . . 00362 \endcode 00363 00364 Thus the edge points are marked where they actually are - in between the pixels. 00365 An important property of the resulting edge image is that it conforms to the notion 00366 of well-composedness as defined by Latecki et al., i.e. connected regions and edges 00367 obtained by a subsequent \ref Labeling do not depend on 00368 whether 4- or 8-connectivity is used. 00369 The non-edge pixels (<TT>.</TT>) in the destination image remain unchanged. 00370 The result conformes to the requirements of a \ref CrackEdgeImage. It can be further 00371 improved by the post-processing operations \ref removeShortEdges() and 00372 \ref closeGapsInCrackEdgeImage(). 00373 00374 The source value type (<TT>SrcAccessor::value_type</TT>) must be a linear algebra, i.e. addition, 00375 subtraction and multiplication of the type with itself, and multiplication 00376 with double and 00377 \ref NumericTraits "NumericTraits" must 00378 be defined. In addition, this type must be less-comparable. 00379 00380 <b> Declarations:</b> 00381 00382 pass arguments explicitly: 00383 \code 00384 namespace vigra { 00385 template <class SrcIterator, class SrcAccessor, 00386 class DestIterator, class DestAccessor, 00387 class GradValue, 00388 class DestValue = DestAccessor::value_type> 00389 void differenceOfExponentialCrackEdgeImage( 00390 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 00391 DestIterator dul, DestAccessor da, 00392 double scale, GradValue gradient_threshold, 00393 DestValue edge_marker = NumericTraits<DestValue>::one()) 00394 } 00395 \endcode 00396 00397 use argument objects in conjunction with \ref ArgumentObjectFactories : 00398 \code 00399 namespace vigra { 00400 template <class SrcIterator, class SrcAccessor, 00401 class DestIterator, class DestAccessor, 00402 class GradValue, 00403 class DestValue = DestAccessor::value_type> 00404 void differenceOfExponentialCrackEdgeImage( 00405 triple<SrcIterator, SrcIterator, SrcAccessor> src, 00406 pair<DestIterator, DestAccessor> dest, 00407 double scale, GradValue gradient_threshold, 00408 DestValue edge_marker = NumericTraits<DestValue>::one()) 00409 } 00410 \endcode 00411 00412 <b> Usage:</b> 00413 00414 <b>\#include</b> <<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>><br> 00415 Namespace: vigra 00416 00417 \code 00418 vigra::BImage src(w,h), edges(2*w-1,2*h-1); 00419 00420 // empty edge image 00421 edges = 0; 00422 ... 00423 00424 // find edges at scale 0.8 with gradient larger than 4.0, mark with 1 00425 vigra::differenceOfExponentialCrackEdgeImage(srcImageRange(src), destImage(edges), 00426 0.8, 4.0, 1); 00427 \endcode 00428 00429 <b> Required Interface:</b> 00430 00431 \code 00432 SrcImageIterator src_upperleft, src_lowerright; 00433 DestImageIterator dest_upperleft; 00434 00435 SrcAccessor src_accessor; 00436 DestAccessor dest_accessor; 00437 00438 SrcAccessor::value_type u = src_accessor(src_upperleft); 00439 double d; 00440 GradValue gradient_threshold; 00441 00442 u = u + u 00443 u = u - u 00444 u = u * u 00445 u = d * u 00446 u < gradient_threshold 00447 00448 DestValue edge_marker; 00449 dest_accessor.set(edge_marker, dest_upperleft); 00450 \endcode 00451 00452 <b> Preconditions:</b> 00453 00454 \code 00455 scale > 0 00456 gradient_threshold > 0 00457 \endcode 00458 00459 The destination image must have twice the size of the source: 00460 \code 00461 w_dest = 2 * w_src - 1 00462 h_dest = 2 * h_src - 1 00463 \endcode 00464 */ 00465 doxygen_overloaded_function(template <...> void differenceOfExponentialCrackEdgeImage) 00466 00467 template <class SrcIterator, class SrcAccessor, 00468 class DestIterator, class DestAccessor, 00469 class GradValue, class DestValue> 00470 void differenceOfExponentialCrackEdgeImage( 00471 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 00472 DestIterator dul, DestAccessor da, 00473 double scale, GradValue gradient_threshold, 00474 DestValue edge_marker) 00475 { 00476 vigra_precondition(scale > 0, 00477 "differenceOfExponentialCrackEdgeImage(): scale > 0 required."); 00478 00479 vigra_precondition(gradient_threshold > 0, 00480 "differenceOfExponentialCrackEdgeImage(): " 00481 "gradient_threshold > 0 required."); 00482 00483 int w = slr.x - sul.x; 00484 int h = slr.y - sul.y; 00485 int x, y; 00486 00487 typedef typename 00488 NumericTraits<typename SrcAccessor::value_type>::RealPromote 00489 TMPTYPE; 00490 typedef BasicImage<TMPTYPE> TMPIMG; 00491 00492 TMPIMG tmp(w,h); 00493 TMPIMG smooth(w,h); 00494 00495 TMPTYPE zero = NumericTraits<TMPTYPE>::zero(); 00496 00497 static const Diff2D right(1,0); 00498 static const Diff2D bottom(0,1); 00499 static const Diff2D left(-1,0); 00500 static const Diff2D top(0,-1); 00501 00502 recursiveSmoothX(srcIterRange(sul, slr, sa), destImage(tmp), scale / 2.0); 00503 recursiveSmoothY(srcImageRange(tmp), destImage(tmp), scale / 2.0); 00504 00505 recursiveSmoothX(srcImageRange(tmp), destImage(smooth), scale); 00506 recursiveSmoothY(srcImageRange(smooth), destImage(smooth), scale); 00507 00508 typename TMPIMG::Iterator iy = smooth.upperLeft(); 00509 typename TMPIMG::Iterator ty = tmp.upperLeft(); 00510 DestIterator dy = dul; 00511 00512 TMPTYPE thresh = (gradient_threshold * gradient_threshold) * 00513 NumericTraits<TMPTYPE>::one(); 00514 00515 // find zero crossings above threshold 00516 for(y=0; y<h-1; ++y, ++iy.y, ++ty.y, dy.y+=2) 00517 { 00518 typename TMPIMG::Iterator ix = iy; 00519 typename TMPIMG::Iterator tx = ty; 00520 DestIterator dx = dy; 00521 00522 for(int x=0; x<w-1; ++x, ++ix.x, ++tx.x, dx.x+=2) 00523 { 00524 TMPTYPE diff = *tx - *ix; 00525 TMPTYPE gx = tx[right] - *tx; 00526 TMPTYPE gy = tx[bottom] - *tx; 00527 00528 if((gx * gx > thresh) && 00529 (diff * (tx[right] - ix[right]) < zero)) 00530 { 00531 da.set(edge_marker, dx, right); 00532 } 00533 if((gy * gy > thresh) && 00534 (diff * (tx[bottom] - ix[bottom]) < zero)) 00535 { 00536 da.set(edge_marker, dx, bottom); 00537 } 00538 } 00539 00540 TMPTYPE diff = *tx - *ix; 00541 TMPTYPE gy = tx[bottom] - *tx; 00542 00543 if((gy * gy > thresh) && 00544 (diff * (tx[bottom] - ix[bottom]) < zero)) 00545 { 00546 da.set(edge_marker, dx, bottom); 00547 } 00548 } 00549 00550 typename TMPIMG::Iterator ix = iy; 00551 typename TMPIMG::Iterator tx = ty; 00552 DestIterator dx = dy; 00553 00554 for(x=0; x<w-1; ++x, ++ix.x, ++tx.x, dx.x+=2) 00555 { 00556 TMPTYPE diff = *tx - *ix; 00557 TMPTYPE gx = tx[right] - *tx; 00558 00559 if((gx * gx > thresh) && 00560 (diff * (tx[right] - ix[right]) < zero)) 00561 { 00562 da.set(edge_marker, dx, right); 00563 } 00564 } 00565 00566 iy = smooth.upperLeft() + Diff2D(0,1); 00567 ty = tmp.upperLeft() + Diff2D(0,1); 00568 dy = dul + Diff2D(1,2); 00569 00570 static const Diff2D topleft(-1,-1); 00571 static const Diff2D topright(1,-1); 00572 static const Diff2D bottomleft(-1,1); 00573 static const Diff2D bottomright(1,1); 00574 00575 // find missing 1-cells below threshold (x-direction) 00576 for(y=0; y<h-2; ++y, ++iy.y, ++ty.y, dy.y+=2) 00577 { 00578 typename TMPIMG::Iterator ix = iy; 00579 typename TMPIMG::Iterator tx = ty; 00580 DestIterator dx = dy; 00581 00582 for(int x=0; x<w-2; ++x, ++ix.x, ++tx.x, dx.x+=2) 00583 { 00584 if(da(dx) == edge_marker) continue; 00585 00586 TMPTYPE diff = *tx - *ix; 00587 00588 if((diff * (tx[right] - ix[right]) < zero) && 00589 (((da(dx, bottomright) == edge_marker) && 00590 (da(dx, topleft) == edge_marker)) || 00591 ((da(dx, bottomleft) == edge_marker) && 00592 (da(dx, topright) == edge_marker)))) 00593 00594 { 00595 da.set(edge_marker, dx); 00596 } 00597 } 00598 } 00599 00600 iy = smooth.upperLeft() + Diff2D(1,0); 00601 ty = tmp.upperLeft() + Diff2D(1,0); 00602 dy = dul + Diff2D(2,1); 00603 00604 // find missing 1-cells below threshold (y-direction) 00605 for(y=0; y<h-2; ++y, ++iy.y, ++ty.y, dy.y+=2) 00606 { 00607 typename TMPIMG::Iterator ix = iy; 00608 typename TMPIMG::Iterator tx = ty; 00609 DestIterator dx = dy; 00610 00611 for(int x=0; x<w-2; ++x, ++ix.x, ++tx.x, dx.x+=2) 00612 { 00613 if(da(dx) == edge_marker) continue; 00614 00615 TMPTYPE diff = *tx - *ix; 00616 00617 if((diff * (tx[bottom] - ix[bottom]) < zero) && 00618 (((da(dx, bottomright) == edge_marker) && 00619 (da(dx, topleft) == edge_marker)) || 00620 ((da(dx, bottomleft) == edge_marker) && 00621 (da(dx, topright) == edge_marker)))) 00622 00623 { 00624 da.set(edge_marker, dx); 00625 } 00626 } 00627 } 00628 00629 dy = dul + Diff2D(1,1); 00630 00631 // find missing 0-cells 00632 for(y=0; y<h-1; ++y, dy.y+=2) 00633 { 00634 DestIterator dx = dy; 00635 00636 for(int x=0; x<w-1; ++x, dx.x+=2) 00637 { 00638 static const Diff2D dist[] = {right, top, left, bottom }; 00639 00640 int i; 00641 for(i=0; i<4; ++i) 00642 { 00643 if(da(dx, dist[i]) == edge_marker) break; 00644 } 00645 00646 if(i < 4) da.set(edge_marker, dx); 00647 } 00648 } 00649 } 00650 00651 template <class SrcIterator, class SrcAccessor, 00652 class DestIterator, class DestAccessor, 00653 class GradValue, class DestValue> 00654 inline 00655 void differenceOfExponentialCrackEdgeImage( 00656 triple<SrcIterator, SrcIterator, SrcAccessor> src, 00657 pair<DestIterator, DestAccessor> dest, 00658 double scale, GradValue gradient_threshold, 00659 DestValue edge_marker) 00660 { 00661 differenceOfExponentialCrackEdgeImage(src.first, src.second, src.third, 00662 dest.first, dest.second, 00663 scale, gradient_threshold, 00664 edge_marker); 00665 } 00666 00667 /********************************************************/ 00668 /* */ 00669 /* removeShortEdges */ 00670 /* */ 00671 /********************************************************/ 00672 00673 /** \brief Remove short edges from an edge image. 00674 00675 This algorithm can be applied as a post-processing operation of 00676 \ref differenceOfExponentialEdgeImage() and \ref differenceOfExponentialCrackEdgeImage(). 00677 It removes all edges that are shorter than <TT>min_edge_length</TT>. The corresponding 00678 pixels are set to the <TT>non_edge_marker</TT>. The idea behind this algorithms is 00679 that very short edges are probably caused by noise and don't represent interesting 00680 image structure. Technically, the algorithms executes a connected components labeling, 00681 so the image's value type must be equality comparable. 00682 00683 If the source image fulfills the requirements of a \ref CrackEdgeImage, 00684 it will still do so after application of this algorithm. 00685 00686 Note that this algorithm, unlike most other algorithms in VIGRA, operates in-place, 00687 i.e. on only one image. Also, the algorithm assumes that all non-edges pixels are already 00688 marked with the given <TT>non_edge_marker</TT> value. 00689 00690 <b> Declarations:</b> 00691 00692 pass arguments explicitly: 00693 \code 00694 namespace vigra { 00695 template <class Iterator, class Accessor, class SrcValue> 00696 void removeShortEdges( 00697 Iterator sul, Iterator slr, Accessor sa, 00698 int min_edge_length, SrcValue non_edge_marker) 00699 } 00700 \endcode 00701 00702 use argument objects in conjunction with \ref ArgumentObjectFactories : 00703 \code 00704 namespace vigra { 00705 template <class Iterator, class Accessor, class SrcValue> 00706 void removeShortEdges( 00707 triple<Iterator, Iterator, Accessor> src, 00708 int min_edge_length, SrcValue non_edge_marker) 00709 } 00710 \endcode 00711 00712 <b> Usage:</b> 00713 00714 <b>\#include</b> <<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>><br> 00715 Namespace: vigra 00716 00717 \code 00718 vigra::BImage src(w,h), edges(w,h); 00719 00720 // empty edge image 00721 edges = 0; 00722 ... 00723 00724 // find edges at scale 0.8 with gradient larger than 4.0, mark with 1 00725 vigra::differenceOfExponentialEdgeImage(srcImageRange(src), destImage(edges), 00726 0.8, 4.0, 1); 00727 00728 // zero edges shorter than 10 pixels 00729 vigra::removeShortEdges(srcImageRange(edges), 10, 0); 00730 \endcode 00731 00732 <b> Required Interface:</b> 00733 00734 \code 00735 SrcImageIterator src_upperleft, src_lowerright; 00736 DestImageIterator dest_upperleft; 00737 00738 SrcAccessor src_accessor; 00739 DestAccessor dest_accessor; 00740 00741 SrcAccessor::value_type u = src_accessor(src_upperleft); 00742 00743 u == u 00744 00745 SrcValue non_edge_marker; 00746 src_accessor.set(non_edge_marker, src_upperleft); 00747 \endcode 00748 */ 00749 doxygen_overloaded_function(template <...> void removeShortEdges) 00750 00751 template <class Iterator, class Accessor, class Value> 00752 void removeShortEdges( 00753 Iterator sul, Iterator slr, Accessor sa, 00754 unsigned int min_edge_length, Value non_edge_marker) 00755 { 00756 int w = slr.x - sul.x; 00757 int h = slr.y - sul.y; 00758 int x,y; 00759 00760 IImage labels(w, h); 00761 labels = 0; 00762 00763 int number_of_regions = 00764 labelImageWithBackground(srcIterRange(sul,slr,sa), 00765 destImage(labels), true, non_edge_marker); 00766 00767 ArrayOfRegionStatistics<FindROISize<int> > 00768 region_stats(number_of_regions); 00769 00770 inspectTwoImages(srcImageRange(labels), srcImage(labels), region_stats); 00771 00772 IImage::Iterator ly = labels.upperLeft(); 00773 Iterator oy = sul; 00774 00775 for(y=0; y<h; ++y, ++oy.y, ++ly.y) 00776 { 00777 Iterator ox(oy); 00778 IImage::Iterator lx(ly); 00779 00780 for(x=0; x<w; ++x, ++ox.x, ++lx.x) 00781 { 00782 if(sa(ox) == non_edge_marker) continue; 00783 if((region_stats[*lx].count) < min_edge_length) 00784 { 00785 sa.set(non_edge_marker, ox); 00786 } 00787 } 00788 } 00789 } 00790 00791 template <class Iterator, class Accessor, class Value> 00792 inline 00793 void removeShortEdges( 00794 triple<Iterator, Iterator, Accessor> src, 00795 unsigned int min_edge_length, Value non_edge_marker) 00796 { 00797 removeShortEdges(src.first, src.second, src.third, 00798 min_edge_length, non_edge_marker); 00799 } 00800 00801 /********************************************************/ 00802 /* */ 00803 /* closeGapsInCrackEdgeImage */ 00804 /* */ 00805 /********************************************************/ 00806 00807 /** \brief Close one-pixel wide gaps in a cell grid edge image. 00808 00809 This algorithm is typically applied as a post-processing operation of 00810 \ref differenceOfExponentialCrackEdgeImage(). The source image must fulfill 00811 the requirements of a \ref CrackEdgeImage, and will still do so after 00812 application of this algorithm. 00813 00814 It closes one pixel wide gaps in the edges resulting from this algorithm. 00815 Since these gaps are usually caused by zero crossing slightly below the gradient 00816 threshold used in edge detection, this algorithms acts like a weak hysteresis 00817 thresholding. The newly found edge pixels are marked with the given <TT>edge_marker</TT>. 00818 The image's value type must be equality comparable. 00819 00820 Note that this algorithm, unlike most other algorithms in VIGRA, operates in-place, 00821 i.e. on only one image. 00822 00823 <b> Declarations:</b> 00824 00825 pass arguments explicitly: 00826 \code 00827 namespace vigra { 00828 template <class SrcIterator, class SrcAccessor, class SrcValue> 00829 void closeGapsInCrackEdgeImage( 00830 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 00831 SrcValue edge_marker) 00832 } 00833 \endcode 00834 00835 use argument objects in conjunction with \ref ArgumentObjectFactories : 00836 \code 00837 namespace vigra { 00838 template <class SrcIterator, class SrcAccessor, class SrcValue> 00839 void closeGapsInCrackEdgeImage( 00840 triple<SrcIterator, SrcIterator, SrcAccessor> src, 00841 SrcValue edge_marker) 00842 } 00843 \endcode 00844 00845 <b> Usage:</b> 00846 00847 <b>\#include</b> <<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>><br> 00848 Namespace: vigra 00849 00850 \code 00851 vigra::BImage src(w,h), edges(2*w-1, 2*h-1); 00852 00853 // empty edge image 00854 edges = 0; 00855 ... 00856 00857 // find edges at scale 0.8 with gradient larger than 4.0, mark with 1 00858 vigra::differenceOfExponentialCrackEdgeImage(srcImageRange(src), destImage(edges), 00859 0.8, 4.0, 1); 00860 00861 // close gaps, mark with 1 00862 vigra::closeGapsInCrackEdgeImage(srcImageRange(edges), 1); 00863 00864 // zero edges shorter than 20 pixels 00865 vigra::removeShortEdges(srcImageRange(edges), 10, 0); 00866 \endcode 00867 00868 <b> Required Interface:</b> 00869 00870 \code 00871 SrcImageIterator src_upperleft, src_lowerright; 00872 00873 SrcAccessor src_accessor; 00874 DestAccessor dest_accessor; 00875 00876 SrcAccessor::value_type u = src_accessor(src_upperleft); 00877 00878 u == u 00879 u != u 00880 00881 SrcValue edge_marker; 00882 src_accessor.set(edge_marker, src_upperleft); 00883 \endcode 00884 */ 00885 doxygen_overloaded_function(template <...> void closeGapsInCrackEdgeImage) 00886 00887 template <class SrcIterator, class SrcAccessor, class SrcValue> 00888 void closeGapsInCrackEdgeImage( 00889 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 00890 SrcValue edge_marker) 00891 { 00892 int w = (slr.x - sul.x) / 2; 00893 int h = (slr.y - sul.y) / 2; 00894 int x, y; 00895 00896 int count1, count2, count3; 00897 00898 static const Diff2D right(1,0); 00899 static const Diff2D bottom(0,1); 00900 static const Diff2D left(-1,0); 00901 static const Diff2D top(0,-1); 00902 00903 static const Diff2D leftdist[] = { 00904 Diff2D(0, 0), Diff2D(-1, 1), Diff2D(-2, 0), Diff2D(-1, -1)}; 00905 static const Diff2D rightdist[] = { 00906 Diff2D(2, 0), Diff2D(1, 1), Diff2D(0, 0), Diff2D(1, -1)}; 00907 static const Diff2D topdist[] = { 00908 Diff2D(1, -1), Diff2D(0, 0), Diff2D(-1, -1), Diff2D(0, -2)}; 00909 static const Diff2D bottomdist[] = { 00910 Diff2D(1, 1), Diff2D(0, 2), Diff2D(-1, 1), Diff2D(0, 0)}; 00911 00912 int i; 00913 00914 SrcIterator sy = sul + Diff2D(0,1); 00915 SrcIterator sx; 00916 00917 // close 1-pixel wide gaps (x-direction) 00918 for(y=0; y<h; ++y, sy.y+=2) 00919 { 00920 sx = sy + Diff2D(2,0); 00921 00922 for(x=2; x<w; ++x, sx.x+=2) 00923 { 00924 if(sa(sx) == edge_marker) continue; 00925 00926 if(sa(sx, left) != edge_marker) continue; 00927 if(sa(sx, right) != edge_marker) continue; 00928 00929 count1 = 0; 00930 count2 = 0; 00931 count3 = 0; 00932 00933 for(i=0; i<4; ++i) 00934 { 00935 if(sa(sx, leftdist[i]) == edge_marker) 00936 { 00937 ++count1; 00938 count3 ^= 1 << i; 00939 } 00940 if(sa(sx, rightdist[i]) == edge_marker) 00941 { 00942 ++count2; 00943 count3 ^= 1 << i; 00944 } 00945 } 00946 00947 if(count1 <= 1 || count2 <= 1 || count3 == 15) 00948 { 00949 sa.set(edge_marker, sx); 00950 } 00951 } 00952 } 00953 00954 sy = sul + Diff2D(1,2); 00955 00956 // close 1-pixel wide gaps (y-direction) 00957 for(y=2; y<h; ++y, sy.y+=2) 00958 { 00959 sx = sy; 00960 00961 for(x=0; x<w; ++x, sx.x+=2) 00962 { 00963 if(sa(sx) == edge_marker) continue; 00964 00965 if(sa(sx, top) != edge_marker) continue; 00966 if(sa(sx, bottom) != edge_marker) continue; 00967 00968 count1 = 0; 00969 count2 = 0; 00970 count3 = 0; 00971 00972 for(i=0; i<4; ++i) 00973 { 00974 if(sa(sx, topdist[i]) == edge_marker) 00975 { 00976 ++count1; 00977 count3 ^= 1 << i; 00978 } 00979 if(sa(sx, bottomdist[i]) == edge_marker) 00980 { 00981 ++count2; 00982 count3 ^= 1 << i; 00983 } 00984 } 00985 00986 if(count1 <= 1 || count2 <= 1 || count3 == 15) 00987 { 00988 sa.set(edge_marker, sx); 00989 } 00990 } 00991 } 00992 } 00993 00994 template <class SrcIterator, class SrcAccessor, class SrcValue> 00995 inline 00996 void closeGapsInCrackEdgeImage( 00997 triple<SrcIterator, SrcIterator, SrcAccessor> src, 00998 SrcValue edge_marker) 00999 { 01000 closeGapsInCrackEdgeImage(src.first, src.second, src.third, 01001 edge_marker); 01002 } 01003 01004 /********************************************************/ 01005 /* */ 01006 /* beautifyCrackEdgeImage */ 01007 /* */ 01008 /********************************************************/ 01009 01010 /** \brief Beautify crack edge image for visualization. 01011 01012 This algorithm is applied as a post-processing operation of 01013 \ref differenceOfExponentialCrackEdgeImage(). The source image must fulfill 01014 the requirements of a \ref CrackEdgeImage, but will <b> not</b> do so after 01015 application of this algorithm. In particular, the algorithm removes zero-cells 01016 marked as edges to avoid staircase effects on diagonal lines like this: 01017 01018 \code 01019 original edge points (*) resulting edge points 01020 01021 . * . . . . * . . . 01022 . * * * . . . * . . 01023 . . . * . => . . . * . 01024 . . . * * . . . . * 01025 . . . . . . . . . . 01026 \endcode 01027 01028 Therfore, this algorithm should only be applied as a vizualization aid, i.e. 01029 for human inspection. The algorithm assumes that edges are marked with <TT>edge_marker</TT>, 01030 and background pixels with <TT>background_marker</TT>. The image's value type must be 01031 equality comparable. 01032 01033 Note that this algorithm, unlike most other algorithms in VIGRA, operates in-place, 01034 i.e. on only one image. 01035 01036 <b> Declarations:</b> 01037 01038 pass arguments explicitly: 01039 \code 01040 namespace vigra { 01041 template <class SrcIterator, class SrcAccessor, class SrcValue> 01042 void beautifyCrackEdgeImage( 01043 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01044 SrcValue edge_marker, SrcValue background_marker) 01045 } 01046 \endcode 01047 01048 use argument objects in conjunction with \ref ArgumentObjectFactories : 01049 \code 01050 namespace vigra { 01051 template <class SrcIterator, class SrcAccessor, class SrcValue> 01052 void beautifyCrackEdgeImage( 01053 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01054 SrcValue edge_marker, SrcValue background_marker) 01055 } 01056 \endcode 01057 01058 <b> Usage:</b> 01059 01060 <b>\#include</b> <<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>><br> 01061 Namespace: vigra 01062 01063 \code 01064 vigra::BImage src(w,h), edges(2*w-1, 2*h-1); 01065 01066 // empty edge image 01067 edges = 0; 01068 ... 01069 01070 // find edges at scale 0.8 with gradient larger than 4.0, mark with 1 01071 vigra::differenceOfExponentialCrackEdgeImage(srcImageRange(src), destImage(edges), 01072 0.8, 4.0, 1); 01073 01074 // beautify edge image for visualization 01075 vigra::beautifyCrackEdgeImage(destImageRange(edges), 1, 0); 01076 01077 // show to the user 01078 window.open(edges); 01079 \endcode 01080 01081 <b> Required Interface:</b> 01082 01083 \code 01084 SrcImageIterator src_upperleft, src_lowerright; 01085 01086 SrcAccessor src_accessor; 01087 DestAccessor dest_accessor; 01088 01089 SrcAccessor::value_type u = src_accessor(src_upperleft); 01090 01091 u == u 01092 u != u 01093 01094 SrcValue background_marker; 01095 src_accessor.set(background_marker, src_upperleft); 01096 \endcode 01097 */ 01098 doxygen_overloaded_function(template <...> void beautifyCrackEdgeImage) 01099 01100 template <class SrcIterator, class SrcAccessor, class SrcValue> 01101 void beautifyCrackEdgeImage( 01102 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01103 SrcValue edge_marker, SrcValue background_marker) 01104 { 01105 int w = (slr.x - sul.x) / 2; 01106 int h = (slr.y - sul.y) / 2; 01107 int x, y; 01108 01109 SrcIterator sy = sul + Diff2D(1,1); 01110 SrcIterator sx; 01111 01112 static const Diff2D right(1,0); 01113 static const Diff2D bottom(0,1); 01114 static const Diff2D left(-1,0); 01115 static const Diff2D top(0,-1); 01116 01117 // delete 0-cells at corners 01118 for(y=0; y<h; ++y, sy.y+=2) 01119 { 01120 sx = sy; 01121 01122 for(x=0; x<w; ++x, sx.x+=2) 01123 { 01124 if(sa(sx) != edge_marker) continue; 01125 01126 if(sa(sx, right) == edge_marker && sa(sx, left) == edge_marker) continue; 01127 if(sa(sx, bottom) == edge_marker && sa(sx, top) == edge_marker) continue; 01128 01129 sa.set(background_marker, sx); 01130 } 01131 } 01132 } 01133 01134 template <class SrcIterator, class SrcAccessor, class SrcValue> 01135 inline 01136 void beautifyCrackEdgeImage( 01137 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01138 SrcValue edge_marker, SrcValue background_marker) 01139 { 01140 beautifyCrackEdgeImage(src.first, src.second, src.third, 01141 edge_marker, background_marker); 01142 } 01143 01144 01145 /** Helper class that stores edgel attributes. 01146 */ 01147 class Edgel 01148 { 01149 public: 01150 /** The edgel's sub-pixel x coordinate. 01151 */ 01152 float x; 01153 01154 /** The edgel's sub-pixel y coordinate. 01155 */ 01156 float y; 01157 01158 /** The edgel's strength (magnitude of the gradient vector). 01159 */ 01160 float strength; 01161 01162 /** 01163 The edgel's orientation. This is the angle 01164 between the x-axis and the edge, so that the bright side of the 01165 edge is on the right. The angle is measured 01166 counter-clockwise in radians like this: 01167 01168 01169 \code 01170 01171 edgel axis 01172 \ (bright side) 01173 (dark \ 01174 side) \ /__ 01175 \\ \ orientation angle 01176 \ | 01177 +------------> x-axis 01178 | 01179 | 01180 | 01181 | 01182 y-axis V 01183 \endcode 01184 01185 So, for example a vertical edge with its dark side on the left 01186 has orientation PI/2, and a horizontal edge with dark side on top 01187 has orientation 0. Obviously, the edge's orientation changes 01188 by PI if the contrast is reversed. 01189 01190 */ 01191 float orientation; 01192 01193 Edgel() 01194 : x(0.0f), y(0.0f), strength(0.0f), orientation(0.0f) 01195 {} 01196 01197 Edgel(float ix, float iy, float is, float io) 01198 : x(ix), y(iy), strength(is), orientation(io) 01199 {} 01200 }; 01201 01202 template <class Image1, class Image2, class BackInsertable> 01203 void internalCannyFindEdgels(Image1 const & gx, 01204 Image1 const & gy, 01205 Image2 const & magnitude, 01206 BackInsertable & edgels) 01207 { 01208 typedef typename Image1::value_type PixelType; 01209 double t = 0.5 / VIGRA_CSTD::sin(M_PI/8.0); 01210 01211 for(int y=1; y<gx.height()-1; ++y) 01212 { 01213 for(int x=1; x<gx.width()-1; ++x) 01214 { 01215 PixelType gradx = gx(x,y); 01216 PixelType grady = gy(x,y); 01217 double mag = magnitude(x, y); 01218 if(mag == 0.0) 01219 continue; 01220 01221 int dx = (int)VIGRA_CSTD::floor(gradx*t/mag + 0.5); 01222 int dy = (int)VIGRA_CSTD::floor(grady*t/mag + 0.5); 01223 01224 int x1 = x - dx, 01225 x2 = x + dx, 01226 y1 = y - dy, 01227 y2 = y + dy; 01228 01229 PixelType m1 = magnitude(x1, y1); 01230 PixelType m3 = magnitude(x2, y2); 01231 01232 if(m1 < mag && m3 <= mag) 01233 { 01234 Edgel edgel; 01235 01236 // local maximum => quadratic interpolation of sub-pixel location 01237 PixelType del = (m1 - m3) / 2.0 / (m1 + m3 - 2.0*mag); 01238 edgel.x = x + dx*del; 01239 edgel.y = y + dy*del; 01240 edgel.strength = mag; 01241 double orientation = VIGRA_CSTD::atan2(-grady, gradx) - M_PI * 1.5; 01242 if(orientation < 0.0) 01243 orientation += 2.0*M_PI; 01244 edgel.orientation = orientation; 01245 edgels.push_back(edgel); 01246 } 01247 } 01248 } 01249 } 01250 01251 /********************************************************/ 01252 /* */ 01253 /* cannyEdgelList */ 01254 /* */ 01255 /********************************************************/ 01256 01257 /** \brief Simple implementation of Canny's edge detector. 01258 01259 This operator first calculates the gradient vector for each 01260 pixel of the image using first derivatives of a Gaussian at the 01261 given scale. Then a very simple non-maxima supression is performed: 01262 for each 3x3 neighborhood, it is determined whether the center pixel has 01263 larger gradient magnitude than its two neighbors in gradient direction 01264 (where the direction is rounded into octands). If this is the case, 01265 a new \ref Edgel is appended to the given vector of <TT>edgels</TT>. The subpixel 01266 edgel position is determined by fitting a parabola 01267 to the three gradient magnitude values 01268 mentioned above. The sub-pixel location of the parabola's tip 01269 and the gradient magnitude and direction (from the pixel center) 01270 are written in the newly created edgel. 01271 01272 <b> Declarations:</b> 01273 01274 pass arguments explicitly: 01275 \code 01276 namespace vigra { 01277 template <class SrcIterator, class SrcAccessor, class BackInsertable> 01278 void cannyEdgelList(SrcIterator ul, SrcIterator lr, SrcAccessor src, 01279 BackInsertable & edgels, double scale); 01280 } 01281 \endcode 01282 01283 use argument objects in conjunction with \ref ArgumentObjectFactories : 01284 \code 01285 namespace vigra { 01286 template <class SrcIterator, class SrcAccessor, class BackInsertable> 01287 void 01288 cannyEdgelList(triple<SrcIterator, SrcIterator, SrcAccessor> src, 01289 BackInsertable & edgels, double scale); 01290 } 01291 \endcode 01292 01293 <b> Usage:</b> 01294 01295 <b>\#include</b> <<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>><br> 01296 Namespace: vigra 01297 01298 \code 01299 vigra::BImage src(w,h); 01300 01301 // empty edgel list 01302 std::vector<vigra::Edgel> edgels; 01303 ... 01304 01305 // find edgels at scale 0.8 01306 vigra::cannyEdgelList(srcImageRange(src), edgels, 0.8); 01307 \endcode 01308 01309 <b> Required Interface:</b> 01310 01311 \code 01312 SrcImageIterator src_upperleft; 01313 SrcAccessor src_accessor; 01314 01315 src_accessor(src_upperleft); 01316 01317 BackInsertable edgels; 01318 edgels.push_back(Edgel()); 01319 \endcode 01320 01321 SrcAccessor::value_type must be a type convertible to float 01322 01323 <b> Preconditions:</b> 01324 01325 \code 01326 scale > 0 01327 \endcode 01328 */ 01329 doxygen_overloaded_function(template <...> void cannyEdgelList) 01330 01331 template <class SrcIterator, class SrcAccessor, class BackInsertable> 01332 void cannyEdgelList(SrcIterator ul, SrcIterator lr, SrcAccessor src, 01333 BackInsertable & edgels, double scale) 01334 { 01335 int w = lr.x - ul.x; 01336 int h = lr.y - ul.y; 01337 01338 // calculate image gradients 01339 typedef typename 01340 NumericTraits<typename SrcAccessor::value_type>::RealPromote 01341 TmpType; 01342 01343 BasicImage<TmpType> tmp(w,h), dx(w,h), dy(w,h); 01344 01345 Kernel1D<double> smooth, grad; 01346 01347 smooth.initGaussian(scale); 01348 grad.initGaussianDerivative(scale, 1); 01349 01350 separableConvolveX(srcIterRange(ul, lr, src), destImage(tmp), kernel1d(grad)); 01351 separableConvolveY(srcImageRange(tmp), destImage(dx), kernel1d(smooth)); 01352 01353 separableConvolveY(srcIterRange(ul, lr, src), destImage(tmp), kernel1d(grad)); 01354 separableConvolveX(srcImageRange(tmp), destImage(dy), kernel1d(smooth)); 01355 01356 combineTwoImages(srcImageRange(dx), srcImage(dy), destImage(tmp), 01357 MagnitudeFunctor<TmpType>()); 01358 01359 01360 // find edgels 01361 internalCannyFindEdgels(dx, dy, tmp, edgels); 01362 } 01363 01364 template <class SrcIterator, class SrcAccessor, class BackInsertable> 01365 inline void 01366 cannyEdgelList(triple<SrcIterator, SrcIterator, SrcAccessor> src, 01367 BackInsertable & edgels, double scale) 01368 { 01369 cannyEdgelList(src.first, src.second, src.third, edgels, scale); 01370 } 01371 01372 /********************************************************/ 01373 /* */ 01374 /* cannyEdgeImage */ 01375 /* */ 01376 /********************************************************/ 01377 01378 /** \brief Detect and mark edges in an edge image using Canny's algorithm. 01379 01380 This operator first calls \ref cannyEdgelList() to generate an 01381 edgel list for the given image. Then it scans this list and selects edgels 01382 whose strength is above the given <TT>gradient_threshold</TT>. For each of these 01383 edgels, the edgel's location is rounded to the nearest pixel, and that 01384 pixel marked with the given <TT>edge_marker</TT>. 01385 01386 <b> Declarations:</b> 01387 01388 pass arguments explicitly: 01389 \code 01390 namespace vigra { 01391 template <class SrcIterator, class SrcAccessor, 01392 class DestIterator, class DestAccessor, 01393 class GradValue, class DestValue> 01394 void cannyEdgeImage( 01395 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01396 DestIterator dul, DestAccessor da, 01397 double scale, GradValue gradient_threshold, DestValue edge_marker); 01398 } 01399 \endcode 01400 01401 use argument objects in conjunction with \ref ArgumentObjectFactories : 01402 \code 01403 namespace vigra { 01404 template <class SrcIterator, class SrcAccessor, 01405 class DestIterator, class DestAccessor, 01406 class GradValue, class DestValue> 01407 void cannyEdgeImage( 01408 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01409 pair<DestIterator, DestAccessor> dest, 01410 double scale, GradValue gradient_threshold, DestValue edge_marker); 01411 } 01412 \endcode 01413 01414 <b> Usage:</b> 01415 01416 <b>\#include</b> <<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>><br> 01417 Namespace: vigra 01418 01419 \code 01420 vigra::BImage src(w,h), edges(w,h); 01421 01422 // empty edge image 01423 edges = 0; 01424 ... 01425 01426 // find edges at scale 0.8 with gradient larger than 4.0, mark with 1 01427 vigra::cannyEdgeImage(srcImageRange(src), destImage(edges), 01428 0.8, 4.0, 1); 01429 \endcode 01430 01431 <b> Required Interface:</b> 01432 01433 see also: \ref cannyEdgelList(). 01434 01435 \code 01436 DestImageIterator dest_upperleft; 01437 DestAccessor dest_accessor; 01438 DestValue edge_marker; 01439 01440 dest_accessor.set(edge_marker, dest_upperleft, vigra::Diff2D(1,1)); 01441 \endcode 01442 01443 <b> Preconditions:</b> 01444 01445 \code 01446 scale > 0 01447 gradient_threshold > 0 01448 \endcode 01449 */ 01450 doxygen_overloaded_function(template <...> void cannyEdgeImage) 01451 01452 template <class SrcIterator, class SrcAccessor, 01453 class DestIterator, class DestAccessor, 01454 class GradValue, class DestValue> 01455 void cannyEdgeImage( 01456 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01457 DestIterator dul, DestAccessor da, 01458 double scale, GradValue gradient_threshold, DestValue edge_marker) 01459 { 01460 std::vector<Edgel> edgels; 01461 01462 cannyEdgelList(sul, slr, sa, edgels, scale); 01463 01464 for(unsigned int i=0; i<edgels.size(); ++i) 01465 { 01466 if(gradient_threshold < edgels[i].strength) 01467 { 01468 Diff2D pix((int)(edgels[i].x + 0.5), (int)(edgels[i].y + 0.5)); 01469 01470 da.set(edge_marker, dul, pix); 01471 } 01472 } 01473 } 01474 01475 template <class SrcIterator, class SrcAccessor, 01476 class DestIterator, class DestAccessor, 01477 class GradValue, class DestValue> 01478 inline void cannyEdgeImage( 01479 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01480 pair<DestIterator, DestAccessor> dest, 01481 double scale, GradValue gradient_threshold, DestValue edge_marker) 01482 { 01483 cannyEdgeImage(src.first, src.second, src.third, 01484 dest.first, dest.second, 01485 scale, gradient_threshold, edge_marker); 01486 } 01487 01488 /********************************************************/ 01489 01490 namespace detail { 01491 01492 template <class DestIterator> 01493 int neighborhoodConfiguration(DestIterator dul) 01494 { 01495 int v = 0; 01496 NeighborhoodCirculator<DestIterator, EightNeighborCode> c(dul, EightNeighborCode::SouthEast); 01497 for(int i=0; i<8; ++i, --c) 01498 { 01499 v = (v << 1) | ((*c != 0) ? 1 : 0); 01500 } 01501 01502 return v; 01503 } 01504 01505 template <class GradValue> 01506 struct SimplePoint 01507 { 01508 Diff2D point; 01509 GradValue grad; 01510 01511 SimplePoint(Diff2D const & p, GradValue g) 01512 : point(p), grad(g) 01513 {} 01514 01515 bool operator<(SimplePoint const & o) const 01516 { 01517 return grad < o.grad; 01518 } 01519 01520 bool operator>(SimplePoint const & o) const 01521 { 01522 return grad > o.grad; 01523 } 01524 }; 01525 01526 template <class SrcIterator, class SrcAccessor, 01527 class DestIterator, class DestAccessor, 01528 class GradValue, class DestValue> 01529 void cannyEdgeImageFromGrad( 01530 SrcIterator sul, SrcIterator slr, SrcAccessor grad, 01531 DestIterator dul, DestAccessor da, 01532 GradValue gradient_threshold, DestValue edge_marker) 01533 { 01534 typedef typename SrcAccessor::value_type PixelType; 01535 typedef typename NormTraits<PixelType>::SquaredNormType NormType; 01536 01537 NormType zero = NumericTraits<NormType>::zero(); 01538 double tan22_5 = M_SQRT2 - 1.0; 01539 typename NormTraits<GradValue>::SquaredNormType g2thresh = squaredNorm(gradient_threshold); 01540 01541 int w = slr.x - sul.x; 01542 int h = slr.y - sul.y; 01543 01544 sul += Diff2D(1,1); 01545 dul += Diff2D(1,1); 01546 Diff2D p(0,0); 01547 01548 for(int y = 1; y < h-1; ++y, ++sul.y, ++dul.y) 01549 { 01550 SrcIterator sx = sul; 01551 DestIterator dx = dul; 01552 for(int x = 1; x < w-1; ++x, ++sx.x, ++dx.x) 01553 { 01554 PixelType g = grad(sx); 01555 NormType g2n = squaredNorm(g); 01556 if(g2n < g2thresh) 01557 continue; 01558 01559 NormType g2n1, g2n3; 01560 // find out quadrant 01561 if(abs(g[1]) < tan22_5*abs(g[0])) 01562 { 01563 // north-south edge 01564 g2n1 = squaredNorm(grad(sx, Diff2D(-1, 0))); 01565 g2n3 = squaredNorm(grad(sx, Diff2D(1, 0))); 01566 } 01567 else if(abs(g[0]) < tan22_5*abs(g[1])) 01568 { 01569 // west-east edge 01570 g2n1 = squaredNorm(grad(sx, Diff2D(0, -1))); 01571 g2n3 = squaredNorm(grad(sx, Diff2D(0, 1))); 01572 } 01573 else if(g[0]*g[1] < zero) 01574 { 01575 // north-west-south-east edge 01576 g2n1 = squaredNorm(grad(sx, Diff2D(1, -1))); 01577 g2n3 = squaredNorm(grad(sx, Diff2D(-1, 1))); 01578 } 01579 else 01580 { 01581 // north-east-south-west edge 01582 g2n1 = squaredNorm(grad(sx, Diff2D(-1, -1))); 01583 g2n3 = squaredNorm(grad(sx, Diff2D(1, 1))); 01584 } 01585 01586 if(g2n1 < g2n && g2n3 <= g2n) 01587 { 01588 da.set(edge_marker, dx); 01589 } 01590 } 01591 } 01592 } 01593 01594 } // namespace detail 01595 01596 /********************************************************/ 01597 /* */ 01598 /* cannyEdgeImageWithThinning */ 01599 /* */ 01600 /********************************************************/ 01601 01602 /** \brief Detect and mark edges in an edge image using Canny's algorithm. 01603 01604 The input pixels of this algorithms must be vectors of length 2 (see Required Interface below). 01605 It first searches for all pixels whose gradient magnitude is larger 01606 than the given <tt>gradient_threshold</tt> and larger than the magnitude of its two neighbors 01607 in gradient direction (where these neighbors are determined by nearest neighbor 01608 interpolation, i.e. according to the octant where the gradient points into). 01609 The resulting edge pixel candidates are then subjected to topological thinning 01610 so that the remaining edge pixels can be linked into edgel chains with a provable, 01611 non-heuristic algorithm. Thinning is performed so that the pixels with highest gradient 01612 magnitude survive. Optionally, the outermost pixels are marked as edge pixels 01613 as well when <tt>addBorder</tt> is true. The remaining pixels will be marked in the destination 01614 image with the value of <tt>edge_marker</tt> (all non-edge pixels remain untouched). 01615 01616 <b> Declarations:</b> 01617 01618 pass arguments explicitly: 01619 \code 01620 namespace vigra { 01621 template <class SrcIterator, class SrcAccessor, 01622 class DestIterator, class DestAccessor, 01623 class GradValue, class DestValue> 01624 void cannyEdgeImageFromGradWithThinning( 01625 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01626 DestIterator dul, DestAccessor da, 01627 GradValue gradient_threshold, 01628 DestValue edge_marker, bool addBorder = true); 01629 } 01630 \endcode 01631 01632 use argument objects in conjunction with \ref ArgumentObjectFactories : 01633 \code 01634 namespace vigra { 01635 template <class SrcIterator, class SrcAccessor, 01636 class DestIterator, class DestAccessor, 01637 class GradValue, class DestValue> 01638 void cannyEdgeImageFromGradWithThinning( 01639 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01640 pair<DestIterator, DestAccessor> dest, 01641 GradValue gradient_threshold, 01642 DestValue edge_marker, bool addBorder = true); 01643 } 01644 \endcode 01645 01646 <b> Usage:</b> 01647 01648 <b>\#include</b> <<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>><br> 01649 Namespace: vigra 01650 01651 \code 01652 vigra::BImage src(w,h), edges(w,h); 01653 01654 vigra::FVector2Image grad(w,h); 01655 // compute the image gradient at scale 0.8 01656 vigra::gaussianGradient(srcImageRange(src), destImage(grad), 0.8); 01657 01658 // empty edge image 01659 edges = 0; 01660 // find edges gradient larger than 4.0, mark with 1, and add border 01661 vigra::cannyEdgeImageFromGradWithThinning(srcImageRange(grad), destImage(edges), 01662 4.0, 1, true); 01663 \endcode 01664 01665 <b> Required Interface:</b> 01666 01667 \code 01668 // the input pixel type must be a vector with two elements 01669 SrcImageIterator src_upperleft; 01670 SrcAccessor src_accessor; 01671 typedef SrcAccessor::value_type SrcPixel; 01672 typedef NormTraits<SrcPixel>::SquaredNormType SrcSquaredNormType; 01673 01674 SrcPixel g = src_accessor(src_upperleft); 01675 SrcPixel::value_type g0 = g[0]; 01676 SrcSquaredNormType gn = squaredNorm(g); 01677 01678 DestImageIterator dest_upperleft; 01679 DestAccessor dest_accessor; 01680 DestValue edge_marker; 01681 01682 dest_accessor.set(edge_marker, dest_upperleft, vigra::Diff2D(1,1)); 01683 \endcode 01684 01685 <b> Preconditions:</b> 01686 01687 \code 01688 gradient_threshold > 0 01689 \endcode 01690 */ 01691 doxygen_overloaded_function(template <...> void cannyEdgeImageFromGradWithThinning) 01692 01693 template <class SrcIterator, class SrcAccessor, 01694 class DestIterator, class DestAccessor, 01695 class GradValue, class DestValue> 01696 void cannyEdgeImageFromGradWithThinning( 01697 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01698 DestIterator dul, DestAccessor da, 01699 GradValue gradient_threshold, 01700 DestValue edge_marker, bool addBorder) 01701 { 01702 int w = slr.x - sul.x; 01703 int h = slr.y - sul.y; 01704 01705 BImage edgeImage(w, h, BImage::value_type(0)); 01706 BImage::traverser eul = edgeImage.upperLeft(); 01707 BImage::Accessor ea = edgeImage.accessor(); 01708 if(addBorder) 01709 initImageBorder(destImageRange(edgeImage), 1, 1); 01710 detail::cannyEdgeImageFromGrad(sul, slr, sa, eul, ea, gradient_threshold, 1); 01711 01712 static bool isSimplePoint[256] = { 01713 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 01714 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 01715 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 01716 1, 1, 1, 0, 0, 0, 1, 1, 1, 1, 0, 1, 0, 1, 0, 1, 0, 1, 01717 0, 0, 0, 0, 0, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1, 1, 01718 0, 1, 1, 0, 1, 0, 0, 1, 0, 1, 0, 1, 0, 1, 0, 0, 0, 0, 01719 0, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 0, 01720 1, 0, 0, 0, 0, 1, 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 0, 1, 01721 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 01722 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 01723 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 01724 0, 1, 0, 0, 0, 0, 0, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 0, 01725 1, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 0, 1, 0, 1, 0, 0, 01726 0, 0, 0, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 01727 1, 0, 1, 0 }; 01728 01729 eul += Diff2D(1,1); 01730 sul += Diff2D(1,1); 01731 int w2 = w-2; 01732 int h2 = h-2; 01733 01734 typedef detail::SimplePoint<GradValue> SP; 01735 // use std::greater becaus we need the smallest gradients at the top of the queue 01736 std::priority_queue<SP, std::vector<SP>, std::greater<SP> > pqueue; 01737 01738 Diff2D p(0,0); 01739 for(; p.y < h2; ++p.y) 01740 { 01741 for(p.x = 0; p.x < w2; ++p.x) 01742 { 01743 BImage::traverser e = eul + p; 01744 if(*e == 0) 01745 continue; 01746 int v = detail::neighborhoodConfiguration(e); 01747 if(isSimplePoint[v]) 01748 { 01749 pqueue.push(SP(p, norm(sa(sul+p)))); 01750 *e = 2; // remember that it is already in queue 01751 } 01752 } 01753 } 01754 01755 static const Diff2D dist[] = { Diff2D(-1,0), Diff2D(0,-1), 01756 Diff2D(1,0), Diff2D(0,1) }; 01757 01758 while(pqueue.size()) 01759 { 01760 p = pqueue.top().point; 01761 pqueue.pop(); 01762 01763 BImage::traverser e = eul + p; 01764 int v = detail::neighborhoodConfiguration(e); 01765 if(!isSimplePoint[v]) 01766 continue; // point may no longer be simple because its neighbors changed 01767 01768 *e = 0; // delete simple point 01769 01770 for(int i=0; i<4; ++i) 01771 { 01772 Diff2D pneu = p + dist[i]; 01773 if(pneu.x == -1 || pneu.y == -1 || pneu.x == w2 || pneu.y == h2) 01774 continue; // do not remove points at the border 01775 01776 BImage::traverser eneu = eul + pneu; 01777 if(*eneu == 1) // point is boundary and not yet in the queue 01778 { 01779 int v = detail::neighborhoodConfiguration(eneu); 01780 if(isSimplePoint[v]) 01781 { 01782 pqueue.push(SP(pneu, norm(sa(sul+pneu)))); 01783 *eneu = 2; // remember that it is already in queue 01784 } 01785 } 01786 } 01787 } 01788 01789 initImageIf(destIterRange(dul, dul+Diff2D(w,h), da), 01790 maskImage(edgeImage), edge_marker); 01791 } 01792 01793 template <class SrcIterator, class SrcAccessor, 01794 class DestIterator, class DestAccessor, 01795 class GradValue, class DestValue> 01796 inline void cannyEdgeImageFromGradWithThinning( 01797 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01798 pair<DestIterator, DestAccessor> dest, 01799 GradValue gradient_threshold, 01800 DestValue edge_marker, bool addBorder) 01801 { 01802 cannyEdgeImageFromGradWithThinning(src.first, src.second, src.third, 01803 dest.first, dest.second, 01804 gradient_threshold, edge_marker, addBorder); 01805 } 01806 01807 template <class SrcIterator, class SrcAccessor, 01808 class DestIterator, class DestAccessor, 01809 class GradValue, class DestValue> 01810 inline void cannyEdgeImageFromGradWithThinning( 01811 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01812 DestIterator dul, DestAccessor da, 01813 GradValue gradient_threshold, DestValue edge_marker) 01814 { 01815 cannyEdgeImageFromGradWithThinning(sul, slr, sa, 01816 dul, da, 01817 gradient_threshold, edge_marker, true); 01818 } 01819 01820 template <class SrcIterator, class SrcAccessor, 01821 class DestIterator, class DestAccessor, 01822 class GradValue, class DestValue> 01823 inline void cannyEdgeImageFromGradWithThinning( 01824 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01825 pair<DestIterator, DestAccessor> dest, 01826 GradValue gradient_threshold, DestValue edge_marker) 01827 { 01828 cannyEdgeImageFromGradWithThinning(src.first, src.second, src.third, 01829 dest.first, dest.second, 01830 gradient_threshold, edge_marker, true); 01831 } 01832 01833 /********************************************************/ 01834 /* */ 01835 /* cannyEdgeImageWithThinning */ 01836 /* */ 01837 /********************************************************/ 01838 01839 /** \brief Detect and mark edges in an edge image using Canny's algorithm. 01840 01841 This operator first calls \ref gaussianGradient() to compute the gradient of the input 01842 image, ad then \ref cannyEdgeImageFromGradWithThinning() to generate an 01843 edge image. See there for more detailed documentation. 01844 01845 <b> Declarations:</b> 01846 01847 pass arguments explicitly: 01848 \code 01849 namespace vigra { 01850 template <class SrcIterator, class SrcAccessor, 01851 class DestIterator, class DestAccessor, 01852 class GradValue, class DestValue> 01853 void cannyEdgeImageWithThinning( 01854 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01855 DestIterator dul, DestAccessor da, 01856 double scale, GradValue gradient_threshold, 01857 DestValue edge_marker, bool addBorder = true); 01858 } 01859 \endcode 01860 01861 use argument objects in conjunction with \ref ArgumentObjectFactories : 01862 \code 01863 namespace vigra { 01864 template <class SrcIterator, class SrcAccessor, 01865 class DestIterator, class DestAccessor, 01866 class GradValue, class DestValue> 01867 void cannyEdgeImageWithThinning( 01868 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01869 pair<DestIterator, DestAccessor> dest, 01870 double scale, GradValue gradient_threshold, 01871 DestValue edge_marker, bool addBorder = true); 01872 } 01873 \endcode 01874 01875 <b> Usage:</b> 01876 01877 <b>\#include</b> <<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>><br> 01878 Namespace: vigra 01879 01880 \code 01881 vigra::BImage src(w,h), edges(w,h); 01882 01883 // empty edge image 01884 edges = 0; 01885 ... 01886 01887 // find edges at scale 0.8 with gradient larger than 4.0, mark with 1, annd add border 01888 vigra::cannyEdgeImageWithThinning(srcImageRange(src), destImage(edges), 01889 0.8, 4.0, 1, true); 01890 \endcode 01891 01892 <b> Required Interface:</b> 01893 01894 see also: \ref cannyEdgelList(). 01895 01896 \code 01897 DestImageIterator dest_upperleft; 01898 DestAccessor dest_accessor; 01899 DestValue edge_marker; 01900 01901 dest_accessor.set(edge_marker, dest_upperleft, vigra::Diff2D(1,1)); 01902 \endcode 01903 01904 <b> Preconditions:</b> 01905 01906 \code 01907 scale > 0 01908 gradient_threshold > 0 01909 \endcode 01910 */ 01911 doxygen_overloaded_function(template <...> void cannyEdgeImageWithThinning) 01912 01913 template <class SrcIterator, class SrcAccessor, 01914 class DestIterator, class DestAccessor, 01915 class GradValue, class DestValue> 01916 void cannyEdgeImageWithThinning( 01917 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01918 DestIterator dul, DestAccessor da, 01919 double scale, GradValue gradient_threshold, 01920 DestValue edge_marker, bool addBorder) 01921 { 01922 // mark pixels that are higher than their neighbors in gradient direction 01923 typedef typename NumericTraits<typename SrcAccessor::value_type>::RealPromote TmpType; 01924 BasicImage<TinyVector<TmpType, 2> > grad(slr-sul); 01925 gaussianGradient(srcIterRange(sul, slr, sa), destImage(grad), scale); 01926 cannyEdgeImageFromGradWithThinning(srcImageRange(grad), destIter(dul, da), 01927 gradient_threshold, edge_marker, addBorder); 01928 } 01929 01930 template <class SrcIterator, class SrcAccessor, 01931 class DestIterator, class DestAccessor, 01932 class GradValue, class DestValue> 01933 inline void cannyEdgeImageWithThinning( 01934 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01935 pair<DestIterator, DestAccessor> dest, 01936 double scale, GradValue gradient_threshold, 01937 DestValue edge_marker, bool addBorder) 01938 { 01939 cannyEdgeImageWithThinning(src.first, src.second, src.third, 01940 dest.first, dest.second, 01941 scale, gradient_threshold, edge_marker, addBorder); 01942 } 01943 01944 template <class SrcIterator, class SrcAccessor, 01945 class DestIterator, class DestAccessor, 01946 class GradValue, class DestValue> 01947 inline void cannyEdgeImageWithThinning( 01948 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01949 DestIterator dul, DestAccessor da, 01950 double scale, GradValue gradient_threshold, DestValue edge_marker) 01951 { 01952 cannyEdgeImageWithThinning(sul, slr, sa, 01953 dul, da, 01954 scale, gradient_threshold, edge_marker, true); 01955 } 01956 01957 template <class SrcIterator, class SrcAccessor, 01958 class DestIterator, class DestAccessor, 01959 class GradValue, class DestValue> 01960 inline void cannyEdgeImageWithThinning( 01961 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01962 pair<DestIterator, DestAccessor> dest, 01963 double scale, GradValue gradient_threshold, DestValue edge_marker) 01964 { 01965 cannyEdgeImageWithThinning(src.first, src.second, src.third, 01966 dest.first, dest.second, 01967 scale, gradient_threshold, edge_marker, true); 01968 } 01969 01970 /********************************************************/ 01971 01972 template <class Image1, class Image2, class BackInsertable> 01973 void internalCannyFindEdgels3x3(Image1 const & grad, 01974 Image2 const & mask, 01975 BackInsertable & edgels) 01976 { 01977 typedef typename Image1::value_type PixelType; 01978 typedef typename PixelType::value_type ValueType; 01979 01980 for(int y=1; y<grad.height()-1; ++y) 01981 { 01982 for(int x=1; x<grad.width()-1; ++x) 01983 { 01984 if(!mask(x,y)) 01985 continue; 01986 01987 ValueType gradx = grad(x,y)[0]; 01988 ValueType grady = grad(x,y)[1]; 01989 double mag = hypot(gradx, grady); 01990 if(mag == 0.0) 01991 continue; 01992 double c = gradx / mag, 01993 s = grady / mag; 01994 01995 Matrix<double> ml(3,3), mr(3,1), l(3,1), r(3,1); 01996 l(0,0) = 1.0; 01997 01998 for(int yy = -1; yy <= 1; ++yy) 01999 { 02000 for(int xx = -1; xx <= 1; ++xx) 02001 { 02002 double u = c*xx + s*yy; 02003 double v = norm(grad(x+xx, y+yy)); 02004 l(1,0) = u; 02005 l(2,0) = u*u; 02006 ml += outer(l); 02007 mr += v*l; 02008 } 02009 } 02010 02011 linearSolve(ml, mr, r); 02012 02013 Edgel edgel; 02014 02015 // local maximum => quadratic interpolation of sub-pixel location 02016 double del = -r(1,0) / 2.0 / r(2,0); 02017 if(std::fabs(del) > 1.5) // don't move by more than about a pixel diameter 02018 del = 0.0; 02019 edgel.x = x + c*del; 02020 edgel.y = y + s*del; 02021 edgel.strength = mag; 02022 double orientation = VIGRA_CSTD::atan2(-grady, gradx) - M_PI * 1.5; 02023 if(orientation < 0.0) 02024 orientation += 2.0*M_PI; 02025 edgel.orientation = orientation; 02026 edgels.push_back(edgel); 02027 } 02028 } 02029 } 02030 02031 02032 /********************************************************/ 02033 /* */ 02034 /* cannyEdgelList3x3 */ 02035 /* */ 02036 /********************************************************/ 02037 02038 /** \brief Improved implementation of Canny's edge detector. 02039 02040 This operator first computes pixels which are crossed by the edge using 02041 cannyEdgeImageWithThinning(). The gradient magnitude in the 3x3 neighborhood of these 02042 pixels are then projected onto the normal of the edge (as determined 02043 by the gradient direction). The edgel's subpixel location is found by fitting a 02044 parabola through the 9 gradient values and determining the parabola's tip. 02045 A new \ref Edgel is appended to the given vector of <TT>edgels</TT>. Since the parabola 02046 is fitted to 9 points rather than 3 points as in cannyEdgelList(), the accuracy is higher. 02047 02048 <b> Declarations:</b> 02049 02050 pass arguments explicitly: 02051 \code 02052 namespace vigra { 02053 template <class SrcIterator, class SrcAccessor, class BackInsertable> 02054 void cannyEdgelList3x3(SrcIterator ul, SrcIterator lr, SrcAccessor src, 02055 BackInsertable & edgels, double scale); 02056 } 02057 \endcode 02058 02059 use argument objects in conjunction with \ref ArgumentObjectFactories : 02060 \code 02061 namespace vigra { 02062 template <class SrcIterator, class SrcAccessor, class BackInsertable> 02063 void 02064 cannyEdgelList3x3(triple<SrcIterator, SrcIterator, SrcAccessor> src, 02065 BackInsertable & edgels, double scale); 02066 } 02067 \endcode 02068 02069 <b> Usage:</b> 02070 02071 <b>\#include</b> <<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>><br> 02072 Namespace: vigra 02073 02074 \code 02075 vigra::BImage src(w,h); 02076 02077 // empty edgel list 02078 std::vector<vigra::Edgel> edgels; 02079 ... 02080 02081 // find edgels at scale 0.8 02082 vigra::cannyEdgelList3x3(srcImageRange(src), edgels, 0.8); 02083 \endcode 02084 02085 <b> Required Interface:</b> 02086 02087 \code 02088 SrcImageIterator src_upperleft; 02089 SrcAccessor src_accessor; 02090 02091 src_accessor(src_upperleft); 02092 02093 BackInsertable edgels; 02094 edgels.push_back(Edgel()); 02095 \endcode 02096 02097 SrcAccessor::value_type must be a type convertible to float 02098 02099 <b> Preconditions:</b> 02100 02101 \code 02102 scale > 0 02103 \endcode 02104 */ 02105 doxygen_overloaded_function(template <...> void cannyEdgelList3x3) 02106 02107 template <class SrcIterator, class SrcAccessor, class BackInsertable> 02108 void cannyEdgelList3x3(SrcIterator ul, SrcIterator lr, SrcAccessor src, 02109 BackInsertable & edgels, double scale) 02110 { 02111 int w = lr.x - ul.x; 02112 int h = lr.y - ul.y; 02113 02114 typedef typename NumericTraits<typename SrcAccessor::value_type>::RealPromote TmpType; 02115 BasicImage<TinyVector<TmpType, 2> > grad(lr-ul); 02116 gaussianGradient(srcIterRange(ul, lr, src), destImage(grad), scale); 02117 02118 UInt8Image edges(lr-ul); 02119 cannyEdgeImageFromGradWithThinning(srcImageRange(grad), destImage(edges), 02120 0.0, 1, false); 02121 02122 // find edgels 02123 internalCannyFindEdgels3x3(grad, edges, edgels); 02124 } 02125 02126 template <class SrcIterator, class SrcAccessor, class BackInsertable> 02127 inline void 02128 cannyEdgelList3x3(triple<SrcIterator, SrcIterator, SrcAccessor> src, 02129 BackInsertable & edgels, double scale) 02130 { 02131 cannyEdgelList3x3(src.first, src.second, src.third, edgels, scale); 02132 } 02133 02134 02135 02136 //@} 02137 02138 /** \page CrackEdgeImage Crack Edge Image 02139 02140 Crack edges are marked <i>between</i> the pixels of an image. 02141 A Crack Edge Image is an image that represents these edges. In order 02142 to accomodate the cracks, the Crack Edge Image must be twice as large 02143 as the original image (precisely (2*w - 1) by (2*h - 1)). A Crack Edge Image 02144 can easily be derived from a binary image or from the signs of the 02145 response of a Laplacean filter. Consider the following sketch, where 02146 <TT>+</TT> encodes the foreground, <TT>-</TT> the background, and 02147 <TT>*</TT> the resulting crack edges. 02148 02149 \code 02150 sign of difference image insert cracks resulting CrackEdgeImage 02151 02152 + . - . - . * . . . 02153 + - - . . . . . . * * * . 02154 + + - => + . + . - => . . . * . 02155 + + + . . . . . . . . * * 02156 + . + . + . . . . . 02157 \endcode 02158 02159 Starting from the original binary image (left), we insert crack pixels 02160 to get to the double-sized image (center). Finally, we mark all 02161 crack pixels whose non-crack neighbors have different signs as 02162 crack edge points, while all other pixels (crack and non-crack) become 02163 region pixels. 02164 02165 <b>Requirements on a Crack Edge Image:</b> 02166 02167 <ul> 02168 <li>Crack Edge Images have odd width and height. 02169 <li>Crack pixels have at least one odd coordinate. 02170 <li>Only crack pixels may be marked as edge points. 02171 <li>Crack pixels with two odd coordinates must be marked as edge points 02172 whenever any of their neighboring crack pixels was marked. 02173 </ul> 02174 02175 The last two requirements ensure that both edges and regions are 4-connected. 02176 Thus, 4-connectivity and 8-connectivity yield identical connected 02177 components in a Crack Edge Image (so called <i>well-composedness</i>). 02178 This ensures that Crack Edge Images have nice topological properties 02179 (cf. L. J. Latecki: "Well-Composed Sets", Academic Press, 2000). 02180 */ 02181 02182 02183 } // namespace vigra 02184 02185 #endif // VIGRA_EDGEDETECTION_HXX
© Ullrich Köthe (ullrich.koethe@iwr.uni-heidelberg.de) |
html generated using doxygen and Python
|