OpenCv学习笔记2（图像变换:几何变换，重映射，直方图，二维离散傅立叶变换，用双目视觉图像计算深度信息）

1.几何变换

首先明确：二者的应用场景相同，都是针对二维图片的变换。仿射变换是透视变换的子集。
仿射变换后平行四边形的各边仍操持平行，透视变换结果允许是梯形等四边形，所以仿射变换是透视变换的子集。
仿射变换在图形中的变换包括：平移、缩放、旋转、斜切及它们的组合形式。这些变换的特点是：平行关系和线段的长度比例保持不变。

OpenCV中使用warpAffine()对图像进行仿射变换，调用参数如下：

warpAffine(src, M, dsize[, dst[, flags[, borderMode[, borderValue]]]])

 src参数是变换的原始图像
 dsize参数为返回图像的大小，返回图像的像素类型和src的相同
 M参数是仿射变换的矩阵，它是一个形状为(2,3)的数组
 flags参数是内插方式
 borderMode是外插方式
 borderValue为背景颜色

仿射变换数学表示： $$\left[\begin{array}{ccc}{a}_{00}& a_{01}& b_0\\ a_{10}& a_{11}& b_1\end{array}\right]\ast \left[\begin{array}{c}x\\ y\\ 1\end{array}\right]=\left[\begin{array}{c}{x}^{\prime }\\ y^{\prime }\end{array}\right]$$
OpenCV提供了getAffineTransform(src, dst)来快速计算仿射变换矩阵。src和dst参数是变换前后的三个点的坐标，它们都是形状为(3, 2)的单精度浮点数数组。

warpPerspective()和warpAffine()类似（用法和参数），也对图像进行几何变换，不过它是在三维空间中进行透视变换，因此它的变换矩阵是3×3的矩阵。这个变换矩阵可以通过getPerspectiveTransform(src, dst)计算。src和dst参数是变换前后的4个点的坐标，它们都是形状为(4, 2)的单精度浮点数数组。

import cv2
import matplotlib.pyplot as pl
pl.rcParams["font.family"] = "SimHei"#直接修改配置字典，设置默认字体
img = cv2.imread("E:/ruanjianDM/jupyternoerbookDM/LYF.jpg")
h, w = img.shape[:2]
src = np.array([[0, 0], [w - 1, 0], [0, h - 1]], dtype=np.float32)# ❶
dst = np.array([[300, 300], [873, 78], [161, 923]], dtype=np.float32) #❷
fig, axes = pl.subplots(1, 1, figsize=(9, 4))
m = cv2.getAffineTransform(src, dst) #❸
result = cv2.warpAffine(img, m, (2* w, 1 * h), borderValue=(255, 255, 255, 255)) #❹
axes.set_title("仿射变换")
axes.imshow(result)

❶src为图9-8中三角形的三个顶点坐标，这三个点分别为图像的左上、右上和左下三个顶点。
❷dst为这三个顶点经过仿射变换之后的坐标，图中用三个箭头连接仿射变换前后的坐标点。
❸调用getAffineTransform()得到仿射变换矩阵m，
❹调用warpAffine()对图像img进行仿射变换，结果图像的大小为原始图像的两倍，背景采用白色填充。

2.重映射-remap

对于图像的各种变换都有一个共同特点：它们从原始图像上的某个位置取出一个像素点，并把它绘制到目标图像上的另外一个位置。从原始坐标到目标坐标的映射不一定是一对一的关系.

OpenCV提供了一个通用的图像映射函数remap()来完成这种计算，其调用参数如下：
remap(src, map1, map2, interpolation[, dst[, borderMode[, borderValue]]])

map1和map2参数是两个大小与原始图像src相同的数组，它们的元素值是图像dst中对应下标的像素点在图像src中的坐标值，其元素可以是整数或单精度浮点数。map1中存储映射的X轴坐标，而map2中存储映射的Y轴坐标。下面的数学公式表示了这种映射关系，其中x和y是目标图像中每个像素的坐标，通过map1和map2分别获得它们在src中的坐标。

