[Happy DSA] 图像的contour trace算法C++实现

问题提出：

给定一个二维图像，基于某个threshold，来提取contours。

在图形图像学中，这个问题有比较好的解决方案，google "coutour trace"，可以得到以下2个比较好的参考文献：

1. http://en.wikipedia.org/wiki/Moore_neighborhood

2. http://www.imageprocessingplace.com/downloads_V3/root_downloads/tutorials/contour_tracing_Abeer_George_Ghuneim/index.html

算法描述：

The following is a formal description of the Moore-Neighbor tracing algorithm:

Input: A square tessellation, T, containing a connected component P of black cells.

Output: A sequence B (b₁, b₂ ,..., b_k) of boundary pixels i.e. the contour.

Define M(a) to be the Moore neighborhood of pixel a. 
Let p denote the current boundary pixel. 
Let c denote the current pixel under consideration i.e. c is in M(p).

Begin

• Set B to be empty.
• From bottom to top and left to right scan the cells of T until a black pixel, s, of P is found.
• Insert s in B.
• Set the current boundary point p to s i.e. p=s
• Backtrack i.e. move to the pixel from which s was entered.
• Set c to be the next clockwise pixel in M(p).
• While c not equal to s do
     If c is black
  • insert c in B
  • set p=c
  • backtrack (move the current pixel c to the pixel from which p was entered)
     else
  • advance the current pixel c to the next clockwise pixel in M(p)
  end While

End 

需要注意的是，

1. 按照先列后行，每列从下到上，每行从左到右的顺序（当然也可以按照其它顺序，只要整个过程一致就可以了），搜索第一个BLACK PIXEL；

2. 有2种情况可以作为trace one contour的条件：

a) 搜索到第一个BLACK PIXEL时，如果它位于某列的第一行，将它的backtrace point定义该列的-1行；

b) 否则进入这个BLACK PIXEL的backtrace point，一定要是WHITE PIXEL。

以下代码段是基于Moore Neighborhood算法的C++实现：

void TraceOneContour(Image<char>* bp, const std::pair<int,int>& startpos,const std::pair<int,int>& btpos, std::vector<int>& contour){int xSize = bp->x_size, ySize = bp->y_size;char BLACK = 1;enum { LKN = 8 };static int cwmat[LKN][2] = {// clock-wise matrix{-1, +1},// left-upper{0, +1},// upper{+1,+1},{+1, 0},{+1, -1},{0, -1},{-1, -1},{-1, 0}// left};contour.clear();int px = startpos.first, py = startpos.second;// current center pointint bx = btpos.first, by = btpos.second ;// bt pointcontour.push_back(py*xSize + px);int cx = bx, cy = by;bool bstop = false;int nVisitStart = 1;static int nStopTimes = 5;do {// get next clockwise pixel for (cx, cy) around (px, py)int cur = -1;for (int i = 0; i < LKN; ++i){if (cwmat[i][0] == (cx - px) &&cwmat[i][1] == (cy - py)){cur = i;break;}}if (cur < 0) break;int last = cur;cur = (cur+1) % LKN;while (cur != last){// loop 8-connected cells around pint prev = (cur - 1 + LKN)%LKN;cx = px + cwmat[cur][0], cy = py + cwmat[cur][1];bx = px + cwmat[prev][0], by = py + cwmat[prev][1];if ((cx == startpos.first && cy == startpos.second)&& (bx == btpos.first && by == btpos.second)){// Jacob's stopping criterionbstop = true;break;}if (cx >= 0 && cx < xSize && cy >=0 && cy < ySize&& bp->pixRef(cx, cy) == BLACK){if (cx == startpos.first && cy == startpos.second){  // In some class, only Jacob's stopping cannot stop algorithm// so add visit time stopping criterionif (++nVisitStart >= nStopTimes){bstop = true;break;}}contour.push_back(cy * xSize + cx);px = cx, py = cy;cx = bx, cy = by;break;}else{cur = (cur+1) % LKN;}}// handle isolated pixel pointif (cur == last) break;if (bstop) break;}  while (1);}template <typename Type, typename SpecFunc>void TraceImageContours(const Image<Type>* in_img, const SpecFunc& specFunc,std::vector<std::vector<int> >& vPolyOut){assert(in_img && in_img->x_size > 0 && in_img->y_size > 0);vPolyOut.clear();int xSize = in_img->x_size, ySize = in_img->y_size;int npix = xSize * ySize;char WHITE = 0, BLACK = 1, GREEN = 2;Image<char> bp(xSize, ySize);std::transform(in_img->pix, in_img->pix + npix, bp.pix, specFunc);std::vector<int> onecontour;std::pair<int,int> bt(0, -1);// firstly column from bottom to top, secondly row from left to rightfor (int i = 0; i < xSize; i++){for (int j = 0; j < ySize; j++){if (bp.pixRef(i, j) == WHITE){bt = std::make_pair(i, j);}else if (bp.pixRef(i, j) == BLACK){if (j == 0) bt = std::make_pair(i, -1);if ((j > 0 && bp.pixRef(i, j-1) == WHITE)|| j == 0){ // start to trace one contourTraceOneContour(&bp, std::make_pair(i, j), bt, onecontour); // mark contour point green(visited)for (BOOST_AUTO(itc, onecontour.begin());itc != onecontour.end(); ++itc){bp.pix[*itc] = GREEN;}if (!onecontour.empty()) vPolyOut.push_back(onecontour);}}}}}

关于Image<T>，可以参考下面的C++类：

template <typename T>struct Image {        int x_size;    int y_size;    T* pix;public:    Image(int w, int h) : x_size(w), y_size(h) {        pix = new T[w*h];    }    Image(int w, int h, T* p) : x_size(w), y_size(h), pix(p) {    }    ~Image() {        delete [] pix;    }         T* detach() {      T* result = pix;      pix = 0;      return result;    }        T& pixRef(int x, int y) {      return pix[y * x_size + x];    }    const T& pixRef(int x, int y) const {      return pix[y * x_size + x];    }};