边缘方向直方图

1. //
   //输入：
   //    image: 灰度图像
   //    dst：边缘图像
   //    edge_hist_feature：边缘直方图特征
   //     #define EDGE_HIST_SIZE 16
   //    double edge_feature[EDGE_HIST_SIZE];
   int calEdgeHistogram(IplImage *image,IplImage * dst,double * edge_hist_feature)
   {
       if(NULL ==  image || NULL == edge_hist_feature)
           return 1;
    
       CvHistogram *hist = 0; // 直方图
       IplImage* canny;//边缘图像
       IplImage* gradient_im;
       IplImage* dx; // x方向的sobel差分
       IplImage* dy; // y方向的sobel差分 
    
       CvMat* canny_m;
       CvMat* gradient; // 梯度值
       CvMat* gradient_dir; //梯度的方向
       CvMat* dx_m; // 
       CvMat* dy_m;
       CvMat* mask;
       CvSize  size;
    
       int i,j;
       float theta;
       float max_val;
    
       int hdims = EDGE_HIST_SIZE;     // 划分HIST的个数，越高越精确
       float hranges_arr[] = {-PI/2,PI/2}; // 直方图的上界和下界
       float* hranges = hranges_arr;
    
       size = cvGetSize(image);
       canny = cvCreateImage(cvGetSize(image),8,1);//边缘图像
       dx = cvCreateImage(cvGetSize(image),32,1);//x方向上的差分
       dy = cvCreateImage(cvGetSize(image),32,1);
       gradient_im = cvCreateImage(cvGetSize(image),32,1);//梯度图像
    
       if(NULL == canny || NULL == dx || NULL == dy || NULL == gradient_im)
           return 2;
    
       canny_m = cvCreateMat(size.height,size.width,CV_32FC1);//边缘矩阵
       dx_m = cvCreateMat(size.height,size.width,CV_32FC1);
       dy_m = cvCreateMat(size.height,size.width,CV_32FC1);
    
       gradient = cvCreateMat(size.height,size.width,CV_32FC1);//梯度矩阵
       gradient_dir = cvCreateMat(size.height,size.width,CV_32FC1);//梯度方向矩阵
       mask = cvCreateMat(size.height,size.width,CV_32FC1);//掩码
       if(NULL == canny_m || NULL == dx_m || NULL == dy_m || NULL == gradient || NULL == gradient_dir || NULL == mask)
           return 2;
    
       cvCanny(image,canny,60,180,3);//边缘检测
       cvCopy(canny,dst);
    
       cvConvert(canny,canny_m);//把图像转换为矩阵
       cvSobel(image,dx,1,0,3);// 一阶X方向的图像差分 :dx
       cvSobel(image,dy,0,1,3);// 一阶Y方向的图像差分 :dy
       cvConvert(dx,dx_m);
       cvConvert(dy,dy_m);
    
       //用cvAdd近似计算梯度
       cvAdd(dx_m,dy_m,gradient); // 梯度值
       cvDiv(dx_m,dy_m,gradient_dir); // 梯度方向
       for(i = 0;i< size.height;i++){
           for(j = 0;j < size.width;j++){
               if(cvmGet(canny_m,i,j)!=0 && cvmGet(dx_m,i,j)!=0){
                   theta = cvmGet(gradient_dir,i,j);
                   theta = atan(theta);
                   cvmSet(gradient_dir,i,j,theta);  
               }else{
                   cvmSet(gradient_dir,i,j,0);
               }
           }
       }
    
       hist = cvCreateHist( 1, &hdims, CV_HIST_ARRAY, &hranges, 1 );  // 创建一个指定尺寸的直方图，并返回创建的直方图指针
       //cvZero(hist); // 清0；
       cvConvert(gradient_dir,gradient_im);//把梯度方向矩阵转化为图像
       cvCalcHist( &gradient_im, hist, 0, canny ); // 计算直方图
       cvGetMinMaxHistValue( hist, 0, &max_val, 0, 0 );  // 只找最大值
       cvConvertScale( hist->bins, hist->bins, max_val ? 255. / max_val : 0., 0 ); // 缩放 bin 到区间 [0,255] ，比例系数
    
       for(int i = 0;i < hdims; i++){
           edge_hist_feature[i] = cvGetReal1D(hist->bins,i) / 255.0;
       }
    
       cvReleaseHist (&hist);
       cvReleaseImage(&gradient_im);
       cvReleaseImage(&canny);
       cvReleaseImage(&dx);    
       cvReleaseImage(&dy);
    
       cvReleaseMat(&canny_m);
       cvReleaseMat(&gradient);
       cvReleaseMat(&gradient_dir);
       cvReleaseMat(&dx_m);
       cvReleaseMat(&dy_m);
       cvReleaseMat(&mask);
    
       return 0;
   }