计算机视觉概念原理与代码实践之图像相关操作代码及分析

系列文章目录

实验一的内容是图像去噪及直方图均衡化的完整代码，其作者是李忆如，该代码位于 Gitee.com 上。

本系列博客着重于计算机视觉的概念原理以及代码实践，其中不涵盖繁杂的数学推导。如果有问题，欢迎在评论区进行讨论并指出，也可以直接私信联系我。

第一章计算机视觉——图像去噪及直方图均衡化（图像增强）

梗概

- 与 opencv 库的效果进行比较分析，并附带 python 代码。

一、实验内容与方法

实验内容：选择一到两幅图像，不借助库函数，通过自己动手编写程序来实现图像的去噪操作，包括添加高斯噪声和椒盐噪声，同时还要进行直方图均衡化处理，并且将自己实现的效果与 OpenCV 的库函数实现的效果进行对比分析。

实验图片：实验选取图片如图1与图2所示：

二、图像的加噪

图像噪声指的是在图像数据里存在的那些并非必要或者多余的干扰信息。因为其概率分布以及产生方式的不同，所以可以划分成不同的类型。在本次实验当中，是以高斯噪声和椒盐噪声作为例子的。

1.高斯噪声

高斯噪声指的是服从高斯分布也就是正态分布的一类噪声。这类噪声通常是由不良照明以及高温所引起的传感器噪声。

加入高斯噪声的代码实现（以灰度图的加噪为例）如下：

定义一个函数 gaussian_noise，该函数接受三个参数，分别是 img、mean 和 sigma，用于给图像添加高斯噪声。
    """
    此函数用将产生的高斯噪声加到图片上
原图被传入，mean 为均值，sigma 为标准差
    返回: gaussian_out : 噪声处理后的图片
    """
    # 将图片灰度标准化
    img = img / 255
    # 产生高斯 noise
噪声等于使用 np.random.normal 函数生成的，其均值为 mean，标准差为 sigma，形状与 img 相同的随机数数组。
    # 将噪声和图片叠加
    gaussian_out = img + noise
    # 将超过 1 的置 1，低于 0 的置 0
gaussian_out 被限制在 0 到 1 之间，即 gaussian_out 被 np.clip 函数处理，其下限为 0，上限为 1
    # 将图片灰度范围的恢复为 0-255
使用 numpy 库中的 np.uint8 函数，将 gaussian_out 乘以 255 后转换为无符号 8 位整数类型，结果存储在 gaussian_out 中。
    # 将噪声范围搞为 0-255
    # noise = np.uint8(noise*255)
会返回噪声，注意返回值，这里也会返回 gaussian_out 。

图 1 中的待实验图片，在加入以 mean = 0 、sigma = 0.1 为例的高斯噪声前后，其情况如图所示。

2.椒盐噪声

椒盐噪声是给图片添加黑白噪点的一种噪声。椒代表黑色的噪点，其颜色为(0,0,0)；盐代表白色的噪点，其颜色为(255,255,255)。可以通过设置 amount 来控制添加噪声的比例，amount 值越大，添加的噪声就越多，图像损坏也就越严重。

加入椒盐噪声的代码实现（以灰度图的加噪为例）如下：

定义一个名为 sp_noise 的函数，该函数接受两个参数，一个是包含噪声的图像 noise_img，另一个是比例 proportion
    """
    添加椒盐噪声
proportion 的值代表加入噪声的量，能够依据需求自行进行调整
    return: img_noise
    """
height 等于 noise_img 的形状中的第一个元素，width 等于 noise_img 的形状中的第二个元素。
获取高度宽度像素值，即 noise_img.shape[1] 。
num 等于 height 乘以 width 再乘以 proportion 的结果  # 一个准备加入多少噪声小点
    for i in range(num):
w 是从 0 到 width 减 1 这个范围内随机生成的整数。
h 是通过 random.randint(0, height - 1) 这个函数生成的一个随机整数，其取值范围在 0 到 height - 1 之间。
        if random.randint(0, 1) == 0:
            noise_img[h, w] = 0
        else:
            noise_img[h, w] = 255
    return noise_img

