matlab练习程序（直方图反向投影）

￡神魔★判官ぃ 2021-09-17 14:40 468阅读 0赞

做meanshift物体跟踪的时候中间有一步叫做直方图反向投影，所以我就先实现了这样一个步骤。

直方图反向投影说白了就是模板匹配，给定一个较小的目标模板，然后再逐个遍历原图像和模板图像相同的图像块的，对比图像块和模板的直方图，然后把比较结果存入一个新的图像中，新图像中的全局极值就是模板在原图像中所在的位置。这里主要麻烦的是怎么比较两个图像块的直方图，Opencv中实现了5种对比的方法，所以我在这里也对应的实现了5种方法。

5种方法分别是correl（相关）、chisqr（卡方）、intersect（相交）、bhattacharyya（名字很长--!!）、emd（earth mover's distance）。

下面是这5种方法的公式，H1是第一个图像块的直方图序列，H2是第二的图像块的直方图序列：

相关法：  ![2012112212330797.png][]

其中：  ![2012112212340330.png][]

卡方法：![2012112212345418.png][]

相交法：![2012112212364983.png][]

bhattacharyya:![2012112212374734.png][]

emd:

![2012112212485680.gif][]

前4种方法比较简单，emd就比较麻烦了，下面是我自己的理解。这个算法使用原序列生成f和d这两个矩阵，然后对这两个矩阵进行处理。d矩阵按我的理解就是H1中各个元素和H2中各个元素在位置上的距离，比如H1和H2都是只有三个元素的序列，H1(1)和H2(1)的距离就是0，所以d(1,1)=0，H1(3)和H2(2)的距离是1，所以d(3,2)=1，以此d=\[0 1 2;1 0 1;2 1 0\]，这里d和H1，H2中的值没关系。f的确定就比较麻烦了，看看约束条件就让人头疼了。按照约束条件是没办法写程序的，不过我想到一个巧妙的方法求的的结果正好满足那四个约束条件：首先不管前三个条件，只使用第四个条件求得矩阵f，然后对f从最右列开始依次减去相邻的左边的列，然后对f从最下行开始依次减去相邻的上边的行，求得的结果正好满足条件。至于为什么这样我还不得其解，原因就指望哪位数学大牛去探索吧。

下面是代码：

main.m

close all;
    clear all;
    clc;
    
    img=imread('lena.jpg');
    imshow(img);
    [m n]=size(img);
    
    w=imcrop();      %这里把要裁剪的图像框出来
    [H W]=size(w);
    hist1=histcount(w);
    
    HH=floor(H/2);
    WW=floor(W/2);
    
    imgn=zeros(m+2*HH+1,n+2*WW+1);
    imgn(HH+1:m+HH,WW+1:n+WW)=img;
    imgn(1:HH,WW+1:n+WW)=img(1:HH,1:n); 
    imgn(1:m+HH,n+WW+1:n+2*WW+1)=imgn(1:m+HH,n:n+WW);
    imgn(m+HH+1:m+2*HH+1,WW+1:n+2*WW+1)=imgn(m:m+HH,WW+1:n+2*WW+1);
    imgn(1:m+2*HH+1,1:WW)=imgn(1:m+2*HH+1,WW+1:2*WW);
    
    re1=imgn;
    re2=imgn;
    re3=imgn;
    re4=imgn;
    re5=imgn;
    for i=HH+1:m+HH
        for j=WW+1:n+WW
            s=imgn(i-HH:i+HH,j-WW:j+WW);
            hist2=histcount(s);
         
            re1(i,j)=correl(hist1,hist2);           %相关法
            re2(i,j)=chisqr(hist1,hist2);           %卡方法
            re3(i,j)=intersect(hist1,hist2);        %相交法
            re4(i,j)=bhattacharyya(hist1,hist2);    %名字很长的法
            re5(i,j)=emd(hist1,hist2);     %由于没有优化，速度实在太慢了，至少运行一晚上，慎用！！
        end
    end
    figure;
    re1=re1(HH+1:m+HH,WW+1:n+WW);
    imshow(mat2gray(re1));
    
    figure;
    re2=re2(HH+1:m+HH,WW+1:n+WW);
    imshow(mat2gray(re2));
    
    figure;
    re3=re3(HH+1:m+HH,WW+1:n+WW);
    imshow(mat2gray(re3));
    
    figure;
    re4=re4(HH+1:m+HH,WW+1:n+WW);
    imshow(mat2gray(re4));
    
    figure;
    re5=re5(HH+1:m+HH,WW+1:n+WW);
    imshow(mat2gray(re5));

histcount.m 统计直方图

function hist=histcount(w)
        [H W]=size(w);
        w=uint8(w);
        hist=zeros(1,256);
        for i=1:H
            for j=1:W
                hist(w(i,j)+1)=hist(w(i,j)+1)+1;
            end
        end
    
        hist=hist/(H*W);
    
    end

correl.m 相关法

