OpenCV-Part4

OpenCV-Part4——图像处理（下）

[TOC]

直方图

定义

• 直方图，是在X轴上具有像素值（不总是从0到255的范围），在Y轴上具有图像中相应像素数的图。

• 直方图，可以用来总体了解图像的强度分布。

• 直方图只是理解图像的另一种方式。通过查看图像的直方图，可以直观地了解该图像的对比度，亮度，强度分布等。

• 

  如图，此直方图是针对灰度图像而非彩色图像绘制的。其左侧区域显示图像中较暗像素的数量，而右侧区域则显示明亮像素的数量。

  从该直方图中，可以看到暗区域多于亮区域，而中间灰度值的数量就非常少。

• 有关直方图的术语：

  • BINS：像素值的区间段数。将整个直方图在 x 轴上分成 n 个子部分，每个子部分的值就是其中所有像素数的总和，每个子部分都称为一个 BIN 。当 x 轴采用 255 个值来展示直方图时，BINS=255。可用于对具体像素精度不在乎，而意图找固定区间段划分的像素总量的情况。在 OpenCV 中也称作 histSize 。
  • DIMS：收集数据的参数量。当仅收集关于强度值的一件事的数据时，这里是1。
  • RANGE：测量的强度值范围。通常为[0, 256]，即所有强度值。

查找并绘制直方图

OpenCV计算直方图

• cv2.calcHist(images，channels，mask，histSize，ranges，hist，accumulate)：查找直方图。
  • images：uint8或float32类型的源图像。需要由列表[img]包裹。
  • channels：计算直方图的通道索引。对于灰度图像，则其值为[0]。对于彩色图像，可以传递[0]、[1]或[2]分别计算蓝色，绿色或红色通道的直方图。需要由列表[channels]包裹。
  • mask：图像掩码。None：查找完整图像的直方图。否则为查找特定区域mask掩码的直方图。
  • histSize：像素值的区间段数。需要放在方括号中。对于全尺寸，我们通过[256]。需要由列表[histSize]包裹。
  • ranges：测量的强度值范围。通常为[0, 256]。
• Numpy 也有计算直方图的函数，但是 cv2.calcHist() 比 np.histogram() 快大约40倍。因此，尽可能使用OpenCV函数。

img = cv2.imread('home.jpg', 0)
hist = cv2.calcHist([img], channels=[0], mask=None, histSize=[256], ranges=[0, 256])

Matplotlib绘制直方图

1. matplotlib.pyplot.hist() ：Matplotlib自带直方图绘图功能，能够直接找到直方图并将其绘制。而无需使用 cv2.calcHist() 查找直方图。

   • import numpy as np
     import cv2
     from matplotlib import pyplot as plt
     
     image = cv2.imread('home.jpg')
     # numpy的ravel函数功能是将多维数组降为一维数组
     plt.hist(image.ravel(), bins=256, range=[0, 256], histtype="step", color='black')  # 对整图的，需要注意和灰度是有区别的！
     plt.hist(image[:, :, 0].ravel(), bins=256, range=[0, 256], histtype="step", color='blue')
     plt.hist(image[:, :, 1].ravel(), bins=256, range=[0, 256], histtype="step", color='green')
     plt.hist(image[:, :, 2].ravel(), bins=256, range=[0, 256], histtype="step", color='red')
     plt.show()
        
     2. `cv2.calcHist()`+`plt.plot()`：先找到直方图数据，再使用matplotlib的法线图。分离了数据和画图，更符合程序解耦设计。
     
        - ```python
          import numpy as np
          import cv2
          from matplotlib import pyplot as plt
          
          image = cv2.imread('home.jpg')
          color = ["blue", "green", "red"]
          ranges = [0, 256]
          # 同时列出数据下标和数据
          for i, color in enumerate(color):
              hist = cv2.calcHist([image], channels=[i], mask=None, histSize=[256], ranges=ranges)
              print(type(hist), hist.shape)  # <class 'numpy.ndarray'>  (256, 1)
              plt.plot(hist, color=color)
              plt.xlim(ranges)
          plt.show()
     ![img](OpenCV-Part4/histogram_rgb_plot.jpg)
     

直方图均衡化

定义

• 直方图均衡化，可以使限制在某一区间的色彩亮度均衡到整个色彩区间上，使像素在直方图上的分布更均匀。
• 直方图均衡化，是用来改善图像的全局亮度和对比度的。例如，在对人脸数据进行训练之前，对人脸图像进行直方图均衡化处理，使其具有相同的光照条件。
• 

全局均衡化

• cv2.equalizeHist()：只能输入灰度图像，输出直方图均衡化后的图像。
• np.hstack((img1, img2))：横向拼接。

equ = cv2.equalizeHist(img)
cv2.imshow('equalization', np.hstack((img, equ)))  # 并排显示
cv2.waitKey(0)

自适应直方图均衡化(CLAHE)

• 当一张图同时拥有较多暗部和亮部时，全局均衡化极其容易丢失信息。
• 自适应直方图均衡，可以将图像分成称为 tiles 的小块（在 OpenCV 中，tileSize 默认为 8x8 ）。并在每个较小的区域中，进行直方图均衡。
• 如果在 tiles 中有噪音，则应用了对比度限制：任何直方图 bin 超出指定的对比度限制（在OpenCV中默认为40），则在应用直方图均衡之前，将这些像素裁剪并均匀地分布到其他bin。
• cv2.createCLAHE(clipLimit, tileGridSize)：CLAHE 自适应直方图均衡化。
  • clipLimit：颜色对比度的阈值。
  • titleGridSize：进行像素均衡化的网格大小。

cl1 = cv2.createCLAHE(clipLimit=2.0, tileGridSize=(8, 8)).apply(img)
cv2.imshow('equalization', np.hstack((img, cl1)))  # 并排显示
cv2.waitKey(0)

二维直方图

定义

• 上述都是在计算一维直方图，因为仅考虑一个特征，即像素的灰度强度值。
• 但是在二维直方图中，需要考虑两个特征：每个像素的 色相(Hue) 和 饱和度(Saturation) 值。
• 二维直方图，通常用于查找颜色直方图。

查找并绘制直方图

OpenCV计算二维直方图

• 对于颜色直方图，我们需要将图像从BGR转换为HSV。对于一维直方图，我们从BGR转换为灰度

• cv2.calcHist()：对于二维直方图，使用相同的函数，但其参数将进行如下修改：

  • channels=[0, 1]：因为需要同时处理H和S平面。
  • histSize=[180, 256]：对于H平面为180，对于S平面为256。
  • ranges=[0, 180, 0, 256]：色相值介于0和180之间，饱和度介于0和256之间。

import cv2

img = cv2.imread('home.jpg')
hsv = cv2.cvtColor(img, cv2.COLOR_BGR2HSV)
hist = cv2.calcHist([hsv], channels=[0, 1], mask=None, histSize=[180, 256], ranges=[0, 180, 0, 256])

Matplotlib绘制二维直方图

• matplotlib.pyplot.imshow()：绘制2D直方图。
• 使用此功能时，插值法应采用最近邻以获得更好的结果。

import cv2
from matplotlib import pyplot as plt

img = cv2.imread("test.jpg")
hsv = cv2.cvtColor(img, cv2.COLOR_BGR2HSV)
hist = cv2.calcHist([hsv], channels=[0, 1], mask=None, histSize=[180, 256], ranges=[0, 180, 0, 256])
plt.imshow(hist, interpolation='nearest')
plt.show()

• 下面是输入图像及其颜色直方图。X轴显示S值，Y轴显示色相。

  

• 在该直方图中，可以在H = 100 和 S = 200 附近看到一些较高的值，对应于天空的蓝色。同样，在 H = 25 和 S = 100 附近可以看到另一个峰值，对应于宫殿的黄色。

YUV色彩空间对彩色图像做直方图均衡化

• YUV色彩空间是把亮度（Luma）与色度（Chroma）分离。

  • “Y”表示亮度，也就是灰度值。
  • “U”表示蓝色通道与亮度的差值。
  • “V”表示红色通道与亮度的差值。

• 对彩色图像进行直方图均衡化时，先将图像从RGB空间转到YUV空间，然后对亮度Y通道进行直方图均衡化得到通道Y”，然后将Y”UV通道进行合并。

import cv2
import numpy as np

img = cv2.imread('input.jpg')

img_yuv = cv2.cvtColor(img, cv2.COLOR_BGR2YUV)
# equalize the histogram of the Y channel
img_yuv[:, :, 0] = cv2.createCLAHE(clipLimit=2.0, tileGridSize=(8, 8)).apply(img_yuv[:, :, 0])
# convert the YUV image back to RGB format
img_output = cv2.cvtColor(img_yuv, cv2.COLOR_YUV2BGR)

cv2.imshow('Color input image', img)
cv2.imshow('Histogram equalized', img_output)
cv2.waitKey(0)

直方图比较

• 对输入的两张图像进行直方图均衡化后，可以对两个图像的直方图进行对比，从对比的结果得到一些的结论。
• cv2.compareHist(H1, H2, method)：对比两个图像的直方图。method有三中可用：
  • cv2.HISTCMP_BHATTACHARYYA：巴氏距离，越小越相似。
  • cv2.HISTCMP_CORREL：相关性，越大越相似。
  • cv2.HISTCMP_CHISQR：卡方，越大越不相似。

import cv2
import numpy as np


def create_rgb_hist(image):
    h, w, c = image.shape
    rgbHist = np.zeros([16 * 16 * 16, 1], np.float32)
    bsize = 256 / 16
    for row in range(h):
        for col in range(w):
            b = image[row, col, 0]
            g = image[row, col, 1]
            r = image[row, col, 2]
            index = np.int(b / bsize) * 16 * 16 + np.int(g / bsize) * 16 + np.int(r / bsize)
            rgbHist[np.int(index), 0] = rgbHist[np.int(index), 0] + 1
    return rgbHist


def hist_compare(image1, image2):
    hist1 = create_rgb_hist(image1)
    hist2 = create_rgb_hist(image2)
    match1 = cv2.compareHist(hist1, hist2, cv2.HISTCMP_BHATTACHARYYA)
    match2 = cv2.compareHist(hist1, hist2, cv2.HISTCMP_CORREL)
    match3 = cv2.compareHist(hist1, hist2, cv2.HISTCMP_CHISQR)
    # 巴氏距离越小越相似，相关性越大越相似，卡方越大越不相似
    print("巴氏距离", match1)
    print("相关性", match2)
    print("卡方", match3)


if __name__ == "__main__":
    img1 = cv2.imread("img5.jpg")
    img2 = cv2.imread("img6.jpg")
    cv2.imshow("img1", img1)
    cv2.imshow("img2", img2)
    hist_compare(img1, img2)
    cv2.waitKey(0)
    cv2.destroyAllWindows()

直方图反投影

定义

• 直方图反投影，常用于图像分割或在图像中查找感兴趣的对象。本质是定位模板图像在输入图像的位置。
• 直方图反投影，一般使用的是颜色直方图，颜色直方图比灰度直方图更容易定义对象。
• 流程：1. 选定一个目标图像，这个图像应包含并尽可能填充我们的检测目标对象；2. 创建该图像的直方图；3. 将该直方图反投影到需要检测的测试图像上，计算出每个像素属于该目标的概率。
• 直方图反投影可以与 camshift 算法等配合使用。
• 我选择 U-Net。

OpenCV的反投影

• cv2.calcBackProject(images, channels, hist, ranges, scale, dst)：
  • hist：目标进行归一化后的直方图。

import numpy as np
import cv2

roi = cv2.imread('rose_red.png')
hsv = cv2.cvtColor(roi, cv2.COLOR_BGR2HSV)
target = cv2.imread('rose.png')
hsvt = cv2.cvtColor(target, cv2.COLOR_BGR2HSV)
# 1. 计算检测目标的直方图
roihist = cv2.calcHist([hsv], [0, 1], None, [180, 256], [0, 180, 0, 256])
# 2. 归一化目标直方图
cv2.normalize(roihist, roihist, 0, 255, cv2.NORM_MINMAX)
# 3. 反传算法
dst = cv2.calcBackProject([hsvt], channels=[0, 1], hist=roihist, ranges=[0, 180, 0, 256], scale=1)
# 4. 用圆盘卷积核对概率图像进行均值滤波
disc = cv2.getStructuringElement(cv2.MORPH_ELLIPSE, (5, 5))
cv2.filter2D(dst, -1, disc, dst)
# 5. 应用阈值作，获得mask
ret, thresh = cv2.threshold(dst, 50, 255, 0)
thresh = cv2.merge((thresh, thresh, thresh))
res = cv2.bitwise_and(target, thresh)
res = np.vstack((target, thresh, res))  # 上下拼接
cv2.imwrite('res.jpg', res)

• 如下示例中，蓝色矩形的区域为目标对象，使用直方图反投影提取整个地面。

  

傅里叶变换

定义

• 傅立叶变换，可以用于分析各种滤波器的频率特性。

• 对于正弦信号：如果幅度在短时间内变化快，则为高频信号；如果变化慢，则为低频信号。

  可以将相同的想法扩展到图像：图像中的振幅在边缘点或噪声急剧变化，边缘和噪声就是图像中的高频内容；如果图像中幅度没有太大变化，则就是低频分量。

• 傅里叶变换，可以将一幅图片分解为正弦和余弦两个分量，即可以将一幅图像从其空间域（spatial domain）转换为频域（frequency domain）。

详细讲解傅里叶的文章：

• https://zhuanlan.zhihu.com/p/19759362
• https://zhuanlan.zhihu.com/p/28478034

• 

• 正弦波就是一个圆周运动在一条直线上的投影。所以频域的基本单元也可以理解为一个始终在旋转的圆：

  

• 时域与频域之间的联系：

  

效果

1. 灰度图 与 傅里叶变换：

   • 可以看到白色区域大多在中心，显示低频率的内容比较多。

   • 

2. 傅里叶变换删去低频内容：

   • 删除图像中的低频内容，即将HPF应用于图像，本质是提取边缘。

   • 

3. 傅里叶变换删去高频内容：

   • 删除图像中的高频内容，即将LPF应用于图像，本质是模糊图像。

   • 

4. 各滤波器是 HPF（High Pass Filter）还是 LPF（Low Pass Filter）：

   • 白点在中间的就是 LPF ，白点在四周的就是 HPF。

   • 

实现

• 对于图像，可以使用2D离散傅里叶变换（DFT, 2D Discrete Fourier Transform）查找频域。同时，可以使用快速傅立叶变换（FFT, Fast Fourier Transform）的快速算法用于DFT的计算。
• OpenCV使用 cv2.dft()、cv2.idft() 实现傅里叶变换，效率更高。
• Numpy使用 np.ifft2()、np.fft.ifftshift() 实现傅里叶变换，使用更友好。
• DFT的性能优化：在一定的阵列尺寸下，DFT计算的性能较好。当数组大小为2的幂时，速度最快。大小为2、3和5的乘积的数组也可以非常有效地处理。
  为达到最佳性能，可以通过OpenCV提供的函数 cv2.getOptimalDFTSize() 寻找最佳尺寸。
  然后将图像填充成最佳性能大小的阵列：1. 对于OpenCV，必须手动填充零；2. 对于Numpy，可以指定FFT计算的新大小，自动填充零。
• 通过使用最优阵列，大概能提升4倍的效率。而OpenCV本身也比Numpy效率快近3倍。

Numpy中的傅里叶变换

import cv2
import numpy as np
from matplotlib import pyplot as plt

img = cv2.imread('test.jpg', 0)

# 1. 使用Numpy实现傅里叶变换 fft.fft2()
f = np.fft.fft2(img)

# 2. 一旦得到结果，零频率分量（直流分量）将出现在左上角。
# 如果要将其置于中心，则需要使用 np.fft.fftshift() 将结果在两个方向上移动。
# 一旦找到了频率变换，就能找到幅度谱。
fshift = np.fft.fftshift(f)
magnitude_spectrum = 20 * np.log(np.abs(fshift))

plt.subplot(121), plt.imshow(img, cmap='gray')
plt.title('Input Image'), plt.xticks([]), plt.yticks([])
plt.subplot(122), plt.imshow(magnitude_spectrum, cmap='gray')
plt.title('Magnitude Spectrum'), plt.xticks([]), plt.yticks([])
plt.show()

# 3. 找到了频率变换，就可以进行高通滤波和重建图像，也就是求逆DFT
rows, cols = img.shape
crow, ccol = rows // 2, cols // 2
fshift[crow - 30:crow + 30, ccol - 30:ccol + 30] = 0
f_ishift = np.fft.ifftshift(fshift)
img_back = np.fft.ifft2(f_ishift)
img_back = np.abs(img_back)

# 4. 众所周知，高通滤波是一种边缘检测操作。这也表明大部分图像数据存在于频谱的低频区域。
# 仔细观察结果可以看到最后一张用JET颜色显示的图像，有一些瑕疵，即显示了一些波纹状的结构，这是振铃效应。
# 这是由于用矩形窗口mask造成的，掩码mask被转换为sinc形状，从而导致此问题。所以矩形窗口不用于过滤，更好的选择是高斯mask。
plt.subplot(131), plt.imshow(img, cmap='gray')
plt.title('Input Image'), plt.xticks([]), plt.yticks([])
plt.subplot(132), plt.imshow(img_back, cmap='gray')
plt.title('Image after HPF'), plt.xticks([]), plt.yticks([])
plt.subplot(133), plt.imshow(img_back)
plt.title('Result in JET'), plt.xticks([]), plt.yticks([])

plt.show()

OpenCV中的傅里叶变换

import cv2
import numpy as np
from matplotlib import pyplot as plt

img = cv2.imread('test.jpg', 0)
rows, cols = img.shape

# 1. 计算DFT效率最佳的尺寸
nrows = cv2.getOptimalDFTSize(rows)
ncols = cv2.getOptimalDFTSize(cols)
print(f"rows: {rows}->{nrows}; cols: {cols}->{ncols}")

# 2. 将大小扩充成2、3和5的乘积
nimg = np.zeros((nrows, ncols))
nimg[:rows, :cols] = img
img = nimg

# 3. OpenCV计算快速傅里叶变换: 输入图像应首先转换为 np.float32 ，然后使用函数 cv2.dft() 和 cv2.idft() 。
# 返回结果与Numpy相同，但有两个通道。第一个通道为有结果的实部，第二个通道为有结果的虚部。
dft = cv2.dft(np.float32(img), flags=cv2.DFT_COMPLEX_OUTPUT)
dft_shift = np.fft.fftshift(dft)
magnitude_spectrum = 20 * np.log(cv2.magnitude(dft_shift[:, :, 0], dft_shift[:, :, 1]))

plt.subplot(121), plt.imshow(img, cmap='gray')
plt.title('Input Image'), plt.xticks([]), plt.yticks([])
plt.subplot(122), plt.imshow(magnitude_spectrum, cmap='gray')
plt.title('Magnitude Spectrum'), plt.xticks([]), plt.yticks([])
plt.show()

# 4. 删除图像中的高频内容(将LPF应用于图像/模糊图像)
# 首先创建一个在低频时具有高值的mask掩码，中心正方形为1，其他均为0。即传递LF内容，在HF区域为0。
rows, cols = img.shape
crow, ccol = rows // 2, cols // 2
mask = np.zeros((rows, cols, 2), np.uint8)
mask[crow - 30:crow + 30, ccol - 30:ccol + 30] = 1

# 5. 应用掩码Mask和求逆DTF
fshift = dft_shift * mask
f_ishift = np.fft.ifftshift(fshift)
img_back = cv2.idft(f_ishift)
img_back = cv2.magnitude(img_back[:, :, 0], img_back[:, :, 1])

plt.subplot(121), plt.imshow(img, cmap='gray')
plt.title('Input Image'), plt.xticks([]), plt.yticks([])
plt.subplot(122), plt.imshow(img_back, cmap='gray')
plt.title('Magnitude Spectrum'), plt.xticks([]), plt.yticks([])
plt.show()

HPF or LPF

import cv2
import numpy as np
from matplotlib import pyplot as plt

# 简单的均值滤波
mean_filter = np.ones((3, 3))
# 构建高斯滤波
x = cv2.getGaussianKernel(5, 10)
gaussian = x * x.T
# 不同的边缘检测算法Scharr-x方向
scharr = np.array([[-3, 0, 3],
                   [-10, 0, 10],
                   [-3, 0, 3]])
# Sobel_x
sobel_x = np.array([[-1, 0, 1],
                    [-2, 0, 2],
                    [-1, 0, 1]])
# Sobel_y
sobel_y = np.array([[-1, -2, -1],
                    [0, 0, 0],
                    [1, 2, 1]])
# 拉普拉斯
laplacian = np.array([[0, 1, 0],
                      [1, -4, 1],
                      [0, 1, 0]])

filters = [mean_filter, gaussian, laplacian, sobel_x, sobel_y, scharr]
filter_name = ['mean_filter', 'gaussian', 'laplacian', 'sobel_x', 'sobel_y', 'scharr_x']
fft_filters = [np.fft.fft2(x) for x in filters]
fft_shift = [np.fft.fftshift(y) for y in fft_filters]
mag_spectrum = [np.log(np.abs(z) + 1) for z in fft_shift]

for i in range(6):
    plt.subplot(2, 3, i + 1), plt.imshow(mag_spectrum[i], cmap='gray')
    plt.title(filter_name[i]), plt.xticks([]), plt.yticks([])
plt.show()

模板匹配

• 模板匹配，在较大图像中查找模板图像位置的方法。
• 我选择深度学习。

OpenCV中的模板匹配

import cv2
from matplotlib import pyplot as plt

img = cv2.imread('test.jpg', 0)
template = cv2.imread('test2.png', 0)
w, h = template.shape[::-1]
# 列表中所有的6种比较方法
methods = ['cv2.TM_CCOEFF', 'cv2.TM_CCOEFF_NORMED', 'cv2.TM_CCORR',
           'cv2.TM_CCORR_NORMED', 'cv2.TM_SQDIFF', 'cv2.TM_SQDIFF_NORMED']
for meth in methods:
    img_ = img.copy()
    method = eval(meth)
    # 应用模板匹配
    res = cv2.matchTemplate(img_, template, method)
    min_val, max_val, min_loc, max_loc = cv2.minMaxLoc(res)
    # 如果方法是TM_SQDIFF或TM_SQDIFF_NORMED，则取最小值
    if method in [cv2.TM_SQDIFF, cv2.TM_SQDIFF_NORMED]:
        top_left = min_loc
    else:
        top_left = max_loc
    bottom_right = (top_left[0] + w, top_left[1] + h)
    cv2.rectangle(img_, top_left, bottom_right, 255, 2)
    
    plt.subplot(121), plt.imshow(res, cmap='gray')
    plt.title('Matching Result'), plt.xticks([]), plt.yticks([])
    plt.subplot(122), plt.imshow(img_, cmap='gray')
    plt.title('Detected Point'), plt.xticks([]), plt.yticks([])
    plt.suptitle(meth)
    plt.show()

多对象模板匹配

import cv2
import numpy as np

img_rgb = cv2.imread('mario.png')
img_gray = cv2.cvtColor(img_rgb, cv2.COLOR_BGR2GRAY)
template = cv2.imread('mario_coin.png', 0)
w, h = template.shape[::-1]
res = cv2.matchTemplate(img_gray, template, cv2.TM_CCOEFF_NORMED)
threshold = 0.8
loc = np.where(res >= threshold)
for pt in zip(*loc[::-1]):
    cv2.rectangle(img_rgb, pt, (pt[0] + w, pt[1] + h), (0, 0, 255), 2)
cv2.imwrite('res.png', img_rgb)