Skip to content

Add cached Sun/Moon calculations in Observer class #578

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
jak574 wants to merge 1 commit into
base: main
Choose a base branch
from
+108 −55
Open

Add cached Sun/Moon calculations in Observer class #578

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
50 changes: 5 additions & 45 deletions astroplan/constraints.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -14,13 +14,14 @@
import astropy.units as u
import numpy as np
from astropy import table
from astropy.coordinates import Galactic, SkyCoord
from astropy.time import Time
from astropy.coordinates import get_body, get_sun, Galactic, SkyCoord
from numpy.lib.stride_tricks import as_strided

# Package
from .exceptions import MissingConstraintWarning
from .moon import moon_illumination
from .observer import _make_cache_key
from .target import get_skycoord
from .utils import time_grid_from_range

Expand All @@ -40,47 +41,6 @@
)


def _make_cache_key(times, targets):
"""
Make a unique key to reference this combination of ``times`` and ``targets``.

Often, we wish to store expensive calculations for a combination of
``targets`` and ``times`` in a cache on an ``observer``` object. This
routine will provide an appropriate, hashable, key to store these
calculations in a dictionary.

Parameters
----------
times : `~astropy.time.Time`
Array of times on which to test the constraint.
targets : `~astropy.coordinates.SkyCoord`
Target or list of targets.

Returns
-------
cache_key : tuple
A hashable tuple for use as a cache key
"""
# make a tuple from times
try:
timekey = tuple(times.jd) + times.shape
except BaseException: # must be scalar
timekey = (times.jd,)
# make hashable thing from targets coords
try:
if hasattr(targets, 'frame'):
# treat as a SkyCoord object. Accessing the longitude
# attribute of the frame data should be unique and is
# quicker than accessing the ra attribute.
targkey = tuple(targets.frame.data.lon.value.ravel()) + targets.shape
else:
# assume targets is a string.
targkey = (targets,)
except BaseException:
targkey = (targets.frame.data.lon,)
return timekey + targkey


def _get_altaz(times, observer, targets, force_zero_pressure=False):
"""
Calculate alt/az for ``target`` at times linearly spaced between
Expand Down Expand Up @@ -466,7 +426,7 @@ def _get_solar_altitudes(self, times, observer, targets):
observer.pressure = 0

# find solar altitude at these times
altaz = observer.altaz(times, get_sun(times))
altaz = observer.altaz(times, observer.get_body("sun", times))
altitude = altaz.alt
# cache the altitude
observer._altaz_cache[aakey] = dict(times=times,
Expand Down Expand Up @@ -542,7 +502,7 @@ def compute_constraint(self, times, observer, targets):
# centred frame, so the separation is as-seen
# by the observer.
# 'get_sun' returns ICRS coords.
sun = get_body('sun', times, location=observer.location)
sun = observer.get_body('sun', times)
targets = get_skycoord(targets)
solar_separation = sun.separation(targets)

Expand Down Expand Up @@ -583,7 +543,7 @@ def __init__(self, min=None, max=None, ephemeris=None):
self.ephemeris = ephemeris

def compute_constraint(self, times, observer, targets):
moon = get_body("moon", times, location=observer.location, ephemeris=self.ephemeris)
moon = observer.get_body("moon", times, ephemeris=self.ephemeris)
# note to future editors - the order matters here
# moon.separation(targets) is NOT the same as targets.separation(moon)
# the former calculates the separation in the frame of the moon coord
Expand Down
113 changes: 103 additions & 10 deletions astroplan/observer.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@
import astropy.units as u
import numpy as np
import pytz
from astropy.coordinates import (EarthLocation, SkyCoord, AltAz, get_sun,
from astropy.coordinates import (EarthLocation, SkyCoord, AltAz,
get_body, Angle, Longitude)
from astropy.time import Time
from astropy.utils.exceptions import AstropyDeprecationWarning
Expand All @@ -26,6 +26,47 @@
MAGIC_TIME = Time(-999, format='jd')


def _make_cache_key(times, targets):
"""
Make a unique key to reference this combination of ``times`` and ``targets``.