3.直方图

3.1 绘制直方图

在NumPy中有三个直方图统计函数：histogram()、histogram2d()histogramdd()，分别对应一维数据、二维数据以及多维数的情况。
下面的程序用histogram()和histogram2d()对图像的颜色分布进行统计。一维直方图是因为只考虑一个灰度值的特征。在2D直方图中要考虑两个图像特征。

import cv2
import matplotlib.pyplot as pl
pl.rcParams["font.family"] = "SimHei"#直接修改配置字典，设置默认字体
img = cv2.imread("E:/ruanjianDM/jupyternoerbookDM/LYF.jpg")
fig, ax = pl.subplots(1, 2, figsize=(12, 5))
colors = ["blue", "green", "red"]

for i in range(3):
    #hist是每个小区间出现的次数，x是区间端点值，0.5 * (x[:-1] + x[1:])为每个区间的中值
    hist, x = np.histogram(img[:,:, i].ravel(), bins=256, range=(0, 256)) #❶
    ax[0].plot(0.5 * (x[:-1] + x[1:]), hist, label=colors[i], color=colors[i])
ax[0].legend(loc="upper left")
ax[0].set_xlim(0, 256)

hist2, x2, y2 = np.histogram2d( #❷
                img[:,:, 0].ravel(), img[:,:, 2].ravel(),
                bins=(100, 100), range=[(0, 256), (0, 256)])
ax[1].imshow(hist2, extent=(0, 256, 0, 256), origin="lower", cmap="gray")
ax[1].set_ylabel("blue")
ax[1].set_xlabel("red")

❶通过histogram()对图像img的三个通道分别进行一维直方图统计，由于被统计的数组必须是一维的，因此这里调用数组的ravel()方法将二维数组转换为一维数组。通过range参数指定统计区间为0～256，bin参数指定将统计区间等分为256份。histogram()返回两个数组hist和x，其中hist为统计结果，长度为bin。而x为统计区间，长度为bin+1。hist[i]的值为数组中满足x[i] <= v < x[i+1]的元素v的个数。

❷用histogram2()对通道0和通道2进行二维直方图统计。被统计的数组是两个一维数组，因此也需要用ravel()方法进行转换。它们分别为图像的通道0和通道2的数据。bins和range参数都变成了有两个元素的序列，分别与两个数组相对应。返回的统计结果hist2是一个二维数组，其形状由bins决定。第0轴与第一个数组相对应，第1轴与第二个数组相对应。它是由两个一维数组的对应元素所构成的二维矢量的分布统计结果。

观察图上图(左)可知红色通道的值普遍较大，因此整个图片呈现暖色调；而从右图不但可以得出红色通道的值比蓝色通道较大的结论，还可以看到几处分布比较密集的领域。例如其中一块的中心坐标大约为(207, 125)，这说明图像中红色值在207附近、蓝色值在125附近的像素点较多。

Opencv给我们提供的函数是cv2.calcHist()，它支持对多幅图像进行N维直方图统计，因此其第一个参数为数组列表。
该函数有5个参数：

image:输入图像，传入时应该用中括号[]括起，
channels:传入图像的通道，如果是灰度图像，那就不用说了，只有一个通道，值为0，如果是彩色图像（有3个通道），
                那么值为0,1,2,中选择一个，对应着  BGR各个通道。这个值也得用[]传入。
mask:掩膜图像。如果统计整幅图，那么为none。主要是如果要统计部分图的直方图，就得构造相应的炎掩膜来计算。
histSize:灰度级的个数，需要中括号，比如[256]
ranges:像素值的范围，通常[0,256]，有的图像如果不是0-256，比如说你来回各种变换导致像素值负值、很大，则需要调整后才可以。

3.2直方图反向映射

