Ellipse rotating in space perimeter function

104 Views Asked by At

what is the function of ellipse perimeter depending on the angle of ellipse rotation in space? Let's start with a very simple one-axis (be it major semiaxis plane) rotation using the approximated perimeter equation:

l ≈ π(3/2 (a+b)-√ab).

My aim is to show how the perimeter value changes for an elliptical defect in heart when seen from different angles.

Please note that I am not a mathematician (I am MD), therefore I would greatly appreciate to keep things as simple as possible.

Thanks for your help, Jan

1

There are 1 best solutions below

9
On

The underlying problem is that when the ellipse is projected to a plane, the projected semiaxes are not the semiaxes of the projected ellipse. For example: Projected ellipse and actual and apparent semi-major axes The two ellipses are exactly the same. The blue lines on the left are the axes of the original ellipse. The red lines on the right are the axes of the projected ellipse.

If the projected semiaxes ("half-axes") are $\vec{a}$ and $\vec{b}$, the points $\vec{p}(\varphi)$ on the projected ellipse are $$\vec{p}(\varphi) = \vec{a}\cos\varphi + \vec{b}\sin\varphi$$

If we have a function $L(a, b)$ that yields the (approximate) perimeter of an ellipse with semiaxes $a$ and $b$, we need functions $u(\vec{a}, \vec{b}) \ge 0$ and $v(\vec{a}, \vec{b}) \ge 0$ that yield the actual lengths of the semiaxes of the projected ellipse. Then, the length of the perimeter of the projected ellipse is $L\left(u(\vec{a}, \vec{b}), v(\vec{a}, \vec{b})\right)$.

By finding the extrema of $\vec{p}(\varphi)$, I found the solutions r1 = $u(\vec{a}, \vec{b})$ and r2 = $v(\vec{a}, \vec{b})$, using Sagemath:

r1 = sqrt( ( ax^4*bx^2 + 2*ax^2*ay^2*bx^2 + ay^4*bx^2 + 2*ax^2*bx^4 - 2*ay^2*bx^4 + bx^6 + 8*ax*ay*bx^3*by + ax^4*by^2 + 2*ax^2*ay^2*by^2 + ay^4*by^2 + 3*bx^4*by^2 + 8*ax*ay*bx*by^3 - 2*ax^2*by^4 + 2*ay^2*by^4 + 3*bx^2*by^4 + by^6 + sqrt(ax^4 + 2*ax^2*ay^2 + ay^4 + bx^4 + 8*ax*ay*bx*by + by^4 + 2*(ax^2 - ay^2)*bx^2 - 2*(ax^2 - ay^2 - bx^2)*by^2)*(ax^2*bx^2 - ay^2*bx^2 + bx^4 + 4*ax*ay*bx*by - ax^2*by^2 + ay^2*by^2 + 2*bx^2*by^2 + by^4)
           ) / ( ax^4 + 2*ax^2*ay^2 + ay^4 + 2*ax^2*bx^2 - 2*ay^2*bx^2 + bx^4 + 8*ax*ay*bx*by - 2*ax^2*by^2 + 2*ay^2*by^2 + 2*bx^2*by^2 + by^4 - sqrt(ax^4 + 2*ax^2*ay^2 + ay^4 + bx^4 + 8*ax*ay*bx*by + by^4 + 2*(ax^2 - ay^2)*bx^2 - 2*(ax^2 - ay^2 - bx^2)*by^2)*(ax^2 + ay^2 - bx^2 - by^2)
           ) )
