pyephem恒星时给出了意想不到的结果

Question

我使用ephem首次，且无法理解oberver.sidereal_time的输出（）

我已经写了几个脚本，以确定从小时角度太阳时。第一个使用ephem来计算赤经和Meeus的天文算法得到一个公式，得到格林威治的平均恒星时间，它可以被转换为与经度一起的局部平均恒星时间。

import sys 
from datetime import datetime, time, timedelta 
import ephem 

def hour_angle(dt, longit, latit, elev): 
    obs = ephem.Observer() 
    obs.date = dt.strftime('%Y/%m/%d %H:%M:%S') 
    obs.lon = longit 
    obs.lat = latit 
    obs.elevation = elev 
    sun = ephem.Sun() 
    sun.compute(obs) 
    # get right ascention 
    ra = ephem.degrees(sun.g_ra) 

    # get sidereal time at greenwich (AA ch12) 
    jd = ephem.julian_date(dt) 
    t = (jd - 2451545.0)/36525 
    theta = 280.46061837 + 360.98564736629 * (jd - 2451545) \ 
      + .000387933 * t**2 - t**3/38710000 

    # hour angle (AA ch13) 
    ha = (theta + longit - ra * 180/ephem.pi) % 360 
    return ha 

def main(): 
    if len(sys.argv) != 6: 
     print 'Usage: hour_angle.py [YYYY/MM/DD] [HH:MM:SS] [longitude] [latitude] [elev]' 
     sys.exit() 
    else: 
     dt = datetime.strptime(sys.argv[1] + ' ' + sys.argv[2], '%Y/%m/%d %H:%M:%S') 
     longit = float(sys.argv[3]) 
     latit = float(sys.argv[4]) 
     elev = float(sys.argv[5]) 

    # get hour angle 
    ha = hour_angle(dt, longit, latit, elev) 

    # convert hour angle to timedelta from noon 
    days = ha/360 
    if days > 0.5: 
     days -= 0.5 
    td = timedelta(days=days) 

    # make solar time 
    solar_time = datetime.combine(dt.date(), time(12)) + td 
    print solar_time 

if __name__ == '__main__': 
    main() 

这使得输出，当我插上一些数据，我希望：

> python hour_angle_ephem.py 2012/11/16 20:34:56 -122.2697 37.8044 3.0 
2012-11-16 12:40:54.697115 

我写的第二个脚本计算赤经相同的方式，但使用ephem的sidereal_time（）获取本地明显恒星时间。

import sys 
from datetime import datetime, time, timedelta 
import math 
import ephem 

def solartime(observer, sun=ephem.Sun()): 
    sun.compute(observer) 
    # sidereal time == ra (right ascension) is the highest point (noon) 
    t = observer.sidereal_time() - sun.ra 
    return ephem.hours(t + ephem.hours('12:00')).norm # .norm for 0..24 

def main(): 
    if len(sys.argv) != 6: 
     print 'Usage: hour_angle.py [YYYY/MM/DD] [HH:MM:SS] [longitude] [latitude] [elev]' 
     sys.exit() 
    else: 
     dt = datetime.strptime(sys.argv[1] + ' ' + sys.argv[2], '%Y/%m/%d %H:%M:%S') 
     longit = float(sys.argv[3]) 
     latit = float(sys.argv[4]) 
     elev = float(sys.argv[5]) 

    obs = ephem.Observer() 
    obs.date = dt.strftime('%Y/%m/%d %H:%M:%S') 
    obs.lon = longit 
    obs.lat = latit 
    obs.elevation = elev 
    solar_time = solartime(obs) 
    print solar_time 

if __name__ == '__main__': 
    main() 

这不会让我输出我期望的。

python hour_angle_ephem2.py 2012/11/16 20:34:56 -122.2697 37.8044 3.0 
9:47:50.83 

AFAIK，两个脚本之间的唯一区别是，第一碱基时角上的本地平均恒星时，而第二个碱基时角上局部明显恒星时，其考虑了地球的章动，我认为这应该是一个非常小的因素。相反，我看到大约三个小时的差异。任何人都可以向我解释发生了什么事？

Answer 1

当您向PyEphem提供一个原始浮点数，并且期望一个角度时，它相信您已经将角度转换为弧度 - 因为它始终将浮点角度视为弧度，以保持一致。但在你的第二个剧本中，你得到的是以度为单位表示的经度和纬度，并将它们提供给PyEphem，就好像它们是弧度一样。你可以看到的结果，如果你添加一个print陈述或两个，看看你Observer看的.lon和.lat属性，如：

print observer.lon #--> -7005:32:16.0 
print observer.lat #--> 2166:01:57.2 

我认为，你不是想要做的仅仅是提供您的原始经度和纬度字符串转换为PyEphem，以便通过在第二个脚本中删除和argv[4]周围的float()调用，将它们解释为人类可读的度数而不是机器可读的弧度。然后你会发现，它返回一个值越接近你期待：

$ python tmp11.py 2012/11/16 20:34:56 -122.2697 37.8044 3.0 
12:40:55.59