影像處理之Connected Component Labeling(CCL) algorithm in C++

在二值化的影像(Binary Image)上，Connected Component Labeling(CCL)是一種想要將一張影像上的物體，切割成多個單一個物件的方法。
假設有一張圖如下(黑背景白前景):

因為影像資料只有0~255(8-bits per pixel)，
會無法得知Pixel之間的連結關係，
所以想要分析這三個物體(不論是物體的寬或高或者是面積)，
須要定義如何取得這三個物體的Pixel連結。

一般來說，Pixel之間的連結有兩種類型，分別為4與8連通，
4連通:

	N2
N3	Pixel	N1
	N4

8連通:

N4	N3	N2
N5	Pixel	N1
N6	N7	N8

例如有一張二值影像(黑底白前景):
則經過4與8連通運算之後結果如下:
4連通:

0	0	0	0	1	0	0	0	0	0
0	0	0	0	0	2	2	2	2	2
0	3	0	0	0	0	0	0	0	2
0	3	0	0	0	0	0	0	2	2
0	3	3	3	0	2	2	2	2	0
0	3	0	3	3	0	0	0	0	0
0	3	0	3	3	0	0	0	0	0
0	3	3	3	3	0	0	0	0	0

8連通:

0	0	0	0	1	0	0	0	0	0
0	0	0	0	0	1	1	1	1	1
0	1	0	0	0	0	0	0	0	1
0	1	0	0	0	0	0	0	1	1
0	1	1	1	0	1	1	1	1	0
0	1	0	1	1	0	0	0	0	0
0	1	0	1	1	0	0	0	0	0
0	1	1	1	1	0	0	0	0	0

那要怎麼取得影像上4連通以及8連通的關係，目前有很多方法可以實作[1]，
我覺得比較容易理解的方法，是先找到起始位置(seed)後，
依照4/8連通的需求，使用DFS(depth-first searching)的演算法[2]。

DFS大致上的pseudo code如下:
輸入: 在一張二值化的影像(I)以及一個起始點(seed，這個起始點必須為前景)
參數: 4/8連通(N)

第一種實作方法－使用遞迴
1. for each N in I(seed) as I(x)
2.   if  I(x) is Foreground
        then Set I(x) as traversed pixel (or labeling) and x as seed, recursively do 1
第二種實作方法－使用Stack，用來暫存開拓後的位置
暫存: Stack(S)
1. PUSH(S, seed), Set I(seed) as traversed pixel
2. while ( not IS_EMPTY(S) )
2.1.    curSeed = POP(S)
2.2.    for each N in I(curSeed) as I(x)
2.3.        if I(x) is Foreground
                   then Set I(x) as traversed pixel (or labeling) and PUSH(S, x)
2.4 end while

然而，我在觀察二值影像資料之後發現，其實可以先將二值影像表示成以下形式[3]:

				1
					2	2	2	2	2
	3								4
	5							6	6
	7	7	7		8	8	8	8
	9		10	10
	11		12	12
	13	13	13	13

把連續的Pixel定義成一個Bound，並且記錄Bound的位置資訊(row & left & right)。
就會變成以下形式:

Bound	row	left	right
1	0	4	5
2	1	5	10
3	2	1	2
4	2	9	10
5	3	1	2
6	3	8	10
7	4	1	4
8	4	5	9
9	5	1	2
10	5	3	5
11	6	1	2
12	6	3	5
13	7	1	5

接下來就是再透過DFS，將每個Bound依照4/8連通的規則連結。
由於Bound的表示法已經將紀錄每個row的連續Pixel位置，所以進行連通運算的時候，
不需要比較同一個row的Bound。

範例如下:

	7	7	7		8	8	8	8
	9		10	10

4連通規則:

If ( Bound(7).right > Bound(9).left && Bound(7).left < Bound(9).right )

 Then Bound(7) & Bound(9) is connected

else
 7 & 9 not connected

8連通規則:

If ( Bound(8).left <= Bound(10).right && (Bound(8).right + 1) >= Bound(10).left )

 Then Bound(8) & Bound(10) is connected

else
  8 & 10 not connected

完整的程式碼如下:

#include 

#include 
#include 

#define PATH "C:\\OpenCV3.2\\build\\x64\\vc11\\lib\\"
#define VER "320"

