I am trying to map coordinates from a quadrilateral to a square. The coordinates are
square: $(500,900)(599,900)(599,999)(500,999)$
quad: $(454,945)(558,951)(598,999)(499,999)$
where the $i^{th}$ entry are corresponding corners. It looks like this:
I am following the method 2 from this article http://www.fmwconcepts.com/imagemagick/bilinearwarp/FourCornerImageWarp2.pdf
Based on its instructions, I can make these two matrices:
$\begin{bmatrix} 1 & 500 & 900 & 450000 \\ 1 & 599 & 900 & 539100 \\ 1 & 599 & 999 & 598401 \\ 1 & 500 & 999 & 499500 \end{bmatrix} \begin{bmatrix} a_0 \\ a_1 \\ a_2 \\ a_3 \end{bmatrix} = \begin{bmatrix} 454 \\ 558 \\ 598 \\ 499 \end{bmatrix}$
$\begin{bmatrix} 1 & 500 & 900 & 450000 \\ 1 & 599 & 900 & 539100 \\ 1 & 599 & 999 & 598401 \\ 1 & 500 & 999 & 499500 \end{bmatrix} \begin{bmatrix} b_0 \\ b_1 \\ b_2 \\ b_3 \end{bmatrix} = \begin{bmatrix} 945 \\ 951 \\ 999 \\ 999 \end{bmatrix}$
If I solve for them I get
$a0=-709.911845730028$
$a1=1.50964187327824$
$a2=0.709621467197225$
$a3=-0.000510152025303541$
$b0=148.305785123967$
$b1=0.611570247933884$
$b2=0.85154576063667$
$b3=-0.000612182430364249$
But then when I run this python script (I want to simulate manually, converting the point (454, 945) on the quad to the point (500,900) on the square), I get a different answer. I get (442.90822654, 1024.0)...
X = 454
Y = 945
a0=-709.911845730028
a1=1.50964187327824
a2=0.709621467197225
a3=-0.000510152025303541
b0=148.305785123967
b1=0.611570247933884
b2=0.85154576063667
b3=-0.000612182430364249
A = b2*a3 - b3*a2
C_one = (b0*a1 - b1*a0)
C = C_one + (b1*X - a1*Y)
B_one = (b0*a3 - b3*a0) + (b2*a1 - b1*a2)
B = B_one + (b3*X - a3*Y)
V = (-B + (B*B - 4*A*C)**0.5 ) / (2*A)
U = (X - a0 - a2*V) / (a1 + a3*V)
print U,V
I tried a test where I warp an image, filled with different quads. It mostly worked except for two areas, ie.. the left eye and right leg. The above example is for the case of the right leg.
original image:
warped image:
warp template:
The quad that I am using for the above example is in the last row in the middle, to the right of the red/maroon quad.
Does anyone know what is going on? Is it possibly a numerical/rounding error? But it did work with all the quads except two of them.
Thanks
Confusion about the direction of the mapping seems to be deeply embedded in the question. The article describes a way of mapping a square to a (curved) quadrilateral using what is often called the "bilinear quadrilateral element" in the literature.
This is not a linear transformation, and because it contains a quadratic term, it is not trivial to invert the mapping.
Added:
The OP has apparently diagnosed the problem. The Python code mapping in the opposite direction, from the quadrilateral to the square, needs to solve a quadratic. When two roots are nearly equal, care ought to be taken in evaluating the roots via the usual quadratic formula to avoid "catastrophic cancellation" by having terms of like sign in the numerator. The OP reports (see Comment) that after taking this measure the mapping works as expected.