Lately I've been doing some image processing work using the most amazing AForge.NET open source library written in C#.

This library contains a whole bunch of useful filters you can play with but unfortunately it doesn't contain any filter for segmentation.

This is my implementation of Grayscale Immersion Watershed Segmentation in C# based on the

Vincent-Soille watershed algorithm which can be found in The Watershed Transform: Defnitions, Algorithms and Parallelization Strategies

paper by Jos B.T.M. Roerdink and Arnold Meijster on page 15. For watershed to work properly the areas you're trying to segment

have to be darker than the background. Blurring is also an important part of the image pre-processing since the algorithm

over segments non-blurred images. An example picture of watershed segmentation:

Original picture	Segmented picture

using System;
using System.Collections.Generic;
using System.Text;
namespace Segmentation
{
public class WatershedPixel
{
public int X;
public int Y;
public int Height;
// labels the pixel as belonging to a unique basin or as a part of the watershed line
public int Label;
// Distance is a work image of distances
public int Distance;
    <span class="kwrd">public</span> WatershedPixel()
    {
        <span class="kwrd">this</span>.X = -1;
        <span class="kwrd">this</span>.Y = -1;
        <span class="kwrd">this</span>.Height = -1000;
        <span class="kwrd">this</span>.Label = -1000;
        <span class="kwrd">this</span>.Distance = -1000;
    }

    <span class="kwrd">public</span> WatershedPixel(<span class="kwrd">int</span> x, <span class="kwrd">int</span> y)
    {
        <span class="kwrd">this</span>.X = x;
        <span class="kwrd">this</span>.Y = y;
        <span class="kwrd">this</span>.Height = -1000;
        <span class="kwrd">this</span>.Label = WatershedCommon.INIT;
        <span class="kwrd">this</span>.Distance = 0;
    }

    <span class="kwrd">public</span> WatershedPixel(<span class="kwrd">int</span> x, <span class="kwrd">int</span> y, <span class="kwrd">int</span> height)
    {
        <span class="kwrd">this</span>.X = x;
        <span class="kwrd">this</span>.Y = y;
        <span class="kwrd">this</span>.Height = height;
        <span class="kwrd">this</span>.Label = WatershedCommon.INIT;
        <span class="kwrd">this</span>.Distance = 0;
    }

    <span class="kwrd">public</span> <span class="kwrd">override</span> <span class="kwrd">bool</span> Equals(Object obj)
    {
        WatershedPixel p = (WatershedPixel)obj;
        <span class="kwrd">return</span> (X == p.X &amp;&amp; X == p.Y);
    }