#ifdef _DEBUG
#pragma comment( lib , PATH"opencv_core"VER"d.lib") 
#pragma comment( lib , PATH"opencv_highgui"VER"d.lib") 
#pragma comment( lib , PATH"opencv_imgcodecs"VER"d.lib") 
#pragma comment( lib , PATH"opencv_imgproc"VER"d.lib") 
#else
#pragma comment( lib , PATH"opencv_core"VER".lib") 
#pragma comment( lib , PATH"opencv_highgui"VER".lib") 
#pragma comment( lib , PATH"opencv_imgcodecs"VER".lib") 
#pragma comment( lib , PATH"opencv_imgproc"VER".lib") 
#endif

using namespace cv ;

///recording the pixel boundary
typedef struct BoundStructTAG {
 //0 based 
 int left , right  ;
 int _8Left , _8Right ;
 int row ;
 //1 based
 int indice;
 bool l_used , r_used ;

 int insertIdx ;

 int objectID ;
 BoundStructTAG(int &l , int &r , int &i , int &rMax) {
  left = l ;
  right = r ;
  indice = i ;
  l_used = r_used = false ;
  _8Left = MAX(0 , l - 1) ;
  _8Right = r ;
 } ;

 inline int isIntersect(BoundStructTAG &n) {
  if ( n.right < this->left ) {
   return -1 ;
  } else  {
   if ( n.left <= this->_8Right ) {
    return 0 ;
   } else {
    return 1 ;
   }
  }
 } ;

 inline int isIntersect4(BoundStructTAG &n) {
  if ( n.right <= this->left ) {
   return -1 ;
  } else {
   if ( n.left < this->right ) {
    return 0 ;
   } else {
    return 1 ;
   }
  }
 } ;

 inline bool isContain(int &x) {
  if ( x >= this->left && x < this->right ) {
   return true ;
  }
  return false ;
 } ;
} BoundStruct;

typedef struct BinaryImageStruct {
 int iWidth , iHeight ;

 BoundStruct *boundaries ;
 int boundaryCount ;
 int *rowIndice ; //indexing the boundary(row)

 BinaryImageStruct() {
  iWidth = iHeight = 0 ;
  boundaryCount = 0 ;
  boundaries = NULL ;
  rowIndice = NULL ;
 } ;

} BinaryImage ;

