2021-10-10

直方图均衡化的python实现

直方图均衡化是图像处理中非常经典的算法，作为图像处理的小白，我对算法进行了python语言的实现。代码中使用了大量的for循环，效率不是很高，但是笔者作为编程小白，还没有能力优化，所以请大家明白原理之后自己动手。

数学原理

其实本来不想讲数学原理，但是对于理解程序还是有帮助的，我是按照数学的推导步骤一步一步向下进行的。

首先考虑连续化灰度值，连续化灰度值便于数学表达式的推导

r 表示原图像的灰度值，其值域为 [0,L-1] r = 0 时像素为黑色，r = L-1 时像素为白色。
s 表示均衡化后图像的灰度值，在r 与 s 之间存在某种映射关系，将此种关系表达为T，有下式：

s=T( r ),0≪r≪L-1

这里有两个要求：

1): T( r ) 在区间 0≪r≪L-1 上为单调递增函数。
2): 0≪T( r )≪L-1。

条件1是为了保证输出灰度值不低于相应的输入值，防止灰度反变换时产生人为缺陷（个人理解：是为了保证不会出现图像灰度反转的情况，例如：映射是单调减函数，低灰度级会映射出高灰度级，这样将会造成部分或者整体灰度的反转，造成不良影响）。

条件2是为了保证输出灰度的范围与输入灰度的范围相同（个人理解:保证输出的图像位数是等于输入的图像位数，例如：256灰度级的图像是8位的图像，输入也应该是256灰度级8位图像）。

一幅图像的灰度级可视为区间固定在 [0,L-1] 的一个随机变量（连续性），对连续性的随机变量我们可以用概率密度函数来描述其数学表达，令 p_r ( r ) 表示随机变量 r 的概率密度函数，p_s (s) 同理。

p_r ( r ) 是原图像的概率密度函数（我们已经知道它的数学表达式），如果我们再知晓 T® 的表达形式，那么可以通过概率论里的：《连续性随机变量的函数的概率密度函数》的内容来得到 p_s (s) 因此我们有：p_s (s)=p_r ( r )*dr/ds 。
对于 T( r ) 我们有重要的变换函数（课本上给的，我也没找见出处）：

这个变换满足条件1与条件2，所以我们可以用于图像的处理。

最终课本给出 p_s (s) 的概率密度函数（公式太长，有兴趣参考冈萨雷斯的数字图像处理课本：p_s (s)=1/(L-1)，0≪r≪L-1
这说明：均衡化之后的随机变量s 由一个均匀PDF表征，它与 p_r ® 的表示形式无关。（这一点可以用在直方图规定化之中）。

再考虑离散值，因为实际的应用我们对图像的处理都是离散的数字图像

对于离散值，我们处理其概率（直方图值）与求和来替代处理概率密度函数与积分。
一幅数字图像中灰度级 r_k出现的概率近似为：p_r (r_k )=n_k/MN ，k=0，1，2，3， …….，L-1 MN是图像的总像素，n_k 是灰度级 r_k 的像素的共个数，L是图像中可能的灰度级的数量(8位图像就是256)
最终：

这个公式是指导我们的工程实现的算法。

技术路线

使用的软件为vs code 。
使用的python库有： CV2，math，matplotlib.pyplot，numpy。不清楚的可以百度一下这些库是做什么的。

代码示例

import numpy as np
import matplotlib.pyplot as plt
import cv2
import math

# 定义函数，功能：展示图像灰度直方图
def Showhist(imgname):
    plt.hist(imgname.ravel(), bins=255, rwidth=0.8, range=(0, 255),density=True)     #设定坐标值域
    plt.xlabel('Gray value')                                                         #设定x轴名字
    plt.ylabel('number')                                                             #设定y轴名字
    plt.title(r'灰度图')                                                             #设定直方图名字
    plt.show()
# 定义函数，功能：展示图像
def Showimg(imgname):
    cv2.imshow('intput', imgname)                                                    #展示图片-imgname，窗口名为intput
    # 设定窗口等待任意按键之后被关闭，不然窗口自动关闭根本看不清
    cv2.waitKey(0)
    cv2.destroyAllWindows()

# 读取照片
img = cv2.imread(r"C:\Users\li_lei\Desktop\picture\lena.jpg", cv2.IMREAD_GRAYSCALE)  #图片地址不能含中文

# 展示原图像
Showimg(img)

# 展示原图像灰度的直方图
Showhist(img)

img1=img.copy()                         #预设的输出图片
# 生成列表用来存储相关数据因为这个demo处理的是8位灰度图像，所以值都是256
hist1 = np.zeros(256, dtype=np.float32) #直方图统计结果-频次
Ph = np.zeros(256, dtype=np.float32)    #灰度频率
s = np.zeros(256, dtype=np.float32)     #均衡化之后的频率 
S = np.zeros(256, dtype=np.float32)     #均衡化之后的灰度级是浮点数
ABS = np.zeros(256, dtype=np.int)       #四舍五入之后的灰度级 
value=255                               #value的值取决于图像的位数，8位就是（256-1）
# 直方图统计
row, col = np.shape(img)
    # 灰度频次计算
for i in range(row):
        for j in range(col):
            hist1[img[i][j]] += 1
    # 灰度频率计算
for i in range(0,256):  
    Ph[i] = hist1[i]/(row*col) 
    # 均衡化之后的频率 
for i in range(1,256):
    s[0] = Ph[0]
    s[i] = (s[i-1]+Ph[i])
    # 均衡化之后的灰度级
for i in range(0,256):
    S[i]=value*s[i]
    # 四舍五入算法  对 ‘均衡化之后的灰度级’ 取整
for i in range(0,256):
    ABS[i]=math.ceil(S[i]-0.5)
    # 进行对文件的灰度值重写，这个功能比较重要 
for i in range(row):
    for j in range(col):
        img1[i,j]=ABS[img[i,j]]        #原图像某位置的灰度值经ABS[]的映射为新的灰度值，这个新的灰度值就是均衡化后的新图像的灰度值
                                       #本质是一个映射列表，原图像的灰度值是输入，是ABS[]的下标，ABS[]内存储的是‘均衡化之后的灰度级’，对于每个下标，都给出相应的值
Showimg(img1)
Showhist(img1)
# 保存图片
cv2.imwrite('ec_lena.jpg', img1)

我们使用lean图像

直方图为：

均衡化之后：

直方图为：

可以发现，虽然直方图不是完全的均衡，但是其效果已经展示出来了，笔者尝试对第二章lena图进行再次均衡化，

直方图为：

发现虽然直方图上更加均衡，但是效果不是特别明显。个人感觉进行二次均衡化的意义不是很大。

文章内照片素材来自网络，若有侵权，联系作者删除，对此给您造成的困扰道歉。
文章内代码内容为原创，仅提供给学习使用！文章写的比较仓促，如有错误，望诸位海涵，如能斧正