Solar System Orbits

This notebook allows you to plot the orbits of solar system objects using standard orbital elements provided by NASA/JPL ans associated organizations.  Since the orbital elements of planets and comets are provided in slightly different formats, two functions are defined, one for each case. For cometary bodies, the following are used:

e = eccentricity
q = perihelion distance
i = inclination angle
omega = argument of pericenter
nodes = longitude of the ascending node

For planets, the following are used:

e = eccentricity
a = length of semimajor axis in astronomical units
i = inclination angle
omega = argument of pericenter
nodes = longitude of the ascending node

The following package is useful for its 3D rotation matrix capabilities.

<< Geometry`Rotations`

Here, we define the two functions used to calculate orbits.

Orbit[{e_, q_, i_, omega_, nodes_}, opts___] := 
      Module[{a, b, atmp, btmp}, 
a = 
          q/(1 - e); 
b = Abs[
              btmp /. Solve[e == Sqrt[1 - btmp^2/a^2], 
                  btmp]][[1]]; 
ParametricPlot3D[
          Append[(RotationMatrix3D[-omega, -i, -nodes] . {a Cos[t] - (a - q), b Sin[t], 0}), PlotStyle /. opts], {t, 0, 2 π}, 
          PlotPoints -> 200]]

PlanetOrbit[{e_, a_, i_, omega_, nodes_}, opts___] := 
      Module[{q, b, atmp, 
          btmp}, 
b = Abs[
              btmp /. Solve[e == Sqrt[1 - btmp^2/a^2], 
                  btmp]][[1]]; 
q = 
          a (1 - e); 
ParametricPlot3D[
          Append[(RotationMatrix3D[-omega, -i, -nodes] . {a Cos[t] - (a - q), b Sin[t], 0}), PlotStyle /. opts], {t, 0, 2 π}, 
          PlotPoints -> 200]]

Here, we define the orbit of 2003 (Subscript[UB, 313]).  This object was recently announced to be 1.5 times the size of Pluto and is considered by some to be the 10th planet. Orbital data obtained from http://en.wikipedia.org/wiki/2003_UB313

ub313 = 
        Orbit[{ .4416129, 37.808, 44.177 Degree, 151.3115 Degree, 
            35.8750 Degree}, 
          PlotStyle -> {RGBColor[Clip[2 Cos[t/2]^2 +  .25], 0, 
                0]}];

Here, we define the orbit of comet Hale-Bopp.  Orbit data obtained from http://eiger.physics.uiowa.edu/~vis/hale-bopp-ephem.html

HaleBopp = 
        Orbit[{ .995103273, 0.914117707, 89.42965142 Degree, 
            130.5909555 Degree, 282.4708087 Degree}, 
          PlotStyle -> {Red}];

Here, we define the orbit of Sedna.  Currently, this is the most distant solar system object known since its most distant point is about 942 AUs from the sun.  It may be 3/4 the size of Pluto.  Sedna data obtained from http://en.wikipedia.org/wiki/90377_Sedna

sedna = 
        PlanetOrbit[{ .849, 502.040, 11.932 Degree, 311.468 Degree, 
            144.544 Degree}, 
          PlotStyle -> {RGBColor[Clip[ .7  Cos[t/2]^2 +  .1], 
                Clip[ .7 Cos[t/2]^2 +  .1], 
                Clip[1 Cos[t/2]^2 +  .1]]}];

Here, we define the orbits of the classic planets.

mercury = 
        PlanetOrbit[{0.20563069,  .38709893 , 7.00487 Degree, 
            48.33167 Degree, 252.25084 Degree}, 
          PlotStyle -> {Gray}];

venus = 
        PlanetOrbit[{ .00677323,  .72333199 , 3.39471 Degree, 
            76.68069 Degree, 181.97973 Degree}, 
          PlotStyle -> {White}];

earth = 
        PlanetOrbit[{ .01671022, 1.00000011 , 
            0.00005 Degree, -11.26064 Degree, 100.46435 Degree}, 
          PlotStyle -> {Blue}];

mars = 
        PlanetOrbit[{ .09341233, 1.52366231 , 1.85061 Degree, 
            49.57854 Degree, 355.45332 Degree}, 
          PlotStyle -> {Red}];

jupiter = 
        PlanetOrbit[{ .04839266, 5.20336301, 1.30530 Degree, 
            100.55615 Degree, 34.40438 Degree}, 
          PlotStyle -> {Orange}];

saturn = 
        PlanetOrbit[{ .05415060, 9.53707032 , 2.48446 Degree, 
            113.71504 Degree, 49.94432 Degree}, 
          PlotStyle -> {Yellow}];

uranus = 
        PlanetOrbit[{ .04716771, 19.19126393, 0.76986 Degree, 
            74.22988 Degree, 313.23218 Degree}, 
          PlotStyle -> {Green}];

neptune = 
        PlanetOrbit[{ .00858587, 30.06896348, 1.76917 Degree, 
            131.72169 Degree, 304.88003 Degree}, 
          PlotStyle -> {Purple}];

pluto = 
        PlanetOrbit[{ .24880766, 39.48168677 , 17.14175 Degree, 
            110.30347 Degree, 238.92881 Degree}, 
          PlotStyle -> {White}];

We will represent the sun as a Sphere and so we load the shapes package for this purpose.

<< Graphics`Shapes`

Now, we combine the orbits and use some options that allow us to control the PlotRange and zoom factor (via PlotRegion).

mydisplay[r_, zoom_] := 
      Show[ub313, mercury, venus, earth, mars, jupiter, saturn, uranus, 
        neptune, pluto, sedna, 
        Graphics3D[{SurfaceColor[Yellow, White, 3], EdgeForm[], 
            Sphere[0.0047, 50, 50]}], 
        PlotRange -> {{-r, r}, {-r, r}, {-r, r}}, 
        Background -> Black, ImageSize -> 600, 
        SphericalRegion -> True, 
        LightSources -> {{{0, 0, 1}, White}}, Boxed -> False, 
        Axes -> False, 
        PlotRegion -> {{-zoom, 1 + zoom}, {-zoom, 1 + zoom}}, 
        ViewPoint -> {1.3, -2.4,  .5}, 
        DisplayFunction -> Identity]

Here, we view the solar system from above out to a distance of about 100 AUs.

Show[mydisplay[1000, 1], DisplayFunction -> $DisplayFunction, 
      ViewPoint -> {0, 0, 2}]

Here, we circle the solar system at a distance of about 1000 AUs and we zoom in about 7 times.

anim = 
        Table[Show[{mydisplay[1000, 7], 
              Graphics3D[{Gray, Line[{{0, 0, 0}, {50, 0, 0}}], 
                  Line[{{0, 0, 0}, {0, 50, 0}}], 
                  Line[{{0, 0, 0}, {0, 0, 50}}]}]}, 
            ViewPoint -> {2 Cos[t], 2 Sin[t], - .5}, 
            DisplayFunction -> $DisplayFunction], {t, 0, 
            2 π,  .05}];

Here, we start at the center of the solar system and zoom out to about 250 AUs.

zoomanim = 
        Table[Show[mydisplay[250, z], 
            DisplayFunction -> $DisplayFunction], {z, 100, 
            1, - .5}];