往期文章回顾:
计算机视觉入门之<零>
计算机视觉入门之图像处理<一>:图像处理基础概念
计算机视觉入门之图像处理<二>:图像处理基础概念
计算机视觉入门之图像处理<三>:图像插值方法
前言
掌握了图像处理基础知识,该系列文章进一步了解图像处理基本方法,本篇文章主要包括图像直方图理论、直方图均衡化原理以及代码实现等方面。
灰度直方图基本概念
灰度范围为
[
0
,
L
−
1
]
\left[ {0,L - 1} \right]
[0,L−1]的数字图像中,将其中每一个灰度级看成一个随机变量,则随机变量的取值情况就表示了该图像的统计特性,改用直方图来描述该特性,则称之为灰度直方图。
灰度直方图是灰度级的离散化函数,表示为:
h
(
r
k
)
=
n
k
h({r_k}) = {n_k}
h(rk)=nk其中
r
k
r_k
rk表示第
k
k
k级灰度值,
n
k
n_k
nk表示像素为
r
k
r_k
rk的像素个数。
灰度直方图的归一化:
p
(
r
k
)
=
n
k
M
N
p({r_k}) = \frac {n_k} {MN}
p(rk)=MNnk式中
M
N
MN
MN表示图像的行数和列数的乘积,归一化后表示的意义是灰度值
r
k
r_k
rk在图像中出现概率的估计。归一化后的直方图所有分量满足下列关系式:
∑
k
=
0
L
−
1
p
(
r
k
)
=
1
\sum_{k=0}^{L-1}p(r_k)=1
k=0∑L−1p(rk)=1
如图所示的四个基本灰度级的特征花粉图像,分别表现为暗、亮、低对比度和高对比度,灰度图横轴表示灰度值
r
k
r_k
rk, 纵轴表示
n
k
n_k
nk或归一化后的值
p
(
r
k
)
p(r_k)
p(rk)。由图可以得出,直方图越靠近灰度值小的一端,图像越暗,反之;直方图集中在灰度级中部的图像具有较低的对比度,而高对比度的图像直方图几乎覆盖了所有灰度级,而且各个灰度值的分布较为均匀。由上可知:一幅图像具有较广的灰度级并且分布均匀,那么该图像具有较高的对比度和较大的灰度变化范围,这种图像可以提供更多的图像细节,更有利于进一步的图像处理。
直方图均衡化理论
![](https://i-blog.csdnimg.cn/blog_migrate/43922a0192a000ab8268f99ef8fab1ab.png#pic_center)
图像直方图均衡化:均表示平均,即在图像中平均到每个灰度级,衡表示平衡、一样多,均衡表明直方图中每个灰度级上的直方图值一样大,直方图一样高,上图很形象地描述了均衡化,当然实际均衡化的效果不可能这样。
假设原始图像的灰度为
r
r
r,目标灰度为
s
s
s,因为灰度值是一个实数,则从数学的角度来看,实数的变化可以通过一个映射函数来实现,即:
s
=
T
(
r
)
,
0
≤
r
≤
L
−
1
s=T(r),0\leq r\leq L-1
s=T(r),0≤r≤L−1该灰度变换函数应当满足以下条件:
- 0 ≤ s = T ( r ) ≤ L − 1 0\leq s=T(r)\leq L-1 0≤s=T(r)≤L−1. 因为灰度变换后还是灰度,所以必须满足灰度值范围。
- T ( r ) T(r) T(r)在 0 ≤ r ≤ L − 1 0\leq r\leq L-1 0≤r≤L−1上是一个严格单调递增函数。①在函数映射变化之后至少要求像素点的亮暗关系不变,比如说, r 1 r_1 r1=30, r 2 r_2 r2=20( r 1 > r 2 r_1>r_2 r1>r2)在灰度图中灰度 r 1 r_1 r1对应的像素点比灰度 r 2 r_2 r2对应的点亮,那么变换之后还应当满足 T ( r 1 ) T(r_1) T(r1)对应的像素点比 T ( r 2 ) T(r_2) T(r2)对应的像素点亮,即 T ( r 1 ) > T ( r 2 ) T(r_1)>T(r_2) T(r1)>T(r2),因此函数 T ( r ) T(r) T(r)应该是增函数;②为什么需要严格的单调?在后面计算中需要用到 r = T − 1 ( s ) r=T^{-1}(s) r=T−1(s),该公式表示的是函数 T ( r ) T(r) T(r)的反函数,计算反函数的前提就是要求 T ( r ) T(r) T(r)单调,不单调则不存在反函数,同时也是保证 s s s到 r r r的反映射是一对一的(不严格单调则反映射存在一对多的情况)因此 T ( r ) T(r) T(r)必须是单调函数。
![](https://i-blog.csdnimg.cn/blog_migrate/2b6002e5677b0ebc52c94a68ef877154.png#pic_center)
上图左边为概率密度函数
p
r
(
r
)
p_r(r)
pr(r),在图像均衡化后最理想的效果就是图像的每个灰度级上的值一样大,再加上要求
p
s
(
s
)
p_s(s)
ps(s)在
[
0
,
L
−
1
]
[0, L-1]
[0,L−1]上的积分为1,所以最直观的
p
s
(
s
)
=
1
L
−
1
p_s(s)=\frac{1}{L-1}
ps(s)=L−11
因为已知
r
r
r和
s
s
s存在的映射关系
T
(
r
)
T(r)
T(r)和
p
r
(
r
)
p_r(r)
pr(r),根据概率学知识就可以使用
p
r
(
r
)
p_r(r)
pr(r)和
T
(
r
)
T(r)
T(r)来并表示
p
s
(
s
)
p_s(s)
ps(s),即:
p
s
(
s
)
=
p
r
(
r
)
∣
d
r
d
s
∣
(1)
p_s(s) = p_r(r) \lvert \frac {dr} {ds} \rvert \tag{1}
ps(s)=pr(r)∣dsdr∣(1)上式中
r
r
r为
s
s
s的函数,即
r
=
T
−
1
(
s
)
r=T^{-1}(s)
r=T−1(s)
将
p
s
(
s
)
=
1
L
−
1
p_s(s)=\frac{1}{L-1}
ps(s)=L−11代入式(1)可得灰度变换函数
T
(
r
)
T(r)
T(r):
1
L
−
1
=
p
r
(
r
)
∣
d
r
d
s
∣
\frac{1}{L-1} = p_r(r) \lvert \frac {dr} {ds} \rvert
L−11=pr(r)∣dsdr∣
1
L
−
1
d
s
=
p
r
(
r
)
d
r
\frac{1}{L-1}ds =p_r(r)dr
L−11ds=pr(r)dr
∫
0
s
1
L
−
1
d
s
=
∫
0
r
p
r
(
r
)
d
r
\int_0^s\frac{1}{L-1}ds =\int_0^rp_r(r)dr
∫0sL−11ds=∫0rpr(r)dr
1
L
−
1
s
=
∫
0
r
p
r
(
r
)
d
r
\frac{1}{L-1}s =\int_0^rp_r(r)dr
L−11s=∫0rpr(r)dr
s
=
(
L
−
1
)
∫
0
r
p
r
(
r
)
d
r
s=(L-1)\int_0^rp_r(r)dr
s=(L−1)∫0rpr(r)dr
T
(
r
)
=
(
L
−
1
)
∫
0
r
p
r
(
r
)
d
r
T(r)=(L-1)\int_0^rp_r(r)dr
T(r)=(L−1)∫0rpr(r)dr
T
(
r
)
=
(
L
−
1
)
∫
0
r
p
r
(
w
)
d
w
(2)
T(r)=(L-1)\int_0^rp_r(w)dw\tag{2}
T(r)=(L−1)∫0rpr(w)dw(2)
验证构造函数是否符合条件:
①
∫
0
r
p
r
(
w
)
d
w
\int_0^r p_r(w)dw
∫0rpr(w)dw表示的是关于随机变量
r
r
r的累计分布函数,其中
p
r
(
w
)
p_r(w)
pr(w)D大于0,因此该积分项随着随机变量
r
r
r的增大而增大,因此符合灰度变换函数的单调递增条件;
②当
r
r
r=
L
−
1
L-1
L−1时,
T
(
r
)
T(r)
T(r)增加到最大,为
L
−
1
L-1
L−1,即
∫
0
L
−
1
p
r
(
w
)
d
w
\int_0^{L-1} p_r(w)dw
∫0L−1pr(w)dw表示灰度直方图中灰度级0~
L
−
1
L-1
L−1各自出现概率之和,为1,符合灰度变换函数的灰度范围条件。
由上述分析可知,
p
s
(
s
)
=
1
L
−
1
p_s(s)=\frac{1}{L-1}
ps(s)=L−11是完全可行的概率密度函数,构造出来的灰度变换函数符合条件。由
p
s
(
s
)
=
1
L
−
1
p_s(s)=\frac{1}{L-1}
ps(s)=L−11可知,它是一个均匀PDF,与
p
r
(
r
)
p_r(r)
pr(r)的具体形式无关。
因为将数字图像的的各个灰度级视为一个离散随机变量,因此使用各个灰度级出现的概率(归一化后直方图的值)与求和来替代概率密度函数与积分,则等式(2)改写为:
s
k
=
T
(
r
k
)
=
(
L
−
1
)
∑
j
=
0
k
p
r
(
r
j
)
=
(
L
−
1
)
M
N
∑
j
=
0
k
n
j
,
k
=
0
,
1
,
2
,
⋯
,
L
−
1
s_k=T(r_k)=(L-1)\sum_{j=0}^kp_r(r_j)= \frac{(L-1)}{MN}\sum_{j=0}^kn_j, k=0, 1, 2,\cdots,L-1
sk=T(rk)=(L−1)j=0∑kpr(rj)=MN(L−1)j=0∑knj,k=0,1,2,⋯,L−1
等式中
∑
j
=
0
k
p
r
(
r
j
)
\sum_{j=0}^kp_r(r_j)
∑j=0kpr(rj)表示累计直方图值
图像直方图均衡化步骤:
①计算灰度值
r
k
r_k
rk出现的概率;
②计算累计直方图值,即
∑
j
=
0
k
p
r
(
r
j
)
\sum_{j=0}^kp_r(r_j)
∑j=0kpr(rj);
③将累计直方图值代入公式计算;
④将计算结果取整;
⑤将新的灰度值放入与原灰度值对应的位置,即可生成均衡化后的图像。
实际例子计算:
如上表格所示,该图像为64*64的3比特数字图像,第一列为灰度级
r
k
r_k
rk,第二列为对应的灰度值个数
n
k
n_k
nk,第三列为直方图值
p
r
(
r
k
)
p_r(r_k)
pr(rk),第四列为累计直方图值(归一化后的直方图和
∑
j
L
−
1
p
r
(
r
k
)
\sum_j^{L-1} p_r(r_k)
∑jL−1pr(rk)),第五列为均衡化后的灰度值
s
k
s_k
sk。第六列为取整后的灰度值
直方图均衡化实现
直方图生成
- 生成灰度直方图
src = cv.imread('lenna.png', 0)#读取原始图像
src = cv.resize(src, (300, 300))#缩小图像尺寸为300*300
array = np.zeros([256])#新建空矩阵用于存放灰度值个数,初始化为0,256表示灰度级为256,
for i in range(src.shape[0]):
for j in range(src.shape[1]):
array[src[i, j]] += 1#遍历整张图像,若灰度值相同则就在矩阵中对应位置加1
# array = array /(300*300)#将数据矩阵中的每一个值都除以300*300,进行归一化处理
plt.title('histogram')
plt.plot(array, c='r')#用实线绘制
plt.hist(src.ravel(), 256)#利用函数绘制灰度图
plt.show()
- 实验结果:
![](https://i-blog.csdnimg.cn/blog_migrate/014418e11b28adf209a8f45fef7d200a.png#pic_center)
图中红色部分(该部分使用实线绘制,故表现为连续,实际应是离散值)为自编代码实现,而蓝色部分为cv内部函数实现.
- 彩色直方图
src = cv.imread('lenna.png', 1)
array = np.zeros([256, 3])
for k in range(src.shape[2]):
for i in range(src.shape[0]):
for j in range(src.shape[1]):
array[src[i, j, k],k] += 1
hist_b1 = cv.calcHist([src], [0], None, [256], [0, 255])#使用cv内部函数绘制B通道灰度图
hist_g1 = cv.calcHist([src], [1], None, [256], [0, 255])#使用cv内部函数绘制G通道灰度图
hist_r1 = cv.calcHist([src], [2], None, [256], [0, 255])#使用cv内部函数绘制R通道灰度图
plt.plot(hist_b1, c='b')
plt.plot(hist_g1, c='g')
plt.plot(hist_r1, c='r')
plt.plot(array[:, 0], 'y*')
plt.plot(array[:, 1], 'c*')
plt.plot(array[:, 2], 'k*')
plt.show()
效果对比:
(同样,实线部分并非真正直方图,只是为了更形象对比两种实现效果)
直方图均衡化
灰度直方图均衡化
def histogram_equalization_fnuc(src):
src_array = np.zeros([256, 2])#建立零矩阵存储灰度值个数
for i in range(src.shape[0]):
for j in range(src.shape[1]):
src_array[src[i, j], 0] = src[i, j]#矩阵的第一列用于存储灰度值
src_array[src[i, j], 1] += 1#矩阵的第二列用于存储灰度值对应的数目
src_array[:, 1] /= (300 * 300)#直方图归一化处理,则矩阵的第二列用于存储灰度值出现的概率
accum = 0
result_img = np.zeros_like(src)#建立一个矩阵用于存储均衡化后的数据
for m in range(256):
accum += src_array[m, 1]#计算累加直方图值(概率)
temp = np.round(accum * 255)#将计算后的像素值取整(灰度图均衡化关键一步)
for i in range(src.shape[0]):
for j in range(src.shape[1]):
if src[i, j] == src_array[m, 0]:#根据之前存储的灰度值,遍历原始图像,判断图像中与当前灰度相同的位置,
# 若相同则在新的图像位置存储均衡化后的灰度值
result_img[i, j] = temp
return result_img
主函数:
def main():
src = cv.imread('lenna.png', 0)#读取灰度图像
src = cv.resize(src, (300, 300))#缩小尺寸
plt.hist(src.ravel(), 256, color='b')#查看原灰度图直方图
result = histogram_equalization_fnuc(src)#对原灰度图进行均衡化
plt.hist(result.ravel(), 256, color='r')#查看原灰度图均衡化后的直方图
plt.show()
cv.imshow("not_func", np.hstack([src, result]))#对比均衡化前后图像效果
cv.waitKey(0)
实现效果:
![]() |
![]() |
上图中蓝色部分为原灰度图的直方图,红色部分为均衡化后的直方图,可以明显看出,均衡化之后图像的灰度级增多,几乎覆盖了8比特图全部灰度级,提高了图像对比度。
彩色直方图均衡化
def main():
b, g, r = cv.split(src)#彩色图像均衡化需要分通道进行
b = histogram_equalization_fnuc(b)
g = histogram_equalization_fnuc(g)
r = histogram_equalization_fnuc(r)
result = cv.merge([b, g, r])#最后将均衡化后的三通道进行合并
cv.imshow("Histogram Equalization", np.hstack([src, result]))#均衡化前后对比,
cv.waitKey(0)
cv.destroyAllWindows()
if __name__ == '__main__':
main()
实现效果:
因为每个通道的灰度图都发生了均衡化,所以三个通道组合的彩色图颜色也发生了变化。
参考文献:[1]数字图像处理:第3版/(美)冈萨雷斯(Gonzalez,R.C),(美)伍兹(Woods, R.E.)著.阮秋琦等译
[2]直方图匹配的数学原理
[Until:2021年1月15日17时,转眼又是半个月……]