Calculating a 2D Transformation Matrix from an Initial and Resultant 2D Matrix

Question

I by no means profess to be a genius when it comes to programming and my current problem has me stumped.

I have found this question Trying to derive a 2D transformation matrix using only the images that seems to at least partially answer my question but the image that should show the solution is no longer available :S

I'm working in C# and not using WPF as neither my input or output needs to be displayed graphically.

In my program I have 2 quadrilaterals, lets call them an input and an output quadrilateral.

The input quad has the co-ords of (2,1),(2,3),(4,4),(3,1) from bottom left clockwise.

The output quad can have any co-ords and will be again listed in order from bottom left clockwise.

Given these 8 co-ordinate pairs, is it possible to calculate a transformation matrix that I could apply to any single co-ordinate pair?

I'm not too hot on Matrices but I am willing to learn if pointed in the right direction.

Many Thanks

Josh

Answer 1

The title of your question is misleading, because it implies that you have the initial and final matrices. What you really have, though, is two sets of points: the beginning points and the ending points.

First, as you mentioned, it's possible to compute a transformation matrix, but not necessarily the transformation matrix. For any particular transformation, there are multiple (infinite?) ways to accomplish it.

Note that the following is just some thoughts. I haven't actually done this. I assume that you have two sets of points, A and B, and you know for certain that B is the result of applying some transformation to A.

The problem is trivial if the only transformation allowed is translation. In that case, you can just take the distance between the bottom left points. That is, if the original points are A[0] through A[3], and the new points are B[0] through B[3], then the transformation is just the X,Y translation: ((B[0].X - A[0].X), (B[0].Y - A[0].Y)).

If scaling is allowed as well, then you can figure the translation and then the scaling. Although to simplify, first you'll want to translate to the origin. In fact, most of this gets simpler if you translate to the origin first.

If rotation is allowed, things get a bit more difficult. Given quadrilaterals, first you have to rotate the thing so that the points are in the same orientation as the original. This will require that you compute the distances between the points and use the ratios to figure out which point is the lower left. Then rotate into the proper place.

The case of rotation, scaling, and translation should be straightforward to solve.

Shearing is a bit more complicated. You should be able to detect shearing, but I don't have a good idea of how to detect the amount of shearing. I think, though, that if you solve the case of rotate/scale/translate, then the shearing solution will become at least a little bit more obvious.

Answer 2

A quick google or a hop, skip and a google found this. I think it will definitely solve your problem.

As I mentioned in comment, you are asking for an isomorphic function to project a single point within a quadrilateral to a point in a second quadrilateral. Rather than work it out manually I've found the below algorithm.

Posting the algorithm here for posterity:

Let p00, p10, p11,and p01 be the vertices of the first quadrilateral listed in counterclockwise order.

Let q00, q10, q11,and q01 be the vertices of the second quadrilateral listed in counterclockwise order.

Define P10 = p10-p00, P11 = p11-p00, P01 = p01-p00, Q10 = q10-q00, Q11 = q11-q00, and Q01 = q01-q00.

Compute a and b so that Q11 = aQ10 + bQ01.

This is a set of two linear equations in two unknowns.

Similarly, compute c and d so that P11 = cP10 + dP01.

It turns out that a >= 0, b >= 0, a + b> 1, c >= 0, d >= 0,and c + d> 1.

Any point p in the quadrilateral can be written as p =(1-x-y)p00 + xp10 + yp01 where x = 0, y = 0,(1 - d)x + c(y - 1) = 0, and d(x - 1)+(1 - c)y = 0.

Any point q in the quadrilateral can be written as q =(1-u-v)q00 + uq10 + vq01 where u = 0, v = 0,(1 - b)u + a(v -1) = 0, and b(u - 1) + (1 - a)v = 0.

The perspective mapping relating (u,v) to (x,y) is u = m0x n0 + n1x + n2y and v = m1y n0 + n1x + n2y where m0 = ad(1 - c - d), m1 = bc(1 - c - d), n0 = cd(1 - a - b), n1 = d(a - c + bc - ad),and n2 = c(b - d - bc + ad).

Effectively the p-quadrilateral is mapped to a “standard” one, <(0, 0),(1, 0),(0, 1),(c,d)>

and the q-quadrilateral is mapped to <(0, 0),(1, 0), (0, 1), (a,b)>.

The (x,y) to (u,v) mapping relates these two.

You can verify that

• (x,y)=(0, 0) is mapped to (u,v)=(0, 0)

• (x,y)=(1, 0) is mapped to (u,v)=(1, 0)

• (x,y)=(0, 1) is mapped to (u,v)=(0, 1)

• (x,y)=(c,d) is mapped to (u,v)=(a,b)

I'm going to give another answer that describes how I'd solve the example in the comments - this answer is getting too long.

Answer 3

To solve an example case - transposing a point (5,5) from a square (0,0)-(10,10) onto a square (0,0)-(20,20) - you are on the right track but either stopped too soon or got lost. The algorithm isn't the simplest but it is quite usable once you understand where you're going.

When dealing with 'squares' which are perpendicular to x & y axes, then yes a = b = c = d = 1.

Let's call the points of your quad p1, p2, p3 and p4. Now think of a as describing a relationship between p2 and p3 (relative to p1), and b as describing a similar relationship between p4 and p3 (relative to p1). Very simplistically, you want to work out a and b where

• the x coordinate of p1 + a * length of side p1-p2 = x coordinate of p3, and
• the y coordinate of p1 + b * length of side p1-p4 = y coordinate of p3.

i.e. starting from the bottom left, how many many bottom edges do I need to add to get to the x coordinate of the top right corner, and similarly for the y coordinate. Although that's a mouthful, when you are dealing with a square the correct result here is definitely 1 and 1. Once you move away from perpendicular squares the it's not quite as simple but this is easy to visualise.

When you solve this

• p =(1-x-y)p00 + xp10 + yp01, where
• x >= 0
• y >= 0
• (1 - d)x + c(y - 1) = 0
• d(x - 1)+(1 - c)y = 0

using a point (5,5) you will get a value of x = 1/2, y = 1/2. Note that these are not the x and y coordinates of your point. Instead they represent where your point would sit if the quadrilateral was projected to a (1,1) size square. In this case they means that your point is positioned 1/2 way across (horizontally) and 1/2 way up (vertically) your polygon - i.e. it's in the middle.

Plug in your a,b,c,d,x and y values into the large equation (following the link rather than trying to read it here!), you also get (u,v) = (1/2, 1/2), where u & v represent the same concept above but in the second square. This is correct because the resulting point is expected to be in the middle of the second polygon.

Finally, when you plug in your u, v into the equation

• q =(1-u-v)q00 + uq10 + vq01

you'll get point q = (10,10), which is what you expected.

I don't know & haven't thought about solving simultaneous equations in code, I'm sure there's a way but it might not be simple but unfortunately I'll have to leave it up to you. I've just done this on a piece of paper and it all worked out; however my maths is a bit rusty to do anything more than the more trivial example.