图 1 的待实验图片加入均值为 0 、标准差为 0.1 的高斯噪声前后如图所示，同时该待实验图片加入比例为 0.1 的椒盐噪声前后也如图所示。

三、图像的去/降噪

噪声的存在对遥感图像的质量产生了严重影响。所以在进行图像增强处理和分类处理之前，需要使用图像的去/降噪方法来进行纠正。本实验选取中值滤波与均值滤波作为例子。

1.均值滤波

滤波器领域内的像素平均值，是通过对目标像素及其周边像素取平均值，然后将该平均值填回目标像素，以此来实现滤波目的的方法，加权平均公式如图所示。

最后，得到经过均值滤波处理后的灰度图。

# 均值降噪转换器
def convert_2d(r):
    n = 3
    # 3*3 滤波器, 每个系数都是 1/9
窗口等于一个由 n 行 n 列组成的数组，该数组的每个元素都为 1 除以 n 的平方，即 np.ones((n, n)) 除以 n ** 2 的结果。
    # 使用滤波器卷积图像
    # mode = same 表示输出尺寸等于输入尺寸
    # boundary 表示采用对称边界条件处理图像边缘
s 等于 scipy 库中的 signal 模块的 convolve2d 函数对 r 进行卷积操作的结果。
窗口，模式为“相同”，边界为“对称”。
    return s.astype(np.uint8)

图 1 中加入了以 mean = 0 、sigma = 0.1 为例的高斯噪声以及以 proportion = 0.1 为例的椒盐噪声。待实验图片降噪前后的情况如图所示。

从图中可以看出，线性非加权的均值滤波器在降噪方面效果不佳，并且滤波之后的图像比较模糊，会丢失较多的信息。所以在实际使用中，通常会使用加权均值滤波器或者其他类型的滤波器。

2.中值滤波

类似于均值滤波的这种方式，是用像素所在领域内的中间值来替代该像素，其公式与原理如同 0 所展示的那样。

使用中值滤波器的代码实现（以灰度图的降噪为例）如下：

# 自定义中值降噪滤波器
定义一个名为 medium_filter 的函数，该函数接受三个参数，分别是图像 im、坐标 x 和 y 以及步长 step。
    sum_s = []
对于 k 而言，k 的取值范围是从负的（步长除以 2）到正的（步长除以 2）再加上 1 。
对于 m 而言，它的取值范围是从负的（步长除以 2）到正的（步长除以 2）再加上 1 。
sum_s 会添加 im 数组中索引为 x + k 和 y + m 位置的元素
    sum_s.sort()
计算总和 sum_s 中索引为 (int(step * step / 2) + 1) 的元素值并返回

加入高斯噪声（以mean = 0 ， sigma = 0.1为例）与椒盐噪声（proportion = 0.1为例）的图1的待实验图片降噪前后如图所示。

从图中可以看出，图片里的噪点明显变少了。这证明中值滤波器的降噪效果比较好，相较于均值滤波器而言。并且，滤波在去除噪音的过程中，能够较好地保留边的锐度以及图像的细节。所以，本实验的去/降噪方法最终选用了中值滤波器。

四、直方图与均衡化 1.直方图

在图像处理中，经常需要对图像的亮度进行分析，也就是分析像素级的分布情况，这时就需要用到直方图，比如颜色直方图、灰度直方图等。

直方图是对灰度级在范围[0，L-1]的数字图像（彩图分解为 RGB）进行构造的。它使用了如图的离散函数，并进行了归一化处理。

直方图构造（以灰度图为例）核心代码如下：

对图像进行灰度图像直方图统计的是 histogram 函数，其参数为 Image 。
    im = Image
    a = [0] * 256  # 构造储存像素数的一维数组
    w = im.shape[0]  # 得到图像的宽高
    h = im.shape[1]
对于 i 在范围 w 中进行遍历，目的是统计所有像素点
        for j in range(h):
            gray = im[i, j]
            a[gray] += 1
    y = a  # y表示灰度值的像素个数
