I will show this using python and numpy.
import numpy as np
First, you have to understand the projection model
def apply_homography(H, p1):
p = H @ p1.T
return (p[:2] / p[2]).T
With some algebraic manipulation you can determine the points at the plane z=1 that produced the given points.
def revert_homography(H, p2):
Hb = np.linalg.inv(H)
# 1 figure out which z coordinate should be added to p2
# order to get z=1 for p1
z = 1/(Hb[2,2] + (Hb[2,0] * p2[:,0] + Hb[2,1]*p2[:,1]))
p2 = np.hstack([p2[:,:2] * z[:,None], z[:, None]])
return p2 @ Hb.T
The projection is not invertible, but under the complanarity assumption it may be inverted successfully.
Now, let's see how to determine the H matrix from the given points (assuming they are coplanar).
If you have the four corners in order in order you can simply specify the (x,y) coordinates of the cornder, then you can use the projection equations to determine the homography matrix like here, or here.
This requires at least 5 points to be determined as there is 9 coefficients, but we can fix one element of the matrix and make it an inhomogeneous equation.
def find_homography(p1, p2):
A = np.zeros((8, 2*len(p1)))
# x2'*(H[2,0]*x1+H[2,1]*x2)
A[6,0::2] = p1[:,0] * p2[:,0]
A[7,0::2] = p1[:,1] * p2[:,0]
# - (H[0,0]*x1+H[0,1]*y1+H[0,2])
A[0,0::2] = -p1[:,0]
A[1,0::2] = -p1[:,1]
A[2,0::2] = -1
# y2'*(H[2,0]*x1+H[2,1]*x2)
A[6,1::2] = p1[:,0] * p2[:,1]
A[7,1::2] = p1[:,1] * p2[:,1]
# - (H[1,0]*x1+H[1,1]*y1+H[1,2])
A[3,1::2] = -p1[:,0]
A[4,1::2] = -p1[:,1]
A[5,1::2] = -1
# assuming H[2,2] = 1 we can pass its coefficient
# to the independent term making an inhomogeneous
# equation
b = np.zeros(2*len(p2))
b[0::2] = -p2[:,0]
b[1::2] = -p2[:,1]
h = np.ones(9)
h[:8] = np.linalg.lstsq(A.T, b, rcond=None)[0]
return h.reshape(3,3)
Here a complete usage example. I pick a random H and transform four random points, this is what you have, I show how to find the transformation matrix H_. Next I create a test set of points, and I show how to find the world coordinates from the image coordinates.
# Pick a random Homography
H = np.random.rand(3,3)
H[2,2] = 1
# Pick a set of random points
p1 = np.random.randn(4, 3);
p1[:,2] = 1;
# The coordinates of the points in the image
p2 = apply_homography(H, p1)
# testing
# Create a set of random points
p_test = np.random.randn(20, 3)
p_test[:,2] = 1;
p_test2 = apply_homography(H, p_test)
# Now using only the corners find the homography
# Find a homography transform
H_ = find_homography(p1, p2)
assert np.allclose(H, H_)
# Predict the plane points for the test points
p_test_predicted = revert_homography(H_, p_test2)
assert np.allclose(p_test_predicted, p_test)