基于陆地移动距离（EMD）的彩色图像直方图距离计算

直方图匹配策略

直方图是用来比较两幅图像相似度的重要途径，常用的距离标准包括：相关、卡方、直方图相交以及Bhattacharyya距离。

光线变化能使得图像颜色值的飘逸，尽管颜色直方图形状没有改变（整体位移），这样的变化使得以上匹配策略失效【我的理解为，此时匹配仍可以得到正确结论，即具有越相似直方图统计特性，会在四种评判中得到更好的结果，但是可能已经丧失其原本的意义了，数值大小不能反应图像之间的差异，即不再具有判断图像相似度的能力】。

举个例子图像A,B,C，我们想从B和C中找到与A更相似的图片（视觉上），当这两者根据某一标准的结果都很差时，即使相对较好的结果我们也不能说它和图像A更相似（虽然在统计上是这样的）。所以需要提出一种比较靠谱的度量方法。

陆地移动距离

更靠谱的距离度量方法就是陆地移动距离（EMD），下面是The Earth Mover’s Distance as a Metric for Image Retrieval（EMD）这篇文献摘要的翻译，看了可能对EMD有更好的认识:

我们研究了一种基于图像内容检索方法的直方图分布特性，陆地移动距离（EMD）。 EMD的意义是将一种直方图分布转化为另一种分布所要付出的最小成本。对于图像检索来说，我们将EMD和矢量量化结合一起，该方法相比以往方法在视觉相似性表现出了更好的方法。EMD相比其他直方图匹配算法，具有更佳的鲁棒性。当用来比较相同质量分布时，EMD度量更加真实。在本文中，我们专注于应用于颜色和纹理，比较EMD了的检索性能与其他距离的检索性能。

因此，EMD距离模型应该还是挺有用途的。

PS：

OPENCV中cvCalcEMD2（）函数能计算的最大直方图签名的行数为15440。
（1）对于1维直方图，区间数目不能超过15440；
（2）对于2维直方图，2个区间数目的乘积不能超过15440，区间相同的话每个区间不能超过124；
（3）对于3维直方图，3个区间数目的乘积不能超过15440，区间相同的话每个区间不能超过24.9；
（4）实际使用的时候，因为同一个程序中的其他部分还要占用内存，所以会比计算的值更少。
具体参见：http://www.cnblogs.com/xrwang/archive/2010/02/02/cvCalcEMD2Exception.html

---

实验：card1、card2分别是光照较强时的图像，card2相对card1有旋转，card3是光照较弱时的图像（由于我是手机拍的，可能不太清楚）
计算了相互间直方图相交以及EMD度量。

#include<cv.h>#include<highgui.h>//输入bgr彩色图像，进行分割，返回计算的直方图//由于彩色图像，维数默认3，size为bin尺寸，range为bin范围//如果 norm =1，返回归一化的直方图CvHistogram* cvCalcbgrHist(IplImage *src,int* size,float** ranges,int norm=1) {    if (src->nChannels==3)    {          //分割单通道        IplImage *b_plane = cvCreateImage(cvGetSize(src),8,1);        IplImage *g_plane = cvCreateImage(cvGetSize(src),8,1);        IplImage *r_plane = cvCreateImage(cvGetSize(src),8,1);        IplImage *plane[] = {b_plane,g_plane,r_plane};        cvCvtPixToPlane(src,b_plane,g_plane,r_plane,0);        CvHistogram *hist;        hist = cvCreateHist(3,size,CV_HIST_ARRAY,ranges,1);        cvCalcHist(plane,hist,0,0);        if(norm ==1)        cvNormalizeHist(hist,1.0);        return hist;    }    else    {      printf("it must be color image");    }}//计算输入两个彩色图像直方图的陆地移动距离（EMD）//输入两个直方图，bins_size相同float cvCalbgrEMD(CvHistogram* hist1,CvHistogram* hist2,int *size){    //首先创建signature，用于cvCalcEMD2的输入    CvMat *sig1;    CvMat *sig2;    int b=0,g=0,r=0;    int b_bins = size[0];    int g_bins = size[1];    int r_bins = size[2];    int numrows = b_bins*g_bins*r_bins;    sig1 = cvCreateMat(numrows,4,CV_32FC1);  //1count + 3 coords = 4    sig2 = cvCreateMat(numrows,4,CV_32FC1);    for( b=0;b<b_bins;b++)        {            for( g=0;g<g_bins;g++)                {                    for( r=0;r<r_bins;r++)                    {                        float bin_val = cvQueryHistValue_3D(hist1,b,g,r);                        int pos = b*g_bins*r_bins+g*r_bins+r;                        cvSet2D(sig1,pos,0,cvScalar(bin_val));                        cvSet2D(sig1,pos,1,cvScalar(b));                        cvSet2D(sig1,pos,2,cvScalar(g));                        cvSet2D(sig1,pos,3,cvScalar(r));                        bin_val = cvQueryHistValue_3D(hist2,b,g,r);                        cvSet2D(sig2,pos,0,cvScalar(bin_val));                        cvSet2D(sig2,pos,1,cvScalar(b));                        cvSet2D(sig2,pos,2,cvScalar(g));                        cvSet2D(sig2,pos,3,cvScalar(r));                    }                 }        }    float emd = cvCalcEMD2(sig1,sig2,CV_DIST_L2);    return emd;}int main(){    IplImage *card1;    IplImage *card2;    IplImage *card3;    card1 = cvLoadImage("card1.jpg");    card2 = cvLoadImage("card2.jpg");    card3 = cvLoadImage("card3.jpg");    //定义基本参数    int b_bins = 20, g_bins =20, r_bins =20;    int hist_size[]  = {b_bins,g_bins,r_bins};    float b_range[] = {0,255};    float g_range[] = {0,255};    float r_range[] = {0,255};    float *range[] = {b_range,g_range,r_range};    //计算直方图    CvHistogram *hist_card1=cvCalcbgrHist(card1,hist_size,range,1);         CvHistogram *hist_card2=cvCalcbgrHist(card2,hist_size,range,1);     CvHistogram *hist_card3=cvCalcbgrHist(card3,hist_size,range,1);     //直方图相交法INTERSECT    double value1 = cvCompareHist(hist_card1,hist_card2,CV_COMP_INTERSECT);    double value2 = cvCompareHist(hist_card1,hist_card3,CV_COMP_INTERSECT);    //EMD    double value3 = cvCalbgrEMD(hist_card1,hist_card2,hist_size);    double value4 = cvCalbgrEMD(hist_card1,hist_card3,hist_size);    cvNamedWindow("card1",1);    cvShowImage("card1",card1);    cvNamedWindow("card2",1);    cvShowImage("card2",card2);    cvNamedWindow("card3",1);    cvShowImage("card3",card3);    printf("INTERSECT Distance : card1&card2 %f\n",value1);    printf("INTERSECT Distance : card1&card3 %f\n",value2);    printf("EMD Distance : card1&card2 %f\n",value3);    printf("EMD Distance : card1&card3 %f\n",value4);    cvWaitKey(0);    return 0;}