function d=correl(H1,H2)
        d=sum(H1-mean(H1).*(H2-mean(H2)))/sqrt(sum((H1-mean(H1)).^2)*sum((H2-mean(H2)).^2));
    end

chisqr.m 卡方法

function d=chisqr(H1,H2)
        d=sum(((H1-H2).^2)/(H1+H2));
    end

intersect.m 相交法

function d=intersect(H1,H2)    
        d=sum(min(H1,H2));
    end

bhattacharyya.m

function d=bhattacharyya(H1,H2)
        d=sqrt(1-sum(sqrt(H1.*H2))/sqrt(sum(H1)*sum(H2)));
    end

emd.m 这里没优化，慎用！

function re=emd(H1,H2)
        m=length(H1);
        n=length(H2);
        f=zeros(m,n);
        d=zeros(m,n);
        
        for i=1:m
            for j=1:n
                if i==j
                    d(i,j)=0;
                end       
                if j>i
                    d(i,j)=j-i;
                end       
                if j<i
                    d(i,j)=i-j;
                end    
                f(i,j)=min(sum(H1(1:i)),sum(H2(1:j)));              
            end   
        end
        
        for i=m:-1:2
            f(:,i)=f(:,i)-f(:,i-1);
        end
        
        for j=n:-1:2
            f(j,:)=f(j,:)-f(j-1,:);
        end
    
        re=(sum(sum(f.*d)))/sum(sum(f));
        
    end

下面是运行效果：

![2012112213134563.jpg][]原图

![2012112213135837.jpg][]使用的模板，就是lena的右眼。

![2012112213143428.jpg][]相关法

![2012112213145541.jpg][]卡方法

![2012112213150793.jpg][]相交法

![2012112213152097.jpg][]bhattacharyya法

![2012112213154566.jpg][]emd法

emd时间最长，效果还不错。卡方法最快，效果最不好。相交法和bhattacharyya效果都挺不错的，时间也不算慢。

参考：

1.[http://www.cnblogs.com/xrwang/archive/2010/02/04/HowToUseHistogram.html][http_www.cnblogs.com_xrwang_archive_2010_02_04_HowToUseHistogram.html]

2.[http://blog.163.com/woshitony111@126/blog/static/71379539201262202820650/][http_blog.163.com_woshitony111_126_blog_static_71379539201262202820650]

3.[http://homepages.inf.ed.ac.uk/rbf/CVonline/LOCAL\_COPIES/RUBNER/emd.htm][http_homepages.inf.ed.ac.uk_rbf_CVonline_LOCAL_COPIES_RUBNER_emd.htm]

4.[http://vision.stanford.edu/~rubner/papers/rubnerIccv98.pdf][http_vision.stanford.edu_rubner_papers_rubnerIccv98.pdf]

[2012112212330797.png]: /images/20210726/1eb4850359314d8d8a7b71bb28c5ae3e.png
[2012112212340330.png]: /images/20210726/2798bb26de9d422796fc1c6d00a1d0ea.png
[2012112212345418.png]: /images/20210726/4fe0188c7d0543cdafd4c1eb63271a73.png
[2012112212364983.png]: /images/20210726/3602600d73084d30b0e2bfad57e16839.png
[2012112212374734.png]: /images/20210726/cf1d1f67663843359505c36c97fc4d23.png
[2012112212485680.gif]: /images/20210726/1cc02629dbe14721987ac6c1b493b057.png
[2012112213134563.jpg]: /images/20210726/372a1c4ee72948ecb3783eb52cde54cb.png
[2012112213135837.jpg]: /images/20210726/f76eff14ebc34fe18b531a07d4cec018.png
[2012112213143428.jpg]: /images/20210726/e7ab4ac77c0e4f01b4bb7b516f5af55b.png
[2012112213145541.jpg]: /images/20210726/57103808bbd94693aa7f9e41aa0b8692.png
[2012112213150793.jpg]: /images/20210726/90c301082ba94530b667d3996467ea67.png
[2012112213152097.jpg]: /images/20210726/719c072642324835a8927f2f403d209f.png
[2012112213154566.jpg]: /images/20210726/3ea787084322486db210629b79fbaa74.png
[http_www.cnblogs.com_xrwang_archive_2010_02_04_HowToUseHistogram.html]: http://www.cnblogs.com/xrwang/archive/2010/02/04/HowToUseHistogram.html
[http_blog.163.com_woshitony111_126_blog_static_71379539201262202820650]: http://blog.163.com/woshitony111@126/blog/static/71379539201262202820650/
[http_homepages.inf.ed.ac.uk_rbf_CVonline_LOCAL_COPIES_RUBNER_emd.htm]: http://homepages.inf.ed.ac.uk/rbf/CVonline/LOCAL_COPIES/RUBNER/emd.htm
[http_vision.stanford.edu_rubner_papers_rubnerIccv98.pdf]: http://vision.stanford.edu/~rubner/papers/rubnerIccv98.pdf