计算机视觉-图像到图像的映射

一、单应性变换

单应性变换是将一个平面内的点映射到另一个平面内的二维投影变换。这里的平面是指图像或者三维中的平面表示。单应性变换具有很强的实用性，比如图像配准，图像纠正和纹理扭曲，以及创建全景图像，我们将频繁的使用单应性变换。本质上，单应性变换H，按照下面的方程映射二维中的点（齐次坐标意义下）：

点的齐次坐标是依赖于其尺度定义的，所以，x=[x,y,w]=[ax,ay,aw]=[x/w,y/w,1]都表示同一个二维点。因此，单应性矩阵H也仅依赖尺度定义，所以，单应性矩阵具有8个独立的自由度。我们通常使用w=1来归一化点，这样，点具有唯一的图像坐标x和y。这个额外的坐标是的我们可以简单地使用一个矩阵来表示变换。

进行点和变换的处理时。我们会按照列优先的原则存储这些点。因此，n个二维点集将会存储为齐次坐标意义下的一个3×n数组。对于聚类和分类的特征，我们一般使用典型的行数组来存储数据。

 原理：

一个完全射影变换具有8个自由度。根据对应点约束，每个对应点对可以写出两个方程，分别对应于x和y坐标。因此，计算单应性矩阵H需要4个对应点对。DLT（Direct Linear Transformation，直接线性变换）是给定4个点或者更多对应点对矩阵，来计算单应性矩阵H的算法。

 代码：

def H_from_points(fp, tp):
    """使用线性DLT方法，计算单应性矩阵H，使fp映射到tp。点自动进行归一化"""

    if fp.shape != tp.shape:
        raise RuntimeError('number of points do not match')

    # 对点进行归一化（对数值计算很重要）
    # --- 映射起始点 ---
    m = mean(fp[:2], axis=1)
    maxstd = max(std(fp[:2], axis=1)) + 1e-9
    C1 = diag([1/maxstd, 1/maxstd, 1])
    C1[0][2] = -m[0]/maxstd
    C1[1][2] = -m[1]/maxstd
    fp = dot(C1,fp)
    
    # --- 映射对应点 ---
    m = mean(tp[:2], axis=1)
    maxstd = max(std(tp[:2], axis=1)) + 1e-9
    C2 = diag([1 / maxstd, 1 / maxstd, 1])
    C2[0][2] = -m[0] / maxstd
    C2[1][2] = -m[1] / maxstd
    tp = dot(C2, tp)
    
    # 创建用于线性方法的矩阵，对于每个对应对，在矩阵中会出现两行数值
    nbr_correspondences = fp.shape[1]
    A = zeros((2 * nbr_correspondences, 9))
    for i in range(nbr_correspondences):
        A[2*i] = [-fp[0][i], -fp[1][i],-1,0,0,0,
                  tp[0][i]*fp[0][i],tp[0][i]*fp[1][i],tp[0][i]]
        A[2*i+1] = [0,0,0,-fp[0][i],-fp[1][i],-1,
                    tp[1][i]*fp[0][i],tp[1][i]*fp[1][i],tp[1][i]]
        
    U,S,V = linalg.svd(A)
    H = V[8].reshape((3,3))
    
    #反归一化
    H = dot(linalg.inv(C2),dot(H,C1))
    
    #归一化，然后返回
    return H / H[2,2]

结果分析：

先对这些点进行归一化操作，使其均值为0，方差为1。使用对应点对来构造矩阵A。最小二乘解即为矩阵SVD分解后所得矩阵V的最后一行。该行经过变换后得到矩阵H，然后对这个矩阵进行处理和归一化，返回输出。

 

二、仿射变换

仿射变换具有6个自由度，因此我们需要三个对应点来估计矩阵H。通过将最后两个元素设置为0，即h~7~=h~8~=0，DLT算法与单应性变换同。

代码：

