Functions that map a quadrilateral to the unit square?

Question

Given some quadrilateral $Q \subset \mathbb R^2$ defined by the vertices $P_i = (x_i,y_i), i=1,2,3,4$ (you can assume they are in positive orientation), is there a function $f: \mathbb R^2 \to \mathbb R^2$ that is particularly easy to compute which satisfies

$$f(Q) \subseteq U \quad \text{ and } \quad f(\mathbb R^2 \setminus Q) \subseteq \mathbb R^2 \setminus U?$$

Here $U = [0,1]\times[0,1]$ denotes the unit square (or $U=[-1,1]\times[-1,1]$ if you prefer)

My first attempt was using a function $f(x,y) := (a + bx + cy +dxy, a' + b'x+c'y +d'xy)$ (known as the 2d perspecitve transformation bilinear interpolation), but determining the coefficients $a,b,c,\ldots, d'$ requires inveriting $4\times 4$ matrix to solve two linear systems of equations.

EDIT: The actual 2d perspective transformation as described here does only produce the desired result if $Q$ is convex, which is not necessarily the case.

Answer 1

The first attempt was bound to fail, as the bilinear transformation suggested only perserves axis aligned lines when mapping a rectangle to some general quadrilateral. The inverse does not have the same form but it is not too difficult to compute, but still not particularly easy. Just for completeness sake here is how this works: Let us first focus on the forward mapping:

Forward mapping: We first define $g: \mathbb R^2 \to \mathbb R^2$ with $g|_U : U \to Q$ with the bilinear ansatz. It turns out we can easily write this as

$$\begin{align} (u,v) := g(x,y) &= (1-x)(1-y)P_0 + x(1-y)P_1 + x y P_2 + (1-x)yP_3 \\ &= \underbrace{P_0}_{=(a,a')} + \underbrace{[-P_0+P_1]}_{=(b,b')}x + \underbrace{[-P_0+P_3]}_{=(c,c')}y + \underbrace{[P_0-P_1+P_2-P_3]}_{=(d,d')}xy \end{align}$$

and we immediately find the coefficients $a,b,\ldots,d'$.

Inverse mapping: The inverse $f:\mathbb R^2 \to \mathbb R^2$ with $f|_Q : Q \to U$ of $g$ can then be computed by isolating $x$ from the equation for $u$:

$$x = \frac{u-a-cy}{b+dy}$$

Pluggin this back into the quation for $v$ and multiplying by the denominator we get

$$v(b+dy) = a' (b+dy) + b'(u-a-cy) + c' y(b+dy) + d' (u-a-cy)$$

which we can rearrange as

$$Ay^2+By+C = 0$$

with

$$\begin{align*} A &= b(a'-v)+b'(u-a) \\ B &= d(a'-v)+d'(u-a) + bc' - cb' \\ C &= dc'-cd' \end{align*}.$$

This lets us compute $y$ and therefore $x$ for some given $(u,v)$.

EDIT: As @amd commented, this doesn't work for nonconvex quadrilaterals, here the plot of such an example:

Answer 2

I expect that we have rather different thresholds for “easy to compute.” A small matrix inversion or two might get a bit tedious, but it’s not very difficult. Both a bilinear map and a planar homography are fairly easy to compute when one of the quadrilaterals is the unit square. Unfortunately, they only meet the criteria of the question when $Q$ is convex.

To handle the general case, I think it’s simplest to split the source quad $Q$ into a pair of triangles, and correspondingly split $U$ along the diagonal from the origin to $(1,1)$. It’s then a simple matter to construct a pair of affine transformations of these $Q$-halves onto the $U$-halves such that they agree along the diagonal.

W.l.o.g. let’s say that we split $Q$ along $P_1P_3$, so that we have $L=\triangle{P_1P_2P_3}$ and $R=\triangle{P_1P_3P_4}$. Let $H$ be the half-plane, including the boundary line $P_1P_3$, in which $P_2$ lies and construct $$\mathtt M_L = \begin{bmatrix}0&1&1 \\ 0&0&1 \\ 1&1&1\end{bmatrix}\begin{bmatrix}P_1&P_2&P_3\\1&1&1\end{bmatrix}^{-1} \\ \mathtt M_R = \begin{bmatrix}0&1&0 \\ 0&1&1 \\ 1&1&1\end{bmatrix}\begin{bmatrix}P_1&P_3&P_4\\1&1&1\end{bmatrix}^{-1}.$$ The matrix multiplications here are almost trivial to compute, as the resulting columns are simple sums of the columns of the right-hand matrix. Working in homogeneous coordinates, the map is then defined piecewise: $$f(P) = \begin{cases}\mathtt M_L P &\text{if }P\in H \\[2ex] \mathtt M_R P &\text{otherwise.}\end{cases}$$ Use whatever method is convenient to determine the half-plane in which $P$ lies, such as checking the sign of $P\cdot(P_3\times P_1)$ (again in homogeneous coordinates).

Or, if you prefer to work with algebraic expressions, construct the inverse map: $$f^{-1}(u,v) = \begin{cases} P_1+u(P_2-P_1)+v(P_3-P_2) &\text{if }u\le v \\[2ex] P_1+u(P_3-P_4)+v(P_4-P_1) &\text{otherwise}\end{cases}$$ and invert those formulas.