反向投影的概念：
反向投影是一种记录给定图像中的像素点如何适应直方图模型像素分布的方式，简单来讲，*反向投影就是首先计算某一特征的直方图模型，然后使用模型去寻找图像中存在的特征。*反向投影在某一位置的值就是原图对应位置像素值在原图像中的总数目。

反向投影作用——目标检测
计算出直方图之后，可以用calcBackProject()将图像中的每点替换为它在直方图中所对应的值。于是在直方图中出现次数越高，图像中对应的像素就越亮。可以用这种方法找出图像中和直方图相匹配的区域。下面用一个实际的例子加以说明.

注意：opencv中0-白色，255-黑色,蓝绿红，与matplotlib中0-黑，255-白，红绿蓝

import cv2
import numpy as np
import matplotlib.pyplot as pl
pl.rcParams["font.family"] = "SimHei"#直接修改配置字典，设置默认字体
#首先载入颜色匹配的模板图像
img = cv2.imread("E:/ruanjianDM/jupyternoerbookDM/fruits_section.jpg")
fig, axes = pl.subplots(2, 3, figsize=(9, 4))
a0,a1,a2,a3,a4,a5=axes[0,0],axes[0,1],axes[0,2],axes[1,0],axes[1,1],axes[1,2]
a0.set_title("模板图像")
a0.imshow(img[:,:,::-1])

#将GRB模型转化为HSV模型，H是色彩;S是深浅,S = 0时，只有灰度;V是明暗，表示色彩的明亮程度
img_hsv = cv2.cvtColor(img, cv2.COLOR_BGR2HSV)
#调用calcHist()对模板图像的色相与饱和度进行二维直方图计算。在色相与饱和度空间进行颜色匹配，能够得到较好的匹配结果。
result = cv2.calcHist([img_hsv], [0, 1], None, 
               [40, 40], [0, 256, 0, 256])
#为了后续的calcBackProject()计算不越界，这里将直方图的最大值缩小到255
result /= np.max(result) / 255
a2.set_title("直方图")
a2.imshow(result,cmap='binary')

#载入目标图像
img2 = cv2.imread("E:/ruanjianDM/jupyternoerbookDM/fruits.jpg") 
a1.set_title("目标图像")
a1.imshow(img2[:,:,::-1])

img_hsv2 = cv2.cvtColor(img2, cv2.COLOR_BGR2HSV)
#调用calcBackProject()将目标图像中的每个像素的颜色变换为其在直方图中所对应的值。它的第一个参数是一个图像
#列表，hist参数指定直方图，返回值是一幅单通道的形状和数值类型与输入图像相同的图像。
#channels、ranges参数和calcHist()的参数含义相同
img_bp = cv2.calcBackProject([img_hsv2], 
                        channels=[0, 1],
                        hist=result,
                        ranges=[0, 256, 0, 256],
                        scale=1)
a3.set_title("匹配结果")
a3.imshow(img_bp,cmap='binary')

#调用threshold()对calcBackProject()的输出图像进行二值化处理
#图像中值小于等于180的点都设置为0，将大于180的设置为255
_, img_th = cv2.threshold(img_bp, 180, 255, cv2.THRESH_BINARY)
#二值化之后的图像进行形态学图像处理
struct = np.ones((3, 3), np.uint8)
img_mp = cv2.morphologyEx(img_th, cv2.MORPH_CLOSE, struct, iterations=5) 
a4.set_title("匹配结果二值化")
a4.imshow(img_th,cmap='binary')
a5.set_title("二值化后的形态学处理")
a5.imshow(img_mp,cmap='binary')

3.3直方图匹配

直方图均衡化的思想是把原图的直方图变换为均匀分布的形式，这样就增加了像素灰度级值的动态范围，从而达到增强图像整体对比度的效果。但是对比度扩大到什么程度是不能进行控制的，给定衣服图像，直方图均衡化的结果是唯一的，它根据图像的灰度信息自动增强整个图像的对比度。简单操作不可控