def Haffine_from_points(fp, tp):
    """计算H仿射变换，使得tp是fp经过仿射变换H得到的"""

    if fp.shape != tp.shape:
        raise RuntimeError('number of points do not match')
        
    # 对点进行归一化（对数值计算很重要）
    # --- 映射起始点 ---
    m = mean(fp[:2], axis=1)
    maxstd = max(std(fp[:2], axis=1)) + 1e-9
    C1 = diag([1/maxstd, 1/maxstd, 1]) 
    C1[0][2] = -m[0]/maxstd
    C1[1][2] = -m[1]/maxstd
    fp_cond = dot(C1,fp)
    
    # --- 映射对应点 ---
    m = mean(tp[:2], axis=1)
    C2 = C1.copy() # 两个点集，必须都进行相同的缩放
    C2[0][2] = -m[0]/maxstd
    C2[1][2] = -m[1]/maxstd
    tp_cond = dot(C2,tp)
    
    # 因为归一化后点的均值为0，所以平移量为0
    A = concatenate((fp_cond[:2],tp_cond[:2]), axis=0)
    U,S,V = linalg.svd(A.T)
    
    # 如Hartley和Zisserman著的Multiplr View Geometry In Computer,Scond Edition所示，
    # 创建矩阵B和C
    tmp = V[:2].T
    B = tmp[:2]
    C = tmp[2:4]
    
    tmp2 = concatenate((dot(C,linalg.pinv(B)),zeros((2,1))), axis=1) 
    H = vstack((tmp2,[0,0,1]))
    
    # 反归一化
    H = dot(linalg.inv(C2),dot(H,C1))
    
    return H / H[2,2]

 

三、图片扭曲（仿射扭曲）

原理：

仿射变换的原理就是进行空间的点坐标的变换，即从二维坐标到二维坐标之间的线性变换，且保持二维图形的“平直性”和“平行性”。其参数公式为：

由于仿射变换具有6个自由度，所以需要三个对应点来估计矩阵H。其中X和Y为原坐标，X’和Y’为变换后的坐标，a~f为6个自由度，其中c、f参数分别为x方向和y方向的平移向量。

阿尔法通道是一个8位的灰度通道，该通道用256级灰度来记录图像中的透明度信息，定义透明、不透明和半透明区域。如果一个像素的alpha通道数值为0，那它就是完全透明的（也就是看不见的），而数值为1则意味着一个完全不透明的像素（传统的数字图像）。我们常说的RGBA图像就是普通的RGB加上这个α通道。说白了就是控制图像透明度，在下面的代码段中：
tp = array([[264,538,540,264],[40,36,605,605],[1,1,1,1]])的最后四个1就表示四个角点的透明度为不透明，以此实现图像的完全覆盖。

代码：

from array import array

from numpy import *
from matplotlib.pyplot import *
from scipy import ndimage
from PIL import Image
im = array(Image.open(r'C:\Users\yifan\Pictures\cat.jpg').convert('L'))
H = array([[1.4,0.05,-100],[0.05,1.5,-100],[0,0,1]])
im2 = ndimage.affine_transform(im, H[:2,:2],(H[0,2],H[1,2]))

gray()
subplot(121)
imshow(im)
axis('off')
subplot(122)
imshow(im2)
axis('off')
show()

结果：

 

四、图像的映射融合(图像中的图像)

图像扭曲，就是把一个图像的角点映射在另一个图像的某个位置，将其覆盖的像素值替换，合成一幅图像。

代码：

from array import array
from mailcap import show
from tkinter import Image

from PCV.geometry.homography import Haffine_from_points
from matplotlib.pyplot import axis, imshow, subplot, gray
from scipy import ndimage


def image_in_image(im1, im2, tp):
    """使用仿射变换将im1放置在im2上，使im1图像的角和tp尽可能的靠近
        tp是齐次表示的，并且是按照从左上角逆时针计算的"""

    # 扭曲的点
    m, n = im1.shape[:2]
    fp = array([[0, m, m, 0], [0, 0, n, n], [1, 1, 1, 1]])

    # 计算仿射变换，并且将其应用于图像im1中
    H = Haffine_from_points(tp, fp)
    im1_t = ndimage.affine_transform(im1, H[:2, :2],
                                     (H[0, 2], H[1, 2]), im2.shape[:2])
    alpha = (im1_t > 0)

    return (1 - alpha) * im2 + alpha * im1_t


im1 = array(Image.open(r'C:\Users\yifan\Pictures\flower.png').convert('L'))
im2 = array(Image.open(r'C:\Users\yifan\Pictures\cat.jpg').convert('L'))

gray()
subplot(131)
imshow(im1)
axis('equal')
axis('off')
subplot(132)
imshow(im2)
axis('equal')
axis('off')

# 选定一些目标点
tp = array([[264, 538, 540, 264], [40, 36, 605, 605], [1, 1, 1, 1]])

im3 = image_in_image(im1, im2, tp)
subplot(133)
imshow(im3)
axis('equal')
axis('off')
show()

结果：

 

结果分析：

图一为flower的灰度图像，图二为cat的灰度图像，通过修改tp = array([[264,538,540,264],[40,36,605,605],[1,1,1,1]])这一行的坐标点可以改变图一映射在图二的位置，其中前四个数字代表四个角点的纵坐标，中间四个数字代表四个角点的横坐标，四个角点的顺序为从左上角开始按照逆时针方向排序，最后四个数字代表四个角点的α通道。