Skip to main content
 首页 » 编程设计

python之用skyfield寻找暮光之城

2025年05月04日14bjzhanghao

我有一个函数可以返回给定位置的日出、日落、正午和黄昏时间。此函数使用 pyephem由于它已被弃用,我想我会重新编写函数以使用 skyfield。但是,skyfield 没有.previous_rising、.next_setting 或.next_transit 函数(至少我找不到),这就是我与 pyephem 一起使用的。

skyfield 确实有函数 find_discrete ,它将在给定时间之间搜索以查找函数何时更改,因此我编写了以下内容来测试:

from datetime import datetime, timedelta 
from skyfield import api, almanac 
import pytz 
 
ts = api.load.timescale() 
planets = api.load('de421.bsp') 
sun, earth = planets['sun'], planets['earth'] 
 
lat, lon =  '44.59028 N', '104.71528 W'  # UFO Mooring Site 
tzn, elv = 'US/Mountain', 1559 
 
tz = pytz.timezone(tzn) 
loc = api.Topos(lat, lon, elevation_m=elv) 
earth_loc = earth + loc 
 
def civil_twil(time): 
    "returns true/false if sun is above -6 degrees" 
    alt, az, dis = earth_loc.at(time).observe(sun).apparent().altaz() 
    return alt > -6 
 
t0 = ts.utc(datetime.now(tz)) 
t1 = ts.utc(tz.normalize(datetime.now(tz) + timedelta(1))) 
civil_dawn = almanac.find_discrete(t0, t1, civil_twil) 

但这只是给我一个错误,即该函数缺少“rough_period”属性,并且文档没有提及可能是什么。我猜它可能只对 almanac 类中定义的函数有效;但同样,它没有被提及。

那么,如何通过天空场找到/确定黄昏时间?

请您参考如下方法:

有一个example在文档中使用 find_discrete 函数查找日出和日落的值。请参阅 almanac.py 中的来源对于这个功能。这里唯一的问题是它仅在太阳顶部明显与地平线齐平时才计算。在下面的示例中,我将 sunrise_sunset 函数更改为 daylength 函数。有了这个,您可以为 daylength 函数提供所需的角度。

from skyfield import api, almanac 
from datetime import datetime, timedelta 
import pytz 
from skyfield.nutationlib import iau2000b 
 
DAYLENGTH_CENTER_HORIZON = 0.0 
DAYLENGTH_TOP_HORIZON = 0.26667 
DAYLENGTH_TOP_HORIZON_APPARENTLY = 0.8333 
DAYLENGTH_CIVIL_TWILIGHT = 6.0 
DAYLENGTH_NAUTICAL_TWILIGHT = 12.0 
DAYLENGTH_ASTRONOMICAL_TWILIGHT = 18.0 
 
def daylength(ephemeris, topos, degrees): 
    """Build a function of time that returns the daylength. 
 
    The function that this returns will expect a single argument that is a  
    :class:`~skyfield.timelib.Time` and will return ``True`` if the sun is up 
    or twilight has started, else ``False``. 
    """ 
    sun = ephemeris['sun'] 
    topos_at = (ephemeris['earth'] + topos).at 
 
    def is_sun_up_at(t): 
        """Return `True` if the sun has risen by time `t`.""" 
        t._nutation_angles = iau2000b(t.tt) 
        return topos_at(t).observe(sun).apparent().altaz()[0].degrees > -degrees 
 
    is_sun_up_at.rough_period = 0.5  # twice a day 
    return is_sun_up_at 
 
ts = api.load.timescale() 
planets = api.load('de421.bsp') 
sun = planets['sun'] 
earth = planets['earth'] 
 
lat, lon = '44.59028 N', '104.71528 W'  # UFO Mooring Site 
tzn, elv = 'US/Mountain', 1559 
 
tz = pytz.timezone(tzn) 
loc = api.Topos(lat, lon, elevation_m=elv) 
 
t0 = ts.utc(datetime.now(tz)) 
t1 = ts.utc(tz.normalize(datetime.now(tz) + timedelta(1))) 
center_time, center_up = almanac.find_discrete(t0, t1, daylength(planets, loc, 
                                                    DAYLENGTH_CENTER_HORIZON)) 
print('Sunrise Sunset center of sun is even with horizon:') 
print(center_time.utc_iso(), center_up) 
 
apparent_top_time, apparent_top_up = almanac.find_discrete(t0, t1, 
        daylength(planets, loc, DAYLENGTH_TOP_HORIZON_APPARENTLY)) 
print('Sunrise Sunset top of sun is apparently even with horizon:') 
print(apparent_top_time.utc_iso(), apparent_top_up) 
 
 
civil_twil_time, civil_twil_up = almanac.find_discrete(t0, t1, 
        daylength(planets, loc, DAYLENGTH_CIVIL_TWILIGHT)) 
print('Civil twilight:') 
print(civil_twil_time.utc_iso(), civil_twil_up) 

这将打印以下结果:

Sunrise Sunset center of sun is even with horizon:
['2019-02-03T14:20:33Z', '2019-02-04T00:05:20Z'] [ True False]

Sunrise Sunset top of sun is apparently even with horizon:
['2019-02-03T14:15:28Z', '2019-02-04T00:10:25Z'] [ True False]

Civil twilight:
['2019-02-03T13:44:36Z', '2019-02-04T00:41:18Z'] [ True False]

第一个列表显示发现更改的时间,第二个列表显示 True 如果太阳升起(或暮光开始),False 当太阳升起时设置。

您的案例中缺少的 rough_period 应该是一个浮点值,表示一天发生的次数。太阳每天升起和落下一次,因此此函数中的事件每天发生两次。这意味着 rough_period 是 0.5。例如,当您想要计算月相时,可以将 rough_period 设置为 7.0(满月轨道为 27.3 天,每个月相为 6.825 天)。请参阅 almanac.py 中有关季节或月相计算的其他示例。源代码。