使用OpenCV和C++實現的分水嶺演算法(Watershed)
阿新 • • 發佈：2019-02-08
分水嶺演算法(watershed)是一種比較基本的數學形態學分割演算法，其基本思想是將灰度影象轉換為梯度影象，將梯度值看作高低起伏的山嶺，將區域性極小值及其鄰域看作一個“集水盆”。設想一個個“集水盆”中存在積水，且水位不斷升高，淹沒梯度較低的地方，當水漫過程停止後，影象就可以被分割成幾塊連通區域。
分水嶺演算法有不同的實現方法。本文要實現的是通過人為標註一些種子點，將這些種子點看作集水盆的底部，利用區域增長的方法，完成影象的分割。試圖實現OpenCV中cv::watershed函式的功能，經過測試，與OpenCV相比分割結果相似，但效能差很多。（前者32ms左右，後者8ms左右，原因可能是迴圈中使用了cv::mat來訪問影象中的元素，改用指標速度可能會提高很多）。
程式執行結果：
參考：
#include <opencv2/imgproc/imgproc.hpp> #include <opencv2/objdetect/objdetect.hpp> #include <opencv2/highgui/highgui.hpp> #include<vector> #include<iostream> #include<queue> #include<fstream> cv::Mat marker_mask; cv::Mat g_markers; cv::Mat img0, img, img_gray,wshed; cv::Point_<int> prev_pt(-1,-1); using std::vector; using std::queue; static void my_watershed(cv::Mat img,cv::Mat& markers,int comp_count); static void mouse_event(int event,int x, int y,int flags, void*) { if(img.rows==0) return; if(event==CV_EVENT_LBUTTONUP||!(flags&CV_EVENT_FLAG_LBUTTON)) prev_pt=cv::Point_<int>(-1,-1); else if(event==CV_EVENT_LBUTTONDOWN) prev_pt=cv::Point2i(x,y); else if(event==CV_EVENT_MOUSEMOVE&&(flags&CV_EVENT_FLAG_LBUTTON)) { cv::Point2i pt(x,y); if(prev_pt.x<0) prev_pt=pt; cv::line(marker_mask,prev_pt,pt,cv::Scalar(255,255,255),1,8,0); cv::line(img,prev_pt,pt,cv::Scalar(255,255,255),1,8,0); prev_pt=pt; cv::imshow("image",img); } } int main() { img0=cv::imread("Lenna.png",1); img=img0.clone(); CvRNG rng = cvRNG(-1); img_gray=img0.clone(); wshed=img0.clone(); marker_mask=cv::Mat(cv::Size(img0.cols,img0.rows),8,1); g_markers=cv::Mat(cv::Size(img0.cols,img0.rows),CV_32S,1); cv::cvtColor(img,marker_mask,CV_BGR2GRAY); cv::cvtColor(marker_mask,img_gray,CV_GRAY2BGR); for(int i=0;i<marker_mask.rows;i++) for(int j=0;j<marker_mask.cols;j++) marker_mask.at<unsigned char>(i,j)=0; for(int i=0;i<g_markers.rows;i++) for(int j=0;j<g_markers.cols;j++) g_markers.at<int>(i,j)=0; cv::imshow("image",img); cv::imshow("watershed transform",wshed); cv::setMouseCallback("image",mouse_event,0); for(;;) { int c=cv::waitKey(0); if((char)c==27) break; if((char)c=='r') { for(int i=0;i<marker_mask.rows;i++) for(int j=0;j<marker_mask.cols;j++) marker_mask.at<unsigned char>(i,j)=0; img0.copyTo(img); cv::imshow("image",img); } if((char)c=='w'||(char)c==' ') { vector<vector<cv::Point>> contours; CvMat* color_tab=0; int comp_count=0; cv::findContours(marker_mask,contours,CV_RETR_CCOMP,CV_CHAIN_APPROX_SIMPLE,cv::Point(0,0)); for(int i=0;i<g_markers.rows;i++) for(int j=0;j<g_markers.cols;j++) g_markers.at<int>(i,j)=0; vector<vector<cv::Point> >::iterator iter=contours.begin(); for(int i=0;i<(int)contours.size();i++) { cv::drawContours(g_markers,contours,i,cv::Scalar::all(comp_count+1), 1,8,vector<cv::Vec4i>()); comp_count++; } if(comp_count==0) continue; color_tab=cvCreateMat(1,comp_count,CV_8UC3); for(int i=0;i<comp_count;i++) { uchar* ptr=color_tab->data.ptr+i*3; ptr[0]=(uchar)(cvRandInt(&rng)%180+50); ptr[1]=(uchar)(cvRandInt(&rng)%180+50); ptr[2]=(uchar)(cvRandInt(&rng)%180+50); } cv::Mat temp=g_markers.clone(); double t=(double)cvGetTickCount(); //my_watershed(img0,g_markers,comp_count); cv::watershed(img0,g_markers); t=(double)cvGetTickCount()-t; std::cout<<"exec time= "<<t/(cvGetTickFrequency()*1000.)<<std::endl; for(int i=0;i<g_markers.rows;i++) for(int j=0;j<g_markers.cols;j++) { int idx=g_markers.at<int>(i,j); uchar* dst=&wshed.at<uchar>(i,j*3); if(idx==-1) dst[0]=dst[1]=dst[2]=(uchar)255; else if(idx<=0||idx>comp_count) dst[0]=dst[1]=dst[2]=(uchar)8; else{ uchar* ptr=color_tab->data.ptr+(idx-1)*3; dst[0]=ptr[0];dst[1]=ptr[1];dst[2]=ptr[2]; } } cv::addWeighted(wshed,0.5,img_gray,0.5,0,wshed); cv::imshow("watershed transform",wshed); cvReleaseMat(&color_tab); } } return 0; } static void my_watershed(cv::Mat img0,cv::Mat& markers,int comp_count) { cv::Mat gray=cv::Mat(cv::Size(img0.rows,img0.cols),8,1); cv::cvtColor(img0,gray,CV_BGR2GRAY); cv::Mat imge=cv::Mat(cv::Size(img0.rows,img0.cols),8,1); cv::Sobel(gray,imge,CV_8U,1,1); vector<queue<cv::Point2i>*>Labeleddata;//影象中各連通區域的點 queue<cv::Point2i>* pque;//某連通區域已包含的點 queue<cv::Point2i> quetem; //用於提取某連通區域中輸入種子點中的初始種子點 vector<int*> SeedCounts; int* Array; cv:: Point2i temp; int row=imge.rows,col=imge.cols; cv::Mat marker_saved=markers.clone(); bool up,down,right,left,uplef,uprig,downlef,downrig; int m,n; for(int i=0;i<comp_count;i++) { Array=new int[256]; SeedCounts.push_back(Array);//統計某waterlevel的各個連通區域中種子點的個數 pque=new queue<cv::Point2i>[256]; Labeleddata.push_back(pque);//儲存該連通區域中種子生長所得的點 } for(int i=0;i<row;i++) for(int j=0;j<col;j++) { if(markers.at<int>(i,j)>0) { temp.x=i; temp.y=j; quetem.push(temp); int num=markers.at<int>(i,j); markers.at<int>(i,j)=-1;//該點已處理,其他種子點生長時將繞過該點 while(!quetem.empty()) { up=down=right=left=uplef=uprig=downlef=downrig=false; temp=quetem.front(); //提取出一個點，在該點的八連通區域內尋找可生長點 m=temp.x; n=temp.y; quetem.pop(); if(m-1>=0)//若上方可生長則新增為新種子 { if(markers.at<int>(m-1,n)==num) { temp.x=m-1; temp.y=n; quetem.push(temp); markers.at<int>(m-1,n)=-1; } else{ up=true; } } if(m-1>=0&&n-1>=0) { if(markers.at<int>(m-1,n-1)==num) { temp.x=m-1; temp.y=n-1; quetem.push(temp); markers.at<int>(m-1,n-1)=-1; } else{ uplef=true; } } if(m+1<=row-1) { if(markers.at<int>(m+1,n)==num) { temp.x=m+1; temp.y=n; quetem.push(temp); markers.at<int>(m+1,n)=-1; } else{ down=true; } } if(m+1<=row-1&&n+1<=col-1) { if(markers.at<int>(m+1,n+1)==num) { temp.x=m+1; temp.y=n+1; quetem.push(temp); markers.at<int>(m+1,n+1)=-1; } else{ downrig=true; } } if(n+1<=col-1) { if(markers.at<int>(m,n+1)==num) { temp.x=m; temp.y=n+1; quetem.push(temp); markers.at<int>(m,n+1)=-1; } else{ right=true; } } if(m-1>=0&&n+1<=col-1) { if(markers.at<int>(m-1,n+1)==num) { temp.x=m-1; temp.y=n+1; quetem.push(temp); markers.at<int>(m-1,n+1)=-1; } else{ uprig=true; } } if(n-1>=0) { if(markers.at<int>(m,n-1)==num) { temp.x=m; temp.y=n-1; quetem.push(temp); markers.at<int>(m,n-1)=-1; } else{ left=true; } } if(m+1<=row-1&&n-1>=0) { if(markers.at<int>(m+1,n-1)==num) { temp.x=m+1; temp.y=n-1; quetem.push(temp); markers.at<int>(m+1,n-1)=-1; } else{ downlef=true; } } //八連通區域中有未標記點，則該點屬於初始種子點 if(up||down||right||left||uplef||downlef||uprig||downrig) { temp.x=m; temp.y=n; Labeleddata[comp_count-1][imge.at<uchar>(m,n)].push(temp); SeedCounts[comp_count-1][imge.at<uchar>(m,n)]++; } } } } bool active; int waterlevel; for(waterlevel=0;waterlevel<180;waterlevel++) { active=true; while(active) //當1-count_com個連通區域都無可生長點時結束迴圈 { active=false; for(int i=0;i<comp_count;i++)//將區域i中將waterlevel梯度以下的點用於區域增長 { if(!Labeleddata[i][waterlevel].empty())//區域增長，經過多次迭代，直至該區域，該waterlevel無可生長點。 { active=true; while(SeedCounts[i][waterlevel]>0) //SeedCount中保留了前一輪生長後各區域，各waterlevel的種子點個數，本輪生長結束後，將根據Labeleddata中的元素個數更新 { SeedCounts[i][waterlevel]--; temp=Labeleddata[i][waterlevel].front(); Labeleddata[i][waterlevel].pop(); m=temp.x; n=temp.y; int num=marker_saved.at<int>(m,n); if(m-1>=0) { if(!marker_saved.at<int>(m-1,n))//上方點未處理過 { temp.x=m-1; temp.y=n; marker_saved.at<int>(m-1,n)=num; if(imge.at<uchar>(m-1,n)<=waterlevel) Labeleddata[i][waterlevel].push(temp); else{ Labeleddata[i][imge.at<uchar>(m-1,n)].push(temp); //本次生長不處理，可能在waterlevel變化到某值時再用於生長 SeedCounts[i][imge.at<uchar>(m-1,n)]++; } } } if(m+1<=row-1) { if(!marker_saved.at<int>(m+1,n)) { temp.x=m+1; temp.y=n; marker_saved.at<int>(m+1,n)=num; if(imge.at<uchar>(m+1,n)<=waterlevel) Labeleddata[i][waterlevel].push(temp); else{ Labeleddata[i][imge.at<uchar>(m+1,n)].push(temp); SeedCounts[i][imge.at<uchar>(m+1,n)]++; } } } if(n+1<=col-1) { if(!marker_saved.at<int>(m,n+1)) { temp.x=m; temp.y=n+1; marker_saved.at<int>(m,n+1)=num; if(imge.at<uchar>(m,n+1)<=waterlevel) Labeleddata[i][waterlevel].push(temp); else{ Labeleddata[i][imge.at<uchar>(m,n+1)].push(temp); SeedCounts[i][imge.at<uchar>(m,n+1)]++; } } } if(n-1>=0) { if(!marker_saved.at<int>(m,n-1)) { temp.x=m; temp.y=n-1; marker_saved.at<int>(m,n-1)=num; if(imge.at<uchar>(m,n-1)<=waterlevel) Labeleddata[i][waterlevel].push(temp); else { Labeleddata[i][imge.at<uchar>(m,n-1)].push(temp); SeedCounts[i][imge.at<uchar>(m,n-1)]++; } } } } SeedCounts[i][waterlevel]=Labeleddata[i][waterlevel].size(); } } } } markers=marker_saved.clone(); }