HOG算子

转自：点击打开链接 

梯度直方图特征(HOG) 是一种对图像局部重叠区域的密集型描述符, 它通过计算局部区域的梯度方向直方图来构成特征。Hog特征结合SVM分类器已经被广泛应用于图像识别中，尤其在行人检测中获得了极大的成功。需要提醒的是，HOG+SVM进行行人检测的方法是法国研究人员Dalal在2005的CVPR上提出的，而如今虽然有很多行人检测算法不断提出，但基本都是以HOG+SVM的思路为主。

        HOG特征是一种局部区域描述符,它通过计算局部区域上的梯度方向直方图来构成人体特征,能够很好地描述人体的边缘。它对光照变化和小量的偏移不敏感。

图像中像素点(x,y)的梯度为

    Dalal提出的Hog特征提取的过程：把样本图像分割为若干个像素的单元（cell），把梯度方向平均划分为9个区间（bin），在每个单元里面对所有像素的梯度方向在各个方向区间进行直方图统计，得到一个9维的特征向量，每相邻的4个单元构成一个块（block），把一个块内的特征向量联起来得到36维的特征向量，用块对样本图像进行扫描，扫描步长为一个单元。最后将所有块的特征串联起来，就得到了人体的特征。例如，对于64*128的图像而言，每2*2的单元（16*16的像素）构成一个块，每个块内有4*9=36个特征，以8个像素为步长，那么，水平方向将有7个扫描窗口，垂直方向将有15个扫描窗口。也就是说，64*128的图片，总共有36*7*15=3780个特征。

在行人检测过程中，除了上面提到的HOG特征提取过程，还包括彩图转灰度，亮度校正等步骤。总结一下，在行人检测中，HOG特征计算的步骤：

（1）将输入的彩图转换为灰度图；

（2）采用Gamma校正法对输入图像进行颜色空间的标准化（归一化）；   目的是调节图像的对比度，降低图像局部的阴影和光照变化所造成的影响，同时可以抑制噪音的干扰；

（3）计算梯度；主要是为了捕获轮廓信息，同时进一步弱化光照的干扰。

（4）将梯度投影到单元的梯度方向；目的是为局部图像区域提供一个编码，

（5）将所有单元格在块上进行归一化；归一化能够更进一步对光照、阴影和边缘进行压缩，通常，每个单元格由多个不同的块共享，但它的归一化是基于不同块的，所以计算结果也不一样。因此，一个单元格的特征会以不同的结果多次出现在最后的向量中。我们将归一化之后的块描述符就称之为HOG描述符。

（6）收集得到检测空间所有块的HOG特征；该步骤就是将检测窗口中所有重叠的块进行HOG特征的收集，并将它们结合成最终的特征向量供分类使用。

function hogf = getHOG(inm)
inm = imresize(inm,[100,100]);

img = double(inm);

%img=imread('f:\111.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 - step          %如果处理的m/step不是整数，最好是i=1:step:m-step
    ii=1;
    for j=1:step:n - step      %注释同上
        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
                if isnan(tmped(p,q))==1  
                    tmped(p,q) = 0;
                end
                Hist(ceil(ang/jiao))=Hist(ceil(ang/jiao))+tmped(p,q);   %ceil向上取整，使用边缘强度加权
            end
        end
        Hist=Hist/sum(Hist);    %方向直方图归一化
        Hist = nan2zero(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 
hogf = f;