1 /************************************************************************/
2 /* */
3 /* Copyright 2012-2013 by Ullrich Koethe */
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_MULTI_LABELING_HXX
37 #define VIGRA_MULTI_LABELING_HXX
38
39 #include "multi_array.hxx"
40 #include "multi_gridgraph.hxx"
41 #include "union_find.hxx"
42 #include "any.hxx"
43
44 namespace vigra{
45
46 namespace labeling_equality{
47
48 struct Yes
49 {
50 char d[100];
51 };
52 struct No
53 {
54 char d[1];
55 };
56
57 template <size_t>
58 struct SizeToType;
59 template <>
60 struct SizeToType<sizeof(Yes)>
61 {
62 typedef VigraTrueType type;
63 };
64 template <>
65 struct SizeToType<sizeof(No)>
66 {
67 typedef VigraFalseType type;
68 };
69
70 template <class Equal>
71 class TakesThreeArguments
72 {
73 public:
74 template <class T>
75 static Yes check(typename T::WithDiffTag*);
76 template <class T>
77 static No check(...);
78
79 typedef typename SizeToType<sizeof(check<Equal>(0))>::type type;
80 static const unsigned int value = type::asBool;
81 };
82
83 template <class Equal, class Data, class Shape>
callEqualImpl(Equal & equal,const Data & u_data,const Data & v_data,const Shape & diff,VigraTrueType)84 bool callEqualImpl(Equal& equal, const Data& u_data, const Data& v_data, const Shape& diff, VigraTrueType)
85 {
86 return equal(u_data, v_data, diff);
87 }
88 template <class Equal, class Data, class Shape>
callEqualImpl(Equal & equal,const Data & u_data,const Data & v_data,const Shape &,VigraFalseType)89 bool callEqualImpl(Equal& equal, const Data& u_data, const Data& v_data, const Shape&, VigraFalseType)
90 {
91 return equal(u_data, v_data);
92 }
93
94 template< class Equal, class Data, class Shape>
callEqual(Equal & equal,const Data & u_data,const Data & v_data,const Shape & diff)95 bool callEqual(Equal& equal, const Data& u_data, const Data& v_data, const Shape& diff)
96 {
97 return callEqualImpl(equal, u_data, v_data, diff, typename TakesThreeArguments<Equal>::type());
98 }
99
100 }
101
102
103 /** \addtogroup Labeling
104 */
105 //@{
106
107 namespace lemon_graph {
108
109 template <class Graph, class T1Map, class T2Map, class Equal>
110 typename T2Map::value_type
labelGraph(Graph const & g,T1Map const & data,T2Map & labels,Equal const & equal)111 labelGraph(Graph const & g,
112 T1Map const & data,
113 T2Map & labels,
114 Equal const & equal)
115 {
116 typedef typename Graph::NodeIt graph_scanner;
117 typedef typename Graph::OutBackArcIt neighbor_iterator;
118 typedef typename T2Map::value_type LabelType;
119
120 vigra::UnionFindArray<LabelType> regions;
121
122 // pass 1: find connected components
123 for (graph_scanner node(g); node != INVALID; ++node)
124 {
125 typename T1Map::value_type center = data[*node];
126
127 // define tentative label for current node
128 LabelType currentIndex = regions.nextFreeIndex();
129
130 for (neighbor_iterator arc(g, node); arc != INVALID; ++arc)
131 {
132 // merge regions if colors are equal
133 if(equal(center, data[g.target(*arc)]))
134 {
135 currentIndex = regions.makeUnion(labels[g.target(*arc)], currentIndex);
136 }
137 }
138 // set label of current node
139 labels[*node] = regions.finalizeIndex(currentIndex);
140 }
141
142 LabelType count = regions.makeContiguous();
143
144 // pass 2: make component labels contiguous
145 for (graph_scanner node(g); node != INVALID; ++node)
146 {
147 labels[*node] = regions.findLabel(labels[*node]);
148 }
149 return count;
150 }
151
152 template <unsigned int N, class DirectedTag, class T1Map, class T2Map, class Equal>
153 typename T2Map::value_type
labelGraph(GridGraph<N,DirectedTag> const & g,T1Map const & data,T2Map & labels,Equal const & equal)154 labelGraph(GridGraph<N, DirectedTag> const & g,
155 T1Map const & data,
156 T2Map & labels,
157 Equal const & equal)
158 {
159 typedef GridGraph<N, DirectedTag> Graph;
160 typedef typename Graph::NodeIt graph_scanner;
161 typedef typename Graph::OutBackArcIt neighbor_iterator;
162 typedef typename T2Map::value_type LabelType;
163 typedef typename Graph::shape_type Shape;
164
165 vigra::UnionFindArray<LabelType> regions;
166
167 // pass 1: find connected components
168 for (graph_scanner node(g); node != INVALID; ++node)
169 {
170 typename T1Map::value_type center = data[*node];
171
172 // define tentative label for current node
173 LabelType currentIndex = regions.nextFreeIndex();
174
175 for (neighbor_iterator arc(g, node); arc != INVALID; ++arc)
176 {
177 Shape diff = g.neighborOffset(arc.neighborIndex());
178 // merge regions if colors are equal
179 if(labeling_equality::callEqual(equal, center, data[g.target(*arc)], diff))
180 {
181 currentIndex = regions.makeUnion(labels[g.target(*arc)], currentIndex);
182 }
183 }
184 // set label of current node
185 labels[*node] = regions.finalizeIndex(currentIndex);
186 }
187
188 LabelType count = regions.makeContiguous();
189
190 // pass 2: make component labels contiguous
191 for (graph_scanner node(g); node != INVALID; ++node)
192 {
193 labels[*node] = regions.findLabel(labels[*node]);
194 }
195 return count;
196 }
197
198
199 template <class Graph, class T1Map, class T2Map, class Equal>
200 typename T2Map::value_type
labelGraphWithBackground(Graph const & g,T1Map const & data,T2Map & labels,typename T1Map::value_type backgroundValue,Equal const & equal)201 labelGraphWithBackground(Graph const & g,
202 T1Map const & data,
203 T2Map & labels,
204 typename T1Map::value_type backgroundValue,
205 Equal const & equal)
206 {
207 typedef typename Graph::NodeIt graph_scanner;
208 typedef typename Graph::OutBackArcIt neighbor_iterator;
209 typedef typename T2Map::value_type LabelType;
210
211 vigra::UnionFindArray<LabelType> regions;
212
213 // pass 1: find connected components
214 for (graph_scanner node(g); node != INVALID; ++node)
215 {
216 typename T1Map::value_type center = data[*node];
217
218 // background always gets label zero
219 if(equal(center, backgroundValue))
220 {
221 labels[*node] = 0;
222 continue;
223 }
224
225 // define tentative label for current node
226 LabelType currentIndex = regions.nextFreeIndex();
227
228 for (neighbor_iterator arc(g, node); arc != INVALID; ++arc)
229 {
230 // merge regions if colors are equal
231 if(equal(center, data[g.target(*arc)]))
232 {
233 currentIndex = regions.makeUnion(labels[g.target(*arc)], currentIndex);
234 }
235 }
236 // set label of current node
237 labels[*node] = regions.finalizeIndex(currentIndex);
238 }
239
240 LabelType count = regions.makeContiguous();
241
242 // pass 2: make component labels contiguous
243 for (graph_scanner node(g); node != INVALID; ++node)
244 {
245 labels[*node] = regions.findLabel(labels[*node]);
246 }
247 return count;
248 }
249
250 template <unsigned int N, class DirectedTag, class T1Map, class T2Map, class Equal>
251 typename T2Map::value_type
labelGraphWithBackground(GridGraph<N,DirectedTag> const & g,T1Map const & data,T2Map & labels,typename T1Map::value_type backgroundValue,Equal const & equal)252 labelGraphWithBackground(GridGraph<N, DirectedTag> const & g,
253 T1Map const & data,
254 T2Map & labels,
255 typename T1Map::value_type backgroundValue,
256 Equal const & equal)
257 {
258 typedef GridGraph<N, DirectedTag> Graph;
259 typedef typename Graph::NodeIt graph_scanner;
260 typedef typename Graph::OutBackArcIt neighbor_iterator;
261 typedef typename T2Map::value_type LabelType;
262 typedef typename Graph::shape_type Shape;
263
264 vigra::UnionFindArray<LabelType> regions;
265
266 // pass 1: find connected components
267 for (graph_scanner node(g); node != INVALID; ++node)
268 {
269 typename T1Map::value_type center = data[*node];
270
271 // background always gets label zero
272 if(labeling_equality::callEqual(equal, center, backgroundValue, Shape()))
273 {
274 labels[*node] = 0;
275 continue;
276 }
277
278 // define tentative label for current node
279 LabelType currentIndex = regions.nextFreeIndex();
280
281 for (neighbor_iterator arc(g, node); arc != INVALID; ++arc)
282 {
283 // merge regions if colors are equal
284 Shape diff = g.neighborOffset(arc.neighborIndex());
285 if(labeling_equality::callEqual(equal, center, data[g.target(*arc)], diff))
286 {
287 currentIndex = regions.makeUnion(labels[g.target(*arc)], currentIndex);
288 }
289 }
290 // set label of current node
291 labels[*node] = regions.finalizeIndex(currentIndex);
292 }
293
294 LabelType count = regions.makeContiguous();
295
296 // pass 2: make component labels contiguous
297 for (graph_scanner node(g); node != INVALID; ++node)
298 {
299 labels[*node] = regions.findLabel(labels[*node]);
300 }
301 return count;
302 }
303
304
305 } // namespace lemon_graph
306
307 /** \brief Option object for labelMultiArray().
308 */
309 class LabelOptions
310 {
311 Any background_value_;
312 NeighborhoodType neighborhood_;
313
314 public:
315
316 /** \brief Set default options.
317 */
LabelOptions()318 LabelOptions()
319 : neighborhood_(DirectNeighborhood)
320 {}
321
322 /** \brief Choose direct or indirect neighborhood.
323
324 Default: <tt>DirectNeighborhood</tt>
325 */
neighborhood(NeighborhoodType n)326 LabelOptions & neighborhood(NeighborhoodType n)
327 {
328 neighborhood_ = n;
329 return *this;
330 }
331
332 /** \brief Query the neighborhood type (direct or indirect).
333 */
getNeighborhood() const334 NeighborhoodType getNeighborhood() const
335 {
336 return neighborhood_;
337 }
338
339 /** \brief Set background value.
340
341 If specified, labelMultiArray() will internally call
342 labelMultiArrayWithBackground() with the given value
343 considered as background and thus ignored. If no
344 background value is set, the array gets labeled completely.
345 Note that the type <tt>T</tt> must correspond to the element
346 type of the data array to be labeled.
347
348 Default: don't ignore any value.
349 */
350 template <class T>
ignoreBackgroundValue(T const & t)351 LabelOptions & ignoreBackgroundValue(T const & t)
352 {
353 background_value_ = t;
354 return *this;
355 }
356
357 /** \brief Check if some background value is to be ignored.
358 */
hasBackgroundValue() const359 bool hasBackgroundValue() const
360 {
361 return bool(background_value_);
362 }
363
364 /** \brief Get the background value to be ignored.
365
366 Throws an exception if the stored background value type
367 is incompatible to the data array's value type.
368 */
369 template <class T>
getBackgroundValue() const370 T getBackgroundValue() const
371 {
372 if(background_value_.empty())
373 return T();
374 vigra_precondition(background_value_.template is_readable<T>(),
375 "LabelOptions::getBackgroundValue<T>(): stored background value is not convertible to T.");
376 return background_value_.template read<T>();
377 }
378 };
379
380 /********************************************************/
381 /* */
382 /* labelMultiArray */
383 /* */
384 /********************************************************/
385
386 /** \brief Find the connected components of a MultiArray with arbitrary many dimensions.
387
388 See also \ref labelMultiArrayBlockwise() for a parallel version of this algorithm.
389
390 By specifying a background value in the \ref vigra::LabelOptions, this function
391 can also realize the behavior of \ref labelMultiArrayWithBackground().
392
393 <b> Declaration:</b>
394
395 \code
396 namespace vigra {
397
398 // specify parameters directly
399 template <unsigned int N, class T, class S1,
400 class Label, class S2,
401 class EqualityFunctor = std::equal_to<T> >
402 Label
403 labelMultiArray(MultiArrayView<N, T, S1> const & data,
404 MultiArrayView<N, Label, S2> labels,
405 NeighborhoodType neighborhood = DirectNeighborhood,
406 EqualityFunctor equal = std::equal_to<T>())
407
408 // specify parameters via LabelOptions
409 template <unsigned int N, class T, class S1,
410 class Label, class S2,
411 class Equal = std::equal<T> >
412 Label
413 labelMultiArray(MultiArrayView<N, T, S1> const & data,
414 MultiArrayView<N, Label, S2> labels,
415 LabelOptions const & options,
416 Equal equal = std::equal<T>());
417
418 }
419 \endcode
420
421 Connected components are defined as regions with uniform values.
422 Thus, the value type <tt>T</tt> of the input array \a data either
423 must be equality comparable, or an EqualityFunctor must be
424 provided that realizes the desired predicate. The destination
425 array's value type <tt>Label</tt> should be large enough to hold
426 the labels without overflow. Region numbers will form a consecutive
427 sequence starting at <b>one</b> and ending with the region number
428 returned by the function (inclusive).
429
430 Argument \a neighborhood specifies the type of connectivity used. It can
431 take the values <tt>DirectNeighborhood</tt> (which corresponds to
432 4-neighborhood in 2D and 6-neighborhood in 3D, default) or
433 <tt>IndirectNeighborhood</tt> (which corresponds to
434 8-neighborhood in 2D and 26-neighborhood in 3D).
435
436 Return: the highest region label used
437
438 <b> Usage:</b>
439
440 <b>\#include</b> \<vigra/multi_labeling.hxx\><br>
441 Namespace: vigra
442
443 \code
444 MultiArray<3,int> src(Shape3(w,h,d));
445 MultiArray<3,int> dest(Shape3(w,h,d));
446
447 // find 6-connected regions
448 int max_region_label = labelMultiArray(src, dest);
449
450 // find 26-connected regions
451 max_region_label = labelMultiArray(src, dest, IndirectNeighborhood);
452
453 // find 6-connected regions, ignore the background value 0
454 // (i.e. call labelMultiArrayWithBackground() internally)
455 max_region_label = labelMultiArray(src, dest,
456 LabelOptions().neighborhood(DirectNeighborhood)
457 .ignoreBackgroundValue(0));
458 \endcode
459
460 <b> Required Interface:</b>
461
462 \code
463 T t;
464 t == t // T is equality comparable
465
466 EqualityFunctor equal; // or a suitable functor is supplied
467 equal(t, t)
468 \endcode
469 */
doxygen_overloaded_function(template<...> unsigned int labelMultiArray)470 doxygen_overloaded_function(template <...> unsigned int labelMultiArray)
471
472 template <unsigned int N, class T, class S1,
473 class Label, class S2,
474 class Equal>
475 inline Label
476 labelMultiArray(MultiArrayView<N, T, S1> const & data,
477 MultiArrayView<N, Label, S2> labels,
478 NeighborhoodType neighborhood,
479 Equal equal)
480 {
481 vigra_precondition(data.shape() == labels.shape(),
482 "labelMultiArray(): shape mismatch between input and output.");
483
484 GridGraph<N, undirected_tag> graph(data.shape(), neighborhood);
485 return lemon_graph::labelGraph(graph, data, labels, equal);
486 }
487
488 template <unsigned int N, class T, class S1,
489 class Label, class S2>
490 inline Label
labelMultiArray(MultiArrayView<N,T,S1> const & data,MultiArrayView<N,Label,S2> labels,NeighborhoodType neighborhood=DirectNeighborhood)491 labelMultiArray(MultiArrayView<N, T, S1> const & data,
492 MultiArrayView<N, Label, S2> labels,
493 NeighborhoodType neighborhood = DirectNeighborhood)
494 {
495 return labelMultiArray(data, labels, neighborhood, std::equal_to<T>());
496 }
497
498 template <unsigned int N, class T, class S1,
499 class Label, class S2>
500 inline Label
labelMultiArray(MultiArrayView<N,T,S1> const & data,MultiArrayView<N,Label,S2> labels,LabelOptions const & options)501 labelMultiArray(MultiArrayView<N, T, S1> const & data,
502 MultiArrayView<N, Label, S2> labels,
503 LabelOptions const & options)
504 {
505 if(options.hasBackgroundValue())
506 return labelMultiArrayWithBackground(data, labels, options.getNeighborhood(),
507 options.template getBackgroundValue<T>());
508 else
509 return labelMultiArray(data, labels, options.getNeighborhood());
510 }
511
512 template <unsigned int N, class T, class S1,
513 class Label, class S2,
514 class Equal>
515 inline Label
labelMultiArray(MultiArrayView<N,T,S1> const & data,MultiArrayView<N,Label,S2> labels,LabelOptions const & options,Equal equal)516 labelMultiArray(MultiArrayView<N, T, S1> const & data,
517 MultiArrayView<N, Label, S2> labels,
518 LabelOptions const & options,
519 Equal equal)
520 {
521 if(options.hasBackgroundValue())
522 return labelMultiArrayWithBackground(data, labels, options.getNeighborhood(),
523 options.template getBackgroundValue<T>(),
524 equal);
525 else
526 return labelMultiArray(data, labels, options.getNeighborhood(), equal);
527 }
528
529 /********************************************************/
530 /* */
531 /* labelMultiArrayWithBackground */
532 /* */
533 /********************************************************/
534
535 /** \brief Find the connected components of a MultiArray with arbitrary many dimensions,
536 excluding the background from labeling.
537
538 From a user's point of view, this function is no longer needed because
539 \ref labelMultiArray() can realizes the same behavior when an appropriate
540 background value is specified in its \ref vigra::LabelOptions. Similarly,
541 \ref labelMultiArrayBlockwise() implements a parallel version of this algorithm.
542
543 <b> Declaration:</b>
544
545 \code
546 namespace vigra {
547
548 template <unsigned int N, class T, class S1,
549 class Label, class S2
550 class Equal = std::equal<T> >
551 Label
552 labelMultiArrayWithBackground(MultiArrayView<N, T, S1> const & data,
553 MultiArrayView<N, Label, S2> labels,
554 NeighborhoodType neighborhood = DirectNeighborhood,
555 T backgroundValue = T(),
556 Equal equal = std::equal<T>());
557
558 }
559 \endcode
560
561 This function is the same as \ref labelMultiArray(), except for
562 the additional parameter \a backgroundValue. Points in the input
563 array \a data with this value (which default to zero) are ignored
564 during labeling, and their output label is automatically set to
565 zero. Region numbers will be a consecutive sequence starting at
566 zero (when background was present) or at one (when no background
567 was present) and ending with the region number returned by the
568 function (inclusive).
569
570 Return: the number of non-background regions found (= highest region label,
571 because background has label 0)
572
573 <b> Usage:</b>
574
575 <b>\#include</b> \<vigra/multi_labeling.hxx\><br>
576 Namespace: vigra
577
578 \code
579 MultiArray<3, int> src(Shape3(w,h,d));
580 MultiArray<3, int> dest(Shape3(w,h,d));
581
582 // find 6-connected regions, ignoring background value zero (the default)
583 int max_region_label = labelMultiArrayWithBackground(src, dest);
584
585 // find 26-connected regions, using -1 as the background value
586 max_region_label = labelMultiArrayWithBackground(src, dest, IndirectNeighborhood, -1);
587 \endcode
588
589 <b> Required Interface:</b>
590
591 \code
592 T t, backgroundValue;
593 t == backgroundValue // T is equality comparable
594
595 EqualityFunctor equal; // or a suitable functor is supplied
596 equal(t, backgroundValue)
597 \endcode
598 */
doxygen_overloaded_function(template<...> unsigned int labelMultiArrayWithBackground)599 doxygen_overloaded_function(template <...> unsigned int labelMultiArrayWithBackground)
600
601 template <unsigned int N, class T, class S1,
602 class Label, class S2,
603 class Equal>
604 inline Label
605 labelMultiArrayWithBackground(MultiArrayView<N, T, S1> const & data,
606 MultiArrayView<N, Label, S2> labels,
607 NeighborhoodType neighborhood,
608 T backgroundValue,
609 Equal equal)
610 {
611 vigra_precondition(data.shape() == labels.shape(),
612 "labelMultiArrayWithBackground(): shape mismatch between input and output.");
613
614 GridGraph<N, undirected_tag> graph(data.shape(), neighborhood);
615 return lemon_graph::labelGraphWithBackground(graph, data, labels, backgroundValue, equal);
616 }
617
618 template <unsigned int N, class T, class S1,
619 class Label, class S2>
620 inline Label
labelMultiArrayWithBackground(MultiArrayView<N,T,S1> const & data,MultiArrayView<N,Label,S2> labels,NeighborhoodType neighborhood=DirectNeighborhood,T backgroundValue=T ())621 labelMultiArrayWithBackground(MultiArrayView<N, T, S1> const & data,
622 MultiArrayView<N, Label, S2> labels,
623 NeighborhoodType neighborhood = DirectNeighborhood,
624 T backgroundValue = T())
625 {
626 return labelMultiArrayWithBackground(data, labels, neighborhood, backgroundValue, std::equal_to<T>());
627 }
628
629 //@}
630
631 } // namespace vigra
632
633 #endif //VIGRA_MULTI_LABELING_HXX
634