Implement support for targets defined by TLE (Two-Line Element set)

astroplan/constraints.py

@@ -15,7 +15,7 @@
 import numpy as np
 from astropy import table
 from astropy.time import Time
-from astropy.coordinates import get_body, get_sun, Galactic, SkyCoord
+from astropy.coordinates import get_body, get_sun, Galactic
 from numpy.lib.stride_tricks import as_strided
 
 # Package
@@ -223,48 +223,45 @@ def __call__(self, observer, targets, times=None,
                  time_range=None, time_grid_resolution=0.5*u.hour,
                  grid_times_targets=False):
         """
-        Compute the constraint for this class
+        Compute the constraint for this class.
 
         Parameters
         ----------
         observer : `~astroplan.Observer`
-            the observation location from which to apply the constraints
-        targets : sequence of `~astroplan.Target`
+            The observation location from which to apply the constraints.
+
+        targets : sequence of `~astroplan.Target` or `~astropy.coordinates.SkyCoord`
             The targets on which to apply the constraints.
-        times : `~astropy.time.Time`
+        times : `~astropy.time.Time` (optional)
             The times to compute the constraint.
-            WHAT HAPPENS WHEN BOTH TIMES AND TIME_RANGE ARE SET?
-        time_range : `~astropy.time.Time` (length = 2)
+        time_range : `~astropy.time.Time` (length = 2) (optional)
             Lower and upper bounds on time sequence.
+            Only used when ``times`` is not provided.
         time_grid_resolution : `~astropy.units.Quantity`
-            Time-grid spacing
+            Time-grid spacing.
         grid_times_targets : bool
-            if True, grids the constraint result with targets along the first
+            If True, grids the constraint result with targets along the first
             index and times along the second. Otherwise, we rely on broadcasting
             the shapes together using standard numpy rules.
+
         Returns
         -------
         constraint_result : 1D or 2D array of float or bool
-            The constraints. If 2D with targets along the first index and times along
+            The constraint values. If 2D with targets along the first index and times along
             the second.
         """
-
         if times is None and time_range is not None:
             times = time_grid_from_range(time_range,
                                          time_resolution=time_grid_resolution)
 
-        if grid_times_targets:
-            targets = get_skycoord(targets)
-            # TODO: these broadcasting operations are relatively slow
-            # but there is potential for huge speedup if the end user
-            # disables gridding and re-shapes the coords themselves
-            # prior to evaluating multiple constraints.
-            if targets.isscalar:
-                # ensure we have a (1, 1) shape coord
-                targets = SkyCoord(np.tile(targets, 1))[:, np.newaxis]
-            else:
-                targets = targets[..., np.newaxis]
-        times, targets = observer._preprocess_inputs(times, targets, grid_times_targets=False)
+        # TODO: broadcasting operations in _preprocess_inputs are relatively slow
+        # but there is potential for huge speedup if the end user
+        # disables gridding and re-shapes the coords themselves
+        # prior to evaluating multiple constraints.
+        times, targets = observer._preprocess_inputs(
+            times, targets, grid_times_targets=grid_times_targets
+        )
+
         result = self.compute_constraint(times, observer, targets)
 
         # make sure the output has the same shape as would result from

astroplan/exceptions.py

@@ -40,3 +40,8 @@ class PlotBelowHorizonWarning(PlotWarning):
 class MissingConstraintWarning(AstroplanWarning):
     """Triggered when a constraint is expected but not supplied"""
     pass
+
+
+class InvalidTLEDataWarning(AstroplanWarning):
+    """TLE data invalid for the requested time"""
+    pass

astroplan/observer.py

