Projection matrix from quadliteral

proj.py

import numpy as np

# 'f'rom is a set of (x, y)
# 't'to  is a set of (u, v)
# H * (x, y, 1) = k*(u, v, 1) for every point
# combine all equtaions and solve: Ah = b
# h is (h0, .., h7), add h8 = 1, zeros for z
# and make a 4x4 projection matrix
if __name__ == '__main__':
    f = [(0.0, 0.0), (0.0, 100.0), (100.0, 100.0), (100.0, 0.0)]
    t = [(0.0, 0.0), (0.0, 100.0), (100.0, 55.0), (100.0, 45.0)]
    A = np.array(
        [[f[0][0], f[0][1], 1, 0, 0, 0, -t[0][0] * f[0][0], -t[0][0] * f[0][1]]
        ,[0, 0, 0, f[0][0], f[0][1], 1, -t[0][1] * f[0][0], -t[0][1] * f[0][1]]
        ,[f[1][0], f[1][1], 1, 0, 0, 0, -t[1][0] * f[1][0], -t[1][0] * f[1][1]]
        ,[0, 0, 0, f[1][0], f[1][1], 1, -t[1][1] * f[1][0], -t[1][1] * f[1][1]]
        ,[f[2][0], f[2][1], 1, 0, 0, 0, -t[2][0] * f[2][0], -t[2][0] * f[2][1]]
        ,[0, 0, 0, f[2][0], f[2][1], 1, -t[2][1] * f[2][0], -t[2][1] * f[2][1]]
        ,[f[3][0], f[3][1], 1, 0, 0, 0, -t[3][0] * f[3][0], -t[3][0] * f[3][1]]
        ,[0, 0, 0, f[3][0], f[3][1], 1, -t[3][1] * f[3][0], -t[3][1] * f[3][1]]
        ])
    b = np.array(
        [t[0][0], t[0][1], t[1][0]
        ,t[1][1], t[2][0], t[2][1]
        ,t[3][0], t[3][1]]);
    print('from', f)
    print('to  ', t)
    print('A=')
    print(A)
    print('B=')
    print(b)
    h = np.linalg.solve(A, b)
    assert(np.allclose(np.dot(A, h), b))
    print('h=')
    print(h)
    print('proj=')
    p = np.array(
        [[h[0], h[1], 0, h[2]]
        ,[h[3], h[4], 0, h[5]]
        ,[0,       0, 1, 0]
        ,[h[6], h[7], 0, 1]
        ])
    print(p)