In [97]:
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")

Project 5A: Image Warping & Mosaicing

Shoot the Pictures

In [241]:
# 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)

Recovering Homographies

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$

In [720]:
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)
In [733]:
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'))
In [734]:
H = computeH(left_p, mid_p)
H_inv = np.linalg.inv(H)
H_inv = H_inv / np.sqrt(np.sum(H_inv**2))

Correspondence Points

I'm using grayscale images to simplify things for now, but I plan to return to the colored images for part (b)

In [561]:
matplotlib.use('module://ipykernel.pylab.backend_inline')

plt.imshow(left, cmap='gray')
plt.scatter(left_p[:, 0], left_p[:, 1], color='r')
Out[561]:
<matplotlib.collections.PathCollection at 0x27891b7d940>
In [562]:
plt.imshow(mid, cmap='gray')
plt.scatter(mid_p[:, 0], mid_p[:, 1], color='r')
Out[562]:
<matplotlib.collections.PathCollection at 0x27892169d90>

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.

In [736]:
from skimage.transform import ProjectiveTransform
warped = transform.warp(left, ProjectiveTransform(matrix=H_inv))

Left image warped to correspond to center image:

In [737]:
plt.imshow(warped, cmap='gray')
Out[737]:
<matplotlib.image.AxesImage at 0x2789e0d7d00>
In [738]:
# 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:

In [739]:
plt.imshow(dst, cmap='gray')
Out[739]:
<matplotlib.image.AxesImage at 0x2789e132e20>
In [761]:
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:

In [743]:
frame = color.rgb2gray(plt.imread('frame.jpg'))
plt.imshow(frame, cmap='gray')
Out[743]:
<matplotlib.image.AxesImage at 0x2789e238e50>
In [746]:
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.

In [753]:
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')
Out[753]:
<matplotlib.collections.PathCollection at 0x2789e508100>
In [764]:
H = computeH(pts, new_pts)
plt.imshow(warpImage(frame, H), cmap='gray')
Out[764]:
<matplotlib.image.AxesImage at 0x2789e6c5ac0>

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.

In [ ]:
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()