[ VIGRA Homepage | Class Index | Function Index | File Index | Main Page ]
![]() |
vigra/edgedetection.hxx | ![]() |
---|
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.4.0, Dec 21 2005 ) */ 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 /* koethe@informatik.uni-hamburg.de or */ 00012 /* vigra@kogs1.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 "vigra/utilities.hxx" 00046 #include "vigra/numerictraits.hxx" 00047 #include "vigra/stdimage.hxx" 00048 #include "vigra/stdimagefunctions.hxx" 00049 #include "vigra/recursiveconvolution.hxx" 00050 #include "vigra/separableconvolution.hxx" 00051 #include "vigra/labelimage.hxx" 00052 #include "vigra/mathutil.hxx" 00053 #include "vigra/pixelneighborhood.hxx" 00054 00055 00056 namespace vigra { 00057 00058 /** \addtogroup EdgeDetection Edge Detection 00059 Edge detectors based on first and second derivatives, 00060 and related post-processing. 00061 */ 00062 //@{ 00063 00064 /********************************************************/ 00065 /* */ 00066 /* differenceOfExponentialEdgeImage */ 00067 /* */ 00068 /********************************************************/ 00069 00070 /** \brief Detect and mark edges in an edge image using the Shen/Castan zero-crossing detector. 00071 00072 This operator applies an exponential filter to the source image 00073 at the given <TT>scale</TT> and subtracts the result from the original image. 00074 Zero crossings are detected in the resulting difference image. Whenever the 00075 gradient at a zero crossing is greater than the given <TT>gradient_threshold</TT>, 00076 an edge point is marked (using <TT>edge_marker</TT>) in the destination image on 00077 the darker side of the zero crossing (note that zero crossings occur 00078 <i>between</i> pixels). For example: 00079 00080 \code 00081 sign of difference image resulting edge points (*) 00082 00083 + - - * * . 00084 + + - => . * * 00085 + + + . . . 00086 \endcode 00087 00088 Non-edge pixels (<TT>.</TT>) remain untouched in the destination image. 00089 The result can be improved by the post-processing operation \ref removeShortEdges(). 00090 A more accurate edge placement can be achieved with the function 00091 \ref differenceOfExponentialCrackEdgeImage(). 00092 00093 The source value type 00094 (<TT>SrcAccessor::value_type</TT>) must be a linear algebra, i.e. addition, 00095 subtraction and multiplication of the type with itself, and multiplication 00096 with double and 00097 \ref NumericTraits "NumericTraits" must 00098 be defined. In addition, this type must be less-comparable. 00099 00100 <b> Declarations:</b> 00101 00102 pass arguments explicitly: 00103 \code 00104 namespace vigra { 00105 template <class SrcIterator, class SrcAccessor, 00106 class DestIterator, class DestAccessor, 00107 class GradValue, 00108 class DestValue = DestAccessor::value_type> 00109 void differenceOfExponentialEdgeImage( 00110 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 00111 DestIterator dul, DestAccessor da, 00112 double scale, GradValue gradient_threshold, 00113 DestValue edge_marker = NumericTraits<DestValue>::one()) 00114 } 00115 \endcode 00116 00117 use argument objects in conjunction with \ref ArgumentObjectFactories: 00118 \code 00119 namespace vigra { 00120 template <class SrcIterator, class SrcAccessor, 00121 class DestIterator, class DestAccessor, 00122 class GradValue, 00123 class DestValue = DestAccessor::value_type> 00124 inline 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 template <class SrcIterator, class SrcAccessor, 00181 class DestIterator, class DestAccessor, 00182 class GradValue, class DestValue> 00183 void differenceOfExponentialEdgeImage( 00184 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 00185 DestIterator dul, DestAccessor da, 00186 double scale, GradValue gradient_threshold, DestValue edge_marker) 00187 { 00188 vigra_precondition(scale > 0, 00189 "differenceOfExponentialEdgeImage(): scale > 0 required."); 00190 00191 vigra_precondition(gradient_threshold > 0, 00192 "differenceOfExponentialEdgeImage(): " 00193 "gradient_threshold > 0 required."); 00194 00195 int w = slr.x - sul.x; 00196 int h = slr.y - sul.y; 00197 int x,y; 00198 00199 typedef typename 00200 NumericTraits<typename SrcAccessor::value_type>::RealPromote 00201 TMPTYPE; 00202 typedef BasicImage<TMPTYPE> TMPIMG; 00203 00204 TMPIMG tmp(w,h); 00205 TMPIMG smooth(w,h); 00206 00207 recursiveSmoothX(srcIterRange(sul, slr, sa), destImage(tmp), scale / 2.0); 00208 recursiveSmoothY(srcImageRange(tmp), destImage(tmp), scale / 2.0); 00209 00210 recursiveSmoothX(srcImageRange(tmp), destImage(smooth), scale); 00211 recursiveSmoothY(srcImageRange(smooth), destImage(smooth), scale); 00212 00213 typename TMPIMG::Iterator iy = smooth.upperLeft(); 00214 typename TMPIMG::Iterator ty = tmp.upperLeft(); 00215 DestIterator dy = dul; 00216 00217 static const Diff2D right(1, 0); 00218 static const Diff2D bottom(0, 1); 00219 00220 00221 TMPTYPE thresh = (gradient_threshold * gradient_threshold) * 00222 NumericTraits<TMPTYPE>::one(); 00223 TMPTYPE zero = NumericTraits<TMPTYPE>::zero(); 00224 00225 for(y=0; y<h-1; ++y, ++iy.y, ++ty.y, ++dy.y) 00226 { 00227 typename TMPIMG::Iterator ix = iy; 00228 typename TMPIMG::Iterator tx = ty; 00229 DestIterator dx = dy; 00230 00231 for(x=0; x<w-1; ++x, ++ix.x, ++tx.x, ++dx.x) 00232 { 00233 TMPTYPE diff = *tx - *ix; 00234 TMPTYPE gx = tx[right] - *tx; 00235 TMPTYPE gy = tx[bottom] - *tx; 00236 00237 if((gx * gx > thresh) && 00238 (diff * (tx[right] - ix[right]) < zero)) 00239 { 00240 if(gx < zero) 00241 { 00242 da.set(edge_marker, dx, right); 00243 } 00244 else 00245 { 00246 da.set(edge_marker, dx); 00247 } 00248 } 00249 if(((gy * gy > thresh) && 00250 (diff * (tx[bottom] - ix[bottom]) < zero))) 00251 { 00252 if(gy < zero) 00253 { 00254 da.set(edge_marker, dx, bottom); 00255 } 00256 else 00257 { 00258 da.set(edge_marker, dx); 00259 } 00260 } 00261 } 00262 } 00263 00264 typename TMPIMG::Iterator ix = iy; 00265 typename TMPIMG::Iterator tx = ty; 00266 DestIterator dx = dy; 00267 00268 for(x=0; x<w-1; ++x, ++ix.x, ++tx.x, ++dx.x) 00269 { 00270 TMPTYPE diff = *tx - *ix; 00271 TMPTYPE gx = tx[right] - *tx; 00272 00273 if((gx * gx > thresh) && 00274 (diff * (tx[right] - ix[right]) < zero)) 00275 { 00276 if(gx < zero) 00277 { 00278 da.set(edge_marker, dx, right); 00279 } 00280 else 00281 { 00282 da.set(edge_marker, dx); 00283 } 00284 } 00285 } 00286 } 00287 00288 template <class SrcIterator, class SrcAccessor, 00289 class DestIterator, class DestAccessor, 00290 class GradValue> 00291 inline 00292 void differenceOfExponentialEdgeImage( 00293 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 00294 DestIterator dul, DestAccessor da, 00295 double scale, GradValue gradient_threshold) 00296 { 00297 differenceOfExponentialEdgeImage(sul, slr, sa, dul, da, 00298 scale, gradient_threshold, 1); 00299 } 00300 00301 template <class SrcIterator, class SrcAccessor, 00302 class DestIterator, class DestAccessor, 00303 class GradValue, class DestValue> 00304 inline 00305 void differenceOfExponentialEdgeImage( 00306 triple<SrcIterator, SrcIterator, SrcAccessor> src, 00307 pair<DestIterator, DestAccessor> dest, 00308 double scale, GradValue gradient_threshold, 00309 DestValue edge_marker) 00310 { 00311 differenceOfExponentialEdgeImage(src.first, src.second, src.third, 00312 dest.first, dest.second, 00313 scale, gradient_threshold, 00314 edge_marker); 00315 } 00316 00317 template <class SrcIterator, class SrcAccessor, 00318 class DestIterator, class DestAccessor, 00319 class GradValue> 00320 inline 00321 void differenceOfExponentialEdgeImage( 00322 triple<SrcIterator, SrcIterator, SrcAccessor> src, 00323 pair<DestIterator, DestAccessor> dest, 00324 double scale, GradValue gradient_threshold) 00325 { 00326 differenceOfExponentialEdgeImage(src.first, src.second, src.third, 00327 dest.first, dest.second, 00328 scale, gradient_threshold, 1); 00329 } 00330 00331 /********************************************************/ 00332 /* */ 00333 /* differenceOfExponentialCrackEdgeImage */ 00334 /* */ 00335 /********************************************************/ 00336 00337 /** \brief Detect and mark edges in a crack edge image using the Shen/Castan zero-crossing detector. 00338 00339 This operator applies an exponential filter to the source image 00340 at the given <TT>scale</TT> and subtracts the result from the original image. 00341 Zero crossings are detected in the resulting difference image. Whenever the 00342 gradient at a zero crossing is greater than the given <TT>gradient_threshold</TT>, 00343 an edge point is marked (using <TT>edge_marker</TT>) in the destination image 00344 <i>between</i> the corresponding original pixels. Topologically, this means we 00345 must insert additional pixels between the original ones to represent the 00346 boundaries between the pixels (the so called zero- and one-cells, with the original 00347 pixels being two-cells). Within VIGRA, such an image is called \ref CrackEdgeImage. 00348 To allow insertion of the zero- and one-cells, the destination image must have twice the 00349 size of the original (precisely, <TT>(2*w-1)</TT> by <TT>(2*h-1)</TT> pixels). Then the algorithm 00350 proceeds as follows: 00351 00352 \code 00353 sign of difference image insert zero- and one-cells resulting edge points (*) 00354 00355 + . - . - . * . . . 00356 + - - . . . . . . * * * . 00357 + + - => + . + . - => . . . * . 00358 + + + . . . . . . . . * * 00359 + . + . + . . . . . 00360 \endcode 00361 00362 Thus the edge points are marked where they actually are - in between the pixels. 00363 An important property of the resulting edge image is that it conforms to the notion 00364 of well-composedness as defined by Latecki et al., i.e. connected regions and edges 00365 obtained by a subsequent \ref Labeling do not depend on 00366 whether 4- or 8-connectivity is used. 00367 The non-edge pixels (<TT>.</TT>) in the destination image remain unchanged. 00368 The result conformes to the requirements of a \ref CrackEdgeImage. It can be further 00369 improved by the post-processing operations \ref removeShortEdges() and 00370 \ref closeGapsInCrackEdgeImage(). 00371 00372 The source value type (<TT>SrcAccessor::value_type</TT>) must be a linear algebra, i.e. addition, 00373 subtraction and multiplication of the type with itself, and multiplication 00374 with double and 00375 \ref NumericTraits "NumericTraits" must 00376 be defined. In addition, this type must be less-comparable. 00377 00378 <b> Declarations:</b> 00379 00380 pass arguments explicitly: 00381 \code 00382 namespace vigra { 00383 template <class SrcIterator, class SrcAccessor, 00384 class DestIterator, class DestAccessor, 00385 class GradValue, 00386 class DestValue = DestAccessor::value_type> 00387 void differenceOfExponentialCrackEdgeImage( 00388 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 00389 DestIterator dul, DestAccessor da, 00390 double scale, GradValue gradient_threshold, 00391 DestValue edge_marker = NumericTraits<DestValue>::one()) 00392 } 00393 \endcode 00394 00395 use argument objects in conjunction with \ref ArgumentObjectFactories: 00396 \code 00397 namespace vigra { 00398 template <class SrcIterator, class SrcAccessor, 00399 class DestIterator, class DestAccessor, 00400 class GradValue, 00401 class DestValue = DestAccessor::value_type> 00402 inline 00403 void differenceOfExponentialCrackEdgeImage( 00404 triple<SrcIterator, SrcIterator, SrcAccessor> src, 00405 pair<DestIterator, DestAccessor> dest, 00406 double scale, GradValue gradient_threshold, 00407 DestValue edge_marker = NumericTraits<DestValue>::one()) 00408 } 00409 \endcode 00410 00411 <b> Usage:</b> 00412 00413 <b>\#include</b> "<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>"<br> 00414 Namespace: vigra 00415 00416 \code 00417 vigra::BImage src(w,h), edges(2*w-1,2*h-1); 00418 00419 // empty edge image 00420 edges = 0; 00421 ... 00422 00423 // find edges at scale 0.8 with gradient larger than 4.0, mark with 1 00424 vigra::differenceOfExponentialCrackEdgeImage(srcImageRange(src), destImage(edges), 00425 0.8, 4.0, 1); 00426 \endcode 00427 00428 <b> Required Interface:</b> 00429 00430 \code 00431 SrcImageIterator src_upperleft, src_lowerright; 00432 DestImageIterator dest_upperleft; 00433 00434 SrcAccessor src_accessor; 00435 DestAccessor dest_accessor; 00436 00437 SrcAccessor::value_type u = src_accessor(src_upperleft); 00438 double d; 00439 GradValue gradient_threshold; 00440 00441 u = u + u 00442 u = u - u 00443 u = u * u 00444 u = d * u 00445 u < gradient_threshold 00446 00447 DestValue edge_marker; 00448 dest_accessor.set(edge_marker, dest_upperleft); 00449 \endcode 00450 00451 <b> Preconditions:</b> 00452 00453 \code 00454 scale > 0 00455 gradient_threshold > 0 00456 \endcode 00457 00458 The destination image must have twice the size of the source: 00459 \code 00460 w_dest = 2 * w_src - 1 00461 h_dest = 2 * h_src - 1 00462 \endcode 00463 */ 00464 template <class SrcIterator, class SrcAccessor, 00465 class DestIterator, class DestAccessor, 00466 class GradValue, class DestValue> 00467 void differenceOfExponentialCrackEdgeImage( 00468 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 00469 DestIterator dul, DestAccessor da, 00470 double scale, GradValue gradient_threshold, 00471 DestValue edge_marker) 00472 { 00473 vigra_precondition(scale > 0, 00474 "differenceOfExponentialCrackEdgeImage(): scale > 0 required."); 00475 00476 vigra_precondition(gradient_threshold > 0, 00477 "differenceOfExponentialCrackEdgeImage(): " 00478 "gradient_threshold > 0 required."); 00479 00480 int w = slr.x - sul.x; 00481 int h = slr.y - sul.y; 00482 int x, y; 00483 00484 typedef typename 00485 NumericTraits<typename SrcAccessor::value_type>::RealPromote 00486 TMPTYPE; 00487 typedef BasicImage<TMPTYPE> TMPIMG; 00488 00489 TMPIMG tmp(w,h); 00490 TMPIMG smooth(w,h); 00491 00492 TMPTYPE zero = NumericTraits<TMPTYPE>::zero(); 00493 00494 static const Diff2D right(1,0); 00495 static const Diff2D bottom(0,1); 00496 static const Diff2D left(-1,0); 00497 static const Diff2D top(0,-1); 00498 00499 recursiveSmoothX(srcIterRange(sul, slr, sa), destImage(tmp), scale / 2.0); 00500 recursiveSmoothY(srcImageRange(tmp), destImage(tmp), scale / 2.0); 00501 00502 recursiveSmoothX(srcImageRange(tmp), destImage(smooth), scale); 00503 recursiveSmoothY(srcImageRange(smooth), destImage(smooth), scale); 00504 00505 typename TMPIMG::Iterator iy = smooth.upperLeft(); 00506 typename TMPIMG::Iterator ty = tmp.upperLeft(); 00507 DestIterator dy = dul; 00508 00509 TMPTYPE thresh = (gradient_threshold * gradient_threshold) * 00510 NumericTraits<TMPTYPE>::one(); 00511 00512 // find zero crossings above threshold 00513 for(y=0; y<h-1; ++y, ++iy.y, ++ty.y, dy.y+=2) 00514 { 00515 typename TMPIMG::Iterator ix = iy; 00516 typename TMPIMG::Iterator tx = ty; 00517 DestIterator dx = dy; 00518 00519 for(int x=0; x<w-1; ++x, ++ix.x, ++tx.x, dx.x+=2) 00520 { 00521 TMPTYPE diff = *tx - *ix; 00522 TMPTYPE gx = tx[right] - *tx; 00523 TMPTYPE gy = tx[bottom] - *tx; 00524 00525 if((gx * gx > thresh) && 00526 (diff * (tx[right] - ix[right]) < zero)) 00527 { 00528 da.set(edge_marker, dx, right); 00529 } 00530 if((gy * gy > thresh) && 00531 (diff * (tx[bottom] - ix[bottom]) < zero)) 00532 { 00533 da.set(edge_marker, dx, bottom); 00534 } 00535 } 00536 00537 TMPTYPE diff = *tx - *ix; 00538 TMPTYPE gy = tx[bottom] - *tx; 00539 00540 if((gy * gy > thresh) && 00541 (diff * (tx[bottom] - ix[bottom]) < zero)) 00542 { 00543 da.set(edge_marker, dx, bottom); 00544 } 00545 } 00546 00547 typename TMPIMG::Iterator ix = iy; 00548 typename TMPIMG::Iterator tx = ty; 00549 DestIterator dx = dy; 00550 00551 for(x=0; x<w-1; ++x, ++ix.x, ++tx.x, dx.x+=2) 00552 { 00553 TMPTYPE diff = *tx - *ix; 00554 TMPTYPE gx = tx[right] - *tx; 00555 00556 if((gx * gx > thresh) && 00557 (diff * (tx[right] - ix[right]) < zero)) 00558 { 00559 da.set(edge_marker, dx, right); 00560 } 00561 } 00562 00563 iy = smooth.upperLeft() + Diff2D(0,1); 00564 ty = tmp.upperLeft() + Diff2D(0,1); 00565 dy = dul + Diff2D(1,2); 00566 00567 static const Diff2D topleft(-1,-1); 00568 static const Diff2D topright(1,-1); 00569 static const Diff2D bottomleft(-1,1); 00570 static const Diff2D bottomright(1,1); 00571 00572 // find missing 1-cells below threshold (x-direction) 00573 for(y=0; y<h-2; ++y, ++iy.y, ++ty.y, dy.y+=2) 00574 { 00575 typename TMPIMG::Iterator ix = iy; 00576 typename TMPIMG::Iterator tx = ty; 00577 DestIterator dx = dy; 00578 00579 for(int x=0; x<w-2; ++x, ++ix.x, ++tx.x, dx.x+=2) 00580 { 00581 if(da(dx) == edge_marker) continue; 00582 00583 TMPTYPE diff = *tx - *ix; 00584 00585 if((diff * (tx[right] - ix[right]) < zero) && 00586 (((da(dx, bottomright) == edge_marker) && 00587 (da(dx, topleft) == edge_marker)) || 00588 ((da(dx, bottomleft) == edge_marker) && 00589 (da(dx, topright) == edge_marker)))) 00590 00591 { 00592 da.set(edge_marker, dx); 00593 } 00594 } 00595 } 00596 00597 iy = smooth.upperLeft() + Diff2D(1,0); 00598 ty = tmp.upperLeft() + Diff2D(1,0); 00599 dy = dul + Diff2D(2,1); 00600 00601 // find missing 1-cells below threshold (y-direction) 00602 for(y=0; y<h-2; ++y, ++iy.y, ++ty.y, dy.y+=2) 00603 { 00604 typename TMPIMG::Iterator ix = iy; 00605 typename TMPIMG::Iterator tx = ty; 00606 DestIterator dx = dy; 00607 00608 for(int x=0; x<w-2; ++x, ++ix.x, ++tx.x, dx.x+=2) 00609 { 00610 if(da(dx) == edge_marker) continue; 00611 00612 TMPTYPE diff = *tx - *ix; 00613 00614 if((diff * (tx[bottom] - ix[bottom]) < zero) && 00615 (((da(dx, bottomright) == edge_marker) && 00616 (da(dx, topleft) == edge_marker)) || 00617 ((da(dx, bottomleft) == edge_marker) && 00618 (da(dx, topright) == edge_marker)))) 00619 00620 { 00621 da.set(edge_marker, dx); 00622 } 00623 } 00624 } 00625 00626 dy = dul + Diff2D(1,1); 00627 00628 // find missing 0-cells 00629 for(y=0; y<h-1; ++y, dy.y+=2) 00630 { 00631 DestIterator dx = dy; 00632 00633 for(int x=0; x<w-1; ++x, dx.x+=2) 00634 { 00635 static const Diff2D dist[] = {right, top, left, bottom }; 00636 00637 int i; 00638 for(i=0; i<4; ++i) 00639 { 00640 if(da(dx, dist[i]) == edge_marker) break; 00641 } 00642 00643 if(i < 4) da.set(edge_marker, dx); 00644 } 00645 } 00646 } 00647 00648 template <class SrcIterator, class SrcAccessor, 00649 class DestIterator, class DestAccessor, 00650 class GradValue, class DestValue> 00651 inline 00652 void differenceOfExponentialCrackEdgeImage( 00653 triple<SrcIterator, SrcIterator, SrcAccessor> src, 00654 pair<DestIterator, DestAccessor> dest, 00655 double scale, GradValue gradient_threshold, 00656 DestValue edge_marker) 00657 { 00658 differenceOfExponentialCrackEdgeImage(src.first, src.second, src.third, 00659 dest.first, dest.second, 00660 scale, gradient_threshold, 00661 edge_marker); 00662 } 00663 00664 /********************************************************/ 00665 /* */ 00666 /* removeShortEdges */ 00667 /* */ 00668 /********************************************************/ 00669 00670 /** \brief Remove short edges from an edge image. 00671 00672 This algorithm can be applied as a post-processing operation of 00673 \ref differenceOfExponentialEdgeImage() and \ref differenceOfExponentialCrackEdgeImage(). 00674 It removes all edges that are shorter than <TT>min_edge_length</TT>. The corresponding 00675 pixels are set to the <TT>non_edge_marker</TT>. The idea behind this algorithms is 00676 that very short edges are probably caused by noise and don't represent interesting 00677 image structure. Technically, the algorithms executes a connected components labeling, 00678 so the image's value type must be equality comparable. 00679 00680 If the source image fulfills the requirements of a \ref CrackEdgeImage, 00681 it will still do so after application of this algorithm. 00682 00683 Note that this algorithm, unlike most other algorithms in VIGRA, operates in-place, 00684 i.e. on only one image. Also, the algorithm assumes that all non-edges pixels are already 00685 marked with the given <TT>non_edge_marker</TT> value. 00686 00687 <b> Declarations:</b> 00688 00689 pass arguments explicitly: 00690 \code 00691 namespace vigra { 00692 template <class Iterator, class Accessor, class SrcValue> 00693 void removeShortEdges( 00694 Iterator sul, Iterator slr, Accessor sa, 00695 int min_edge_length, SrcValue non_edge_marker) 00696 } 00697 \endcode 00698 00699 use argument objects in conjunction with \ref ArgumentObjectFactories: 00700 \code 00701 namespace vigra { 00702 template <class Iterator, class Accessor, class SrcValue> 00703 inline 00704 void removeShortEdges( 00705 triple<Iterator, Iterator, Accessor> src, 00706 int min_edge_length, SrcValue non_edge_marker) 00707 } 00708 \endcode 00709 00710 <b> Usage:</b> 00711 00712 <b>\#include</b> "<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>"<br> 00713 Namespace: vigra 00714 00715 \code 00716 vigra::BImage src(w,h), edges(w,h); 00717 00718 // empty edge image 00719 edges = 0; 00720 ... 00721 00722 // find edges at scale 0.8 with gradient larger than 4.0, mark with 1 00723 vigra::differenceOfExponentialEdgeImage(srcImageRange(src), destImage(edges), 00724 0.8, 4.0, 1); 00725 00726 // zero edges shorter than 10 pixels 00727 vigra::removeShortEdges(srcImageRange(edges), 10, 0); 00728 \endcode 00729 00730 <b> Required Interface:</b> 00731 00732 \code 00733 SrcImageIterator src_upperleft, src_lowerright; 00734 DestImageIterator dest_upperleft; 00735 00736 SrcAccessor src_accessor; 00737 DestAccessor dest_accessor; 00738 00739 SrcAccessor::value_type u = src_accessor(src_upperleft); 00740 00741 u == u 00742 00743 SrcValue non_edge_marker; 00744 src_accessor.set(non_edge_marker, src_upperleft); 00745 \endcode 00746 */ 00747 template <class Iterator, class Accessor, class Value> 00748 void removeShortEdges( 00749 Iterator sul, Iterator slr, Accessor sa, 00750 unsigned int min_edge_length, Value non_edge_marker) 00751 { 00752 int w = slr.x - sul.x; 00753 int h = slr.y - sul.y; 00754 int x,y; 00755 00756 IImage labels(w, h); 00757 labels = 0; 00758 00759 int number_of_regions = 00760 labelImageWithBackground(srcIterRange(sul,slr,sa), 00761 destImage(labels), true, non_edge_marker); 00762 00763 ArrayOfRegionStatistics<FindROISize<int> > 00764 region_stats(number_of_regions); 00765 00766 inspectTwoImages(srcImageRange(labels), srcImage(labels), region_stats); 00767 00768 IImage::Iterator ly = labels.upperLeft(); 00769 Iterator oy = sul; 00770 00771 for(y=0; y<h; ++y, ++oy.y, ++ly.y) 00772 { 00773 Iterator ox(oy); 00774 IImage::Iterator lx(ly); 00775 00776 for(x=0; x<w; ++x, ++ox.x, ++lx.x) 00777 { 00778 if(sa(ox) == non_edge_marker) continue; 00779 if((region_stats[*lx].count) < min_edge_length) 00780 { 00781 sa.set(non_edge_marker, ox); 00782 } 00783 } 00784 } 00785 } 00786 00787 template <class Iterator, class Accessor, class Value> 00788 inline 00789 void removeShortEdges( 00790 triple<Iterator, Iterator, Accessor> src, 00791 unsigned int min_edge_length, Value non_edge_marker) 00792 { 00793 removeShortEdges(src.first, src.second, src.third, 00794 min_edge_length, non_edge_marker); 00795 } 00796 00797 /********************************************************/ 00798 /* */ 00799 /* closeGapsInCrackEdgeImage */ 00800 /* */ 00801 /********************************************************/ 00802 00803 /** \brief Close one-pixel wide gaps in a cell grid edge image. 00804 00805 This algorithm is typically applied as a post-processing operation of 00806 \ref differenceOfExponentialCrackEdgeImage(). The source image must fulfill 00807 the requirements of a \ref CrackEdgeImage, and will still do so after 00808 application of this algorithm. 00809 00810 It closes one pixel wide gaps in the edges resulting from this algorithm. 00811 Since these gaps are usually caused by zero crossing slightly below the gradient 00812 threshold used in edge detection, this algorithms acts like a weak hysteresis 00813 thresholding. The newly found edge pixels are marked with the given <TT>edge_marker</TT>. 00814 The image's value type must be equality comparable. 00815 00816 Note that this algorithm, unlike most other algorithms in VIGRA, operates in-place, 00817 i.e. on only one image. 00818 00819 <b> Declarations:</b> 00820 00821 pass arguments explicitly: 00822 \code 00823 namespace vigra { 00824 template <class SrcIterator, class SrcAccessor, class SrcValue> 00825 void closeGapsInCrackEdgeImage( 00826 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 00827 SrcValue edge_marker) 00828 } 00829 \endcode 00830 00831 use argument objects in conjunction with \ref ArgumentObjectFactories: 00832 \code 00833 namespace vigra { 00834 template <class SrcIterator, class SrcAccessor, class SrcValue> 00835 inline 00836 void closeGapsInCrackEdgeImage( 00837 triple<SrcIterator, SrcIterator, SrcAccessor> src, 00838 SrcValue edge_marker) 00839 } 00840 \endcode 00841 00842 <b> Usage:</b> 00843 00844 <b>\#include</b> "<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>"<br> 00845 Namespace: vigra 00846 00847 \code 00848 vigra::BImage src(w,h), edges(2*w-1, 2*h-1); 00849 00850 // empty edge image 00851 edges = 0; 00852 ... 00853 00854 // find edges at scale 0.8 with gradient larger than 4.0, mark with 1 00855 vigra::differenceOfExponentialCrackEdgeImage(srcImageRange(src), destImage(edges), 00856 0.8, 4.0, 1); 00857 00858 // close gaps, mark with 1 00859 vigra::closeGapsInCrackEdgeImage(srcImageRange(edges), 1); 00860 00861 // zero edges shorter than 20 pixels 00862 vigra::removeShortEdges(srcImageRange(edges), 10, 0); 00863 \endcode 00864 00865 <b> Required Interface:</b> 00866 00867 \code 00868 SrcImageIterator src_upperleft, src_lowerright; 00869 00870 SrcAccessor src_accessor; 00871 DestAccessor dest_accessor; 00872 00873 SrcAccessor::value_type u = src_accessor(src_upperleft); 00874 00875 u == u 00876 u != u 00877 00878 SrcValue edge_marker; 00879 src_accessor.set(edge_marker, src_upperleft); 00880 \endcode 00881 */ 00882 template <class SrcIterator, class SrcAccessor, class SrcValue> 00883 void closeGapsInCrackEdgeImage( 00884 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 00885 SrcValue edge_marker) 00886 { 00887 int w = (slr.x - sul.x) / 2; 00888 int h = (slr.y - sul.y) / 2; 00889 int x, y; 00890 00891 int count1, count2, count3; 00892 00893 static const Diff2D right(1,0); 00894 static const Diff2D bottom(0,1); 00895 static const Diff2D left(-1,0); 00896 static const Diff2D top(0,-1); 00897 00898 static const Diff2D leftdist[] = { 00899 Diff2D(0, 0), Diff2D(-1, 1), Diff2D(-2, 0), Diff2D(-1, -1)}; 00900 static const Diff2D rightdist[] = { 00901 Diff2D(2, 0), Diff2D(1, 1), Diff2D(0, 0), Diff2D(1, -1)}; 00902 static const Diff2D topdist[] = { 00903 Diff2D(1, -1), Diff2D(0, 0), Diff2D(-1, -1), Diff2D(0, -2)}; 00904 static const Diff2D bottomdist[] = { 00905 Diff2D(1, 1), Diff2D(0, 2), Diff2D(-1, 1), Diff2D(0, 0)}; 00906 00907 int i; 00908 00909 SrcIterator sy = sul + Diff2D(0,1); 00910 SrcIterator sx; 00911 00912 // close 1-pixel wide gaps (x-direction) 00913 for(y=0; y<h; ++y, sy.y+=2) 00914 { 00915 sx = sy + Diff2D(2,0); 00916 00917 for(x=2; x<w; ++x, sx.x+=2) 00918 { 00919 if(sa(sx) == edge_marker) continue; 00920 00921 if(sa(sx, left) != edge_marker) continue; 00922 if(sa(sx, right) != edge_marker) continue; 00923 00924 count1 = 0; 00925 count2 = 0; 00926 count3 = 0; 00927 00928 for(i=0; i<4; ++i) 00929 { 00930 if(sa(sx, leftdist[i]) == edge_marker) 00931 { 00932 ++count1; 00933 count3 ^= 1 << i; 00934 } 00935 if(sa(sx, rightdist[i]) == edge_marker) 00936 { 00937 ++count2; 00938 count3 ^= 1 << i; 00939 } 00940 } 00941 00942 if(count1 <= 1 || count2 <= 1 || count3 == 15) 00943 { 00944 sa.set(edge_marker, sx); 00945 } 00946 } 00947 } 00948 00949 sy = sul + Diff2D(1,2); 00950 00951 // close 1-pixel wide gaps (y-direction) 00952 for(y=2; y<h; ++y, sy.y+=2) 00953 { 00954 sx = sy; 00955 00956 for(x=0; x<w; ++x, sx.x+=2) 00957 { 00958 if(sa(sx) == edge_marker) continue; 00959 00960 if(sa(sx, top) != edge_marker) continue; 00961 if(sa(sx, bottom) != edge_marker) continue; 00962 00963 count1 = 0; 00964 count2 = 0; 00965 count3 = 0; 00966 00967 for(i=0; i<4; ++i) 00968 { 00969 if(sa(sx, topdist[i]) == edge_marker) 00970 { 00971 ++count1; 00972 count3 ^= 1 << i; 00973 } 00974 if(sa(sx, bottomdist[i]) == edge_marker) 00975 { 00976 ++count2; 00977 count3 ^= 1 << i; 00978 } 00979 } 00980 00981 if(count1 <= 1 || count2 <= 1 || count3 == 15) 00982 { 00983 sa.set(edge_marker, sx); 00984 } 00985 } 00986 } 00987 } 00988 00989 template <class SrcIterator, class SrcAccessor, class SrcValue> 00990 inline 00991 void closeGapsInCrackEdgeImage( 00992 triple<SrcIterator, SrcIterator, SrcAccessor> src, 00993 SrcValue edge_marker) 00994 { 00995 closeGapsInCrackEdgeImage(src.first, src.second, src.third, 00996 edge_marker); 00997 } 00998 00999 /********************************************************/ 01000 /* */ 01001 /* beautifyCrackEdgeImage */ 01002 /* */ 01003 /********************************************************/ 01004 01005 /** \brief Beautify crack edge image for visualization. 01006 01007 This algorithm is applied as a post-processing operation of 01008 \ref differenceOfExponentialCrackEdgeImage(). The source image must fulfill 01009 the requirements of a \ref CrackEdgeImage, but will <b> not</b> do so after 01010 application of this algorithm. In particular, the algorithm removes zero-cells 01011 marked as edges to avoid staircase effects on diagonal lines like this: 01012 01013 \code 01014 original edge points (*) resulting edge points 01015 01016 . * . . . . * . . . 01017 . * * * . . . * . . 01018 . . . * . => . . . * . 01019 . . . * * . . . . * 01020 . . . . . . . . . . 01021 \endcode 01022 01023 Therfore, this algorithm should only be applied as a vizualization aid, i.e. 01024 for human inspection. The algorithm assumes that edges are marked with <TT>edge_marker</TT>, 01025 and background pixels with <TT>background_marker</TT>. The image's value type must be 01026 equality comparable. 01027 01028 Note that this algorithm, unlike most other algorithms in VIGRA, operates in-place, 01029 i.e. on only one image. 01030 01031 <b> Declarations:</b> 01032 01033 pass arguments explicitly: 01034 \code 01035 namespace vigra { 01036 template <class SrcIterator, class SrcAccessor, class SrcValue> 01037 void beautifyCrackEdgeImage( 01038 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01039 SrcValue edge_marker, SrcValue background_marker) 01040 } 01041 \endcode 01042 01043 use argument objects in conjunction with \ref ArgumentObjectFactories: 01044 \code 01045 namespace vigra { 01046 template <class SrcIterator, class SrcAccessor, class SrcValue> 01047 inline 01048 void beautifyCrackEdgeImage( 01049 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01050 SrcValue edge_marker, SrcValue background_marker) 01051 } 01052 \endcode 01053 01054 <b> Usage:</b> 01055 01056 <b>\#include</b> "<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>"<br> 01057 Namespace: vigra 01058 01059 \code 01060 vigra::BImage src(w,h), edges(2*w-1, 2*h-1); 01061 01062 // empty edge image 01063 edges = 0; 01064 ... 01065 01066 // find edges at scale 0.8 with gradient larger than 4.0, mark with 1 01067 vigra::differenceOfExponentialCrackEdgeImage(srcImageRange(src), destImage(edges), 01068 0.8, 4.0, 1); 01069 01070 // beautify edge image for visualization 01071 vigra::beautifyCrackEdgeImage(destImageRange(edges), 1, 0); 01072 01073 // show to the user 01074 window.open(edges); 01075 \endcode 01076 01077 <b> Required Interface:</b> 01078 01079 \code 01080 SrcImageIterator src_upperleft, src_lowerright; 01081 01082 SrcAccessor src_accessor; 01083 DestAccessor dest_accessor; 01084 01085 SrcAccessor::value_type u = src_accessor(src_upperleft); 01086 01087 u == u 01088 u != u 01089 01090 SrcValue background_marker; 01091 src_accessor.set(background_marker, src_upperleft); 01092 \endcode 01093 */ 01094 template <class SrcIterator, class SrcAccessor, class SrcValue> 01095 void beautifyCrackEdgeImage( 01096 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01097 SrcValue edge_marker, SrcValue background_marker) 01098 { 01099 int w = (slr.x - sul.x) / 2; 01100 int h = (slr.y - sul.y) / 2; 01101 int x, y; 01102 01103 SrcIterator sy = sul + Diff2D(1,1); 01104 SrcIterator sx; 01105 01106 static const Diff2D right(1,0); 01107 static const Diff2D bottom(0,1); 01108 static const Diff2D left(-1,0); 01109 static const Diff2D top(0,-1); 01110 01111 // delete 0-cells at corners 01112 for(y=0; y<h; ++y, sy.y+=2) 01113 { 01114 sx = sy; 01115 01116 for(x=0; x<w; ++x, sx.x+=2) 01117 { 01118 if(sa(sx) != edge_marker) continue; 01119 01120 if(sa(sx, right) == edge_marker && sa(sx, left) == edge_marker) continue; 01121 if(sa(sx, bottom) == edge_marker && sa(sx, top) == edge_marker) continue; 01122 01123 sa.set(background_marker, sx); 01124 } 01125 } 01126 } 01127 01128 template <class SrcIterator, class SrcAccessor, class SrcValue> 01129 inline 01130 void beautifyCrackEdgeImage( 01131 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01132 SrcValue edge_marker, SrcValue background_marker) 01133 { 01134 beautifyCrackEdgeImage(src.first, src.second, src.third, 01135 edge_marker, background_marker); 01136 } 01137 01138 01139 /** Helper class that stores edgel attributes. 01140 */ 01141 class Edgel 01142 { 01143 public: 01144 /** The edgel's sub-pixel x coordinate. 01145 */ 01146 float x; 01147 01148 /** The edgel's sub-pixel y coordinate. 01149 */ 01150 float y; 01151 01152 /** The edgel's strength (magnitude of the gradient vector). 01153 */ 01154 float strength; 01155 01156 /** 01157 The edgel's orientation. This is the angle 01158 between the x-axis and the edge, so that the bright side of the 01159 edge is on the right. The angle is measured 01160 counter-clockwise in radians like this: 01161 01162 01163 \code 01164 01165 edgel axis 01166 \ (bright side) 01167 (dark \ 01168 side) \ /__ 01169 \\ \ orientation angle 01170 \ | 01171 +------------> x-axis 01172 | 01173 | 01174 | 01175 | 01176 y-axis V 01177 \endcode 01178 01179 So, for example a vertical edge with its dark side on the left 01180 has orientation PI/2, and a horizontal edge with dark side on top 01181 has orientation 0. Obviously, the edge's orientation changes 01182 by PI if the contrast is reversed. 01183 01184 */ 01185 float orientation; 01186 01187 Edgel() 01188 : x(0.0f), y(0.0f), strength(0.0f), orientation(0.0f) 01189 {} 01190 01191 Edgel(float ix, float iy, float is, float io) 01192 : x(ix), y(iy), strength(is), orientation(io) 01193 {} 01194 }; 01195 01196 template <class Image1, class Image2> 01197 void internalCannyFindEdgels(Image1 const & dx, 01198 Image1 const & dy, 01199 Image2 const & magnitude, 01200 std::vector<Edgel> & edgels) 01201 { 01202 typedef typename Image1::value_type PixelType; 01203 01204 PixelType zero = NumericTraits<PixelType>::zero(); 01205 double tan22_5 = M_SQRT2 - 1.0; 01206 01207 for(int y=1; y<dx.height()-1; ++y) 01208 { 01209 for(int x=1; x<dx.width()-1; ++x) 01210 { 01211 bool maximum_found = false; 01212 Edgel edgel; 01213 01214 PixelType gradx = dx(x,y); 01215 PixelType grady = dy(x,y); 01216 01217 // find out quadrant 01218 if(abs(grady) < tan22_5*abs(gradx)) 01219 { 01220 // north-south edge 01221 PixelType m1 = magnitude(x-1, y); 01222 PixelType m2 = magnitude(x, y); 01223 PixelType m3 = magnitude(x+1, y); 01224 01225 if(m1 < m2 && m3 <= m2) 01226 { 01227 edgel.y = y; 01228 01229 // local maximum => quadratic interpolation of sub-pixel location 01230 PixelType del = (m1 - m3) / 2.0; 01231 del /= (m1 + m3 - 2.0*m2); 01232 edgel.x = x + del; 01233 edgel.strength = m2; 01234 01235 maximum_found = true; 01236 } 01237 } 01238 else if(abs(gradx) < tan22_5*abs(grady)) 01239 { 01240 // west-east edge 01241 PixelType m1 = magnitude(x, y-1); 01242 PixelType m2 = magnitude(x, y); 01243 PixelType m3 = magnitude(x, y+1); 01244 01245 if(m1 < m2 && m3 <= m2) 01246 { 01247 edgel.x = x; 01248 01249 // local maximum => quadratic interpolation of sub-pixel location 01250 PixelType del = (m1 - m3) / 2.0; 01251 del /= (m1 + m3 - 2.0*m2); 01252 edgel.y = y + del; 01253 edgel.strength = m2; 01254 01255 maximum_found = true; 01256 } 01257 } 01258 else if(gradx*grady < zero) 01259 { 01260 // north-west-south-east edge 01261 PixelType m1 = magnitude(x+1, y-1); 01262 PixelType m2 = magnitude(x, y); 01263 PixelType m3 = magnitude(x-1, y+1); 01264 01265 if(m1 < m2 && m3 <= m2) 01266 { 01267 // local maximum => quadratic interpolation of sub-pixel location 01268 PixelType del = (m1 - m3) / 2.0; 01269 del /= (m1 + m3 - 2.0*m2); 01270 edgel.x = x - del; 01271 edgel.y = y + del; 01272 edgel.strength = m2; 01273 01274 maximum_found = true; 01275 } 01276 } 01277 else 01278 { 01279 // north-east-south-west edge 01280 PixelType m1 = magnitude(x-1, y-1); 01281 PixelType m2 = magnitude(x, y); 01282 PixelType m3 = magnitude(x+1, y+1); 01283 01284 if(m1 < m2 && m3 <= m2) 01285 { 01286 // local maximum => quadratic interpolation of sub-pixel location 01287 PixelType del = (m1 - m3) / 2.0; 01288 del /= (m1 + m3 - 2.0*m2); 01289 edgel.x = x + del; 01290 edgel.y = y + del; 01291 edgel.strength = m2; 01292 01293 maximum_found = true; 01294 } 01295 } 01296 01297 if(maximum_found) 01298 { 01299 double orientation = VIGRA_CSTD::atan2(-grady, gradx) - M_PI * 1.5; 01300 if(orientation < 0.0) 01301 orientation += 2.0*M_PI; 01302 edgel.orientation = orientation; 01303 edgels.push_back(edgel); 01304 } 01305 } 01306 } 01307 } 01308 01309 /********************************************************/ 01310 /* */ 01311 /* cannyEdgelList */ 01312 /* */ 01313 /********************************************************/ 01314 01315 /** \brief Simple implementation of Canny's edge detector. 01316 01317 This operator first calculates the gradient vector for each 01318 pixel of the image using first derivatives of a Gaussian at the 01319 given scale. Then a very simple non-maxima supression is performed: 01320 for each 3x3 neighborhood, it is determined whether the center pixel has 01321 larger gradient magnitude than its two neighbors in gradient direction 01322 (where the direction is rounded into octands). If this is the case, 01323 a new \ref Edgel is appended to the given vector of <TT>edgels</TT>. The subpixel 01324 edgel position is determined by fitting a parabola 01325 to the three gradient magnitude values 01326 mentioned above. The sub-pixel location of the parabola's tip 01327 and the gradient magnitude and direction are written in the newly created edgel. 01328 01329 <b> Declarations:</b> 01330 01331 pass arguments explicitly: 01332 \code 01333 namespace vigra { 01334 template <class SrcIterator, class SrcAccessor> 01335 void cannyEdgelList(SrcIterator ul, SrcIterator lr, SrcAccessor src, 01336 std::vector<Edgel> & edgels, double scale); 01337 } 01338 \endcode 01339 01340 use argument objects in conjunction with \ref ArgumentObjectFactories: 01341 \code 01342 namespace vigra { 01343 template <class SrcIterator, class SrcAccessor> 01344 void 01345 cannyEdgelList(triple<SrcIterator, SrcIterator, SrcAccessor> src, 01346 std::vector<Edgel> & edgels, double scale); 01347 } 01348 \endcode 01349 01350 <b> Usage:</b> 01351 01352 <b>\#include</b> "<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>"<br> 01353 Namespace: vigra 01354 01355 \code 01356 vigra::BImage src(w,h); 01357 01358 // empty edgel list 01359 std::vector<vigra::Edgel> edgels; 01360 ... 01361 01362 // find edgels at scale 0.8 01363 vigra::cannyEdgelList(srcImageRange(src), edgels, 0.8); 01364 \endcode 01365 01366 <b> Required Interface:</b> 01367 01368 \code 01369 SrcImageIterator src_upperleft; 01370 SrcAccessor src_accessor; 01371 01372 src_accessor(src_upperleft); 01373 \endcode 01374 01375 SrcAccessor::value_type must be a type convertible to float 01376 01377 <b> Preconditions:</b> 01378 01379 \code 01380 scale > 0 01381 \endcode 01382 */ 01383 template <class SrcIterator, class SrcAccessor> 01384 void cannyEdgelList(SrcIterator ul, SrcIterator lr, SrcAccessor src, 01385 std::vector<Edgel> & edgels, double scale) 01386 { 01387 int w = lr.x - ul.x; 01388 int h = lr.y - ul.y; 01389 01390 // calculate image gradients 01391 typedef typename 01392 NumericTraits<typename SrcAccessor::value_type>::RealPromote 01393 TmpType; 01394 01395 BasicImage<TmpType> tmp(w,h), dx(w,h), dy(w,h); 01396 01397 Kernel1D<double> smooth, grad; 01398 01399 smooth.initGaussian(scale); 01400 grad.initGaussianDerivative(scale, 1); 01401 01402 separableConvolveX(srcIterRange(ul, lr, src), destImage(tmp), kernel1d(grad)); 01403 separableConvolveY(srcImageRange(tmp), destImage(dx), kernel1d(smooth)); 01404 01405 separableConvolveY(srcIterRange(ul, lr, src), destImage(tmp), kernel1d(grad)); 01406 separableConvolveX(srcImageRange(tmp), destImage(dy), kernel1d(smooth)); 01407 01408 combineTwoImages(srcImageRange(dx), srcImage(dy), destImage(tmp), 01409 MagnitudeFunctor<TmpType>()); 01410 01411 01412 // find edgels 01413 internalCannyFindEdgels(dx, dy, tmp, edgels); 01414 } 01415 01416 template <class SrcIterator, class SrcAccessor> 01417 inline void 01418 cannyEdgelList(triple<SrcIterator, SrcIterator, SrcAccessor> src, 01419 std::vector<Edgel> & edgels, double scale) 01420 { 01421 cannyEdgelList(src.first, src.second, src.third, edgels, scale); 01422 } 01423 01424 /********************************************************/ 01425 /* */ 01426 /* cannyEdgeImage */ 01427 /* */ 01428 /********************************************************/ 01429 01430 /** \brief Detect and mark edges in an edge image using Canny's algorithm. 01431 01432 This operator first calls \ref cannyEdgelList() to generate an 01433 edgel list for the given image. Then it scans this list and selects edgels 01434 whose strength is above the given <TT>gradient_threshold</TT>. For each of these 01435 edgels, the edgel's location is rounded to the nearest pixel, and that 01436 pixel marked with the given <TT>edge_marker</TT>. 01437 01438 <b> Declarations:</b> 01439 01440 pass arguments explicitly: 01441 \code 01442 namespace vigra { 01443 template <class SrcIterator, class SrcAccessor, 01444 class DestIterator, class DestAccessor, 01445 class GradValue, class DestValue> 01446 void cannyEdgeImage( 01447 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01448 DestIterator dul, DestAccessor da, 01449 double scale, GradValue gradient_threshold, DestValue edge_marker); 01450 } 01451 \endcode 01452 01453 use argument objects in conjunction with \ref ArgumentObjectFactories: 01454 \code 01455 namespace vigra { 01456 template <class SrcIterator, class SrcAccessor, 01457 class DestIterator, class DestAccessor, 01458 class GradValue, class DestValue> 01459 inline void cannyEdgeImage( 01460 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01461 pair<DestIterator, DestAccessor> dest, 01462 double scale, GradValue gradient_threshold, DestValue edge_marker); 01463 } 01464 \endcode 01465 01466 <b> Usage:</b> 01467 01468 <b>\#include</b> "<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>"<br> 01469 Namespace: vigra 01470 01471 \code 01472 vigra::BImage src(w,h), edges(w,h); 01473 01474 // empty edge image 01475 edges = 0; 01476 ... 01477 01478 // find edges at scale 0.8 with gradient larger than 4.0, mark with 1 01479 vigra::cannyEdgeImage(srcImageRange(src), destImage(edges), 01480 0.8, 4.0, 1); 01481 \endcode 01482 01483 <b> Required Interface:</b> 01484 01485 see also: \ref cannyEdgelList(). 01486 01487 \code 01488 DestImageIterator dest_upperleft; 01489 DestAccessor dest_accessor; 01490 DestValue edge_marker; 01491 01492 dest_accessor.set(edge_marker, dest_upperleft, vigra::Diff2D(1,1)); 01493 \endcode 01494 01495 <b> Preconditions:</b> 01496 01497 \code 01498 scale > 0 01499 gradient_threshold > 0 01500 \endcode 01501 */ 01502 template <class SrcIterator, class SrcAccessor, 01503 class DestIterator, class DestAccessor, 01504 class GradValue, class DestValue> 01505 void cannyEdgeImage( 01506 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01507 DestIterator dul, DestAccessor da, 01508 double scale, GradValue gradient_threshold, DestValue edge_marker) 01509 { 01510 std::vector<Edgel> edgels; 01511 01512 cannyEdgelList(sul, slr, sa, edgels, scale); 01513 01514 for(unsigned int i=0; i<edgels.size(); ++i) 01515 { 01516 if(gradient_threshold < edgels[i].strength) 01517 { 01518 Diff2D pix((int)(edgels[i].x + 0.5), (int)(edgels[i].y + 0.5)); 01519 01520 da.set(edge_marker, dul, pix); 01521 } 01522 } 01523 } 01524 01525 template <class SrcIterator, class SrcAccessor, 01526 class DestIterator, class DestAccessor, 01527 class GradValue, class DestValue> 01528 inline void cannyEdgeImage( 01529 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01530 pair<DestIterator, DestAccessor> dest, 01531 double scale, GradValue gradient_threshold, DestValue edge_marker) 01532 { 01533 cannyEdgeImage(src.first, src.second, src.third, 01534 dest.first, dest.second, 01535 scale, gradient_threshold, edge_marker); 01536 } 01537 01538 /********************************************************/ 01539 01540 namespace detail { 01541 01542 template <class DestIterator> 01543 int neighborhoodConfiguration(DestIterator dul) 01544 { 01545 int v = 0; 01546 NeighborhoodCirculator<DestIterator, EightNeighborCode> c(dul, EightNeighborCode::SouthEast); 01547 for(int i=0; i<8; ++i, --c) 01548 { 01549 v = (v << 1) | ((*c != 0) ? 1 : 0); 01550 } 01551 01552 return v; 01553 } 01554 01555 template <class GradValue> 01556 struct SimplePoint 01557 { 01558 Diff2D point; 01559 GradValue grad; 01560 01561 SimplePoint(Diff2D const & p, GradValue g) 01562 : point(p), grad(g) 01563 {} 01564 01565 bool operator<(SimplePoint const & o) const 01566 { 01567 return grad < o.grad; 01568 } 01569 }; 01570 01571 } // namespace detail 01572 01573 /********************************************************/ 01574 /* */ 01575 /* cannyEdgeImageWithThinning */ 01576 /* */ 01577 /********************************************************/ 01578 01579 /** \brief Detect and mark edges in an edge image using Canny's algorithm. 01580 01581 This operator first calls \ref cannyEdgeImage() to generate an 01582 edge image. The resulting edge pixels are then subjected to topological thinning 01583 so that the remaining edge pixels can be linked into edgel chains with a provable, 01584 non-heuristic algorithm. Optionally, the outermost pixels are marked as edge pixels 01585 as well when <tt>addBorder</tt> is true. 01586 01587 <b> Declarations:</b> 01588 01589 pass arguments explicitly: 01590 \code 01591 namespace vigra { 01592 template <class SrcIterator, class SrcAccessor, 01593 class DestIterator, class DestAccessor, 01594 class GradValue, class DestValue> 01595 void cannyEdgeImageWithThinning( 01596 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01597 DestIterator dul, DestAccessor da, 01598 double scale, GradValue gradient_threshold, 01599 DestValue edge_marker, bool addBorder = true); 01600 } 01601 \endcode 01602 01603 use argument objects in conjunction with \ref ArgumentObjectFactories: 01604 \code 01605 namespace vigra { 01606 template <class SrcIterator, class SrcAccessor, 01607 class DestIterator, class DestAccessor, 01608 class GradValue, class DestValue> 01609 void cannyEdgeImageWithThinning( 01610 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01611 pair<DestIterator, DestAccessor> dest, 01612 double scale, GradValue gradient_threshold, 01613 DestValue edge_marker, bool addBorder = true); 01614 } 01615 \endcode 01616 01617 <b> Usage:</b> 01618 01619 <b>\#include</b> "<a href="edgedetection_8hxx-source.html">vigra/edgedetection.hxx</a>"<br> 01620 Namespace: vigra 01621 01622 \code 01623 vigra::BImage src(w,h), edges(w,h); 01624 01625 // empty edge image 01626 edges = 0; 01627 ... 01628 01629 // find edges at scale 0.8 with gradient larger than 4.0, mark with 1, annd add border 01630 vigra::cannyEdgeImageWithThinning(srcImageRange(src), destImage(edges), 01631 0.8, 4.0, 1, true); 01632 \endcode 01633 01634 <b> Required Interface:</b> 01635 01636 see also: \ref cannyEdgelList(). 01637 01638 \code 01639 DestImageIterator dest_upperleft; 01640 DestAccessor dest_accessor; 01641 DestValue edge_marker; 01642 01643 dest_accessor.set(edge_marker, dest_upperleft, vigra::Diff2D(1,1)); 01644 \endcode 01645 01646 <b> Preconditions:</b> 01647 01648 \code 01649 scale > 0 01650 gradient_threshold > 0 01651 \endcode 01652 */ 01653 template <class SrcIterator, class SrcAccessor, 01654 class DestIterator, class DestAccessor, 01655 class GradValue, class DestValue> 01656 void cannyEdgeImageWithThinning( 01657 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01658 DestIterator dul, DestAccessor da, 01659 double scale, GradValue gradient_threshold, 01660 DestValue edge_marker, bool addBorder) 01661 { 01662 int w = slr.x - sul.x; 01663 int h = slr.y - sul.y; 01664 01665 BImage edgeImage(w, h, BImage::value_type(0)); 01666 BImage::traverser eul = edgeImage.upperLeft(); 01667 BImage::Accessor ea = edgeImage.accessor(); 01668 if(addBorder) 01669 initImageBorder(destImageRange(edgeImage), 1, 1); 01670 cannyEdgeImage(sul, slr, sa, eul, ea, 01671 scale, gradient_threshold, 1); 01672 01673 static bool isSimplePoint[256] = { 01674 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 01675 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 01676 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 01677 1, 1, 1, 0, 0, 0, 1, 1, 1, 1, 0, 1, 0, 1, 0, 1, 0, 1, 01678 0, 0, 0, 0, 0, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1, 1, 01679 0, 1, 1, 0, 1, 0, 0, 1, 0, 1, 0, 1, 0, 1, 0, 0, 0, 0, 01680 0, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 0, 01681 1, 0, 0, 0, 0, 1, 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 0, 1, 01682 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 01683 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 01684 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 01685 0, 1, 0, 0, 0, 0, 0, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 0, 01686 1, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 0, 1, 0, 1, 0, 0, 01687 0, 0, 0, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 01688 1, 0, 1, 0 }; 01689 01690 eul += Diff2D(1,1); 01691 sul += Diff2D(1,1); 01692 int w2 = w-2; 01693 int h2 = h-2; 01694 01695 typedef detail::SimplePoint<GradValue> SP; 01696 std::priority_queue<SP, std::vector<SP> > pqueue; 01697 01698 Diff2D p(0,0); 01699 for(; p.y < h2; ++p.y) 01700 { 01701 for(p.x = 0; p.x < w2; ++p.x) 01702 { 01703 BImage::traverser e = eul + p; 01704 if(*e == 0) 01705 continue; 01706 int v = detail::neighborhoodConfiguration(e); 01707 if(isSimplePoint[v]) 01708 { 01709 pqueue.push(SP(p, norm(sa(sul+p)))); 01710 *e = 2; // remember that it is already in queue 01711 } 01712 } 01713 } 01714 01715 static const Diff2D dist[] = { Diff2D(-1,0), Diff2D(0,-1), 01716 Diff2D(1,0), Diff2D(0,1) }; 01717 01718 while(pqueue.size()) 01719 { 01720 p = pqueue.top().point; 01721 pqueue.pop(); 01722 01723 BImage::traverser e = eul + p; 01724 int v = detail::neighborhoodConfiguration(e); 01725 if(!isSimplePoint[v]) 01726 continue; // point may no longer be simple because its neighbors changed 01727 01728 *e = 0; // delete simple point 01729 01730 for(int i=0; i<4; ++i) 01731 { 01732 Diff2D pneu = p + dist[i]; 01733 if(pneu.x == -1 || pneu.y == -1 || pneu.x == w2 || pneu.y == h2) 01734 continue; // do not remove points at the border 01735 01736 BImage::traverser eneu = eul + pneu; 01737 if(*eneu == 1) // point is boundary and not yet in the queue 01738 { 01739 int v = detail::neighborhoodConfiguration(eneu); 01740 if(isSimplePoint[v]) 01741 { 01742 pqueue.push(SP(pneu, norm(sa(sul+pneu)))); 01743 *eneu = 2; // remember that it is already in queue 01744 } 01745 } 01746 } 01747 } 01748 01749 initImageIf(destIterRange(dul, dul+Diff2D(w,h), da), 01750 maskImage(edgeImage), edge_marker); 01751 } 01752 01753 template <class SrcIterator, class SrcAccessor, 01754 class DestIterator, class DestAccessor, 01755 class GradValue, class DestValue> 01756 inline void cannyEdgeImageWithThinning( 01757 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01758 pair<DestIterator, DestAccessor> dest, 01759 double scale, GradValue gradient_threshold, 01760 DestValue edge_marker, bool addBorder) 01761 { 01762 cannyEdgeImageWithThinning(src.first, src.second, src.third, 01763 dest.first, dest.second, 01764 scale, gradient_threshold, edge_marker, addBorder); 01765 } 01766 01767 template <class SrcIterator, class SrcAccessor, 01768 class DestIterator, class DestAccessor, 01769 class GradValue, class DestValue> 01770 inline void cannyEdgeImageWithThinning( 01771 SrcIterator sul, SrcIterator slr, SrcAccessor sa, 01772 DestIterator dul, DestAccessor da, 01773 double scale, GradValue gradient_threshold, DestValue edge_marker) 01774 { 01775 cannyEdgeImageWithThinning(sul, slr, sa, 01776 dul, da, 01777 scale, gradient_threshold, edge_marker, true); 01778 } 01779 01780 template <class SrcIterator, class SrcAccessor, 01781 class DestIterator, class DestAccessor, 01782 class GradValue, class DestValue> 01783 inline void cannyEdgeImageWithThinning( 01784 triple<SrcIterator, SrcIterator, SrcAccessor> src, 01785 pair<DestIterator, DestAccessor> dest, 01786 double scale, GradValue gradient_threshold, DestValue edge_marker) 01787 { 01788 cannyEdgeImageWithThinning(src.first, src.second, src.third, 01789 dest.first, dest.second, 01790 scale, gradient_threshold, edge_marker, true); 01791 } 01792 01793 //@} 01794 01795 /** \page CrackEdgeImage Crack Edge Image 01796 01797 Crack edges are marked <i>between</i> the pixels of an image. 01798 A Crack Edge Image is an image that represents these edges. In order 01799 to accomodate the cracks, the Crack Edge Image must be twice as large 01800 as the original image (precisely (2*w - 1) by (2*h - 1)). A Crack Edge Image 01801 can easily be derived from a binary image or from the signs of the 01802 response of a Laplacean filter. Consider the following sketch, where 01803 <TT>+</TT> encodes the foreground, <TT>-</TT> the background, and 01804 <TT>*</TT> the resulting crack edges. 01805 01806 \code 01807 sign of difference image insert cracks resulting CrackEdgeImage 01808 01809 + . - . - . * . . . 01810 + - - . . . . . . * * * . 01811 + + - => + . + . - => . . . * . 01812 + + + . . . . . . . . * * 01813 + . + . + . . . . . 01814 \endcode 01815 01816 Starting from the original binary image (left), we insert crack pixels 01817 to get to the double-sized image (center). Finally, we mark all 01818 crack pixels whose non-crack neighbors have different signs as 01819 crack edge points, while all other pixels (crack and non-crack) become 01820 region pixels. 01821 01822 <b>Requirements on a Crack Edge Image:</b> 01823 01824 <ul> 01825 <li>Crack Edge Images have odd width and height. 01826 <li>Crack pixels have at least one odd coordinate. 01827 <li>Only crack pixels may be marked as edge points. 01828 <li>Crack pixels with two odd coordinates must be marked as edge points 01829 whenever any of their neighboring crack pixels was marked. 01830 </ul> 01831 01832 The last two requirements ensure that both edges and regions are 4-connected. 01833 Thus, 4-connectivity and 8-connectivity yield identical connected 01834 components in a Crack Edge Image (so called <i>well-composedness</i>). 01835 This ensures that Crack Edge Images have nice topological properties 01836 (cf. L. J. Latecki: "Well-Composed Sets", Academic Press, 2000). 01837 */ 01838 01839 01840 } // namespace vigra 01841 01842 #endif // VIGRA_EDGEDETECTION_HXX
© Ullrich Köthe (koethe@informatik.uni-hamburg.de) |
html generated using doxygen and Python
|