直方图均衡化（OpenCV4）

关键词　　

　　直方图（histogram）：直方图是图像的灰度——像素数统计图，即对于每个灰度值，统计在图像中具有该灰度值的像素个数，并绘制成图形，称为灰度直方图（简称直方图）。

　　直方图模型（Gray-level histogram）表示图像中不同灰度级出现的相对频率。

详细介绍

　　直方图均衡化通常用来扩大图像的动态范围。容易想到的是，当一幅图像的灰度值大部分集中在一道灰度区间里时，因为像素的亮度都相近，会使图像的细节不明显，整幅图像看起来很模糊。简而言之，直方图均衡化是通过拉伸像素强度分布范围来增强图像对比度的一种方法。

　　数字图像是离散化的数值矩阵，它的直方图可被视为一个离散函数，表示数字图像中每个灰度级与其出现概率见的统计关系。假设一副数字图像f(x,y)的像素总数为N，r_k表示第k个灰度级对应的灰度，n_k表示灰度为r_k的像素个数（频数），若用横坐标表示灰度级，用纵坐标表示频数，则直方图可被定义为P(r_k)=n_k/N，其中，P(r_k)表示灰度r_k出现的相对频数即概率。直方图在一定的程度上能够反映数字图像的概貌性描述，包括图像的灰度范围、灰度分布、整幅图像的亮度均值和明暗对比度等，并可以此为基础进行分析来得出对图像进一步处理的重要依据。直方图均衡化也叫做直方图均匀化，就是把给定图像的直方图变成均匀分布的直方图，是一种较为常用的灰度增强算法。直方图均衡化通常包括以下三个主要步骤。

　　（1）预处理。输入图像，计算图像的直方图。

　　（2）灰度变换表。根据输入图像的直方图计算灰度值变换表。

　　（3）查表变换。执行变换x'=H(x)，表示堆在步骤1中得到的直方图使用步骤2得到的灰度值变换表进行查表的操作，通过遍历整幅图像的每一个像素元，将原始图像灰度值x放入变换表H(x)，可得到变换后得到的新灰度值x'。

　　根据信息论，图像在经过直方图均衡化后，将会包含更多的信息量进而能突出某些图像特征。假设图像具有n级灰度，其第i级灰度出现的概率为p，则该级灰度所含的信息量为：

 　　整幅图像的信息量为：

　　信息论已经证明，具有均匀分布直方图的图像，其信息量H最大。即当p₀=p₁=...=p_n-1=1/n时，上式有最大值。

 直方图模型

　　记r_k∈[0,1]为灰度级，灰度级为r_k的像素个数。则直方图：h(r_k)=n_k。则归一化的直方图：

　　n：像素总数；p_k(r_k)：原始图像灰度分布的概率密度函数。如果将r_k归一化到[0,1]之间，则r_k便可以视作区间[0,1]的随机变量。直方图均衡化处理：假设原图的灰度化值变量为r，变换后新图的灰度值变量为s，我们希望寻找一个灰度变换函数T：s=T(r)，使得概率密度函数p_r(r)变换成希望的概率密度函数p_s(s)。灰度变换函数应该满足：

　　（1）T(r)在区间[0,1]中单调递增且单值；

　　（2）∀r∈[0,1]，有T(r)∈[0,1]。

　　满足以上条件的一个重要的直方图均衡化灰度变换函数（均匀分布的随机变量）为，（原始图像灰度r的累积分布函数（CDF））

　　根据该方程可以由原图像的各像素直接得到直方图均衡化给灰度所占百分比。例子：

　　统计原始图像的直方图：

　　8个灰度级，总计64*64=4096点。注意：离散均衡不可能拉平。

 　　仅存5个灰级，宏观拉平，微观不可能平，层次减少，对比度提高。均衡化的直方图：

