RHT_VC6.rar_RHT_RHT_VC6_Randomized Hough_随机霍夫变换_霍夫变换

/*-----------本程序仅demo-------------------*/ /*---------------如有参阅，欢迎----*/ #include <cv.h> #include <highgui.h> #include <math.h> #include<stdlib.h> #include<time.h> typedef struct tagD { struct tagD* D; CvPoint v; unsigned char data; }D;//定义边缘数据的结构 typedef struct tagCell { struct tagCell* cell; int a; int b; int r; int score; }cell; //定义参数结构 typedef struct tagCircle { int a; int b; int r; }Circle; //定义参数结构 typedef struct tagP { struct tagCell* pCell; struct tagCell* cCell; struct tagCell* eCell; }P;//链表头、尾 typedef struct tagDaa { const char* filename; unsigned int Max_Lim; unsigned int Mim_Lim; unsigned int T1; unsigned int Block1; unsigned int Hlock1; int score; double A; double K; double SF; }Daa; //改进的RHT检测圆函数 int RHT_Y(Daa* G_c); int main(int argc, char** argv) { const char* filename = "cell.BMP"; IplImage* img=cvLoadImage(filename); IplImage* gray = cvCreateImage( cvGetSize(img), 8, 1 ); CvMemStorage* storage = cvCreateMemStorage(0); cvCvtColor( img, gray, CV_BGR2GRAY ); cvSmooth( gray, gray, CV_GAUSSIAN, 9, 9 ); // smooth it, otherwise a lot of false circles may be detected CvSeq* circles = cvHoughCircles( gray, storage, CV_HOUGH_GRADIENT, 1, gray->height/8, 50, 70 ); int i; for( i = 0;i < circles->total; i++ ) { float* p = (float*)cvGetSeqElem( circles, i ); cvCircle( img, cvPoint(cvRound(p[0]),cvRound(p[1])), cvRound(p[2]), CV_RGB(255,0,0), 3, 8, 0 ); } cvNamedWindow( "circl", 1 ); cvShowImage( "circl", img ); Daa* G_c = new Daa; G_c->filename = filename; //图片路径 //切割参数>=1; //当切割参数为1，1时为标准的随机霍夫变换 //当切割参数不为1，1时为改进的随机霍夫变换 G_c->Block1 = 5; //对图象横向切割,仅RHT_Y使用 G_c->Hlock1 = 5; //对图象纵线切割,仅RHT_Y G_c->Max_Lim = 200; //限制检测的最大半径 G_c->Mim_Lim = 10; //限制检测的最小圆半径 G_c->T1 = 500; //越大越好，但是速度越慢, G_c->K = 2; //确定真实圆的计算值，一般取1.5~2，过小易产生过多无效圆 G_c->SF = 0.8; //判断真实圆的累加条件,0.5~0.8越小真实圆判定越准确； G_c->A = 3; //比较假象圆r,a,b的阕值 G_c->score = 3; //假象圆判断条件 //本程序会比较慢，请耐心等待结果； //由于程序仅是Demo版本，改进的RHT注意调节合适的分割范围； RHT_Y(G_c); return 0; } //改进的RHT检测圆函数 int RHT_Y(Daa* G_c) { Daa* GG = G_c; IplImage* img = cvLoadImage( GG->filename, 0); IplImage* im = cvLoadImage( GG->filename, 1); if( img == 0) return -1; //读取图象文件,设置为灰度图象 IplImage* iml = NULL; iml = cvCreateImage(cvGetSize(img), IPL_DEPTH_8U, 1); cvSmooth( img, iml, CV_BLUR, 7,7); IplImage* dst = NULL; dst = cvCreateImage(cvGetSize(img), IPL_DEPTH_8U, 1); cvCanny( iml, dst, MAX(50/2,1), 50, 3 ); // cvCanny(img, dst, 50, 150, 3); cvNamedWindow( "dst", 1 ); cvShowImage( "dst", dst ); Circle RealC[1000]; int realC = 0; CvMat stub, *dimg = (CvMat*)dst; dimg = cvGetMat( dimg, &stub); const uchar* image; int step, width, height; image = dimg->data.ptr; step = dimg->step; width = dimg->cols; height = dimg->rows; D* PRData = new D; unsigned int Hlock = GG->Block1; unsigned int Block = GG->Hlock1; srand((int)time(0)); uchar ImgData[500*500]; { for(int ai = 0; ai < height; ai++ ) { for(int aj = 0; aj < width ; aj++ ) { ImgData[ai * step + aj] = image[ai * step + aj]; } } } for (unsigned int hr = 0; hr < Hlock; hr++) { for (unsigned int br = 0; br < Block; br++) { unsigned int Size = 0; D* RData = new D[height*width]; { for(unsigned int ai = height / Block * br; ai < height / Block * (br + 1); ai++ ) { for(unsigned int aj = width / Hlock * hr; aj < width / Hlock * (hr + 1); aj++ ) { if( ImgData[ai * step + aj] != 0 ) { RData[Size].data = (unsigned char)ImgData[ai * step + aj]; RData[Size].v.x=aj; RData[Size].v.y=ai; //测试用 //cvLine( im, cvPoint(aj, ai), cvPoint(aj, ai), CV_RGB(0, 0, 255), 1,8, 0); Size++; } } } } //测试用 /* cvNamedWindow( "cirl", 1 ); cvShowImage( "cirl", im ); cvWaitKey(0); */ unsigned int Tx = GG->T1; int nSize = Size; P* pCell = new P; unsigned int Index[3]; double a,b,M,r[3]; int X1,X2,X3,X4; CvPoint P1,P2,P3; cell* Tep = NULL; pCell->pCell = NULL; PRData->D = RData; while(Tx--) { if (nSize >= 10) { //生成随机数 do { //随机数种子 Index[0]=2 + (int)((nSize-1) * rand() / (RAND_MAX + 1.0)); nSize--; Index[1]=2 + (int)((nSize-1) * rand() / (RAND_MAX + 1.0)); nSize--; Index[2]=2 + (int)((nSize-1) * rand() / (RAND_MAX + 1.0)); nSize--; nSize +=3; }while(Index[0] == Index[1] || Index[0] == Index[2]); } else break; P1 = PRData->D[Index[0]].v; P2 = PRData->D[Index[1]].v; P3 = PRData->D[Index[2]].v; X1 = P2.x*P2.x+P2.y*P2.y-P1.x*P1.x-P1.y*P1.y; X2 = 2*(P2.y-P1.y); X3 = P3.x*P3.x+P3.y*P3.y-P1.x*P1.x-P1.y*P1.y; X4 = 2*(P3.y-P1.y); M = 4.0*((P2.x-P1.x)*(P3.y-P1.y)-(P3.x-P1.x)*(P2.y-P1.y)); a = fabs((X1*X4-X2*X3) / M); X2 = 2*(P2.x-P1.x); X4 = 2*(P3.x-P1.x); b = fabs((X2*X3 - X1*X4)/M); r[0] = sqrt((P1.x-a)*(P1.x-a)+(P1.y-b)*(P1.y-b)); r[1] = sqrt((P2.x-a)*(P2.x-a)+(P2.y-b)*(P2.y-b)); r[2] = sqrt((P3.x-a)*(P3.x-a)+(P3.y-b)*(P3.y-b)); if (!cvIsInf(r[0]) && abs(int(r[0] - r[1])) <= 1 && abs(int(r[0] - r[2])) <= 1 && (unsigned int)r[0] <= GG->Max_Lim&& (unsigned int)r[0] >= GG->Mim_Lim) { //测试用 /* if (abs((int)r[0] - 194) < 5) if (abs((int)a - 205) < 3) if (abs((int)b - 205) <3) SS++; */ cell* ca = new cell; ca->score=1; ca->a = (int)a; ca->b = (int)b; ca->r = (int)r[0]; BOOL Flag = FALSE; if (realC > 0) { for (int mk = 0; mk < realC; mk++) { if (sqrt(pow(RealC[mk].r - ca->r, 2) + pow(RealC[mk].a - ca->a, 2) + pow(RealC[mk].b - ca->b, 2)) < GG->A) { Flag = TRUE; delete ca; break; } } } if (Flag == TRUE) continue; if (pCell->pCell == NULL) { pCell->pCell = ca; pCell->eCell = ca; pCell->eCell->cell=NULL; } else { pCell->cCell = pCell->pCell; Tep = NULL; for (int i =0; pCell->cCell != NULL; i++) { if (sqrt(pow(pCell->cCell->r - ca->r, 2) + pow(pCell->cCell->a - ca->a, 2) + pow(pCell->cCell->b - ca->b, 2)) < GG->A) { Flag = TRUE; pCell->cCell->score +=1; if (pCell->cCell->score >= GG->score) { int cc = 0; { BOOL* Buf = new BOOL[nSize]; for (int sss = 0; sss < nSize; sss++) { if (fabs(pCell->cCell->r - sqrt(pow(PRData->D[sss].v.x - pCell->cCell->a, 2) + pow(PRData->D[sss].v.y - pCell->cCell->b, 2))) < GG->SF) { Buf[sss] = TRUE; cc++; } else Buf[sss] = FALSE; } if (cc > GG->K * ca->r) { RealC[realC].a = ca->a; RealC[realC].b = ca->b; RealC[realC].r = ca->r; realC++; Size = nSize - cc; { for(int ai = 0; ai < height; ai++ ) { for(int aj = 0; aj < width ; aj++ ) { if (ImgData[ai * step + aj] != 0) { if (fabs(pCell->cCell->r - sqrt(pow(aj - pCell->cCell->a, 2) + pow(ai - pCell->cCell->b, 2))) < GG->SF) { ImgData[ai * step + aj] =0;