HOG（Histograms of Oriented Gradients）用于目标定位检测

Halcom

HOG（Histograms of Oriented Gradients）用于目标定位检测：
视频链接：
视频：http://pan.baidu.com/s/1mi1LFFM
录制的视频是算法底层原理讲解，底层代码实现，方便大家真正掌握算法实质，开发出更加出色的算法。录制视频的初衷是：避免读者朋友利用大把时间学习已有常见算法，本系列视频旨在让读者朋友快速深入了解这些常见算法原理，有更多的时间去研究更加高大上算法（价值）。

具体链接在halcom.cn论坛，联系人QQ：3283892722
该论坛是一个学习交流平台，我会逐一的和大家分享学习。
欢迎大家录制视频，你可在论坛进行打赏分享。
视频专用播放器：http://halcom.cn/forum.php?mod=viewthread&tid=258&extra=page%3D1

       方向梯度直方图（Histogram of Oriented Gradient, HOG）特征通过计算和统计图像局部区域的梯度方向直方图来构成特征。Hog特征结合SVM分类器已经被广泛应用于图像识别中，尤其在行人检测中获得了极大的成功。
       HOG特征提取算法的实现过程（图像目标定位检测）：
1）灰度化，如果是彩色图像，采用img = rgb2gray(img)灰度化；Gray = 0.3*R+0.59*G+0.11*B;
2）采用Gamma校正法对输入图像进行颜色空间的标准化（归一化）；目的是调节图像的对比度，降低图像局部的阴影和光照变化所造成的影响，同时可以抑制噪音的干扰；也可以简便地直接采用 img = sqrt(img);注意：当构建模板时，没有进行Gamma校准，则检测时，也不应该进行Gamma校准。
3）计算图像每个像素的梯度（包括赋值和相位）；主要是为了捕获轮廓信息，相位保证图像信息的特征唯一性。   
   滤波掩膜如[-1,0,1]，具体看sobel、canny、prewitt等（MATLAB图像滤波去噪分析及其应用）
   mag = sqrt( dx^2 + dy^2 )
   phase = atan2(dy, dx)
   dx为水平梯度，dy为垂直梯度；
4）将图像划分成num_bins个cells（在这里，选择每个cell的大小为图像大小，每个cell图像由不同角度成分构成）；

num_bins=9     % 平分360度

5）统计每个cell的梯度方向直方图（一般cell设为9个，将360度9等分，求解每等分角度下的梯度方向直方图），即可形成每个cell的descriptor；

All_cond = mag_cond & phase_cond & (phase >= min_orientation) & (phase < max_orientation);

mag_cond、phase_cond分别为幅值mag、相位phase的二值化图像；min_orientation、max_orientation角度区间，9等分，就有9个区间；All_cond大小与原图像等大。

6）将每个cell分割为N个block（例如8*8的block，一般为图像的1/8），每个block为cell图像的扫描块（互不重叠），对所有block的特征分别求和，再reshape即得到cell图像的HOG特征descriptor。

All_cond_blocks = im2col(All_cond, [8,8], 'distinct');

All_cond_blocks = reshape(sum(All_cond_blocks), hog_height, hog_width);

7）归一化统计的HoG特征；
        利用构建的模板图像的HoG特征，就可以在测试图像中进行查找匹配了。
1）读取测试图像；
2）计算测试图像的HoG特征；
3）将每个梯度方向的HoG特征分别于模板的HoG特征分别卷积，再相加，得到response；

1. response = zeros(hog_height, hog_width);
   
2. for i=1:num_bins
   
3.     response = response + conv2(hog_test(:,:,i), template(:,:,i), 'same');
   
4. end

复制代码

4）求解response的最大值对应的坐标信息x,y，就是测试图像最相似的目标区域；
5）标记测试图像的区域即可。

主程序为：

1. clc,clear,close all
   
2. warning off
   
3. load('template_images_pos.mat');  % 正样本模板图像
   
4. load('template_images_neg.mat');  % 负样本模板图像
   
