直方图均衡

1 直方图均衡

如图3.16所示，是4个基本灰度级为特征的花粉图像：暗图像、亮图像、低对比度图像和高对比度图像，右侧显示了与这些图像对应的直方图。

由四张图的对比我们可知，如一幅图像的像素倾向于占据整个可能的灰度级并分布均匀，则该图像会有高对比度的外观并展示灰色调的变换，最终效果将会是一幅灰度细节丰富且动态范围较大的图像。图像均衡化即是实现这种效果。

2 直方图均衡化推导过程

连续灰度值的情况下，一幅图像的灰度级可以看成【0,L-1】内的随机变量。随机变量的基本描绘子是其概率密度函数（PDF）。r表示输入图像的灰度值，s表示均衡化之后的图像的灰度值，pr(r)和ps(s)表示随机变量r和s的PDF，因为r为输入图像的像素，所有pr(r)可以求得，则变换之后的s的PDF：

令

由上式可以看出s的PDF始终是均匀的。

对于离散值的灰度值，我们处理其概率与求和来代替PDF和积分。一幅图像中rk出现的概率近似为：

M表示行 N表示列 MN表示图像中像素的总数，nk表示灰度为rk的个数

3 直方图均衡使用举例

其他与s0，s1求法相同，把他们近似为最接近的整数。

这些是均衡化之后的直方图的值，r0映射为1，r1映射为3，r2映射为5，r3映射为6，r4映射为6，r5，r6,r7映射为7，根据r中对应的像素个数，

1对应的像素个数为790,概率为0.19，

3对应的像素个数为1023，对应的概率为0.25,

5对应的像素个数为850,对应的概率为0.21，

6对应的像素个数为985，对应的概率为0.24，

7对应的像素个数为448，对应的概率为0.11

下面第一张是原始的图像的直方图，第二张是均衡化之后的直方图

4 C++实现

1. 求输入图像中每一个灰度级别的元素个数
2. 根据每一个灰度级别的个数求该灰度级别所占的概率
3. 根据下式求新的灰度级别

 可以简化为（L-1）与前k个概率之和的乘积

    4.找到r和s的对应关系

    5.根据对应关系为图像赋新的灰度值

#include <iostream>
#include<opencv2/opencv.hpp>

int main() {
    cv::Mat src=cv::imread("../1.jpg",0);
    int gray[256]={0};//记录灰度级别下的像素个数
    int mn=src.cols*src.rows;//像素总数
    double r_prob[256]={0};//输入图像的概率
    double pdf[256]={0};//累积概率密度
    double r[256]={0};//r与s的映射
    double s_prob[256]={0};//均衡化之后的概率
    //统计每个灰度级别下的像素个数
    cv::Mat dst=src.clone();
    for(int row=0;row<src.rows;row++){
        for(int col=0;col<src.cols;col++){
            int g=src.at<uchar>(row,col);
            gray[g]++;
        }
    }

    //计算每一个像素级别的概率
    for(int i=0;i<256;i++){
        r_prob[i]=((double)gray[i])/mn;
    }

    pdf[0]=r_prob[0];
    for(int i=1;i<256;i++){
        pdf[i]=pdf[i-1]+r_prob[i];
    }

    //每个输入的灰度级别对应的输出的灰度级别
    for(int i=0;i<256;i++)
    {
        int c=cvRound(255*(pdf[i]));
        r[i]=c;
    }

    for(int row=0;row<src.rows;row++){
        for(int col=0;col<src.cols;col++){
            int g=src.at<uchar>(row,col);
            dst.at<uchar>(row,col)=r[g];//重新赋值新的灰度级别
        }
    }

    cv::imshow("src",src);
    cv::imshow("dst",dst);

    int bins=256;
    int hist_size[]={bins};
    float range[]={0,256};
    const float* ranges[]={range};
    cv::Mat histsrc,histdst;
    int channels[]={0};


    cv::calcHist(&src,1,channels,cv::Mat(),histsrc,1,hist_size,ranges,true,false);
    cv::calcHist(&dst,1,channels,cv::Mat(),histdst,1,hist_size,ranges);

    std::cout<<histsrc.size()<<std::endl;

    int scale=1;
    cv::Mat srcImage(256*scale,256,CV_8U,cv::Scalar(0));
    cv::Mat dstImage(256*scale,256,CV_8U,cv::Scalar(0));

    double maxValue=0;
    double minValue=0;

    cv::minMaxLoc(histsrc,&minValue,&maxValue,0,0);
    int hpt=cv::saturate_cast<int>(0.9*256);
    for(int i=0;i<256;i++){
        float binValue=histsrc.at<float>(i);
        int realValue=cv::saturate_cast<int>(binValue*hpt/maxValue);
        cv::rectangle(srcImage,cv::Point(i*scale,255),cv::Point((i+1)*scale-1,256-realValue),
                      cv::Scalar(255));
    }

    cv::minMaxLoc(histdst,&minValue,&maxValue,0,0);
    hpt=cv::saturate_cast<int>(0.9*256);
    for(int i=0;i<256;i++){
        float binValue=histdst.at<float>(i);
        int realValue=cv::saturate_cast<int>(binValue*hpt/maxValue);
        cv::rectangle(dstImage,cv::Point(i*scale,255),cv::Point((i+1)*scale-1,256-realValue),
                      cv::Scalar(255));
    }

    cv::imshow("srcImage",srcImage);
    cv::imshow("dstImage",dstImage);
    cv::waitKey(0);
    return 0;
}