from skimage import color, transform
import skimage.io as skio
import numpy as np
import os
import numpy as np
import matplotlib
import matplotlib.pyplot as plt
import seaborn as sns
import pickle
import cv2
# Ignore warnings
import warnings
warnings.filterwarnings("ignore")
# matplotlib.use('module://ipykernel.pylab.backend_inline')
# matplotlib.use('TkAgg')
# img = plt.imread('car_mid.jpg')
# plt.imshow(img)
# pts = plt.ginput(n=8, timeout=0)
with open('car_mid.p', 'wb') as fp:
pickle.dump(pts, fp)
The homography transform between a set of corresponding points between two images is defined by:
$p' = Hp$, where $H$ is a $3x3$ matrix with 8 unknowns. We only need four equations (pairs of correspondence points) to recover the 8 unknowns, but the recovery may be unstable and prone to noise. Alternatively, we can solve an overdetermined system by providing more than 4 correspondences and rearranging $H$ into a $9x1$ vector. Our new system is $PH = \vec{0}$, where $P$ is the $Nx9$ matrix composed by stacking the following submatrix $p_i$ for each correspondence pair (source: http://6.869.csail.mit.edu/fa12/lectures/lecture13ransac/lecture13ransac.pdf)
\begin{equation*} p_i = \begin{bmatrix} -x_i & -y_i & -1 & 0 & 0 & 0 & x_ix_i' & y_ix_i' & x_i' \\ 0 & 0 & 0 & -x_i & -y_i & -1 & x_iy_i' & y_iy_i' & y_i' \end{bmatrix} \end{equation*}and
\begin{equation*} H = \begin{bmatrix} h_1 \\ h_2 \\ h_3 \\ h_4 \\ h_5 \\ h_6 \\ h_7 \\ h_8 \\ h_9 \end{bmatrix} \end{equation*}where the least squares solution to $H$ is the right singular vector corresponding to the smallest eigenvalue of the Singular Value Decomposition of $P$
def computeH(im1_pts, im2_pts):
P = np.empty((0, 9))
for i in range(len(im1_pts)):
im1 = im1_pts[i]
im2 = im2_pts[i]
p_1 = np.array([-im1[0], -im1[1], -1, 0, 0, 0, im1[0]*im2[0], im1[1]*im2[0], im2[0]])
p_2 = np.array([0, 0, 0, -im1[0], -im1[1], -1, im1[0]*im2[1], im1[1]*im2[1], im2[1]])
P = np.vstack((P, p_1.reshape(1, 9)))
P = np.vstack((P, p_2.reshape(1, 9)))
return np.linalg.svd(P)[2][-1].reshape(3, 3)
with open('car_left.p', 'rb') as fp:
left_p = pickle.load(fp)
with open('car_mid.p', 'rb') as fp:
mid_p = pickle.load(fp)
left = color.rgb2gray(plt.imread('car_left.jpg'))
mid = color.rgb2gray(plt.imread('car_mid.jpg'))
H = computeH(left_p, mid_p)
H_inv = np.linalg.inv(H)
H_inv = H_inv / np.sqrt(np.sum(H_inv**2))
I'm using grayscale images to simplify things for now, but I plan to return to the colored images for part (b)
matplotlib.use('module://ipykernel.pylab.backend_inline')
plt.imshow(left, cmap='gray')
plt.scatter(left_p[:, 0], left_p[:, 1], color='r')
plt.imshow(mid, cmap='gray')
plt.scatter(mid_p[:, 0], mid_p[:, 1], color='r')
Here, i'm using skimage.transform
as a sanity check that my homography matrix is correct. My goal is to transform the left-most image to correspond with the middle image.
from skimage.transform import ProjectiveTransform
warped = transform.warp(left, ProjectiveTransform(matrix=H_inv))
Left image warped to correspond to center image:
plt.imshow(warped, cmap='gray')
# Code borrowed from: https://stackoverflow.com/questions/46520123/how-do-i-use-opencvs-remap-function
h, w = mid.shape[:2]
indy, indx = np.indices((h, w), dtype=np.float32)
lin_homg_ind = np.array([indx.ravel(), indy.ravel(), np.ones_like(indx).ravel()])
# warp the coordinates of src to those of true_dst
map_ind = H_inv @ lin_homg_ind
map_x, map_y = map_ind[:-1]/map_ind[-1] # ens|ure homogeneity
map_x = map_x.reshape(h, w).astype(np.float32)
map_y = map_y.reshape(h, w).astype(np.float32)
# remap!
dst = cv2.remap(left, map_x, map_y, cv2.INTER_LINEAR)
Rectified (left) image produced with linear interpolation using cv2.map()
and the inverse Homography transform:
plt.imshow(dst, cmap='gray')
def warpImage(im, H):
H_inv = np.linalg.inv(H)
H_inv = H_inv / np.sqrt(np.sum(H_inv**2))
h, w = im.shape[:2]
indy, indx = np.indices((h, w), dtype=np.float32)
lin_homg_ind = np.array([indx.ravel(), indy.ravel(), np.ones_like(indx).ravel()])
# warp the coordinates of src to those of true_dst
map_ind = H_inv @ lin_homg_ind
map_x, map_y = map_ind[:-1]/map_ind[-1] # ens|ure homogeneity
map_x = map_x.reshape(h, w).astype(np.float32)
map_y = map_y.reshape(h, w).astype(np.float32)
# remap!
dst = cv2.remap(im, map_x, map_y, cv2.INTER_LINEAR)
return dst
Let's try to rectify this crooked frame to be front parallel:
frame = color.rgb2gray(plt.imread('frame.jpg'))
plt.imshow(frame, cmap='gray')
matplotlib.use('module://ipykernel.pylab.backend_inline')
matplotlib.use('TkAgg')
plt.imshow(frame)
pts = plt.ginput(n=4, timeout=0)
with open('frame.p', 'wb') as fp:
pickle.dump(pts, fp)
True correspondence points shown in blue, front-aligned correspondence points shown in red.
matplotlib.use('module://ipykernel.pylab.backend_inline')
with open('frame.p', 'rb') as fp:
frame_p = np.array(pickle.load(fp))
new_pts = np.array([(500, 500), (2500, 500), (500, 3500), (2500, 3500)])
plt.imshow(frame, cmap='gray')
plt.scatter(frame_p[:, 0], frame_p[:, 1], s=50)
plt.scatter(new_pts[:, 0], new_pts[:, 1], s=50, color='r')
H = computeH(pts, new_pts)
plt.imshow(warpImage(frame, H), cmap='gray')
It seems to have flipped the image around, but the image rectification algorithm seems functional other than that. This should easily be fixable by defining more correspondence points.
FILE = "main.html"
with open(FILE, 'r') as html_file:
content = html_file.read()
# Get rid off prompts and source code
content = content.replace("div.input_area {
display: none;","div.input_area {
display: none;\n\tdisplay: none;")
content = content.replace(".prompt {
display: none;",".prompt {
display: none;\n\tdisplay: none;")
f = open(FILE, 'w')
f.write(content)
f.close()