I have a level surface in $(x,y,z)$ space. For concreteness, let's say
$\frac{1}{r^2}x^2 + \left(\frac{\sin \theta}{r}\right)^2y^2 + \left(\frac{\cos \theta}{r}\right)^2z^2 + \frac{2 \cos \theta \sin \theta}{r^2} xy = 1$
which is describes a circle with radius $r$ initially in the $xz$ plane that is tilted by angle $\theta$. How can I get the level surface of the corresponding projected ellipse in the $xy$ plane?
In other words, if I have an observer at $z = \infty$, how can I describe the equation of the ellipse they will see?
So as I stated in my comment, your set $S$ defined by your equation equation (that I will rewrite for more convenience) $$x^2 + \sin^2(\theta) y^2 + \cos^2(\theta) z^2 + 2\cos(\theta) \sin (\theta) xy = r^2$$ Defines a $2$-dimensional smooth manifold in $\mathbb{R}^3$ (that will very probably look like a deformed sphere) so its projection parallel to the $z$ axis on the $x,y$ plane will be a filled ellipse.
Now the projection $p$ is defined to be $p((x,y,z)) = (x,y)$, thus if $\cos(\theta)\neq 0$ $$\begin{align*} p(S) &= \bigcup_{z \in \mathbb{R}}\left\{(x,y) \in \mathbb{R}^2\mid x^2 + \sin^2(\theta) y^2 + \cos^2(\theta) z^2 + 2\cos(\theta) \sin (\theta) xy = r^2 \right\}\\ &=\bigcup_{z' \geqslant 0}\left\{(x,y) \in \mathbb{R}^2\mid x^2 + \sin^2(\theta) y^2 + z' + 2\cos(\theta) \sin (\theta) xy = r^2 \right\}\\ &=\bigcup_{z' \geqslant 0}\left\{(x,y) \in \mathbb{R}^2\mid x^2 + \sin^2(\theta) y^2 + 2\cos(\theta) \sin (\theta) xy = r^2 - z' \right\}\\ &=\left\{(x,y) \in \mathbb{R}^2\mid x^2 + \sin^2(\theta) y^2 + 2\cos(\theta) \sin (\theta) xy \leqslant r^2 \right\} \end{align*}$$
Now let's rewrite this latter inéquation.
$$\begin{align*} &x^2 + \sin^2(\theta) y^2 + 2\cos(\theta) \sin (\theta) xy \leqslant r^2 \\ \iff& (1-\cos^2(\theta))x^2 + \cos^2(\theta) x^2 + \sin^2(\theta) y^2 + 2\cos(\theta) \sin (\theta) xy \leqslant r^2 \\ \iff& (\sin \theta \, x)^2 + (\cos \theta \, x + \sin \theta \, y)^2 \leqslant r^2 \\ \end{align*}$$
Now we have two cases, first let us suppose that $\sin \theta \neq 0$. This equation is just the equation of a disk in the basis $((\sin \theta ,0), (\cos \theta , \sin \theta)$. Basically the basis change matrix is $$\begin{pmatrix} \sin \theta & 0 \\ \cos \theta & \sin \theta \end{pmatrix}$$
Its inverse is $$A=\begin{pmatrix} \sin^{-1} \theta & 0 \\ -\frac{\cos \theta}{\sin^2 \theta} & \sin^{-1} \theta \end{pmatrix} = \frac{1}{sin^2 \theta} \begin{pmatrix} \sin \theta & 0 \\ -\cos \theta & \sin \theta \end{pmatrix}$$
Thus $p(S) = A\cdot D_r$ where $D$ denotes the disk of radius $r$ in $\mathbb{R}^2$, so a filled rotated ellipse.
Now if $\sin \theta = 0$, we have $(\cos\theta \, x)^2 \leqslant r^2$ i.e. $|x| \leqslant \frac{r}{|\cos\theta|}$ so your solution space is $\left[-\frac{r}{|\cos\theta|},\frac{r}{|\cos\theta|} \right] \times \mathbb{R}$, so basically a large vertical line.
Now if $\cos \theta = 0$,
$$\begin{align*} p(S) &= \left\{(x,y) \in \mathbb{R}^2\mid x^2 + \sin^2(\theta) y^2 = r^2 \right\}\\ &= \left\{(x,y) \in \mathbb{R}^2\mid x^2 + y^2 = r^2 \right\} \end{align*}$$
So it is a circle of radius r.
Let us get back to the study of the ellipse defined by $A(\theta) D_r$, to know it better it suffices to compute the polar decomposition of $A = OS$ with $O = O_2(\mathbb{R})$ and $s \in S_n^{++}(\mathbb{R})$. More particularly what is interesting to study is the eigenvalues of $S$ which means the square roots of the eigenvalues of $A^T A$, this will give you the distorsion of the circle to make your ellipse.
Let us compute $$A^T A = \frac{1}{\sin^4 \theta} \begin{pmatrix} 1 & -sin \theta \cos \theta \\ -sin \theta \cos \theta & \sin^{2} \theta \end{pmatrix} = \frac{1}{\sin^4 \theta }R$$ This will of course be diagonalizable and the characteristic polynomial of $R$ is $$\begin{align*} \chi_R(X) &= (X-1)(X-\sin^2\theta) - \sin^2\theta \cos^2 \theta \\ &= X^2 -(1 + \sin^2\theta)X +sin^2 \theta (1-\cos^2\theta) \\ &= X^2 -(1 + \sin^2\theta)X +sin^4 \theta \end{align*}$$
Now search for its roots, $$\begin{align*} \Delta &= (1+\sin^2 \theta)^2 - 4 \sin^4 \theta \\ &= 1 + 2 \sin^2 \theta-3\sin^4 \theta \\ &= (1-\sin^2\theta)(1 + 3 \sin^2\theta) \\ &= \cos^2\theta (1 + 3 \sin^2 \theta) \end{align*}$$
This gives you two distinct eigenvalues $$r_\pm = \frac{1 + \sin^2\theta \pm |\cos \theta| \sqrt{1 + 3 \sin^2 \theta}}{2}$$
So the eigenvalues of $S$ are $$s_\pm = \frac{1}{\sin^2\theta}\sqrt{\frac{1 + \sin^2\theta \pm |\cos \theta| \sqrt{1 + 3 \sin^2 \theta}}{2}}$$
Which means the distortion (length of big axis / length of small axis) will be $$\begin{align*} d &=\sqrt{\frac{1 + \sin^2\theta + |\cos \theta| \sqrt{1 + 3 \sin^2 \theta}}{1 + \sin^2\theta - |\cos \theta| \sqrt{1 + 3 \sin^2 \theta}}}\\ &=\sqrt{\frac{ \left( 1 + \sin^2\theta + |\cos \theta| \sqrt{1 + 3 \sin^2 \theta }\right)^2 }{r_+ r-}}\\ &=\frac{ 1 + \sin^2\theta + |\cos \theta| \sqrt{1 + 3 \sin^2 \theta } }{\sqrt{4\sin^4 \theta}}\\ &= \frac{ 1 + \sin^2\theta + |\cos \theta| \sqrt{1 + 3 \sin^2 \theta } }{2\sin^2 \theta} \end{align*} $$
So basically your ellipse will have $r_+$ for big axis and $r_-$ for small axis.
The computation of the angle is possible but very fastidious.