1
$\begingroup$

Given an ellipse in a 3-dimenion system defined by the center at $C$ with the major semi-axis $\vec{a}$ and minor semi-axis $\vec{b}$, I want to find the minimum and maximum points with respect to each axis such that I can evaluate the dimensions of the smallest cuboid (perpendicular to the Cartesian coordinate system) inscribing the ellipse.

I thought of setting a moving point $P$ such that the ellipse is the locus of the $P$, so by changing the angle $\theta$ between $\vec{a}$ and $\vec{CP}$, the coordinates of $P$ can be evaluated in terms of $\theta$ (as $i(\theta)$, $j(\theta)$, $k(\theta)$),and hence the maximum and minimum points can be found with $\frac{d}{d\theta}i(\theta)$, $\frac{d}{d\theta}j(\theta)$, $\frac{d}{d\theta}k(\theta)$. However, $i(\theta)$ has to be found with the polar ellipse equation with $$\vec{CP}=\frac{|\vec{a}||\vec{b}|}{\sqrt{|\vec{a}\cos{\theta}|^2+|\vec{b}\sin\theta|^2}}\times \frac{\vec{a}\tan\theta+\vec{b}}{|\vec{a}\tan\theta+\vec{b}|}$$Expanding $\vec{a}$ and $\vec{b}$ into the three dimensions, the equation becomes very complicated.

Would there be a simpler method doing so?

$\endgroup$
2
$\begingroup$

Yes, there is a much simpler method.

If you have semi-axes $\vec{a} = ( x_a , y_a , z_a )$ and $\vec{b} = ( x_b , y_b , z_b )$, you can parametrize the ellipse using $t \in \mathbb{R}$, $-1 \le t \le 1$: $$\vec{p}(t) = \vec{a} t \pm \vec{b} \sqrt{1 - t^2} \tag{1}\label{1}$$ The $+$ sign refers to the half of the ellipse that is on the same side as $\vec{b}$ with respect to the center of the ellipse, and $-$ to the other half.

We can find the minimum and maximum $x$ by solving $$\frac{\partial \vec{p}(t)}{\partial x} = 0 \tag{X}\label{X}$$ for $t$; for $y$ by solving $$\frac{\partial \vec{p}(t)}{\partial y} = 0 \tag{Y}\label{Y}$$ for $t$; and for $z$ by solving $$\frac{\partial \vec{p}(t)}{\partial z} = 0 \tag{Z}\label{Z}$$ for $t$.

Equation $\eqref{X}$ yields $$t = \pm \frac{x_a}{\sqrt{x_a^2 + x_b^2}}$$ equation $\eqref{Y}$ yields $$t = \pm \frac{y_a}{\sqrt{y_a^2 + y_b^2}}$$ and equation $\eqref{Z}$ yields $$t = \pm \frac{z_a}{\sqrt{z_a^2 + z_b^2}}$$ Note that in all three cases we can see that $-1 \le t \le +1$; in other words, all the values of $t$ we thus obtain are valid.

This leaves you four candidate values per coordinate you need to choose from. For example, if we use $$\begin{array}{l} t_{+} = + \frac{x_a}{\sqrt{x_a^2 + x_b^2}} \\ t_{-} = - \frac{x_a}{\sqrt{x_a^2 + x_b^2}} \end{array}$$ then $$\begin{cases} x_1 = x_a t_{+} + x_b \sqrt{1 - t_{+}^2} \\ x_2 = x_a t_{+} - x_b \sqrt{1 - t_{+}^2} \\ x_3 = x_a t_{-} + x_b \sqrt{1 - t_{-}^2} \\ x_4 = x_a t_{-} - x_b \sqrt{1 - t_{-}^2} \end{cases} \qquad \iff \qquad x_i = \frac{ \pm x_a^2 \pm x_b^2 }{\sqrt{x_a^2 + x_b^2}}$$ Note that by our definition in equation $\eqref{1}$, and because $-1 \le t_{+} \le +1$ and $-1 \le t_{-} \le +1$, all four refer to valid $x$ coordinates on the ellipsoid. One is the minimum, and another the maximum: $$\begin{cases} x_{min} = \min \left ( x_1 , x_2 , x_3 , x_4 \right ) = \frac{- x_a^2 - x_b^2}{\sqrt{x_a^2 + x_b^2}} = -\sqrt{x_a^2 + x_b^2} \\ x_{max} = \max \left ( x_1 , x_2 , x_3 , x_4 \right ) = \frac{x_a^2 + x_b^2}{\sqrt{x_a^2 + x_b^2}} = \sqrt{x_a^2 + x_b^2} \end{cases}$$

The same applies to $y$ and $z$ axes, too.

