数字图像处理（三）直方图均衡

发表于 更新于 分类于 阅读次数： Waline： 本文字数： 6.9k 阅读时长 ≈ 6 分钟

在数字图像处理中，直方图均衡是调整图像亮度，对比度，图像增晰等操作中常用的做法。对于图像中各个不同的颜色进行直方图统计，采取统计数据对于颜色进行重新映射，从而达到调整对比度，图像增晰的目的。

1.数字图像的直方图

直方图是指对于整幅图像所有的像素点颜色进行统计，从而得到不同颜色的统计数据，将该数据以直方图的形式画出，即为图像的直方图。

例如，对于下面的这两幅图像：

其所对应的直方图为：

由图可见，高对比度的图像颜色分布较广，低对比度图像颜色集中分布在某一区域。第二张图像中的像素值集中分布在一个颜色上，造成图像模糊。

使用图像直方图可以清晰地看出图像的颜色分布情况，从而对于如何进行图形的增晰，提高对比度有了下一步的指导作用

2.直方图均衡

直方图均衡，顾名思义，就是在图像直方图上进行图像颜色的均衡。具体来说就是将之前分布过于集中的颜色在直方图上分散开来，让其映射到应当占据的像素值上。具体方法如下(以3位8色为例)：

假设一张64*64像素(M=64,NJ=64,MN=4096)，3位(L=8)的图像来说，其颜色空间为[0,L-1]，即[0,7]。对整张图像的像素值进行统计，得到统计数据如下：

\(r_k\) \(n_k\) \(p_r(r_k)=n_k/MN\)
\(r_0=0\) 790 0.19
\(r_1=1\) 1023 0.25
\(r_2=2\) 850 0.21
\(r_3=3\) 656 0.16
\(r_4=4\) 329 0.08
\(r_5=5\) 245 0.06
\(r_6=6\) 122 0.03
\(r_7=7\) 81 0.02

直方图均衡变换公式为：\(s_k=T(r_k)=(L-1)\sum_{j=0}^kp_r(r_j)=\frac{L-1}{MN}\sum_{j=0}^{k}n_j\)

对于该直方图：\(s_0=T(r_0)=7*\sum_{j=0}^0p_r(r_j)=7p_r(r_0)=1.33\)

类似的，有：\(s_1=T(r_1)=7\sum_{j=0}^1p_r(r_j)=7p_r(r_0)+7p_r(r_1)=3.08\)

以及\(s_1=4.55,s_3=5.67,s_4=6.23,s_5=6.65,s_6=6.86,s_7=7.00\)

将所有的\(s\)值近似为最接近的整数(四舍五入) \[ s_0=1.33\rightarrow1~~~~~~~s_1=3.08\rightarrow3~~~~~~~s_2=4.55\rightarrow5~~~~s_3=5.67\rightarrow6 \]

\[ s_4=6.23\rightarrow6~~~~~~~s_5=6.65\rightarrow7~~~~~~~s_6=6.86\rightarrow7~~~~s_7=7.00\rightarrow7 \]

将原有图像通过\(s_k=T(r_k)\)映射到新图像中，得到的结果即为直方图均衡的结果。

以上图中的长城为例，对于第二张图进行直方图均衡，得到的结果与原图对比如下：

原图像整体偏暗，经过直方图均衡之后亮度提升,对比度提升，长城的细节纹理更加清晰。

3.目标直方图均衡

以一幅图像的直方图为标准，对待均衡的图像按照该标准进行直方图均衡。方法如下：

首先对于待配准图像计算其均衡变换：\(s_k=T(r_k)=(L-1)\sum_{j=0}^kp_r(r_j)~~~~k=0,1,2...L-1\)

计算出标准图像的直方图均衡变换：\(G(z_p)=(L-1)\sum_{i=0}^qp_z(z_i)\)

对于每一个待配准图像颜色\(r_k\)，有\(s_k\)与其对应；对于该\(s_k\)，找到与其最近的\(G(z_q)\)，再将其反变换到颜色\(z_q\)，可以建立从\(r_k\rightarrow z_q\)的映射，从而将待配准的图像以标准图像为模板进行直方图均衡。

还是以上面的长城图像为例，以原图为目标进行直方图均衡，得到结果对比如下：

用来匹配的原图为：

可以看出，相比于上面的直方图均衡，目标直方图均衡在整体的色彩上更加偏向于用于做为目标的标准图，观感上与原图更加接近。

4.局部直方图均衡

顾名思义，局部直方图是在图像的每一个小图像块例如\(3*3,5*5 7*7\)大小的局部进行直方图均衡的操作，算法部分类似于直方图均衡，在此不再赘述，得到的效果如下：

附录

参考文献：

[1]数字图像处理[M]：第三版/（美）拉斐尔·C·冈萨雷斯（Rafael C.Gonzalez），（美）理查德·E·伍兹（Richard E. Woods）著；阮秋琦等译，—北京：电子工业出版社，2017.5

源码：

1.直方图统计：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
import cv2 as cv
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import mlab
from matplotlib import rcParams
def graybar(img,win_num,img_name,flag):
    num=np.zeros((256))
    img_shape=img.shape
    if flag:
        for i in range(img_shape[0]):
            for j in range(img_shape[1]):
                num[img[i][j][0]]+=1
    else:
        for i in range(img_shape[0]):
            for j in range(img_shape[1]):
                if img[i][j][0]==0:
                    img[i][j][0]=255
                num[img[i][j][0]]+=1
    fig=plt.figure(win_num)
    rects=plt.bar(range(256),num,0.2)
    plt.title('gray_bar_'+img_name)
    plt.show()
def flag(str):
    if '1' in str:
        return 0
    return 1