x 是通过在 range(256)中遍历每个元素 i 得到的列表  # x 表示 0 到 256 的灰度

图1（原图）与图2（原图）的直方图如图所示。

2.直方图均衡化

直方图均衡化：已知一个图像具有灰度概率密度分布。对这个图像进行一种灰度变换，其变换函数的离散形式如图所示。通过这种变换，使图像演变为具有均匀灰度概率密度分布的新图像，并且尽量让新图像占满全部灰度。这样做增加了像素灰度值的动态范围，从而达到增强图像整体对比度的效果。

直方图均衡化的实现步骤如图所示

五、灰度图均衡化 1.非调库均衡化

依据直方图以及均衡化的定义，对经过降噪处理后的图 1 进行直方图的构建以及均衡化操作。非调库灰度图均衡化的代码实现情况如下：

定义一个均衡化函数，该函数接收两个参数，一个是灰度像素数 a，另一个是图片 Image。
    im = array(Image)
    b = [0] * 256  # 储存灰度像素所占的比例
    c = [0] * 256  # 储存累计分布概率
    w = im.shape[0]  # 图像的宽和高
    h = im.shape[1]
    mn = w * h  # mn表示像素点总数
对于 i 在 a 的长度范围内进行遍历：  # 用于计算灰度分布密度
        b[i] = a[i] / mn
    进行计算累计直方图数据 a 的操作。
        if i == 1:
            c[i] = b[i]
        else:
            c[i] = c[i - 1] + b[i]
            a[i] = int(255 * c[i])
    在范围 w 里的每一个 i 都用于更新像素点的灰度值。
        for j in range(h):
            im[i, j] = a[Image[i, j]]
    return im

均衡化前后直方图如图所示：

均衡化前后图片如图所示：

分析：从图中能够看出，经过均衡化处理后，直方图的灰度分布变得均匀了，并且几乎占满了所有可能的灰度范围，这符合相关的定义。观察该图可以发现，均衡化之后的图片对比度有所增强，同时能够看到更多的细节，这也满足了预期，说明均衡化是成功的。

2. 调库均衡化

Python 的 Opencv（cv2）库中，有一个灰度均衡化函数是 cv2.equalizeHist(img)，它返回的就是经过均衡化处理后的灰度图像。

六、彩图均衡化

彩图本质上是 RGB 三通道图。在进行均衡化时，可以根据不同的考虑使用不同的方法。本实验以拆分合并（复现库函数）、HSI 模型与 RGB 模型互换以及库函数这三种方式为例，来实现彩图的均衡化。

1.拆分合并均衡化

核心原理：先把彩图的 RGB 三通道进行拆分。接着对拆分后的每一通道采用类似灰度均衡化的方式进行均衡化操作。最后再将经过均衡化处理后的各通道进行合并，从而得到均衡化后的彩图。核心代码如下：

def my_equalizehist(single):
使用 cv 计算单个图像的直方图，参数为 [single]，通道为 [0]，不使用掩码 None，直方图的 bin 数量为 256，范围为 [0, 256]，计算得到的直方图赋值给 hist 。
    num_of_pixels = single.size
    ratio = np.zeros(256)
    transf_map = np.zeros(256)
    result = single.copy()
    j = 0
    for i in hist:
        if j > 0:
ratio[j]等于 i 除以 num_of_pixels 的结果再加上 ratio[j - 1]
        else:
            ratio[j] = i / num_of_pixels
transf_map[j] 等于 ratio[j] 乘以 255 的结果进行四舍五入
        j = j + 1
对于 i 而言，其在 range(single.shape[0])的范围内进行循环。
j 在 single 的形状的第二个维度的范围内进行遍历。
result 的第 i 行第 j 列等于 transf_map 中 single 的第 i 行第 j 列所对应的值
        # result[result==j]=k[j]
    return result

均衡化前后图片如图：

2.库函数均衡化

在 Python 的 Opencv（cv2）库中，用于彩图均衡化的函数是 cv.equalizeHist(img)。它可以用来均衡化各个通道，并且能够使用 cv.merge 进行合并。其具体用法如下：