Often, we wish to store expensive calculations for a combination of
``targets`` and ``times`` in a cache on an ``observer``` object. This
routine will provide an appropriate, hashable, key to store these
calculations in a dictionary.

Parameters
----------
times : `~astropy.time.Time`
Array of times on which to test the constraint.
targets : `~astropy.coordinates.SkyCoord`
Target or list of targets.

Returns
-------
cache_key : tuple
A hashable tuple for use as a cache key
"""
# make a tuple from times
try:
timekey = tuple(times.ravel().jd) + times.shape
except BaseException: # must be scalar
timekey = (times.jd,)
# make hashable thing from targets coords
try:
if hasattr(targets, 'frame'):
# treat as a SkyCoord object. Accessing the longitude
# attribute of the frame data should be unique and is
# quicker than accessing the ra attribute.
targkey = tuple(targets.frame.data.lon.value.ravel()) + targets.shape
else:
# assume targets is a string.
targkey = (targets,)
except BaseException:
targkey = (targets.frame.data.lon,)
return timekey + targkey


# Handle deprecated MAGIC_TIME variable
def deprecation_wrap_module(mod, deprecated):
"""Return a wrapped object that warns about deprecated accesses"""
Expand Down Expand Up @@ -530,6 +571,58 @@ def _preprocess_inputs(self, time, target=None, grid_times_targets=False):
.format(time.shape, target.shape))
return time, target

def get_body(self, body, time, ephemeris=None):
"""
Get a solar system body from `~astropy.coordinates.get_body` for the
current observer location. Results are cached to speed up multiple
calls for the same body and time.

Parameters
----------
body : str
Name of solar system body.

time : `~astropy.time.Time` or other (see below)
The time at which the observation is taking place. Will be used as
the ``obstime`` attribute in the resulting frame or coordinate. This
will be passed in as the first argument to the `~astropy.time.Time`
initializer, so it can be anything that `~astropy.time.Time` will
accept (including a `~astropy.time.Time` object)

ephemeris : str, optional
Ephemeris to use. If not given, use the one set with
``astropy.coordinates.solar_system_ephemeris.set`` (which is
set to 'builtin' by default).

Returns
-------
`~astropy.coordinates.SkyCoord`
The position of the solar system body at the requested time.

Examples
--------
>>> from astroplan import Observer
>>> from astropy.time import Time
>>> from astropy.coordinates import SkyCoord
>>> apo = Observer.at_site("APO")
>>> time = Time('2001-02-03 04:05:06')
>>> target = apo.get_body("sun", time)
"""
if not isinstance(time, Time):
time = Time(time)

_cache_key = _make_cache_key(time, body)

if not hasattr(self, "_body_cache"):
self._body_cache = {}

if _cache_key not in self._body_cache:
self._body_cache[_cache_key] = get_body(
body, time, location=self.location, ephemeris=ephemeris
)

return self._body_cache[_cache_key]

def altaz(self, time, target=None, obswl=None, grid_times_targets=False):
"""
Get an `~astropy.coordinates.AltAz` frame or coordinate.
Expand Down Expand Up @@ -591,9 +684,9 @@ def altaz(self, time, target=None, obswl=None, grid_times_targets=False):
"""
if target is not None:
if target is MoonFlag:
target = get_body("moon", time, location=self.location)
target = self.get_body("moon", time)
elif target is SunFlag:
target = get_sun(time)
target = self.get_body("sun", time)

time, target = self._preprocess_inputs(time, target, grid_times_targets)

Expand Down Expand Up @@ -1346,7 +1439,7 @@ def sun_rise_time(self, time, which='nearest', horizon=0*u.degree, n_grid_points
>>> print("ISO: {0.iso}, JD: {0.jd}".format(sun_rise)) # doctest: +SKIP
ISO: 2001-02-02 14:02:50.554, JD: 2451943.08531
"""
return self.target_rise_time(time, get_sun(time), which, horizon,
return self.target_rise_time(time, self.get_body("sun", time), which, horizon,
n_grid_points=n_grid_points)

@u.quantity_input(horizon=u.deg)
Expand Down Expand Up @@ -1397,7 +1490,7 @@ def sun_set_time(self, time, which='nearest', horizon=0*u.degree, n_grid_points=
>>> print("ISO: {0.iso}, JD: {0.jd}".format(sun_set)) # doctest: +SKIP
ISO: 2001-02-04 00:35:42.102, JD: 2451944.52479
"""
return self.target_set_time(time, get_sun(time), which, horizon,
return self.target_set_time(time, self.get_body("sun", time), which, horizon,
n_grid_points=n_grid_points)

def noon(self, time, which='nearest', n_grid_points=150):
Expand Down Expand Up @@ -1426,7 +1519,7 @@ def noon(self, time, which='nearest', n_grid_points=150):
`~astropy.time.Time`
Time at solar noon
"""
return self.target_meridian_transit_time(time, get_sun(time), which,
return self.target_meridian_transit_time(time, self.get_body("sun", time), which,
n_grid_points=n_grid_points)

def midnight(self, time, which='nearest', n_grid_points=150):
Expand Down Expand Up @@ -1455,7 +1548,7 @@ def midnight(self, time, which='nearest', n_grid_points=150):
`~astropy.time.Time`
Time at solar midnight
"""
return self.target_meridian_antitransit_time(time, get_sun(time), which,
return self.target_meridian_antitransit_time(time, self.get_body("sun", time), which,
n_grid_points=n_grid_points)

# Twilight convenience functions
Expand Down Expand Up @@ -1812,7 +1905,7 @@ def moon_altaz(self, time, ephemeris=None):
if not isinstance(time, Time):
time = Time(time)

moon = get_body("moon", time, location=self.location, ephemeris=ephemeris)
moon = self.get_body("moon", time, ephemeris=ephemeris)
return self.altaz(time, moon)

def sun_altaz(self, time):
Expand Down Expand Up @@ -1841,7 +1934,7 @@ def sun_altaz(self, time):
if not isinstance(time, Time):
time = Time(time)

sun = get_sun(time)
sun = self.get_body("sun", time)
return self.altaz(time, sun)

@u.quantity_input(horizon=u.deg)
Expand Down Expand Up @@ -1951,7 +2044,7 @@ def is_night(self, time, horizon=0*u.deg, obswl=None):
if not isinstance(time, Time):
time = Time(time)

solar_altitude = self.altaz(time, target=get_sun(time), obswl=obswl).alt
solar_altitude = self.altaz(time, target=self.get_body("sun", time), obswl=obswl).alt

if solar_altitude.isscalar:
return bool(solar_altitude < horizon)
Expand Down