图像直方图比较

本文主要是介绍图像直方图比较，希望对大家解决编程问题提供一定的参考价值，需要的开发者们随着小编来一起学习吧！

对于直方图的比较，我们可以使用 OpenCV 提供的函数 compareHist() 进行比较，从而得到一个数值，表示两个直方图的匹配程度（相似性）。

原理

对于两个直方图（ $H_{1}$ 和 $H_{2}$ ），我们首先选用一个指标（ $d(H_{1}, H_{2})$ ），表示两个直方图的匹配度。

根据 Histogram Comparison 参考资料^[1]，OpenCV的 compareHist 函数计算直方图匹配度所提供的指标有以下4种：

HISTCMP_CORREL：相关性
HISTCMP_CHISQR：卡方
HISTCMP_INTERSECT：交集
HISTCMP_BHATTACHARYYA：巴氏距离

而根据 HistCompMethods^[2]，OpenCV的 compareHist 函数计算直方图匹配度所提供的指标有以下7种（实际为6种，巴氏和海林格算一种）：

HISTCMP_CORREL：相关性（Correlation）

$d(H_1,H_2) = \frac{\sum_I (H_1(I) - \bar{H_1}) (H_2(I) - \bar{H_2})}{\sqrt{\sum_I(H_1(I) - \bar{H_1})^2 \sum_I(H_2(I) - \bar{H_2})^2}}$

其中，

$\bar{H_1} = \frac{1}{N}\sum_I H_1(I)$

且 $N$ 是直方图桶的数量。
HISTCMP_CHISQR：卡方（Chi-Square）

$d(H_1,H_2) = \sum _I \frac{\left(H_1(I)-H_2(I)\right)^2}{H_1(I)}$
HISTCMP_INTERSECT：交集（Intersection）

$d(H_1,H_2) = \sum _I \min (H_1(I), H_2(I))$
HISTCMP_BHATTACHARYYA：巴氏距离（Bhattacharyya distance），事实上 OpenCV 计算的是与巴氏距离系数相关的海林格距离（Hellinger distance）。

$d(H_1,H_2) = \sqrt{1 - \frac{1}{\sqrt{\bar{H_1} \bar{H_2} N^2}} \sum_I \sqrt{H_1(I) \cdot H_2(I)}}$
HISTCMP_HELLINGER：海林格距离（Hellinger distance），同巴氏距离。
HISTCMP_CHISQR_ALT：修正卡方（Alternative Chi-Square）

$d(H_1,H_2) = 2 × \sum _I \frac{\left(H_1(I)-H_2(I)\right)^2}{H_1(I)+H_2(I)}$

常被用于材质比较。
HISTCMP_KL_DIV：库尔巴克-莱布勒散度（Kullback-Leibler divergence）

$d(H_1,H_2) = \sum _I H_1(I) \log \left(\frac{H_1(I)}{H_2(I)}\right)$

OpenCV 中的声明与定义（以 4.10 版为例）

OpenCV 中 compareHist函数声明：

  /** @brief Compares two histograms.The function cv::compareHist compares two dense or two sparse histograms using the specified method.The function returns \f$d(H_1, H_2)\f$ .While the function works well with 1-, 2-, 3-dimensional dense histograms, it may not be suitablefor high-dimensional sparse histograms. In such histograms, because of aliasing and samplingproblems, the coordinates of non-zero histogram bins can slightly shift. To compare such histogramsor more general sparse configurations of weighted points, consider using the #EMD function.@param H1 First compared histogram.@param H2 Second compared histogram of the same size as H1 .@param method Comparison method, see #HistCompMethods*/CV_EXPORTS_W double compareHist( InputArray H1, InputArray H2, int method );/** @overload */CV_EXPORTS double compareHist( const SparseMat& H1, const SparseMat& H2, int method );

