CODEWITHSUNDEEP
numpytransformleast-squares
Felipe Moser
Felipe Moser

Reputation: 365

Numpy least-squares solution not accurate results

I'm calculating the affine transformation I need from a a few points in a 3D space and using numpy.linalg.lstsq to do so. However, the results I'm getting, while not terribly far off, are not accurate enough, even in trivially simple examples:

m = 100
xy = np.array([[0, 0, 0],
           [m, 0, 0],
           [m, m, 0],
           [0, m, 0],
           [0, 0, m],
           [m, 0, m],
           [m, m, m],
           [0, m, m]])
uv = np.array([[0.5, 0, 0],
           [m + 0.5, 0, 0],
           [m+ 0.5, m, 0],
           [0.5, m, 0],
           [0.5, 0, m],
           [m+ 0.5, 0, m],
           [m+ 0.5, m, m],
           [0.5, m, m]])

pts_a = np.hstack([uv, np.ones((uv.shape[0], 1))])
pts_b = np.hstack([xy, np.ones((xy.shape[0], 1))])
solution_1 = np.linalg.lstsq(pts_a, pts_b, rcond=None)[0]

The result I'm expecting from the above code is:

[[1, 0, 0, -0.5],
 [0, 1, 0, 0],
 [0, 0, 1, 0],
 [0, 0, 0, 1]])

The result I'm getting: 

 [[ 1.00000000e+00  3.49047642e-16  3.60109527e-16 -5.00000000e-01]
 [ 1.77081442e-16  1.00000000e+00 -3.93150475e-16  1.80460546e-15]
 [ 2.21351803e-16 -3.11848610e-16  1.00000000e+00 -6.28251374e-15]
 [ 2.76689754e-18  1.06035619e-17 -1.19061095e-17  1.00000000e+00]]

Those small differences make a considerable difference in my results. Any ideas how to solve it? NOTE: I can ONLY use numpy and math for my project, so using a different library is sadly not possible! Thanks!

Upvotes: 2

Views: 1014

Answers (1)

ead
ead

Reputation: 34367

Actually, the difference isn't small but quite big - you have the wrong sign for solution[0,3].

The problem is thath you didn't calculate the desired transformation T but the inverse of this transformation, i.e. T^-1.

Let's do some math:

T*X=U, with X - original vectors
            U - transformed vectors

transposing it =>
X^t * T^t = U^t
 |     |     |
\|/   \|/   \|/
 A  *  x  =  b

In your program A=pts_b and b=pts_a, that means the transformation T is (you have to swap pts_b and pts_b and to transpose the result to get the right matrix):

T = np.linalg.lstsq(pts_b, pts_a)[0].T

and voila:

>>> T
array([[  1.00000000e+00,  -8.15320034e-17,  -6.59194921e-17,  5.00000000e-01],
       [ -4.97379910e-16,   1.00000000e+00,   7.77156117e-16,  -1.02678283e-14],
       [ -2.13162819e-16,   4.44089210e-16,   1.00000000e+00,   1.91513472e-15],
       [ -4.44089205e-18,  -8.84708973e-17,   9.88792381e-17,   1.00000000e+00]])

PS: You have solved the equation:

X^t  =  U^t * (T^t)^(-1)
 |       |     |
\|/     \|/   \|/
 b   =   A  *  x

Upvotes: 1

Related Questions