3
\$\begingroup\$

The 3D ellipse is represented by its center and the vectors for its major and minor axes and the 2D ellipse by its center, the length of the major and minor axes and the angle between the major axis and the x-axis. I need the general projection on an axis-aligned plane, discarding a coordinate.

The problem is that the axes after projection are no longer are aligned to the true axes of the 2D ellipse.

\$\endgroup\$
0

2 Answers 2

3
\$\begingroup\$

Let's express a point \$\vec p\$ on your 3D ellipse as a parametric equation, using the vectors \$ \vec c\$ for the center, \$ \vec M \$ for the semi-major axis and \$ \vec m \$ for the semi-minor axis, and a parameter \$ 0 \le t \le 2 \pi\$:

$$\vec p(t) = \cos(t) \vec M + \sin(t) \vec m + \vec c$$

From this we can derive three points on your 2D ellipse, at the tip of each of the 3D ellipse's major axes, and one "diagonal" point between them:

$$ (M_x + c_x, M_y + c_y)\\ (m_x + c_x, m_y + c_y)\\ (d_x, d_y) = (\frac {M_x + m_x} {\sqrt{2}} + c_x, \frac {M_y + m_y} {\sqrt{2}} + c_y)$$

These will be three solutions to the general equation of an ellipse centered at a point \$\vec c\$:

$$ A (x - c_x)^2 + B (x - c_x) (y - c_y) + C (y - c_y)^2 = 1$$

We can clear out all instances of \$\vec c\$ to shift our problem to the origin for now (and add it back later). We can then express this equation and our three known solutions as a system of linear equations:

$$\begin{bmatrix} M_x^2 & M_x M_y & M_y^2\\ m_x^2 & m_x m_y & m_y^2\\ d_x^2 & d_x d_y & d_y^2\\ \end{bmatrix} \begin{bmatrix} A\\B\\C \end{bmatrix} = \begin{bmatrix} 1\\1\\1 \end{bmatrix} $$

