I am given the locations of three points:
p1 = [1.0, 1.0, 1.0]
p2 = [1.0, 2.0, 1.0]
p3 = [1.0, 1.0, 2.0]
and their transformed counterparts:
p1_prime = [2.414213562373094, 5.732050807568877, 0.7320508075688767]
p2_prime = [2.7677669529663684, 6.665063509461097, 0.6650635094610956]
p3_prime = [2.7677669529663675, 5.665063509461096, 1.6650635094610962]
The affine transformation matrix is of the form
trans_mat = np.array([[…, …, …, …],
[…, …, …, …],
[…, …, …, …],
[…, …, …, …]])
such that with
import numpy as np
def transform_pt(point, trans_mat):
a = np.array([point[0], point[1], point[2], 1])
ap = np.dot(a, trans_mat)[:3]
return [ap[0], ap[1], ap[2]]
you would get:
transform_pt(p1, trans_mat) == p1_prime
transform_pt(p2, trans_mat) == p2_prime
transform_pt(p3, trans_mat) == p3_prime
Assuming the transformation is homogeneous (consists of only rotations and translations), how can I determine this transformation matrix?
From a CAD program, I know the matrix is:
trans_mat = np.array([[0.866025403784, -0.353553390593, -0.353553390593, 0],
[0.353553390593, 0.933012701892, -0.066987298108, 0],
[0.353553390593, -0.066987298108, 0.933012701892, 0],
[0.841081377402, 5.219578794378, 0.219578794378, 1]])
I'd like to know how this can be found.
Six points alone is not enough to uniquely determine the affine transformation. However, based on what you had asked in a question earlier (shortly before it was deleted) as well as your comment, it would seem that you are not merely looking for an affine transformation, but a homogeneous affine transformation.
This answer by robjohn provides the solution to the problem. Although it solves a more general problem with many points, the solution for 6 points can be found at the very bottom of the answer. I shall transcribe it here in a more programmer-friendly format:
import numpy as np
def recover_homogenous_affine_transformation(p, p_prime):
'''
Find the unique homogeneous affine transformation that
maps a set of 3 points to another set of 3 points in 3D
space:
p_prime == np.dot(p, R) + t
where `R` is an unknown rotation matrix, `t` is an unknown
translation vector, and `p` and `p_prime` are the original
and transformed set of points stored as row vectors:
p = np.array((p1, p2, p3))
p_prime = np.array((p1_prime, p2_prime, p3_prime))
The result of this function is an augmented 4-by-4
matrix `A` that represents this affine transformation:
np.column_stack((p_prime, (1, 1, 1))) == \
np.dot(np.column_stack((p, (1, 1, 1))), A)
Source: https://math.stackexchange.com/a/222170 (robjohn)
'''
# construct intermediate matrix
Q = p[1:] - p[0]
Q_prime = p_prime[1:] - p_prime[0]
# calculate rotation matrix
R = np.dot(np.linalg.inv(np.row_stack((Q, np.cross(*Q)))),
np.row_stack((Q_prime, np.cross(*Q_prime))))
# calculate translation vector
t = p_prime[0] - np.dot(p[0], R)
# calculate affine transformation matrix
return np.column_stack((np.row_stack((R, t)),
(0, 0, 0, 1)))
For your sample inputs, this recovers the exact same matrix as what you had obtained from the CAD program:
>>> recover_homogenous_affine_transformation(
np.array(((1.0,1.0,1.0),
(1.0,2.0,1.0),
(1.0,1.0,2.0))),
np.array(((2.4142135623730940, 5.732050807568877, 0.7320508075688767),
(2.7677669529663684, 6.665063509461097, 0.6650635094610956),
(2.7677669529663675, 5.665063509461096, 1.6650635094610962))))
array([[ 0.8660254 , -0.35355339, -0.35355339, 0. ],
[ 0.35355339, 0.9330127 , -0.0669873 , 0. ],
[ 0.35355339, -0.0669873 , 0.9330127 , 0. ],
[ 0.84108138, 5.21957879, 0.21957879, 1. ]])