HOG（Histogram of Oriented Gradient）方向梯度直方图

思路来源：HOG（Histogram of Oriented Gradient）方向梯度直方图特征描述图像局部梯度或边缘方向的分布，选择HOG而不是其他低级特征，如像素的亮度、梯度信息或它们的组合特征，原因在于像素亮度只能表征两个图像块之间亮度差异，而梯度特征对噪声比较敏感，对于梯度方向直方图特征，它不仅能很好的描述图像块局部结构信息，而且对噪声或异常样本具有很好的鲁棒性。

核心思想是一幅图像中的物体的表象和形状可以被梯度或边缘的方向密度分布很好地描述。其实现方法是先将图像分成小的叫做方格单元连通区域；然后采集方格单元中各像素点的梯度方向或边缘方向直方图；最后把这些直方图组合起来就可以构成特征描述子。为了提高精确度，还可以把这些局部直方图在图像的更大的区间(block)中进行对比度归一化，此方法通过先计算各直方图在这个区间(block)中的密度，然后根据这个密度值对区间中的各个方格单元做归一化。通过这个归一化后，能对光照变化和阴影获得更好的稳定性。

优点：首先，由于HOG是在图像的局部方格单元上操作，所以它对图像几何的和光学的形变都能保持很好的不变性，这两种形变只会出现在方向梯度直方图特征的核心思想是一幅图像中的物体的表象和形状可以被梯度或边缘的方向密度分布很好地描述。其实现方法是先将图像分成小的叫做方格单元连通区域；然后采集方格单元中各像素点的梯度方向或边缘方向直方图；最后把这些直方图组合起来就可以构成特征描述子。为了提高精确度，还可以把这些局部直方图在图像的更大的区间(block)中进行对比度归一化，此方法通过先计算各直方图在这个区间(block)中的密度，然后根据这个密度值对区间中的各个方格单元做归一化。通过这个归一化后，能对光照变化和阴影获得更好的稳定性。

算法流程图：



简单步骤：

1.对原图像gamma校正，img=sqrt(img);

2.求图像竖直边缘，水平边缘，边缘强度，边缘斜率。

3.将图像每16*16（取其他也可以）个像素分到一个cell中。对于256*256的lena来说，就分成了16*16个cell了。

4.对于每个cell求其梯度方向直方图。通常取9（取其他也可以）个方向（特征），也就是每360/9=40度分到一个方向，方向大小按像素边缘强度加权。最后归一化直方图。

5.每2*2（取其他也可以）个cell合成一个block，所以这里就有（16-1）*（16-1）=225个block。

6.所以每个block中都有2*2*9个特征，一共有225个block，所以总的特征有225*36个。

当然一般HOG特征都不是对整幅图像取的，而是对图像中的一个滑动窗口取的。
输入lena图像：



HOG特征：



左边坐标为block的个数，右边左边为每一个block中特征数。

Matlab代码：

clear all; close all; clc;

img=double(imread('lena.jpg'));
imshow(img,[]);
[m n]=size(img);

img=sqrt(img);      %伽马校正

%下面是求边缘
fy=[-1 0 1];        %定义竖直模板
fx=fy';             %定义水平模板
Iy=imfilter(img,fy,'replicate');    %竖直边缘
Ix=imfilter(img,fx,'replicate');    %水平边缘
Ied=sqrt(Ix.^2+Iy.^2);              %边缘强度
Iphase=Iy./Ix;              %边缘斜率，有些为inf,-inf,nan，其中nan需要再处理一下

%下面是求cell
step=16;                %step*step个像素作为一个单元
orient=9;               %方向直方图的方向个数
jiao=360/orient;        %每个方向包含的角度数
Cell=cell(1,1);              %所有的角度直方图,cell是可以动态增加的，所以先设了一个
ii=1;
jj=1;
for i=1:step:m          %如果处理的m/step不是整数，最好是i=1:step:m-step
ii=1;
for j=1:step:n      %注释同上
tmpx=Ix(i:i+step-1,j:j+step-1);
tmped=Ied(i:i+step-1,j:j+step-1);
tmped=tmped/sum(sum(tmped));        %局部边缘强度归一化
tmpphase=Iphase(i:i+step-1,j:j+step-1);
Hist=zeros(1,orient);               %当前step*step像素块统计角度直方图,就是cell
for p=1:step
for q=1:step
if isnan(tmpphase(p,q))==1  %0/0会得到nan，如果像素是nan，重设为0
tmpphase(p,q)=0;
end
ang=atan(tmpphase(p,q));    %atan求的是[-90 90]度之间
ang=mod(ang*180/pi,360);    %全部变正，-90变270
if tmpx(p,q)<0              %根据x方向确定真正的角度
if ang<90               %如果是第一象限
ang=ang+180;        %移到第三象限
end
if ang>270              %如果是第四象限
ang=ang-180;        %移到第二象限
end
end
ang=ang+0.0000001;          %防止ang为0
Hist(ceil(ang/jiao))=Hist(ceil(ang/jiao))+tmped(p,q);   %ceil向上取整，使用边缘强度加权
end
end
Hist=Hist/sum(Hist);    %方向直方图归一化
Cell{ii,jj}=Hist;       %放入Cell中
ii=ii+1;                %针对Cell的y坐标循环变量
end
jj=jj+1;                    %针对Cell的x坐标循环变量
end

%下面是求feature,2*2个cell合成一个block,没有显式的求block
[m n]=size(Cell);
feature=cell(1,(m-1)*(n-1));
for i=1:m-1
for j=1:n-1
f=[];
f=[f Cell{i,j}(:)' Cell{i,j+1}(:)' Cell{i+1,j}(:)' Cell{i+1,j+1}(:)'];
feature{(i-1)*(n-1)+j}=f;
end
end

%到此结束，feature即为所求
%下面是为了显示而写的
l=length(feature);
f=[];
for i=1:l
f=[f;feature{i}(:)'];
end
figure
mesh(f)

Reference:《NavneetDalalThesis》

所有关于HOG代码：方向梯度直方图HOG