filename=['citywall.bmp','citywall1.bmp','citywall2.bmp','elain.bmp','elain1.bmp','elain2.bmp','elain3.bmp','lena.bmp','lena1.bmp','lena2.bmp','lena3.bmp','woman.BMP','woman1.bmp','woman2.bmp']
for i in filename:
    img=cv.imread(i)
    graybar(img,0,i,flag(i))

2.直方图均衡：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
import cv2 as cv
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import mlab
from matplotlib import rcParams
def graybar_trans(img,flag):
    num=np.zeros((256))
    img_shape=img.shape
    if flag:
        for i in range(img_shape[0]):
            for j in range(img_shape[1]):
                num[img[i][j][0]]+=1
    else:
        for i in range(img_shape[0]):
            for j in range(img_shape[1]):
                if img[i][j][0]==0:
                    img[i][j][0]=255
                num[img[i][j][0]]+=1
    fre=num/np.sum(num)
    trans=np.zeros((256))
    trans[0]=fre[0]*256
    for i in range(256):
        if i:
            trans[i]=trans[i-1]+fre[i]*255
    trans=trans.astype(np.uint8)
    new_img=np.zeros((img_shape[0],img_shape[1]))
    for i in range(img_shape[0]):
        for j in range(img_shape[1]):
            new_img[i][j]=trans[img[i][j][0]]
    return new_img.astype(np.uint8)
def flag(str):
    if '1' in str:
        return 0
    return 1
filename=['citywall.bmp','citywall1.bmp','citywall2.bmp','elain.bmp','elain1.bmp','elain2.bmp','elain3.bmp','lena.bmp','lena1.bmp','lena2.bmp','lena3.bmp','woman.BMP','woman1.bmp','woman2.bmp']
for i in filename:
    img=cv.imread(i)
    result=graybar_trans(img,flag(i))
    # cv.namedWindow(i,cv.WINDOW_FREERATIO)
    # cv.imshow(i,result)
    # cv.waitKey(0)
    cv.imwrite('bar_'+i,result)

3.目标直方图均衡：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
import cv2 as cv
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import mlab
from matplotlib import rcParams
def find(a,x):
    shape_tmp=a.shape
    for i in range(shape_tmp[0]):
        if x-a[i]<=0:
            return i
    return shape_tmp[0]-1
def graybar_trans(img,cnt,flag):
    file_tem=['citywall.bmp','elain.bmp','lena.bmp','woman.bmp']
    img_tem=cv.imread(file_tem[cnt])
    num,num_temp=np.zeros((256)),np.zeros((256))
    img_shape=img.shape
    if flag:
        for i in range(img_shape[0]):
            for j in range(img_shape[1]):
                num[img[i][j][0]]+=1
                num_temp[img_tem[i][j][0]]+=1
    else:
        for i in range(img_shape[0]):
            for j in range(img_shape[1]):
                if img[i][j][0]==0:
                    img[i][j][0]=255
                num[img[i][j][0]]+=1
                num_temp[img_tem[i][j][0]]+=1
    fre=num/np.sum(num)
    fre_temp=num_temp/np.sum(num_temp)
    z,z_temp=np.zeros((256)),np.zeros((256))
    z[0],z_temp[0]=fre[0]*256,fre_temp[0]*256
    for i in range(256):
        if i:
            z[i]=z[i-1]+fre[i]*255
            z_temp[i]=z_temp[i-1]+fre_temp[i]*255
    trans=np.zeros((256))
    for i in range(256):
        trans[i]=find(z_temp,z[i])
    trans=trans.astype(np.uint8)
    new_img=np.zeros((img_shape[0],img_shape[1]))
    for i in range(img_shape[0]):
        for j in range(img_shape[1]):
            new_img[i][j]=trans[img[i][j][0]]
    return new_img.astype(np.uint8)
def cnt(str):
    if 'citywall' in str:
        return 0
    if 'elain' in str:
        return 1
    if 'lena' in str:
        return 2
    if 'woman' in str:
        return 3
def flag(str):
    if '1' in str:
        return 0
    return 1
filename=['citywall.bmp','citywall1.bmp','citywall2.bmp','elain.bmp','elain1.bmp','elain2.bmp','elain3.bmp','lena.bmp','lena1.bmp','lena2.bmp','lena3.bmp','woman.BMP','woman1.bmp','woman2.bmp']
for i in filename:
    img=cv.imread(i)
    result=graybar_trans(img,cnt(i),flag(i))
    # cv.namedWindow(i,cv.WINDOW_FREERATIO)
    # cv.imshow(i,result)
    # cv.waitKey(0)
    cv.imwrite("aim_bar_"+i,result)

4.局部直方图均衡:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
import cv2 as cv
import numpy as np

def trans(imname):
    img=cv.imread(imname)
    shape=img.shape
    cnt_i=shape[0]//7+1
    cnt_j=shape[0]//7+1
    result=np.zeros((shape[0],shape[1]))
    for i in range(cnt_i):
        for j in range(cnt_j):
            num=np.zeros(256)
            for point_i in range(i*7,i*7+7):
                for point_j in range(j*7,j*7+7):
                    if point_i<shape[0] and point_j<shape[1]:
                        num[img[point_i][point_j][0]]+=1
            fre=num/np.sum(num)
            T=np.zeros((256))
            T[0]=fre[0]*256
            for k in range(256):
                if k:
                    T[k]=T[k-1]+fre[k]*255
            for point_i in range(i*7,i*7+7):
                for point_j in range(j*7,j*7+7):
                    if point_i<shape[0] and point_j<shape[1]:
                        result[point_i][point_j]=T[img[point_i][point_j][0]]
    result=result.astype(np.uint8)
    cv.imwrite("trans_"+imname,result)

trans('lena.bmp')
trans('elain.bmp')