В приведенном ниже коде вычисление расстояния от Земли до Солнца с помощью пифема хорошо согласуется с вычислением, полученным по небу, при использовании sun.compute (дата), но не при выполнении sun.compute (obs) с наблюдателем на Земле в точкета же дата.
Это потому, что в случае "ephem2" sun.earth_distance - это расстояние от Солнца до наблюдателя, а не до центра Земли?
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)
Вывод:
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
Редактировать:
Как предложил Брэндон, добавление объекта топос с полем неба дает тот же результат, что и sun.compute (наблюдатель).
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
Выход:
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
Оставшаяся разница в ~ 12 км зависит от времени / местоположения