直方图匹配又称直方图规定化，即变换原图的直方图为规定的某种形式的直方图，从而使两幅图像具有类似的色调和反差。直方图匹配属于非线性点运算。
原理：对两个直方图都做均衡化，变成相同的归一化的均匀直方图，以此均匀直方图为媒介，再对参考图像做均衡化的逆运算。可以控制达到预期的目标

#直方图匹配算法（目标图像，原始图像）
def histogram_match(src, dst):
    res = np.zeros_like(dst)
    #预存各个图像对象的三个通道的数据
    cdf_src = np.zeros((3, 256))
    cdf_dst = np.zeros((3, 256))
    cdf_res = np.zeros((3, 256))
    kw = dict(bins=256, range=(0, 256))

    for ch in (0, 1, 2):
        #计算原图像src与目标图像dst的归一化之后的直方图统计，得到的结果为概率密度分布
        hist_src, x = np.histogram(src[:, :, ch], **kw) 
        hist_dst, _ = np.histogram(dst[:, :, ch], **kw)
        hist_src=hist_src/np.sum(hist_src)
        hist_dst=hist_dst/np.sum(hist_dst)
        #用cumsum()对概率密度分布进行累加，得到累计分布
        cdf_src[ch] = np.cumsum(hist_src) 
        cdf_dst[ch] = np.cumsum(hist_dst)
        #在原图像的累计分布中搜索目标图像的累计分布所对应的下标index
        index = np.searchsorted(cdf_src[ch], cdf_dst[ch], side="left") 
        #调用clip()将index的取值范围限制在0～255之间
        np.clip(index, 0, 255, out=index) 
        #然后用它对目标图像进行映射，即将目标图像中每个像素值v替换为index[v]。
        res[:, :, ch] = index[dst[:, :, ch]] 
        hist_res, _ = np.histogram(res[:, :, ch], **kw)
        hist_res=hist_res/np.sum(hist_res)
        cdf_res[ch] = np.cumsum(hist_res)
        #返回原始图像，目标图像，匹配后图像的累计分布  
    return res, (cdf_src, cdf_dst, cdf_res),x

import cv2
import numpy as np
import matplotlib.pyplot as pl
pl.rcParams["font.family"] = "SimHei"#直接修改配置字典，设置默认字体
fig, axes = pl.subplots(2, 3, figsize=(10, 6))
a0,a1,a2,a3,a4,a5=axes[0,0],axes[0,1],axes[0,2],axes[1,0],axes[1,1],axes[1,2]
dst = cv2.imread("E:/ruanjianDM/jupyternoerbookDM/picture/summer.jpg")
src = cv2.imread("E:/ruanjianDM/jupyternoerbookDM/picture/autumn.jpg")
a0.set_title("(原始)秋季图像")
a0.imshow(src[:,:,::-1])
a0.get_xaxis().set_visible(False)
a0.get_yaxis().set_visible(False)
a1.set_title("(目标)夏季图像")
a1.imshow(dst[:,:,::-1])
a1.get_xaxis().set_visible(False)
a1.get_yaxis().set_visible(False)
#一幅秋景的直方图复制给夏景图像
res, cdfs,x = histogram_match(src, dst)
a2.set_title("(匹配)处理后的夏季图像")
a2.imshow(res[:,:,::-1])
a2.get_xaxis().set_visible(False)
a2.get_yaxis().set_visible(False)
a3.set_title("蓝通道累计分布")
a3.plot(0.5 * (x[:-1] + x[1:]),cdfs[1][0],color='b', linestyle=':')
a3.plot(0.5 * (x[:-1] + x[1:]),cdfs[2][0],color='y', linestyle='-')
a3.plot(0.5 * (x[:-1] + x[1:]),cdfs[0][0],color='k', linestyle='--')
a3.legend(['夏季','匹配后','秋季'],loc="lower right")
a4.set_title("绿通道累计分布")
a4.plot(0.5 * (x[:-1] + x[1:]),cdfs[1][1],color='b', linestyle=':')
a4.plot(0.5 * (x[:-1] + x[1:]),cdfs[2][1],color='y', linestyle='-')
a4.plot(0.5 * (x[:-1] + x[1:]),cdfs[0][1],color='k', linestyle='--')
a4.legend(['夏季','匹配后','秋季'],loc="lower right")
a5.set_title("红通道累计分布")
a5.plot(0.5 * (x[:-1] + x[1:]),cdfs[1][2],color='b', linestyle=':')
a5.plot(0.5 * (x[:-1] + x[1:]),cdfs[2][2],color='y', linestyle='-')
a5.plot(0.5 * (x[:-1] + x[1:]),cdfs[0][2],color='k', linestyle='--')
a5.legend(['夏季','匹配后','秋季'],loc="lower right")

