Home · Pages · Index · Overviews

Watershed Module Reference

Routines to produce a watershed partition of an image. More...

 #include <water.shed.h>

Routines

Partition * Build_WatershedG (Pixel_APart *image I, boolean iscon2n, boolean color)

Label_Array * Label_Watershed (Pixel_APart *image, Integer_Array *labels RO,
                                int *nbasins O, boolean iscon2n)

Partition * Build_Seeded_WatershedG (Pixel_APart *image I, boolean iscon2n, Vector *seeds)

void Average_Watershed_Heights (Partition *shed, int *num O, int *den O, int64 *sqr O, int null)

Detailed Description

The idea of a watershed partition of an image is well established and informally appeals to exactly the idea of a watershed in nature, i.e. a region for which all rain drops landing in the region flow to the same minimum basin. Proceeding more formally, a minimum is a connected set of voxels of the same value for which every adjacent voxel not in the set has a greater value. In a watershed partition, there is a catchment basin enclosing each minimum, that is the connected set of all voxels for which there is a connected path to its minimum along which values do not increase (a rain drop flows downhill to the minimum). For continuous domains this definition would suffice, but for integer valued domains there is the problem of ties which is solved in the algorithm of Vincent and Soile (IEEE Trans. on Pat. Analysis & Mach. Learning 13, 6 (1991)), by having the water flow to the geodesically closer basin. A watershed partition is a set of non-overlapping catchment basins that cover the entire image.

A watershed partition of a UINT8_ or UINT16_TYPE array or slice image, using a variant of the algorithm of Vincent and Soille, is generated by Build_Watershed. The partition is returned as a Partition object whose label array is a coloring of the watershed in the smallest type (UINT8, UINT16, or UINT32) that accommodates the number of colors in the labeling. A current important limitation is that the routine restricts the size of image to 231-1. Note carefully that as a consequence all indices on this page are simply of type int, as opposed to the int64 required in the rest of the library. The routine consumes 4 bytes per voxel in working storage along with two integer stacks whose size is generally a small fraction of the image size. In the future the 32-bit limitation on size will be removed as large memories continue to become less expensive and more common place.

Build_Watershed can actually be broken down into two distinct parts. In the first part of the process, a label field of the watershed partition of the image is built in a UINT32_TYPE Label_Array of the same shape as the image or its underlying array wherein the pixels of catchment basin i are labeled with i+1 for i ∈ [0,nbasins-1]. This is accomplished with a call to Label_Watershed. In the second phase, a Partition object is produced by calling Make_Partition (described in the Partition class page) with the watershed labeling produced in the first step. In fact, the entire Partition class was initially developed for watersheds and was later identified as an independently useful concept. The depth of a P_Vertex is the value of its catchment basin minimum, and the height of an P_Edge is the height of the dam or ridge that separates two catchment basins of the watershed. Thus one can get a basic representation of a watershed partition in the label field array if that is all one desires, or one can further build a Partition graph for the watershed decomposition. We illustrate this by giving the code for Build_Watershed.

 Partition *Build_Watershed(Pixel_APart *image, boolean iscon2n, boolean color)
 { Label_Array *h;
   Partition   *s;
   int          n;

   h = Make_Array_With_Shape(PLAIN_KIND,INT32_TYPE,AForm_Shape(image));
   Label_Watershed(image,h,&n,iscon2n);
   s = Make_Partition(image,h,n,iscon2n,color);
   if (Get_Partition_Color_Count(s) < 256)
     Convert_Array_Inplae(h,PLAIN_KIND,UINT8_TYPE,8,0);
   else if (Get_Partition_Color_Count(s) < 0x10000)
     Convert_Array_Inplae(h,PLAIN_KIND,UINT16_TYPE,16,0);
   Pack_Array(h);
   return (s);
 }

