Simple implementation of watershed Algorithms

Without looking at the watershed, the VC image processing book rarely introduces this algorithm, and cximage also finds related algorithms. Let's look at matlab and opencv. MATLAB encapsulates the watershed algorithm and does not look at the source code. opencv is too high to write, and the level difference is not added to opencv. I found the source code of several algorithms on the Internet. One is written by a foreigner using a template, and the other by a Chinese. Foreigners write Nb, but they cannot read it. Chinese people read it in the fog. Finally, they can find a way to write it.

 

First of all, we declare that this is a simple implementation and only applicable to some simple situations. The effect is acceptable, but there are deficiencies. In some cases, the split may split the original single prospect. Take a look at cainiao like me. High pass can guide.

 

What is watershed first,

Definition:

Image segmentation refers to the separation of different regions with special meanings in the image. These regions are mutually independent, and each region meets the consistency of a specific region.

This code is not fully implemented and can only meet the current work needs.

There are many principles of watershed on the Internet. I just want to explain my simple understanding of this algorithm: this algorithm can solve the problem. You can divide it into prospects, but there are some interference or other reasons, the split foreground can be connected to one another, and the foreground cannot be differentiated. At this time, we can use this watershed algorithm to help.

(Let's take a look at a few pictures of the principles of watershed. This is too troublesome for typing. Let's look at the pictures directly)

 

