Skip to content
Snippets Groups Projects
Commit 1775e87b authored by vjshah3's avatar vjshah3
Browse files

add MP3

parent 0ae697f1
No related branches found
No related tags found
No related merge requests found
Showing
with 674 additions and 0 deletions
File added
%% Cell type:markdown id: tags:
# Image Stitching (Python)
## Usage
This code snippet provides an overall code structure and some interactive plot interfaces for the Stitching Pairs of Images section of MP 3. In main function, we outline the required functionalities step by step. Feel free to make modifications on the starter code if it's necessary.
## Package installation
- `opencv`
- `numpy`
- `skimage`
- `scipy`
%% Cell type:markdown id: tags:
# Common imports
%% Cell type:code id: tags:
``` python
import numpy as np
import skimage
import cv2
import matplotlib.pyplot as plt
from scipy.spatial import distance
import scipy
```
%% Cell type:markdown id: tags:
# Helper functions
%% Cell type:code id: tags:
``` python
def imread(fname):
"""
read image into np array from file
"""
return skimage.io.imread(fname)
def imread_bw(fname):
"""
read image as gray scale format
"""
return cv2.cvtColor(imread(fname), cv2.COLOR_BGR2GRAY)
def imshow(img):
"""
show image
"""
skimage.io.imshow(img)
def get_sift_data(img):
"""
detect the keypoints and compute their SIFT descriptors with opencv library
"""
sift = cv2.SIFT_create()
kp, des = sift.detectAndCompute(img, None)
return kp, des
def plot_inlier_matches(ax, img1, img2, inliers):
"""
plot the match between two image according to the matched keypoints
:param ax: plot handle
:param img1: left image
:param img2: right image
:inliers: x,y in the first image and x,y in the second image (Nx4)
"""
res = np.hstack([img1, img2])
ax.set_aspect('equal')
ax.imshow(res, cmap='gray')
ax.plot(inliers[:,0], inliers[:,1], '+r')
ax.plot(inliers[:,2] + img1.shape[1], inliers[:,3], '+r')
ax.plot([inliers[:,0], inliers[:,2] + img1.shape[1]],
[inliers[:,1], inliers[:,3]], 'r', linewidth=0.4)
ax.axis('off')
```
%% Cell type:markdown id: tags:
# Your implementations
%% Cell type:code id: tags:
``` python
def get_best_matches(img1, img2, num_matches):
kp1, des1 = get_sift_data(img1)
kp2, des2 = get_sift_data(img2)
kp1, kp2 = np.array(kp1), np.array(kp2)
# Find distance between descriptors in images
dist = scipy.spatial.distance.cdist(des1, des2, 'sqeuclidean')
# Write your code to get the matches according to dist
# <YOUR CODE>
pass
def ransac(...):
"""
write your ransac code to find the best model, inliers, and residuals
"""
# <YOUR CODE>
pass
def compute_homography(...):
"""
write your code to compute homography according to the matches
"""
# <YOUR CODE>
pass
def warp_images(...):
"""
write your code to stitch images together according to the homography
"""
# <YOUR CODE>
pass
```
%% Cell type:markdown id: tags:
# Main functions
%% Cell type:markdown id: tags:
#### Load images
%% Cell type:code id: tags:
``` python
img1 = imread('./stitch/left.jpg')
img2 = imread('./stitch/right.jpg')
```
%% Cell type:markdown id: tags:
#### Part (3) compute and display the initial SIFT matching result
%% Cell type:code id: tags:
``` python
data = get_best_matches(img1, img2, 300)
fig, ax = plt.subplots(figsize=(20,10))
plot_inlier_matches(ax, img1, img2, data)
fig.savefig('sift_match.pdf', bbox_inches='tight')
```
%% Cell type:markdown id: tags:
#### Part (4) performn RANSAC to get the homography and inliers
%% Cell type:code id: tags:
``` python
# display the inlier matching, report the average residual
# <YOUR CODE>
# print("Average residual:", np.average(best_model_errors))
# print("Inliers:", max_inliers)
# fig.savefig('ransac_match.pdf', bbox_inches='tight')
```
%% Cell type:markdown id: tags:
#### Part (5) warp images to stitch them together
%% Cell type:code id: tags:
``` python
# display and report the stitching results
# <YOUR CODE>
# cv2.imwrite('stitched_images.jpg', im[:,:,::-1]*255.,
# [int(cv2.IMWRITE_JPEG_QUALITY), 90])
```
MP3/Q1/extra_credits/opt_01/park_center.jpg