　　小结：1）因为直方图是近似的概率密度函数，所以用离散灰度级进行变换时很少得到完全平坦的结果； 2） 变换后灰度级减少，即出现灰度“简并”现象，造成一些灰度层次的损失。直方图均衡化是一种非线性变换。直方图均衡的特点：增加像素灰度值的动态范围，提高图像对比度。

　　均衡化优点 能自动增强整个图像的对比度，但具体的增强效果不易控制，处理的结果是全局均衡的直方图，实际中需特定形状的直方图，从而有选择的增强某个灰度值范围内的对比度。分别对原始直方图和规定化处理后的直方图进行均衡化处理。

　　直方图均衡化的缺陷：不能用于交互方式的图象增强应用，因为直方图均衡化只能产生唯一一个结果。 恒定值直方图近似 希望通过一个指定的函数（如高斯函数）或用交互图形产生一个特定的直方图。根据这个直方图确定一个灰度级变换T(r),使由T产生的新图象的直方图符合指定的直方图。

　　

示例程序

#pragma once
#include "opencv2/highgui/highgui.hpp"
#include "opencv2/imgproc/imgproc.hpp"
#include <iostream>

using namespace std;
using namespace cv;

#define WINDOW_NAME1 "【原始图】" 

void drawHis(Mat srcImage, Mat dstImage)
{
    //将色调量化为30个等级，将饱和度量化为32个等级
    int hueBinNum = 30;//色调的直方图直条数量
    int saturationBinNum = 32;//饱和度的直方图直条数量
    int histSize[] = { hueBinNum, saturationBinNum };
    // 定义色调的变化范围为0到179
    float hueRanges[] = { 0, 180 };
    //定义饱和度的变化范围为0（黑、白、灰）到255（纯光谱颜色）
    float saturationRanges[] = { 0, 256 };
    const float* ranges[] = { hueRanges, saturationRanges };
    MatND dstHist;
    //参数准备，calcHist函数中将计算第0通道和第1通道的直方图
    int channels[] = { 0, 1 };

    //【3】正式调用calcHist，进行直方图计算
    calcHist(&dstImage,//输入的数组
        1, //数组个数为1
        channels,//通道索引
        Mat(), //不使用掩膜
        dstHist, //输出的目标直方图
        2, //需要计算的直方图的维度为2
        histSize, //存放每个维度的直方图尺寸的数组
        ranges,//每一维数值的取值范围数组
        true, // 指示直方图是否均匀的标识符，true表示均匀的直方图
        false);//累计标识符，false表示直方图在配置阶段会被清零

    //【4】为绘制直方图准备参数
    double maxValue = 0;//最大值
    minMaxLoc(dstHist, 0, &maxValue, 0, 0);//查找数组和子数组的全局最小值和最大值存入maxValue中
    int scale = 10;
    Mat histImg = Mat::zeros(saturationBinNum * scale, hueBinNum * 10, CV_8UC3);

    //【5】双层循环，进行直方图绘制
    for (int hue = 0; hue < hueBinNum; hue++)
        for (int saturation = 0; saturation < saturationBinNum; saturation++)
        {
            float binValue = dstHist.at<float>(hue, saturation);//直方图组距的值
            int intensity = cvRound(binValue * 255 / maxValue);//强度

            //正式进行绘制
            rectangle(histImg, Point(hue * scale, saturation * scale),
                Point((hue + 1) * scale - 1, (saturation + 1) * scale - 1),
                Scalar::all(intensity), FILLED);
        }

    //【6】显示效果图
    imshow("素材图", srcImage);
    imshow("H-S 直方图", histImg);
}

void drawHis1(Mat srcImage, string name)
{
    MatND dstHist;       // 在cv中用CvHistogram *hist = cvCreateHist
    int dims = 1;         // 一维
    float hranges[] = { 0, 255 };    // 灰度区间
    const float* ranges[] = { hranges };   // 这里需要为const类型
    int size = 256;  // 灰度级数
    int channels = 0; 

    calcHist(&srcImage, 1, &channels, Mat(), dstHist, dims, &size, ranges);    // cv 中是cvCalcHist
    int scale = 1;

    Mat dstImage(size * scale, size, CV_8U, Scalar(0)); // 全0矩阵

    double minValue = 0;
    double maxValue = 0;
    minMaxLoc(dstHist, &minValue, &maxValue, 0, 0);  //  在cv中用的是cvGetMinMaxHistValue

    int hpt = saturate_cast<int>(0.9 * size);
    for (int i = 0; i < 256; i++)
    {
        float binValue = dstHist.at<float>(i);           //   注意hist中是float类型    而在OpenCV1.0版中用cvQueryHistValue_1D
        int realValue = saturate_cast<int>(binValue * hpt / maxValue);
        rectangle(dstImage, Point(i * scale, size - 1), Point((i + 1) * scale - 1, size - realValue), Scalar(255));
    }
    imshow(name, dstImage);
}

void drawHisrgb(Mat srcImage)
{
    int bins = 256;
    int hist_size[] = { bins };
    float range[] = { 0, 256 };
    const float* ranges[] = { range };
    MatND redHist, grayHist, blueHist;
    int channels_r[] = { 0 };

    //【3】进行直方图的计算（红色分量部分）
    calcHist(&srcImage, 1, channels_r, Mat(), //不使用掩膜
        redHist, 1, hist_size, ranges,
        true, false);

    //【4】进行直方图的计算（绿色分量部分）
    int channels_g[] = { 1 };
    calcHist(&srcImage, 1, channels_g, Mat(), // do not use mask
        grayHist, 1, hist_size, ranges,
        true, // the histogram is uniform
        false);

    //【5】进行直方图的计算（蓝色分量部分）
    int channels_b[] = { 2 };
    calcHist(&srcImage, 1, channels_b, Mat(), // do not use mask
        blueHist, 1, hist_size, ranges,
        true, // the histogram is uniform
        false);

    //-----------------------绘制出三色直方图------------------------
    //参数准备
    double maxValue_red, maxValue_green, maxValue_blue;
    minMaxLoc(redHist, 0, &maxValue_red, 0, 0);
    minMaxLoc(grayHist, 0, &maxValue_green, 0, 0);
    minMaxLoc(blueHist, 0, &maxValue_blue, 0, 0);
    int scale = 1;
    int histHeight = 256;
    Mat histImage = Mat::zeros(histHeight, bins * 3, CV_8UC3);

    //正式开始绘制
    for (int i = 0; i < bins; i++)
    {
        //参数准备
        float binValue_red = redHist.at<float>(i);
        float binValue_green = grayHist.at<float>(i);
        float binValue_blue = blueHist.at<float>(i);
        int intensity_red = cvRound(binValue_red * histHeight / maxValue_red);  //要绘制的高度
        int intensity_green = cvRound(binValue_green * histHeight / maxValue_green);  //要绘制的高度
        int intensity_blue = cvRound(binValue_blue * histHeight / maxValue_blue);  //要绘制的高度

        //绘制红色分量的直方图
        rectangle(histImage, Point(i * scale, histHeight - 1),
            Point((i + 1) * scale - 1, histHeight - intensity_red),
            Scalar(255, 0, 0));

        //绘制绿色分量的直方图
        rectangle(histImage, Point((i + bins) * scale, histHeight - 1),
            Point((i + bins + 1) * scale - 1, histHeight - intensity_green),
            Scalar(0, 255, 0));

        //绘制蓝色分量的直方图
        rectangle(histImage, Point((i + bins * 2) * scale, histHeight - 1),
            Point((i + bins * 2 + 1) * scale - 1, histHeight - intensity_blue),
            Scalar(0, 0, 255));

    }

    //在窗口中显示出绘制好的直方图
    imshow("图像的RGB直方图", histImage);
}

