Finding the missing centre or vertex of a projected rectangle.

Question

In the attached image BCDE is a rectangular object that is projected through focal point J to produce a perspective image with vertices B'C'D'E' on the vertical blue plane. The plane of the rectangle is not coplanar to the blue screen. The 2D coordinates of points B'C'D' on the blue plane are known. The final vertex of the projected image is not know and two dashed sides D'Q and B'Q are completely obscured or missing so we cannot simply extend those sides to find the missing corner(Q).

The 3D Geogebra construction of the setup is here.

The orientation of the original rectangle, can be rotated around any axis using the sliders in the left hand menu. There are also sliders that can be used to alter the location of the focal point. (You might need to scroll the list on the left to find the sliders).

The problem is to find the 2D coordinates on the blue plane, of vertex Q or the projected centre $Z_1$. Finding either of these points would provide a 4th 2D coordinate that is the minimum required to carry out a standard homographic matrix transformation. The preferred solution is one that does not require knowledge of the 3D focal point location, if that is possible.

This question is closely related to and inspired by this question.

Answer 1

In addition to the $2D$ coordinates of $B', C', D'$, we need the $3D$ coordinates of the focal point (point $J$), AND the ratio of the side lengths $r =\dfrac{ | BC | }{ | C D | } $.

With this, points $B, C, D$ can be expressed as