///
//input :
//   img -> Run-length based binary image
//output :
//   stats -> object information. each row represent an object(CV_32S), the information ordered by x1, y1, width, height, area, number of run-length of an object
//   centroid -> object's centroid. each row represent an object(CV_64F)
//   objectLinkedList -> object's run-length linked list index
//   objStartLinkedList -> object's run-length start index
//parameter :
//   conn -> conn = 4 or 8
//return : ammount of detected objects
inline int LabelingImage(BinaryImage* img , cv::Mat &stats , cv::Mat &centroid , int* &objLinkedList , int* &objStartLinkedList, int conn) {
 uchar *boundUsed = (uchar *)std::malloc(sizeof(uchar) * img->boundaryCount) ;
 int bp = 0 ;
 //stack
 int *openList = (int *)std::malloc(sizeof(int) * img->boundaryCount) ;
 int sp = 0 ;

#define IS_USED(x) (boundUsed[x])
#define PUSH(x) { openList[sp++] = x; }
#define POP() ( openList[--sp] )
#define CLEAR() { sp = 0 ; }
#define IS_EMPTY() ( sp == 0 ) 
#define SET_USED(x) { boundUsed[x] = 0xff ;}
 std::memset( boundUsed , 0 , sizeof(uchar) * img->boundaryCount ) ;
 int *labels = (int *)std::malloc( (img->boundaryCount) * sizeof(int) ) ;
 
 //find connectivity using DFS
 int objs = 0 ;
 for ( int i = 1 ; i < img->boundaryCount ; ++i ) {
  if ( !IS_USED(i) ) {
   CLEAR() ;
   PUSH(i) ;
   objs++ ;
   
   while ( !IS_EMPTY() ) {
    int curIdx = POP() ;
    SET_USED(curIdx) ;
    labels[curIdx] = objs ;
    BoundStruct *curBound = img->boundaries + curIdx ;
    int row = curBound->row ;
    uchar breakDown = 0 ;
    int nxtRow = row + 1 ;
    int prevRow = row - 1 ;
    if ( prevRow >= 0 ) {
     int rc = img->rowIndice[prevRow + 1] - img->rowIndice[prevRow] ;
     BoundStruct* prevB = img->boundaries + img->rowIndice[prevRow] ;
     switch ( conn ) {
     case 4:
      for ( int j = 0 ; j < rc && !breakDown ; ++j ) {
       BoundStruct* bound = prevB + j ;
       if ( !IS_USED(bound->indice) ) {
        int inter = curBound->isIntersect4(*bound);
        switch (inter)
        {
        case 0:
         PUSH(bound->indice) ;
         break ;
        case 1:
         breakDown = 0xff ;
         break ;
        default:
         break;
        }
       }
      }
      break ;
     case 8:
      for ( int j = 0 ; j < rc && !breakDown ; ++j ) {
       BoundStruct* bound = prevB + j ;
       if ( !IS_USED(bound->indice) ) {
        int inter = curBound->isIntersect(*bound);

        switch (inter)
        {
        case 0:
         PUSH(bound->indice) ;
         break ;
        case 1:
         breakDown = 0xff ;
         break ;
        default:
         break;
        }
       }
      }
      break ;
     } //end of switch(conn)
    } // end of prev searching
    if ( nxtRow < img->iHeight ) {
     breakDown = 0 ;
     int rc = img->rowIndice[nxtRow + 1] - img->rowIndice[nxtRow] ;
     BoundStruct* nextB = img->boundaries + img->rowIndice[nxtRow] ;
     switch ( conn ) {
     case 4:
      for ( int j = 0 ; j < rc && !breakDown ; ++j ) {
       BoundStruct* bound = nextB + j ;
       if ( !IS_USED(bound->indice) ) {
        int inter = curBound->isIntersect4(*bound);
        switch (inter)
        {
        case 0:
         PUSH(bound->indice) ;
         break ;
        case 1:
         breakDown = 0xff ;
         break ;
        default:
         break;
        }
       }
      }
      break ;
     case 8:
      for ( int j = 0 ; j < rc && !breakDown ; ++j ) {
       BoundStruct* bound = nextB + j ;
       if ( !IS_USED(bound->indice) ) {
        int inter = curBound->isIntersect(*bound);

        switch (inter)
        {
        case 0:
         PUSH(bound->indice) ;
         break ;
        case 1:
         breakDown = 0xff ;
         break ;
        default:
         break;
        }
       }
      }
      break ;
     } //end of switch(conn)
     
    }
   } ;
  }
 } // end of for
 //-----------------------------------------------------------------
 //collect object info
 stats = cv::Mat(objs , 6 , CV_32SC1) ;
 centroid = cv::Mat(objs , 2 , CV_64FC1) ;

 int *objss = (int *)stats.data;
 double *tmp2Ptr = (double *)centroid.data;


 objLinkedList = (int *)std::malloc( sizeof(int) * img->boundaryCount ) ;
 objStartLinkedList = (int *)std::malloc( sizeof(int) * objs ) ;

 int *objLastIndice = (int *)std::malloc( sizeof(int) * objs ) ;

 std::memset(objLastIndice , 0 , sizeof(int) * objs) ;
 objLinkedList[0] = 0 ;

 for ( int i = 0 ; i < objs ; ++i , objss += 6 , tmp2Ptr += 2) {
  objss[0] = INT_MAX ;  //x1
  objss[1] = INT_MAX ;  //y1
  objss[2] = 0 ;  //x2
  objss[3] = 0 ;  //y2
  objss[4] = 0 ;  //area
  objss[5] = 0 ;  //bound counts
  tmp2Ptr[0] = 0 ;
  tmp2Ptr[1] = 0 ;
 }



 for ( int i = 0 ; i < img->iHeight ; ++i ) {

  const int curRowSize = img->rowIndice[i + 1] - img->rowIndice[i] ;
  if ( curRowSize == 0 ) continue ;
  BoundStruct *curRow = img->boundaries + img->rowIndice[i] ;
  
  for ( int j = 0 ; j < curRowSize ; ++j ) {
   int left = curRow[j].left , right = curRow[j].right ;
   int indice = labels[curRow[j].indice];
   indice = indice - 1;


   int lastIdx = objLastIndice[indice];
   objLastIndice[indice] = curRow[j].indice ;
   if ( lastIdx != 0 ) {
    objLinkedList[lastIdx] = curRow[j].indice;
   } else {
    objStartLinkedList[indice] = curRow[j].indice;
   }
   curRow[j].objectID = indice ;
   objLinkedList[curRow[j].indice] = 0 ;
   left = curRow[j].left ;
   int numOfPixels = (right - left);

   objss = (int *)(stats.data + (indice * stats.step[0])) ;
   objss[0] = MIN(objss[0] , left) ;
   objss[1] = MIN(objss[1] , i ) ;
   objss[2] = MAX(objss[2] , right) ;
   objss[3] = MAX(objss[3] , i ) ;

   objss[4] = objss[4] + numOfPixels ;
   objss[5] = objss[5] + 1 ;

   tmp2Ptr = (double *)(centroid.data + (indice * centroid.step[0])) ;
   tmp2Ptr[0] += ((right + left - 1) * numOfPixels / 2) ;
   tmp2Ptr[1] += ((right - left) * i) ;
  }
 }
 
 int *statPtr = (int *)stats.data ;
 double *cPtr = (double *)centroid.data ;
 for ( int i = 0 ; i < objs ; ++i , statPtr += 6 , cPtr += 2) {
  statPtr[2] = MIN( img->iWidth , statPtr[2] - statPtr[0] );  //transfer to width
  statPtr[3] = MIN( img->iHeight , statPtr[3] - statPtr[1] + 1 ); //transfer to height

  cPtr[0] = cPtr[0] / (double)statPtr[4] ;
  cPtr[1] = cPtr[1] / (double)statPtr[4] ;
 }

#undef PUSH
#undef POP
#undef CLEAR
#undef IS_EMPTY
#undef IS_USED
#undef SET_USED
 std::free(objLastIndice) ;
 std::free(boundUsed) ;
 std::free(labels) ;
 std::free(openList) ;
 return objs ;
} ;


