watershed

Watershed transform

Syntax

L = watershed(A)
L = watershed(A, conn)

Description

L = watershed(A) computes a label matrix identifying the watershed regions of the input matrix A, which can have any dimension. The elements of L are integer values greater than or equal to 0. The elements labeled 0 do not belong to a unique watershed region. These are called watershed pixels. The elements labeled 1 belong to the first watershed region, the elements labeled 2 belong to the second watershed region, and so on.

By default, watershed uses 8-connected neighborhoods for 2-D inputs and 26-connected neighborhoods for 3-D inputs. For higher dimensions, watershed uses the connectivity given by conndef(ndims(A),'maximal').

L = watershed(A, conn) specifies the connectivity to be used in the watershed computation. conn can have any of the following scalar values.

Value

Meaning

Two-dimensional connectivities

4

4-connected neighborhood

8

8-connected neighborhood

Three-dimensional connectivities

6

6-connected neighborhood

18

18-connected neighborhood

26

26-connected neighborhood

Connectivity can be defined in a more general way for any dimension by using for conn a 3-by-3-by- ...-by-3 matrix of 0's and 1's. The 1-valued elements define neighborhood locations relative to the center element of conn. Note that conn must be symmetric about its center element.

Remarks

The watershed transform algorithm used by this function changed in version 5.4 (R2007a) of the Image Processing Toolbox™ software. The previous algorithm occasionally produced labeled watershed basins that were not contiguous. If you need to obtain the same results as the previous algorithm, use the function watershed_old.

Class Support

A can be a numeric or logical array of any dimension, and it must be nonsparse. The output array L is an unsigned integer type.

Examples

2-D Example

  1. Make a binary image containing two overlapping circular objects.

    center1 = -10;
    center2 = -center1;
    dist = sqrt(2*(2*center1)^2);
    radius = dist/2 * 1.4;
    lims = [floor(center1-1.2*radius) ceil(center2+1.2*radius)];
    [x,y] = meshgrid(lims(1):lims(2));
    bw1 = sqrt((x-center1).^2 + (y-center1).^2) <= radius;
    bw2 = sqrt((x-center2).^2 + (y-center2).^2) <= radius;
    bw = bw1 | bw2;
    figure, imshow(bw,'InitialMagnification','fit'), title('bw')
    
  2. Compute the distance transform of the complement of the binary image.

    D = bwdist(~bw);
    figure, imshow(D,[],'InitialMagnification','fit')
    title('Distance transform of ~bw')
  3. Complement the distance transform, and force pixels that don't belong to the objects to be at -Inf.

    D = -D;
    D(~bw) = -Inf;
  4. Compute the watershed transform and display the resulting label matrix as an RGB images.

    L = watershed(D);
    rgb = label2rgb(L,'jet',[.5 .5 .5]);
    figure, imshow(rgb,'InitialMagnification','fit')
    title('Watershed transform of D')

3-D Example

  1. Make a 3-D binary image containing two overlapping spheres.

    center1 = -10;
    center2 = -center1;
    dist = sqrt(3*(2*center1)^2);
    radius = dist/2 * 1.4;
    lims = [floor(center1-1.2*radius) ceil(center2+1.2*radius)];
    [x,y,z] = meshgrid(lims(1):lims(2));
    bw1 = sqrt((x-center1).^2 + (y-center1).^2 + ...
               (z-center1).^2) <= radius;
    bw2 = sqrt((x-center2).^2 + (y-center2).^2 + ...
               (z-center2).^2) <= radius;
    bw = bw1 | bw2;
    figure, isosurface(x,y,z,bw,0.5), axis equal, title('BW')
    xlabel x, ylabel y, zlabel z
    xlim(lims), ylim(lims), zlim(lims)
    view(3), camlight, lighting gouraud
    
  2. Compute the distance transform.

    D = bwdist(~bw);
    figure, isosurface(x,y,z,D,radius/2), axis equal
    title('Isosurface of distance transform')
    xlabel x, ylabel y, zlabel z
    xlim(lims), ylim(lims), zlim(lims)
    view(3), camlight, lighting gouraud
  3. Complement the distance transform, force nonobject pixels to be -Inf, and then compute the watershed transform.

    D = -D;
    D(~bw) = -Inf;
    L = watershed(D);
    figure
    isosurface(x,y,z,L==2,0.5)
    isosurface(x,y,z,L==3,0.5)
    axis equal
    title('Segmented objects')
    xlabel x, ylabel y, zlabel z
    xlim(lims), ylim(lims), zlim(lims)
    view(3), camlight, lighting gouraud

More About

expand all

Algorithms

watershed uses the Fernand Meyer algorithm [1].

References

[1] Meyer, Fernand, "Topographic distance and watershed lines," Signal Processing , Vol. 38, July 1994, pp. 113-125.

Was this topic helpful?