$ B = J + \lambda_1 (JB') $

$ C = J + \lambda_2 (JC') $

$ D = J + \lambda_3 (JD') $

where $\lambda_1, \lambda_2, \lambda_3 \lt 0 $

Now impose the orthogonality condition between the segment $BC$ and the segment $CD$, then

$ ( \lambda_2 (JC') - \lambda_1 (JB') ) \cdot ( \lambda_3 (JD') - \lambda_2 (JC') ) = 0 $

The second equation comes from the ratio $r$, and we can write

$ ( \lambda_2 (JC') - \lambda_1 (JB') ) = r \| \lambda_3 (JD') - \lambda_2 (JC') \| $

which when squaring becomes:

$ ( \lambda_2 (JC') - \lambda_1 (JB') ) \cdot ( \lambda_2 (JC') - \lambda_1 (JB') ) = r^2 ( \lambda_3 (JD') - \lambda_2 (JC') ) \cdot ( \lambda_3 (JD') - \lambda_2 (JC') ) $

These two equations define $(\lambda_1, \lambda_2, \lambda_3)$ up to a scale. So for our purposes, we take $\lambda_1 = -1 $ and solve for $\lambda_2$ and $\lambda_3$.

This solution can be obtained in the Mathematica environment for example, by using Mathematica script. There can be more than one solution, but only one of them is valid. Identifying the valid solution is possible if you know the direction of the normal vector to the plane of the original rectangle. If not, then any set $\lambda_1 \lt 0, \lambda_2 \lt 0 , \lambda_3 \lt 0 $ is valid. I did not use Mathematica, instead I used VBA script in the Microsoft Excel Macro environment. I've include my VBA script program below.

Once we've obtained a valid set of $\lambda_1, \lambda_2, \lambda_3$, (they must all be negative) then we now have the coordinates of points $B, C, D$. Point $E$ can be constructed as follows

$ E = B + D - C $

Now consider the line

$ L(t) = J + \mu (J E) $

We want to intersect this line with the projection plane that contains $B', C', D'$. We already have the equation of this plane. It is of the form

$ n \cdot ( p - p_0) = 0 $

So we can determine $\mu$ by plugging $L(t)$ into the plane equation

$ n \cdot ( J - p_0 + \mu (JE) ) = 0 $

And this will give $\mu$ (which should be a negative value). Using this value for $\mu$, the coordinates of $E'$ are given by

$ E' = J + \mu (JE) $

Below is a VBA script code listing of the program that I used to find $\lambda_1, \lambda_2, \lambda_3$ in an implementation of the above setup and procedure. This program, when run, demonstrates that the method described above works and finds the fourth unknown vertex in the image of the rectangle. The source code for the function "intersect_three_quadrics" can be found in This Excel file. This function solves a system of three quadratic equations in three unknowns. The file contains other functions that solve a system of two quadratic equations in two unknowns. These are included in this file as macros (VBA script). Click on the link, to open the online file, then click on "Editing" and choose "Open in Desktop App". This will open the file in your desktop Excel program. Click "View" then select "Macros".

Public Sub fourth_vertex_in_the_projection()
Dim p(4, 3), q(4, 3), f(3), n(3), p0(3) As Double
Dim u1(3), u2(3), u3(3), u4(3), u5(3) As Double
Dim r As Double
Dim q1(4, 4, 3), q1c(4, 4, 3), ipoints(8, 3) As Double
Dim cp(4, 3) As Double
Dim v1(3), v2(3), v3(3), v4(3), n1(3), n2(3) As Double

' p's are the vertices of the rectangle,  q's are their projection through the focal point f.

' n is the normal to the projection plane, p0 is a point on that plane.

p(1, 1) = 25
p(1, 2) = 16
p(1, 3) = 200

u1(1) = 4 / 2
u1(2) = -8 / 2
u1(3) = 1 / 2

u2(1) = -4
u2(2) = -1
u2(3) = 8

r = 0.5

For i = 1 To 3
p(2, i) = p(1, i) + u1(i)
p(3, i) = p(2, i) + u2(i)
p(4, i) = p(3, i) - u1(i)
Next i

f(1) = 0
f(2) = 0
f(3) = 0

p0(1) = -30
p0(2) = -30
p0(3) = -10

n(1) = 2
n(2) = 3
n(3) = 4


' first thing, is to find the projections of p's and store them in the q's

' L(t) = f + t (f p)

For i = 1 To 4

For j = 1 To 3
u4(j) = p(i, j) - f(j)
Next j

' n . ( L - p0 ) = 0

t = (dot(n, p0) - dot(n, f)) / dot(n, u4)

For j = 1 To 3
q(i, j) = f(j) + t * u4(j)
Next j

Next i

' Next, we will use q(1,.), q(2,.), q(3, .) to write orthogonality condition

For j = 1 To 3
u1(j) = q(1, j) - f(j)
u2(j) = q(2, j) - f(j)
u3(j) = q(3, j) - f(j)
Next j

' now  ( t2 u2 - t1 u1 ) . ( t3 u3 - t2 u2) = 0

' so

'  t2 t3 (u2.u3) - t2^2 (u2.u2) - t1 t3 (u1.u3) + t1 t2 (u1.u2) = 0

' first entry in Q is  t1 = -1

q1(1, 4, 1) = 0.5
q1(4, 1, 1) = q1(1, 4, 1)

q1(4, 4, 1) = 1

' second entry

q1(2, 3, 2) = 0.5 * dot(u2, u3)
q1(3, 2, 2) = q1(2, 3, 2)

q1(2, 2, 2) = -dot(u2, u2)

q1(1, 3, 2) = -0.5 * dot(u1, u3)
q1(3, 1, 2) = q1(1, 3, 2)

q1(1, 2, 2) = 0.5 * dot(u1, u2)
q1(2, 1, 2) = q1(1, 2, 2)


' third entry

' ( t2 u2 - t1 u1)^2 = r^2 ( t3 u3 - t2 u2)^2

'  so that   t2^2 (u2.u2) + t1^2 (u1.u1) - 2 t1 t2 (u1.u2) =

'              r^2 ( t3^2 (u3.u3) + t2^2 (u2.u2) - 2 t2 t3 (u2.u3) )


q1(2, 2, 3) = (1 - r ^ 2) * dot(u2, u2)

q1(1, 1, 3) = dot(u1, u1)

q1(1, 2, 3) = -dot(u1, u2)

q1(2, 1, 3) = q1(1, 2, 3)

q1(3, 3, 3) = -(r ^ 2) * dot(u3, u3)

q1(2, 3, 3) = (r ^ 2) * dot(u2, u3)

q1(3, 2, 3) = q1(2, 3, 3)


Call intersect_three_quadrics(q1, q1c, ipoints)

MsgBox ("no. of solutions =" + Str(ipoints(0, 0)))

For i = 1 To ipoints(0, 0)

If ipoints(i, 1) > 0 Or ipoints(i, 2) > 0 Or ipoints(i, 3) > 0 Then GoTo lnext_i

t1 = ipoints(i, 1)
t2 = ipoints(i, 2)
t3 = ipoints(i, 3)

' now construct the vertices of the scaled rectangle

For j = 1 To 3

cp(1, j) = f(j) + t1 * u1(j)
cp(2, j) = f(j) + t2 * u2(j)
cp(3, j) = f(j) + t3 * u3(j)

v1(j) = cp(1, j) - cp(2, j)
v2(j) = cp(3, j) - cp(2, j)

v3(j) = p(1, j) - p(2, j)
v4(j) = p(3, j) - p(2, j)

Next j

' check the normal vector

Call cross(v1, v2, n1)

Call cross(v3, v4, n2)

Call normalize1(n1)

Call normalize1(n2)

For j = 1 To 3
ActiveSheet.Cells(j, 4) = n1(j)
ActiveSheet.Cells(j, 5) = n2(j)
Next j

' now construct the fourth vertex

For j = 1 To 3
  cp(4, j) = cp(1, j) + cp(3, j) - cp(2, j)
Next j

' finally , project cp(4) to the projection plane

' n . ( L - p0 ) = 0

For j = 1 To 3
u4(j) = cp(4, j) - f(j)
Next j

t = (dot(n, p0) - dot(n, f)) / dot(n, u4)

For j = 1 To 3
u5(j) = f(j) + t * u4(j)
Next j

For j = 1 To 3
ActiveSheet.Cells(j, 1) = q(4, j)
ActiveSheet.Cells(j, 2) = u5(j)
Next j

MsgBox ("")

lnext_i:
Next i

End Sub