图像处理基础01-直方图均衡化的推导和编程实现

1 简介

直方图均衡化是将图像转化为另一幅图像，转化后图像像素值的分布更接近均匀分布。原本图像的像素值（灰度值）可能集中在某一区域，这样我们看到的图像其实是比较模糊的，灰度没有层次感。直方图均衡化能够增加图像灰度值的动态范围，从而达到增强图像对比度的目的，使图像看起来更清晰。
直方图均衡属于空域图像增强，但并没有考虑图像的空间信息。但是又能得到很好的视觉效果，很有意思。
直方图相同时，图像不一定相同。一个典型例子就是图像翻转后，直方图相同，但图像不同。

2 推导

首先明确我们的目的是得到一个变换关系来对图像进行变换，具体为对所有值为r的像素点，将r改为s：
$s = T(r)$
通过变换，改变了每种灰度级的个数，变换后图像的灰度级符合均匀分布。
冈萨雷斯的书里给了几个公式，但没有推导过程，直接给出了 $T(r)$ 的公式(书中的式3.3-4)，并做了证明，证明使用(3.3-4)变换后的图像，灰度级符合均匀分布。其实根据证明过程也可以知道 $T(r)$ 是如何得出的。
那么 $T(r)$ 是如何得出的？用 $p_r(r)$ 和 $p_s(s)$ 表示r和s的概率密度(PDF)，则有以下性质(3.3-3)：
$p_s(s)ds = p_r(r)dr$
$p_s(s)ds$ 的意义是区间 $(s,s+ds)$ 上的概率， $p_r(r)dr$ 同理。显然两者是相等的。可以通过下图理解一下，下图的 $H_B(r)$ 应为 $p_r(s)$ ， $H_A(r)$ 应为 $p_s(s)$ 。

则有 $\frac{dT(r)}{dr} = \frac{ds}{dr} = \frac{p_r(r)}{p_s(s)}$
因为 $S$ 为均匀分布，即 $p_s(s) = \frac{1}{L-1}$ ， $L$ 为灰度级数(一般为256)，所以
$T(r) = {(L-1)}\int_0^{r}p_r(w)dw$
即得到了(3.3-4)。显然上式就是求r的分布函数，然后乘上 $(L-1)$ 。

对一张图像，其均衡化后的图像是唯一的，可以利用这一点来进行直方图匹配（规定化）。

3 编程实现

实际上图像中的像素值是离散随机变量，那么有
$s_k = T(r_k) = (L-1)\sum_{j=0}^{k}p_r(r_j) = (L-1)\sum_{j=0}^{k} \frac{n_j}{n} \qquad k = 0,1,2, ... ,L-1$
其中 $n$ 为像素点总数， $n_j$ 为像素值为j的像素点数量， $L$ 为灰度级数。注意， $\frac{n_j}{n}$ 得到的可能是小数，所以得到的 $s_k$ 可能为小数，灰度值为小数时也可以显示，但如果要求得到的数据为指定位数（如8位）那就进行取整。
均衡后的直方图看起来并不是均匀分布，我认为主要是因为取整。

3.1 python和opencv

hist, bins = np.histogram(img_gray.flatten(), 256, [0,256],normed=False)  
# 计算直方图 参数依次为 区间数 区间范围(不在该范围的数据会被忽略)

plt.hist(img_gray.flatten(),256,[0,256], normed=False, color = 'r')  # 画直方图
img_new = cv2.equalizeHist(img_gray) # 直方图均衡

3.2 纯C语言

直接调用函数dtype* hist_equalization(dtype *img, int img_size) 即可实现直方图均衡化。
图像处理基础的整个工程详见我的github，之后会陆续增加其他图像处理方法。

#include "core.h"
#include<stdlib.h>
/*
typedef int num_t;
typedef char dtype;
*/


num_t* histogram(dtype *img, int img_size)
{
        int gray_level = 256;
        num_t* hist = (num_t*)malloc(sizeof(num_t)*(gray_level));  // each gray level may have great number
        for(int i=0;i<gray_level;i++)
        {
                for(int j=0;j<img_size;j++)
                {
                        if(img[j]==i)
                        {
                                hist[i] +=1;
                        }
                }
        }
        return hist;
}


dtype* hist_equalization(dtype *img, int img_size)
{       
        int gray_level = 256;
        num_t* hist = histogram(img, img_size);
        num_t* hist_acc = (num_t*)malloc(sizeof(num_t)*gray_level);  // Accumulate histogram data
        dtype* new_img = (dtype*)malloc(sizeof(dtype)*img_size);
        
        for(int i=0;i<gray_level;i++)
        {       
                if(i==0)
                {       
                        hist_acc[i] = hist[i];
                }
                else    
                {       
                        hist_acc[i] = hist_acc[i-1]+hist[i];   
                        hist_acc[i-1] = hist_acc[i-1]*(gray_level-1)/img_size;          // don't need hist_acc[i-1] to store accumulate number any more
                }
        }
        //hist_acc[gray_level-1] =  hist_acc[gray_level-1]*(gray_level-1)/img_size;    //actually hist_acc[gray_level-1] equals to img_size before reassignment
        hist_acc[gray_level-1] = gray_level-1;
        
        for(int i=0;i<img_size;i++)
        {       
                int pixel_value = img[i];
                new_img[i] = hist_acc[pixel_value];
        }
        
        free(hist_acc);
        free(hist);
        return new_img;
}

参考文献：
[1] https://wenku.baidu.com/view/6c55ecbdc77da26925c5b090
[2] 数字图像处理，冈萨雷斯
[3] 数字图像处理学习指导与题解，陈青

图像处理基础01-直方图均衡化的推导和编程实现

1 简介

2 推导

3 编程实现

3.1 python和opencv

3.2 纯C语言

相关推荐