5. Itest = imread('../Images/test5.jpg');
   
6. if size(Itest,3)>1
   
7.     Itest = rgb2gray(Itest);
   
8. end
   
9. if ~isfloat(Itest)
   
10.     Itest = im2double(Itest);
    
11. end
    
12. % Gamma矫正
    
13. Itest = sqrt(Itest+0.00000001);
    
14. % figure,imshow(Itest)
    
15. num_objects = 1;  % 找最近似的一个目标
    
16. % 构建目标模板
    
17. template = hog_pos(template_images_pos);
    
18. % 目标HoG检测
    
19. [x,y,score] = HoG_detect(Itest,template,num_objects);
    
20. % 绘图
    
21. draw_detection(Itest, num_objects, x, y);

复制代码

计算每个目标模板的HoG特征的 平均值和

1. function template = hog_pos(template_images_pos)
   
2. % 计算每个目标模板的HoG特征的 平均值和
   
3. 
   
4. image_size = size(template_images_pos{1});
   
5. % num_bins = 9
   
6. template_size = [image_size/8, 9];
   
7. template = zeros( template_size );
   
8. 
   
9. for i=1:size(template_images_pos,2)  % 正样本模板的个数
   
10.     template = template + HoG_feature( template_images_pos{i} );  % 累加
    
11. end
    
12. 
    
13. template = template./ size(template_images_pos,2);  % 取平均的操作

复制代码

HoG特征计算（该算子类比于LBP算子，区别于HoG特征分类，HoG目标定位采用整幅图直接切块处理，HoG特征分类时，采用16x16块进行8x8blocks处理，其实是相似的）：

1. function phaseHist = HoG_feature(Itest)
   
2. % HoG特征, 对每一幅图像，求解9个角度范围的特征
   
3. if size(Itest,3)>1
   
4.     Itest = rgb2gray(Itest);
   
5. end
   
6. if ~isfloat(Itest)
   
7.     Itest = im2double(Itest);
   
8. end
   
9. % 求幅值和相位
   
10. [mag, phase] = mygradient(Itest);
    
11. [height,width] = size(Itest);
    
12. hog_height = ceil(height/8);
    
13. hog_width = ceil(width/8);
    
14. 
    
15. % 二值化
    
16. mag_cond = im2bw(mag);
    
17. phase_cond = phase ~= 0.0;
    
18. 
    
19. num_bins = 9;   % 平分360度
    
20. phaseHist = zeros(hog_height, hog_width, num_bins );
    
21. % 计算每一个cell
    
22. for i=1:num_bins
    
23.     min_orientation = -pi + (i-1)*2*pi/num_bins;
    
24.     max_orientation = -pi + i*2*pi/num_bins;
    
25.     All_cond = mag_cond & phase_cond & (phase >= min_orientation) & (phase < max_orientation);
    
26.     % 计算8x8的block块，整个All_cond图像进行8x8块不重叠扫描，不够的区域填充0
    
27.     % im2col作用是： If A = [A11 A12; A21 A22], where each Aij is M-by-N, then B = [A11(:) A21(:) A12(:) A22(:)].
    
28.     All_cond_blocks = im2col(All_cond, [8,8], 'distinct');
    
29.     % 对每一个block求和，再reshape回来
    
30.     All_cond_blocks = reshape(sum(All_cond_blocks), hog_height, hog_width);
    
31.     phaseHist(:,:,i) = All_cond_blocks;
    
32. end
    
33. 
    
34. %% 归一化操作
    
35. phaseHist_sum = sum(phaseHist,3);
    
36. phaseHist_sum = phaseHist_sum + (phaseHist_sum==0);
    
37. phaseHist_sum_norm = zeros(hog_height, hog_width, num_bins );
    
38. for i=1:num_bins
    
39.     phaseHist_sum_norm(:,:,i) = phaseHist_sum;
    
40. end
    
41. phaseHist = phaseHist./phaseHist_sum_norm;

复制代码

HoG_detect算子：