The watershed algorithm tends to oversegment an image as it is very sensitive to small fluctuations in signal. One approach is to smooth the image with a Gaussian filter prior to applying the watershed, but a far superior approach is to progressively merge watersheds for which the difference between their depth and the dam height separating them is small. To faciliatate this the Partition class has routines that allow one to progressively merge regions in either the fixed order of dam height (Static_Collapse), or more flexibly in any dynamic order defined by the user (General_Collapse). For example, for a watershed partition w, the routine Average_Watershed_Heights fills in vectors num and den of length Get_Partition_Edge_Count(p) so that for the jth edge, num[j] is the sum of the pixel values along the boundary separating the two basins, den[j] is the number of pixels in along the boundary, and therefore num[j]/den[j] is the average dam height of the boundary separating the two regions of the edge. In addition, if sqr is not NULL, then upon return sqr[j] contains the sum of the square of each pixel value along the boundary separating the two basins, so that from this one can compute the variance and standard deviation of pixels along a boundary if desired. Moreover, pixels of value null are not included in the sums and counts, allowing one to mask certain pixels if desired (otherwise pick null to be out of the range of values in the image). The average dam height between two basins may be a more useful measure of separability in cases where there are 'weaknesses' in the barrier, such as irregularly stained biological boundaries. One may then wish to collapse basins (regions) in order of average dam height, where a little thought reveals that average dam height changes as basins are merged, because this implies dams (boundaries) must also be merged. A complete example of such a merging using the values provided by Average_Watershed_Heights can be found here.


Routine Documentation

Partition * Build_WatershedG (Pixel_APart *image I, boolean iscon2n, boolean color)

Compute the watershed partition of the UINT8_ or UINT16_TYPE array or slice image with respect to the connectivity indicated by iscon2n, using a variant of the algorithm of Vincent and Soille, and generate a Partition object modeling the partition. The Partition object creates a new reference to image and its label array is a coloring iff color is true. A current important limitation of the Partition class is that this routine restricts the size of image to 231-1. Note carefully that as a consequence all indices for this class are simply of type int, as opposed to the int64 required in the rest of the library. The routine takes linear time in the size of the array and consumes 4 bytes per pixel in working storage along with two integer stacks whose size is generally a small fraction of the image size. In the future the 32-bit limitation on size will be removed as large memories continue to become less expensive and more common place.

Label_Array * Label_Watershed (Pixel_APart *image, Integer_Array *labels RO,
                                int *nbasins O, boolean iscon2n)

Produce a label field of the watershed partition of image in the supplied integer array labels with respect to the connectivity indicated by iscon2n. That is, in the integer array labels, the pixels of catchment basin i are labeled with i+1 for i ∈ [0,nbasins-1]. The number of basins is returned in nbasins. This label field is in effect a watershed partitioning. This label field, can be used to produce a Partition graph by Make_Partition.

If image is an array, then labels must have the same shape as image. On the otherhand, if image is a slice of an array, then labels may have either the shape of the slice or the shape of the array of the slice. In the latter case, the routine only accesses the pixels of labels corresponding to the pixels of the slice, thus permitting one to compute the watersheds of disjoint slices that cover an array in different threads that share a single label field. This is about 15% more efficient than having each slice use its own label field (proportional to the size of the slice) due to the indexing issues involved with correlating traversals in arrays of two different sizes.

Partition * Build_Seeded_WatershedG (Pixel_APart *image I, boolean iscon2n, Vector *seeds)

Seeded watershed. Build a watershed graph of image w.r.t iscon2n connectivity, but all basins not involving an index in the vector "seeds" are merged with the seed- containing-basins with the lowest dam height between them. The type of Vector must be of OFFS_TYPE and the integers should be in [0,image->size-1].

void Average_Watershed_Heights (Partition *shed, int *num O, int *den O, int64 *sqr O, int null)

Compute the average watershed height between any two adjacent basins. The caller supplies the arrays num and den which must be of size Get_Partition_Edge(shed) (the number of edges in the watershed graph). For the P_Vertex with index d, the average watershed height between its two regions is num[d]/den[d]. In addition, if sqr is not NULL, then upon return sqr[j] contains the sum of the square of each pixel value along the boundary separating the two basins, so that from this one can compute the variance and standard deviation of pixels along a boundary if desired. Moreover, pixels of value null are not included in the sums and counts, allowing one to mask certain pixels if desired (otherwise pick null to be out of the range of values in the image).