It took me more than half a day to understand the algorithm. My understanding (there may be errors or deficiencies) Is that m is the basin with the smallest gray value (or the largest relative background) as the lowest, in fact, that is, the prospects are the most away from the background (this certainly lacks understanding. I hope you will be guided by the audience ). According to the continuous addition of water in this basin, the water level of the basin is constantly improved (that is, it is constantly approaching the background, there may be conflicts with other regions. At this time, we need to adopt a harmonious attitude to prevent such incidents, that is, to build a dam in a conflicting place. This dam is what we need to find.

There are three possibilities for raising the Basin Water Level

1: The strength above the basin, that is, the gray scale of the region is adjacent to the basin and has nothing to do with other forces. This situation must be swallowed up.

2: adjacent strength of the basin, that is, the gray scale of the region is not only adjacent to the basin, but also related to other basins. This is the dam we need to find.

3: when the water level reaches a certain height, the new basin appears, neither the region nor the original basin. The water level is the lowest level in the region. At this time, a new basin was created, that is, a new evil force emerged.

Based on the above understanding, I wrote a Watershed Based on the image distance change, which is suitable for the Cross of some small areas.

First, we need to binarization the image to split the location areas of interest. Binarization, as long as they can be separated.

The second step is distance transformation. The Code comes from the Principle and Implementation of a VC image algorithm. Let's look at this title.

// Calculate the distance value of the binarization image <br/> bool processing: distance (byte * & pdata) <br/>{</P> <p> // array for storing pixel values <br/> If (pdata! = NULL) <br/>{< br/> return false; <br/>}< br/> int * pnbinary = new int [m_size]; <br/> int * pnstore = new int [m_size]; <br/> pdata = new byte [m_size]; <br/> // binarization the image <br/> for (Int J = 0; j <m_height; ++ J) <br/>{< br/> for (INT I = 0; I <m_biwidth; ++ I) <br/>{< br/> // The White pixel is the background, save as 0 <br/> If (m_pbinary [I + J * m_biwidth]) <br/> {<br/> pnbinary [I + J * m_biwidth] = 0; <br/> pnstore [I + J * m_biwidth] = 0; <br/>}< br/> // black pixels are saved as 1 <br/> else <br/> {<br/> pnbinary [I + J * m_biwidth] = 1; <br/> pnstore [I + J * m_biwidth] = 1; <br/>}</P> <p> int S = 1; </P> <p> while (S = 1) <br/> {<br/> S = 0; </P> <p> // accumulate the distance <br/> for (Int J = 1; j <m_height-1; ++ J) <br/> {<br/> for (INT I = 1; I <m_biwidth-1; ++ I) <br/>{</P> <p> // if it is the background, perform the next loop <br/> If (pnbinary [I + J * m_biwidth] = 0) <br/> continue; </P> <p> // if the current pixel value is the same as the surrounding value, add a pixel value <br/> If (pnbinary [I + J * m_biwidth] = pnbinary [m_biwidth * (J-1) + I] & pnbinary [I + J * m_biwidth] = pnbinary [m_biwidth * (J + 1) + I]) <br/> If (pnbinary [I + J * m_biwidth] = pnbinary [m_biwidth * j + I-1] & pnbinary [I + J * m_biwidth] = pnbinary [m_biwidth * J + I + 1]) <br/> {<br/> pnstore [m_biwidth * j + I] ++; <br/> S = 1; <br/>}</P> <p> // Save the current result before the next round of loop <br/> (Int J = 0; j <m_height; j ++) <br/> for (INT I = 1; I <m_biwidth; I ++) <br/> pnbinary [m_biwidth * j + I] = pnstore [m_biwidth * j + I]; <br/>}</P> <p> for (Int J = 0; j <m_height; ++ J) <br/> {<br/> for (INT I = 0; I <m_biwidth; ++ I) <br/> {<br/> // pdata [I + J * m_biwidth] = max (0, (25-pnstore [J * m_biwidth + I]) * 10 + 5); // <br/> pdata [I + J * m_biwidth] = max (0,255-min (255, pnbinary [I + J * m_biwidth]); <br/>}</P> <p> Delete pnstore; <br/> Delete pnbinary; </P> <p> return true; <br/>}

This function is written in a class, so the basic parameters of some images are members of the class. I will not explain this function.

// Watershed processing binarization Distance Image <br/> bool processing: watershed (byte * data, int * hist) <br/>{< br/> // modify to try mechanism <br/> If (Data = NULL | hist = NULL) <br/>{< br/> return false; <br/>}< br/> queue <point> * statistic = new queue <point> [1, 255]; // record the same level basin queue. The background is not processed <br/> int * wdata = new int [m_size]; <br/> memset (wdata, 255, m_size * sizeof (INT); <br/> bool cue = true; // variable mark <br/> int start_queue = 0; // mark the start queue <br/> int lab_cout = 256; // basin id </P & Gt; <p> point pt; </P> <p> const int wshed = 255; // dam <br/> // all data is pushed into the queue // edge points are not processed <br/> for (INT I = 1; I <m_biwidth-1; ++ I) <br/> {<br/> for (Int J = 1; j <m_height-1; ++ J) <br/> {<br/> wdata [I + J * m_biwidth] = data [I + J * m_biwidth]; <br/> If (data [I + J * m_biwidth] <255) <br/>{< br/> PT. X = I; <br/> PT. y = J; <br/> statistic [DATA [I + J * m_biwidth]. push (PT ); <br/>}< br/> // complete the initial watermark code for the first low position basin. Description 1 <br/> (; start_queue <256 & hist [Start_queue] = 0; ++ start_queue ); </P> <p> // complete searching for all source water level basin codes. Description 2 <br/> while (hist [start_queue]> 0) <br/> {<br/> Pt = statistic [start_queue]. front (); <br/> statistic [start_queue]. pop (); <br/> wdata [PT. X + Pt. y * m_biwidth] = lab_cout; <br/> -- hist [start_queue]; <br/> cue = true; <br/> while (cue) <br/>{< br/> cue = false; <br/> int size = (INT) statistic [start_queue]. size (); <br/> for (INT I = 0; I <size; ++ I) <br/> {<br/> int T; <br/> Pt = statistic [start_queue]. front (); <br/> statistic [start_queue]. pop (); <br/> // search for lab in 8 directions <br/> T = wdata [PT. X + (pt. y-1) * m_biwidth]; <br/> If (T> 255) <br/> {<br/> wdata [PT. X + Pt. y * m_biwidth] = lab_cout; <br/> -- hist [start_queue]; <br/> cue = true; <br/> continue; <br/>}< br/> T = wdata [PT. x-1 + (pt. y-1) * m_biwidth]; <br/> If (T> 255) <br/> {<br/> wdata [PT. X + Pt. y * m_biwidth] = lab_cout; <br/> -- hist [start_queue]; <Br/> cue = true; <br/> continue; <br/>}< br/> T = wdata [PT. x-1 + Pt. y * m_biwidth]; <br/> If (T> 255) <br/> {<br/> wdata [PT. X + Pt. y * m_biwidth] = lab_cout; <br/> -- hist [start_queue]; <br/> cue = true; <br/> continue; <br/>}< br/> T = wdata [PT. x-1 + (pt. Y + 1) * m_biwidth]; <br/> If (T> 255) <br/>{< br/> wdata [PT. X + Pt. y * m_biwidth] = lab_cout; <br/> -- hist [start_queue]; <br/> cue = true; <br/> continue; <br/>}< br/> T = wdata [PT. X + (pt. Y + 1) * m_biwidth]; <br/> If (T> 255) <br/>{< br/> wdata [PT. X + Pt. y * m_biwidth] = lab_cout; <br/> -- hist [start_queue]; <br/> cue = true; <br/> continue; <br/>}< br/> T = wdata [PT. X + 1 + (pt. Y + 1) * m_biwidth]; <br/> If (T> 255) <br/>{< br/> wdata [PT. X + Pt. y * m_biwidth] = lab_cout; <br/> -- hist [start_queue]; <br/> cue = true; <br/> continue; <br/>}< br/> T = wdata [PT. X + 1 + (pt. y) * m_biwidth]; <br/> If (T> 255) <br/>{< br/> wdata [PT. X + Pt. y * m_biwidth] = lab_cout; <br/> -- hist [start_queue]; <br/> cue = true; <br/> continue; <br/>}< br/> T = wdata [PT. X + 1 + (pt. y-1) * m_biwidth]; <br/> If (T> 255) <br/> {<br/> wdata [PT. X + Pt. y * m_biwidth] = lab_cout; <br/> -- hist [start_queue]; <br/> cue = true; <br/> continue; <br/>}< br/> statistic [start_queue]. push (PT); // not processed until it is pushed into the queue <br/>}</P> <p >}< br/> ++ lab_cout; <br/>}< br/> // the process of searching for the basin and raising the water // The background is not processed </P> <p> for (in T K = start_queue + 1; k <255; ++ K) <br/>{</P> <p> while (hist [k]> 0) <br/>{< br/> cue = true; <br/> while (cue) <br/>{< br/> cue = false; <br/> int size = (INT) statistic [K]. size (); <br/> for (INT I = 0; I <size; ++ I) <br/>{< br/> int T; <br/> int lab = 0; <br/> Pt = statistic [K]. front (); <br/> statistic [K]. pop (); <br/> // search for lab in 8 directions <br/> T = wdata [PT. X + (pt. y-1) * m_biwidth]; <br/> If (T> 255) <br/> {<br/> If (lab = 0) lab = T; <br/> Else if (T! = Lab) lab = wshed; </P> <p >}< br/> T = wdata [PT. x-1 + (pt. y-1) * m_biwidth]; <br/> If (T> 255) <br/> {<br/> If (lab = 0) lab = T; <br/> else if (T! = Lab) lab = wshed; </P> <p >}< br/> T = wdata [PT. x-1 + Pt. y * m_biwidth]; <br/> If (T> 255) <br/> {<br/> If (lab = 0) lab = T; <br/> else if (T! = Lab) lab = wshed; </P> <p >}< br/> T = wdata [PT. x-1 + (pt. Y + 1) * m_biwidth]; <br/> If (T> 255) <br/>{< br/> If (lab = 0) lab = T; <br/> else if (T! = Lab) lab = wshed; </P> <p >}< br/> T = wdata [PT. X + (pt. Y + 1) * m_biwidth]; <br/> If (T> 255) <br/>{< br/> If (lab = 0) lab = T; <br/> else if (T! = Lab) lab = wshed; </P> <p >}< br/> T = wdata [PT. X + 1 + (pt. Y + 1) * m_biwidth]; <br/> If (T> 255) <br/>{< br/> If (lab = 0) lab = T; <br/> else if (T! = Lab) lab = wshed; </P> <p >}< br/> T = wdata [PT. X + 1 + (pt. y) * m_biwidth]; <br/> If (T> 255) <br/> {<br/> If (lab = 0) lab = T; <br/> else if (T! = Lab) lab = wshed; </P> <p >}< br/> T = wdata [PT. X + 1 + (pt. y-1) * m_biwidth]; <br/> If (T> 255) <br/> {<br/> If (lab = 0) lab = T; <br/> else if (T! = Lab) lab = wshed; </P> <p >}< br/> If (lab = wshed) <br/>{< br/> -- hist [k]; <br/> wdata [PT. X + Pt. y * m_biwidth] = wshed; <br/> continue; <br/>}< br/> else <br/>{< br/> If (lab = 0) <br/> {<br/> statistic [K]. push (PT); // not processed until it is pushed into the queue <br/>}< br/> If (lab> 255) <br/> {<br/> wdata [PT. X + Pt. y * m_biwidth] = lab; <br/> -- hist [k]; <br/> cue = true; <br/>}</P> <p >}< br/> If (! Statistic [K]. size () <br/>{< br/> Pt = statistic [K]. front (); <br/> statistic [K]. pop (); <br/> wdata [PT. X + Pt. y * m_biwidth] = ++ lab_cout; <br/> -- hist [k]; <br/>}</P> <p> for (INT I = 0; I <m_biwidth; ++ I) <br/>{< br/> for (Int J = 0; j <m_height; ++ J) <br/> {<br/> If (wdata [I + J * m_biwidth]> 255) <br/> {<br/> data [I + J * m_biwidth] = 0; <br/>}< br/> else <br/> data [I + J * m_biwidth] = wdata [I + J * m_biwidth]; <br/>}< br/> Delete [] statistic; <br/> return true; <br/>}

This watershed has two parameters. The first is the distance data worth reaching, and the second is the histogram statistics of the distance data.

Histogram statistics determine the number of cycles. When a gray scale is processed, such as swallowed up by the water level or become a dam, it is considered to be processed.

Note 1: First, find the lowest gray level in the histogram statistics as the starting level.

NOTE 2: Find the Basin Based on the water level. The basin is now in the middle of 2. One independent, two connected to the relevant basin. We do not know whether to connect because it is based on histogram statistics. So there will be the following situation: first, we can find the initial Basin in eight directions, then we will drown this point as the basin, and set the record mark to true. The second 8 direction does not have a basin mark, so it may be the sphere of influence of the basin at the moment, but it has not developed yet, or it has nothing to do with the basin at the moment. When the process is completed, it indicates that the first basin with a low water level has been found. The rest is the other basin under this water level. Now we add a basin marker and start looking for the second basin. How to cycle until all the initial water level basins are found.

... (It's Too Late To write it again tomorrow)