interpolation between arbitrary quadrilaterals

Bilinear interpolation works just as well for non-rectangular quadrilaterals. Consider that you have the coordinates (s,t) each with values in the range 0..1, in the manner that I'm sure you've already read about.



The coordinate s will represent the texel's position between A and B, as well as between D and C. Considered the two points along these lines; we will call theme P and Q. The coordinate t is the relative position between P and Q.

As you can see, this works well for both cases. The only real problem is that the texture will be deformed to fit the quad.

thanks for the reply, but what I can't seem to get is how to find the (s,t) coordinates given the point's absolute position. I know the 3D positions of the 4 vertices of the non-rectangular original quad and I also know the 3D position of the point inside this Quad, I can easily project this into a 2D representation, but I can't see how to get the (s,t) coordinates you show in the second drawing from this.

I tried expressing the positions of the points to find s and t like this:

P = A + (B-A)s
Q = D + (C-D)s

thus the resulting point, which I know (M) is:

M = P + (Q-P)t

this gives me a non-linear system of equations on s and t :(




Thanks in advance


Diego

Well, the math gets a little rough, but the answer is elegant.

P = A+(B-A)s = A(1-s) + Bs
Q = D+(C-D)s = D(1-s) + Cs

M = P+(Q-P)t = P(1-t) + Qt

= A(1-s)(1-t) + B(1-t)s + D(1-s)t + Cst
= A(1+st-s-t) + B(s-st) + Cst + D(t-st)
= A + Ast - As - At + Bs - Bst + Cst + Dt - Dst
= (1-s)A + st(A-D+C-B) - t(A-D)

This might seem like a problem, having two unknowns and one equation, but that's not really the case. What we really have is:

M_x = (1-s)A_x + st(A-D+C-B)_x - t(A-D)_x
M_y = (1-s)A_y + st(A-D+C-B)_y - t(A-D)_y

We can solve for s in each case:

M_x = (s)(t(A-D+C-B)_x - A_x) + A_x - t(A-D)_x
M_y = (s)(t(A-D+C-B)_y - A_y) + A_y - t(A-D)_y

s = (M_x + t(A-D)_x - A_x)/(t(A-D+C-B)_x - A_x)
s = (M_y + t(A-D)_y - A_y)/(t(A-D+C-B)_y - A_y)

Setting these equal and solving for t:
(M_x + t(A-D)_x - A_x)/(t(A-D+C-B)_x - A_x) = (M_y + t(A-D)_y - A_y)/(t(A-D+C-B)_y - A_y)

Also, for my sanity in typing this:
T = M-A
U = A-D
V = A-D+C-B

(T_x+tU_x)(tV_y-A_y) = (T_y+tU_y)(tV_x-A_x)

-A_yT_x + t(T_xV_y - U_xA_y) + t²U_xV_y = -A_yT_x + t(T_xV_y - U_xA_y) + t²U_xV_y

Lastly:
t²(U_xV_y - U_yV_x) + t(T_xV_y - U_xA_y - T_yV_x + U_yA_x) + (A_xT_y - A_yT_x) = 0

Which can be stuffed into the quadratic equation. From here it should be trivial to find a solution, if one exists, where t is between 0 and 1. Then by plugging it back into the above equation, you can find an s, if one exists, that meets the same conditions.

(The solution can probably be simplified, and in fact I see many cross product candidates. Maybe I'll come back later.)

Okay, I couldn't help it:

t = [-b ± √(b²-4ac)]/2a

where:
a = (A-D)⊗(C-B)
b = M⊗(A-D+C-B) - A⊗(C-B)
c = A⊗M

I am too tired to make any assumptions as to whether or not this will work with any vector not of the form

Thank you so much for all your help!, the solution is really elegant :)
As to whether this would work outside the z=0 plane, as originally I have 4 (x,y,z) points, I would actually have 3 formulas instead of 2, but one of them would be unnecessary because all the points lie on a plane. However, due to floating point error, this might not be the case when you plug everything in, so I guess the best way to deal with this would be to project everything into a 2D coordinate system.

It's really great how in the resulting formula for t, most of the values can be precalculated during the creation of the Quad. Does anyone know if this is how OpenGL does texture mapping?

The "standard" answer (which OpenGL uses, AFAIK) is to simply not use quads; use triangles. Then you have a nice, simple, and unique mapping: barycentric interpolation. Using triangles is also no limitation, since obviously you can split quads into triangles.

Here's how it works:

Say you have a triangle with vertices (p1, p2, p3) and a point p. Also, let's say the triangle has edge vectors e2=p2-p1 and e3=p3-p1, and normal n=(p2-p1)x(p3-p1), where 'x' is the cross product. Then, form the 3x3 matrix


M = [e2, e3, n]

whose columns are the vectors introduced above, and compute

q = M^(-1) (p - p1) .

The vector q=(q1,q2,q3) that you get gives the coordinates of p in the "local frame" of the triangle. If p lies on the plane of the triangle, then q3=0. The barycentric coordinates of the projection of p onto the plane of the triangle are given by b=(b1,b2,b3)=(1-q1-q2, q1, q2). (Clearly if p is on the plane then these are simply the barycentric coordinates of p).

Now, suppose that vertex p1 of the triangle has been assigned the texture coordinates (u1, v1); p2 has been assigned (u2, v2), and so on. The u,v coordinates of p in the texture are simply

(u,v) = b1 (u1, v1) + b2 (u2, v2) + b3 (u3, v3) .

Since you can compute M^(-1) ahead of time, you see that this method requires just some multiplications, additions, and subtraction -- no square roots, etc.

Cheers.

[Edited by - Emergent on June 24, 2009 10:45:23 AM]

I see a problem with that approach however, imagine (I can't figure out how to include an image :S) you have a quad shaped like this:

P0=(0,0)
P1=(1,0)
P2=(1,1)
P3(0,10)

Also, assume you're assigning the texture coordinates (0,0), (1,0), (1,1) and (0,1) to these 4 points resp.

If you divide it into 2 triangles (say (P0 P1 P2) and (P0 P2 P3)) and chose to interpolate depending on which triangle the point is in, you would have a very big difference in the stretching of the texture the minute you step beyond the P0 P2 line from one triangle to another. Maybe you can treat this as collateral damage given that the 2-triangle approach is more efficient, but I think this is an unwanted effect

You do not want to use a bilinear interpolation method. Use a perspective transformation instead so that you preserve lines and quadratics. See my Miscellaneous page, files Wm4QuadToQuadTransforms.{h,cpp}. The comments in the header show how to do the bilinear interpolation and what problems arise. The second PDF link on that page shows how to perform the perspective transformation.

Dave Eberly:

That looks promising, but I'll have to devote time to get to understand it, a time I don't have given that the deadline for my work is next week. I plan to continue working on the raytracer though, so I may look into Perspective Transformations in the future.


erissian:

I was checking the formulas you derived and I found you ommited the term Bs when stepping from this line:

= A + Ast - As - At + Bs - Bst + Cst + Dt - Dst

to this one

= (1-s)A + st(A-D+C-B) - t(A-D)



I added it and arrived to the following:

s = [ M - A + t(A-D) ] / [ B - A + t(A-D+C-B) ]


This gives the following when substituting for T, U and V:


(T_x + tU_x)*(tV_y + B_y - A_y) =
(T_y + tU_y)*(tV_x + B_x - A_x)

and:

t²U⊗V + t(T⊗V + U⊗(B-A)) + T⊗(B-A) = 0

so the coefficients of the quadratic formula are:

a = (A-D)⊗(C-B)
b = (M-A)⊗(A-D+C-B) + (A-D)⊗(B-A)
c = (M-A)⊗(B-A)

However, there must be something wrong with my calculations, because the texture is not displaying properly. Substituting A by M yields a correct t (0), and the same happens for B (0), but substituting for C doesn't...

Quote:

Original post by diegovar However, there must be something wrong with my calculations, because the texture is not displaying properly. Substituting A by M yields a correct t (0), and the same happens for B (0), but substituting for C doesn't...

I ran through this again (after sleeping) and have a successfully tested solution. If you still need it, I'll post it. (Otherwise, it's too much typing!)

In the end, the coefficients were right. For the sake of completeness, here's a, b and c of the quadratic.


a = (A-D)⊗(C-B)
b = (M-A)⊗(A-D+C-B) + (A-D)⊗(B-A)
c = (M-A)⊗(B-A)

Here's a pic showing the texture coords at work.


I must've messed up the code previously, because it worked, we presented our work on Thursday and the professor was very pleased :D.

Here's a couple of pictures of our resulting work.





Thanks a lot for all your help!


Diego

Quote:

Original post by diegovar we presented our work on Thursday and the professor was very pleased :D.

Good!

Quote:

Here's a couple of pictures of our resulting work.

Those all look fantastic!

I know this is a VERY old thread, but I encountered this problem as well and found the only solution here.
Unfortunately, I think erissian's formula is flawed. There is a small error which kept me busy for a while, but I solved it and thought I could share the solution:

erissian wrote:
M[sub]x[/sub] = (1-s)A[sub]x[/sub] + st(A-D+C-B)[sub]x[/sub] - t(A-D)[sub]x[/sub]
M[sub]y[/sub] = (1-s)A[sub]y[/sub] + st(A-D+C-B)[sub]y[/sub] - t(A-D)[sub]y[/sub]

But it should be:
M[sub]x[/sub] = [color=#ff0000]sB[sub]x[/sub]+(1-s)A[sub]x[/sub] + st(A-D+C-B)[sub]x[/sub] - t(A-D)[sub]x[/sub]
M[sub]y[/sub] = [color=#ff0000]sB[sub]y[/sub]+(1-s)A[sub]y[/sub] + st(A-D+C-B)[sub]y[/sub] - t(A-D)[sub]y[/sub]

This simplifies to these coefficients of the quadratic equation at[sup]2[/sup] + bt + ct = 0
a = C?D+D?B+A?C+B?A
b = M?C+D?M+B?M+M?A+C?A+B?D+2(A?B)
c = M?B+A?M+B?A

I don't know if these are the same values diegovar calculated in the end, but they work