immersion

Synopsis

Prototype:
y = immersion(f, offsets)
Description:
The immersion function implements the watershed transform algorithm by Vincent and Soille [VincentSoillePAMI1991]. It implements the algorithm on its strict form. This means that there is no guarantee of dividing lines between every catchment basin. Nevertheless, some lines are produced, and receive a label zero.
Definition:
Flooding-WT
Exploration:
Breadth-First

Algorithm

Algoritmo 1: Immersion

References

[VincentSoillePAMI1991]L. Vincent and P. Soille, ``Watersheds in digital spaces: An efficient algorithm based on immersion simulations,'' IEEE Transactions on Pattern Analysis and Machine Intelligence, vol. 13, no. 6, pp. 583–598, 1991.

Source Code

 1 from ipdp.common import *
 2 
 3 #constants
 4 MASK = -2
 5 INIT = -1
 6 WSHED = 0
 7 FICTITIOUS_PIXEL = -1
 8 
 9 
10 def immersion(im, offsets):
11 
12     # initialise variables
13     ws = wsImage(im)
14     N, im, lab, D = ws.begin(offsets)
15 
16     lab[:] = INIT
17     basins = 0
18     cur_dist = 0
19 
20     dist = ws.makeWorkCopy(0)
21 
22     queue = wsQueue()
23 
24     # "sorting"
25     levels = dict()
26     for p in D:
27 
28         if levels.has_key(im[p]):
29             levels[im[p]].append(p)
30         else:
31             levels[im[p]] = [p]
32 
33     # watershed
34     H = sorted(levels.keys())
35     for h in H:
36         for p in levels[h]:
37             lab[p] = MASK
38             for q in N(p):
39 
40                 if lab[q] > 0 or lab[q] == WSHED:
41                     dist[p] = 1
42                     queue.push(p)
43                     break
44 
45         cur_dist = 1
46         queue.push(FICTITIOUS_PIXEL)
47 
48         while True:
49 
50             p = queue.pop()
51             if p == FICTITIOUS_PIXEL:
52                 if queue.empty():
53                     break
54                 else:
55                     queue.push(FICTITIOUS_PIXEL)
56                     cur_dist += 1
57                     p = queue.pop()
58 
59             for q in N(p):
60 
61                 if dist[q] < cur_dist and (lab[q] > 0 or lab[q] == WSHED):
62                     if lab[q] > 0:
63                         if lab[p] == MASK or lab[p] == WSHED:
64                             lab[p] = lab[q]
65                         elif lab[p] != lab[q]:
66                             lab[p] = WSHED
67                     elif lab[p] == MASK:
68                         lab[p] = WSHED
69                 elif lab[q] == MASK and dist[q] == 0:
70                     dist[q] = cur_dist + 1
71                     queue.push(q)
72 
73         for p in levels[h]:
74             dist[p] = 0
75             if lab[p] == MASK:
76                 basins += 1
77                 queue.push(p)
78                 lab[p] = basins
79                 while not queue.empty():
80                     q = queue.pop()
81                     for u in N(q):
82                         if lab[u] != MASK:
83                             continue
84 
85                         queue.push(u)
86                         lab[u] = basins
87 
88     return ws.end()