    <span class="kwrd">public</span> <span class="kwrd">override</span> <span class="kwrd">int</span> GetHashCode()
    {
        <span class="kwrd">return</span> X;
    }
    <span class="kwrd">public</span> <span class="kwrd">override</span> <span class="kwrd">string</span> ToString()
    {
        <span class="kwrd">return</span> <span class="str">&#34;Height = &#34;</span> + Height + <span class="str">&#34;; X = &#34;</span> + X.ToString() + <span class="str">&#34;, Y = &#34;</span> + Y.ToString() + 

                   "; Label = " + Label.ToString() + "; Distance = " + Distance.ToString();
}
}
<span class="kwrd">public</span> <span class="kwrd">class</span> WatershedGrayscale : FilterGrayToGray
{
    <span class="preproc">#region</span> Variables
    <span class="kwrd">private</span> WatershedPixel FictitiousPixel = <span class="kwrd">new</span> WatershedPixel();
    <span class="kwrd">private</span> <span class="kwrd">int</span> _currentLabel = 0;
    <span class="kwrd">private</span> <span class="kwrd">int</span> _currentDistance = 0;
    <span class="kwrd">private</span> FifoQueue _fifoQueue = <span class="kwrd">new</span> FifoQueue();
    <span class="rem">// each pixel can be accesesd from 2 places: a dictionary for faster direct lookup of neighbouring pixels </span>
    <span class="rem">// or from a height ordered list</span>
    <span class="rem">// sorted array of pixels according to height</span>
    <span class="kwrd">private</span> List&lt;List&lt;WatershedPixel&gt;&gt; _heightSortedList;
    <span class="rem">// string in the form &#34;X,Y&#34; is used as a key for the dictionary lookup of a pixel</span>
    <span class="kwrd">private</span> Dictionary&lt;<span class="kwrd">string</span>, WatershedPixel&gt; _pixelMap;
    <span class="kwrd">private</span> <span class="kwrd">int</span> _watershedPixelCount = 0;
    <span class="kwrd">private</span> <span class="kwrd">int</span> _numberOfNeighbours = 8;
    <span class="kwrd">private</span> <span class="kwrd">bool</span> _borderInWhite;
    <span class="kwrd">private</span> <span class="kwrd">int</span> _pictureWidth = 0;
    <span class="kwrd">private</span> <span class="kwrd">int</span> _pictureHeight = 0;
    <span class="preproc">#endregion</span>

    <span class="preproc">#region</span> Constructors
    <span class="kwrd">public</span> WatershedGrayscale()
        : <span class="kwrd">this</span>(8)
    {}

    <span class="kwrd">public</span> WatershedGrayscale(<span class="kwrd">int</span> numberOfNeighbours)
    {
        <span class="kwrd">if</span> (numberOfNeighbours != 8 &amp;&amp; numberOfNeighbours != 4)
            <span class="kwrd">throw</span> <span class="kwrd">new</span> Exception(<span class="str">&#34;Invalid number of neighbour pixels to check. Valid values are 4 and 8.&#34;</span>);
        _borderInWhite = <span class="kwrd">true</span>;
        _numberOfNeighbours = numberOfNeighbours;
        _heightSortedList = <span class="kwrd">new</span> List&lt;List&lt;WatershedPixel&gt;&gt;(256);
        <span class="kwrd">for</span> (<span class="kwrd">int</span> i = 0; i &lt; 256; i++)
            _heightSortedList.Add(<span class="kwrd">new</span> List&lt;WatershedPixel&gt;());
    }
    <span class="preproc">#endregion</span>

    <span class="preproc">#region</span> Properties
    <span class="rem">/// &lt;summary&gt;</span>
    <span class="rem">/// number of neighbours to check for each pixel. valid values are 8 and 4</span>
    <span class="rem">/// &lt;/summary&gt;</span>
    <span class="kwrd">public</span> <span class="kwrd">int</span> NumberOfNeighbours
    {
        get { <span class="kwrd">return</span> _numberOfNeighbours; }
        set
        {
            <span class="kwrd">if</span> (<span class="kwrd">value</span> != 8 &amp;&amp; <span class="kwrd">value</span> != 4)
                <span class="kwrd">throw</span> <span class="kwrd">new</span> Exception(<span class="str">&#34;Invalid number of neighbour pixels to check. Valid values are 4 and 8.&#34;</span>);
            _numberOfNeighbours = <span class="kwrd">value</span>;
        }
    }

    <span class="rem">/// &lt;summary&gt;</span>
    <span class="rem">/// Number of labels/basins found</span>
    <span class="rem">/// &lt;/summary&gt;</span>
    <span class="kwrd">public</span> <span class="kwrd">int</span> LabelCount
    {
        get { <span class="kwrd">return</span> _currentLabel; }
        set { _currentLabel = <span class="kwrd">value</span>; }
    }

    <span class="rem">/// &lt;summary&gt;</span>
    <span class="rem">/// True: border is drawn in white. False: border is drawn in black</span>
    <span class="rem">/// &lt;/summary&gt;</span>
    <span class="rem">/// &lt;value&gt;&lt;/value&gt;</span>
    <span class="kwrd">public</span> <span class="kwrd">bool</span> BorderInWhite
    {
        get { <span class="kwrd">return</span> _borderInWhite; }
        set { _borderInWhite = <span class="kwrd">value</span>; }
    }
    <span class="preproc">#endregion</span>

    <span class="kwrd">private</span> <span class="kwrd">void</span> CreatePixelMapAndHeightSortedArray(BitmapData data)
    {
        _pictureWidth = data.Width;
        _pictureHeight = data.Height;
        <span class="rem">// pixel map holds every pixel thus size of (_pictureWidth * _pictureHeight)</span>
        _pixelMap = <span class="kwrd">new</span> Dictionary&lt;<span class="kwrd">string</span>, WatershedPixel&gt;(_pictureWidth * _pictureHeight);
        <span class="kwrd">unsafe</span>
        {                
            <span class="kwrd">int</span> offset = data.Stride - data.Width;
            <span class="kwrd">byte</span>* ptr = (<span class="kwrd">byte</span>*)(data.Scan0);

            <span class="rem">// get histogram of all values in grey = height</span>
            <span class="kwrd">for</span> (<span class="kwrd">int</span> y = 0; y &lt; data.Height; y++)
            {
                <span class="kwrd">for</span> (<span class="kwrd">int</span> x = 0; x &lt; data.Width; x++, ptr++)
                {
                    WatershedPixel p = <span class="kwrd">new</span> WatershedPixel(x, y, *ptr);
                    <span class="rem">// add every pixel to the pixel map</span>
                    _pixelMap.Add(p.X.ToString() + <span class="str">&#34;,&#34;</span> + p.Y.ToString(), p);
                    _heightSortedList[*ptr].Add(p);
                }
                ptr += offset;
            }
        }
        <span class="kwrd">this</span>._currentLabel = 0;
    }

    <span class="kwrd">private</span> <span class="kwrd">void</span> Segment()
    {            
        <span class="rem">// Geodesic SKIZ (skeleton by influence zones) of each level height</span>
        <span class="kwrd">for</span> (<span class="kwrd">int</span> h = 0; h &lt; _heightSortedList.Count; h++)
        {
            <span class="rem">// get all pixels for current height</span>
            <span class="kwrd">foreach</span> (WatershedPixel heightSortedPixel <span class="kwrd">in</span> _heightSortedList[h])
            {
                heightSortedPixel.Label = WatershedCommon.MASK;                    
                <span class="rem">// for each pixel on current height get neighbouring pixels</span>
                List&lt;WatershedPixel&gt; neighbouringPixels = GetNeighbouringPixels(heightSortedPixel);
                <span class="rem">// initialize queue with neighbours at level h of current basins or watersheds</span>
                <span class="kwrd">foreach</span> (WatershedPixel neighbouringPixel <span class="kwrd">in</span> neighbouringPixels)
                {
                    <span class="kwrd">if</span> (neighbouringPixel.Label &gt; 0 || neighbouringPixel.Label == WatershedCommon.WSHED)
                    {
                        heightSortedPixel.Distance = 1;
                        _fifoQueue.AddToEnd(heightSortedPixel);
                        <span class="kwrd">break</span>;
                    }                        
                }
            }
            _currentDistance = 1;
            _fifoQueue.AddToEnd(FictitiousPixel);
            <span class="rem">// extend basins</span>
            <span class="kwrd">while</span> (<span class="kwrd">true</span>)
            {
                WatershedPixel p = _fifoQueue.RemoveAtFront();
                <span class="kwrd">if</span> (p.Equals(FictitiousPixel))
                {
                    <span class="kwrd">if</span> (_fifoQueue.IsEmpty)
                        <span class="kwrd">break</span>;
                    <span class="kwrd">else</span>
                    {
                        _fifoQueue.AddToEnd(FictitiousPixel);
                        _currentDistance++;
                        p = _fifoQueue.RemoveAtFront();
                    }
                }
                List&lt;WatershedPixel&gt; neighbouringPixels = GetNeighbouringPixels(p);
                <span class="rem">// labelling p by inspecting neighbours</span>
                <span class="kwrd">foreach</span> (WatershedPixel neighbouringPixel <span class="kwrd">in</span> neighbouringPixels)
                {
                    <span class="rem">// neighbouringPixel belongs to an existing basin or to watersheds</span>
                    <span class="rem">// in the original algorithm the condition is:</span>
                    <span class="rem">//   if (neighbouringPixel.Distance &lt; currentDistance &amp;&amp; </span>
                    <span class="rem">//      (neighbouringPixel.Label &gt; 0 || neighbouringPixel.Label == WatershedCommon.WSHED))</span>
                    <span class="rem">//   but this returns incomplete borders so the this one is used                        </span>
                    <span class="kwrd">if</span> (neighbouringPixel.Distance &lt;= _currentDistance &amp;&amp; 
                       (neighbouringPixel.Label &gt; 0 || neighbouringPixel.Label == WatershedCommon.WSHED))
                    {
                        <span class="kwrd">if</span> (neighbouringPixel.Label &gt; 0)
                        {
                            <span class="rem">// the commented condition is also in the original algorithm </span>
                            <span class="rem">// but it also gives incomplete borders</span>
                            <span class="kwrd">if</span> (p.Label == WatershedCommon.MASK <span class="rem">/*|| p.Label == WatershedCommon.WSHED*/</span>) 
                                p.Label = neighbouringPixel.Label;
                            <span class="kwrd">else</span> <span class="kwrd">if</span> (p.Label != neighbouringPixel.Label)
                            {
                                p.Label = WatershedCommon.WSHED;
                                _watershedPixelCount++;
                            }
                        }
                        <span class="kwrd">else</span> <span class="kwrd">if</span> (p.Label == WatershedCommon.MASK)
                        {
                            p.Label = WatershedCommon.WSHED;
                            _watershedPixelCount++;
                        }
                    }
                    <span class="rem">// neighbouringPixel is plateau pixel</span>
                    <span class="kwrd">else</span> <span class="kwrd">if</span> (neighbouringPixel.Label == WatershedCommon.MASK &amp;&amp; neighbouringPixel.Distance == 0)
                    {
                        neighbouringPixel.Distance = _currentDistance + 1;
                        _fifoQueue.AddToEnd(neighbouringPixel);
                    }
                }
            }
            <span class="rem">// detect and process new minima at height level h</span>
            <span class="kwrd">foreach</span> (WatershedPixel p <span class="kwrd">in</span> _heightSortedList[h])
            {
                <span class="rem">// reset distance to zero</span>
                p.Distance = 0;
                <span class="rem">// if true then p is inside a new minimum </span>
                <span class="kwrd">if</span> (p.Label == WatershedCommon.MASK)
                {
                    <span class="rem">// create new label</span>
                    _currentLabel++;
                    p.Label = _currentLabel;
                    _fifoQueue.AddToEnd(p);
                    <span class="kwrd">while</span> (!_fifoQueue.IsEmpty)
                    {
                        WatershedPixel q = _fifoQueue.RemoveAtFront();
                        <span class="rem">// check neighbours of q</span>
                        List&lt;WatershedPixel&gt; neighbouringPixels = GetNeighbouringPixels(q);
                        <span class="kwrd">foreach</span> (WatershedPixel neighbouringPixel <span class="kwrd">in</span> neighbouringPixels)
                        {
                            <span class="kwrd">if</span> (neighbouringPixel.Label == WatershedCommon.MASK)
                            {
                                neighbouringPixel.Label = _currentLabel;
                                _fifoQueue.AddToEnd(neighbouringPixel);
                            }
                        }
                    }
                }
            }
        }
    }

    <span class="kwrd">private</span> List&lt;WatershedPixel&gt; GetNeighbouringPixels(WatershedPixel centerPixel)
    {
        List&lt;WatershedPixel&gt; temp = <span class="kwrd">new</span> List&lt;WatershedPixel&gt;();
        <span class="kwrd">if</span> (_numberOfNeighbours == 8)
        {
            <span class="rem">/*</span>

                CP = Center pixel
                (X,Y) – get all 8 connected 
                |-1,-1|0,-1|1,-1|
                |-1, 0| CP |1, 0|
                |-1,+1|0,+1|1,+1|
                */             

// -1, -1                
if ((centerPixel.X - 1) >= 0 && (centerPixel.Y - 1) >= 0)
temp.Add(_pixelMap[(centerPixel.X - 1).ToString() + "," + (centerPixel.Y - 1).ToString()]);
//  0, -1
if ((centerPixel.Y - 1) >= 0)
temp.Add(_pixelMap[centerPixel.X.ToString() + "," + (centerPixel.Y - 1).ToString()]);
//  1, -1
if ((centerPixel.X + 1) < _pictureWidth && (centerPixel.Y - 1) >= 0)
temp.Add(_pixelMap[(centerPixel.X + 1).ToString() + "," + (centerPixel.Y - 1).ToString()]);
// -1, 0
if ((centerPixel.X - 1) >= 0)
temp.Add(_pixelMap[(centerPixel.X - 1).ToString() + "," + centerPixel.Y.ToString()]);
//  1, 0
if ((centerPixel.X + 1) < _pictureWidth)
temp.Add(_pixelMap[(centerPixel.X + 1).ToString() + "," + centerPixel.Y.ToString()]);
// -1, 1
if ((centerPixel.X - 1) >= 0 && (centerPixel.Y + 1) < _pictureHeight)
temp.Add(_pixelMap[(centerPixel.X - 1).ToString() + "," + (centerPixel.Y + 1).ToString()]);
//  0, 1
if ((centerPixel.Y + 1) < _pictureHeight)
temp.Add(_pixelMap[centerPixel.X.ToString() + "," + (centerPixel.Y + 1).ToString()]);
//  1, 1
if ((centerPixel.X + 1) < _pictureWidth && (centerPixel.Y + 1) < _pictureHeight)
temp.Add(_pixelMap[(centerPixel.X + 1).ToString() + "," + (centerPixel.Y + 1).ToString()]);
}
else
{
/*
                CP = Center pixel, N/A = not used
                (X,Y) – get only 4 connected 
                | N/A |0,-1| N/A |
                |-1, 0| CP |+1, 0|
                | N/A |0,+1| N/A |
                */
//  -1, 0
if ((centerPixel.X - 1) >= 0)
temp.Add(_pixelMap[(centerPixel.X - 1).ToString() + "," + centerPixel.Y.ToString()]);
//  0, -1
if ((centerPixel.Y - 1) >= 0)
temp.Add(_pixelMap[centerPixel.X.ToString() + "," + (centerPixel.Y - 1).ToString()]);
//  1, 0
if ((centerPixel.X + 1) < _pictureWidth)
temp.Add(_pixelMap[(centerPixel.X + 1).ToString() + "," + centerPixel.Y.ToString()]);
//  0, 1
if ((centerPixel.Y + 1) < _pictureHeight)
temp.Add(_pixelMap[centerPixel.X.ToString() + "," + (centerPixel.Y + 1).ToString()]);
}
return temp;
}
    <span class="kwrd">private</span> <span class="kwrd">void</span> DrawWatershedLines(BitmapData data)
    {
        <span class="kwrd">if</span> (_watershedPixelCount == 0)
            <span class="kwrd">return</span>;

        <span class="kwrd">byte</span> watershedColor = 255;
        <span class="kwrd">if</span> (!_borderInWhite)
            watershedColor = 0;
        <span class="kwrd">unsafe</span>
        {
            <span class="kwrd">int</span> offset = data.Stride - data.Width;
            <span class="kwrd">byte</span>* ptr = (<span class="kwrd">byte</span>*)(data.Scan0);

            <span class="kwrd">for</span> (<span class="kwrd">int</span> y = 0; y &lt; data.Height; y++)
            {
                <span class="kwrd">for</span> (<span class="kwrd">int</span> x = 0; x &lt; data.Width; x++, ptr++)
                {
                    <span class="rem">// if the pixel in our map is watershed pixel then draw it</span>
                    <span class="kwrd">if</span> (_pixelMap[x.ToString() + <span class="str">&#34;,&#34;</span> + y.ToString()].Label == WatershedCommon.WSHED)
                        *ptr = watershedColor;
                }
                ptr += offset;
            }
        }
    }

    <span class="kwrd">protected</span> <span class="kwrd">override</span> <span class="kwrd">void</span> ProcessFilter(BitmapData imageData)
    {
        CreatePixelMapAndHeightSortedArray(imageData);
        Segment();
        DrawWatershedLines(imageData);
    }
}

<span class="kwrd">public</span> <span class="kwrd">class</span> FifoQueue
{
    List&lt;WatershedPixel&gt; queue = <span class="kwrd">new</span> List&lt;WatershedPixel&gt;();

    <span class="kwrd">public</span> <span class="kwrd">int</span> Count
    {
        get { <span class="kwrd">return</span> queue.Count; }
    }

    <span class="kwrd">public</span> <span class="kwrd">void</span> AddToEnd(WatershedPixel p)
    {
        queue.Add(p);
    }

    <span class="kwrd">public</span> WatershedPixel RemoveAtFront()
    {
        WatershedPixel temp = queue[0];
        queue.RemoveAt(0);
        <span class="kwrd">return</span> temp;
    }

    <span class="kwrd">public</span> <span class="kwrd">bool</span> IsEmpty
    {
        get { <span class="kwrd">return</span> (queue.Count == 0); }
    }

    <span class="kwrd">public</span> <span class="kwrd">override</span> <span class="kwrd">string</span> ToString()
    {
        <span class="kwrd">return</span> <span class="kwrd">base</span>.ToString() + <span class="str">&#34; Count = &#34;</span> + queue.Count.ToString();
    }
}

<span class="kwrd">public</span> <span class="kwrd">class</span> WatershedCommon
{
    <span class="preproc">#region</span> Constants
    <span class="kwrd">public</span> <span class="kwrd">const</span> <span class="kwrd">int</span> INIT = -1;
    <span class="kwrd">public</span> <span class="kwrd">const</span> <span class="kwrd">int</span> MASK = -2;
    <span class="kwrd">public</span> <span class="kwrd">const</span> <span class="kwrd">int</span> WSHED = 0;
    <span class="preproc">#endregion</span>
}

}

This filter will also become the part of the AForge.NET library in the future so to use this algorithm you have to have it installed.

The performance can be improved as I've gone for the readability factor since performance wasn't very important on my project.