In the code below, computing the earth-sun distance with pyephem agrees well with the one compute from skyfield when using sun.compute(date), but not when doing sun.compute(obs) with obs an observer on earth at the same date.
Is this because in the "ephem2" case the sun.earth_distance is the distance from the sun to the observer,and not to the Earth center ?
from datetime import datetime
import ephem
from skyfield.api import load, utc
from numpy import deg2rad
date = datetime(2012,12,10,20,tzinfo=utc)
# skyfield
planets = load('de421.bsp')
earth,sun = planets['earth'],planets['sun']
ts = load.timescale()
t = ts.utc(date)
astrometric = earth.at(t).observe(sun)
ra, dec, distance = astrometric.radec()
d_skyfield = distance.m
# ephem
sun = ephem.Sun()
sun.compute(date)
d_ephem = sun.earth_distance*ephem.meters_per_au
# ephem2
obs = ephem.Observer()
obs.date = date
obs.lon = deg2rad(-97.4856)
obs.lat = deg2rad(36.604)
obs.elevation = 320
obs.pressure = 0
sun = ephem.Sun()
sun.compute(obs)
d_ephem2 = sun.earth_distance*ephem.meters_per_au
print 'Skyfield, Ephem',d_skyfield,d_ephem
print 'diff',(d_skyfield-d_ephem)/1000.0,'km'
print 'ratio',d_skyfield/(d_ephem)
print '\nSkyfield, Ephem2',d_skyfield,d_ephem2
print 'diff',(d_skyfield-d_ephem2)/1000.0,'km'
print 'ratio',d_skyfield/(d_ephem2)
print '\nEphem, Ephem2',d_ephem,d_ephem2
print 'diff',(d_ephem-d_ephem2)/1000.0,'km'
print 'ratio',d_ephem/(d_ephem2)
Output:
Skyfield, Ephem 147310866157.34006 1.47310854128e+11
diff 12.029628997802734 km
ratio 1.000000081661525
Skyfield, Ephem2 147310866157.34006 1.47308027525e+11
diff 2838.632373474121 km
ratio 1.000019270045368
Ephem, Ephem2 1.47310854128e+11 1.47308027525e+11
diff 2826.60274448 km
ratio 1.00001918838
Edit:
As suggested by Brandon, adding a topos object with skyfield gives the same result as sun.compute(obsever).
from skyfield.toposlib import Topos
# skyfield2
lamont = Topos(latitude_degrees=36.604,longitude_degrees=-97.4856,elevation_m=320)
lamont = earth + lamont
astrometric = lamont.at(t).observe(sun)
ra, dec, distance = astrometric.radec()
d_skyfield2 = distance.m
Output:
Skyfield, Skyfield2 147310866157.34006 147308040439.59784
diff 2825.7177422180175 km
ratio 1.0000191823727598
Skyfield2, Ephem2 147308040439.59784 1.47308027525e+11
diff 12.914631256103515 km
ratio 1.0000000876709265
The remaining difference of ~12 km varies with time/location
Yes, I would expect the difference to be because the measure is made to the observer's position; can you use a Skyfield "Topos" object to do the same measurement with Skyfield and see if they agree?