<< unwrapp SIP - Scilab Image Processing Toolbox xgetpixel >>

SIP - Scilab Image Processing Toolbox >> SIP - Scilab Image Processing Toolbox > watershed


image segmentation with markers

Calling Sequence

w = watershed(img [,markers, nhood])


Grayscale image array. Preferably the gradient of the image to segment.
image with the markers (seeds). Each mark must have a unique label from 1 to N, where N is the number of marks. If markers is '-1' or omitted, then the regional minima of im will be taken as the markers.
A scalar. The connectivity to consider in the algorithm. May be 4 or 8, for 4-neighborhood or 8-neighborhood, resp.
image with the watershed regions (catchment basins), each with a unique number, from 1 to N. The elements labeled 1 belong to the first watershed region, the elements labeled 2 belong to the second basin, and so on. If markers != -1 or is omitted, the regions will have the same label as the corresponding supplied markers.


Function watershed computes the Watershed transform for image segmentation. This operation is also known as morphological sup-reconstruction. It is a region-growing algorithm which partitions the image into regions around each marker. Read the references or search the Internet for more theoretical information.

In the example we show a useful application of this operator: separation of overlapping objects.


// Suppose we have an image of many round or oval objects, such as a 
// microscope image of blood cells.  After thresholding the image, we
// end up with a binary image like this:

xset('auto clear', 'on');
a = gray_imread(SIPDIR + 'images/disks.bmp');

// We want to make the computer count the number of cells in the image, 
// but there are some circles that are overlapping, thus forming a single 
// connected component. Watershed is classically used for separating
// mingled objects like these. 
// First, calculate the distance transform:

a = 1-a;
d = bwdist(a);
d = normal(sqrt(d),255);   // normalize it to 0-255 range
imshow(d+1, hotcolormap(256));   

// The latter command shows the distance transform in shades from
// black to red to yellow to white. The brighter the color, the
// greater the distance of a point to the background.
// If you have the ENRICO toolbox, you can nicely plot the distance
// transform in 3D using sadesurf. ENRICO is not necessary for this
// example, but anyway you may download it at:
//       http://www.weizmann.ac.il/~fesegre/scistuff.html
// Note that the peaks of the distance transform are in the middle of
// each blob. The idea is to run watershed segmentation using these
// peaks as markers. For this, we invert the distance transform so
// that the peaks become the regional minima:

d = 255 - d;
imshow(d+1, hotcolormap(256));

// Now we "and" the distance transform with the original image, so
// that the background remains dark.

d = d .* a;

// Finally, run watershed segmentation. It automatically detects the
// regional minima for us:

w = watershed(d/255);   // input to watershed must be in 0-1 range
imshow(w, rand(256,3));
// 'w' is an image with a unique number for each watershed region.
// The imshow with a random colormap displays each region with a
// unique arbitrary color. Note how the regions were correctly separated by
// watershed, except for the hardest cases. It is extremely easy to
// count the number of regions:

n = max(w) - 1     // 26 regions minus the background

// The computer found 25 regions, but there are 20, an error of about 20%
// Let's improve this result. In the cases with many overlapping
// circles, the result would be perfect if it weren't for the small
// spurious regions. These are much smaller than the circles,
// so we can safely eliminate the regions with less than 100 pixels:

w2 = w;
for i=1:n
   f = find(w==i);        // coordinates of i-th region
   if size(f,'*') < 100
      w2(f) = 26;         // merge small regions with the background

imshow(w2, rand(256,3));   // note how the small regions are gone

// Now we count again, using a different way:

n = size( unique(w2), '*') - 1    // subtract 1 for the background

// Now it's 100% correct! We have an authomatic method wich is surprisingly 
// robust. This is specially useful to deal with bigger images in large 
// ammount. 
// Enjoy!
// TIP #1: Another way to improve the results is to do a median
// filtering of the distance transform. This will remove many spurious
// minima. Use mogrify(img, ['-median', '2']). Slight Gaussian smoothing also
// works well.
// TIP #2: for grayscale image segmentation, calculate the image
// gradient before watershed. Use edge(img,'sobel',-1) for this.
// TIP #3: use xgetpixel in a for loop (or something similar) to select
// the seed pixels to be used with watershed segmentation.

xset('auto clear', 'off');


"The Image Foresting Tranform: Theory, Algorithms, and Applications" A.X. Falcao, R.A. Lotufo, and J.Stolfi, IEEE T. Pattern Analysis and Machine Intelligence, (accepted for publication).

The IFT home page and the GIFT free software:


The original algorithm certainly that of Vincent and Soille, although it differs from the one we used in SIP/Animal:

"Watersheds in digital spaces: an efficient algorithm based on immersion simulations." IEEE T. Pattern Analysis and Machine Intelligence, 1991.



The latest version of the Scilab Image Processing toolbox can be found at


See Also

<< unwrapp SIP - Scilab Image Processing Toolbox xgetpixel >>