Therefore, if the center of the ellipse is at $( x_0 , y_0 , z_0 )$, then $$\begin{cases} x_{min} = x_0 - \sqrt{x_a^2 + x_b^2} \\ y_{min} = y_0 - \sqrt{y_a^2 + y_b^2} \\ z_{min} = z_0 - \sqrt{z_a^2 + z_b^2} \end{cases}, \qquad \begin{cases} x_{max} = x_0 + \sqrt{x_a^2 + x_b^2} \\ y_{max} = y_0 + \sqrt{y_a^2 + y_b^2} \\ z_{max} = z_0 + \sqrt{z_a^2 + z_b^2} \end{cases}$$

Note that the coordinates $\left ( x_{min} , y_{min} , z_{min} \right )$ and $\left ( x_{max} , y_{max} , z_{max} \right )$ specify the axis-aligned cuboid that encloses the ellipsoid, not the cuboid that can fit inside the ellipsoid.


There is no need to take my word for it, obviously. As an example, here is a simple awk script you can use to test it. check.awk:

#!/usr/bin/awk -f

BEGIN {
    if (N < 1)
        N = 1000

    srand()

    while (1) {
        xa = 2.0*rand() - 1.0
        ya = 2.0*rand() - 1.0
        za = 2.0*rand() - 1.0
        aa = xa*xa + ya*ya + za*za
        if (aa > 0.25 && aa < 1.0) {
            a = sqrt(aa)
            xa = xa / a
            ya = ya / a
            za = za / a
            break
        }
    }

    while (1) {
        xb = 2.0*rand() - 1.0
        yb = 2.0*rand() - 1.0
        zb = 2.0*rand() - 1.0
        bb = xb*xb + yb*yb + zb*zb
        if (bb > 0.25 && bb < 1.0) {
            b = sqrt(bb)
            xb = xb / b
            yb = yb / b
            zb = zb / b
            ab = xa*xb + ya*yb + za*zb
            if (ab > -0.25 && ab < 0.25) {
                xb = xb - ab*xa
                yb = yb - ab*ya
                zb = zb - ab*za
                b = sqrt(xb*xb + yb*yb + zb*zb)
                xb = xb / b
                yb = yb / b
                zb = zb / b
                break
            }
        }
    }

    a = 0.5 + 9.0 * rand()
    xa = xa * a
    ya = ya * a
    za = za * a

    b = 0.5 + 9.0 * rand()
    xb = xb * b
    yb = yb * b
    zb = zb * b

    xe = sqrt(xa*xa + xb*xb)
    ye = sqrt(ya*ya + yb*yb)
    ze = sqrt(za*za + zb*zb)

    xmin = 0
    ymin = 0
    zmin = 0

    xmax = 0
    ymax = 0
    zmax = 0

    for (i = -N; i <= N; i++) {
        ta = i/N
        tb = sqrt(1.0 - ta*ta)

        x = ta*xa + tb*xb
        y = ta*ya + tb*yb
        z = ta*za + tb*zb
        if (x < xmin) xmin = x
        if (x > xmax) xmax = x
        if (y < ymin) ymin = y
        if (y > ymax) ymax = y
        if (z < zmin) zmin = z
        if (z > zmax) zmax = z

        printf "%12.9f %12.9f %12.9f\n", x, y, z
    }

    for (i = N-1; i >= -N; i--) {
        ta = i/N
        tb = -sqrt(1.0 - ta*ta) 

        x = ta*xa + tb*xb
        y = ta*ya + tb*yb
        z = ta*za + tb*zb
        if (x < xmin) xmin = x
        if (x > xmax) xmax = x
        if (y < ymin) ymin = y
        if (y > ymax) ymax = y
        if (z < zmin) zmin = z
        if (z > zmax) zmax = z

        printf "%12.9f %12.9f %12.9f\n", x, y, z
    }

    printf "- - -  0 0 0  0 0 0\n"
    printf "- - -  %12.9f %12.9f %12.9f  %12.9f %12.9f %12.9f\n", xa, ya, za, xb, yb, zb

    printf "# a = (%9.6f, %9.6f, %9.6f)\n", xa, ya, za
    printf "# b = (%9.6f, %9.6f, %9.6f)\n", xb, yb, zb
    printf "# a.a = %12.9f\n", xa*xa + ya*ya + za*za
    printf "# a.b = %12.9f\n", xa*xb + ya*yb + za*zb
    printf "# b.b = %12.9f\n", xb*xb + yb*yb + zb*zb

    printf "# Extrema:   %9.6f   %9.6f   %9.6f\n", xe, ye, ze
    printf "# Minimum:   %9.6f   %9.6f   %9.6f\n", xmin, ymin, zmin
    printf "# Maximum:   %9.6f   %9.6f   %9.6f\n", xmax, ymax, zmax

    printf "a = (%9.6f, %9.6f, %9.6f)\n", xa, ya, za > "/dev/stderr"
    printf "b = (%9.6f, %9.6f, %9.6f)\n", xb, yb, zb > "/dev/stderr"
    printf "a.a = %12.9f\n", xa*xa + ya*ya + za*za > "/dev/stderr"
    printf "a.b = %12.9f\n", xa*xb + ya*yb + za*zb > "/dev/stderr"
    printf "b.b = %12.9f\n", xb*xb + yb*yb + zb*zb > "/dev/stderr"

    printf "Extrema:   %9.6f   %9.6f   %9.6f\n", xe, ye, ze > "/dev/stderr"
    printf "Minimum:   %9.6f   %9.6f   %9.6f\n", xmin, ymin, zmin > "/dev/stderr"
    printf "Maximum:   %9.6f   %9.6f   %9.6f\n", xmax, ymax, zmax > "/dev/stderr"

    exit(0)
}