r2 = sqrt( ( ax^4*bx^2 + 2*ax^2*ay^2*bx^2 + ay^4*bx^2 + 2*ax^2*bx^4 - 2*ay^2*bx^4 + bx^6 + 8*ax*ay*bx^3*by + ax^4*by^2 + 2*ax^2*ay^2*by^2 + ay^4*by^2 + 3*bx^4*by^2 + 8*ax*ay*bx*by^3 - 2*ax^2*by^4 + 2*ay^2*by^4 + 3*bx^2*by^4 + by^6 - sqrt(ax^4 + 2*ax^2*ay^2 + ay^4 + bx^4 + 8*ax*ay*bx*by + by^4 + 2*(ax^2 - ay^2)*bx^2 - 2*(ax^2 - ay^2 - bx^2)*by^2)*(ax^2*bx^2 - ay^2*bx^2 + bx^4 + 4*ax*ay*bx*by - ax^2*by^2 + ay^2*by^2 + 2*bx^2*by^2 + by^4)
           ) / ( ax^4 + 2*ax^2*ay^2 + ay^4 + 2*ax^2*bx^2 - 2*ay^2*bx^2 + bx^4 + 8*ax*ay*bx*by - 2*ax^2*by^2 + 2*ay^2*by^2 + 2*bx^2*by^2 + by^4 + sqrt(ax^4 + 2*ax^2*ay^2 + ay^4 + bx^4 + 8*ax*ay*bx*by + by^4 + 2*(ax^2 - ay^2)*bx^2 - 2*(ax^2 - ay^2 - bx^2)*by^2)*(ax^2 + ay^2 - bx^2 - by^2)
           ) )

where $\vec{a} = ( \text{ax} , \text{ay} )$, $\vec{b} = ( \text{bx} , \text{by} )$, and ^ denotes exponentiation. It can be cleaned up some, but that is drudge work I'll leave to whoever wants to use it.

I tested a quick-and-horribly-dirty version of the above in C,

#include <math.h>

void semiaxes(const double ax, const double ay,
              const double bx, const double by,
              double length[2])
{
    const double  ax2 = ax*ax, ay2 = ay*ay, bx2 = bx*bx, by2 = by*by;
    const double  ax4 = ax2*ax2, ay4 = ay2*ay2, bx4 = bx2*bx2, by4 = by2*by2;
    const double  c1 = sqrt(ax4 + 2*ax2*ay2 + ay4 + bx4 + 8*ax*ay*bx*by + by4 + 2*(ax2 - ay2)*bx2 - 2*(ax2 - ay2 - bx2)*by2);
    const double  s1 = ( ax4*bx2 + 2*ax2*ay2*bx2 + ay4*bx2 + 2*ax2*bx4 - 2*ay2*bx4 + bx2*bx4 + 8*ax*ay*bx*bx2*by + ax4*by2 + 2*ax2*ay2*by2 + ay4*by2 + 3*bx4*by2 + 8*ax*ay*bx*by*by2 - 2*ax2*by4 + 2*ay2*by4 + 3*bx2*by4 + by2*by4 + c1*(ax2*bx2 - ay2*bx2 + bx4 + 4*ax*ay*bx*by - ax2*by2 + ay2*by2 + 2*bx2*by2 + by4)
                       ) / ( ax4 + 2*ax2*ay2 + ay4 + 2*ax2*bx2 - 2*ay2*bx2 + bx4 + 8*ax*ay*bx*by - 2*ax2*by2 + 2*ay2*by2 + 2*bx2*by2 + by4 - c1*(ax2 + ay2 - bx2 - by2) );
    const double  s2 = ( ax4*bx2 + 2*ax2*ay2*bx2 + ay4*bx2 + 2*ax2*bx4 - 2*ay2*bx4 + bx2*bx4 + 8*ax*ay*bx*bx2*by + ax4*by2 + 2*ax2*ay2*by2 + ay4*by2 + 3*bx4*by2 + 8*ax*ay*bx*by*by2 - 2*ax2*by4 + 2*ay2*by4 + 3*bx2*by4 + by2*by4 - c1*(ax2*bx2 - ay2*bx2 + bx4 + 4*ax*ay*bx*by - ax2*by2 + ay2*by2 + 2*bx2*by2 + by4)
                       ) / ( ax4 + 2*ax2*ay2 + ay4 + 2*ax2*bx2 - 2*ay2*bx2 + bx4 + 8*ax*ay*bx*by - 2*ax2*by2 + 2*ay2*by2 + 2*bx2*by2 + by4 + c1*(ax2 + ay2 - bx2 - by2) );
    length[0] = sqrt(s1);
    length[1] = sqrt(s2);
}

with a number of random semiaxis vectors. Calculating the perimeter lengths numerically, the results obtained by direct evaluation of the projected ellipse and those by calculating the perimeter of an axis-aligned ellipse where the semiaxis lengths were provided by semiaxes(), were in excellent agreement. In other words, this seems to work.