Calculating right ascension for two planets and their separation using python

Question

I am experimenting with the Python library PyEphem for astronomy, and I am trying to recreate/understand some of the basic calculations there.

There is a function called separation that, given two planets and a date/time, it calculates the separation (angle) between those 2 planets with respect to their "x" projection in a plane.

So for 2018/1/1 for planets Mercury and Mars, we have:

import ephem
import math

mercury = ephem.Mercury('2018/1/1')
mars = ephem.Mars('2018/1/1')

s1 = ephem.separation(mercury, mars)
print(math.degrees(s1))

Which returns

33.792384499568264

But if I wanted to calculate this without the separation function, then I think that the calculation would be as simple as the "right ascension" of 1 minus the "right ascension" of the other:

math.degrees(mercury.ra) - math.degrees(mars.ra)

Which returns

35.114532008671574

Why are the angles different? Since I am not including a latitude and longitude of the observer, all the calculations are supposed to be geo-centric, according to PyEphem.

Is anybody familiar with the calculations going on behind PyEphem, or another library with built-in ephemeris that can produce consistent results for separation?

Answer 1

The issue appears to be in your understanding of the separation function. You state that

it calculates the separation (angle) between those 2 planets with respect to their "x" projection in a plane.

however, the documentation states

The separation() function computes the angle in degrees between two bodies, as measured by their right ascension and declination.

In general separation() can measure the angle between any pair of spherical coordinates, so long as the elements of each coordinate are spherical longitude (angle around the sphere) followed by spherical latitude (angle above or below its equator). Each pair should be provided as a two-item sequence like a tuple or list. Appropriate coordinate pairs include right ascension and declination; heliocentric longitude and latitude; azimuth and altitude; and even the geographic longitude and latitude of two locations on earth.

This means that this function is calculating the angular separation between the two objects, not just the separation in a particular plane. If you want to perform the same calculation (as the documentation says they're performing), it looks like you actually want the general angular distance formula.

$$\theta = \cos^{-1}\left[\sin(\delta_1)sin(\delta_2)+\cos(\delta_1)\cos(\delta_2)\cos(\alpha_1-\alpha_2)\right]$$

where $\delta\in[-\pi/2,\pi/2]$ is the declination, and $\alpha\in[0,2\pi]$ is the Right Ascension.