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

seededregiongrowing.hxx
1/************************************************************************/
2/* */
3/* Copyright 1998-2010 by Ullrich Koethe, Hans Meine */
4/* */
5/* This file is part of the VIGRA computer vision library. */
6/* The VIGRA Website is */
7/* http://hci.iwr.uni-heidelberg.de/vigra/ */
8/* Please direct questions, bug reports, and contributions to */
9/* ullrich.koethe@iwr.uni-heidelberg.de or */
10/* vigra@informatik.uni-hamburg.de */
11/* */
12/* Permission is hereby granted, free of charge, to any person */
13/* obtaining a copy of this software and associated documentation */
14/* files (the "Software"), to deal in the Software without */
15/* restriction, including without limitation the rights to use, */
16/* copy, modify, merge, publish, distribute, sublicense, and/or */
17/* sell copies of the Software, and to permit persons to whom the */
18/* Software is furnished to do so, subject to the following */
19/* conditions: */
20/* */
21/* The above copyright notice and this permission notice shall be */
22/* included in all copies or substantial portions of the */
23/* Software. */
24/* */
25/* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND */
26/* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES */
27/* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND */
28/* NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT */
29/* HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, */
30/* WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING */
31/* FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR */
32/* OTHER DEALINGS IN THE SOFTWARE. */
33/* */
34/************************************************************************/
35
36#ifndef VIGRA_SEEDEDREGIONGROWING_HXX
37#define VIGRA_SEEDEDREGIONGROWING_HXX
38
39#include <vector>
40#include <stack>
41#include <queue>
42#include "utilities.hxx"
43#include "stdimage.hxx"
44#include "stdimagefunctions.hxx"
45#include "pixelneighborhood.hxx"
46#include "bucket_queue.hxx"
47#include "multi_shape.hxx"
48
49namespace vigra {
50
51namespace detail {
52
53template <class COST>
54class SeedRgPixel
55{
56public:
57 Point2D location_, nearest_;
58 COST cost_;
59 int count_;
60 int label_;
61 int dist_;
62
63 SeedRgPixel()
64 : location_(0,0), nearest_(0,0), cost_(0), count_(0), label_(0)
65 {}
66
67 SeedRgPixel(Point2D const & location, Point2D const & nearest,
68 COST const & cost, int const & count, int const & label)
69 : location_(location), nearest_(nearest),
70 cost_(cost), count_(count), label_(label)
71 {
72 int dx = location_.x - nearest_.x;
73 int dy = location_.y - nearest_.y;
74 dist_ = dx * dx + dy * dy;
75 }
76
77 void set(Point2D const & location, Point2D const & nearest,
78 COST const & cost, int const & count, int const & label)
79 {
80 location_ = location;
81 nearest_ = nearest;
82 cost_ = cost;
83 count_ = count;
84 label_ = label;
85
86 int dx = location_.x - nearest_.x;
87 int dy = location_.y - nearest_.y;
88 dist_ = dx * dx + dy * dy;
89 }
90
91 struct Compare
92 {
93 // must implement > since priority_queue looks for largest element
94 bool operator()(SeedRgPixel const & l,
95 SeedRgPixel const & r) const
96 {
97 if(r.cost_ == l.cost_)
98 {
99 if(r.dist_ == l.dist_) return r.count_ < l.count_;
100
101 return r.dist_ < l.dist_;
102 }
103
104 return r.cost_ < l.cost_;
105 }
106 bool operator()(SeedRgPixel const * l,
107 SeedRgPixel const * r) const
108 {
109 if(r->cost_ == l->cost_)
110 {
111 if(r->dist_ == l->dist_) return r->count_ < l->count_;
112
113 return r->dist_ < l->dist_;
114 }
115
116 return r->cost_ < l->cost_;
117 }
118 };
119
120 struct Allocator
121 {
122 ~Allocator()
123 {
124 while(!freelist_.empty())
125 {
126 delete freelist_.top();
127 freelist_.pop();
128 }
129 }
130
131 SeedRgPixel *
132 create(Point2D const & location, Point2D const & nearest,
133 COST const & cost, int const & count, int const & label)
134 {
135 if(!freelist_.empty())
136 {
137 SeedRgPixel * res = freelist_.top();
138 freelist_.pop();
139 res->set(location, nearest, cost, count, label);
140 return res;
141 }
142
143 return new SeedRgPixel(location, nearest, cost, count, label);
144 }
145
146 void dismiss(SeedRgPixel * p)
147 {
148 freelist_.push(p);
149 }
150
151 std::stack<SeedRgPixel<COST> *> freelist_;
152 };
153};
154
155struct UnlabelWatersheds
156{
157 int operator()(int label) const
158 {
159 return label < 0 ? 0 : label;
160 }
161};
162
163} // namespace detail
164
165/** \addtogroup Superpixels
166*/
167//@{
168
169/********************************************************/
170/* */
171/* seededRegionGrowing */
172/* */
173/********************************************************/
174
175/** Choose between different types of Region Growing */
177 CompleteGrow = 0,
178 KeepContours = 1,
179 StopAtThreshold = 2,
180 SRGWatershedLabel = -1
181};
182
183/** \brief Region Segmentation by means of Seeded Region Growing.
184
185 This algorithm implements seeded region growing as described in
186
187 R. Adams, L. Bischof: <em>"Seeded Region Growing"</em>, IEEE Trans. on Pattern
188 Analysis and Maschine Intelligence, vol 16, no 6, 1994, and
189
190 Ullrich K&ouml;the:
191 <em><a href="http://hci.iwr.uni-heidelberg.de/people/ukoethe/papers/index.php#cite_primary_segmentation">Primary Image Segmentation</a></em>,
192 in: G. Sagerer, S.
193 Posch, F. Kummert (eds.): Mustererkennung 1995, Proc. 17. DAGM-Symposium,
194 Springer 1995
195
196 The seed image is a partly segmented image which contains uniquely
197 labeled regions (the seeds) and unlabeled pixels (the candidates, label 0).
198 The highest seed label found in the seed image is returned by the algorithm.
199
200 Seed regions can be as large as you wish and as small as one pixel. If
201 there are no candidates, the algorithm will simply copy the seed image
202 into the output image. Otherwise it will aggregate the candidates into
203 the existing regions so that a cost function is minimized.
204 Candidates are taken from the neighborhood of the already assigned pixels,
205 where the type of neighborhood is determined by parameter <tt>neighborhood</tt>
206 which can take the values <tt>FourNeighborCode()</tt> (the default)
207 or <tt>EightNeighborCode()</tt>. The algorithm basically works as follows
208 (illustrated for 4-neighborhood, but 8-neighborhood works in the same way):
209
210 <ol>
211
212 <li> Find all candidate pixels that are 4-adjacent to a seed region.
213 Calculate the cost for aggregating each candidate into its adjacent region
214 and put the candidates into a priority queue.
215
216 <li> While( priority queue is not empty and termination criterion is not fulfilled)
217
218 <ol>
219
220 <li> Take the candidate with least cost from the queue. If it has not
221 already been merged, merge it with it's adjacent region.
222
223 <li> Put all candidates that are 4-adjacent to the pixel just processed
224 into the priority queue.
225
226 </ol>
227
228 </ol>
229
230 <tt>SRGType</tt> can take the following values:
231
232 <DL>
233 <DT><tt>CompleteGrow</tt> <DD> produce a complete tesselation of the volume (default).
234 <DT><tt>KeepContours</tt> <DD> keep a 1-voxel wide unlabeled contour between all regions.
235 <DT><tt>StopAtThreshold</tt> <DD> stop when the boundary indicator values exceed the
236 threshold given by parameter <tt>max_cost</tt>.
237 <DT><tt>KeepContours | StopAtThreshold</tt> <DD> keep 1-voxel wide contour and stop at given <tt>max_cost</tt>.
238 </DL>
239
240 The cost is determined jointly by the source image and the
241 region statistics functor. The source image contains feature values for each
242 pixel which will be used by the region statistics functor to calculate and
243 update statistics for each region and to calculate the cost for each
244 candidate. The <TT>RegionStatisticsArray</TT> must be compatible to the
245 \ref ArrayOfRegionStatistics functor and contains an <em> array</em> of
246 statistics objects for each region. The indices must correspond to the
247 labels of the seed regions. The statistics for the initial regions must have
248 been calculated prior to calling <TT>seededRegionGrowing()</TT> (for example by
249 means of \ref inspectTwoImagesIf()).
250
251 For each candidate
252 <TT>x</TT> that is adjacent to region <TT>i</TT>, the algorithm will call
253 <TT>stats[i].cost(as(x))</TT> to get the cost (where <TT>x</TT> is a <TT>SrcIterator</TT>
254 and <TT>as</TT> is
255 the SrcAccessor). When a candidate has been merged with a region, the
256 statistics are updated by calling <TT>stats[i].operator()(as(x))</TT>. Since
257 the <TT>RegionStatisticsArray</TT> is passed by reference, this will overwrite
258 the original statistics.
259
260 If a candidate could be merged into more than one regions with identical
261 cost, the algorithm will favour the nearest region. If <tt>StopAtThreshold</tt> is active,
262 and the cost of the current candidate at any point in the algorithm exceeds the optional
263 <tt>max_cost</tt> value (which defaults to <tt>NumericTraits<double>::max()</tt>),
264 region growing is aborted, and all voxels not yet assigned to a region remain unlabeled.
265
266 In some cases, the cost only depends on the feature value of the current
267 pixel. Then the update operation will simply be a no-op, and the <TT>cost()</TT>
268 function returns its argument. This behavior is implemented by the
269 \ref SeedRgDirectValueFunctor. With <tt>SRGType == KeepContours</tt>,
270 this is equivalent to the watershed algorithm.
271
272 <b> Declarations:</b>
273
274 pass 2D array views:
275 \code
276 namespace vigra {
277 template <class T1, class S1,
278 class TS, class AS,
279 class T2, class S2,
280 class RegionStatisticsArray, class Neighborhood>
281 TS
282 seededRegionGrowing(MultiArrayView<2, T1, S1> const & src,
283 MultiArrayView<2, TS, AS> const & seeds,
284 MultiArrayView<2, T2, S2> labels,
285 RegionStatisticsArray & stats,
286 SRGType srgType = CompleteGrow,
287 Neighborhood n = FourNeighborCode(),
288 double max_cost = NumericTraits<double>::max());
289 }
290 \endcode
291
292 \deprecatedAPI{seededRegionGrowing}
293 pass \ref ImageIterators and \ref DataAccessors :
294 \code
295 namespace vigra {
296 template <class SrcIterator, class SrcAccessor,
297 class SeedImageIterator, class SeedAccessor,
298 class DestIterator, class DestAccessor,
299 class RegionStatisticsArray, class Neighborhood>
300 typename SeedAccessor::value_type
301 seededRegionGrowing(SrcIterator srcul, SrcIterator srclr, SrcAccessor as,
302 SeedImageIterator seedsul, SeedAccessor aseeds,
303 DestIterator destul, DestAccessor ad,
304 RegionStatisticsArray & stats,
305 SRGType srgType = CompleteGrow,
306 Neighborhood neighborhood = FourNeighborCode(),
307 double max_cost = NumericTraits<double>::max());
308 }
309 \endcode
310 use argument objects in conjunction with \ref ArgumentObjectFactories :
311 \code
312 namespace vigra {
313 template <class SrcIterator, class SrcAccessor,
314 class SeedImageIterator, class SeedAccessor,
315 class DestIterator, class DestAccessor,
316 class RegionStatisticsArray, class Neighborhood>
317 typename SeedAccessor::value_type
318 seededRegionGrowing(triple<SrcIterator, SrcIterator, SrcAccessor> src,
319 pair<SeedImageIterator, SeedAccessor> seeds,
320 pair<DestIterator, DestAccessor> dest,
321 RegionStatisticsArray & stats,
322 SRGType srgType = CompleteGrow,
323 Neighborhood neighborhood = FourNeighborCode(),
324 double max_cost = NumericTraits<double>::max());
325 }
326 \endcode
327 \deprecatedEnd
328
329 <b> Usage:</b>
330
331 <b>\#include</b> <vigra/seededregiongrowing.hxx><br>
332 Namespace: vigra
333
334 Example: implementation of the voronoi tesselation
335
336 \code
337 MultiArray<2, int> points(w,h);
338 MultiArray<2, float> dist(x,y);
339
340 int max_region_label = 100;
341
342 // throw in some random points:
343 for(int i = 1; i <= max_region_label; ++i)
344 points(w * rand() / RAND_MAX , h * rand() / RAND_MAX) = i;
345
346 // calculate Euclidean distance transform
347 distanceTransform(points, dist, 2);
348
349 // init statistics functor
350 ArrayOfRegionStatistics<SeedRgDirectValueFunctor<float> > stats(max_region_label);
351
352 // find voronoi region of each point (the point image is overwritten with the
353 // voronoi region labels)
354 seededRegionGrowing(dist, points, points, stats);
355 \endcode
356
357 \deprecatedUsage{seededRegionGrowing}
358 \code
359 vigra::BImage points(w,h);
360 vigra::FImage dist(x,y);
361
362 // empty edge image
363 points = 0;
364 dist = 0;
365
366 int max_region_label = 100;
367
368 // throw in some random points:
369 for(int i = 1; i <= max_region_label; ++i)
370 points(w * rand() / RAND_MAX , h * rand() / RAND_MAX) = i;
371
372 // calculate Euclidean distance transform
373 vigra::distanceTransform(srcImageRange(points), destImage(dist), 2);
374
375 // init statistics functor
376 vigra::ArrayOfRegionStatistics<vigra::SeedRgDirectValueFunctor<float> >
377 stats(max_region_label);
378
379 // find voronoi region of each point
380 vigra:: seededRegionGrowing(srcImageRange(dist), srcImage(points),
381 destImage(points), stats);
382 \endcode
383 <b> Required Interface:</b>
384 \code
385 SrcIterator src_upperleft, src_lowerright;
386 SeedImageIterator seed_upperleft;
387 DestIterator dest_upperleft;
388
389 SrcAccessor src_accessor;
390 SeedAccessor seed_accessor;
391 DestAccessor dest_accessor;
392
393 RegionStatisticsArray stats;
394
395 // calculate costs
396 RegionStatisticsArray::value_type::cost_type cost =
397 stats[seed_accessor(seed_upperleft)].cost(src_accessor(src_upperleft));
398
399 // compare costs
400 cost < cost;
401
402 // update statistics
403 stats[seed_accessor(seed_upperleft)](src_accessor(src_upperleft));
404
405 // set result
406 dest_accessor.set(seed_accessor(seed_upperleft), dest_upperleft);
407 \endcode
408 \deprecatedEnd
409
410 Further requirements are determined by the <TT>RegionStatisticsArray</TT>.
411*/
412doxygen_overloaded_function(template <...> void seededRegionGrowing)
413
414template <class SrcIterator, class SrcAccessor,
415 class SeedImageIterator, class SeedAccessor,
416 class DestIterator, class DestAccessor,
418typename SeedAccessor::value_type
423 RegionStatisticsArray & stats,
424 SRGType srgType,
426 double max_cost)
427{
428 int w = srclr.x - srcul.x;
429 int h = srclr.y - srcul.y;
430 int count = 0;
431
432 SrcIterator isy = srcul, isx = srcul; // iterators for the src image
433
434 typedef typename SeedAccessor::value_type LabelType;
435 typedef typename RegionStatisticsArray::value_type RegionStatistics;
436 typedef typename RegionStatistics::cost_type CostType;
437 typedef detail::SeedRgPixel<CostType> Pixel;
438
439 typename Pixel::Allocator allocator;
440
441 typedef std::priority_queue<Pixel *, std::vector<Pixel *>,
442 typename Pixel::Compare> SeedRgPixelHeap;
443
444 // copy seed image in an image with border
445 IImage regions(w+2, h+2);
446 IImage::Iterator ir = regions.upperLeft() + Diff2D(1,1);
447 IImage::Iterator iry, irx;
448
449 initImageBorder(destImageRange(regions), 1, SRGWatershedLabel);
450 copyImage(seedsul, seedsul+Diff2D(w,h), aseeds, ir, regions.accessor());
451
452 // allocate and init memory for the results
453
455 int cneighbor, maxRegionLabel = 0;
456
457 typedef typename Neighborhood::Direction Direction;
458 int directionCount = Neighborhood::DirectionCount;
459
460 Point2D pos(0,0);
461 for(isy=srcul, iry=ir, pos.y=0; pos.y<h;
462 ++pos.y, ++isy.y, ++iry.y)
463 {
464 for(isx=isy, irx=iry, pos.x=0; pos.x<w;
465 ++pos.x, ++isx.x, ++irx.x)
466 {
467 if(*irx == 0)
468 {
469 // find candidate pixels for growing and fill heap
470 for(int i=0; i<directionCount; i++)
471 {
472 // cneighbor = irx[dist[i]];
473 cneighbor = irx[Neighborhood::diff((Direction)i)];
474 if(cneighbor > 0)
475 {
476 CostType cost = stats[cneighbor].cost(as(isx));
477
478 Pixel * pixel =
479 allocator.create(pos, pos+Neighborhood::diff((Direction)i), cost, count++, cneighbor);
480 pheap.push(pixel);
481 }
482 }
483 }
484 else
485 {
486 vigra_precondition((LabelType)*irx <= (LabelType)stats.maxRegionLabel(),
487 "seededRegionGrowing(): Largest label exceeds size of RegionStatisticsArray.");
488 if(maxRegionLabel < *irx)
489 maxRegionLabel = *irx;
490 }
491 }
492 }
493
494 // perform region growing
495 while(pheap.size() != 0)
496 {
497 Pixel * pixel = pheap.top();
498 pheap.pop();
499
500 Point2D pos = pixel->location_;
501 Point2D nearest = pixel->nearest_;
502 int lab = pixel->label_;
503 CostType cost = pixel->cost_;
504
505 allocator.dismiss(pixel);
506
507 if((srgType & StopAtThreshold) != 0 && cost > max_cost)
508 break;
509
510 irx = ir + pos;
511 isx = srcul + pos;
512
513 if(*irx) // already labelled region / watershed?
514 continue;
515
516 if((srgType & KeepContours) != 0)
517 {
518 for(int i=0; i<directionCount; i++)
519 {
520 cneighbor = irx[Neighborhood::diff((Direction)i)];
521 if((cneighbor>0) && (cneighbor != lab))
522 {
523 lab = SRGWatershedLabel;
524 break;
525 }
526 }
527 }
528
529 *irx = lab;
530
531 if((srgType & KeepContours) == 0 || lab > 0)
532 {
533 // update statistics
534 stats[*irx](as(isx));
535
536 // search neighborhood
537 // second pass: find new candidate pixels
538 for(int i=0; i<directionCount; i++)
539 {
540 if(irx[Neighborhood::diff((Direction)i)] == 0)
541 {
542 CostType cost = stats[lab].cost(as(isx, Neighborhood::diff((Direction)i)));
543
544 Pixel * new_pixel =
545 allocator.create(pos+Neighborhood::diff((Direction)i), nearest, cost, count++, lab);
546 pheap.push(new_pixel);
547 }
548 }
549 }
550 }
551
552 // free temporary memory
553 while(pheap.size() != 0)
554 {
555 allocator.dismiss(pheap.top());
556 pheap.pop();
557 }
558
559 // write result
560 transformImage(ir, ir+Point2D(w,h), regions.accessor(), destul, ad,
561 detail::UnlabelWatersheds());
562
563 return (LabelType)maxRegionLabel;
564}
565
566template <class SrcIterator, class SrcAccessor,
567 class SeedImageIterator, class SeedAccessor,
568 class DestIterator, class DestAccessor,
569 class RegionStatisticsArray, class Neighborhood>
570inline typename SeedAccessor::value_type
571seededRegionGrowing(SrcIterator srcul,
572 SrcIterator srclr, SrcAccessor as,
573 SeedImageIterator seedsul, SeedAccessor aseeds,
574 DestIterator destul, DestAccessor ad,
575 RegionStatisticsArray & stats,
576 SRGType srgType,
577 Neighborhood n)
578{
579 return seededRegionGrowing(srcul, srclr, as,
580 seedsul, aseeds,
581 destul, ad,
582 stats, srgType, n, NumericTraits<double>::max());
583}
584
585
586
587template <class SrcIterator, class SrcAccessor,
588 class SeedImageIterator, class SeedAccessor,
589 class DestIterator, class DestAccessor,
590 class RegionStatisticsArray>
591inline typename SeedAccessor::value_type
592seededRegionGrowing(SrcIterator srcul,
593 SrcIterator srclr, SrcAccessor as,
594 SeedImageIterator seedsul, SeedAccessor aseeds,
595 DestIterator destul, DestAccessor ad,
596 RegionStatisticsArray & stats,
597 SRGType srgType)
598{
599 return seededRegionGrowing(srcul, srclr, as,
600 seedsul, aseeds,
601 destul, ad,
602 stats, srgType, FourNeighborCode());
603}
604
605template <class SrcIterator, class SrcAccessor,
606 class SeedImageIterator, class SeedAccessor,
607 class DestIterator, class DestAccessor,
608 class RegionStatisticsArray>
609inline typename SeedAccessor::value_type
610seededRegionGrowing(SrcIterator srcul,
611 SrcIterator srclr, SrcAccessor as,
612 SeedImageIterator seedsul, SeedAccessor aseeds,
613 DestIterator destul, DestAccessor ad,
614 RegionStatisticsArray & stats)
615{
616 return seededRegionGrowing(srcul, srclr, as,
617 seedsul, aseeds,
618 destul, ad,
619 stats, CompleteGrow);
620}
621
622template <class SrcIterator, class SrcAccessor,
623 class SeedImageIterator, class SeedAccessor,
624 class DestIterator, class DestAccessor,
625 class RegionStatisticsArray, class Neighborhood>
626inline typename SeedAccessor::value_type
627seededRegionGrowing(triple<SrcIterator, SrcIterator, SrcAccessor> img1,
628 pair<SeedImageIterator, SeedAccessor> img3,
629 pair<DestIterator, DestAccessor> img4,
630 RegionStatisticsArray & stats,
631 SRGType srgType,
632 Neighborhood n,
633 double max_cost = NumericTraits<double>::max())
634{
635 return seededRegionGrowing(img1.first, img1.second, img1.third,
636 img3.first, img3.second,
637 img4.first, img4.second,
638 stats, srgType, n, max_cost);
639}
640
641template <class SrcIterator, class SrcAccessor,
642 class SeedImageIterator, class SeedAccessor,
643 class DestIterator, class DestAccessor,
644 class RegionStatisticsArray>
645inline typename SeedAccessor::value_type
646seededRegionGrowing(triple<SrcIterator, SrcIterator, SrcAccessor> img1,
647 pair<SeedImageIterator, SeedAccessor> img3,
648 pair<DestIterator, DestAccessor> img4,
649 RegionStatisticsArray & stats,
650 SRGType srgType)
651{
652 return seededRegionGrowing(img1.first, img1.second, img1.third,
653 img3.first, img3.second,
654 img4.first, img4.second,
655 stats, srgType, FourNeighborCode());
656}
657
658template <class SrcIterator, class SrcAccessor,
659 class SeedImageIterator, class SeedAccessor,
660 class DestIterator, class DestAccessor,
661 class RegionStatisticsArray>
662inline typename SeedAccessor::value_type
663seededRegionGrowing(triple<SrcIterator, SrcIterator, SrcAccessor> img1,
664 pair<SeedImageIterator, SeedAccessor> img3,
665 pair<DestIterator, DestAccessor> img4,
666 RegionStatisticsArray & stats)
667{
668 return seededRegionGrowing(img1.first, img1.second, img1.third,
669 img3.first, img3.second,
670 img4.first, img4.second,
671 stats, CompleteGrow);
672}
673
674template <class T1, class S1,
675 class TS, class AS,
676 class T2, class S2,
677 class RegionStatisticsArray, class Neighborhood>
678inline TS
679seededRegionGrowing(MultiArrayView<2, T1, S1> const & img1,
680 MultiArrayView<2, TS, AS> const & img3,
681 MultiArrayView<2, T2, S2> img4,
682 RegionStatisticsArray & stats,
683 SRGType srgType,
684 Neighborhood n,
685 double max_cost = NumericTraits<double>::max())
686{
687 vigra_precondition(img1.shape() == img3.shape(),
688 "seededRegionGrowing(): shape mismatch between input and output.");
689 return seededRegionGrowing(srcImageRange(img1),
690 srcImage(img3),
691 destImage(img4),
692 stats, srgType, n, max_cost);
693}
694
695template <class T1, class S1,
696 class TS, class AS,
697 class T2, class S2,
698 class RegionStatisticsArray>
699inline TS
700seededRegionGrowing(MultiArrayView<2, T1, S1> const & img1,
701 MultiArrayView<2, TS, AS> const & img3,
702 MultiArrayView<2, T2, S2> img4,
703 RegionStatisticsArray & stats,
704 SRGType srgType)
705{
706 vigra_precondition(img1.shape() == img3.shape(),
707 "seededRegionGrowing(): shape mismatch between input and output.");
708 return seededRegionGrowing(srcImageRange(img1),
709 srcImage(img3),
710 destImage(img4),
711 stats, srgType, FourNeighborCode());
712}
713
714template <class T1, class S1,
715 class TS, class AS,
716 class T2, class S2,
717 class RegionStatisticsArray>
718inline TS
719seededRegionGrowing(MultiArrayView<2, T1, S1> const & img1,
720 MultiArrayView<2, TS, AS> const & img3,
721 MultiArrayView<2, T2, S2> img4,
722 RegionStatisticsArray & stats)
723{
724 vigra_precondition(img1.shape() == img3.shape(),
725 "seededRegionGrowing(): shape mismatch between input and output.");
726 return seededRegionGrowing(srcImageRange(img1),
727 srcImage(img3),
728 destImage(img4),
729 stats, CompleteGrow);
730}
731
732/********************************************************/
733/* */
734/* fastSeededRegionGrowing */
735/* */
736/********************************************************/
737
738template <class SrcIterator, class SrcAccessor,
739 class DestIterator, class DestAccessor,
740 class RegionStatisticsArray, class Neighborhood>
741typename DestAccessor::value_type
742fastSeededRegionGrowing(SrcIterator srcul, SrcIterator srclr, SrcAccessor as,
743 DestIterator destul, DestAccessor ad,
744 RegionStatisticsArray & stats,
745 SRGType srgType,
746 Neighborhood,
747 double max_cost,
748 std::ptrdiff_t bucket_count = 256)
749{
750 typedef typename DestAccessor::value_type LabelType;
751
752 vigra_precondition((srgType & KeepContours) == 0,
753 "fastSeededRegionGrowing(): the turbo algorithm doesn't support 'KeepContours', sorry.");
754
755 int w = srclr.x - srcul.x;
756 int h = srclr.y - srcul.y;
757
758 SrcIterator isy = srcul, isx = srcul; // iterators for the src image
759 DestIterator idy = destul, idx = destul; // iterators for the dest image
760
761 BucketQueue<Point2D, true> pqueue(bucket_count);
762 LabelType maxRegionLabel = 0;
763
764 Point2D pos(0,0);
765 for(isy=srcul, idy = destul, pos.y=0; pos.y<h; ++pos.y, ++isy.y, ++idy.y)
766 {
767 for(isx=isy, idx=idy, pos.x=0; pos.x<w; ++pos.x, ++isx.x, ++idx.x)
768 {
769 LabelType label = ad(idx);
770 if(label != 0)
771 {
772 vigra_precondition(label <= stats.maxRegionLabel(),
773 "fastSeededRegionGrowing(): Largest label exceeds size of RegionStatisticsArray.");
774
775 if(maxRegionLabel < label)
776 maxRegionLabel = label;
777
778 AtImageBorder atBorder = isAtImageBorder(pos.x, pos.y, w, h);
779 if(atBorder == NotAtBorder)
780 {
781 NeighborhoodCirculator<DestIterator, Neighborhood> c(idx), cend(c);
782 do
783 {
784 if(ad(c) == 0)
785 {
786 std::ptrdiff_t priority = (std::ptrdiff_t)stats[label].cost(as(isx));
787 pqueue.push(pos, priority);
788 break;
789 }
790 }
791 while(++c != cend);
792 }
793 else
794 {
795 RestrictedNeighborhoodCirculator<DestIterator, Neighborhood>
796 c(idx, atBorder), cend(c);
797 do
798 {
799 if(ad(c) == 0)
800 {
801 std::ptrdiff_t priority = (std::ptrdiff_t)stats[label].cost(as(isx));
802 pqueue.push(pos, priority);
803 break;
804 }
805 }
806 while(++c != cend);
807 }
808 }
809 }
810 }
811
812 // perform region growing
813 while(!pqueue.empty())
814 {
815 Point2D pos = pqueue.top();
816 std::ptrdiff_t cost = pqueue.topPriority();
817 pqueue.pop();
818
819 if((srgType & StopAtThreshold) != 0 && cost > max_cost)
820 break;
821
822 idx = destul + pos;
823 isx = srcul + pos;
824
825 std::ptrdiff_t label = ad(idx);
826
827 AtImageBorder atBorder = isAtImageBorder(pos.x, pos.y, w, h);
828 if(atBorder == NotAtBorder)
829 {
830 NeighborhoodCirculator<DestIterator, Neighborhood> c(idx), cend(c);
831
832 do
833 {
834 std::ptrdiff_t nlabel = ad(c);
835 if(nlabel == 0)
836 {
837 ad.set(label, idx, c.diff());
838 std::ptrdiff_t priority =
839 std::max((std::ptrdiff_t)stats[label].cost(as(isx, c.diff())), cost);
840 pqueue.push(pos+c.diff(), priority);
841 }
842 }
843 while(++c != cend);
844 }
845 else
846 {
847 RestrictedNeighborhoodCirculator<DestIterator, Neighborhood>
848 c(idx, atBorder), cend(c);
849 do
850 {
851 std::ptrdiff_t nlabel = ad(c);
852 if(nlabel == 0)
853 {
854 ad.set(label, idx, c.diff());
855 std::ptrdiff_t priority =
856 std::max((std::ptrdiff_t)stats[label].cost(as(isx, c.diff())), cost);
857 pqueue.push(pos+c.diff(), priority);
858 }
859 }
860 while(++c != cend);
861 }
862 }
863
864 return maxRegionLabel;
865}
866
867template <class SrcIterator, class SrcAccessor,
868 class DestIterator, class DestAccessor,
869 class RegionStatisticsArray, class Neighborhood>
870inline typename DestAccessor::value_type
871fastSeededRegionGrowing(SrcIterator srcul, SrcIterator srclr, SrcAccessor as,
872 DestIterator destul, DestAccessor ad,
873 RegionStatisticsArray & stats,
874 SRGType srgType,
875 Neighborhood n)
876{
877 return fastSeededRegionGrowing(srcul, srclr, as,
878 destul, ad,
879 stats, srgType, n, NumericTraits<double>::max(), 256);
880}
881
882template <class SrcIterator, class SrcAccessor,
883 class DestIterator, class DestAccessor,
884 class RegionStatisticsArray>
885inline typename DestAccessor::value_type
886fastSeededRegionGrowing(SrcIterator srcul, SrcIterator srclr, SrcAccessor as,
887 DestIterator destul, DestAccessor ad,
888 RegionStatisticsArray & stats,
889 SRGType srgType)
890{
891 return fastSeededRegionGrowing(srcul, srclr, as,
892 destul, ad,
893 stats, srgType, FourNeighborCode());
894}
895
896template <class SrcIterator, class SrcAccessor,
897 class DestIterator, class DestAccessor,
898 class RegionStatisticsArray>
899inline typename DestAccessor::value_type
900fastSeededRegionGrowing(SrcIterator srcul, SrcIterator srclr, SrcAccessor as,
901 DestIterator destul, DestAccessor ad,
902 RegionStatisticsArray & stats)
903{
904 return fastSeededRegionGrowing(srcul, srclr, as,
905 destul, ad,
906 stats, CompleteGrow);
907}
908
909template <class SrcIterator, class SrcAccessor,
910 class DestIterator, class DestAccessor,
911 class RegionStatisticsArray, class Neighborhood>
912inline typename DestAccessor::value_type
913fastSeededRegionGrowing(triple<SrcIterator, SrcIterator, SrcAccessor> src,
914 pair<DestIterator, DestAccessor> dest,
915 RegionStatisticsArray & stats,
916 SRGType srgType,
917 Neighborhood n,
918 double max_cost,
919 std::ptrdiff_t bucket_count = 256)
920{
921 return fastSeededRegionGrowing(src.first, src.second, src.third,
922 dest.first, dest.second,
923 stats, srgType, n, max_cost, bucket_count);
924}
925
926template <class T1, class S1,
927 class T2, class S2,
928 class RegionStatisticsArray, class Neighborhood>
929inline T2
930fastSeededRegionGrowing(MultiArrayView<2, T1, S1> const & src,
931 MultiArrayView<2, T2, S2> dest,
932 RegionStatisticsArray & stats,
933 SRGType srgType,
934 Neighborhood n,
935 double max_cost,
936 std::ptrdiff_t bucket_count = 256)
937{
938 vigra_precondition(src.shape() == dest.shape(),
939 "fastSeededRegionGrowing(): shape mismatch between input and output.");
940 return fastSeededRegionGrowing(srcImageRange(src),
941 destImage(dest),
942 stats, srgType, n, max_cost, bucket_count);
943}
944
945/********************************************************/
946/* */
947/* SeedRgDirectValueFunctor */
948/* */
949/********************************************************/
950
951/** \brief Statistics functor to be used for seeded region growing.
952
953 This functor can be used if the cost of a candidate during
954 \ref seededRegionGrowing() is equal to the feature value of that
955 candidate and does not depend on properties of the region it is going to
956 be merged with.
957
958 <b>\#include</b> <vigra/seededregiongrowing.hxx><br>
959 Namespace: vigra
960*/
961template <class Value>
963{
964 public:
965 /** the functor's argument type
966 */
967 typedef Value argument_type;
968
969 /** the functor's result type (unused, only necessary for
970 use of SeedRgDirectValueFunctor in \ref vigra::ArrayOfRegionStatistics
971 */
972 typedef Value result_type;
973
974 /** \deprecated use argument_type
975 */
976 typedef Value value_type;
977
978 /** the return type of the cost() function
979 */
980 typedef Value cost_type;
981
982 /** Do nothing (since we need not update region statistics).
983 */
984 void operator()(argument_type const &) const {}
985
986 /** Return argument (since cost is identical to feature value)
987 */
988 cost_type const & cost(argument_type const & v) const
989 {
990 return v;
991 }
992};
993
994//@}
995
996} // namespace vigra
997
998#endif // VIGRA_SEEDEDREGIONGROWING_HXX
999
Two dimensional difference vector.
Definition diff2d.hxx:186
int y
Definition diff2d.hxx:392
int x
Definition diff2d.hxx:385
Two dimensional point or position.
Definition diff2d.hxx:586
Class for a single RGB value.
Definition rgbvalue.hxx:128
Statistics functor to be used for seeded region growing.
Definition seededregiongrowing.hxx:963
Value value_type
Definition seededregiongrowing.hxx:976
Value result_type
Definition seededregiongrowing.hxx:972
cost_type const & cost(argument_type const &v) const
Definition seededregiongrowing.hxx:988
Value argument_type
Definition seededregiongrowing.hxx:967
Value cost_type
Definition seededregiongrowing.hxx:980
void operator()(argument_type const &) const
Definition seededregiongrowing.hxx:984
size_type size() const
Definition tinyvector.hxx:913
void transformImage(...)
Apply unary point transformation to each pixel.
SRGType
Definition seededregiongrowing.hxx:176
void seededRegionGrowing(...)
Region Segmentation by means of Seeded Region Growing.
FourNeighborhood::NeighborCode FourNeighborCode
Definition pixelneighborhood.hxx:379
AtImageBorder
Encode whether a point is near the image border.
Definition pixelneighborhood.hxx:69
@ NotAtBorder
&#160;
Definition pixelneighborhood.hxx:70
void initImageBorder(...)
Write value to the specified border pixels in the image.
AtImageBorder isAtImageBorder(int x, int y, int width, int height)
Find out whether a point is at the image border.
Definition pixelneighborhood.hxx:111
void copyImage(...)
Copy source image into destination image.

© Ullrich Köthe (ullrich.koethe@iwr.uni-heidelberg.de)
Heidelberg Collaboratory for Image Processing, University of Heidelberg, Germany

html generated using doxygen and Python
vigra 1.12.2