You can solve this system for the unknowns \$A\$, \$B\$, \$C\$ (eg. using Cramer's Rule or other solving techniques of your liking) to get the equation of the 2D ellipse.

Taking the formulae from Wikipedia, that leaves us with...

$$\begin{align} a,b&=\frac {-{\sqrt {2{\Big (}(B^{2}-4AC){\Big )}\left((A+C)\pm {\sqrt {(A-C)^{2}+B^{2}}}\right)}}}{B^{2}-4AC}\\ \\ \theta &= \begin{cases}\arctan \left({\frac {1}{B}}\left(C-A-{\sqrt {(A-C)^{2}+B^{2}}}\right)\right)&{\text{for }}B\neq 0\\0&{\text{for }}B=0,\ A<C\\90^{\circ }&{\text{for }}B=0,\ A>C\end{cases} \end{align}$$

...for the length of the semi-major axis \$a\$, semi-minor axis \$b\$, and the angle \$\theta\$ from the positive x axis to the major axis.

\$\endgroup\$
3
\$\begingroup\$

Christopher Brierley Jones describes the process on his website how an ellipse can be fit perfectly into any convex quadrilateral. The process is very similar to what DMGregory describes in his answer to extract the ellipse parameters from the matrix representation.

The key insight is that instead of projecting the ellipse directly, you can project the bounding rectangle. That will always result in an convex quadrilateral. Then you can use Jones process to fit an ellipse into the resulting quadrilateral.

You already mentioned that the axes are no longer aligned after the projection. Note that not even the center of the resulting ellipse aligns with the projected center of the circle. Jones describes the relation between the original circle and the projected ellipse as "inscrutable".

enter image description here

I implemented function for mapping a convex quad to an {cx,cy,ra,rb,angle} ellipse in JavaScript.

You can see a live demo here or a video here

// W,X,Y,Z are {x:number,y:number} vertices of the convex qudriliteral.
function quadToEllipse(W,X,Y,Z) {
    // Reconstruct matrix that transforms the unit square ((-1,-1), (1,1)) into quad (W,X,Y,Z)
    const m00 =  X.x * Y.x * Z.y - W.x * Y.x * Z.y - X.x * Y.y * Z.x + W.x * Y.y * Z.x - 
                W.x * X.y * Z.x + W.y * X.x * Z.x + W.x * X.y * Y.x - W.y * X.x * Y.x;
    const m01 =  W.x * Y.x * Z.y - W.x * X.x * Z.y - X.x * Y.y * Z.x + X.y * Y.x * Z.x - 
                W.y * Y.x * Z.x + W.y * X.x * Z.x + W.x * X.x * Y.y - W.x * X.y * Y.x;
    const m02 =  X.x * Y.x * Z.y - W.x * X.x * Z.y - W.x * Y.y * Z.x - X.y * Y.x * Z.x + 
                W.y * Y.x * Z.x + W.x * X.y * Z.x + W.x * X.x * Y.y - W.y * X.x * Y.x;
    const m10 =  X.y * Y.x * Z.y - W.y * Y.x * Z.y - W.x * X.y * Z.y + W.y * X.x * Z.y - 
                X.y * Y.y * Z.x + W.y * Y.y * Z.x + W.x * X.y * Y.y - W.y * X.x * Y.y;
    const m11 = -X.x * Y.y * Z.y + W.x * Y.y * Z.y + X.y * Y.x * Z.y - W.x * X.y * Z.y - 
                W.y * Y.y * Z.x + W.y * X.y * Z.x + W.y * X.x * Y.y - W.y * X.y * Y.x;
    const m12 =  X.x * Y.y * Z.y - W.x * Y.y * Z.y + W.y * Y.x * Z.y - W.y * X.x * Z.y - 
                X.y * Y.y * Z.x + W.y * X.y * Z.x + W.x * X.y * Y.y - W.y * X.y * Y.x;
    const m20 =  X.x * Z.y - W.x * Z.y - X.y * Z.x + W.y * Z.x - X.x * Y.y + W.x * Y.y + X.y * Y.x - W.y * Y.x;
    const m21 =  Y.x * Z.y - X.x * Z.y - Y.y * Z.x + X.y * Z.x + W.x * Y.y - W.y * Y.x - W.x * X.y + W.y * X.x;
    const m22 =  Y.x * Z.y - W.x * Z.y - Y.y * Z.x + W.y * Z.x + X.x * Y.y - X.y * Y.x + W.x * X.y - W.y * X.x;

    // invert matrix
    const determinant = +m00*(m11*m22-m21*m12) -m01*(m10*m22-m12*m20) +m02*(m10*m21-m11*m20);
    if(determinant == 0) return null;
    const invdet = 1/determinant; 
    const J =  (m11*m22-m21*m12)*invdet;
    const K = -(m01*m22-m02*m21)*invdet;
    const L =  (m01*m12-m02*m11)*invdet;
    const M = -(m10*m22-m12*m20)*invdet;
    const N =  (m00*m22-m02*m20)*invdet;
    const O = -(m00*m12-m10*m02)*invdet;
    const P =  (m10*m21-m20*m11)*invdet;
    const Q = -(m00*m21-m20*m01)*invdet;
    const R =  (m00*m11-m10*m01)*invdet;

    // extract ellipse coefficients from matrix
    const a = J*J + M*M - P*P;
    const b = J*K + M*N - P*Q;
    const c = K*K + N*N - Q*Q;
    const d = J*L + M*O - P*R;
    const f = K*L + N*O - Q*R;
    const g = L*L + O*O - R*R;

    // deduce ellipse center from coefficients
    const centerX = (c*d - b*f) / (b*b - a*c);
    const centerY = (a*f - b*d) / (b*b - a*c);

    // deduce ellipse radius from coefficients
    const radiusA = Math.sqrt(2*(a*f*f + c*d*d + g*b*b - 2*b*d*f - a*c*g)/((b*b - a*c) * (Math.sqrt((a-c)*(a-c) + 4*b*b) - (a+c))));
    const radiusB = Math.sqrt(2*(a*f*f + c*d*d + g*b*b - 2*b*d*f - a*c*g)/((b*b - a*c) * (-Math.sqrt((a-c)*(a-c) + 4*b*b) - (a+c))));

    // deduce ellipse rotation from coefficients
    let angle = 0;
    if(b==0 && a <= c) {
        angle = 0;
    } else if(b == 0 && a >= c) {
        angle = Math.PI / 2;
    } else if(b != 0 && a > c) {
        angle = Math.PI / 2 + 0.5 * (Math.PI / 2 - Math.atan2((a-c), (2*b)));
    } else if(b != 0 && a <= c) {
        angle = Math.PI / 2 + 0.5 * (Math.PI / 2 - Math.atan2((a-c), (2*b)));
    }
            
    return {
        centerX, centerY, radiusA, radiusB, angle
    }

Update: See this link for an example where this technique is used for drawing ellipses onto the faces of a rotating 3d cube as shown below:

ellipses on the sides of a 3d cube

\$\endgroup\$
0

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .