From that very page you posted, there's a link to the source
code. I'll explain the bilinear transformation in
http://www.antigrain.com/__code/include/agg_trans_bilinear.h.html
The idea here is to find a transformation of the form:
output_x = a * input_x + b * input_x * input_y + c * input_y + d
output_y = e * input_x + f * input_x * input_y + g * input_y + h
The term "bilinear" comes from each of those equations being linear in
either of the input coordinates by themselves. We want to solve for
the right values of a, b, c, and d. Say you have the reference
rectangle r1, r2, r3, r4 which you want to map to (0,0), (1,0), (0,1),
(1,1) (or some image coordinate system).
For a,b,c,d:
0 = a * r1_x + b * r1_x * r1_y + c * r1_y + d
1 = a * r2_x + b * r2_x * r2_y + c * r2_y + d
0 = a * r3_x + b * r3_x * r3_y + c * r3_y + d
1 = a * r4_x + b * r4_x * r4_y + c * r4_y + d
For e,f,g,h:
0 = e * r1_x + f * r1_x * r1_y + g * r1_y + h
0 = e * r2_x + f * r2_x * r2_y + g * r2_y + h
1 = e * r3_x + f * r3_x * r3_y + g * r3_y + h
1 = e * r4_x + f * r4_x * r4_y + g * r4_y + h
You can solve this however you like best. (If you're familiar with
matrix notation, these are two matrix equations for which the matrix
is the same, and then you simply need to find the LU decomposition
once, and solve the two unknown vectors). The coefficients are then
applied to map the interior of the rectangle to the position in the
rectangle.
If by any chance you're looking for the inverse transform, that is,
if you want to know where a given pixel will land, you simply switch
inputs and outputs:
For a,b,c,d:
r1_x = a * 0 + b * 0 * 0 + c * 0 + d
r2_x = a * 1 + b * 1 * 0 + c * 0 + d
r3_x = a * 0 + b * 0 * 1 + c * 1 + d
r4_x = a * 1 + b * 1 * 1 + c * 1 + d
For e,f,g,h:
r1_y = e * 0 + f * 0 * 0 + g * 0 + h
r2_y = e * 1 + f * 1 * 0 + g * 0 + h
r3_y = e * 0 + f * 0 * 1 + g * 1 + h
r4_y = e * 0 + f * 0 * 1 + g * 1 + h