直方图比较方法 EMD算法

这里介绍的是一种直方图的比较方法。EMD算法，这里不同于我们熟知的目标检测算法，这里EMD的全称是Earth Mover's Distance，即陆地移动距离，是一种度量准则。光线变化能引起颜色直方图的漂移，引起颜色值位置发生变化，这时候一些直方图的比较匹配策略失效。EMD算法本身是一个通用的算法，它实际上度量的是怎样讲一个直方图的形状转变为另一个直方图的形状，包括直方图的部分（或全部）到一个新的位置，可以在直方图任何维的直方图上进行这种度量。
有关EMD陆地移动距离的具体内容可参考：Learning openCV 之 陆地移动距离，或参考后面附的论文。

另，这种方法是不是可以改进粒子滤波算法，有待进一步探讨。


实现代码：
// histogram_EMD.cpp : Defines the entry point for the console application.
//

#include "stdafx.h"
#include
#include

int _tmain (int argc, _TCHAR * argv[])
{
IplImage* src , *src1;


src=cvLoadImage ("lena.jpg", 1);
src1=cvLoadImage ("lena1.jpg", 1);

// Compute the HSV image and decompose it into separate planes.
//
IplImage* hsv = cvCreateImage( cvGetSize(src ), 8, 3 );
IplImage* hsv1 = cvCreateImage( cvGetSize(src ), 8, 3 );
cvCvtColor( src , hsv, CV_BGR2HSV );
cvCvtColor( src1 , hsv1, CV_BGR2HSV );
IplImage* h_plane = cvCreateImage( cvGetSize(src ), 8, 1 );
IplImage* s_plane = cvCreateImage( cvGetSize(src ), 8, 1 );
IplImage* v_plane = cvCreateImage( cvGetSize(src ), 8, 1 );
IplImage* h_plane1 = cvCreateImage( cvGetSize(src ), 8, 1 );
IplImage* s_plane1 = cvCreateImage( cvGetSize(src ), 8, 1 );
IplImage* v_plane1 = cvCreateImage( cvGetSize(src ), 8, 1 );
IplImage* planes [] = { h_plane, s_plane,v_plane };
IplImage* planes1 [] = { h_plane1, s_plane1,v_plane1 };
cvCvtPixToPlane( hsv , h_plane, s_plane, v_plane , 0 );
cvCvtPixToPlane( hsv1 , h_plane1, s_plane1, v_plane1 , 0 );
// Build the histogram and compute its contents.
//
int h_bins = 30, s_bins = 32, v_bins = 8;
CvHistogram* hist , *hist1;
{
int hist_size [] = { h_bins, s_bins ,v_bins };
float h_ranges [] = { 0, 180 }; // hue is [0,180]
float s_ranges [] = { 0, 255 };
float v_ranges [] = { 0, 255 };
float* ranges [] = { h_ranges, s_ranges ,v_ranges };
hist = cvCreateHist (
3,
hist_size,
CV_HIST_ARRAY,
ranges,
1
);
hist1 = cvCreateHist (
3,
hist_size,

CV_HIST_ARRAY,
ranges,
1
);
}
cvCalcHist( planes , hist, 0, 0 ); //Compute histogram
cvNormalizeHist( hist , 1.0 ); //Normalize it
cvCalcHist( planes1 , hist1, 0, 0 ); //Compute histogram
cvNormalizeHist( hist1 , 1.0 ); //Normalize it

CvMat* sig1 ,*sig2;
int numrows = h_bins* s_bins*v_bins ;
//Create matrices to store the signature in
//
sig1 = cvCreateMat (numrows, 4, CV_32FC1);
sig2 = cvCreateMat (numrows, 4, CV_32FC1);
int h ,s, v;
for( h = 0; h < h_bins; h ++ )
{
for( s = 0; s < s_bins; s ++ )
{
for( v =0; v < v_bins; v ++)
{
float bin_val = cvQueryHistValue_3D( hist, h , s, v );

cvSet2D(sig1 ,h* s_bins*v_bins + s* v_bins+v ,0,cvScalar( bin_val)); //bin value
cvSet2D(sig1 ,h* s_bins*v_bins + s* v_bins+v ,1,cvScalar( h)); //Coord 1
cvSet2D(sig1 ,h* s_bins*v_bins + s* v_bins+v ,2,cvScalar( s)); //Coord 2
cvSet2D(sig1 ,h* s_bins*v_bins + s* v_bins+v ,3,cvScalar( v)); //Coord 3

bin_val = cvQueryHistValue_3D ( hist1, h, s ,v );
cvSet2D(sig2 ,h* s_bins*v_bins + s* v_bins+v ,0,cvScalar( bin_val)); //bin value
cvSet2D(sig2 ,h* s_bins*v_bins + s* v_bins+v ,1,cvScalar( h)); //Coord 1
cvSet2D(sig2 ,h* s_bins*v_bins + s* v_bins+v ,2,cvScalar( s)); //Coord 2
cvSet2D(sig2 ,h* s_bins*v_bins + s* v_bins+v ,3,cvScalar( v)); //Coord 3

}
}
}

float emd = cvCalcEMD2( sig1,sig2 ,CV_DIST_L2);
printf(" 相似度为： %3.1f%%",(1-emd )*100);

cvNamedWindow( "Source1" , 1 );
cvMoveWindow("Source1" , 0, 0);
cvShowImage( "Source1" , src );
cvNamedWindow( "Source2" , 1 );
cvMoveWindow("Source2" , 550, 0);
cvShowImage( "Source2" , src1 );
cvWaitKey(0);

return 0;
}