# Visualising orbits from different viewpoints in Python

I have a python program that I use to visualise a planetary orbital system. I convert the orbits to xyz coordinates and then have a basic orrery that views the system from 'above', by simply plotting the x and y to produce a kind of cross section, (the z is ignored for the plotting, but is accurately calculated from the elements)

I would like to be able to rotate the view so that, for example the system could be viewed from various angles. I expect the maths is really simple but I am no geometrist.

Any hints would be appreciated. I would like to stick with Python3 and numpy. I have tried googling but the results made my head hurt.

• "I would like to stick with Python 3 and NumPy" well, you will need a visualization library anyway, right? I recommend you to check out docs.poliastro.space/en/stable/examples/Plotting%20in%203D.html – astrojuanlu May 27 at 16:42
• I don't need a visualisation library - I already have a basic orrery coded because most visualisation libraries depend on some kind of heavy gui environment. I agree poliastro is brilliant (thanks for the link), as is matplotlib, I could also use pygame or just plot to a html canvas, however before I can do that I need to be able to do the projection. @tfb's answer is where I am going I think. – deMangler May 28 at 7:53

This answer doesn't talk about how to do it Python at all: rather how to deal with the rotation. I think once you can do that then turning the maths into Python is simple. Initially I'll assume that you are computing positions in terms of three orthogonal axes, and the positions look like $$(x, y, z)$$, and you're just projecting these down onto the $$(x, y)$$ plane (so, no perspective). So the question is given a set of rotated axes, $$(x', y', z')$$, how do you convert from one set to the other. Once you have $$x' = x'(x, y, z)$$ and similarly for the rest, you can compute the coordinates in the $$(x', y', z')$$ coordinates & then project down onto the $$(x', y')$$ plane.

The way to do this is to define rotation matrices: there are three angles you need to know. To see why, consider the new axes: there are two angles which define where the $$z'$$ axis is, and then you can spin the whole coordinate system around that axis, which is another angle.

And now I am going to make a mess of this, because I always get confused between rotation of vectors and rotation of coordinate systems: there are sign differences.

The rotation matrices are simply 3 by 3 matrices $$R$$ such that $$R^TR = RR^T = I$$ and $$\det{R} = 1$$, where $$R^T$$ means the transpose of $$R$$. These in fact are elements of a representation of the special orthogonal group in 3 dimensions, $$SO(3)$$: it's worth looking that up.

We can define these things as products of three basic rotations, about the $$x, y, z$$ axes respectively:

\begin{align} R_x(\theta) &= \begin{bmatrix} 1 & 0 & 0\\ 0 & \cos\theta & -\sin\theta\\ 0 & \sin\theta & \cos\theta \end{bmatrix}\\ R_y(\theta) &= \begin{bmatrix} \cos\theta & 0 & \sin\theta\\ 0 & 1 & 0\\ -\sin\theta & 0 & \cos\theta\\ \end{bmatrix}\\ R_z(\theta) &= \begin{bmatrix} \cos\theta & -\sin\theta & 0\\ \sin\theta & \cos\theta & 0\\ 0 & 0 & 1 \end{bmatrix} \end{align}

Then a general rotation about three angles $$\alpha, \beta, \gamma$$ is

$$R(\alpha,\beta,\gamma) = R_x(\alpha)R_y(\beta)R_z(\gamma)$$

Note that the multiplications here are of course matrix multiplication: the thing that in Numpy is np.matmul, and in particular they're not element-wise multiplication. Then, finally, to compute the new coordinates you do

$$[x',y',z'] = R(\alpha,\beta,\gamma)[x,y,z]^T$$

Where, again, this is matrix multiplication of course.

So here's an example: if we start off with a point at $$(x, 0, 0)$$, then what are its coordinates in a set of axes rotated by $$\theta$$ about $$z$$? Well, you can do the multiplications (or get your tame algebra system to do them for you) and the answer is $$(x', y', z') = (x\cos\theta, x\sin\theta,0)$$, and it's clear that the rotation of the axes is $$\theta$$, clockwise. I think this means I've botched a sign somewhere, but it doesn't really matter. And after a rotation about the $$x$$ axis you'll get $$(x',y',z') = (x, 0, 0)$$ which is obviously geometrically correct.

Combining rotations becomes less intuitive I think, but the maths just works.

I find by far the best way to work out what's going on is to write a program which implements the transformations, and then take some plots you understand and transform them, which will fairly quickly show you what the various angles mean and what mistakes you've made.

• Very helpful indeed. I can understand your answer very clearly. Thanks. :) – deMangler May 28 at 7:54

Start with answers to How can I plot a satellite's orbit in 3D from a TLE using Python and Skyfield?

Plotting in 3D makes my head hurt too, but for some reason I like it when my head hurts. If you like you can paste python into your question; blocks of text that are indented by 4 spaces appear as a "code block". You can have a look at the 3d plotting in answers to that linked question and if you have more specific questions go ahead and edit here. That script is old and probably can be written in a nicer way now. Also see this answer and this and this

And if you want to have even more fun with Python, download Blender and join Blender SE! See this and this. Blender's interface is Python and you can just paste big chunks of Python into it and animate, make movies, etc.

• Thanks - those are some good links. :) I am trying to keep my code as light as possible so I am trying to avoid too may libraries. I pretty much stick to Python3 and numpy at the moment. I will update my question. – deMangler May 27 at 9:57