OpenCV 中 HistCompMethods 枚举：

  enum HistCompMethods {/** Correlation\f[d(H_1,H_2) =  \frac{\sum_I (H_1(I) - \bar{H_1}) (H_2(I) - \bar{H_2})}{\sqrt{\sum_I(H_1(I) - \bar{H_1})^2 \sum_I(H_2(I) - \bar{H_2})^2}}\f]where\f[\bar{H_k} =  \frac{1}{N} \sum _J H_k(J)\f]and \f$N\f$ is a total number of histogram bins. */HISTCMP_CORREL        = 0,/** Chi-Square\f[d(H_1,H_2) =  \sum _I  \frac{\left(H_1(I)-H_2(I)\right)^2}{H_1(I)}\f] */HISTCMP_CHISQR        = 1,/** Intersection\f[d(H_1,H_2) =  \sum _I  \min (H_1(I), H_2(I))\f] */HISTCMP_INTERSECT     = 2,/** Bhattacharyya distance(In fact, OpenCV computes Hellinger distance, which is related to Bhattacharyya coefficient.)\f[d(H_1,H_2) =  \sqrt{1 - \frac{1}{\sqrt{\bar{H_1} \bar{H_2} N^2}} \sum_I \sqrt{H_1(I) \cdot H_2(I)}}\f] */HISTCMP_BHATTACHARYYA = 3,HISTCMP_HELLINGER     = HISTCMP_BHATTACHARYYA, //!< Synonym for HISTCMP_BHATTACHARYYA/** Alternative Chi-Square\f[d(H_1,H_2) =  2 * \sum _I  \frac{\left(H_1(I)-H_2(I)\right)^2}{H_1(I)+H_2(I)}\f]This alternative formula is regularly used for texture comparison. See e.g. @cite Puzicha1997 */HISTCMP_CHISQR_ALT    = 4,/** Kullback-Leibler divergence\f[d(H_1,H_2) = \sum _I H_1(I) \log \left(\frac{H_1(I)}{H_2(I)}\right)\f] */HISTCMP_KL_DIV        = 5};

OpenCV 中 compareHist 函数定义：

// C O M P A R E   H I S T O G R A M S double cv::compareHist( InputArray _H1, InputArray _H2, int method )
{CV_INSTRUMENT_REGION();Mat H1 = _H1.getMat(), H2 = _H2.getMat();const Mat* arrays[] = {&H1, &H2, 0};Mat planes[2];NAryMatIterator it(arrays, planes);double result = 0;int j;CV_Assert( H1.type() == H2.type() && H1.depth() == CV_32F );double s1 = 0, s2 = 0, s11 = 0, s12 = 0, s22 = 0;CV_Assert( it.planes[0].isContinuous() && it.planes[1].isContinuous() );for( size_t i = 0; i < it.nplanes; i++, ++it ){const float* h1 = it.planes[0].ptr<float>();const float* h2 = it.planes[1].ptr<float>();const int len = it.planes[0].rows*it.planes[0].cols*H1.channels();j = 0;if( (method == CV_COMP_CHISQR) || (method == CV_COMP_CHISQR_ALT)){for( ; j < len; j++ ){double a = h1[j] - h2[j];double b = (method == CV_COMP_CHISQR) ? h1[j] : h1[j] + h2[j];if( fabs(b) > DBL_EPSILON )result += a*a/b;}}else if( method == CV_COMP_CORREL ){
#if (CV_SIMD_64F || CV_SIMD_SCALABLE_64F)v_float64 v_s1 = vx_setzero_f64();v_float64 v_s2 = vx_setzero_f64();v_float64 v_s11 = vx_setzero_f64();v_float64 v_s12 = vx_setzero_f64();v_float64 v_s22 = vx_setzero_f64();for ( ; j <= len - VTraits<v_float32>::vlanes(); j += VTraits<v_float32>::vlanes()){v_float32 v_a = vx_load(h1 + j);v_float32 v_b = vx_load(h2 + j);// 0-1v_float64 v_ad = v_cvt_f64(v_a);v_float64 v_bd = v_cvt_f64(v_b);v_s12 = v_muladd(v_ad, v_bd, v_s12);v_s11 = v_muladd(v_ad, v_ad, v_s11);v_s22 = v_muladd(v_bd, v_bd, v_s22);v_s1 = v_add(v_s1, v_ad);v_s2 = v_add(v_s2, v_bd);// 2-3v_ad = v_cvt_f64_high(v_a);v_bd = v_cvt_f64_high(v_b);v_s12 = v_muladd(v_ad, v_bd, v_s12);v_s11 = v_muladd(v_ad, v_ad, v_s11);v_s22 = v_muladd(v_bd, v_bd, v_s22);v_s1 = v_add(v_s1, v_ad);v_s2 = v_add(v_s2, v_bd);}s12 += v_reduce_sum(v_s12);s11 += v_reduce_sum(v_s11);s22 += v_reduce_sum(v_s22);s1 += v_reduce_sum(v_s1);s2 += v_reduce_sum(v_s2);
#elif CV_SIMD && 0 //Disable vectorization for CV_COMP_CORREL if f64 is unsupported due to low precisionv_float32 v_s1 = vx_setzero_f32();v_float32 v_s2 = vx_setzero_f32();v_float32 v_s11 = vx_setzero_f32();v_float32 v_s12 = vx_setzero_f32();v_float32 v_s22 = vx_setzero_f32();for (; j <= len - VTraits<v_float32>::vlanes(); j += VTraits<v_float32>::vlanes()){v_float32 v_a = vx_load(h1 + j);v_float32 v_b = vx_load(h2 + j);v_s12 = v_muladd(v_a, v_b, v_s12);v_s11 = v_muladd(v_a, v_a, v_s11);v_s22 = v_muladd(v_b, v_b, v_s22);v_s1 += v_a;v_s2 += v_b;}s12 += v_reduce_sum(v_s12);s11 += v_reduce_sum(v_s11);s22 += v_reduce_sum(v_s22);s1 += v_reduce_sum(v_s1);s2 += v_reduce_sum(v_s2);
#endiffor( ; j < len; j++ ){double a = h1[j];double b = h2[j];s12 += a*b;s1 += a;s11 += a*a;s2 += b;s22 += b*b;}}else if( method == CV_COMP_INTERSECT ){
#if (CV_SIMD_64F || CV_SIMD_SCALABLE_64F)v_float64 v_result = vx_setzero_f64();for ( ; j <= len - VTraits<v_float32>::vlanes(); j += VTraits<v_float32>::vlanes()){v_float32 v_src = v_min(vx_load(h1 + j), vx_load(h2 + j));v_result = v_add(v_result, v_add(v_cvt_f64(v_src), v_cvt_f64_high(v_src)));}result += v_reduce_sum(v_result);
#elif CV_SIMDv_float32 v_result = vx_setzero_f32();for (; j <= len - VTraits<v_float32>::vlanes(); j += VTraits<v_float32>::vlanes()){v_float32 v_src = v_min(vx_load(h1 + j), vx_load(h2 + j));v_result = v_add(v_result, v_src);}result += v_reduce_sum(v_result);
#endiffor( ; j < len; j++ )result += std::min(h1[j], h2[j]);}else if( method == CV_COMP_BHATTACHARYYA ){
#if (CV_SIMD_64F || CV_SIMD_SCALABLE_64F)v_float64 v_s1 = vx_setzero_f64();v_float64 v_s2 = vx_setzero_f64();v_float64 v_result = vx_setzero_f64();for ( ; j <= len - VTraits<v_float32>::vlanes(); j += VTraits<v_float32>::vlanes()){v_float32 v_a = vx_load(h1 + j);v_float32 v_b = vx_load(h2 + j);v_float64 v_ad = v_cvt_f64(v_a);v_float64 v_bd = v_cvt_f64(v_b);v_s1 = v_add(v_s1, v_ad);v_s2 = v_add(v_s2, v_bd);v_result = v_add(v_result, v_sqrt(v_mul(v_ad, v_bd)));v_ad = v_cvt_f64_high(v_a);v_bd = v_cvt_f64_high(v_b);v_s1 = v_add(v_s1, v_ad);v_s2 = v_add(v_s2, v_bd);v_result = v_add(v_result, v_sqrt(v_mul(v_ad, v_bd)));}s1 += v_reduce_sum(v_s1);s2 += v_reduce_sum(v_s2);result += v_reduce_sum(v_result);
#elif CV_SIMD && 0 //Disable vectorization for CV_COMP_BHATTACHARYYA if f64 is unsupported due to low precisionv_float32 v_s1 = vx_setzero_f32();v_float32 v_s2 = vx_setzero_f32();v_float32 v_result = vx_setzero_f32();for (; j <= len - VTraits<v_float32>::vlanes(); j += VTraits<v_float32>::vlanes()){v_float32 v_a = vx_load(h1 + j);v_float32 v_b = vx_load(h2 + j);v_s1 += v_a;v_s2 += v_b;v_result += v_sqrt(v_a * v_b);}s1 += v_reduce_sum(v_s1);s2 += v_reduce_sum(v_s2);result += v_reduce_sum(v_result);
#endiffor( ; j < len; j++ ){double a = h1[j];double b = h2[j];result += std::sqrt(a*b);s1 += a;s2 += b;}}else if( method == CV_COMP_KL_DIV ){for( ; j < len; j++ ){double p = h1[j];double q = h2[j];if( fabs(p) <= DBL_EPSILON ) {continue;}if(  fabs(q) <= DBL_EPSILON ) {q = 1e-10;}result += p * std::log( p / q );}}elseCV_Error( cv::Error::StsBadArg, "Unknown comparison method" );}if( method == CV_COMP_CHISQR_ALT )result *= 2;else if( method == CV_COMP_CORREL ){size_t total = H1.total();double scale = 1./total;double num = s12 - s1*s2*scale;double denom2 = (s11 - s1*s1*scale)*(s22 - s2*s2*scale);result = std::abs(denom2) > DBL_EPSILON ? num/std::sqrt(denom2) : 1.;}else if( method == CV_COMP_BHATTACHARYYA ){s1 *= s2;s1 = fabs(s1) > FLT_EPSILON ? 1./std::sqrt(s1) : 1.;result = std::sqrt(std::max(1. - result*s1, 0.));}return result;
}

应用示例

比较两个图像的灰度图匹配度。

#include "opencv2/imgcodecs.hpp"
#include "opencv2/highgui.hpp"
#include "opencv2/imgproc.hpp"
#include <iostream>using namespace cv;
using std::cout;
using std::endl;
using std::string;// 命令行参数
const char* keys ="{ help  h| | Print help message. }""{ @input1 |wukong1.jpg | Path to input image 1. }""{ @input2 |wukong2.jpg | Path to input image 2. }";const string CMP_METHODS[]{"Correlation", "Chi-Squared", "Intersection", "Bhattacharyya", "Hellinger", "Alternative Chi-Squared", "Kullback-Leibler Divergence"};int main(int argc, char** argv)
{// 命令行解析CommandLineParser parser(argc, argv, keys);// 加载图像Mat img1 = imread(parser.get<String>("@input1"));Mat img2 = imread(parser.get<String>("@input2"));// 检查图像是否加载成功if (img1.empty() || img2.empty()){cout << "Could not open or find the image!\n" << endl;cout << "Usage: " << argv[0] << " <Input_Image1> <Input_Image2>" << endl;parser.printMessage();return EXIT_FAILURE;}// 转换为灰度图Mat gray1, gray2;cvtColor(img1, gray1, COLOR_BGR2GRAY);cvtColor(img2, gray2, COLOR_BGR2GRAY);// 计算直方图int histSize = 256; // 256个灰度级，桶数// 灰度级范围float range[] = { 0, 256 };const float* histRange = { range };// 其它直方图计算参数bool uniform = true, accumulate = false;Mat hist1, hist2;calcHist(&gray1, 1, 0, Mat(), hist1, 1, &histSize, &histRange, uniform, accumulate);calcHist(&gray2, 1, 0, Mat(), hist2, 1, &histSize, &histRange, uniform, accumulate);// 归一化直方图normalize(hist1, hist1, 0, 1, NORM_MINMAX, -1, Mat());normalize(hist2, hist2, 0, 1, NORM_MINMAX, -1, Mat());// 比较直方图double correlation = compareHist(hist1, hist2, HISTCMP_CORREL);cout << "Similarity between the two images with " << CMP_METHODS[0] << ": " << correlation << endl;double chiSquared = compareHist(hist1, hist2, HISTCMP_CHISQR);cout << "Similarity between the two images with " << CMP_METHODS[1] << ": " << chiSquared << endl;double intersection = compareHist(hist1, hist2, HISTCMP_INTERSECT);cout << "Similarity between the two images with " << CMP_METHODS[2] << ": " << intersection << endl;double bhattacharyya = compareHist(hist1, hist2, HISTCMP_BHATTACHARYYA);cout << "Similarity between the two images with " << CMP_METHODS[3] << ": " << bhattacharyya << endl;double hellinger = compareHist(hist1, hist2, HISTCMP_HELLINGER);cout << "Similarity between the two images with " << CMP_METHODS[4] << ": " << hellinger << endl;double alternativeChiSquared = compareHist(hist1, hist2, HISTCMP_CHISQR_ALT);cout << "Similarity between the two images with " << CMP_METHODS[5] << ": " << alternativeChiSquared << endl;double kullbackLeibler = compareHist(hist1, hist2, HISTCMP_KL_DIV);cout << "Similarity between the two images with " << CMP_METHODS[6] << ": " << kullbackLeibler << endl;// 显示图像imshow("Image 1", img1);imshow("Image 2", img2);waitKey(0);return EXIT_SUCCESS;
}

import cv2
import numpy as np# 读取图像
image1_path = '../data/Histogram_Comparison_Source_0.jpg'
image2_path = '../data/Histogram_Comparison_Source_1.jpg'image1 = cv2.imread(image1_path)
image2 = cv2.imread(image2_path)# 转换为灰度图像
gray_image1 = cv2.cvtColor(image1, cv2.COLOR_BGR2GRAY)
gray_image2 = cv2.cvtColor(image2, cv2.COLOR_BGR2GRAY)# 计算灰度图像的直方图
hist1 = cv2.calcHist([gray_image1], [0], None, [256], [0, 256])
hist2 = cv2.calcHist([gray_image2], [0], None, [256], [0, 256])# 归一化直方图
cv2.normalize(hist1, hist1, alpha=0, beta=1, norm_type=cv2.NORM_MINMAX)
cv2.normalize(hist2, hist2, alpha=0, beta=1, norm_type=cv2.NORM_MINMAX)# 使用compareHist函数比较直方图
comparison_method = cv2.HISTCMP_CORREL  # 相关性比较方法
corr_similarity = cv2.compareHist(hist1, hist2, comparison_method)
chi_squared_similarity = cv2.compareHist(hist1, hist2, cv2.HISTCMP_CHISQR)
intersection_similarity = cv2.compareHist(hist1, hist2, cv2.HISTCMP_INTERSECT)
bhattacharyya_similarity = cv2.compareHist(hist1, hist2, cv2.HISTCMP_BHATTACHARYYA)
hellinger_similarity = cv2.compareHist(hist1, hist2, cv2.HISTCMP_HELLINGER)
alternative_chi_squared_similarity = cv2.compareHist(hist1, hist2, cv2.HISTCMP_CHISQR_ALT)
kbl_divergence_similarity = cv2.compareHist(hist1, hist2, cv2.HISTCMP_KL_DIV)print(f"Similarity between the two images with Correlation: {corr_similarity}")
print(f"Similarity between the two images with Chi-Squared: {chi_squared_similarity}")
print(f"Similarity between the two images with Intersection: {intersection_similarity}")
print(f"Similarity between the two images with Bhattacharyya: {bhattacharyya_similarity}")
print(f"Similarity between the two images with Hellinger: {hellinger_similarity}")
print(f"Similarity between the two images with Alternative Chi-Squared: {alternative_chi_squared_similarity}")
print(f"Similarity between the two images with Kullback-Leibler Divergence: {kbl_divergence_similarity}")# 显示结果
cv2.imshow('Image 1', image1)
cv2.imshow('Image 2', image2)
cv2.waitKey(0)
cv2.destroyAllWindows()

Similarity between the two images with Correlation: 0.5185838942269291
Similarity between the two images with Chi-Squared: 53.894400613305585
Similarity between the two images with Intersection: 25.52862914837897
Similarity between the two images with Bhattacharyya: 0.358873990739656
Similarity between the two images with Hellinger: 0.358873990739656
Similarity between the two images with Alternative Chi-Squared: 34.90277478480459
Similarity between the two images with Kullback-Leibler Divergence: 71.12863163014705