void main() {
 cv::Mat img = cv::imread("GrayImage.bmp" , cv::IMREAD_GRAYSCALE) ;

 BinaryImage *bin = (BinaryImage *)std::malloc(sizeof(BinaryImage));
 int w = (img.cols / 2) + 1 ;
 int h = img.rows ;
 bin->boundaries = (BoundStruct *)std::malloc((sizeof(BoundStruct) * (w * h) + 1) * sizeof(BoundStruct)) ;
 bin->iHeight = h ;
 bin->iWidth = img.cols ;
 bin->rowIndice = (int *)std::malloc(sizeof(int) * ( h + 1 ) ) ;
 bin->rowIndice[0] = 1 ;

 //encode to RLE(run-length encoding)
 int indice = 1 ;
 for ( int i = 0 ; i < img.rows ; ++i ) {
  uchar *ptr = img.data + (i * img.step[0]) ;
  for ( int j = 0 ; j < img.cols ; ++j ) {
   if ( ptr[j] ) {
    int l = j ;
    while ( ptr[j] && j < img.cols) {
     ++j ;
    }

    bin->boundaries[indice] = BoundStruct(l , j , indice , img.cols) ;
    bin->boundaries[indice].row = i ;
    indice++;
   }
  }
  bin->rowIndice[i + 1] = indice ;
 }
 bin->boundaryCount = indice ;
 //----------------------------------------------
 //finding connectivity
 cv::Mat stat_ , centroid ;
 int *objStart , *objLinkedList ;
 int objects = LabelingImage(bin , stat_ , centroid , objLinkedList , objStart , 8 ) ;
 //----------------------------------------------

 //output blob result
 int *objInfo = (int *)stat_.data ;
 double *centroidPtr = (double *)centroid.data ;

 for ( int i = 0 ; i < objects ; ++i , objInfo += 6 , centroidPtr += 2) {
  printf( "object id : %d, BoundingBox : [%d,%d,%d,%d], Area : %d, Centroid=(%f,%f)\r\n" , 
   i , objInfo[0] , objInfo[1] , objInfo[2] , objInfo[3] , objInfo[4] , 
   centroidPtr[0] , centroidPtr[1] ) ;
 }
 //-------------------------------------------------

 std::free(objStart) ;
 std::free(objLinkedList) ;
 std::free(bin->rowIndice) ;
 std::free(bin->boundaries) ;
 std::free(bin) ;
}

[1] : https://en.wikipedia.org/wiki/Connected-component_labeling
[2] : https://stackoverflow.com/questions/14465297/connected-component-labelling
[3] : https://en.wikipedia.org/wiki/Run-length_encoding