void Hisbijiao()
{
    Mat srcImage_base, hsvImage_base;
    Mat srcImage_test1, hsvImage_test1;
    Mat srcImage_test2, hsvImage_test2;
    Mat hsvImage_halfDown;

    //【2】载入基准图像(srcImage_base) 和两张测试图像srcImage_test1、srcImage_test2，并显示
    srcImage_base = imread("1.jpg", 1);
    srcImage_test1 = imread("2.jpg", 1);
    srcImage_test2 = imread("3.jpg", 1);
    //显示载入的3张图像
    imshow("基准图像", srcImage_base);
    imshow("测试图像1", srcImage_test1);
    imshow("测试图像2", srcImage_test2);

    // 【3】将图像由BGR色彩空间转换到 HSV色彩空间
    cvtColor(srcImage_base, hsvImage_base, COLOR_BGR2HSV);
    cvtColor(srcImage_test1, hsvImage_test1, COLOR_BGR2HSV);
    cvtColor(srcImage_test2, hsvImage_test2, COLOR_BGR2HSV);

    //【4】创建包含基准图像下半部的半身图像(HSV格式)
    hsvImage_halfDown = hsvImage_base(Range(hsvImage_base.rows / 2, hsvImage_base.rows - 1), Range(0, hsvImage_base.cols - 1));

    //【5】初始化计算直方图需要的实参
    // 对hue通道使用30个bin,对saturatoin通道使用32个bin
    int h_bins = 50; int s_bins = 60;
    int histSize[] = { h_bins, s_bins };
    // hue的取值范围从0到256, saturation取值范围从0到180
    float h_ranges[] = { 0, 256 };
    float s_ranges[] = { 0, 180 };
    const float* ranges[] = { h_ranges, s_ranges };
    // 使用第0和第1通道
    int channels[] = { 0, 1 };

    // 【6】创建储存直方图的 MatND 类的实例:
    MatND baseHist;
    MatND halfDownHist;
    MatND testHist1;
    MatND testHist2;

    // 【7】计算基准图像，两张测试图像，半身基准图像的HSV直方图:
    calcHist(&hsvImage_base, 1, channels, Mat(), baseHist, 2, histSize, ranges, true, false);
    normalize(baseHist, baseHist, 0, 1, NORM_MINMAX, -1, Mat());

    calcHist(&hsvImage_halfDown, 1, channels, Mat(), halfDownHist, 2, histSize, ranges, true, false);
    normalize(halfDownHist, halfDownHist, 0, 1, NORM_MINMAX, -1, Mat());

    calcHist(&hsvImage_test1, 1, channels, Mat(), testHist1, 2, histSize, ranges, true, false);
    normalize(testHist1, testHist1, 0, 1, NORM_MINMAX, -1, Mat());

    calcHist(&hsvImage_test2, 1, channels, Mat(), testHist2, 2, histSize, ranges, true, false);
    normalize(testHist2, testHist2, 0, 1, NORM_MINMAX, -1, Mat());


    //【8】按顺序使用4种对比标准将基准图像的直方图与其余各直方图进行对比:
    for (int i = 0; i < 4; i++)
    {
        //进行图像直方图的对比
        int compare_method = i;
        double base_base = compareHist(baseHist, baseHist, compare_method);
        double base_half = compareHist(baseHist, halfDownHist, compare_method);
        double base_test1 = compareHist(baseHist, testHist1, compare_method);
        double base_test2 = compareHist(baseHist, testHist2, compare_method);
        //输出结果
        printf(" 方法 [%d] 的匹配结果如下：

 【基准图 - 基准图】：%f, 【基准图 - 半身图】：%f,【基准图 - 测试图1】： %f, 【基准图 - 测试图2】：%f 
-----------------------------------------------------------------
", i, base_base, base_half, base_test1, base_test2);
    }

    printf("检测结束。");
}

Mat g_srcImage; Mat g_hsvImage; Mat g_hueImage;
int g_bins = 30;//直方图组距
void on_BinChange(int, void*)
{
    //【1】参数准备
    MatND hist;
    int histSize = MAX(g_bins, 2);
    float hue_range[] = { 0, 180 };
    const float* ranges = { hue_range };

    //【2】计算直方图并归一化
    calcHist(&g_hueImage, 1, 0, Mat(), hist, 1, &histSize, &ranges, true, false);
    normalize(hist, hist, 0, 255, NORM_MINMAX, -1, Mat());

    //【3】计算反向投影
    MatND backproj;
    calcBackProject(&g_hueImage, 1, 0, hist, backproj, &ranges, 1, true);

    //【4】显示反向投影
    imshow("反向投影图", backproj);

    //【5】绘制直方图的参数准备
    int w = 400; int h = 400;
    int bin_w = cvRound((double)w / histSize);
    Mat histImg = Mat::zeros(w, h, CV_8UC3);

    //【6】绘制直方图
    for (int i = 0; i < g_bins; i++)
    {
        rectangle(histImg, Point(i * bin_w, h), Point((i + 1) * bin_w, h - cvRound(hist.at<float>(i) * h / 255.0)), Scalar(100, 123, 255), -1);
    }

    //【7】显示直方图窗口
    imshow("直方图", histImg);
}

void drawhr()
{
    //【1】读取源图像，并转换到 HSV 空间
    g_srcImage = imread("3.png", 1);
    if (!g_srcImage.data) 
    { 
        printf("读取图片错误，请确定目录下是否有imread函数指定图片存在~！ 
"); 
        // return false; 
    }
    cvtColor(g_srcImage, g_hsvImage, COLOR_BGR2HSV);

    //【2】分离 Hue 色调通道
    g_hueImage.create(g_hsvImage.size(), g_hsvImage.depth());
    int ch[] = { 0, 0 };
    mixChannels(&g_hsvImage, 1, &g_hueImage, 1, ch, 1);

    //【3】创建 Trackbar 来输入bin的数目
    namedWindow(WINDOW_NAME1, WINDOW_AUTOSIZE);
    createTrackbar("色调组距 ", WINDOW_NAME1, &g_bins, 180, on_BinChange);
    on_BinChange(0, 0);//进行一次初始化

    //【4】显示效果图
    imshow(WINDOW_NAME1, g_srcImage);

}

#include "histogram.h"

int main()
{
    Mat srcImage, dstImage;
    srcImage = imread("1.png", 1);
    if (!srcImage.data) 
    { 
        printf("读取图片错误，请确定目录下是否有imread函数指定图片存在~！ 
"); 
        return false; 
    }

    /* 二维 */
    // 转换为 HLS（HIS） 通道图像
    Mat img;
    img = imread("1.jpg");
    Mat hlsimg;
    cvtColor(img, hlsimg, COLOR_RGB2HSV);
    drawHis(img, hlsimg);

    /* rgb直方图 */
    drawHisrgb(img);

    /* 直方图比较 */
    Hisbijiao();

    /* 反向投影 */
    drawhr();

    /* 实验2补 */
    cvtColor(srcImage, srcImage, COLOR_BGR2GRAY); // 灰度化图像
    drawHis1(srcImage, "原始的灰度图的直方图"); // 绘画一维直方图
    imshow("原始图", srcImage);
    equalizeHist(srcImage, dstImage);
    imshow("经过直方图均衡化后的图", dstImage);
    drawHis1(dstImage, "均衡化后图像的直方图"); // 绘画一维直方图

    waitKey();
    destroyAllWindows();
    return 0;
}

 ---------------------------------------------------------continue---------------------------------------------------

参考文献 

[1] 刘衍琦,詹福宇,王建德.计算机视觉与深度学习实战：以MATLAB、Python为工具[M].北京：电子工业出版社.2019.11.

[2]阮秋琦,阮宇智.数字图像处理．电子工业出版社[引用日期2020-06-06].