使用 cv 计算 res 的直方图，参数为[0]，不使用掩码 None，分 256 个区间，区间范围是[0,256]，计算结果存储在 hist_res 中
plt.plot(hist_res);
plt.show()
EB = cv.equalizeHist(B)
EG = cv.equalizeHist(G)
ER = cv.equalizeHist(R)
equal_test 是将 EB、EG、ER 进行合并操作，即将它们合并在一起，形成一个新的元素。具体来说，是通过 cv.merge 函数来实现这个合并操作的，将 EB、EG、ER 作为参数传递给 cv.merge 函数，然后该函数会将它们合并为一个新的元素，并将结果赋值给 equal_test。

3.HSI模型RGB模型互换均衡化

HSI 颜色模型是一个数学模型。这个数学模型高度抽象且模拟。它能够满足计算机数字化颜色管理的需要。HSI 模型是以人的视觉系统为起始点，直接运用颜色的三个要素，即色调（Hue）、饱和度（Saturation）以及亮度（Intensity，有时也被翻译为密度或灰度）去对颜色进行描述。RBG 转换到 HSI 空间的公式呈现于此，HSI 转化为 RGB 的公式也呈现于此。

代码如下：

import cv2
import numpy as np
导入 matplotlib 库的绘图模块 plt 。
import collections
def rgb_hsi(rgb_image):
    # 保存原始图像的行列数
行数等于 rgb_image 的形状中的第一个维度的整数值
cols 等于 rgb_image 的形状中第一个维度的长度，也就是 rgb_image.shape[1] 的值，它是一个整数。
    # 图像复制
    hsi_image = rgb_image.copy()
    # 通道拆分
    b = rgb_image[:, :, 0]
    g = rgb_image[:, :, 1]
    r = rgb_image[:, :, 2]
    # 归一化到[0,1]
    b = b / 255.0
    g = g / 255.0
    r = r / 255.0
    for i in range(rows):
        for j in range(cols):
num 等于 0.5 乘以（（r 数组中第 i 行第 j 列的元素减去 g 数组中第 i 行第 j 列的元素）加上（r 数组中第 i 行第 j 列的元素减去 b 数组中第 i 行第 j 列的元素））
den 等于 r 数组中第 i 行第 j 列元素减去 g 数组中第 i 行第 j 列元素的差的平方，加上 r 数组中第 i 行第 j 列元素减去 b 数组中第 i 行第 j 列元素的差，乘以 g 数组中第 i 行第 j 列元素减去 b 数组中第 i 行第 j 列元素的差，然后对这个结果开平方
            if den == 0:
                hsi_h = 0
            else:
theta 等于 np.arccos(num / den)的结果转换为浮点数
                if b[i, j] <= g[i, j]:
                    hsi_h = theta
                else:
                    hsi_h = 2*np.pi - theta
            min_RGB = min(min(b[i, j], g[i, j]), r[i, j])
            sum = b[i, j]+g[i, j]+r[i, j]
            if sum == 0:
                hsi_s = 0
            else:
                hsi_s = 1 - 3*min_RGB/sum
            hsi_h = hsi_h/(2*np.pi)
            hsi_i = sum/3.0
            # 输出HSI图像，扩充到255以方便显示，一般H分量在[0,2pi]之间，S和I在[0,1]之间
            hsi_image[i, j, 0] = hsi_h*255
            hsi_image[i, j, 1] = hsi_s*255
            hsi_image[i, j, 2] = hsi_i*255
    return hsi_image
def hsi_rgb(hsi_image):
    # 保存原始图像的行列数
    rows = np.shape(hsi_image)[0]
    cols = np.shape(hsi_image)[1]
    # 对原始图像进行复制
    rgb_image = hsi_image.copy()
    # 对图像进行通道拆分
    hsi_h = hsi_image[:, :, 0]
    hsi_s = hsi_image[:, :, 1]
    hsi_i = hsi_image[:, :, 2]
    # 把通道归一化到[0,1]
    hsi_h = hsi_h / 255.0
    hsi_s = hsi_s / 255.0
    hsi_i = hsi_i / 255.0
    B, G, R = hsi_h, hsi_s, hsi_i
    for i in range(rows):
        for j in range(cols):
            hsi_h[i, j] *= 360
            if 0 <= hsi_h[i, j] < 120:
                B = hsi_i[i, j] * (1 - hsi_s[i, j])
                R = hsi_i[i, j] * (1 + (hsi_s[i, j] * np.cos(hsi_h[i, j] * np.pi / 180)) / np.cos(
                    (60 - hsi_h[i, j]) * np.pi / 180))
                G = 3 * hsi_i[i, j] - (R + B)
            elif 120 <= hsi_h[i, j] < 240:
                hsi_h[i, j] = hsi_h[i, j] - 120
                R = hsi_i[i, j] * (1 - hsi_s[i, j])
                G = hsi_i[i, j] * (1 + (hsi_s[i, j] * np.cos(hsi_h[i, j] * np.pi / 180)) / np.cos(
                    (60 - hsi_h[i, j]) * np.pi / 180))
                B = 3 * hsi_i[i, j] - (R + G)
            elif 240 <= hsi_h[i, j] <= 300:
                hsi_h[i, j] = hsi_h[i, j] - 240
                G = hsi_i[i, j] * (1 - hsi_s[i, j])
                B = hsi_i[i, j] * (1 + (hsi_s[i, j] * np.cos(hsi_h[i, j] * np.pi / 180)) / np.cos(
                    (60 - hsi_h[i, j]) * np.pi / 180))
                R = 3 * hsi_i[i, j] - (G + B)
            rgb_image[i, j, 0] = B * 255
            rgb_image[i, j, 1] = G * 255
            rgb_image[i, j, 2] = R * 255
    return rgb_image
# 计算灰度图的直方图
def draw_histogram(grayscale):
    # 对图像进行通道拆分
    hsi_i = grayscale[:, :, 2]
    color_key = []
    color_count = []
    color_result = []
    histogram_color = list(hsi_i.ravel())  # 将多维数组转换成一维数组
    color = dict(collections.Counter(histogram_color))  # 统计图像中每个亮度级出现的次数
    color = sorted(color.items(), key=lambda item: item[0])  # 根据亮度级大小排序
    for element in color:
        key = list(element)[0]
        count = list(element)[1]
        color_key.append(key)
        color_count.append(count)
    for i in range(0, 256):
        if i in color_key:
            num = color_key.index(i)
            color_result.append(color_count[num])
        else:
            color_result.append(0)
    color_result = np.array(color_result)
    return color_result
def histogram_equalization(histogram_e, lut_e, image_e):
    sum_temp = 0
    cf = []
    for i in histogram_e:
        sum_temp += i
        cf.append(sum_temp)
    for i, v in enumerate(lut_e):
        lut_e[i] = int(255.0 * (cf[i] / sum_temp) + 0.5)
    equalization_result = lut_e[image_e]
    return equalization_result
x = []
for i in range(0, 256):  # 横坐标
    x.append(i)
# 原图及其直方图
rgb_image = cv2.imread("lena.jpg")
cv2.imshow('rgb', rgb_image)
histogram = draw_histogram(rgb_image)
plt.bar(x, histogram)  # 绘制原图直方图
plt.savefig('before_histogram.png')
plt.show()
# rgb转hsi
hsi_image = rgb_hsi(rgb_image)
cv2.imshow('hsi_image', hsi_image)
cv2.imwrite('hsi_result.png', hsi_image)
# hsi在亮度分量上均衡化
histogram_1 = draw_histogram(hsi_image)
lut = np.zeros(256, dtype=hsi_image.dtype)  # 创建空的查找表
result = histogram_equalization(histogram_1, lut, hsi_image)  # 均衡化处理
cv2.imshow('his_color_image', result)
cv2.imwrite('his_color.png', result)  # 保存均衡化后图片
# hsi转rgb
image_equ = cv2.imread('his_color.png')  # 读取图像
rgb_result = hsi_rgb(image_equ)
cv2.imshow('rgb_image', rgb_result)
cv2.imwrite('gbr_result.png', rgb_result)
rgb = cv2.imread("gbr_result.png")
histogram_2 = draw_histogram(rgb)
plt.bar(x, histogram_2)
plt.savefig('after_histogram.png')
plt.show()
plt.plot(x, lut)  # 绘制灰度级变换曲线图
plt.savefig('Grayscale_transformation_curve.png')
plt.show()
cv2.waitKey(0)
cv2.destroyAllWindows()

七、图像增强

图像增强是一种能让图像更契合特定应用的图像处理技术。它主要被分成空间域图像增强和频率域图像增强这两类。常见的图像增强方法框架呈现的情况如所示。除了直方图均衡化之外，以下是几种常见的图像增强算法的简介以及它们的原理：

对数变换（即非线性变换）：利用公式 s = c log (1 + r) 来进行灰度变换。这种变换能将窄带低灰度的输入图像值，映射为宽带输出值。灰度变换前后的情况如图 35 所示，图像增强的样例（傅里叶频谱）则如图 36 。

图35 对数变换灰度变化前后

图36 对数变换图像增强样例（图源ppt）

图像的减法运算属于代数运算。它是通过对两幅输入图像进行点对点的加、减、乘、除计算，从而得到输出图像的运算。其中，减法运算通过公式 C(x,y) = f(x,y) - h(x,y)来得到图像。这种运算常用于除去背景噪声，比如缓慢变化的加性非随机噪声。图像增强的样例如图 37 所示。

图37 减法图像增强样例（图源ppt）

均值滤波器和中值滤波器属于平滑滤波器，本实验用到了它们。其原理可在实验步骤与过程中看到，核心是用像素领域内的某些值来代替该像素。

锐化滤波器包含 Laplace 算子和 Sobel 算子。Laplace 算子是微分算子，它主要强调图像中灰度的突变区域，而不强调灰度级缓慢变化的区域。Sobel 算子主要用在工业检测图像的边缘检测。在算子运行时，采取类似卷积的方式，将模板在图像上移动，并在每个位置上计算对应中心像素的梯度值。对一幅灰度图像求梯度所得的结果是一幅梯度图。Sobel 算子的定义如图 38 所示，锐化前后的图像如图 39 所示。

图38 Sobel算子定义

图39 锐化前后图片样例（图源网络）

频率域滤波，具体来说，是先将图像进行傅里叶变换。接着在变换域进行处理。之后再进行傅里叶反变换，从而转换回空间域。其原理在于，用傅里叶变换表示的函数特征，完全能够通过傅里叶反变换来重建，并且不会丢失任何信息，因为任何周期或非周期函数都可以表示为不同频率的正弦函数和余弦函数之和的形式。

高通滤波法：理想高通滤波器（IHPF）的定义如图 40 呈现。滤波前后的样例如图 41 所示。

图40 理想高通滤波器（IHPF）定义

图41 高通滤波处理图像前后

除了图 34 框架和上述提及的方法外，还有一些经典算法可用于图像增强，像灰度世界算法、自动白平衡、自动色彩均衡等。并且近年来，大量深度学习算法被应用于图像增强以及图像处理领域。

八、参考资料

数字图像处理编程方面：用 Python 语言且不调用 OpenCV 函数来实现直方图统计以及直方图均衡，这是 CSDN 博客中的内容。

彩色图像的直方图均衡化（python 实现）

图像处理方面的一个内容是图像增强，这在 QtHalcon 的博客 - CSDN 博客中有相关介绍。

图像增强主要有哪些方法呢？比如对比度增强、亮度调整、色彩变换等方法。对比度增强可以通过提高图像的对比度来使图像更加清晰；亮度调整可以改变图像的整体亮度，使其更亮或更暗；色彩变换可以对图像的色彩进行调整，如改变颜色的饱和度、色调等。这些都是常见的图像增强方法。