目录
图像到图像的映射
1、单应性变换
单应性变换是将一个平面内的点映射到另一个平面内的二维投影变换。单应性变换具有很强的实用性,比如图像配准、图像纠正和纹理扭曲,以及创建全景图像。
单应性矩阵具有H具有8个独立自由度。通常使用w=1来归一化点,这样,点具有唯一的图像坐标x和y。
创建homography.py文件,该函数可以实现兑点进行归一化和转换齐次坐标:
def normalize(points):
for row in points:
row /= points[-1]
return points
def make_homog(points):
return vstack((points,ones((1,points.shape[1]))))
1.1直接线性变换算法
单应性矩阵可以由两幅图像中对应点对计算出来。根据对应点对约束,每个对应点可以写出两个方程,分别对应x和y坐标。因此,计算单应性矩阵H需要4个对应点对。
使用SVD算法找到H的最小二乘解,添加H_from_points(fp,tp)函数到homography.py文件中,代码如下所示:
def H_from_points(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]
1.2仿射变换
仿射变换具有6个自由度,需要三个对应点来估计矩阵H。通过将最后两个元素设置为0,即h21、h22设置为0。仿射变换可以通过上述的DLT算法估算出来。
添加Haffine_from_points(fp,tp)函数到homograph.py文件中,代码如下所示:
def Haffine_from_points(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_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)
A = concatenate((fp_cond[:2], tp_cond[:2]), axis=0)
U, S, V = linalg.svd(A.T)
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]
homography.py
from numpy import *
from scipy import ndimage
def normalize(points):
for row in points:
row /= points[-1]
return points
def make_homog(points):
return vstack((points,ones((1,points.shape[1]))))
def H_from_points(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]
def Haffine_from_points(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_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)
A = concatenate((fp_cond[:2], tp_cond[:2]), axis=0)
U, S, V = linalg.svd(A.T)
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]
2、图像扭曲
对图像块应用仿射变换,将其称为图像扭曲。该操作不仅应用在计算机图形学中,而且经常出现在计算机视觉算法中。扭曲操作可以使用Scipy工具包中的ndimage包来实现。命令如下:
transformed_im = ndimage.affine_transform(im,A,b,size)
代码如下:
from scipy import ndimage
from PIL import Image
from pylab import *
from scipy import ndimage
im = array(Image.open('113.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]))
figure()
gray()
imshow(im2)
show()
运行结果:
2.1图像中的图像
将图片或者图片的一部分放置在另一幅图片中,使得它们能够和指定的区域或者标记物对齐。
将image_in_image()函数添加到warp.py文件中,该函数使用仿射变换将im1放置在im2上,使im1图像的角和tp尽可能的靠近。tp是齐次表示的,并且是按照从左上角逆时针计算的。
warp.py
# -*- coding: utf-8 -*-
# import matplotlib.delaunay as md
from scipy.spatial import Delaunay
from scipy import ndimage
from pylab import *
from numpy import *
from PCV.geometry import homography
def image_in_image(im1, im2, tp):
m, n = im1.shape[:2]
fp = array([[0, m, m, 0], [0, 0, n, n], [1, 1, 1, 1]])
# n是im2投影到im1的图像的长,m是宽
# im1四个顶点的坐标
# 得到单应性矩阵
H = homography.Haffine_from_points(tp, fp)
# 两个图形的变换关系(位移+缩放),只取H的前两行
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
使用仿射将一幅图放置在另一幅图中
import warp
from pylab import *
from numpy import *
from PIL import Image
im1 = array(Image.open('bao.jpg').convert('L'))
im2 = array(Image.open('113.jpg').convert('L'))
#tp = array([[200,460,475,210],[20,16,569,546],[1,1,1,1]])
#tp = array([[369,598,575,326],[512,481,803,782],[1,1,1,1]])
#tp = array([[264,538,540,264],[40,36,605,605],[1,1,1,1]])
tp = array([[564,638,690,464],[36,24,468,455],[1,1,1,1]])
#tp = array([[120,138,340,264],[116,104,268,150],[1,1,1,1]])
im3 = warp.image_in_image(im1,im2,tp)
figure()
gray()
#imshow(im1)
imshow(im2)
#imshow(im3)
axis('equal')
axis('off')
show()