It generates two random orthogonal unit vectors, $\vec{a}$ and $\vec{b}$, scales them to a random length between $0.5$ and $9.5$, and then computes four thousand points (by default) on the ellipse (centered at origin), to find the actual minimum and maximum coordinates. Those points are printed to standard output, so you can verify it is an ellipse. The test results are also printed to standard error.

Basically, I run the above using

./check.awk > out && (printf 'splot "out" u 1:2:3 notitle w lines lc -1, "out" u 4:5:6 notitle w lines lc 1, "out" u 7:8:9 notitle w lines lc 3\n'; read line) | gnuplot

which lets me look at the ellipse interactively in Gnuplot.

Here is a typical result:

a = ( 0.545091, -0.541625,  1.059682)
b = ( 5.159619,  2.694286, -1.276956)
a.a =  1.713407215
a.b =  0.000000000
b.b = 35.511465778
Extrema:    5.188332    2.748188    1.659380
Minimum:   -5.188332   -2.748188   -1.659380
Maximum:    5.188332    2.748188    1.659380

This shows that the extrema, $\pm\sqrt{x_a^2 + x_b^2}$, $\pm\sqrt{y_a^2 + y_b^2}$, $\pm\sqrt{z_a^2 + z_b^2}$, are indeed correct.

$\endgroup$
  • $\begingroup$ Observation: $\frac{x_a^2+x_b^2}{\sqrt{x_a^2+x_b^2}}$ can be simplified into $\sqrt{x_a^2+x_b^2}$. The modulus of the vector from the center $\vec{c} = (x_0,y_0,z_0)$ to the maxiumum vertex of the cuboid $\vec{v_{max}} = (x_{max}, y_{max}, z_{max})$ is $\sqrt{x_a^2+y_a^2+z_a^2+x_b^2+y_b^2+z_b^2}$, which is equivalent to $\sqrt{|\vec{a}|^2 + |\vec{b}|^2}$. Since the major and minor axes are perpendicular by definition, $|\vec{v_{max}}-\vec{c}|$ is equal to $|\vec{a}+\vec{b}|$ (half of the diagonal of the rectangle inscribing the ellipse). Doesn't this feel a bit wrong? $\endgroup$ – SOFe Jun 4 '17 at 15:24
  • $\begingroup$ It's feeling wrong because the cuboid inscribing an obliquely positioned ellipse should increase in diagonal length when it becomes less upright, right? $\endgroup$ – SOFe Jun 4 '17 at 15:37
  • $\begingroup$ @SOFe: Thanks for the observation; I edited it into the answer. $\endgroup$ – Nominal Animal Jun 4 '17 at 19:45
  • 1
    $\begingroup$ @SOFe: $\lvert \vec{v}_{max} - \vec{c} \rvert = \sqrt{ x_a^2 + y_a^2 + z_a^2 + x_b^2 + y_b^2 + z_b^2 } = \sqrt{ \lvert \vec{a} \rvert^2 + \lvert \vec{b} \rvert^2 }$. $\vec{v}_{max} = \vec{c} + \left ( \sqrt{x_a^2 + x_b^2} ,\, \sqrt{y_a^2 + y_b^2} ,\, \sqrt{z_a^2 + z_b^2} \right )$ and $\vec{v}_{min} = \vec{c} - \left ( \sqrt{x_a^2 + x_b^2} ,\, \sqrt{y_a^2 + y_b^2} ,\, \sqrt{z_a^2 + z_b^2} \right )$. The diagonal of the rectangle is $\vec{v}_{max} - \vec{v}_{min} = \left ( 2 \sqrt{x_a^2 + x_b^2} ,\, 2 \sqrt{y_a^2 + y_b^2} ,\, 2 \sqrt{z_a^2 + z_b^2} \right )$. $\endgroup$ – Nominal Animal Jun 4 '17 at 19:55
0
$\begingroup$

With the ellipse in this form:

$e: \begin{cases} x=c_x+a_x\cos\theta+b_x\sin\theta\\ y=c_y+a_y\cos\theta+b_y\sin\theta&0\leq\theta\lt 2\pi\\ z=c_z+a_z\cos\theta+b_z\sin\theta \end{cases}$

You have simply to find the extrema for each coordinate separately. e.g.

$\dfrac{dx}{d\theta}=-a_x\sin\theta_x+b_x\cos\theta_x=0$

$\tan\theta_x=b_x/a_x$ Giving the two values for $x$ that you are looking for the cuboid to be constructed.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.