810 KiB

MP3/Q1/extra_credits/opt_01/park_left.jpg

725 KiB

MP3/Q1/extra_credits/opt_01/park_right.jpg

708 KiB

MP3/Q1/extra_credits/opt_02/state_farm_center.jpg

629 KiB

MP3/Q1/extra_credits/opt_02/state_farm_left.jpg

740 KiB

MP3/Q1/extra_credits/opt_02/state_farm_right.jpg

740 KiB

MP3/Q1/extra_credits/opt_03/hessel_center.jpg

1.14 MiB

MP3/Q1/extra_credits/opt_03/hessel_left.jpg

1.33 MiB

MP3/Q1/extra_credits/opt_03/hessel_right.jpg

1.54 MiB

MP3/Q1/ransac_match.png

1.54 MiB

MP3/Q1/result.png

792 KiB

MP3/Q1/stitch/left.jpg

340 KiB

MP3/Q1/stitch/right.jpg

353 KiB

%% Cell type:markdown id: tags:
# Fundamental Matrix Estimation, Camera Calibration, Triangulation (Python)
## Usage
This code snippet provides an overall code structure and some interactive plot interfaces for the Fundamental Matrix Estimation, Camera Calibration, Triangulation section of MP 3. We outline the required functionalities step by step. Feel free to make modifications on the starter code if it's necessary.
## Package installation
- You will need [GUI backend](https://matplotlib.org/faq/usage_faq.html#what-is-a-backend) to enable interactive plots in `matplotlib`.
- `numpy`
- `PIL`
%% Cell type:markdown id: tags:
#### Common imports
%% Cell type:code id: tags:
``` python
from PIL import Image
import numpy as np
import matplotlib.pyplot as plt
```
%% Cell type:markdown id: tags:
#### Part (1)
%% Cell type:code id: tags:
``` python
def get_residual(F, p1, p2):
"""
Function to compute the residual average residual on frame 2
param: F (3x3): fundamental matrix: (pt in frame 2).T * F * (pt in frame 1) = 0
param: p1 (Nx2): 2d points on frame 1
param: p2 (Nx2): 2d points on frame 2
"""
P1 = np.c_[p1, np.ones((p1.shape[0],1))].transpose()
P2 = np.c_[p2, np.ones((p2.shape[0],1))].transpose()
L2 = np.matmul(F, P1).transpose()
L2_norm = np.sqrt(L2[:,0]**2 + L2[:,1]**2)
L2 = L2 / L2_norm[:,np.newaxis]
pt_line_dist = np.multiply(L2, P2.T).sum(axis = 1)
return np.mean(np.square(pt_line_dist))
def plot_fundamental(ax, F, p1, p2, I):
"""
Function to display epipolar lines and corresponding points
param: F (3x3): fundamental matrix: (pt in frame 2).T * F * (pt in frame 1) = 0
param: p1 (Nx2): 2d points on frame 1
param: p2 (Nx2): 2d points on frame 2
param: I: frame 2
"""
N = p1.shape[0]
P1 = np.c_[p1, np.ones((N,1))].transpose()
P2 = np.c_[p2, np.ones((N,1))].transpose()
L2 = np.matmul(F, P1).transpose() # transform points from
# the first image to get epipolar lines in the second image
L2_norm = np.sqrt(L2[:,0]**2 + L2[:,1]**2)
L2 = L2 / L2_norm[:,np.newaxis]
pt_line_dist = np.multiply(L2, P2.T).sum(axis=1)
closest_pt = p2 - (L2[:,0:2]*pt_line_dist[:,np.newaxis])
# Find endpoints of segment on epipolar line (for display purposes).
# offset from the closest point is 10 pixels
pt1 = closest_pt - np.c_[L2[:,1], -L2[:,0]]*10
pt2 = closest_pt + np.c_[L2[:,1], -L2[:,0]]*10
# Display points and segments of corresponding epipolar lines.
# You will see points in red corsses, epipolar lines in green
# and a short cyan line that denotes the shortest distance between
# the epipolar line and the corresponding point.
ax.set_aspect('equal')
ax.imshow(np.array(I))
ax.plot(p2[:,0],p2[:,1], '+r')
ax.plot([p2[:,0], closest_pt[:,0]],[p2[:,1], closest_pt[:,1]], 'r')
ax.plot([pt1[:,0], pt2[:,0]],[pt1[:,1], pt2[:,1]], 'g')
```
%% Cell type:code id: tags:
``` python
# write your code here for part estimating essential matrices
def fit_fundamental(...):
"""
Solves for the fundamental matrix using the matches with unnormalized method.
"""
# <YOUR CODE>
pass
def fit_fundamental_normalized(...):
"""
Solve for the fundamental matrix using the matches with normalized method.
"""
# <YOUR CODE>
pass
```
%% Cell type:code id: tags:
``` python
# Fundamental matrix estimation
name = 'library'
# You also need to report results for name = 'lab'
# name = 'lab'
I1 = Image.open('./{:s}1.jpg'.format(name))
I2 = Image.open('./{:s}2.jpg'.format(name))
matches = np.loadtxt('./{:s}_matches.txt'.format(name))
N = len(matches);
## Display two images side-by-side with matches
## this code is to help you visualize the matches, you don't need
## to use it to produce the results for the assignment
I3 = np.zeros((I1.size[1],I1.size[0]*2,3))
I3[:,:I1.size[0],:] = I1;
I3[:,I1.size[0]:,:] = I2;
fig, ax = plt.subplots()
ax.set_aspect('equal')
ax.plot(matches[:,0],matches[:,1], '+r')
ax.plot( matches[:,2]+I1.size[0],matches[:,3], '+r')
ax.plot([matches[:,0], matches[:,2]+I1.size[0]],[matches[:,1], matches[:,3]], 'r')
ax.imshow(np.array(I3).astype(np.uint8))
# non-normalized method
F = fit_fundamental(...) # <YOUR CODE>
pt1_2d = matches[:, :2]
pt2_2d = matches[:, 2:]
v2 = get_residual(F, pt1_2d, pt2_2d)
v1 = get_residual(F.T, pt2_2d, pt1_2d)
print('{:s}: residual in frame 2 (non-normalized method) = '.format(name), v2)
print('{:s}: residual in frame 1 (non-normalized method) = '.format(name), v1)
print('{:s}: residual combined (non-normalized method) = '.format(name), (v1+v2)/2)
# Plot epipolar lines in image I2
fig, ax = plt.subplots()
plot_fundamental(ax, F, pt1_2d, pt2_2d, I2)
# Plot epipolar lines in image I1
fig, ax = plt.subplots()
plot_fundamental(ax, F.T, pt2_2d, pt1_2d, I1)
# normalized method
F = fit_fundamental_normalized(...) # <YOUR CODE>
pt1_2d = matches[:, :2]
pt2_2d = matches[:, 2:]
v2 = get_residual(F, pt1_2d, pt2_2d)
v1 = get_residual(F.T, pt2_2d, pt1_2d)
print('{:s}: residual in frame 2 (normalized method) = '.format(name), v2)
print('{:s}: residual in frame 1 (normalized method) = '.format(name), v1)
print('{:s}: residual combined (normalized method) = '.format(name), (v1+v2)/2)
# Plot epipolar lines in image I2
fig, ax = plt.subplots()
plot_fundamental(ax, F, pt1_2d, pt2_2d, I2)
# Plot epipolar lines in image I1
fig, ax = plt.subplots()
plot_fundamental(ax, F.T, pt2_2d, pt1_2d, I1)
```
%% Cell type:markdown id: tags:
#### Part (2)
%% Cell type:code id: tags:
``` python
def evaluate_points(M, points_2d, points_3d):
"""
Visualize the actual 2D points and the projected 2D points calculated from
the projection matrix
You do not need to modify anything in this function, although you can if you
want to
:param M: projection matrix 3 x 4
:param points_2d: 2D points N x 2
:param points_3d: 3D points N x 3
:return:
"""
N = len(points_3d)
points_3d = np.hstack((points_3d, np.ones((N, 1))))
points_3d_proj = np.dot(M, points_3d.T).T
u = points_3d_proj[:, 0] / points_3d_proj[:, 2]
v = points_3d_proj[:, 1] / points_3d_proj[:, 2]
residual = np.sum(np.hypot(u-points_2d[:, 0], v-points_2d[:, 1]))
points_3d_proj = np.hstack((u[:, np.newaxis], v[:, np.newaxis]))
return points_3d_proj, residual
# Write your code here for camera calibration (lab)
def camera_calibration(...):
"""
write your code to compute camera matrix
"""
# <YOUR CODE>
pass
# Load 3D points, and their corresponding locations in
# the two images.
pts_3d = np.loadtxt('./lab_3d.txt')
matches = np.loadtxt('./lab_matches.txt')
# <YOUR CODE> print lab camera projection matrices:
lab1_proj = camera_calibration(...)
lab2_proj = camera_calibration(...)
print('lab 1 camera projection')
print(lab1_proj)
print('')
print('lab 2 camera projection')
print(lab2_proj)
# <YOUR CODE> evaluate the residuals for both estimated cameras
_, lab1_res = evaluate_points(...)
print('residuals between the observed 2D points and the projected 3D points:')
print('residual in lab1:', lab1_res)
_, lab2_res = evaluate_points(...)
print('residual in lab2:', lab2_res)
```
%% Cell type:code id: tags:
``` python
lib1_proj = np.loadtxt('./library1_camera.txt')
lib2_proj = np.loadtxt('./library2_camera.txt')
print('library1 camera projection')
print(lib1_proj)
print('library2 camera projection')
print(lib2_proj)
```
%% Cell type:markdown id: tags:
#### Part (3)
%% Cell type:code id: tags:
``` python
# Write your code here for computing camera centers
def calc_camera_center(...):
"""
write your code to get camera center in the world
from the projection matrix
"""
# <YOUR CODE>
pass
# <YOUR CODE> compute the camera centers using
# the projection matrices
lab1_c = calc_camera_center(...)
lab2_c = calc_camera_center(...)
print('lab1 camera center', lab1_c)
print('lab2 camera center', lab2_c)
# <YOUR CODE> compute the camera centers with the projection matrices
lib1_c = calc_camera_center(...)
lib2_c = calc_camera_center(...)
print('library1 camera center', lib1_c)
print('library2 camera center', lib2_c)
```
%% Cell type:markdown id: tags:
#### Part (4)
%% Cell type:code id: tags:
``` python
# Write your code here for triangulation
from mpl_toolkits.mplot3d import Axes3D
def triangulation(...):
"""
write your code to triangulate the points in 3D
"""
# <YOUR CODE>
pass
def evaluate_points_3d(...):
"""
write your code to evaluate the triangulated 3D points
"""
# <YOUR CODE>
pass
# triangulate the 3D point cloud for the lab data
matches_lab = np.loadtxt('./lab_matches.txt')
points_3d_gt = np.loadtxt('./lab_3d.txt')
points_3d_lab = triangulation(...) # <YOUR CODE>
res_3d_lab = evaluate_points_3d(...) # <YOUR CODE>
print('Mean 3D reconstuction error for the lab data: ', round(np.mean(res_3d_lab), 5))
_, res_2d_lab1 = evaluate_points(lab1_proj, lab_pt1, points_3d_lab)
_, res_2d_lab2 = evaluate_points(lab2_proj, lab_pt2, points_3d_lab)
print('2D reprojection error for the lab 1 data: ', np.mean(res_2d_lab1))
print('2D reprojection error for the lab 2 data: ', np.mean(res_2d_lab2))
# visualization of lab point cloud
camera_centers = np.vstack((lab1_c, lab2_c))
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.scatter(points_3d_lab[:, 0], points_3d_lab[:, 1], points_3d_lab[:, 2], c='b', label='Points')
ax.scatter(camera_centers[:, 0], camera_centers[:, 1], camera_centers[:, 2], c='g', s=50, marker='^', label='Camera Centers')
ax.legend(loc='best')
# triangulate the 3D point cloud for the library data
matches_lib = np.loadtxt('./library_matches.txt')
points_3d_lib = triangulation(...) # <YOUR CODE>
_, res_2d_lib1 = evaluate_points(lib1_proj, lib_pt1, points_3d_lib)
_, res_2d_lib2 = evaluate_points(lib2_proj, lib_pt2, points_3d_lib)
print('2D reprojection error for the library 1 data: ', np.mean(res_2d_lib1))
print('2D reprojection error for the library 2 data: ', np.mean(res_2d_lib2))
# visualization of library point cloud
camera_centers_library = np.vstack((lib1_c, lib2_c))
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.scatter(points_3d_lib[:, 0], points_3d_lib[:, 1], points_3d_lib[:, 2], c='b', label='Points')
ax.scatter(camera_centers_library[:, 0], camera_centers_library[:, 1],
camera_centers_library[:, 2], c='g', s=90,
marker='^', label='Camera Centers')
ax.view_init(azim=-45, elev=45)
ax.legend(loc='best')
```
MP3/Q2/lab1.jpg

91.7 KiB

MP3/Q2/lab2.jpg

105 KiB

312.747 309.140 30.086
305.796 311.649 30.356
307.694 312.358 30.418
310.149 307.186 29.298
311.937 310.105 29.216
311.202 307.572 30.682
307.106 306.876 28.660
309.317 312.490 30.230
307.435 310.151 29.318
308.253 306.300 28.881
306.650 309.301 28.905
308.069 306.831 29.189
309.671 308.834 29.029
308.255 309.955 29.267
307.546 308.613 28.963
311.036 309.206 28.913
307.518 308.175 29.069
309.950 311.262 29.990
312.160 310.772 29.080
311.988 312.709 30.514
8.800000000000000000e+02 2.140000000000000000e+02 7.310000000000000000e+02 2.380000000000000000e+02
4.300000000000000000e+01 2.030000000000000000e+02 2.200000000000000000e+01 2.480000000000000000e+02
2.700000000000000000e+02 1.970000000000000000e+02 2.040000000000000000e+02 2.300000000000000000e+02
8.860000000000000000e+02 3.470000000000000000e+02 9.030000000000000000e+02 3.420000000000000000e+02
7.450000000000000000e+02 3.020000000000000000e+02 6.350000000000000000e+02 3.160000000000000000e+02
9.430000000000000000e+02 1.280000000000000000e+02 8.670000000000000000e+02 1.770000000000000000e+02
4.760000000000000000e+02 5.900000000000000000e+02 9.580000000000000000e+02 5.720000000000000000e+02
4.190000000000000000e+02 2.140000000000000000e+02 3.280000000000000000e+02 2.440000000000000000e+02
3.170000000000000000e+02 3.350000000000000000e+02 4.260000000000000000e+02 3.860000000000000000e+02
7.830000000000000000e+02 5.210000000000000000e+02 1.064000000000000000e+03 4.700000000000000000e+02
2.350000000000000000e+02 4.270000000000000000e+02 4.800000000000000000e+02 4.950000000000000000e+02
6.650000000000000000e+02 4.290000000000000000e+02 9.640000000000000000e+02 4.190000000000000000e+02
6.550000000000000000e+02 3.620000000000000000e+02 6.950000000000000000e+02 3.740000000000000000e+02
4.270000000000000000e+02 3.330000000000000000e+02 5.050000000000000000e+02 3.720000000000000000e+02
4.120000000000000000e+02 4.150000000000000000e+02 6.450000000000000000e+02 4.520000000000000000e+02
7.460000000000000000e+02 3.510000000000000000e+02 6.920000000000000000e+02 3.590000000000000000e+02
4.340000000000000000e+02 4.150000000000000000e+02 7.120000000000000000e+02 4.440000000000000000e+02
5.250000000000000000e+02 2.340000000000000000e+02 4.650000000000000000e+02 2.630000000000000000e+02
7.160000000000000000e+02 3.080000000000000000e+02 5.910000000000000000e+02 3.240000000000000000e+02
6.020000000000000000e+02 1.870000000000000000e+02 4.470000000000000000e+02 2.130000000000000000e+02
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment