Given orthogonal vibrations, how can I find the magnitude and direction of the major axis of the resulting ellipse

conic sectionsgeometryparametrictrigonometry

The Context

In my work, I am using a vibration table to test the resonant properties and survivability of a structure under different shaking regimes. The table is driven by an eccentrically weighted rotary motor which excites vibrations in the Y and Z axes (green and blue respectively in the image below).

vibration table with axes: x-red, y-green, z-blue

Using 3-axis accelerometers on the table and at various points on the structure, I can record the excitation and response vibrations in each axis (example plot at 20Hz below).

example plot with Y and Z vibrations

However, two sinusoids of differing amplitudes and phases define an elliptical polarization (the same data as above is plotted parametrically below).

parametric plot of the Y and Z vibrations

The Question

Based on the answers to other Stack Exchange questions (here and here), I know that the parametric form of an ellipse, centered at the origin and rotated by an angle $\theta$ is
$$ a_y(t) = A\cos\big(2\pi ft\big)\cos\big(\theta\big) – B\sin\big(2\pi ft\big)\sin\big(\theta\big) $$
$$ a_z(t) = A\cos\big(2\pi ft\big)\sin\big(\theta\big) + B\sin\big(2\pi ft\big)\cos\big(\theta\big) $$
where $A$ and $B$ are the major and minor radii respectively.

My question is, how can I convert my axes-aligned sinusoids
$$ a_y(t) = Y\sin\big(2\pi ft\big) $$
$$ a_z(t) = Z\sin\big(2\pi ft + \phi\big) $$
into the above form, such that I can immediately read off the rotation angle and radii?

My goal is to know the direction and magnitude of the major radius of the acceleration ellipse, so that I can better align my structure in the direction of strongest vibration.

What I've Tried So Far

Based on the process discussed in this article, I’ve attempted to simply set the two forms equal to one another and solve for $A$, $B$, and $\theta$ as functions of $Y$, $Z$, and $\phi$. However, that process has yielded inconsistent results, and none of them seem to agree with the ellipse as I’ve plotted it.

I also have a gut feeling that the major and minor radii correspond to the eigenvalues/eigenvectors of some transformation matrix, but I have no idea how to go about constructing said matrix.

Numerically, the sinusoids I've plotted above have parameters
$ Y = 0.7847g $, $ Z = 1.0498g $, and $ \phi = 0.9537 $. The major and minor radii are approximately $R = 1.1805g$ and $r = 0.5696g$ respectively, and the rotation angle is approximately $\theta = 1.0330$.

Best Answer

It’s fairly easy to convert your parametric equations to an implicit Cartesian equation. Once you have that, you can use a variety of methods to extract the ellipse’s parameters. To reduce clutter a bit, I’ll set $\omega = 2\pi f$.

First, expand the second equation using the formula for the sine of a sum and isolate $\cos{\omega t}$ (don’t worry about zero denominators at this point): $$\cos{\omega t} = {a_z/Z - \sin{\omega t} \cos\phi \over \sin\phi}.$$ From the first equation we have $\sin{\omega t} = a_y/Y$. Substitute into the above and then use the trig identity $\cos^2{\omega t} + \sin^2{\omega t} = 1$ to eliminate $\omega t$: $$\left({a_z/Z - a_y/Y \cos\phi \over \sin\phi}\right)^2 + \left(a_y/Y\right)^2 = 1,$$ then rearrange into the more standard-looking $$ {a_y^2 \over Y^2} - 2{\cos\phi \over YZ} a_y a_z + {a_z^2 \over Z^2} = \sin^2\phi.$$ If you like, you can multiply through by $Y^2Z^2$ to eliminate all of the denominators.

From here you could proceed by constructing the corresponding symmetric matrix and then computing its eigensystem, but you can save yourself a quite a bit of work by using fairly well-known formulas. In terms of the general conic equation $Ax^2+Bxy+Cy^2+Dx+Ey+F=0$ we have $$A = \frac1{Y^2} \\ B = -2{\cos\phi \over YZ} \\ C = \frac1{Z^2} \\ D=E=0 \\ F=-\sin^2\phi.$$ Set $\Delta = B^2-4AC$. Using the formulas from the ellipse article on Wikipedia we have $$\tan\Theta = {C-A-\sqrt{(A+C)^2+\Delta} \over B}$$ for the inclination of the ellipse’s major axis from the $x$-axis. With that angle in hand, you can determine the semiaxis lengths by applying the inverse rotation, intersecting the line with this direction and its perpendicular with the ellipse, or just using the formulas from the same source. Most of the terms vanish, leaving $$a,b = -{\sqrt{2\Delta F \left(A+C \pm \sqrt{(A+C)^2+\Delta}\right)} \over \Delta}.$$ I’ll leave plugging in the values of these coefficients and simplifying to you.

Another way to attack this is to use a somewhat different parameterization of an ellipse: $$\mathbf c + \mathbf u \cos\tau + \mathbf v \sin\tau,$$ where $\mathbf u$ and $\mathbf v$ are conjugate radii of the ellipse and $\mathbf c$ is its center (in this case the origin). The two vectors are conjugate if $\mathbf v$ is parallel to $\dot{\mathbf u}$, or vice-versa. Setting $\mathbf u = [Y\sin{\omega t_1},Z\sin(\omega t_1+\phi)]$ and $\mathbf v = [Y\sin{\omega t_2},Z\sin(\omega t_2+\phi)]$, that constraint can be expressed as $$\begin{align} \det\begin{bmatrix}\dot{\mathbf u} \\ \mathbf v\end{bmatrix} &= \det\begin{bmatrix} \omega Y \cos{\omega t_1} & \omega Z \cos(\omega t_1+\phi) \\ Y\sin{\omega t_2} & Z\sin(\omega t_2+\phi) \end{bmatrix} \\ &= \omega YZ \left( \cos{\omega t_1}\sin(\omega t_2+\phi) - \cos(\omega t_1+\phi)\sin{\omega t_2} \right) = 0. \end{align}$$ When $\mathbf u$ and $\mathbf v$ are orthogonal, they coincide with the semiaxes of the ellipse. Thus, we want $$\begin{align} \mathbf u\cdot\mathbf v &= [Y\sin{\omega t_1},Z\sin(\omega t_1+\phi)]\cdot[Y\sin{\omega t_2},Z\sin(\omega t_2+\phi)] \\ &= Y^2 \sin{\omega t_1}\sin{\omega t_2} + Z^2 \sin(\omega t_1+\phi)\sin(\omega t_2+\phi) = 0.\end{align}$$ In principle, you can solve these equations for $t_1$ and $t_2$, which gives you both the principal axes of the ellipse ($\mathbf u$ and $\mathbf v$) and its semiaxis lengths (their norms).