OpenCV-Python|图像处理模块 — cv2.remap()函数的使用

@[TOC](OpenCV-Python|图像处理模块 — cv2.remap()函数的使用)

前言

OpenCV-Python|图像处理模块 — cv2.remap()函数的使用

原理

什么是重映射？

• 重映射是从图像中的一个位置获取像素并将其放置在新图像中的另一位置的过程。
• 为了完成映射过程，可能有必要对非整数像素位置进行一些插值，因为源图像和目标图像之间并不总是一对一的像素对应关系。
• 对每个像素位置$(x,y)$，重映射可以表示为$$g(x,y)=f(h(x,y))$$其中$g()$是重映射图，$f()$是源图，$h(x,y)$是在坐标$(x,y)$上进行操作的映射函数。
• 举一个简单的例子，有一张图片$I$以及一个重映射$$h(x,y)=(I.cols-x,y)$$那么图像将会在$x$方向进行翻转。如下图所示：
  
  可以看出红圈相对于$x$水平方向改变了位置。
  

重映射的简单应用

现在我们来学习重映射的简单应用。
假设我们有两幅图像img1和img2，称为左图和右图，并且通过特征提取与匹配得到了img1到img2的投影变换矩阵H

我们可以通过cv.remap()函数来将img2映射到img1对应位置上并合成：


这个过程是怎么实现的呢？

首先，通过单应性矩阵H可以计算出右图映射在左图上的四个顶点，由于H是左图到右图的单应性，所以右图映射到左图就是H \ (x,y,1)，(x,y,1)表示右图的顶点。

TL = np.linalg.solve(H, np.array([0,0,1]))
TL = np.round(TL/TL[-1])
BL = np.linalg.solve(H, np.array([0,img2.shape[0]-1,1]))
BL = np.round(BL/BL[-1])
TR = np.linalg.solve(H, np.array([img2.shape[1]-1,0,1]))
TR = np.round(TR/TR[-1])
BR = np.linalg.solve(H, np.array([img2.shape[1]-1,img2.shape[0]-1,1]))
BR = np.round(BR/BR[-1])

接着可以求出右图映射到左图的横纵坐标的最小最大值：

u0_im_ = int(min(TL[0], BL[0], TR[0], BR[0]));   u1_im_ = int(max(TL[0], BL[0], TR[0], BR[0]))
v0_im_ = int(min(TL[1], BL[1], TR[1], BR[1]));   v1_im_ = int(max(TL[1], BL[1], TR[1], BR[1]))
print(u0_im_, u1_im_, v0_im_, v1_im_)

# 计算结果：
# 281 938 -14 488

计算出的结果是以img1图像坐标系为参考坐标系的。

然后可以求出拼接画布的尺寸：

u0 = min(0, u0_im_)
u1 = max(img1.shape[1]-1, u1_im_)
ur = np.arange(u0, u1 + 1)
v0 = min(0, v0_im_)
v1 = max(img1.shape[0]-1, v1_im_)
vr = np.arange(v0, v1 + 1)
cw = u1 - u0 + 1
ch = v1 - v0 + 1
print(u0, u1, v0, v1, ch, cw)

# 计算结果：
# 0 938 -14 488 503 939

计算出的结果是同样是以img1图像坐标系为参考坐标系的。
此时，我们得到了在img1图像坐标系下拼接画布的x和y坐标范围，即ur和vr，通过np.meshgrid()函数生成二维坐标网格。

u, v = np.meshgrid(ur, vr)

下面就是cv.remap()函数应用的地方。对于img1，img1到img1的单应性矩阵是单位矩阵，所以映射矩阵map_x和map_y就是u，v：

u = np.float32(u);  v = np.float32(v)   # remap函数要求映射矩阵为CV_32F
warped_img1 = cv2.remap(img1, u, v, cv2.INTER_LINEAR, borderMode=cv2.BORDER_REFLECT_101)
mask1 = np.ones((img1.shape[0],img1.shape[1]))
warped_mask1 = cv2.remap(mask1, u, v, cv2.INTER_LINEAR)

对于img2，img1到img2经过了投影变换H，所以map_x和map_y是用H * u 和 H * v并归一化：

z_ = H[2,0]*u + H[2,1]*v + H[2,2]
map_x = (H[0,0]*u + H[0,1]*v + H[0,2]) / z_
map_y = (H[1,0]*u + H[1,1]*v + H[1,2]) / z_
map_x = np.float32(map_x);  map_y = np.float32(map_y)
warped_img2 = cv2.remap(img2, map_x, map_y, cv2.INTER_LINEAR, borderMode=cv2.BORDER_REFLECT_101)
mask2 = np.ones((img2.shape[0],img2.shape[1]))
warped_mask2 = cv2.remap(mask2, map_x, map_y, cv2.INTER_LINEAR)

最后，通过掩码操作就可以将两幅图像合成。

完整代码：

import numpy as np
import cv2

# read img1 and img2
img1 = cv2.imread('yosemite1.jpg')
img2 = cv2.imread('yosemite2.jpg')
cv2.imshow('img', np.concatenate((img1,img2),axis=1))
cv2.waitKey(1)

# Feature extraction and matching
ft_detector = cv2.SIFT_create()
keyPoints1, descriptors1 = ft_detector.detectAndCompute(img1, None)
keyPoints2, descriptors2 = ft_detector.detectAndCompute(img2, None)

bf = cv2.BFMatcher(crossCheck=False)
matches = bf.match(descriptors1, descriptors2)
matches = sorted(matches, key = lambda x:x.distance)
sourcePoints = np.float32([ keyPoints1[m.queryIdx].pt for m in matches ]).reshape(-1, 1, 2)
destinationPoints = np.float32([ keyPoints2[m.trainIdx].pt for m in matches ]).reshape(-1, 1, 2)
# Obtain homography
H, _ = cv2.findHomography(sourcePoints, destinationPoints, method=cv2.RANSAC, ransacReprojThreshold=5.0)
print(H)

# 映射右图的四个顶点
TL = np.linalg.solve(H, np.array([0,0,1]))
TL = np.round(TL/TL[-1])
BL = np.linalg.solve(H, np.array([0,img2.shape[0]-1,1]))
BL = np.round(BL/BL[-1])
TR = np.linalg.solve(H, np.array([img2.shape[1]-1,0,1]))
TR = np.round(TR/TR[-1])
BR = np.linalg.solve(H, np.array([img2.shape[1]-1,img2.shape[0]-1,1]))
BR = np.round(BR/BR[-1])

# img2映射后的坐标范围
u0_im_ = int(min(TL[0], BL[0], TR[0], BR[0]));   u1_im_ = int(max(TL[0], BL[0], TR[0], BR[0]))
v0_im_ = int(min(TL[1], BL[1], TR[1], BR[1]));   v1_im_ = int(max(TL[1], BL[1], TR[1], BR[1]))
print(u0_im_, u1_im_, v0_im_, v1_im_)

# 拼接画布的尺寸
u0 = min(0, u0_im_)
u1 = max(img1.shape[1]-1, u1_im_)
ur = np.arange(u0, u1 + 1)
v0 = min(0, v0_im_)
v1 = max(img1.shape[0]-1, v1_im_)
vr = np.arange(v0, v1 + 1)
cw = u1 - u0 + 1
ch = v1 - v0 + 1
print(u0, u1, v0, v1, ch, cw)

u, v = np.meshgrid(ur, vr)

u = np.float32(u);  v = np.float32(v)   # remap函数要求映射矩阵为CV_32F
warped_img1 = cv2.remap(img1, u, v, cv2.INTER_LINEAR, borderMode=cv2.BORDER_REFLECT_101)
mask1 = np.ones((img1.shape[0],img1.shape[1]))
warped_mask1 = cv2.remap(mask1, u, v, cv2.INTER_LINEAR)

z_ = H[2,0]*u + H[2,1]*v + H[2,2]
map_x = (H[0,0]*u + H[0,1]*v + H[0,2]) / z_
map_y = (H[1,0]*u + H[1,1]*v + H[1,2]) / z_
map_x = np.float32(map_x);  map_y = np.float32(map_y)
warped_img2 = cv2.remap(img2, map_x, map_y, cv2.INTER_LINEAR, borderMode=cv2.BORDER_REFLECT_101)
mask2 = np.ones((img2.shape[0],img2.shape[1]))
warped_mask2 = cv2.remap(mask2, map_x, map_y, cv2.INTER_LINEAR)

mass = warped_mask1 + warped_mask2
mass[mass==0] = np.nan
output = np.zeros_like(warped_img1)
for c in range(3):
    output[:,:,c] = (warped_img1[:,:,c] * warped_mask1 + warped_img2[:,:,c] * warped_mask2) / mass

cv2.imshow('warped_img1', np.uint8(warped_img1 * warped_mask1[..., np.newaxis].repeat(3, axis=-1)))
cv2.imshow('warped_img2', np.uint8(warped_img2 * warped_mask2[..., np.newaxis].repeat(3, axis=-1)))
cv2.imshow('output_img', output)
cv2.waitKey(0)
cv2.destroyAllWindows()