Vincent and Soille watershed algorithm by immersion [2]

1: procedure Watershed by Immersion

2: INPUT: digital grey scale image G = (D, E, )

3: OUTPUT: labeled watershed image lab on D.

4: #define INIT −1 (*initial value of lab image*)

5: #define MASK −2 (*initial value at each level*)

6: #define WSHED 0 (*label of the watershed pixels*)

7: #define FICTITIOUS (−1, −1) (*fictitious pixel *)

8: curlab ← 0 (*curlab is the current label*)

9: fifo_init (queue)

10: for all p ∈ D do

11: lab [p]← INIT; dist [p]← 0 (* is a work image of distances*)

12: end for

13: SORT pixels in increasing order of grey values (minimum hmin, maximum hmax)

14: (*Start the Flooding*)

15: for h = hmin to hmax do (*Geodesic SKIZ of level h −1 inside level h*)

16: for all with im [p] = do (*mask all pixels at level h*)

17: (*these are directly accessible because of the sorting step*)

18: lab [p] ← mask

19: if has a neighbor q with (lab [q] > 0 or lab [q] = WSHED) then

20: (*Initialize queue with neighbours at level h of current basins or watersheds*)

21: dist [p] ← 1; fifo_add (p, queue)

22: end if

23: end for

24: curdist ← 1; fifo_add (FICTITIOUS, queue)

25: loop (*extend basins*)

26: p ← fifo_remove (queue)

27: if p = FICTITIOUS then

28: if fifo_empty (queue) then

29: BREAK

30: else

31: fifo_add (FICTITIOUS, queue); curdist ← curdist + 1;

32: p ← fifo_remove (queue)

33: end if

34: end if

35: for all do (*labeling p by inspecting neighbors*)

36: if dist[q] < curdist and (lab[q] > 0 or lab[q] = WSHED) then

37: (*q belongs to an existing basin or to watersheds *)

38: if lab[q] > 0 then

39: if = MASK or lab[p] = WSHED then

40: lab[p] ← lab[q]

41: else if lab[p] ≠ lab[q] then

42: lab[p] ← WSHED

43: end if

44: else if lab[p] = MASK then

45: lab[p] WSHED

46: end if

47: else if lab[q] = MASK and dsit[q] = 0 then (*q is plateau pixel*)

58: dsit[q] ← curdist + 1; fifo_add (p, queue)

49: end if

50: end for