class CategoricalData:
    """Container for categorical (i.e. non-numerical) sensor data.

    This container allows simple manipulation and interpolation of a time series
    of non-numerical data represented as discrete events. The data is stored as
    a list of sensor values and two integer arrays:

    * `unique_values` stores one copy of each unique object in the data series
    * `events` stores the time indices (dumps) where each event occurs
    * `indices` stores indices linking each event to the `unique_values` list

    The __getitem__ interface (i.e. `data[dump]`) returns the data associated
    with the last event before the requested dump(s), in effect doing a
    zeroth-order interpolation of the data at each event. Events can be added
    and removed and realigned, and the container can be split along the time
    axis, amongst other functionality.

    Parameters
    ----------
    sensor_values : sequence, length *N*
        Sequence of sensor values (of any type, preferably not None [see Notes])
    events : sequence of non-negative ints, length *N* + 1
        Corresponding monotonic sequence of dump indices where each sensor value
        came into effect. The last event is one past the last dump where the
        final sensor value applied, and therefore equal to the total number of
        dumps for which sensor values were specified.

    Attributes
    ----------
    unique_values : list, length *M*
        List of unique sensor values in order they were found in `sensor_values`
        with any :class:`ComparableArrayWrapper` objects unwrapped
    indices : array of int, shape (*N*,)
        Array of indices into `unique_values`, one per sensor event
    dtype : :class:`numpy.dtype` object
        Sensor data type as NumPy dtype (found on demand from `unique_values`)

    Notes
    -----
    Any object values wrapped in a :class:`ComparableArrayWrapper` will be
    unwrapped before adding it to `unique_values`. When adding, removing and
    comparing values to this container, any object values will be wrapped again
    temporarily to ensure proper comparisons.

    It is discouraged to have a sensor value of None as this value is given
    a special meaning in methods such as :meth:`CategoricalData.add` and
    :meth:`sensor_to_categorical`. On the other hand, it is the most sensible
    dummy object value and any Nones entering through this initialiser will
    probably not cause any issues.

    It is better to make `unique_values` a list instead of an array because an
    array assimilates objects such as tuples, lists and other arrays. The
    alternative is an array of :class:`ComparableArrayWrapper` objects but
    these then need to be unpacked at some later stage which is also tricky.

    """

    def __init__(self, sensor_values, events):
        values, self.indices = unique_in_order(sensor_values, return_inverse=True)
        self.unique_values = [ComparableArrayWrapper.unwrap(v) for v in values]
        self.events = np.asarray(events)

    @property
    def _comparable_values(self):
        """Comparable version of unique values, wrapping any objects."""
        return [ComparableArrayWrapper(value) for value in self.unique_values]

    def _lookup(self, dumps):
        """Look up relevant indices occurring at specified `dumps`.

        Parameters
        ----------
        dumps : int or sequence of int
            Specified dumps

        Returns
        -------
        indices : int or array of int
            Corresponding sensor value indices at specified dumps

        """
        preceding_events = self.events.searchsorted(dumps, side='right') - 1
        if np.any(preceding_events < 0) or np.any(preceding_events >= len(self.indices)):
            raise IndexError('Some dumps in (%s) are outside event range: %d <= dumps < %d' %
                             (dumps, self.events[0], self.events[-1]))
        return self.indices[preceding_events]

    def __getitem__(self, key):
        """Look up sensor value at selected dumps.

        Be aware that the dtype of the returned array may differ from that of
        the :class:`CategoricalData` object if the sensor values are
        multi-dimensional arrays themselves, in effect falling back to the
        underlying dtype. For large multi-dimensional sensor values this
        method may also cause memory issues as it will duplicate these arrays
        into the final output array.

        Parameters
        ----------
        key : int or slice or sequence of int or sequence of bool
            Index or slice specification selecting certain dumps

        Returns
        -------
        val : object or array of objects
            Sensor values at selected dumps, either single value or array of them

        """
        if isinstance(key, slice):
            # Convert slice notation to the corresponding sequence of dump indices
            key = list(range(*key.indices(self.events[-1])))
        # Convert sequence of bools (one per dump) to sequence of indices where key is True
        elif np.asarray(key).dtype == np.bool and len(np.asarray(key)) == self.events[-1]:
            key = np.nonzero(key)[0]
        indices = self._lookup(key)
        # Interpret indices as either a sequence of ints or a single int
        try:
            values = [self.unique_values[index] for index in indices]
        except TypeError:
            return self.unique_values[indices]
        # Handle empty selections specially to ensure proper dtype and shape
        if not values:
            all_possible_values = np.array(self.unique_values)
            dtype = all_possible_values.dtype
            shape = all_possible_values.shape
            return np.empty((0,) + shape[1:], dtype)
        return np.array(values)

    def __repr__(self):
        """Short human-friendly string representation of categorical data object."""
        return "<katdal.CategoricalData events={} values={} type={} at {:#x}>".format(
               len(self.indices), len(self.unique_values), self.dtype, id(self))

    def __str__(self):
        """Long human-friendly string representation of categorical data object."""
        index_width = len(str(self.events[-1] - 1))
        return '\n'.join([('%*d - %*d: %s' % (index_width, segm.start, index_width, segm.stop - 1, val))
                          for segm, val in self.segments()])

    def __len__(self):
        """Length operator indicates number of events produced by sensor."""
        return len(self.indices)

    def _bool_per_dump(self, bool_per_value):
        """Turn list of bools per unique value into an array of bools per dump."""
        bool_per_event = np.atleast_1d(np.array(bool_per_value)[self.indices])
        bool_per_dump = np.empty(self.events[-1], dtype=np.bool)
        for n, (start, end) in enumerate(zip(self.events[:-1], self.events[1:])):
            bool_per_dump[start:end] = bool_per_event[n]
        return bool_per_dump

    def __eq__(self, other):
        """Equality comparison operator."""
        return self._bool_per_dump([value == other for value in self._comparable_values])

    def __ne__(self, other):
        """Inequality comparison operator."""
        return self._bool_per_dump([value != other for value in self._comparable_values])

    def __lt__(self, other):
        """Less-than comparison operator."""
        return self._bool_per_dump([value < other for value in self._comparable_values])

    def __gt__(self, other):
        """Greather-than comparison operator."""
        return self._bool_per_dump([value > other for value in self._comparable_values])

    def __le__(self, other):
        """Less-than-or-equal comparison operator."""
        return self._bool_per_dump([value <= other for value in self._comparable_values])

    def __ge__(self, other):
        """Greater-than-or-equal comparison operator."""
        return self._bool_per_dump([value >= other for value in self._comparable_values])

    @property
    def dtype(self):
        """Sensor value type."""
        return infer_dtype(self.unique_values)

    def segments(self):
        """Generator that iterates through events and returns segment and value.

        Yields
        ------
        segment : slice object
            The slice representing range of dump indices of the current segment
        value : object
            Sensor value associated with segment

        """
        for start, end, ind in zip(self.events[:-1], self.events[1:], self.indices):
            yield slice(start, end), self.unique_values[ind]

    def add(self, event, value=None):
        """Add or override sensor event.

        This adds a new event to the container, with a new value or a duplicate
        of the existing value at that dump. If the new event coincides with an
        existing one, it overrides the value at that dump.

        Parameters
        ----------
        event : int
            Dump of event to add or override
        value : object, optional
            New value for event (duplicate current value at this dump by default)

        """
        # If value has not been seen before, add it to unique_values (and create new index for it)
        if value is not None:
            try:
                value_index = self._comparable_values.index(value)
            except ValueError:
                value_index = len(self.unique_values)
                self.unique_values += [value]
        else:
            value_index = self._lookup(event)
        # If new event coincides with existing event, simply change value of that event, else insert new event
        event_index = self.events.searchsorted(event)
        before, after = event_index, (event_index + 1 if self.events[event_index] == event else event_index)
        self.indices = np.r_[self.indices[:before], [value_index], self.indices[after:]]
        self.events = np.r_[self.events[:before], [event], self.events[after:]]

    def remove(self, value):
        """Remove sensor value, remapping indices and merging segments in process.

        If the sensor value does not exist, do nothing.

        Parameters
        ----------
        value : object
            Sensor value to remove from container

        """
        try:
            index = self._comparable_values.index(value)
        except ValueError:
            pass
        else:
            keep = (self.indices != index)
            remap = np.arange(len(self.unique_values))
            remap[index:] -= 1
            self.indices = remap[self.indices[keep]]
            self.events = np.r_[self.events[:-1][keep], self.events[-1]]
            del self.unique_values[index]

    def add_unmatched(self, segments, match_dist=1):
        """Add duplicate events for segment starts that don't match sensor events.

        Given a sequence of segments, this matches each segment start to the
        nearest sensor event dump (within `match_dist`). Any unmatched segment
        starts are added as duplicate sensor events (or ignored if they fall
        outside the sensor event range).

        Parameters
        ----------
        segments : sequence of int
            Monotonically increasing sequence of segment starts, including an
            extra element at the end that is one past the end of the last segment
        match_dist : int, optional
            Maximum distance in dumps that signify a match between events

        """
        # Identify unmatched segment starts
        segments = np.asarray(segments)
        unmatched = segments[np.abs(self.events[np.newaxis, :] - segments[:, np.newaxis]).min(axis=1) > match_dist]
        # Add these dumps as duplicate events, ignoring those that are out of bounds
        for segm in unmatched:
            try:
                self.add(segm)
            except IndexError:
                pass

    def align(self, segments):
        """Align sensor events with segment starts, possibly discarding events.

        Given a sequence of segments, this moves each sensor event dump onto the
        nearest segment start. If more than one event ends up in the same segment,
        only keep the last event, discarding the rest.

        The end result is that the sensor event dumps become a subset of the
        segment starts and there cannot be more sensor events than segments.

        Parameters
        ----------
        segments : sequence of int
            Monotonically increasing sequence of segment starts, including an
            extra element at the end that is one past the end of the last segment

        """
        # For each event, pick the segment with the closest start to it and then shift event onto segment start
        segments_with_event = np.abs(self.events[np.newaxis, :] - segments[:, np.newaxis]).argmin(axis=0)
        events = segments[segments_with_event]
        # When multiple sensor events are associated with the same segment, only keep the final one
        final = np.nonzero(np.diff(events) > 0)[0]
        subset, self.indices = np.unique(self.indices[final], return_inverse=True)
        self.unique_values = [self.unique_values[index] for index in subset]
        self.events = np.r_[events[final], events[-1]]

    def partition(self, segments):
        """Partition dataset into multiple sets along time axis.

        Given a sequence of segments, split the container into a sequence of
        containers, one per segment. Each container contains only the events
        occurring within its corresponding segment, with event dumps relative to
        the start of the segment, and the containers share the same unique
        values.

        Parameters
        ----------
        segments : sequence of int
            Monotonically increasing sequence of segment starts, including an
            extra element at the end that is one past the end of the last segment

        Returns
        -------
        split_data : sequence of :class:`CategoricalData` objects
            Resulting multiple datasets in chronological order

        """
        # Ignore last element in event list, as it is not a real event but a placeholder for dataset length
        events = self.events[:-1]
        # Find segment starts in event sequence, associating dumps before first event with it, ditto for ones past last
        initial_indices = self.indices[(events.searchsorted(segments[:-1], side='right') - 1).clip(0, len(events) - 1)]
        split_data = []
        for start, end, initial_index in zip(segments[:-1], segments[1:], initial_indices):
            segment_events = (events >= start) & (events < end)
            # Bypass the normal CategoricalData initialiser to ensure that each cat_data has the same unique_values
            cat_data = CategoricalData([], [])
            cat_data.unique_values = self.unique_values
            cat_data.indices = self.indices[segment_events]
            cat_data.events = events[segment_events] - start
            # Insert initial event if it is not there, and pad events with data segment length
            if len(cat_data.events) == 0 or cat_data.events[0] != 0:
                cat_data.indices = np.r_[initial_index, cat_data.indices]
                cat_data.events = np.r_[0, cat_data.events, end - start]
            else:
                cat_data.events = np.r_[cat_data.events, end - start]
            split_data.append(cat_data)
        return split_data

    def remove_repeats(self):
        """Remove repeated events of the same value."""
        changes = np.nonzero([1] + np.diff(self.indices).tolist())[0]
        self.indices = self.indices[changes]
        self.events = np.r_[self.events[changes], self.events[-1]] def sensor_to_categorical(sensor_timestamps, sensor_values, dump_midtimes,
                          dump_period, transform=None, initial_value=None,
                          greedy_values=None, allow_repeats=False, **kwargs):
    """Align categorical sensor events with dumps and clean up spurious events.

    This converts timestamped sensor data into a categorical dataset by
    comparing the sensor timestamps to a series of dump timestamps and assigning
    each sensor event to the dump in which it occurred. When multiple sensor
    events happen in the same dump, only the last one is kept. The first dump is
    guaranteed to have a valid value by either using the supplied `initial_value`
    or extrapolating the first proper value back in time. The sensor data may
    be transformed before events that repeat values are potentially discarded.
    Finally, events with values marked as "greedy" take precedence over normal
    events when both occur within the same dump (either changing from or to the
    greedy value, or if the greedy value occurs completely within a dump).

    XXX Future improvements include picking the event with the longest duration
    within a dump as opposed to the final event, and "snapping" event boundaries
    to dump boundaries with a given tolerance (e.g. 5-10% of dump period).

    Parameters
    ----------
    sensor_timestamps : sequence of float, length *M*
        Sequence of sensor timestamps (typically UTC seconds since Unix epoch)
    sensor_values : sequence, length *M*
        Corresponding sequence of sensor values [potentially wrapped]
    dump_midtimes : sequence of float, length *N*
        Sequence of dump midtimes (same reference as sensor timestamps)
    dump_period : float
        Duration of each dump, in seconds
    transform : callable or None, optional
        Transform [unwrapped] sensor values before fixing initial value,
        mapping dumps to events and discarding repeats
    initial_value : object or None, optional
        Sensor value [transformed, unwrapped] to use for dump = 0 up to first
        proper event (force first proper event to start at dump = 0 by default)
    greedy_values : sequence or None, optional
        List of [transformed, unwrapped] sensor values considered "greedy"
    allow_repeats : {False, True}, optional
        If False, discard sensor events that do not change [transformed] value

    Returns
    -------
    data : :class:`CategoricalData` object
        Constructed categorical dataset [unwraps any wrapped values]

    """
    sensor_timestamps = np.atleast_1d(sensor_timestamps)
    sensor_values = np.atleast_1d(sensor_values)
    dump_endtimes = dump_midtimes + 0.5 * dump_period
    num_dumps = len(dump_endtimes)
    # Insert an extra prior dump to collect sensor values before first dump
    dump_endtimes = np.r_[dump_endtimes[0] - dump_period, dump_endtimes]
    # Check if sensor values are objects wrapped in ComparableArrayWrappers
    wrapped_values = len(sensor_values) and isinstance(sensor_values[0],
                                                       ComparableArrayWrapper)
    # Convert sensor event times to dump indices by picking the dump during
    # which each sensor event occurred (more precisely: closest dump centroid)
    events = dump_endtimes.searchsorted(sensor_timestamps) - 1
    # Get rid of excess events before the first dump and after the last dump
    # The dump index of prior events is -1 and of later events is `num_dumps`
    first_proper_event = events.searchsorted(-1, side='right')
    # Shift the final prior event (if any) to the start of the first dump
    if first_proper_event > 0:
        first_proper_event -= 1
        events[first_proper_event] = 0
    one_past_last_event = events.searchsorted(num_dumps)
    within_dumps = slice(first_proper_event, one_past_last_event)
    sensor_values = sensor_values[within_dumps]
    events = events[within_dumps]
    # Apply optional transform to sensor values
    if transform is not None:
        if wrapped_values:
            orig_transform = transform
            def transform(value):   # noqa: E306
                """Unwrap wrapped value, transform and rewrap."""
                return ComparableArrayWrapper(orig_transform(value.unwrapped))
        sensor_values = np.array([transform(y) for y in sensor_values])
    # Force first dump to have valid sensor value
    # (insert initial value or let the first proper value apply from the start)
    if events[0] != 0 and initial_value is not None:
        if wrapped_values:
            initial_value = ComparableArrayWrapper(initial_value)
        sensor_values = np.r_[[initial_value], sensor_values]
        events = np.r_[0, events]
    events[0] = 0
    # Clean up dump->event mapping, taking into account greedy values
    greedy_values = () if greedy_values is None else greedy_values
    greedy = [value in greedy_values for value in sensor_values]
    # Add one-past-last-dump terminator (will be removed again by `cleaned_up`)
    events = np.r_[events, num_dumps]
    # NB: `events` is mutated by `_single_event_per_dump`
    cleaned_up = list(_single_event_per_dump(events, greedy))
    sensor_values = sensor_values[cleaned_up]
    events = events[cleaned_up]
    # Discard sensor events that do not change the (transformed) sensor value
    # (i.e. that repeat the previous value)
    if not allow_repeats:
        changes_value = [n for n in range(len(sensor_values)) if n == 0 or
                         sensor_values[n] != sensor_values[n - 1]]
        sensor_values = sensor_values[changes_value]
        events = events[changes_value]
    # Last event is fixed at one-past-last-dump to indicate end of last segment
    return CategoricalData(sensor_values, np.r_[events, num_dumps]) DEFAULT_SENSOR_PROPS = {
    '*nd_coupler': {'categorical': True, 'greedy_values': (True,), 'initial_value': '0',
                    'transform': lambda x: x not in ('0', 'False', 0)},
    '*nd_pin': {'categorical': True, 'greedy_values': (True,), 'initial_value': '0',
                'transform': lambda x: x not in ('0', 'False', 0)},
    'Observation/label': {'initial_value': '', 'transform': str, 'allow_repeats': True},
    'Observation/scan_state': {'allow_repeats': True},
} DEFAULT_VIRTUAL_SENSORS = {
    'Timestamps/mjd': _calc_mjd, 'Antennas/{ant}/lst': _calc_lst,
    'Antennas/{ant}/ra': _calc_radec, 'Antennas/{ant}/dec': _calc_radec,
    'Antennas/{ant}/parangle': _calc_parangle,
    'Antennas/{ant}/target_[xy]_{projection}_{coordsys}': _calc_target_coords,
    'Antennas/{ant}/basis_[uvw]': _calc_uvw_basis,
    'Antennas/{ant}/[uvw]': _calc_uvw_per_ant,
} class BrokenFile(Exception):
    """Data set could not be loaded because file is inconsistent or misses critical bits.""" class DataSet:
    """Base class for accessing a visibility data set.

    This provides a simple interface to a generic file (or files) containing
    visibility data (both single-dish and interferometer data supported).
    The data are not loaded into memory on opening the file, but are accessible
    via properties after typically selecting a subset of the data. This allows
    the reading of huge files.

    Parameters
    ----------
    name : string
        Name / identifier of data set
    ref_ant : string, optional
        Name of reference antenna, used to partition data set into scans
        (default is first antenna in use by script)
    time_offset : float, optional
        Offset to add to all correlator timestamps, in seconds
    url : string, optional
        Location of data set (either local filename or full URL accepted)

    Attributes
    ----------
    version : string
        Format version string
    observer : string
        Name of person that recorded the data set
    description : string
        Short description of the purpose of the data set
    experiment_id : string
        Experiment ID, a unique string used to link the data files of an
        experiment together with blog entries, etc.
    obs_params : dict mapping string to string or list of strings
        Observation parameters, typically set in observation script
    obs_script_log : list of strings
        Observation script output log (useful for debugging)

    subarrays : list of :class:`SubArray` objects
        List of all subarrays in data set
    subarray : int
        Index of currently selected subarray
    ants : list of :class:`katpoint.Antenna` objects
        List of selected antennas
    inputs : array of strings
        List of selected correlator input labels ('ant1h')
    corr_products : array of strings, shape (*B*, 2)
        Array of selected correlation products as pairs of input labels
        (e.g. [('ant1h', 'ant1h'), ('ant1h', 'ant2h')])
    receivers : dict mapping string to string or list of strings
        Identifier of the active receiver on each antenna

    spectral_windows : list of :class:`SpectralWindow` objects
        List of all spectral windows in data set
    spw : int
        Index of currently selected spectral window
    channel_width : float
        Channel bandwidth of selected spectral window, in Hz
    freqs / channel_freqs : array of float, shape (*F*,)
        Centre frequency of each selected channel, in Hz
    channels : array of int, shape (*F*,)
        Original channel indices of selected channels

    dump_period : float
        Dump period, in seconds
    sensor : :class:`SensorCache` object
        Sensor cache
    catalogue : :class:`katpoint.Catalogue` object
        Catalogue of all targets / sources / fields in data set
    start_time : :class:`katpoint.Timestamp` object
        Timestamp of start of first sample in file, in UT seconds since Unix epoch
    end_time : :class:`katpoint.Timestamp` object
        Timestamp of end of last sample in file, in UT seconds since Unix epoch
    dumps : array of int, shape (*T*,)
        Original dump indices of selected dumps
    scan_indices : list of int
        List of currently selected scans as indices
    compscan_indices : list of int
        List of currently selected compound scans as indices
    target_indices : list of int
        List of currently selected targets as indices into catalogue
    target_projection : {'ARC', 'SIN', 'TAN', 'STG', 'CAR'}, optional
        Type of spherical projection for target coordinates
    target_coordsys : {'azel', 'radec'}, optional
        Spherical pointing coordinate system for target coordinates
    shape : tuple of 3 ints
        Shape of selected visibility data array, as (*T*, *F*, *B*)
    size : int
        Size of selected visibility data array, in bytes
    applycal_products : list of string
        List of calibration products that will be applied to data

    """

    def __init__(self, name, ref_ant='', time_offset=0.0, url=''):
        self.name = name
        self.ref_ant = ref_ant
        self.time_offset = time_offset
        self.url = parse_url_or_path(url).geturl()

        self.version = ''
        self.observer = ''
        self.description = ''
        self.experiment_id = ''
        self.obs_params = {}
        self.obs_script_log = []

        self.subarrays = []
        self.subarray = -1
        self.ants = []
        self.inputs = []
        self.corr_products = np.empty(shape=(0, 2), dtype='|S5')
        self.receivers = {}

        self.spectral_windows = []
        self.spw = -1
        self.channel_width = 0.0
        self.freqs = self.channel_freqs = np.empty(0)
        self.channels = np.empty(0, dtype=np.int)

        self.dump_period = 0.0
        self.sensor = {}
        self.catalogue = katpoint.Catalogue()
        self.start_time = katpoint.Timestamp(0.0)
        self.end_time = katpoint.Timestamp(0.0)
        self.dumps = np.empty(0, dtype=np.int)
        self.scan_indices = []
        self.compscan_indices = []
        self.target_indices = []
        self.target_projection = 'ARC'
        self.target_coordsys = 'azel'
        self.shape = (0, 0, 0)
        self.size = 0
        self.applycal_products = []

        self._selection = {}
        self._time_keep = []
        self._freq_keep = []
        self._corrprod_keep = []
        self._weights_keep = 'all'
        self._flags_keep = 'all'

    def __repr__(self):
        """Short human-friendly string representation of data set object."""
        class_name = self.__class__.__name__
        return f"<katdal.{class_name} '{self.name}' shape {self.shape} at {id(self):#x}>"

    def __str__(self):
        """Verbose human-friendly string representation of data set."""
        # Start with static file information
        descr = ['===============================================================================',
                 f'Name: {self.name} (version {self.version})',
                 '-------------------------------------------------------------------------------',
                 f'URL: {self.url}',
                 '===============================================================================',
                 'Observer: {}  Experiment ID: {}'.format(self.observer if self.observer else '-',
                                                          self.experiment_id if self.experiment_id else '-'),
                 "Description: '{}'".format(self.description if self.description else 'No description'),
                 f'Observed from {self.start_time.local()} to {self.end_time.local()}',
                 'Dump rate / period: {:.5f} Hz / {:.3f} s'.format(1 / self.dump_period, self.dump_period),
                 'Subarrays: {}'.format(len(self.subarrays)),
                 '  ID  Antennas                            Inputs  Corrprods']
        for n, sub in enumerate(self.subarrays):
            ant_names = ','.join([ant.name for ant in sub.ants])
            descr.append('  %2d  %28s  %2d      %3d' %
                         (n, ant_names.ljust(7 * 4 + 6), len(sub.inputs), len(sub.corr_products)))
        descr += ['Spectral Windows: %d' % (len(self.spectral_windows),),
                  '  ID Band Product  CentreFreq(MHz)  Bandwidth(MHz)  Channels  ChannelWidth(kHz)']
        for n, spw in enumerate(self.spectral_windows):
            descr.append('  %2d %-4s %-9s %9.3f        %8.3f          %5d     %9.3f' %
                         (n, spw.band if spw.band else '-', spw.product if spw.product else '-',
                          spw.centre_freq / 1e6,
                          spw.channel_width / 1e6 * spw.num_chans,
                          spw.num_chans, spw.channel_width / 1e3))
        # Now add dynamic information, which depends on the current selection criteria
        descr += ['-------------------------------------------------------------------------------',
                  'Data selected according to the following criteria:']
        for k, v in sorted(self._selection.items()):
            descr.append('  {}={}'.format(k, f"'{v}'" if isinstance(v, str) else v))
        descr.append('-------------------------------------------------------------------------------')
        descr.append('Shape: (%d dumps, %d channels, %d correlation products) => Size: %s' %
                     tuple(list(self.shape) + ['%.3f %s' % ((self.size / 1e9, 'GB') if self.size > 1e9 else
                                                            (self.size / 1e6, 'MB') if self.size > 1e6 else
                                                            (self.size / 1e3, 'KB'))]))
        autocorrs = np.array([(inpA[:-1] == inpB[:-1]) for inpA, inpB in self.corr_products])
        crosscorrs = np.array([(inpA[:-1] != inpB[:-1]) for inpA, inpB in self.corr_products])
        descr.append('Antennas: %s  Inputs: %d  Autocorr: %s  Crosscorr: %s' %
                     (','.join([('*' + ant.name if ant.name == self.ref_ant else ant.name) for ant in self.ants]),
                      len(self.inputs), 'yes' if np.any(autocorrs) else 'no', 'yes' if np.any(crosscorrs) else 'no'))
        chan_min, chan_max = self.channels.argmin(), self.channels.argmax()
        descr.append('Channels: %d (index %d - %d, %8.3f MHz - %8.3f MHz), each %7.3f kHz wide' %
                     (len(self.channels), self.channels[chan_min], self.channels[chan_max],
                      self.freqs[chan_min] / 1e6, self.freqs[chan_max] / 1e6, self.channel_width / 1e3))
        # Discover maximum name and tag string lengths for targets beforehand
        name_len, tag_len = 4, 4
        for n in self.target_indices:
            target = self.catalogue.targets[n]
            name_len = max(len(target.name), name_len)
            tag_len = max(len(' '.join(target.tags[1:] if target.body_type != 'xephem' else target.tags[2:])), tag_len)
        descr += ['Targets: %d selected out of %d in catalogue' % (len(self.target_indices), len(self.catalogue)),
                  '  ID  %s  Type      RA(J2000)     DEC(J2000)  %s  Dumps  ModelFlux(Jy)' %
                  ('Name'.ljust(name_len), 'Tags'.ljust(tag_len))]
        for n in self.target_indices:
            target = self.catalogue.targets[n]
            target_type = target.body_type if target.body_type != 'xephem' else target.tags[1]
            tags = ' '.join(target.tags[1:] if target.body_type != 'xephem' else target.tags[2:])
            ra, dec = target.radec() if target_type == 'radec' else ('-', '-')
            # Calculate average target flux over selected frequency band
            flux_spectrum = target.flux_density(self.freqs / 1e6)
            flux_valid = ~np.isnan(flux_spectrum)
            flux = '{:9.2f}'.format(flux_spectrum[flux_valid].mean()) if np.any(flux_valid) else ''
            target_dumps = ((self.sensor.get('Observation/target_index') == n) & self._time_keep).sum()
            descr.append('  %2d  %s  %s  %11s  %11s  %s  %5d  %s' %
                         (n, target.name.ljust(name_len), target_type.ljust(8),
                          ra, dec, tags.ljust(tag_len), target_dumps, flux))
        scans = self.sensor.get('Observation/scan_index')
        compscans = self.sensor.get('Observation/compscan_index')
        total_scans, total_compscans = len(scans.unique_values), len(compscans.unique_values)
        descr += ['Scans: %d selected out of %d total       Compscans: %d selected out of %d total' %
                  (len(self.scan_indices), total_scans, len(self.compscan_indices), total_compscans),
                  '  Date        Timerange(UTC)       ScanState  CompScanLabel  Dumps  Target']
        current_date = ''
        for scan, state, target in self.scans():
            start, end = time.gmtime(self.timestamps[0]), time.gmtime(self.timestamps[-1])
            start_date, start_time = time.strftime('%d-%b-%Y/', start), time.strftime('%H:%M:%S', start)
            end_date, end_time = time.strftime('%d-%b-%Y/', end), time.strftime('%H:%M:%S', end)
            timerange = (start_date if start_date != current_date else '') + start_time
            current_date = start_date
            timerange += ' - ' + (end_date if end_date != current_date else '') + end_time
            compscan, label = self.sensor['Observation/compscan_index'][0], self.sensor['Observation/label'][0]
            descr.append('  %31s  %3d:%5s  %3d:%9s  %5d  %3d:%s' %
                         (timerange, scan, state.ljust(5), compscan, label.ljust(9),
                          self.shape[0], self.target_indices[0], target.name))
        return '\n'.join(descr)

    def _fix_flux_freq_range(self):
        """Ensure that target flux density models are valid on data frequencies."""
        # Pad the minimum and maximum channel freqs with a few extra channels
        min_freq = min((spw.channel_freqs.min() - 2 * spw.channel_width) / 1e6
                       for spw in self.spectral_windows)
        max_freq = max((spw.channel_freqs.max() + 2 * spw.channel_width) / 1e6
                       for spw in self.spectral_windows)
        for target in self.catalogue:
            model = target.flux_model
            if not model:
                continue
            if min_freq < model.min_freq_MHz and \
               min_freq > 0.6 * model.min_freq_MHz or \
               max_freq > model.max_freq_MHz and \
               max_freq < 1.6 * model.max_freq_MHz:
                new_min_freq = min(min_freq, model.min_freq_MHz)
                new_max_freq = max(max_freq, model.max_freq_MHz)
                logger.warning('Extending flux density model frequency range '
                               'of %r from %d-%d MHz to %d-%d MHz', target.name,
                               model.min_freq_MHz, model.max_freq_MHz,
                               new_min_freq, new_max_freq)
                model.min_freq_MHz = new_min_freq
                model.max_freq_MHz = new_max_freq

    def _set_keep(self, time_keep=None, freq_keep=None, corrprod_keep=None,
                  weights_keep=None, flags_keep=None):
        """Set time, frequency and/or correlation product selection masks.

        Set the selection masks for those parameters that are present. Also
        include weights and flags selections as options.

        Parameters
        ----------
        time_keep : array of bool, shape (*T*,), optional
            Boolean selection mask with one entry per timestamp
        freq_keep : array of bool, shape (*F*,), optional
            Boolean selection mask with one entry per frequency channel
        corrprod_keep : array of bool, shape (*B*,), optional
            Boolean selection mask with one entry per correlation product
        weights_keep : 'all' or string or sequence of strings, optional
            Names of selected weight types (or 'all' for the lot)
        flags_keep : 'all' or string or sequence of strings, optional
            Names of selected flag types (or 'all' for the lot)

        """
        if time_keep is not None:
            self._time_keep = time_keep
            # Ensure that sensor cache gets updated time selection
            if self.sensor:
                self.sensor._set_keep(self._time_keep)
        if freq_keep is not None:
            self._freq_keep = freq_keep
        if corrprod_keep is not None:
            self._corrprod_keep = corrprod_keep
        if weights_keep is not None:
            self._weights_keep = weights_keep
        if flags_keep is not None:
            self._flags_keep = flags_keep

    def select(self, **kwargs):
        """Select subset of data, based on time / frequency / corrprod filters.

        This applies a set of selection criteria to the data set, which updates
        the data set properties and attributes to match the selection. In other
        words, the :meth:`timestamps` and :meth:`vis` methods will return the
        selected subset of the data, while attributes such as :attr:`ants`,
        :attr:`channel_freqs` and :attr:`shape` are updated. The sensor cache
        will also return the selected subset of sensor data via the __getitem__
        interface. This function returns nothing, but modifies the existing
        data set in-place.

        The selection criteria are divided into groups, based on whether they
        affect the time, frequency or correlation product dimension::

        * Time: `dumps`, `timerange`, `scans`, `compscans`, `targets`
        * Frequency: `channels`, `freqrange`
        * Correlation product: `corrprods`, `ants`, `inputs`, `pol`

        The `subarray` and `spw` criteria are special, as they affect multiple
        dimensions (time + correlation product and time + frequency,
        respectively), are always active and are forced to be a single index.

        If there are multiple criteria on the same dimension within a select()
        call, they are ANDed together, while multiple items within the same
        criterion (e.g. `targets=['Hyd A', 'Vir A']`) are ORed together. When a
        second select() call is done, all new selections replace previous
        selections on the same dimension, while existing selections on other
        dimensions are preserved. The `reset` parameter finetunes this behaviour.

        If :meth:`select` is called without any parameters the selection is
        reset to the original data set.

        In addition, the `weights` and `flags` criteria are lists of names that
        select which weights and flags to include in the corresponding data set
        property.

        Parameters
        ----------
        strict : {True, False}, optional
            True if select() raises TypeError if it encounters an unknown kwarg

        dumps : int or slice or sequence of ints or sequence of bools, optional
            Select dumps by index, slice or boolean mask of length *T*
            (keep dumps where mask is True)
        timerange : sequence of 2 :class:`katpoint.Timestamp` objects
                    or equivalent, optional
            Select range of times between given start and end times
        scans : int or string or sequence, optional
            Select scans by index or state (or negate state by prepending '~')
        compscans : int or string or sequence, optional
            Select compscans by index or label (or negate label by prepending '~')
        targets : int or string or :class:`katpoint.Target` object or sequence,
                  optional
            Select targets by index or name or description or object

        spw : int, optional
            Select spectral window by index (only one may be active)
        channels : int or slice or sequence of ints or sequence of bools, optional
            Select frequency channels by index, slice or boolean mask of length
            *F* (keep channels where mask is True)
        freqrange : sequence of 2 floats, optional
            Select range of frequencies between start and end frequencies, in Hz

        subarray : int, optional
            Select subarray by index (only one may be active)
        corrprods : int or slice or sequence of ints or sequence of bools or
                    sequence of string pairs or {'auto', 'cross'}, optional
            Select correlation products by index, slice or boolean mask of length
            *B* (keep products where mask is True). Alternatively, select by
            value via a sequence of string pairs, or select all autocorrelations
            via 'auto' or all cross-correlations via 'cross'.
        ants : string or :class:`katpoint.Antenna` object or sequence, optional
            Select antennas by name or object. If all antennas specified are
            prefaced by a ~ this is treated as a deselection and these antennas
            are excluded.
        inputs : string or sequence of strings, optional
            Select inputs by label
        pol : string or sequence of strings
              {'H', 'V', 'HH', 'VV', 'HV', 'VH'}, optional
            Select polarisation terms

        weights : 'all' or string or sequence of strings, optional
            List of names of weights to be multiplied together, as a sequence
            or string of comma-separated names (combine all weights by default)
        flags : 'all' or string or sequence of strings, optional
            List of names of flags that will be OR'ed together, as a sequence
            or string of comma-separated names (use all flags by default). An
            empty string or sequence discards all flags.

        reset : {'auto', '', 'T', 'F', 'B', 'TF', 'TB', 'FB', 'TFB'}, optional
            Remove existing selections on specified dimensions before applying
            the new selections. The default 'auto' option clears those dimensions
            that will be modified by the new selections and leaves the selections
            on unaffected dimensions intact except if `select` is called without
            any parameters, in which case all selections are cleared. By setting
            `reset` to '', new selections apply on top of existing selections.

        Raises
        ------
        TypeError
            If a keyword argument is unknown and `strict` is enabled
        IndexError
            If `spw` or `subarray` is out of range

        """
        time_selectors = ['dumps', 'timerange', 'scans', 'compscans', 'targets']
        freq_selectors = ['channels', 'freqrange']
        corrprod_selectors = ['corrprods', 'ants', 'inputs', 'pol']
        # Check if keywords are valid and raise exception only if this is explicitly enabled
        valid_kwargs = time_selectors + freq_selectors + corrprod_selectors + \
            ['spw', 'subarray', 'weights', 'flags', 'reset', 'strict']
        # Check for definition of strict
        strict = kwargs.get('strict', True)
        if strict and set(kwargs.keys()) - set(valid_kwargs):
            raise TypeError("select() got unexpected keyword argument(s) %s, valid ones are %s" %
                            (list(set(kwargs.keys()) - set(valid_kwargs)), valid_kwargs))
        # If select() is called without arguments, reset all selections
        reset = 'TFB' if not kwargs else kwargs.pop('reset', 'auto')
        kwargs['spw'] = spw = kwargs.get('spw', self.spw)
        if spw >= len(self.spectral_windows):
            raise IndexError('Data set has %d spectral window(s): spw should be in range 0..%d, is %d instead' %
                             (len(self.spectral_windows), len(self.spectral_windows) - 1, spw))
        kwargs['subarray'] = subarray = kwargs.get('subarray', self.subarray)
        if subarray >= len(self.subarrays):
            raise IndexError('Data set has %d subarray(s): subarray should be in range 0..%d, is %d instead' %
                             (len(self.subarrays), len(self.subarrays) - 1, subarray))
        # In 'auto' mode, only reset flags for those dimensions that will be affected by selectors
        if reset == 'auto':
            reset = 'T' if set(kwargs.keys()).intersection(time_selectors) else ''
            reset += 'F' if set(kwargs.keys()).intersection(freq_selectors) else ''
            reset += 'B' if set(kwargs.keys()).intersection(corrprod_selectors) else ''
        # Change spectral window and/or subarray
        if spw != self.spw:
            reset += 'TF'
            self.spw = spw
        if subarray != self.subarray:
            reset += 'TB'
            self.subarray = subarray
        # Reset the selection flags on the appropriate dimensions
        if 'T' in reset:
            self._time_keep[:] = True
            self._time_keep &= (self.sensor.get('Observation/spw_index') == spw)
            self._time_keep &= (self.sensor.get('Observation/subarray_index') == subarray)
            for key in time_selectors:
                self._selection.pop(key, None)
        # Since the number of freqs / corrprods may change due to spw / subarray, create these flags afresh
        if 'F' in reset:
            self._freq_keep = np.ones(self.spectral_windows[self.spw].num_chans, dtype=np.bool)
            for key in freq_selectors:
                self._selection.pop(key, None)
        if 'B' in reset:
            self._corrprod_keep = np.ones(len(self.subarrays[self.subarray].corr_products), dtype=np.bool)
            for key in corrprod_selectors:
                self._selection.pop(key, None)
        # Now add the new selection criteria to the list (after the existing ones were kept or culled)
        self._selection.update(kwargs)

        for k, v in self._selection.items():
            # Selections that affect time axis
            if k == 'dumps':
                if np.asarray(v).dtype == np.bool:
                    self._time_keep &= v
                else:
                    dump_keep = np.zeros(len(self._time_keep), dtype=np.bool)
                    dump_keep[v] = True
                    self._time_keep &= dump_keep
            elif k == 'timerange':
                start_time = katpoint.Timestamp(v[0]).secs + 0.5 * self.dump_period
                end_time = katpoint.Timestamp(v[1]).secs - 0.5 * self.dump_period
                self._time_keep &= (self.sensor.timestamps[:] >= start_time)
                self._time_keep &= (self.sensor.timestamps[:] <= end_time)
            elif k in ('scans', 'compscans'):
                scans = _selection_to_list(v)
                scan_keep = np.zeros(len(self._time_keep), dtype=np.bool)
                scan_sensor = self.sensor.get('Observation/scan_state' if k == 'scans' else 'Observation/label')
                scan_index_sensor = self.sensor.get(f'Observation/{k[:-1]}_index')
                for scan in scans:
                    if isinstance(scan, numbers.Integral):
                        scan_keep |= (scan_index_sensor == scan)
                    elif scan[0] == '~':
                        scan_keep |= ~(scan_sensor == scan[1:])
                    else:
                        scan_keep |= (scan_sensor == scan)
                self._time_keep &= scan_keep
            elif k == 'targets':
                targets = v if is_iterable(v) else [v]
                target_keep = np.zeros(len(self._time_keep), dtype=np.bool)
                target_index_sensor = self.sensor.get('Observation/target_index')
                for t in targets:
                    try:
                        if isinstance(t, numbers.Integral):
                            target_index = t
                        elif isinstance(t, katpoint.Target) or isinstance(t, str) and ',' in t:
                            target_index = self.catalogue.targets.index(t)
                        else:
                            target_index = self.catalogue.targets.index(self.catalogue[t])
                    except (KeyError, ValueError):
                        # Warn here, in case the user gets the target subtly wrong and wonders why it is not selected
                        logger.warning("Skipping unknown selected target '%s'", t)
                        continue
                    target_keep |= (target_index_sensor == target_index)
                self._time_keep &= target_keep
            # Selections that affect frequency axis
            elif k == 'channels':
                if np.asarray(v).dtype == np.bool:
                    self._freq_keep &= v
                else:
                    chan_keep = np.zeros(len(self._freq_keep), dtype=np.bool)
                    chan_keep[v] = True
                    self._freq_keep &= chan_keep
            elif k == 'freqrange':
                start_freq = v[0] + 0.5 * self.spectral_windows[self.spw].channel_width
                end_freq = v[1] - 0.5 * self.spectral_windows[self.spw].channel_width
                self._freq_keep &= (self.spectral_windows[self.spw].channel_freqs >= start_freq)
                self._freq_keep &= (self.spectral_windows[self.spw].channel_freqs <= end_freq)
            # Selections that affect corrprod axis
            elif k == 'corrprods':
                if v == 'auto':
                    self._corrprod_keep &= [(inpA[:-1] == inpB[:-1])
                                            for inpA, inpB in self.subarrays[self.subarray].corr_products]
                elif v == 'cross':
                    self._corrprod_keep &= [(inpA[:-1] != inpB[:-1])
                                            for inpA, inpB in self.subarrays[self.subarray].corr_products]
                else:
                    v = np.asarray(v)
                    if v.ndim == 2 and v.shape[1] == 2:
                        all_corrprods = self.subarrays[self.subarray].corr_products
                        v = v.tolist()
                        v = np.array([list(cp) in v for cp in all_corrprods])
                    if np.asarray(v).dtype == np.bool:
                        self._corrprod_keep &= v
                    else:
                        cp_keep = np.zeros(len(self._corrprod_keep), dtype=np.bool)
                        cp_keep[v] = True
                        self._corrprod_keep &= cp_keep
            elif k == 'ants':
                ants = _selection_to_list(v)
                ant_names = [(ant.name if isinstance(ant, katpoint.Antenna) else ant) for ant in ants]
                if _is_deselection(ant_names):
                    ant_names = [ant_name[1:] for ant_name in ant_names]
                    self._corrprod_keep &= [(inpA[:-1] not in ant_names and inpB[:-1] not in ant_names)
                                            for inpA, inpB in self.subarrays[self.subarray].corr_products]
                else:
                    self._corrprod_keep &= [(inpA[:-1] in ant_names and inpB[:-1] in ant_names)
                                            for inpA, inpB in self.subarrays[self.subarray].corr_products]
            elif k == 'inputs':
                inps = _selection_to_list(v)
                self._corrprod_keep &= [(inpA in inps and inpB in inps)
                                        for inpA, inpB in self.subarrays[self.subarray].corr_products]
            elif k == 'pol':
                pols = _selection_to_list(v)
                # Lower case and strip out empty strings
                pols = [i.lower() for i in pols if i]

                # Proceed if we have a selection
                if len(pols) > 0:
                    # If given a selection assume we keep nothing
                    keep = np.zeros(self._corrprod_keep.shape, dtype=np.bool)

                    # or separate polarisation selections together
                    for polAB in pols:
                        polAB = polAB * 2 if polAB in ('h', 'v') else polAB
                        keep |= [(inpA[-1] == polAB[0] and inpB[-1] == polAB[1])
                                 for inpA, inpB in self.subarrays[self.subarray].corr_products]

                    # and into final corrprod selection
                    self._corrprod_keep &= keep

            # Selections that affect weights and flags
            elif k == 'weights':
                self._weights_keep = v
            elif k == 'flags':
                self._flags_keep = v

        # Update the relevant data members based on selection made
        # These would all be more efficient as properties, but at the expense of extra lines of code...
        self.shape = (self._time_keep.sum(), self._freq_keep.sum(), self._corrprod_keep.sum())
        self.size = np.prod(self.shape, dtype=np.int64) * np.dtype('complex64').itemsize
        if not self.size:
            logger.warning('The selection criteria resulted in an empty data set')
        # This is quicker than indexing np.arange()
        self.dumps = self._time_keep.nonzero()[0]
        self.channels = self._freq_keep.nonzero()[0]
        self.freqs = self.channel_freqs = self.spectral_windows[self.spw].channel_freqs[self._freq_keep]
        self.channel_width = self.spectral_windows[self.spw].channel_width
        self.corr_products = self.subarrays[self.subarray].corr_products[self._corrprod_keep]
        self.inputs = sorted(set(np.ravel(self.corr_products)))
        input_ants = {inp[:-1] for inp in self.inputs}
        self.ants = [ant for ant in self.subarrays[self.subarray].ants if ant.name in input_ants]
        # Ensure that updated selections make their way to sensor cache, as
        # well as any underlying datasets and data lazy indexers that need it
        self._set_keep(self._time_keep, self._freq_keep, self._corrprod_keep,
                       self._weights_keep, self._flags_keep)
        # Figure out which scans, compscans and targets are included in selection
        self.scan_indices = sorted(set(self.sensor['Observation/scan_index']))
        self.compscan_indices = sorted(set(self.sensor['Observation/compscan_index']))
        self.target_indices = sorted(set(self.sensor['Observation/target_index']))

    def scans(self):
        """Generator that iterates through scans in data set.

        This iterates through the currently selected list of scans, returning
        the scan index, scan state and associated target object. In addition,
        after each iteration the data set will reflect the scan selection, i.e.
        the timestamps, visibilities, sensor values, etc. will be those of the
        current scan. The scan selection applies on top of any existing
        selection.

        Yields
        ------
        scan : int
            Scan index
        state : string
            Scan state
        target : :class:`katpoint.Target` object
            Target associated with scan

        """
        scans = self.scan_indices[:]
        # This is the active selection onto which scan selection will be added
        preselection = dict(self._selection.items())
        # This will ensure that the original selection is properly restored
        preselection['reset'] = 'T'
        old_timekeep = self._time_keep.copy()
        state_data = self.sensor.get('Observation/scan_state')
        for scan in scans:
            # Add scan selection on top of existing selection
            self.select(scans=scan, reset='')
            state = state_data.unique_values[state_data.indices[scan]]
            target = self.catalogue.targets[self.target_indices[0]]
            yield scan, state, target
            # A quick way to reset the time selection to the original one
            self._set_keep(old_timekeep.copy())
            self._selection.pop('scans', None)
        # Restore original selection more thoroughly
        self.select(**preselection)

    def compscans(self):
        """Generator that iterates through compound scans in data set.

        This iterates through the currently selected list of compound scans,
        returning the compound scan index, label and the first associated target
        object. In addition, after each iteration the data set will reflect the
        compound scan selection, i.e. the timestamps, visibilities, sensor
        values, etc. will be those of the current compound scan. The compound
        scan selection applies on top of any existing selection.

        Yields
        ------
        compscan : int
            Compound scan index
        label : string
            Compound scan label
        target : :class:`katpoint.Target` object
            First target associated with compound scan

        """
        compscans = self.compscan_indices[:]
        # This is the active selection onto which compscan selection will be added
        preselection = dict(list(self._selection.items()))
        # This will ensure that the original selection is properly restored
        preselection['reset'] = 'T'
        old_timekeep = self._time_keep.copy()
        for compscan in compscans:
            # Add scan selection on top of existing selection
            self.select(compscans=compscan, reset='')
            label_data = self.sensor.get('Observation/label')
            label = label_data.unique_values[label_data.indices[compscan]]
            target = self.catalogue.targets[self.target_indices[0]]
            yield compscan, label, target
            # A quick way to reset the time selection to the original one
            self._set_keep(old_timekeep.copy())
            self._selection.pop('compscans', None)
        # Restore original selection more thoroughly
        self.select(**preselection)

    # - - - - - - - - - - - - - - Format-specific properties - - - - - - - - - - - - - - - - - -

    @property
    def timestamps(self):
        """Visibility timestamps in UTC seconds since Unix epoch.

        The timestamps are returned as an array of float64, shape (*T*,), with
        one timestamp per integration aligned with the integration *midpoint*.

        """
        raise NotImplementedError

    @property
    def vis(self):
        r"""Complex visibility data as a function of time, frequency and corrprod.

        The visibility data are returned as an array of complex64, shape
        (*T*, *F*, *B*), with time along the first dimension, frequency along the
        second dimension and correlation product ("baseline") index along the
        third dimension. The array always has all three dimensions, even for
        scalar (single) values. The number of integrations *T* matches the
        length of :meth:`timestamps`, the number of frequency channels *F*
        matches the length of :meth:`freqs` and the number of correlation
        products *B* matches the length of :meth:`corr_products`.

        The sign convention of the imaginary part is consistent with an
        electric field of :math:`e^{i(\omega t - jz)}` i.e. phase that
        increases with time.
        """
        raise NotImplementedError

    @property
    def weights(self):
        """Visibility weights as a function of time, frequency and baseline.

        The weights data are returned as an array indexer of float32, shape
        (*T*, *F*, *B*), with time along the first dimension, frequency along the
        second dimension and correlation product ("baseline") index along the
        third dimension. The number of integrations *T* matches the length of
        :meth:`timestamps`, the number of frequency channels *F* matches the
        length of :meth:`freqs` and the number of correlation products *B*
        matches the length of :meth:`corr_products`.

        """
        raise NotImplementedError

    @property
    def flags(self):
        """Visibility flags as a function of time, frequency and baseline.

        The flags data are returned as an array indexer of bool, shape
        (*T*, *F*, *B*), with time along the first dimension, frequency along the
        second dimension and correlation product ("baseline") index along the
        third dimension. The number of integrations *T* matches the length of
        :meth:`timestamps`, the number of frequency channels *F* matches the
        length of :meth:`freqs` and the number of correlation products *B*
        matches the length of :meth:`corr_products`.

        """
        raise NotImplementedError

    @property
    def temperature(self):
        """Air temperature in degrees Celsius."""
        raise NotImplementedError

    @property
    def pressure(self):
        """Barometric pressure in millibars."""
        raise NotImplementedError

    @property
    def humidity(self):
        """Relative humidity as a percentage."""
        raise NotImplementedError

    @property
    def wind_speed(self):
        """Wind speed in metres per second."""
        raise NotImplementedError

    @property
    def wind_direction(self):
        """Wind direction as an azimuth angle in degrees."""
        raise NotImplementedError

    # - - - - - - - - - - - - Virtual sensors exposed as properties - - - - - - - - - - - - - - -

    @property
    def mjd(self):
        """Visibility timestamps in Modified Julian Days (MJD).

        The timestamps are returned as an array of float64, shape (*T*,), with
        one timestamp per integration aligned with the integration *midpoint*.

        """
        return self.sensor['Timestamps/mjd']

    @property
    def lst(self):
        """Local sidereal time at the reference antenna in hours.

        The sidereal times are returned in an array of float, shape (*T*,).

        """
        return self.sensor[f'Antennas/{self.ref_ant}/lst'] * (12 / np.pi)

    def _sensor_per_ant(self, base_name):
        """Extract a single sensor per antenna and safely stack the results."""
        def sensor_data(ant_name):
            return self.sensor[f'Antennas/{ant_name}/{base_name}']
        return np.column_stack([sensor_data(ant.name) for ant in self.ants]) \
            if self.ants else np.zeros((self.shape[0], 0))

    def _delta_sensor_per_corrprod(self, base_name):
        """Extract a single sensor per corrprod and safely stack the results.

        The sensor is specialised to return the difference between the sensors
        of the two antennas making up the correlation product, since that is
        what (u, v, w) sensors need.
        """
        def difference(antA, antB):
            coord1 = self.sensor[f'Antennas/{antA}/{base_name}']
            coord2 = self.sensor[f'Antennas/{antB}/{base_name}']
            return coord1 - coord2
        return np.column_stack([difference(inpA[:-1], inpB[:-1])
                                for inpA, inpB in self.corr_products]) \
            if len(self.corr_products) else np.zeros((self.shape[0], 0))

    @property
    def az(self):
        """Azimuth angle of each dish in degrees.

        The azimuth angles are returned in an array of float, shape (*T*, *A*).

        """
        return rad2deg(self._sensor_per_ant('az'))

    @property
    def el(self):
        """Elevation angle of each dish in degrees.

        The elevation angles are returned in an array of float, shape (*T*, *A*).

        """
        return rad2deg(self._sensor_per_ant('el'))

    @property
    def ra(self):
        """Right ascension of the actual pointing of each dish in J2000 degrees.

        The right ascensions are returned in an array of float, shape (*T*, *A*).

        """
        return rad2deg(self._sensor_per_ant('ra'))

    @property
    def dec(self):
        """Declination of the actual pointing of each dish in J2000 degrees.

        The declinations are returned in an array of float, shape (*T*, *A*).

        """
        return rad2deg(self._sensor_per_ant('dec'))

    @property
    def parangle(self):
        """Parallactic angle of the actual pointing of each dish in degrees.

        The parallactic angle is the position angle of the observer's vertical
        on the sky, measured from north toward east. This is the angle between
        the great-circle arc connecting the celestial North pole to the dish
        pointing direction, and the great-circle arc connecting the zenith above
        the antenna to the pointing direction, or the angle between the
        *hour circle* and *vertical circle* through the pointing direction.
        It is returned as an array of float, shape (*T*, *A*).

        """
        return rad2deg(self._sensor_per_ant('parangle'))

    @property
    def target_x(self):
        """Target *x* coordinate of each dish in degrees.

        The target coordinates are projections of the spherical coordinates of
        the dish pointing direction to a plane with the target position at the
        origin. The type of projection (e.g. ARC, SIN, etc.) and spherical
        pointing coordinate system (e.g. azel or radec) can be set via the
        :attr:`target_projection` and :attr:`target_coordsys` attributes,
        respectively. The target *x* coordinates are returned as an array of
        float, shape (*T*, *A*).

        """
        name = f'target_x_{self.target_projection}_{self.target_coordsys}'
        return rad2deg(self._sensor_per_ant(name))

    @property
    def target_y(self):
        """Target *y* coordinate of each dish in degrees.

        The target coordinates are projections of the spherical coordinates of
        the dish pointing direction to a plane with the target position at the
        origin. The type of projection (e.g. ARC, SIN, etc.) and spherical
        pointing coordinate system (e.g. azel or radec) can be set via the
        :attr:`target_projection` and :attr:`target_coordsys` attributes,
        respectively. The target *y* coordinates are returned as an array of
        float, shape (*T*, *A*).

        """
        name = f'target_y_{self.target_projection}_{self.target_coordsys}'
        return rad2deg(self._sensor_per_ant(name))

    @property
    def u(self):
        """U coordinate for each correlation product in metres.

        This calculates the *u* coordinate of the baseline vector of each
        correlation product as a function of time while tracking the target.
        It is returned as an array of float, shape (*T*, *B*). The sign
        convention is :math:`u_1 - u_2` for baseline (ant1, ant2).

        """
        return self._delta_sensor_per_corrprod('u')

    @property
    def v(self):
        """V coordinate for each correlation product in metres.

        This calculates the *v* coordinate of the baseline vector of each
        correlation product as a function of time while tracking the target.
        It is returned as an array of float, shape (*T*, *B*). The sign
        convention is :math:`v_1 - v_2` for baseline (ant1, ant2).

        """
        return self._delta_sensor_per_corrprod('v')

    @property
    def w(self):
        """W coordinate for each correlation product in metres.

        This calculates the *w* coordinate of the baseline vector of each
        correlation product as a function of time while tracking the target.
        It is returned as an array of float, shape (*T*, *B*).The sign
        convention is :math:`w_1 - w_2` for baseline (ant1, ant2).

        """
        return self._delta_sensor_per_corrprod('w') class Subarray:
    """Subarray specification.

    A subarray is determined by the specific correlation products produced by the
    correlator and the antenna objects associated with the inputs found in the
    correlation products.

    Parameters
    ----------
    ants : sequence of :class:`katpoint.Antenna` objects
        List of antenna objects, culled to contain only antennas found in
        `corr_products`
    corr_products : sequence of (string, string) pairs, length *B*
        Correlation products as pairs of input labels, e.g. ('ant1h', 'ant2v'),
        exposed as an array of strings with shape (*B*, 2)

    Attributes
    ----------
    inputs : list of strings
        List of correlator input labels found in `corr_products`, e.g. 'ant1h'

    """

    def __init__(self, ants, corr_products):
        self.corr_products = np.array([(inpA.lower(), inpB.lower())
                                       for inpA, inpB in corr_products])
        # Extract all inputs (and associated antennas) from corr product list
        self.inputs = sorted(set(np.ravel(self.corr_products)))
        input_ants = {inp[:-1] for inp in self.inputs}
        # Only keep antennas that are involved in correlation products
        self.ants = [ant for ant in ants if ant.name in input_ants]

    def __repr__(self):
        """Short human-friendly string representation of subarray object."""
        return "<katdal.Subarray antennas={} inputs={} corrprods={} at {:#x}>".format(
               len(self.ants), len(self.inputs), len(self.corr_products), id(self))

    @property
    def _description(self):
        """Complete string representation, used internally for comparisons."""
        ants = '\n'.join(ant.description for ant in self.ants)
        corrprods = ' '.join(f'{inpA},{inpB}' for inpA, inpB in self.corr_products)
        return '\n'.join((ants, corrprods))

    def __eq__(self, other):
        """Equality comparison operator."""
        return self._description == (other._description
                                     if isinstance(other, Subarray) else other)

    def __ne__(self, other):
        """Inequality comparison operator."""
        return not (self == other)

    def __lt__(self, other):
        """Less-than comparison operator (needed for sorting and np.unique)."""
        return self._description < (other._description
                                    if isinstance(other, Subarray) else other)

    def __hash__(self):
        """Base hash on description string, just like equality operator."""
        return hash(self._description) class WrongVersion(Exception):
    """Trying to access data using accessor class with the wrong version.""" def _robust_target(description):
    """Robust build of :class:`katpoint.Target` object from description string."""
    if not description:
        return katpoint.Target('Nothing, special')
    try:
        return katpoint.Target(description)
    except ValueError:
        logger.warning("Invalid target description '%s' - replaced with dummy target", description)
        return katpoint.Target('Nothing, special') def _selection_to_list(names, **groups):
    """Normalise string of comma-separated names or sequence of names / objects.

    Parameters
    ----------
    names : string / object or sequence of strings / objects
        A string of comma-separated names or a sequence of names / objects
    groups : dict, optional
        Each extra keyword is the name of a predefined list of names / objects

    Returns
    -------
    list : list of strings / objects
        List of names / objects
    """
    if isinstance(names, str):
        if not names:
            return []
        elif names in groups:
            return list(groups[names])
        else:
            return [name.strip() for name in names.split(',')]
    elif is_iterable(names):
        return list(names)
    else:
        return [names] DESCRIPTIONS = ('reserved - bit 0',
                'predefined static flag list',
                'flag based on live CAM information',
                'no data was received',
                'RFI detected in ingest',
                'RFI predicted from space based pollutants',
                'RFI detected in calibration',
                'some correction/postprocessing step could not be applied') NAMES = ('reserved0', 'static', 'cam', 'data_lost',
         'ingest_rfi', 'predicted_rfi', 'cal_rfi', 'postproc') class LazyIndexer:
    """Two-stage deferred indexer for objects with expensive __getitem__ calls.

    This class was originally designed to extend and speed up the indexing
    functionality of HDF5 datasets as found in :mod:`h5py`, but works on any
    equivalent object (defined as any object with `shape`, `dtype` and
    `__getitem__` members) where a call to __getitem__ may be very expensive.
    The following discussion focuses on the HDF5 use case as the main example.

    Direct extraction of a subset of an HDF5 dataset via the __getitem__
    interface (i.e. `dataset[index]`) has a few issues:

    1. Data access can be very slow (or impossible) if a very large dataset is
       fully loaded into memory and then indexed again at a later stage
    2. Advanced indexing (via boolean masks or sequences of integer indices) is
       only supported on a single dimension in the current version of h5py (2.0)
    3. Even though advanced indexing has limited support, simple indexing (via
       single integer indices or slices) is frequently much faster.

    This class wraps an :class:`h5py.Dataset` or equivalent object and exposes a
    new __getitem__ interface on it. It efficiently composes two stages of
    indexing: a first stage specified at object instantiation time and a second
    stage that applies on top of the first stage when __getitem__ is called on
    this object. The data are only loaded after the combined index is determined,
    addressing issue 1.

    Furthermore, advanced indexing is allowed on any dimension by decomposing
    the selection as a series of slice selections covering contiguous segments
    of the dimension to alleviate issue 2. Finally, this also allows faster
    data retrieval by extracting a large slice from the HDF5 dataset and then
    performing advanced indexing on the resulting :class:`numpy.ndarray` object
    instead, in response to issue 3.

    The `keep` parameter of the :meth:`__init__` and :meth:`__getitem__` methods
    accepts a generic index or slice specification, i.e. anything that would be
    accepted by the :meth:`__getitem__` method of a :class:`numpy.ndarray` of
    the same shape as the dataset. This could be a single integer index, a
    sequence of integer indices, a slice object (representing the colon operator
    commonly used with __getitem__, e.g. representing `x[1:10:2]` as
    `x[slice(1,10,2)]`), a sequence of booleans as a mask, or a tuple containing
    any number of these (typically one index item per dataset dimension). Any
    missing dimensions will be fully selected, and any extra dimensions will
    be ignored.

    Parameters
    ----------
    dataset : :class:`h5py.Dataset` object or equivalent
        Underlying dataset or array object on which lazy indexing will be done.
        This can be any object with shape, dtype and __getitem__ members.
    keep : NumPy index expression, optional
        First-stage index as a valid index or slice specification
        (supports arbitrary slicing or advanced indexing on any dimension)
    transforms : list of :class:`LazyTransform` objects or None, optional
        Chain of transforms to be applied to data after final indexing. The
        chain as a whole may only add or drop dimensions at the end of data
        shape without changing the preserved dimensions.

    Attributes
    ----------
    name : string
        Name of HDF5 dataset (or empty string for unnamed ndarrays, etc.)

    Raises
    ------
    InvalidTransform
        If transform chain does not obey restrictions on changing the data shape

    """

    def __init__(self, dataset, keep=slice(None), transforms=None):
        self.dataset = dataset
        self.transforms = [] if transforms is None else transforms
        self.name = getattr(self.dataset, 'name', '')
        # Ensure that keep is a tuple (then turn it into a list to simplify further processing)
        keep = list(keep) if isinstance(keep, tuple) else [keep]
        # Ensure that keep is same length as data shape (truncate or pad with blanket slices as necessary)
        keep = keep[:len(dataset.shape)] + [slice(None)] * (len(dataset.shape) - len(keep))
        # Ensure that each index in lookup is an array of integer indices, or None
        self._lookup = []
        for dim_keep, dim_len in zip(keep, dataset.shape):
            if isinstance(dim_keep, slice):
                # Turn slice into array of integer indices (or None if it contains all indices)
                dim_keep = dim_keep.indices(dim_len)
                dim_keep = np.arange(*dim_keep) if dim_keep != slice(None).indices(dim_len) else None
            else:
                dim_keep = np.atleast_1d(dim_keep)
                # Turn boolean mask into integer indices (True means keep that index), or None if all is True
                if dim_keep.dtype == np.bool and len(dim_keep) == dim_len:
                    dim_keep = np.nonzero(dim_keep)[0] if not dim_keep.all() else None
            self._lookup.append(dim_keep)
        # Shape of data array after first-stage indexing and before transformation
        self._initial_shape = tuple([(len(dim_keep) if dim_keep is not None else dim_len)
                                     for dim_keep, dim_len in zip(self._lookup, self.dataset.shape)])
        # Type of data array before transformation
        self._initial_dtype = self.dataset.dtype
        # Test validity of shape and dtype
        self.shape, self.dtype

    def __repr__(self):
        """Short human-friendly string representation of lazy indexer object."""
        return "<katdal.{} '{}': shape {}, type {} at {:#x}>".format(
               self.__class__.__name__, self.name, self.shape, self.dtype, id(self))

    def _name_shape_dtype(self, name, shape, dtype):
        """Helper function to create strings for display (limits dtype length)."""
        dtype_str = (str(dtype)[:50] + '...') if len(str(dtype)) > 50 else str(dtype)
        return f"{name} -> {shape} {dtype_str}"

    def __str__(self):
        """Verbose human-friendly string representation of lazy indexer object."""
        shape, dtype = self._initial_shape, self._initial_dtype
        descr = [self._name_shape_dtype(self.name, shape, dtype)]
        for transform in self.transforms:
            shape, dtype = transform.new_shape(shape), transform.dtype if transform.dtype is not None else dtype
            descr += ['-> ' + self._name_shape_dtype(transform.name, shape, dtype)]
        return '\n'.join(descr)

    def __len__(self):
        """Length operator."""
        return self.shape[0]

    def __iter__(self):
        """Iterator."""
        for index in range(len(self)):
            yield self[index]

    def __getitem__(self, keep):
        """Extract a selected array from the underlying dataset.

        This applies the given second-stage index on top of the first-stage index
        and retrieves the relevant data from the dataset as an array, optionally
        transforming it afterwards.

        Parameters
        ----------
        keep : NumPy index expression
            Second-stage index as a valid index or slice specification
            (supports arbitrary slicing or advanced indexing on any dimension)

        Returns
        -------
        data : array
            Extracted output array

        """
        ndim = len(self.dataset.shape)
        # Ensure that keep is a tuple (then turn it into a list to simplify further processing)
        keep = list(keep) if isinstance(keep, tuple) else [keep]
        # The original keep tuple will be passed to data transform chain
        original_keep = tuple(keep)
        # Ensure that keep is same length as data dimension (truncate or pad with blanket slices as necessary)
        keep = keep[:ndim] + [slice(None)] * (ndim - len(keep))
        # Map current selection to original data indices based on any existing initial selection, per data dimension
        keep = [(dkeep if dlookup is None else dlookup[dkeep]) for dkeep, dlookup in zip(keep, self._lookup)]
        # Iterate over dimensions of dataset, storing information on selection on each dimension:
        # `selection` is a list with one element per dimension; each element is a list of contiguous segments along
        # the dimension, and each segment is represented by a tuple of 3 elements:
        # (dataset selection, post-selection, output array selection)
        # Similarly, `segment_sizes` is a list of lists of segment lengths (empty lists for scalar-selected dimensions)
        selection, segment_sizes = [], []
        for dim_keep, dim_len in zip(keep, self.dataset.shape):
            if np.isscalar(dim_keep):
                # If selection is a scalar, pass directly to dataset selector and remove dimension from output
                selection.append([(dim_keep, None, None)])
                segment_sizes.append([])
            elif isinstance(dim_keep, slice):
                # If selection is a slice, pass directly to dataset selector without post-selection
                start, stop, stride = dim_keep.indices(dim_len)
                segm_size = len(range(start, stop, stride))
                selection.append([(slice(start, stop, stride), slice(None), slice(0, segm_size, 1))])
                segment_sizes.append([segm_size])
            elif len(dim_keep) == 0:
                # If selection is empty, pass to post-selector, as HDF5 datasets do not support zero-length selection
                selection.append([(slice(0, 1, 1), slice(0, 0, 1), slice(0, 0, 1))])
                segment_sizes.append([0])
            else:
                # Anything else is advanced indexing via bool or integer sequences
                dim_keep = np.atleast_1d(dim_keep)
                # Turn boolean mask into integer indices (True means keep that index)
                if dim_keep.dtype == np.bool and len(dim_keep) == dim_len:
                    dim_keep = np.nonzero(dim_keep)[0]
                elif not np.all(dim_keep == np.unique(dim_keep)):
                    raise TypeError('LazyIndexer cannot handle duplicate or unsorted advanced integer indices')
                # Split indices into multiple contiguous segments (specified by first and one-past-last data indices)
                jumps = np.nonzero(np.diff(dim_keep) > 1)[0]
                first = [dim_keep[0]] + dim_keep[jumps + 1].tolist()
                last = dim_keep[jumps].tolist() + [dim_keep[-1]]
                segments = np.c_[first, np.array(last) + 1]
                if len(dim_keep) > 0.2 * dim_len and len(segments) > 1:
                    # If more than 20% of data are selected in 2 or more separate segments (the Ratcliffian benchmark),
                    # select data at dataset level with a single slice spanning segments and then postselect the ndarray
                    selection.append([(slice(segments[0, 0], segments[-1, 1], 1),
                                       dim_keep - dim_keep[0], slice(0, len(dim_keep), 1))])
                    segment_sizes.append([len(dim_keep)])
                else:
                    # Turn each segment into a separate slice at dataset level without post-selection,
                    # and construct contiguous output slices of the same segment sizes
                    segm_sizes = [end - start for start, end in segments]
                    segm_starts = np.cumsum([0] + segm_sizes)
                    selection.append([(slice(start, end, 1), slice(None), slice(segm_starts[n], segm_starts[n + 1], 1))
                                      for n, (start, end) in enumerate(segments)])
                    segment_sizes.append(segm_sizes)
        # Short-circuit the selection if all dimensions are selected with scalars (resulting in a scalar output)
        if segment_sizes == [[]] * ndim:
            out_data = self.dataset[tuple([select[0][0] for select in selection])]
        else:
            # Use dense N-dimensional meshgrid to slice data set into chunks, based on segments along each dimension
            chunk_indices = np.mgrid[[slice(0, len(select), 1) for select in selection]]
            # Pre-allocate output ndarray to have the correct shape and dtype (will be at least 1-dimensional)
            out_data = np.empty([np.sum(segments) for segments in segment_sizes if segments], dtype=self.dataset.dtype)
            # Iterate over chunks, extracting them from dataset and inserting them into the right spot in output array
            for chunk_index in chunk_indices.reshape(ndim, -1).T:
                # Extract chunk from dataset (don't use any advanced indexing here, only scalars and slices)
                dataset_select = tuple([select[segment][0] for select, segment in zip(selection, chunk_index)])
                chunk = self.dataset[dataset_select]
                # Perform post-selection on chunk (can be fancier / advanced indexing because chunk is now an ndarray)
                post_select = [select[segment][1] for select, segment in zip(selection, chunk_index)]
                # If any dimensions were dropped due to scalar indexing, drop them from post_select/out_select tuples
                post_select = tuple([select for select in post_select if select is not None])
                # Do post-selection one dimension at a time, as ndarray does not allow simultaneous advanced indexing
                # on more than one dimension. This caters for the scenario where more than one dimension satisfies
                # the Ratcliffian benchmark (the only way to get advanced post-selection).
                for dim in range(len(chunk.shape)):
                    # Only do post-selection on this dimension if non-trivial (otherwise an unnecessary copy happens)
                    if not (isinstance(post_select[dim], slice) and post_select[dim] == slice(None)):
                        # Prepend the appropriate number of colons to the selection to place it at correct dimension
                        chunk = chunk[[slice(None)] * dim + [post_select[dim]]]
                # Determine appropriate output selection and insert chunk into output array
                out_select = [select[segment][2] for select, segment in zip(selection, chunk_index)]
                out_select = tuple([select for select in out_select if select is not None])
                out_data[out_select] = chunk
        # Apply transform chain to output data, if any
        return reduce(lambda data, transform: transform(data, original_keep), self.transforms, out_data)

    @property
    def shape(self):
        """Shape of data array after first-stage indexing and transformation, i.e. ``self[:].shape``."""
        new_shape = reduce(lambda shape, transform: transform.new_shape(shape), self.transforms, self._initial_shape)
        # Do a quick test of shape transformation as verification of the transform chain
        allowed_shapes = [self._initial_shape[:(n + 1)] for n in range(len(self._initial_shape))]
        if new_shape[:len(self._initial_shape)] not in allowed_shapes:
            raise InvalidTransform('Transform chain may only add or drop dimensions at the end of data shape: '
                                   f'final shape is {new_shape}, expected one of {allowed_shapes}')
        return new_shape

    @property
    def dtype(self):
        """Type of data array after transformation, i.e. ``self[:].dtype``."""
        return reduce(lambda dtype, transform: transform.dtype if transform.dtype is not None else dtype,
                      self.transforms, self._initial_dtype) class LazyTransform:
    """Transformation to be applied by LazyIndexer after final indexing.

    A :class:`LazyIndexer` potentially applies a chain of transforms to the
    data after the final second-stage indexing is done. These transforms are
    restricted in their capabilities to simplify the indexing process.
    Specifically, when it comes to the data shape, transforms may only::

      - add dimensions at the end of the data shape, or
      - drop dimensions at the end of the data shape.

    The preserved dimensions are not allowed to change their shape or
    interpretation so that the second-stage indexing matches the first-stage
    indexing on these dimensions. The data type (aka `dtype`) is allowed to
    change.

    Parameters
    ----------
    name : string or None, optional
        Name of transform
    transform : function, signature ``data = f(data, keep)``, optional
        Transform to apply to data (`keep` is user-specified second-stage index)
    new_shape : function, signature ``new_shape = f(old_shape)``, optional
        Function that predicts data array shape tuple after first-stage indexing
        and transformation, given its original shape tuple as input.
        Restrictions apply as described above.
    dtype : :class:`numpy.dtype` object or equivalent or None, optional
        Type of output array after transformation (None if same as input array)

    """

    def __init__(self, name=None, transform=lambda d, k: d, new_shape=lambda s: tuple(s), dtype=None):
        self.name = 'unnamed' if name is None else name
        self.transform = transform
        self.new_shape = new_shape
        self.dtype = np.dtype(dtype) if dtype is not None else None

    def __repr__(self):
        """Short human-friendly string representation of lazy transform object."""
        class_name = self.__class__.__name__
        dtype = 'unchanged' if self.dtype is None else self.dtype
        return f"<katdal.{class_name} '{self.name}': type '{dtype}' at {id(self):#x}>"

    def __call__(self, data, keep):
        """Transform data (`keep` is user-specified second-stage index)."""
        return self.transform(data, keep) class RecordSensorGetter(SensorGetter):
    """Raw (uninterpolated) sensor data in record array form.

    This is a wrapper for uninterpolated sensor data which resembles a record
    array with fields 'timestamp', 'value' and optionally 'status'. This is
    also the typical format of HDF5 datasets used to store sensor data.

    Technically, the data is interpreted as a NumPy "structured" array, which
    is a simpler version of a recarray that only provides item-style access to
    fields and not attribute-style access.

    Object-valued sensors are not treated specially in this class, as it is
    assumed that any wrapping already occurred in the construction of the
    recarray-like `data` input and will be reflected in its dtype. The original
    HDF5 sensor datasets also did not contain any objects as they only support
    standard KATCP types, so there was no need for wrapping there.

    Parameters
    ----------
    data : recarray-like, with fields 'timestamp', 'value' and optionally 'status'
        Uninterpolated sensor data as structured array or equivalent (such as
        an :class:`h5py.Dataset`)
    name : string or None, optional
        Sensor name (assumed to be data.name by default, if it exists)

    """

    def __init__(self, data, name=None):
        name = name if name is not None else getattr(data, 'name', '')
        super().__init__(name)
        self._data = data

    def get(self):
        """Extract timestamp, value and status of each sensor data point.

        Values are passed through :func:`to_str`.
        """
        timestamp = np.asarray(self._data['timestamp'])
        value = to_str(np.asarray(self._data['value']))
        try:
            status = self._data['status']
        except ValueError:
            status = None
        return SensorData(self.name, timestamp, value, status)

    def __repr__(self):
        """Short human-friendly string representation of sensor data object."""
        return "<katdal.{} '{}' len={} type={} at {:#x}>".format(
               self.__class__.__name__, self.name,
               len(self._data), self._data['value'].dtype, id(self)) class SensorCache(MutableMapping):
    """Container for sensor data providing name lookup, interpolation and caching.

    *Sensor data* is defined as a one-dimensional time series of values. The
    values may be numerical or non-numerical (*categorical*), and the timestamps
    are monotonically increasing but not necessarily regularly spaced.

    A *sensor cache* stores sensor data with dictionary-like lookup based on
    the sensor name. Since the extraction of sensor data from e.g. HDF5 files
    may be costly, the data is first represented in uncached (raw) form as
    :class:`SensorGetter` objects, which typically wrap the underlying sensor
    HDF5 datasets. After extraction, the sensor data are stored either as
    a NumPy array (for numerical data) or as a :class:`CategoricalData` object
    (for non-numerical data).

    The sensor cache stores a timestamp array (or indexer) onto which the sensor
    data will be interpolated, together with a boolean selection mask that
    selects a subset of the interpolated data as the final output. Interpolation
    is linear for numerical data and zeroth-order for non-numerical data. Both
    extraction and selection may be enabled or disabled through the appropriate
    use of the two main interfaces that retrieve sensor data:

    * The __getitem__ interface (i.e. `cache[sensor]`) presents a simple
      high-level interface to the end user that always extracts the sensor data
      and selects the requested subset from it. In addition, the return type is
      always a NumPy array.

    * The get() interface (i.e. `cache.get(sensor)`) is an advanced interface
      for library builders that provides full control of the extraction process
      via *sensor properties*. It does not apply selection by default, as this
      is more convenient for library routines.

    In addition, the sensor cache may contain *virtual sensors* which calculate
    their values based on the values of other sensors. They are identified by
    pattern templates that potentially match multiple sensor names.

    Parameters
    ----------
    cache : mapping from string to :class:`SensorGetter` objects
        Initial sensor cache mapping sensor names to raw (uncached) sensor data
    timestamps : array of float
        Correlator data timestamps onto which sensor values will be interpolated,
        as UTC seconds since Unix epoch
    dump_period : float
        Dump period, in seconds
    keep : int or slice or sequence of int or sequence of bool, optional
        Default time selection specification that will be applied to sensor data
        (this can be disabled on data retrieval)
    props : dict, optional
        Default properties that govern how sensor data are interpreted and
        interpolated (this can be overridden on data retrieval). Can use ``*``
        as a wildcard anywhere in the key.
    virtual : dict mapping string to function, optional
        Virtual sensors, specified as a pattern matching the virtual sensor name
        and a corresponding function that will create the sensor (together with
        any associated virtual sensors)
    aliases : dict mapping string to string, optional
        Alternate names for sensors, as a dictionary mapping each alias to the
        original sensor name suffix. This will create additional sensors with
        the aliased names and the data of the original sensors.
    store : string, optional
        Hostname / endpoint of katstore webserver to access additional sensors
    """

    def __init__(self, cache, timestamps, dump_period, keep=slice(None),
                 props=None, virtual={}, aliases={}, store=None):
        super().__init__()
        # This needs to be an RLock because instantiating a virtual sensor
        # may require further sensor lookups (hopefully without a loop, which
        # would really cause problems).
        self._lock = threading.RLock()
        # Store internals of the cache in a regular dict
        self._raw = dict(cache)
        self.timestamps = timestamps
        self.dump_period = dump_period
        self.keep = keep
        self.props = props if props is not None else {}
        # Add virtual sensor templates
        self.virtual = virtual
        # Add sensor aliases
        for alias, original in aliases.items():
            self.add_aliases(alias, original)
        self.store = store

    def __str__(self):
        """Verbose human-friendly string representation of sensor cache object."""
        with self._lock:
            names = sorted([key for key in self.keys()])
            maxlen = max([len(name) for name in names])
            objects = [self.get(name, extract=False) for name in names]
        obj_reprs = [(f"<numpy.ndarray shape={obj.shape} type={obj.dtype} at {id(obj):#x}>"
                      if isinstance(obj, np.ndarray) else repr(obj)) for obj in objects]
        actual = [f'{name!s:{maxlen}} : {obj_repr}'
                  for name, obj_repr in zip(names, obj_reprs)]
        virtual = [f'{pat!s:{maxlen}} : <function {func.__module__}.{func.__name__}>'
                   for pat, func in self.virtual.items()]
        return '\n'.join(['Actual sensors', '--------------'] + actual +
                         ['\nVirtual sensors', '---------------'] + virtual)

    def __repr__(self):
        """Short human-friendly string representation of sensor cache object."""
        with self._lock:
            sensors = [self.get(name, extract=False) for name in self.keys()]
        return "<katdal.{} sensors={} cached={} virtual={} at {:#x}>".format(
               self.__class__.__name__, len(sensors),
               np.sum([not isinstance(s, SensorGetter) for s in sensors]),
               len(self.virtual), id(self))

    def __getitem__(self, name):
        """Sensor values interpolated to correlator data timestamps.

        Time selection is enforced as this is a user-facing sensor data
        extraction method, and end users expect their selections to apply.
        This has the added benefit of a consistent array return type.

        Parameters
        ----------
        name : string
            Sensor name

        Returns
        -------
        sensor_data : array
           Interpolated sensor data as 1-D array, one value per selected timestamp

        Raises
        ------
        KeyError
            If sensor name was not found in cache

        """
        return self.get(name, select=True)

    def _set_keep(self, keep=None):
        """Set time selection for sensor values."""
        if keep is not None:
            self.keep = keep

    @staticmethod
    def _get_props(name, prop_map, **kwargs):
        """Retrieve properties for a sensor.

        Sensor names in `prop_map` may contain ``*`` wildcard characters. All
        matching entries are merged, and `prop_map` is updated in place with
        the merged result.

        Parameters
        ----------
        name : str
            Sensor name
        prop_map : dict
            Maps sensor names to mappings of properties.
        kwargs
            Extra properties to apply, overriding those in `prop_map`.
        """
        # Look up properties associated with this specific sensor
        props = prop_map.setdefault(name, {})
        # Look up properties associated with this class of sensor
        for key, val in prop_map.items():
            if '*' in key:
                regex = '.*'.join(re.escape(part) for part in key.split('*'))
                if re.match('^' + regex + '$', name):
                    props.update(val)
        # Any properties passed directly to this method takes precedence
        props.update(kwargs)
        return props

    @staticmethod
    def _extract(sensor_getter, timestamps, dump_period, **props):
        sensor_data = sensor_getter.get()
        # Clean up sensor data if non-empty
        if sensor_data:
            time_offset = props.get('time_offset', 0)
            sensor_data.timestamp += time_offset
            # Sort sensor events in chronological order and discard duplicates and unreadable sensor values
            sensor_data = remove_duplicates_and_invalid_values(sensor_data)
        if not sensor_data:
            sensor_data = dummy_sensor_getter(sensor_data.name, value=props.get('initial_value'),
                                              dtype=sensor_data.value.dtype).get()
            logger.warning("No usable data found for sensor '%s' - replaced with dummy data (%r)" %
                           (sensor_data.name, sensor_data.value[0]))
        # Determine if sensor produces categorical or numerical data
        # (float data are non-categorical, by default)
        categ = props.get('categorical', not np.issubdtype(sensor_data.value.dtype, np.floating))
        props['categorical'] = categ
        if categ:
            sensor_data = sensor_to_categorical(sensor_data.timestamp, sensor_data.value,
                                                timestamps, dump_period, **props)
        else:
            # Interpolate numerical data onto data timestamps
            sensor_timestamps = sensor_data.timestamp
            # Warn if sensor data will be extrapolated to start or end
            # of data set with potentially bogus results
            if len(sensor_timestamps) > 1:
                if sensor_timestamps[0] > timestamps[0]:
                    logger.warning("First data point for sensor '%s' only arrives %g seconds into data set" %
                                   (sensor_data.name, sensor_timestamps[0] - timestamps[0]))
                if sensor_timestamps[-1] < timestamps[-1]:
                    logger.warning("Last data point for sensor '%s' arrives %g seconds "
                                   "before end of data set" %
                                   (sensor_data.name, timestamps[-1] - sensor_timestamps[-1]))
            sensor_data = np.interp(timestamps, sensor_timestamps, sensor_data.value)
        return sensor_data

    def add_aliases(self, alias, original):
        """Add alternate names / aliases for sensors.

        Search for sensors with names ending in the `original` suffix and form
        a corresponding alternate name by replacing `original` with `alias`.
        The new aliased sensors will re-use the data of the original sensors.

        Parameters
        ----------
        alias : string
            The new sensor name suffix that replaces `original`
        original : string
            Sensors with names that end in this will get aliases

        """
        for name, data in list(self._raw.items()):
            if name.endswith(original):
                self._raw[name.replace(original, alias)] = data

    def get(self, name, select=False, extract=True, **kwargs):
        """Sensor values interpolated to correlator data timestamps.

        Time selection is disabled by default, as this is a more advanced data
        extraction method typically called by library routines that want to
        operate on the full array of sensor values. For additional allowed
        parameters when extracting categorical data, see the docstring for
        :func:`sensor_to_categorical`.

        Parameters
        ----------
        name : string
            Sensor name
        select : {False, True}, optional
            True if preset time selection will be applied to interpolated data
        extract : {True, False}, optional
            True if sensor data should be extracted, interpolated and cached
        categorical : {None, True, False}, optional
            Interpret sensor data as categorical or numerical (by default, data
            of type float is numerical and of any other type is categorical)
        kwargs : dict, optional
            Additional parameters are passed to :func:`sensor_to_categorical`

        Returns
        -------
        data : array or :class:`CategoricalData` or :class:`SensorGetter` object
            If extraction is disabled, this will be a :class:`SensorGetter` object
            for uncached sensors. If selection is enabled, this will be a 1-D
            array of values, one per selected timestamp. If selection is
            disabled, this will be a 1-D array of values (of the same length as
            the :attr:`timestamps` attribute) for numerical data, and a
            :class:`CategoricalData` object for categorical data.

        Raises
        ------
        ValueError
            If select=True and extract=False, as select requires interpolation
        KeyError
            If sensor name was not found in cache and did not match virtual template

        """
        if select and not extract:
            raise ValueError('Cannot apply selection on raw sensor data')
        with self._lock:
            try:
                # First try to load the actual sensor data from cache
                sensor_data = self._raw[name]
            except KeyError:
                # Otherwise, iterate through virtual sensor templates and look for a match
                for pattern, create_sensor in self.virtual.items():
                    # Expand variable names enclosed in braces to the relevant regular expression
                    # (match anything but slashes, which are preferred delimiters in virtual sensor names)
                    pattern = re.sub(r'(\{[a-zA-Z_]\w*\})',
                                     lambda m: '(?P<{}>[^/]+)'.format(m.group(0)[1:-1]), pattern)
                    match = re.match(pattern, name)
                    if match:
                        # Call sensor creation function with extracted variables from sensor name
                        sensor_data = create_sensor(self, name, **match.groupdict())
                        break
                else:
                    if self.store:
                        # Katstore samples sensors at least once every 10 minutes
                        # Go that far back to support sporadic discrete sensors
                        start_time = self.timestamps[0] - self.dump_period - 600
                        end_time = self.timestamps[-1] + self.dump_period + 60
                        sensor_data = get_sensor_from_katstore(
                            self.store, name, start_time, end_time)
                    else:
                        raise KeyError(f"Unknown sensor '{name}' (does not match actual name or "
                                       "virtual template and no sensor store provided)")
            # If this is the first time this sensor is accessed, extract its data and store it in cache, if enabled
            if isinstance(sensor_data, SensorGetter) and extract:
                props = self._get_props(name, self.props, **kwargs)
                # If this is the first time any sensor is accessed, obtain all data timestamps via indexer
                self.timestamps = self.timestamps[:] if not isinstance(self.timestamps, np.ndarray) else self.timestamps
                sensor_data = self._extract(sensor_data, self.timestamps, self.dump_period, **props)
                self._raw[name] = sensor_data
        return sensor_data[self.keep] if select else sensor_data

    def get_with_fallback(self, sensor_type, names):
        """Sensor values interpolated to correlator data timestamps.

        Get data for a type of sensor that may have one of several names.
        Try each name in turn until something works, or crash sensibly.

        Parameters
        ----------
        sensor_type : string
            Name of sensor class / type, used for informational purposes only
        names : sequence of strings
            Sensor names to try until one of them provides data

        Returns
        -------
        sensor_data : array
           Interpolated sensor data as 1-D array, one value per selected timestamp

        Raises
        ------
        KeyError
            If none of the sensor names were found in the cache

        """
        for name in names:
            try:
                return self.get(name, select=True)
            except KeyError:
                logger.debug('Could not find %s sensor with name %r, trying next option', sensor_type, name)
        raise KeyError(f'Could not find any {sensor_type} sensor, tried {names}')

    # MutableMapping abstract methods

    def __setitem__(self, key, item):
        with self._lock:
            self._raw[key] = item

    def __delitem__(self, key):
        with self._lock:
            del self._raw[key]

    def __iter__(self):
        return iter(self._raw)

    def __len__(self):
        return len(self._raw)

    def __contains__(self, key):
        # __contains__ is implemented by MutableMapping via __getitem__, but
        # that does unnecessary extraction. This approach is cheaper but only
        # reflects keys that have been explicitly created or cached.
        with self._lock:
            return key in self._raw def to_str(value):
    """Convert string-likes to the native string type.

    Bytes are decoded to str, with surrogateencoding error handler.

    Tuples, lists, dicts and numpy arrays are processed recursively, with the
    exception that numpy structured types with string or object fields won't
    be handled.
    """
    if isinstance(value, np.ndarray) and value.dtype.kind == 'S':
        return np.char.decode(value, 'utf-8', 'surrogateescape')
    elif isinstance(value, bytes):
        return value.decode('utf-8', 'surrogateescape')

    # We use type(value) so that subclasses are reconstructed correctly
    if isinstance(value, (list, tuple)):
        return type(value)(to_str(item) for item in value)
    elif isinstance(value, dict):
        return type(value)((to_str(key), to_str(val)) for key, val in value.items())
    elif isinstance(value, np.ndarray) and value.dtype == 'O':
        return np.vectorize(to_str, otypes='O')(value)
    else:
        return value class SpectralWindow:
    """Spectral window specification.

    A spectral window is determined by the number of frequency channels produced
    by the correlator and their corresponding centre frequencies, as well as the
    channel width. The channels are assumed to be regularly spaced and to be the
    result of either lower-sideband downconversion (channel frequencies
    decreasing with channel index) or upper-sideband downconversion (frequencies
    increasing with index). For further information the receiver band and
    correlator product names are also available.

    .. warning::

        Instances should be treated as immutable. Changing the attributes will
        lead to inconsistencies between them.

    Parameters
    ----------
    centre_freq : float
        Centre frequency of spectral window, in Hz
    channel_width : float
        Bandwidth of each frequency channel, in Hz
    num_chans : int
        Number of frequency channels
    product : string, optional
        Name of data product / correlator mode
    sideband : {-1, +1}, optional
        Type of downconversion (-1 => lower sideband, +1 => upper sideband)
    band : {'L', 'UHF', 'S', 'X', 'Ku'}, optional
        Name of receiver / band
    bandwidth : float, optional
        The bandwidth of the whole spectral window, in Hz. If specified,
        `channel_width` is ignored and computed from the bandwidth. If not
        specified, bandwidth is computed from the channel width. Specifying
        this is a good idea if the channel width cannot be exactly represented
        in floating point.

    Attributes
    ----------
    channel_freqs : array of float, shape (*F*,)
        Centre frequency of each frequency channel (assuming LSB mixing), in Hz
    """

    def __init__(self, centre_freq, channel_width, num_chans, product=None,
                 sideband=-1, band='L', bandwidth=None):
        if bandwidth is None:
            bandwidth = channel_width * num_chans
        else:
            channel_width = bandwidth / num_chans
        self.centre_freq = centre_freq
        self.channel_width = channel_width
        self.bandwidth = bandwidth
        self.num_chans = num_chans
        self.product = product if product is not None else ''
        self.sideband = sideband
        self.band = band
        # channel_freqs is computed on demand
        self._channel_freqs_lock = threading.Lock()
        self._channel_freqs = None

    @property
    def channel_freqs(self):
        with self._channel_freqs_lock:
            if self._channel_freqs is None:
                # Don't subtract half a channel width as channel 0 is centred on 0 Hz in baseband
                # We use self.bandwidth and self.num_chans to avoid rounding
                # errors that might accumulate if channel_width is inexact.
                self._channel_freqs = self.centre_freq + self.sideband * self.bandwidth * (
                    np.arange(self.num_chans) - self.num_chans // 2) / self.num_chans
            return self._channel_freqs

    def __repr__(self):
        """Short human-friendly string representation of spectral window object."""
        band = self.band if self.band else 'unknown',
        product = repr(self.product) if self.product else 'unknown'
        return (f"<katdal.SpectralWindow {band}-band product={product} "
                f"centre={self.centre_freq/1e6:.3f} MHz bandwidth={self.bandwidth/1e6:.3f} MHz "
                f"channels={self.num_chans} at {id(self):#x}>")

    @property
    def _description(self):
        """Complete hashable representation, used internally for comparisons."""
        # Pick values that enable a sensible ordering of spectral windows
        # Using self.bandwidth is generally redundant but may play a role in
        # obscure rounding cases.
        return (self.centre_freq,
                -self.channel_width, self.num_chans, self.sideband,
                self.band, self.product, -self.bandwidth)

    def __eq__(self, other):
        """Equality comparison operator."""
        return self._description == (
            other._description if isinstance(other, SpectralWindow) else other)

    def __ne__(self, other):
        """Inequality comparison operator."""
        return not (self == other)

    def __lt__(self, other):
        """Less-than comparison operator (needed for sorting and np.unique)."""
        return self._description < (
            other._description if isinstance(other, SpectralWindow) else other)

    def __hash__(self):
        """Base hash on description tuple, just like equality operator."""
        return hash(self._description)

    def subrange(self, first, last):
        """Get a new :class:`SpectralWindow` representing a subset of the channels.

        The returned :class:`SpectralWindow` covers the same frequencies as
        channels [first, last) of the original.

        Raises
        ------
        IndexError
            If [first, last) is not a (non-empty) subinterval of the channels
        """
        if not (0 <= first < last <= self.num_chans):
            raise IndexError('channel indices out of range')
        channel_shift = (first + last) // 2 - self.num_chans // 2
        num_chans = last - first
        # We use self.bandwidth and self.num_chans to avoid rounding errors
        # that might accumulate if channel_width is inexact.
        centre_freq = self.centre_freq \
            + channel_shift * self.bandwidth * self.sideband / self.num_chans
        return SpectralWindow(
            centre_freq, self.channel_width, num_chans,
            self.product, self.sideband, self.band,
            bandwidth=self.bandwidth * num_chans / self.num_chans)

    def rechannelise(self, num_chans):
        """Get a new :class:`SpectralWindow` with a different number of channels.

        The returned :class:`SpectralWindow` covers the same frequencies as the
        original, but dividing the bandwidth into a different number of
        channels.
        """
        if num_chans == self.num_chans:
            return self
        # Find the centre of the bandwidth (whereas centre_freq is the centre
        # of the middle channel)
        centre_freq = self.centre_freq
        if self.num_chans % 2 == 0:
            centre_freq -= self.sideband * 0.5 * self.channel_width
        channel_width = self.bandwidth / num_chans
        # Now convert to the centre of the new middle channel
        if num_chans % 2 == 0:
            centre_freq += self.sideband * 0.5 * channel_width
        return SpectralWindow(
            centre_freq, channel_width, num_chans,
            self.product, self.sideband, self.band,
            bandwidth=self.bandwidth)
import logging
import pathlib
import secrets
import h5py
import katpoint
import numpy as np
from .categorical import CategoricalData, sensor_to_categorical
from .dataset import (DEFAULT_SENSOR_PROPS, DEFAULT_VIRTUAL_SENSORS,
                      BrokenFile, DataSet, Subarray, WrongVersion,
                      _robust_target, _selection_to_list)
from .flags import DESCRIPTIONS as FLAG_DESCRIPTIONS
from .flags import NAMES as FLAG_NAMES
from .lazy_indexer import LazyIndexer, LazyTransform
from .sensordata import RecordSensorGetter, SensorCache, to_str
from .spectral_window import SpectralWindow
        irregular = abs(expected_dumps - num_dumps) >= 0.01
        if irregular:
            # Warn the user, as this is anomalous
            logger.warning("Irregular timestamps detected in file '%s': expected %.3f dumps "
                           "based on dump period and start/end times, got %d instead",
                           filename, expected_dumps, num_dumps)
            if quicklook:
                logger.warning("Quicklook option selected - partitioning data based on synthesised timestamps instead")
        if not irregular or quicklook:
            # Estimate timestamps by assuming they are uniformly spaced (much quicker than loading them from file).
            # This is useful for the purpose of segmenting data set, where accurate timestamps are not that crucial.
            # The real timestamps are still loaded when the user explicitly asks for them.
            data_timestamps = self._timestamps[0] + self.dump_period * np.arange(num_dumps)
        else:
            # Load the real timestamps instead (could take several seconds on a large data set)
            data_timestamps = self._timestamps[:num_dumps]
        # Move timestamps from start of each dump to the middle of the dump
        data_timestamps += 0.5 * self.dump_period + self.time_offset
        if data_timestamps[0] < 1e9:
            logger.warning("File '%s' has invalid first correlator timestamp (%f)", filename, data_timestamps[0])
        self._time_keep = np.ones(num_dumps, dtype=np.bool)
        self.start_time = katpoint.Timestamp(data_timestamps[0] - 0.5 * self.dump_period)
        self.end_time = katpoint.Timestamp(data_timestamps[-1] + 0.5 * self.dump_period)
        self._keepdims = keepdims

        # ------ Extract flags ------

        # Check if flag group is present, else use dummy flag data
        self._flags = markup_group['flags'] if 'flags' in markup_group else \
            dummy_dataset('dummy_flags', shape=self._vis.shape[:-1], dtype=np.uint8, value=0)
        # Obtain flag descriptions from file or recreate default flag description table
        self._flags_description = to_str(markup_group['flags_description'][:]) \
            if 'flags_description' in markup_group else np.array(list(zip(FLAG_NAMES, FLAG_DESCRIPTIONS)))
        self._flags_select = np.array([0], dtype=np.uint8)
        self._flags_keep = 'all'

        # ------ Extract weights ------

        # Check if weight group present, else use dummy weight data
        self._weights = markup_group['weights'] if 'weights' in markup_group else \
            dummy_dataset('dummy_weights', shape=self._vis.shape[:-1], dtype=np.float32, value=1.0)
        # Obtain weight descriptions from file or recreate default weight description table
        self._weights_description = to_str(markup_group['weights_description'][:]) \
            if 'weights_description' in markup_group else np.array(list(zip(WEIGHT_NAMES, WEIGHT_DESCRIPTIONS)))
        self._weights_select = []
        self._weights_keep = 'all'

        # ------ Extract sensors ------

        # Populate sensor cache with all HDF5 datasets below sensor group that fit the description of a sensor
        cache = {}

        def register_sensor(name, obj):
            """A sensor is defined as a non-empty dataset with expected dtype."""
            if isinstance(obj, h5py.Dataset) and obj.shape != () and \
               obj.dtype.names == ('timestamp', 'value', 'status'):
                # Rename pedestal sensors from the old regime to become sensors of the corresponding antenna
                name = ('Antennas/ant' + name[13:]) if name.startswith('Pedestals/ped') else name
                cache[name] = RecordSensorGetter(obj, name)
        sensors_group.visititems(register_sensor)
        # Use estimated data timestamps for now, to speed up data segmentation
        self.sensor = SensorCache(cache, data_timestamps, self.dump_period, keep=self._time_keep,
                                  props=SENSOR_PROPS, virtual=VIRTUAL_SENSORS, aliases=SENSOR_ALIASES)

        # ------ Extract subarrays ------

        # By default, only pick antennas that were in use by the script
        script_ants = to_str(config_group['Observation'].attrs['script_ants']).split(',')
        self.ref_ant = script_ants[0] if not ref_ant else ref_ant
        # Original list of correlation products as pairs of input labels
        corrprods = get_single_value(config_group['Correlator'], 'bls_ordering')
        if len(corrprods) != self._vis.shape[2]:
            # Apply k7_capture baseline mask after the fact, in the hope that it fixes correlation product mislabelling
            corrprods = np.array([cp for cp in corrprods if cp[0][:-1] in script_ants and cp[1][:-1] in script_ants])
            # If there is still a mismatch between labels and data shape, file is considered broken (maybe bad labels?)
            if len(corrprods) != self._vis.shape[2]:
                raise BrokenFile('Number of baseline labels (containing expected antenna names) '
                                 'received from correlator (%d) differs from number of baselines in data (%d)' %
                                 (len(corrprods), self._vis.shape[2]))
            else:
                logger.warning('Reapplied k7_capture baseline mask to fix unexpected number of baseline labels')
        # All antennas in configuration as katpoint Antenna objects
        ants = [katpoint.Antenna(to_str(config_group['Antennas'][name].attrs['description']))
                for name in config_group['Antennas']]
        self.subarrays = [Subarray(ants, corrprods)]
        self.sensor['Observation/subarray'] = CategoricalData(self.subarrays, [0, len(data_timestamps)])
        self.sensor['Observation/subarray_index'] = CategoricalData([0], [0, len(data_timestamps)])
        # Store antenna objects in sensor cache too, for use in virtual sensor calculations
        for ant in ants:
            self.sensor[f'Antennas/{ant.name}/antenna'] = CategoricalData([ant], [0, len(data_timestamps)])
        # Extract array reference from first antenna (first 5 fields of description)
        array_ant_fields = ['array'] + ants[0].description.split(',')[1:5]
        array_ant = katpoint.Antenna(','.join(array_ant_fields))
        self.sensor['Antennas/array/antenna'] = CategoricalData([array_ant], [0, len(data_timestamps)])

        # ------ Extract spectral windows / frequencies ------

        # Ideally we would like to use calculated center-frequency-hz sensor produced by k7_capture (better for nband)
        if self.version >= '2.1':
            centre_freq = self.sensor.get('RFE/center-frequency-hz')
        else:
            # Fall back to basic RFE7 LO frequency, as this supported multiple spectral windows before k7_capture did
            # This assumes WBC mode, though (NBC modes only fully supported since HDF5 v2.1)
            centre_freq = self.sensor.get('RFE/rfe7.lo1.frequency')
            centre_freq.unique_values = [freq - 4200e6 for freq in centre_freq.unique_values]
        num_chans = get_single_value(config_group['Correlator'], 'n_chans')
        if num_chans != self._vis.shape[1]:
            raise BrokenFile(f'Number of channels received from correlator ({num_chans}) '
                             f'differs from number of channels in data ({self._vis.shape[1]})')
        bandwidth = get_single_value(config_group['Correlator'], 'bandwidth')
        channel_width = bandwidth / num_chans
        try:
            mode = self.sensor.get('DBE/dbe.mode').unique_values[0]
        except (KeyError, IndexError):
            # Guess the mode for version 2.0 files that haven't been re-augmented
            mode = 'wbc' if num_chans <= 1024 else 'wbc8k' if bandwidth > 200e6 else 'nbc'
        self.spectral_windows = [SpectralWindow(spw_centre, channel_width, num_chans, mode)
                                 for spw_centre in centre_freq.unique_values]
        self.sensor['Observation/spw'] = CategoricalData([self.spectral_windows[idx] for idx in centre_freq.indices],
                                                         centre_freq.events)

        self.sensor['Observation/spw_index'] = CategoricalData(centre_freq.indices, centre_freq.events)