将秋季的直方图复制给夏季图像，匹配后的夏季图像的色调跟秋季类似。通过下面的红绿蓝通道可知，匹配后的夏季图与秋季的直方图重合。

4.二维离散傅立叶变换

图像数据可以看作二维离散信号，对其进行二维离散傅立叶变换，能将其转换为频域信号，将原始图像分解为众多二维正弦波的叠加。由于NumPy已经提供了二维离散傅立叶变换的函数，因此本部分主要使用NumPy的相关函数进行说明。

对于一个NN的二维实数信号x进行二维快速傅立叶变换之后，得到表示频域信号的NN个复数元素的数组X。其中X[i,j]和X[N-i,N-j]共轭，并且X[0, 0]、X[0, N/2]、X[N/2, 0]、X[N/2, N/2]这4个元素的虚部为0。所以图像包含的数据量并没有增加。

from numpy import fft
x = np.random.rand(4, 4)
X = fft.fft2(x)
#allclose比较两个array是不是每一元素都相等
print(x)
print(X)
print(np.allclose(X[1:, 1:], X[3:0:-1, 3:0:-1].conj()))# 共轭复数
print(X[0,0],X[0,2],X[2,0],X[2,2])#虚部为0

x2 = fft.ifft2(X) # 将频域信号转换回空域信号,x2是一个复数数组，只是虚部接近于0
np.allclose(x, x2) # 和原始信号进行比较，返回结果为True

频域信号中的每个元素都对应空域信号中的一个二维正弦波，如果只选择频域信号中的一部分转换回空域信号，就相当于对空域信号进行了滤波处理。下面演示将频域信号中的不同区域转换回空域信号之后的滤波效果。

N=256
import cv2
import numpy as np
import matplotlib.pyplot as pl

pl.rcParams["font.family"] = "SimHei"#直接修改配置字典，设置默认字体
fig, axes = pl.subplots(2, 3, figsize=(10, 6))
a0,a1,a2,a3,a4,a5=axes[0,0],axes[0,1],axes[0,2],axes[1,0],axes[1,1],axes[1,2]

#首先载入一幅彩色图像，并将其转换为灰度图像。由于FFT运算的最佳大小为2的整
#数次幂，因此使用resize()将图像的大小改为256*256
img= cv2.imread("E:/ruanjianDM/jupyternoerbookDM/picture/LYF.jpg",0)
img = cv2.resize(img, (N, N))
a0.set_title("灰度图像")
a0.imshow(img,cmap='gray')

img_freq = fft.fft2(img)#得到频域信号
#为了能将频域信号作为图像显示，计算它的每个元素的模值，并取对数，得到数组img_mag
img_mag = np.log10(np.abs(img_freq))
a1.set_title("频域图像")
a1.imshow(img_mag,cmap='gray')

#fftshift()将两个对角线上的方块对调，即1、3象限对调，2、4象限对调。
#这样图像的中部与低频对应，而4角与高频信号对应。
img_mag_shift = fft.fftshift(img_mag)#将零频移到数组中间
a2.set_title("移位后频域图像")
a2.imshow(img_mag_shift,cmap='gray')

"下面是选择频域的一部分将其转换为空域信号"
rects = [(80, 125, 85, 130), (90, 90, 95, 95),
        (150, 10, 250, 250), (110, 110, 146, 146)]
filtered_results = []
for i, (x0, y0, x1, y1) in enumerate(rects):
    #mask是一个布尔数组，其形状和频域信号数组一样
    mask = np.zeros((N, N), dtype=np.bool) 
    #将其中坐标在指定的矩形范围之内的元素设置为True
    mask[x0:x1,y0:y1 ] = True 
    #同时选择共轭对称的部分，否则通过ifft2()转换回空域信号时虚部将不会为0
    mask[N - x1:N - x0 , N - y1:N - y0] = True 
    #通过fftshift()对mask数组进行移位，使得它和频域信号img_freq匹配
    mask = fft.fftshift(mask)#将正方形区域的True移到四周
    #对原始图像频域进行滤波之后的频域信号
    img_freq2 = img_freq * mask 
    #转换为空域信号
    img_filtered = fft.ifft2(img_freq2).real 
    filtered_results.append(img_filtered)
a3.set_title("过滤图像1")
a3.imshow(filtered_results[0],cmap='gray')
a4.set_title("过滤图像2")
a4.imshow(filtered_results[1],cmap='gray')
a5.set_title("过滤图像4")
a5.imshow(filtered_results[3],cmap='gray')

结论
1.从频域图像看：四个角较亮，对应低频信号（别忘了opencv中0-白色，255-黑色呦），中心对应高频信号
2.移位后频域图像：是将低频信号移到中心
3.过滤图像分别是过滤图像1,2,4是rects[0],rects[1],rects[4]过滤的频域图像对应的空间域图像

5.用双目视觉图像计算深度信息

import cv2
import numpy as np
import matplotlib.pyplot as pl
pl.rcParams["font.family"] = "SimHei"#直接修改配置字典，设置默认字体
img_left = cv2.pyrDown(cv2.imread('E:/ruanjianDM/jupyternoerbookDM/picture/aloeL.jpg'))
img_right = cv2.pyrDown(cv2.imread('E:/ruanjianDM/jupyternoerbookDM/picture/aloeR.jpg'))

img_left = cv2.cvtColor(img_left, cv2.COLOR_BGR2RGB)
img_right = cv2.cvtColor(img_right, cv2.COLOR_BGR2RGB)

stereo_parameters = dict(
        SADWindowSize = 5,
        numDisparities = 192,
        preFilterCap = 4,
        minDisparity = -24,
        uniquenessRatio = 1,
        speckleWindowSize = 150,
        speckleRange = 2,
        disp12MaxDiff = 10,
        fullDP = False,
        P1 = 600,
        P2 = 2400)
stereo = cv2.StereoSGBM(**stereo_parameters)
#StereoSGBM.compute()计算视差信息,视差越大，目标点到照相机的距离越小
disparity = stereo.compute(img_left, img_right).astype(np.float32) / 16

h, w = img_left.shape[:2]
ygrid, xgrid = np.mgrid[:h, :w]
ygrid = ygrid.astype(np.float32)
xgrid = xgrid.astype(np.float32)
res = cv2.remap(img_right, xgrid - disparity, ygrid, cv2.INTER_LINEAR)

fig, axes = pl.subplots(1, 3, figsize=(9, 3))
axes[0].imshow(img_left)
axes[0].imshow(img_right, alpha=0.5)
axes[1].imshow(disparity, cmap="gray")
axes[2].imshow(img_left)
axes[2].imshow(res, alpha=0.5)
for ax in axes:
    ax.axis("off")
fig.subplots_adjust(0, 0, 1, 1, 0, 0)

上面程序没有运行出来，先放在上面