[알고리즘]

harris corner detector code

Neo Park 2012. 10. 25. 11:30

 

CornerPoints DibHarrisCorner(CDib& dib, double th)
{
 register int i, j, x, y;
 
 int w = dib.GetWidth();
 int h = dib.GetHeight();
 
 BYTE** ptr = dib.GetPtr();
 
 //-------------------------------------------------------------------------
 // (fx)*(fx), (fx)*(fy), (fy)*(fy) 계산
 //-------------------------------------------------------------------------
 
 double** dx2 = new double*[h];
 double** dy2 = new double*[h];
 double** dxy = new double*[h];
 
 for( i = 0 ; i < h ; i++ ) {
      dx2[i] = new double[w];
      dy2[i] = new double[w];
      dxy[i] = new double[w];
      memset(dx2[i], 0, sizeof(int)*w);
      memset(dy2[i], 0, sizeof(int)*w);
      memset(dxy[i], 0, sizeof(int)*w);
 }
 
 double tx, ty;
 for( j = 1 ; j < h-1 ; j++ ) {
      for( i = 1 ; i < w-1 ; i++ )  {
           tx = (ptr[j-1][i+1] + ptr[j][i+1] + ptr[j+1][i+1] 
                   - ptr[j-1][i-1] - ptr[j][i-1] - ptr[j+1][i-1]) / 6.0;
           ty = (ptr[j+1][i-1] + ptr[j+1][i] + ptr[j+1][i+1] 
                   - ptr[j-1][i-1] - ptr[j-1][i] - ptr[j-1][i+1]) / 6.0;
 
       dx2[j][i] = tx * tx;
       dy2[j][i] = ty * ty;
       dxy[j][i] = tx * ty;
      }
 }
 
 //-------------------------------------------------------------------------
 // 가우시안 필터링
 //-------------------------------------------------------------------------
 
 double** gdx2 = new double*[h]; // 가우시언 필터링 수행후.
 double** gdy2 = new double*[h];
 double** gdxy = new double*[h];
 
 for( i = 0 ; i < h ; i++ ) {
      gdx2[i] = new double[w];
      gdy2[i] = new double[w];
      gdxy[i] = new double[w];


      memset(gdx2[i], 0, sizeof(double)*w);
      memset(gdy2[i], 0, sizeof(double)*w);
      memset(gdxy[i], 0, sizeof(double)*w);
 }
 
 double g[5][5] = { {1, 4, 6, 4, 1}, {4, 16, 24, 16, 24}, 
                          {6, 24, 36, 24, 6}, {4, 16, 24, 16, 24}, {1, 4, 6, 4, 1} };
 
 for( y = 0 ; y < 5 ; y++ ) {
      for( x = 0 ; x < 5 ; x++ )  {
            g[y][x] /= 256.;
      }
 }
 
 double tx2, ty2, txy;
 for( j = 2 ; j < h-2 ; j++ ) {
      for( i = 2 ; i < w-2 ; i++ )  {
           tx2 = ty2 = txy = 0;


           for( y = 0 ; y < 5 ; y++ )  {
                for( x = 0 ; x < 5 ; x++ )   {
                     tx2 += ( dx2[j+y-2][i+x-2] * g[y][x] );
                     ty2 += ( dy2[j+y-2][i+x-2] * g[y][x] );
                     txy += ( dxy[j+y-2][i+x-2] * g[y][x] );
                }

           }
           gdx2[j][i] = tx2;
           gdy2[j][i] = ty2;
           gdxy[j][i] = txy;
      }
 }
 
 //-------------------------------------------------------------------------
 // 코너 응답 함수 생성
 //-------------------------------------------------------------------------
 
 double** crf = new double*[h];
 for( i = 0 ; i < h ; i++ ) {
      crf[i] = new double[w];
      memset(crf[i], 0, sizeof(double)*w);
 }
 
 double k = 0.04;
 for( j = 2 ; j < h-2 ; j++ ) {
      for( i = 2 ; i < w-2 ; i++ )  {
           crf[j][i] = (gdx2[j][i]*gdy2[j][i] - gdxy[j][i]*gdxy[j][i])
                         - k*(gdx2[j][i] + gdy2[j][i])*(gdx2[j][i] + gdy2[j][i]);
      }
 }
 
 //-------------------------------------------------------------------------
 // 임계값보다 큰 국지적 최대값을 찾아 코너 포인트로 결정
 //-------------------------------------------------------------------------
 
 CornerPoints cp;
 cp.num = 0;
 
 for( j = 2 ; j < h-2 ; j++ ) {
      for( i = 2 ; i < w-2 ; i++ )  {
           if( crf[j][i] > th )   {
                if(  crf[j][i] > crf[j-1][i] && crf[j][i] > crf[j-1][i+1] &&
                     crf[j][i] > crf[j][i+1] && crf[j][i] > crf[j+1][i+1] &&
                     crf[j][i] > crf[j+1][i] && crf[j][i] > crf[j+1][i-1] &&
                     crf[j][i] > crf[j][i-1] && crf[j][i] > crf[j-1][i-1] )   {


                         if( cp.num < MAX_CORNER )    {
                              cp.x[cp.num] = i;
                              cp.y[cp.num] = j;
                              cp.num++;
                         }
                }
           }
      }
 }
 
 //-------------------------------------------------------------------------
 // 동적 할당한 메모리 해제
 //-------------------------------------------------------------------------
 
 for( i = 0 ; i < h ; i++ ) {
      delete [] dx2[i];
      delete [] dy2[i];
      delete [] dxy[i];
      delete [] gdx2[i];
      delete [] gdy2[i];
      delete [] gdxy[i];
      delete [] crf[i];
 }
 
 delete [] dx2;
 delete [] dy2;
 delete [] dxy;
 delete [] gdx2;
 delete [] gdy2;
 delete [] gdxy;
 delete [] crf;
 
 return cp;
}

 

 

 

참조 : http://blog.naver.com/PostView.nhn?blogId=adoda&logNo=10091847142