1. function [x,y,score] = HoG_detect(Itest,template,num_objects)
   
2. % Itest:输入的待检测图像（原始图像）
   
3. % template：模板HoG特征图像
   
4. % num_objects  找最近似的目标个数
   
5. % x,y 输出的目标检测位置坐标
   
6. % score：相似度得分
   
7. hog_test = HoG_feature( Itest );
   
8. [hog_height, hog_width, num_bins] = size(hog_test);
   
9. 
   
10. template = rot90(template,2);  % 逆时针旋转2*90度
    
11. 
    
12. response = zeros(hog_height, hog_width);
    
13. for i=1:num_bins
    
14.     response = response + conv2(hog_test(:,:,i), template(:,:,i), 'same');
    
15. end
    
16. 
    
17. % 最大得分的一个目标
    
18. % [x, y] = find( response == max(max( response )) );
    
19. 
    
20. % 网格化划分，灰度图可以显示为一个3D图
    
21. [X,Y] = meshgrid(1:hog_width,1:hog_height);
    
22. columnVector_response = [Y(:), X(:), response(:)];
    
23. % -3表示按照第3列，降序排列
    
24. sorted_response = sortrows(columnVector_response, -3);
    
25. 
    
26. % 输出结果的初始化
    
27. x = zeros(1,num_objects);
    
28. y = zeros(1,num_objects);
    
29. score = zeros(1,num_objects);
    
30. num_objects_current = 0;
    
31. for i=1:size(sorted_response,1)
    
32.     y_ = 8*sorted_response(i,1);
    
33.     x_ = 8*sorted_response(i,2);
    
34.     score_ = sorted_response(i,3);
    
35.     distances = sqrt((y-y_).^2 + (x-x_).^2);
    
36.     num_overlaps = sum(distances < 80*sqrt(2));   % 阈值比较
    
37.    
    
38.     if (num_overlaps == 0)
    
39.         num_objects_current = num_objects_current + 1;
    
40.         y(num_objects_current) = y_ - 4;
    
41.         x(num_objects_current) = x_ - 4;
    
42.         score(num_objects_current) = score_;
    
43.     end
    
44.    
    
45.     if (num_objects_current >= num_objects)  % 最多只取num_objects个目标
    
46.         break
    
47.     end
    
48. end

复制代码

绘图显示结果：

1. function draw_detection(Itest, ndet, x, y, scale)
   
2.     %display top ndet detections
   
3.     figure; clf; imshow(Itest);
   
4.     for i = 1:ndet
   
5.       % draw a rectangle.  use color to encode confidence of detection
   
6.       %  top scoring are green, fading to red
   
7.       hold on;
   
8.       if nargin < 5
   
9.           scale = ones(size(x,2),1);
   
10.       end
    
11.       h = rectangle('Position',[x(i)-round(100/scale(i)) y(i)-round(100/scale(i)) round(200/scale(i)) round(200/scale(i))],'EdgeColor',[(i/ndet) ((ndet-i)/ndet)  0],'LineWidth',3,'Curvature',[0.3 0.3]);
    
12.       hold off;
    
13.     end
    
14. end

复制代码

参考：
【1】INRIA Person Dataset：http://pascal.inrialpes.fr/data/human/
【2】百度网盘INRIA Person数据库链接：http://pan.baidu.com/s/1c2AivU4
【3】http://blog.csdn.net/hujingshuang/article/details/47337707/

Halcom

1. function [mag,phase] = mygradient(I)
   
2. % 求幅值和相位
   
3. if size(I,3)>1
   
4.     I = rgb2gray(I);
   
5. end
   
6. if ~isfloat(I)
   
7.     I = im2double(I);
   
8. end
   
9. 
   
10. gradFilter = fspecial('prewitt');
    
11. dx = imfilter(I, gradFilter, 'replicate');
    
12. dy = imfilter(I, gradFilter', 'replicate');
    
13. mag = sqrt(dx.*dx + dy.*dy);
    
14. phase = atan2(dy, dx);

复制代码