@@ -507,22 +507,44 @@ def _preprocess_inputs(self, time, target=None, grid_times_targets=False):
             the shapes together using standard numpy rules. Useful for grid
             searches for rise/set times etc.
         """
-        # make sure we have a non-scalar time
         if not isinstance(time, Time):
             time = Time(time)
 
+        # In grid mode, scalar time should still behave like a length-1 time axis
+        if grid_times_targets and time.isscalar:
+            time = time[None]  # shape (1,)
+
         if target is None:
             return time, None
 
+        # Remember whether target is a single time-dependent target
+        is_multiple_targets = isinstance(target, (list, tuple))
+        is_target_time_dependent = (
+            callable(getattr(target, "get_skycoord", None)) and not
+            hasattr(target, "coord")
+        )
+        is_single_time_dependent_target = (not is_multiple_targets) and is_target_time_dependent
+
         # convert any kind of target argument to non-scalar SkyCoord
-        target = get_skycoord(target)
+        target = get_skycoord(target, times=time)
+
         if grid_times_targets:
+            # Only ambiguous case: a single time-dependent target produces shape == time.shape
+            # but grid mode requires a leading target axis (1, ...).
+            if (
+                is_single_time_dependent_target
+                and (not target.isscalar)
+                and (target.shape == time.shape)
+            ):
+                target = target[np.newaxis, ...]
+
+            # Ensure at least one targets axis for scalar targets
             if target.isscalar:
-                # ensure we have a (1, 1) shape coord
-                target = SkyCoord(np.tile(target, 1))[:, np.newaxis]
-            else:
-                while target.ndim <= time.ndim:
-                    target = target[:, np.newaxis]
+                target = SkyCoord(np.tile(target, 1))  # shape (1,)
+
+            # Make target have one more dim than time (first targets axis, then time axes)
+            while target.ndim < 1 + time.ndim:
+                target = target[..., np.newaxis]
 
         elif not self._is_broadcastable(target.shape, time.shape):
             raise ValueError('Time and Target arguments cannot be broadcast '

astroplan/scheduling.py

@@ -11,6 +11,7 @@
 from astropy import units as u
 from astropy.time import Time
 from astropy.table import Table
+from astropy.coordinates import ConvertError
 
 from .utils import time_grid_from_range, stride_array
 from .constraints import AltitudeConstraint
@@ -118,7 +119,7 @@ def __init__(self, blocks, observer, schedule, global_constraints=[]):
         self.observer = observer
         self.schedule = schedule
         self.global_constraints = global_constraints
-        self.targets = get_skycoord([block.target for block in self.blocks])
+        self.targets = [block.target for block in self.blocks]
 
     def create_score_array(self, time_resolution=1*u.minute):
         """
@@ -147,8 +148,9 @@ def create_score_array(self, time_resolution=1*u.minute):
                     applied_score = constraint(self.observer, block.target,
                                                times=times)
                     score_array[i] *= applied_score
+        targets = get_skycoord(self.targets, times=times)
         for constraint in self.global_constraints:
-            score_array *= constraint(self.observer, self.targets, times,
+            score_array *= constraint(self.observer, targets, times,
                                       grid_times_targets=True)
         return score_array
 
@@ -265,45 +267,108 @@ def open_slots(self):
         return [slot for slot in self.slots if not slot.occupied]
 
     def to_table(self, show_transitions=True, show_unused=False):
-        # TODO: allow different coordinate types
+        def _format_target_info(target):
+            if hasattr(target, "coord"):
+                try:
+                    return target.coord.icrs.to_string("hmsdms")
+                except ConvertError:
+                    return repr(target.coord)
+            if hasattr(target, "alt") and hasattr(target, "az"):
+                parts = [
+                    "alt={:.6f} deg".format(u.Quantity(target.alt).to_value(u.deg)),
+                    "az={:.6f} deg".format(u.Quantity(target.az).to_value(u.deg)),
+                ]
+                if hasattr(target, "pressure"):
+                    try:
+                        p_hpa = u.Quantity(target.pressure).to_value(u.hPa)
+                    except (TypeError, ValueError, u.UnitConversionError):
+                        pass
+                    else:
+                        if p_hpa != 0.0:
+                            parts.append("pressure={:.3f} hPa".format(p_hpa))
+                return ", ".join(parts)
+            if hasattr(target, "satellite"):
+                return (
+                    f"#{target.satellite.model.satnum} "
+                    f"epoch {target.satellite.epoch.utc_strftime(format='%Y-%m-%d %H:%M:%S')}"
+                )
+            return ""
+
         target_names = []
         start_times = []
         end_times = []
         durations = []
+        target_types = []
+        target_info = []
         ra = []
         dec = []
         config = []
+
         for slot in self.slots:
-            if hasattr(slot.block, 'target'):
+            if hasattr(slot.block, "target"):
                 start_times.append(slot.start.iso)
                 end_times.append(slot.end.iso)
                 durations.append(slot.duration.to(u.minute).value)
                 target_names.append(slot.block.target.name)
-                ra.append(u.Quantity(slot.block.target.ra))
-                dec.append(u.Quantity(slot.block.target.dec))
+                target_types.append(slot.block.target.__class__.__name__)
+                try:
+                    ra.append(u.Quantity(slot.block.target.ra))
+                except (AttributeError, NotImplementedError):
+                    ra.append("")
+                try:
+                    dec.append(u.Quantity(slot.block.target.dec))
+                except (AttributeError, NotImplementedError):
+                    dec.append("")
+                target_info.append(_format_target_info(slot.block.target))
                 config.append(slot.block.configuration)
             elif show_transitions and slot.block:
                 start_times.append(slot.start.iso)
                 end_times.append(slot.end.iso)
                 durations.append(slot.duration.to(u.minute).value)
-                target_names.append('TransitionBlock')
-                ra.append('')
-                dec.append('')
+                target_names.append("TransitionBlock")
+                target_types.append("TransitionBlock")
+                ra.append("")
+                dec.append("")
+                target_info.append("")
                 changes = list(slot.block.components.keys())
-                if 'slew_time' in changes:
-                    changes.remove('slew_time')
+                if "slew_time" in changes:
+                    changes.remove("slew_time")
                 config.append(changes)
             elif slot.block is None and show_unused:
                 start_times.append(slot.start.iso)
                 end_times.append(slot.end.iso)
                 durations.append(slot.duration.to(u.minute).value)
-                target_names.append('Unused Time')
-                ra.append('')
-                dec.append('')
-                config.append('')
-        return Table([target_names, start_times, end_times, durations, ra, dec, config],
-                     names=('target', 'start time (UTC)', 'end time (UTC)',
-                            'duration (minutes)', 'ra', 'dec', 'configuration'))
+                target_names.append("Unused Time")
+                target_types.append("")
+                ra.append("")
+                dec.append("")
+                target_info.append("")
+                config.append("")
+
+        return Table(
+            [
+                target_names,
+                start_times,
+                end_times,
+                durations,
+                ra,
+                dec,
+                target_types,
+                target_info,
+                config,
+            ],
+            names=(
+                "target",
+                "start time (UTC)",
+                "end time (UTC)",
+                "duration (minutes)",
+                "ra",
+                "dec",
+                "target type",
+                "target info",
+                "configuration",
+            ),
+        )
 
     def new_slots(self, slot_index, start_time, end_time):
         """
@@ -996,9 +1061,8 @@ def __call__(self, oldblock, newblock, start_time, observer):
             # use the constraints cache for now, but should move that machinery
             # to observer
             from .constraints import _get_altaz
-            from .target import get_skycoord
             if oldblock.target != newblock.target:
-                targets = get_skycoord([oldblock.target, newblock.target])
+                targets = get_skycoord([oldblock.target, newblock.target], times=start_time)
                 aaz = _get_altaz(start_time, observer, targets)['altaz']
                 sep = aaz[0].separation(aaz[1])
                 if sep/self.slew_rate > 1 * u.second: