From b4cd2c8e86675015dee72708e761dc6ec3f9a797 Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Fri, 13 Mar 2020 18:08:00 -0400 Subject: [PATCH 01/19] Add code to support and test computed MSIDs This includes the cleaned MUPS valve temperature computed MSIDs. --- Ska/engarchive/derived/base.py | 4 +- Ska/engarchive/derived/comps.py | 51 ++++++ Ska/engarchive/derived/mups_valve.py | 250 +++++++++++++++++++++++++++ Ska/engarchive/fetch.py | 37 ++++ Ska/engarchive/tests/test_comps.py | 132 ++++++++++++++ 5 files changed, 473 insertions(+), 1 deletion(-) create mode 100644 Ska/engarchive/derived/comps.py create mode 100644 Ska/engarchive/derived/mups_valve.py create mode 100644 Ska/engarchive/tests/test_comps.py diff --git a/Ska/engarchive/derived/base.py b/Ska/engarchive/derived/base.py index 3f65c566..c523d0de 100644 --- a/Ska/engarchive/derived/base.py +++ b/Ska/engarchive/derived/base.py @@ -1,7 +1,6 @@ # Licensed under a 3-clause BSD style license - see LICENSE.rst from Chandra.Time import DateTime -from .. import fetch, fetch_eng import Ska.Numpy import numpy as np from .. import cache @@ -35,6 +34,8 @@ def calc(self, data): raise NotImplementedError def fetch(self, start, stop): + from .. import fetch + unit_system = fetch.get_units() # cache current units and restore after fetch fetch.set_units(self.unit_system) dataset = fetch.MSIDset(self.rootparams, start, stop) @@ -89,6 +90,7 @@ def fetch(self, start, stop): return dataset def __call__(self, start, stop): + from .. import fetch_eng dataset = fetch_eng.MSIDset(self.rootparams, start, stop, filter_bad=True) # Translate state codes "ON" and "OFF" to 1 and 0, respectively. diff --git a/Ska/engarchive/derived/comps.py b/Ska/engarchive/derived/comps.py new file mode 100644 index 00000000..109448e7 --- /dev/null +++ b/Ska/engarchive/derived/comps.py @@ -0,0 +1,51 @@ +class ComputedMsid: + # Global dict of registered computed MSIDs + msid_classes = {} + + # Default MSID attributes that are provided + msid_attrs = ('times', 'vals', 'bads') + + def __init_subclass__(cls, **kwargs): + super().__init_subclass__(**kwargs) + if not cls.__name__.startswith('Comp_'): + raise ValueError(f'comp class name {cls.__name__} must start with Comp_') + cls.msid_classes[cls.__name__.upper()] = cls + + @property + def fetch_eng(self): + from .. import fetch_eng + return fetch_eng + + @property + def fetch_sci(self): + from .. import fetch_sci + return fetch_sci + + @property + def fetch_cxc(self): + from .. import fetch_cxc + return fetch_cxc + + def __call__(self, start, stop): + dat = self.get_MSID(start, stop) + out = {attr: getattr(dat, attr) for attr in self.msid_attrs} + return out + + +class Comp_MUPS_Clean: + msid_attrs = ComputedMsid.msid_attrs + ('vals_raw', 'vals_nan', 'vals_corr', + 'vals_model', 'source') + + def get_MSID(self, start, stop): + from .mups_valve import fetch_clean_msid + dat = fetch_clean_msid(self.msid.lower(), start, stop, + dt_thresh=5.0, median=7, model_spec=None) + return dat + + +class Comp_PM2THV1T(Comp_MUPS_Clean, ComputedMsid): + msid = 'pm2thv1t' + + +class Comp_PM1THV2T(Comp_MUPS_Clean, ComputedMsid): + msid = 'pm1thv2t' diff --git a/Ska/engarchive/derived/mups_valve.py b/Ska/engarchive/derived/mups_valve.py new file mode 100644 index 00000000..51602295 --- /dev/null +++ b/Ska/engarchive/derived/mups_valve.py @@ -0,0 +1,250 @@ +# coding: utf-8 + +""" + Fetch clean MUPS valve temperature telemetry + + This makes use of the temperature correction code provided by Scott Blanchard (to correct + for a resistor dropout in the thermistor data) and xija thermal model developed by + Matt Dahmer. + + The basic cleaning algorithm is very simple: + + - Fetch raw telemetry from the cheta archive. + - Compute a xija thermal model prediction for the same timespan (actually starting a few + days in advance to burn out uncertainty in the pseudo-node value). + - Accept either raw telemetry or corrected data which are within a tolerance (5 degF) of + the model. + - In the gaps where the model diverges from both raw and corrected temperatures, + "repropagate" the model starting from the last accepted temperature value. + This effectively takes out much of the systematic model error, which can be up to + 15 degF after a transition to a new attitude. + - In some cases this allows recovery of additional data, while in others the data are + not recoverable: either due to partial disconnect of the parallel resistor or full + disconnects where the output voltage exceeds 5.12 V. + + The output of this function is a `fetch.Msid` object with some bonus attributes, + documented below. In particular the output cleaned data are labeled so one knows + exactly where each data point came from. + + The function is fast, and one can get 5 years of cleaned telemetry in a few seconds + on a modern laptop with SSD drive. + + This cleaning technique recovers on average about 90% of data for PM2THV1T. Since + 2015, about 60% of telemetry is good (no dropout) while 30% is in a recoverable + fully-dropped state (and 10% is not recoverable). + + ``` + def fetch_clean_msid(msid, start, stop=None, dt_thresh=5.0, median=7, model_spec=None): + + Fetch a cleaned version of telemetry for ``msid``. + + If not supplied the model spec will be downloaded from github using this URL: + https://raw.githubusercontent.com/sot/chandra_models/master/chandra_models/xija/ + mups_valve/{msid}_spec.json + + This function returns a `fetch.Msid` object like a normal fetch but with extra attributes: + + - vals: cleaned telemetry (either original or corrected telemetry, or xija model prediction) + - source: label for each vals data point + - 0: unrecoverable, so use xija model value + - 1: original telemetry + - 2: corrected telemetry + - vals_raw: raw (uncleaned) telemetry + - vals_nan: cleaned telem but with np.nan at points where data are unrecoverable (this is + for plotting) + - vals_corr: telemetry with the MUPS correction applied + - vals_model: xija model prediction + + :param start: start time + :param stop: stop time (default=NOW) + :param dt_thresh: tolerance for matching model to data in degF (default=5 degF) + :param median: length of median filter (default=7, use 0 to disable) + :param model_spec: file name or URL containing relevant xija model spec + + :returns: fetch.Msid object + ``` +""" + +import os + +import numpy as np +import numba +from astropy.utils.data import download_file +from scipy.interpolate import interp1d +from scipy.ndimage import median_filter + +from cheta import fetch_eng +from xija import XijaModel + + +def c2f(degc): + degf = 32 + degc * 1.8 + return degf + + +# Define MUPS valve thermistor point pair calibration table (from TDB). +pp_counts = np.array([0, 27, 36, 44, 55, 70, 90, 118, 175, 195, 210, 219, 226, 231, 235, 255]) +pp_temps = np.array([369.5305, 263.32577, 239.03652, 222.30608, 203.6944, 183.2642, 161.0796, + 134.93818, 85.65725, 65.6537, 47.3176, 33.50622, 19.9373, 7.42435, -5.79635, + -111.77265]) + +count_to_degf = interp1d(pp_counts, pp_temps) +degf_to_counts = interp1d(pp_temps, pp_counts) + +# Define MUPS valve thermistor voltage point pair calibration table. +count_to_volts = lambda counts: counts / 256 * 5.12 +volts_to_counts = lambda volts: volts / 5.12 * 256 + +# Voltage and Temperature, with and without resistor (see flight note 447 for +# the source of these numbers). + +volt_with_resistor = [4.153325779, 3.676396578, 3.175100371, 2.587948965, 2.435, 2.025223702, + 1.538506813, 1.148359251, 0.63128179, 0.354868907, 0.208375569] +volt_without_resistor = [28.223, 15, 9.1231, 5.5228, 4.87, 3.467, 2.249, + 1.5027, 0.7253, 0.38276, 0.21769] +temp_without_resistor = [50, 77, 100, 125, 130, 150, 175, 200, 250, 300, 350] + +volt_without_resistor_to_temp_without_resistor = interp1d(volt_without_resistor, + temp_without_resistor) +temp_without_resistor_to_volt_without_resistor = interp1d(temp_without_resistor, + volt_without_resistor) +volt_with_resistor_to_volt_without_resistor = interp1d(volt_with_resistor, volt_without_resistor) +volt_without_resistor_to_volt_with_resistor = interp1d(volt_without_resistor, volt_with_resistor) + + +def get_corr_mups_temp(temp): + """ Calculate a MUPS valve thermistor corrected temperature. + Args: + temp (float, int): Temperature in Fahrenheit to which a correction will be applied. + + Returns: + (float): Corrected temperaure + + """ + # Convert observed temperature to voltage + count = degf_to_counts(temp) + volt = count_to_volts(count) + + # Convert voltage as read, assuming no resistor, to what it would be with the resistor + new_volt = volt_without_resistor_to_volt_with_resistor(volt) + + # Convert this new voltage to counts and, in turn, to a new temperature + new_count = volts_to_counts(new_volt) + new_temp = count_to_degf(new_count) + + return new_temp + + +@numba.autojit() +def select_using_model(data1, data2, model, dt_thresh, out, source): + """Select from either ``data1`` or ``data2`` using ``model`` to get clean data. + + This does a 2-step process: + + 1. Select either `data1` or `data2` where either is within `dt_thresh` of `model`. + 2. Fill in gaps where possible by propagating the model from previous good data, + thus reducing the impact of model errors. + + :param data1: ndarray of data (e.g. observed temperature telemetry) + :param data2: ndarray of data (e.g. corrected temperature telemetry) + :param model: ndarray with model prediction for data + :param dt_thresh: model error tolerance for selecting data1 or data2 + :param out: ndarray where output cleaned data is stored (must be same length as data) + :param source: ndarray where output source labeling is stored (same length as data) + """ + # Select all data where either data1 or data2 are within dt_thresh of model + for data, source_id in [(data1, 1), (data2, 2)]: + ok = np.abs(data - model) < dt_thresh + out[ok] = data[ok] + source[ok] = source_id + + # Fill in gaps where possible by propagating model from last good data value + for ii in range(len(data1)): + if source[ii] != 0: + # Already got a good value so move on. + continue + + if ii == 0: + model_prop_ii = model[0] + else: + # Propagate model using last point + delta-model + model_prop_ii = out[ii-1] + (model[ii] - model[ii-1]) # noqa + + if abs(data1[ii] - model_prop_ii) < dt_thresh: + out[ii] = data1[ii] + source[ii] = 1 + elif abs(data2[ii] - model_prop_ii) < dt_thresh: + out[ii] = data2[ii] + source[ii] = 2 + else: + # No good match, use propagated model value and mark as such + out[ii] = model_prop_ii + source[ii] = 0 + + +def fetch_clean_msid(msid, start, stop=None, dt_thresh=5.0, median=7, model_spec=None): + """Fetch a cleaned version of telemetry for ``msid``. + + If not supplied the model spec will be downloaded from github using this URL: + https://raw.githubusercontent.com/sot/chandra_models/master/chandra_models/xija/ + mups_valve/{msid}_spec.json + + This function returns a `fetch.Msid` object like a normal fetch but with extra attributes: + + - vals: cleaned telemetry (either original or corrected telemetry, or xija model prediction) + - source: label for each vals data point + - 0: unrecoverable, so use xija model value + - 1: original telemetry + - 2: corrected telemetry + - vals_raw: raw (uncleaned) telemetry + - vals_nan: cleaned telem but with np.nan at points where data are unrecoverable (this is + for plotting) + - vals_corr: telemetry with the MUPS correction applied + - vals_model: xija model prediction + + :param start: start time + :param stop: stop time (default=NOW) + :param dt_thresh: tolerance for matching model to data in degF (default=5 degF) + :param median: length of median filter (default=7, use 0 to disable) + :param model_spec: file name or URL containing relevant xija model spec + + :returns: fetch.Msid object + """ + if model_spec is None or not os.path.exists(model_spec): + url = (f'https://raw.githubusercontent.com/sot/chandra_models/master/chandra_models/xija/' + f'mups_valve/{msid}_spec.json') + model_spec = download_file(url, cache=True) + + dat = fetch_eng.Msid(msid, start, stop) + + t_obs = dat.vals + if median: + t_obs = median_filter(t_obs, size=median) + t_corr = get_corr_mups_temp(t_obs) + + mdl = XijaModel(model_spec=model_spec, + start=dat.times[0] - 400000, + stop=dat.times[-1]) + mdl.comp['mups0'].set_data(75) + mdl.make() + mdl.calc() + + t_model = np.interp(xp=mdl.times, fp=c2f(mdl.comp[msid].mvals), x=dat.times) + + out = np.zeros_like(dat.vals) + source = np.zeros(len(dat.vals), dtype=int) + select_using_model(t_obs, t_corr, t_model, dt_thresh=dt_thresh, out=out, source=source) + + dat.vals_raw = dat.vals + + dat.vals = out + dat.vals_nan = out.copy() + + dat.bads = (source == 0) + dat.vals_nan[dat.bads] = np.nan + + dat.source = source + dat.vals_corr = t_corr + dat.vals_model = t_model + + return dat diff --git a/Ska/engarchive/fetch.py b/Ska/engarchive/fetch.py index 6d5004ca..9c036703 100644 --- a/Ska/engarchive/fetch.py +++ b/Ska/engarchive/fetch.py @@ -29,6 +29,7 @@ from . import cache from . import remote_access from .remote_access import ENG_ARCHIVE +from .derived.comps import ComputedMsid from . import __version__ # noqa from Chandra.Time import DateTime @@ -150,6 +151,8 @@ def get_msids(cls, source): elif source == 'maude': import maude out = list(maude.MSIDS.keys()) + elif source == 'comp': + out = list(key for key in content if key.startswith('COMP_')) else: raise ValueError('source must be "cxc" or "msid"') @@ -295,6 +298,9 @@ def load_msid_names(all_msid_names_files): # Load the MSID names all_colnames = load_msid_names(all_msid_names_files) +# Add computed MSIDs +all_colnames['computed'] = list(ComputedMsid.msid_classes) + # Save the names for k, colnames in six.iteritems(all_colnames): content.update((x, k) for x in sorted(colnames) @@ -401,6 +407,8 @@ def _msid_glob(msid, source): :param msid: input MSID glob :returns: tuple (msids, MSIDs) """ + if msid.lower().startswith('comp_'): + source = 'comp' source_msids = data_source.get_msids(source) @@ -578,6 +586,12 @@ def _get_data(self): logger.info('Getting data for %s between %s to %s', self.msid, self.datestart, self.datestop) + if self.content == 'computed': + if self.stat: + raise ValueError('stats are not supported for computed MSIDs') + self._get_comp_data() + return + # Avoid stomping on caller's filetype 'ft' values with _cache_ft() with _cache_ft(): ft['content'] = self.content @@ -620,6 +634,29 @@ def _get_data(self): # telemetry to existing CXC values. self._get_msid_data_from_maude(*args) + def _get_comp_data(self): + logger.info(f'Getting comp data for {self.msid}') + + # Get the ComputedMsid subclass corresponding to MSID and do the + # computation. This returns a dict of MSID attribute values. + comp_cls = ComputedMsid.msid_classes[self.MSID] + dat = comp_cls()(self.tstart, self.tstop) + + # Allow upstream class to be a bit sloppy on times and include samples + # outside the time range. This can happen with classes that inherit + # from DerivedParameter. + ok = (dat['times'] >= self.tstart) & (dat['times'] <= self.tstop) + all_ok = np.all(ok) + + # Apply attributes to self + self.colnames = list(dat) + for attr, val in dat.items(): + if (not all_ok + and isinstance(val, np.ndarray) + and len(val) == len(dat['times'])): + val = val[ok] + setattr(self, attr, val) + def _get_stat_data(self): """Do the actual work of getting stats values for an MSID from HDF5 files""" diff --git a/Ska/engarchive/tests/test_comps.py b/Ska/engarchive/tests/test_comps.py new file mode 100644 index 00000000..c5502486 --- /dev/null +++ b/Ska/engarchive/tests/test_comps.py @@ -0,0 +1,132 @@ +# Licensed under a 3-clause BSD style license - see LICENSE.rst + +"""Test that computed MSIDs work as expected.""" + +import numpy as np +import pytest + +from ..derived.base import DerivedParameter +from ..derived.comps import ComputedMsid + +try: + import maude + date1 = '2016:001:00:00:00.1' + date2 = '2016:001:00:00:02.0' + maude.get_msids(msids='ccsdsid', start=date1, stop=date2) +except Exception: + HAS_MAUDE = False +else: + HAS_MAUDE = True + + +class Comp_Plus_Five: + """Silly base comp to add 5 to the value""" + def get_MSID(self, start, stop): + dat = self.fetch_eng.MSID(self.msid, start, stop) + dat.vals = dat.vals + 5 + return dat + + +class Comp_TEPHIN_Plus_Five(Comp_Plus_Five, ComputedMsid): + msid = 'tephin' + + +class Comp_CSS1_NPM_SUN(DerivedParameter, ComputedMsid): + """Coarse Sun Sensor Counts 1 filtered for NPM and SA Illuminated + + Defined as CSS-1 current converted back into counts + (AOCSSI1 * 4095 / 5.49549) when in NPM (AOPCADMD==1) and SA is illuminated + (AOSAILLM==1). Otherwise, "Bads" flag is set equal to one. + + """ + rootparams = ['aocssi1', 'aopcadmd', 'aosaillm'] + time_step = 1.025 + max_gap = 10.0 + + def __call__(self, start, stop): + # Get an interpolated MSIDset for rootparams + msids = self.fetch(start, stop) + + # Do the computation and set bad values + npm_sun = ((msids['aopcadmd'].vals == 'NPNT') & + (msids['aosaillm'].vals == 'ILLM')) + msids.bads = msids.bads | ~npm_sun + css1_npm_sun = msids['aocssi1'].vals * 4095 / 5.49549 + + out = {'vals': css1_npm_sun, + 'times': msids.times, + 'bads': msids.bads} + return out + + +# Need to add classes before importing fetch +from .. import fetch_eng, fetch + + +def test_comp_from_derived_parameter(): + """Test that on-the-fly comp gives same result as derived parameter from + same code. + """ + dat1 = fetch_eng.Msid('comp_css1_npm_sun', '2020:001', '2020:010') + dat2 = fetch_eng.Msid('dp_css1_npm_sun', '2020:001', '2020:010') + + for attr in ('vals', 'times', 'bads'): + assert np.all(getattr(dat1, attr) == getattr(dat2, attr)) + + +def test_simple_comp(): + dat1 = fetch_eng.Msid('tephin', '2020:001', '2020:010') + dat2 = fetch_eng.Msid('comp_tephin_plus_five', '2020:001', '2020:010') + assert np.all(dat1.vals + 5 == dat2.vals) + assert np.all(dat1.times == dat2.times) + assert np.all(dat1.bads == dat2.bads) + + +@pytest.mark.skipif("not HAS_MAUDE") +def test_simple_comp_with_maude(): + with fetch.data_source('maude'): + dat1 = fetch_eng.Msid('tephin', '2020:001', '2020:003') + dat2 = fetch_eng.Msid('comp_tephin_plus_five', '2020:001', '2020:003') + assert np.all(dat1.vals + 5 == dat2.vals) + assert np.all(dat1.times == dat2.times) + assert np.all(dat1.bads == dat2.bads) + + +def test_comp_with_stat(): + with pytest.raises(ValueError, + match='stats are not supported for computed MSIDs'): + fetch_eng.Msid('comp_tephin_plus_five', + '2020:001', '2020:010', stat='5min') + + +def test_mups_valve(): + colnames = ['times', 'vals', 'bads', 'vals_raw', + 'vals_nan', 'vals_corr', 'vals_model', 'source'] + + dat = fetch.MSID('comp_PM2THV1T', '2020:001', '2020:010') + assert len(dat.vals) == 36661 + assert np.count_nonzero(dat.source != 0) == 34499 + assert dat.colnames == colnames + for attr in colnames: + assert len(dat.vals) == len(getattr(dat, attr)) + + dat = fetch.Msid('comp_PM2THV1T', '2020:001', '2020:010') + assert len(dat.vals) == 34499 # Some bad values + assert dat.colnames == colnames + for attr in colnames: + if attr != 'bads': + assert len(dat.vals) == len(getattr(dat, attr)) + + dat = fetch.MSID('comp_PM1THV2T', '2020:001', '2020:010') + assert len(dat.vals) == 36661 # Same as PM2THV1T + assert dat.colnames == colnames + for attr in colnames: + assert len(dat.vals) == len(getattr(dat, attr)) + + dat = fetch.Msid('comp_PM1THV2T', '2020:001', '2020:010') + assert len(dat.vals) == 36240 # Some bad values + assert len(dat.source) == 36240 # Filtering applies to sources + assert dat.colnames == colnames + for attr in colnames: + if attr != 'bads': + assert len(dat.vals) == len(getattr(dat, attr)) From 95aab7ba92602bf2f44e5d861d24b23f7b980df2 Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Sun, 15 Mar 2020 10:55:01 -0400 Subject: [PATCH 02/19] Much-improved implementation and API Comps have a generalized MSID match regex to allow a comp that has arguments from the regex. Also check for comps dynamically instead of using pre-computed content hacks. This also removes issues with import order. --- Ska/engarchive/derived/comps.py | 76 ++++++++++++++++++++++-------- Ska/engarchive/fetch.py | 26 +++++----- Ska/engarchive/tests/test_comps.py | 25 +++++----- 3 files changed, 80 insertions(+), 47 deletions(-) diff --git a/Ska/engarchive/derived/comps.py b/Ska/engarchive/derived/comps.py index 109448e7..34f04d32 100644 --- a/Ska/engarchive/derived/comps.py +++ b/Ska/engarchive/derived/comps.py @@ -1,15 +1,41 @@ +import re + + class ComputedMsid: # Global dict of registered computed MSIDs - msid_classes = {} + msid_classes = [] - # Default MSID attributes that are provided + # Standard base MSID attributes that must be provided msid_attrs = ('times', 'vals', 'bads') + # Extra MSID attributes that are provided beyond times, vals, bads + extra_msid_attrs = () + def __init_subclass__(cls, **kwargs): super().__init_subclass__(**kwargs) + + # Validate class name, msid_attrs, msid_match if not cls.__name__.startswith('Comp_'): raise ValueError(f'comp class name {cls.__name__} must start with Comp_') - cls.msid_classes[cls.__name__.upper()] = cls + + cls.msid_attrs = ComputedMsid.msid_attrs + cls.extra_msid_attrs + + if not hasattr(cls, 'msid_match'): + raise ValueError(f'comp {cls.__name__} must define msid_match') + + # Force match to include entire line + cls.msid_match = cls.msid_match + '$' + + cls.msid_classes.append(cls) + + @classmethod + def get_matching_comp_cls(cls, msid): + for comp_cls in ComputedMsid.msid_classes: + match = re.match(comp_cls.msid_match + '$', msid, re.IGNORECASE) + if match: + return comp_cls + + return None @property def fetch_eng(self): @@ -26,26 +52,38 @@ def fetch_cxc(self): from .. import fetch_cxc return fetch_cxc - def __call__(self, start, stop): - dat = self.get_MSID(start, stop) - out = {attr: getattr(dat, attr) for attr in self.msid_attrs} - return out + def __call__(self, start, stop, msid): + match = re.match(self.msid_match, msid, re.IGNORECASE) + if not match: + raise RuntimeError(f'unexpected mismatch of {msid} with {self.msid_match}') + match_args = [arg.lower() for arg in match.groups()] + msid_attrs = self.get_msid_attrs(start, stop, msid.lower(), match_args) + if set(msid_attrs) != set(self.msid_attrs): + raise ValueError(f'computed class did not return expected attributes') -class Comp_MUPS_Clean: - msid_attrs = ComputedMsid.msid_attrs + ('vals_raw', 'vals_nan', 'vals_corr', - 'vals_model', 'source') + return msid_attrs - def get_MSID(self, start, stop): - from .mups_valve import fetch_clean_msid - dat = fetch_clean_msid(self.msid.lower(), start, stop, - dt_thresh=5.0, median=7, model_spec=None) - return dat + def get_msid_attrs(self, start, stop, msid, msid_args): + """Get the attributes required for this MSID. + + TODO: detailed docs here since this is the main user-defined method + """ + raise NotImplementedError() -class Comp_PM2THV1T(Comp_MUPS_Clean, ComputedMsid): - msid = 'pm2thv1t' +class Comp_MUPS_Valve_Temp_Clean(ComputedMsid): + msid_match = r'comp_(pm2thv1t|pm1thv2t)' + extra_msid_attrs = ('vals_raw', 'vals_nan', 'vals_corr', 'vals_model', 'source') + + def get_msid_attrs(self, start, stop, msid, msid_args): + from .mups_valve import fetch_clean_msid + + # Get cleaned MUPS valve temperature data as an MSID object + dat = fetch_clean_msid(msid_args[0], start, stop, + dt_thresh=5.0, median=7, model_spec=None) + # Convert to dict as required by the get_msids_attrs API + msid_attrs = {attr: getattr(dat, attr) for attr in self.msid_attrs} -class Comp_PM1THV2T(Comp_MUPS_Clean, ComputedMsid): - msid = 'pm1thv2t' + return msid_attrs diff --git a/Ska/engarchive/fetch.py b/Ska/engarchive/fetch.py index 9c036703..eb21b589 100644 --- a/Ska/engarchive/fetch.py +++ b/Ska/engarchive/fetch.py @@ -151,8 +151,6 @@ def get_msids(cls, source): elif source == 'maude': import maude out = list(maude.MSIDS.keys()) - elif source == 'comp': - out = list(key for key in content if key.startswith('COMP_')) else: raise ValueError('source must be "cxc" or "msid"') @@ -298,9 +296,6 @@ def load_msid_names(all_msid_names_files): # Load the MSID names all_colnames = load_msid_names(all_msid_names_files) -# Add computed MSIDs -all_colnames['computed'] = list(ComputedMsid.msid_classes) - # Save the names for k, colnames in six.iteritems(all_colnames): content.update((x, k) for x in sorted(colnames) @@ -383,6 +378,11 @@ def msid_glob(msid): msids = collections.OrderedDict() MSIDS = collections.OrderedDict() + # First check if `msid` matches a computed class. This does not allow + # for globs, and here the output MSIDs is the single computed class. + if ComputedMsid.get_matching_comp_cls(msid) is not None: + return [msid], [msid.upper()] + sources = data_source.sources(include_test=False) for source in sources: ms, MS = _msid_glob(msid, source) @@ -407,9 +407,6 @@ def _msid_glob(msid, source): :param msid: input MSID glob :returns: tuple (msids, MSIDs) """ - if msid.lower().startswith('comp_'): - source = 'comp' - source_msids = data_source.get_msids(source) MSID = msid.upper() @@ -586,10 +583,11 @@ def _get_data(self): logger.info('Getting data for %s between %s to %s', self.msid, self.datestart, self.datestop) - if self.content == 'computed': + comp_cls = ComputedMsid.get_matching_comp_cls(self.msid) + if comp_cls: if self.stat: raise ValueError('stats are not supported for computed MSIDs') - self._get_comp_data() + self._get_comp_data(comp_cls) return # Avoid stomping on caller's filetype 'ft' values with _cache_ft() @@ -634,13 +632,11 @@ def _get_data(self): # telemetry to existing CXC values. self._get_msid_data_from_maude(*args) - def _get_comp_data(self): + def _get_comp_data(self, comp_cls): logger.info(f'Getting comp data for {self.msid}') - # Get the ComputedMsid subclass corresponding to MSID and do the - # computation. This returns a dict of MSID attribute values. - comp_cls = ComputedMsid.msid_classes[self.MSID] - dat = comp_cls()(self.tstart, self.tstop) + # Do computation. This returns a dict of MSID attribute values. + dat = comp_cls()(self.tstart, self.tstop, self.msid) # Allow upstream class to be a bit sloppy on times and include samples # outside the time range. This can happen with classes that inherit diff --git a/Ska/engarchive/tests/test_comps.py b/Ska/engarchive/tests/test_comps.py index c5502486..c1752e66 100644 --- a/Ska/engarchive/tests/test_comps.py +++ b/Ska/engarchive/tests/test_comps.py @@ -5,6 +5,7 @@ import numpy as np import pytest +from .. import fetch_eng, fetch from ..derived.base import DerivedParameter from ..derived.comps import ComputedMsid @@ -19,19 +20,20 @@ HAS_MAUDE = True -class Comp_Plus_Five: +class Comp_Val_Plus_Five(ComputedMsid): """Silly base comp to add 5 to the value""" - def get_MSID(self, start, stop): - dat = self.fetch_eng.MSID(self.msid, start, stop) - dat.vals = dat.vals + 5 - return dat + msid_match = r'comp_(\w+)_plus_five' + def get_msid_attrs(self, start, stop, msid, msid_args): + dat = self.fetch_eng.MSID(msid_args[0], start, stop) -class Comp_TEPHIN_Plus_Five(Comp_Plus_Five, ComputedMsid): - msid = 'tephin' + out = {'vals': dat.vals + 5, + 'bads': dat.bads, + 'times': dat.times} + return out -class Comp_CSS1_NPM_SUN(DerivedParameter, ComputedMsid): +class Comp_CSS1_NPM_SUN(ComputedMsid, DerivedParameter): """Coarse Sun Sensor Counts 1 filtered for NPM and SA Illuminated Defined as CSS-1 current converted back into counts @@ -42,8 +44,9 @@ class Comp_CSS1_NPM_SUN(DerivedParameter, ComputedMsid): rootparams = ['aocssi1', 'aopcadmd', 'aosaillm'] time_step = 1.025 max_gap = 10.0 + msid_match = 'comp_css1_npm_sun' - def __call__(self, start, stop): + def get_msid_attrs(self, start, stop, msid, msid_args): # Get an interpolated MSIDset for rootparams msids = self.fetch(start, stop) @@ -59,10 +62,6 @@ def __call__(self, start, stop): return out -# Need to add classes before importing fetch -from .. import fetch_eng, fetch - - def test_comp_from_derived_parameter(): """Test that on-the-fly comp gives same result as derived parameter from same code. From 3d1b2c20e01d8c4705e5a34a031b67b00130fbe7 Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Sun, 15 Mar 2020 10:59:40 -0400 Subject: [PATCH 03/19] Change built-in MUPS valve temp comp name --- Ska/engarchive/derived/comps.py | 2 +- Ska/engarchive/tests/test_comps.py | 8 ++++---- 2 files changed, 5 insertions(+), 5 deletions(-) diff --git a/Ska/engarchive/derived/comps.py b/Ska/engarchive/derived/comps.py index 34f04d32..d19887af 100644 --- a/Ska/engarchive/derived/comps.py +++ b/Ska/engarchive/derived/comps.py @@ -73,7 +73,7 @@ def get_msid_attrs(self, start, stop, msid, msid_args): class Comp_MUPS_Valve_Temp_Clean(ComputedMsid): - msid_match = r'comp_(pm2thv1t|pm1thv2t)' + msid_match = r'(pm2thv1t|pm1thv2t)_clean' extra_msid_attrs = ('vals_raw', 'vals_nan', 'vals_corr', 'vals_model', 'source') def get_msid_attrs(self, start, stop, msid, msid_args): diff --git a/Ska/engarchive/tests/test_comps.py b/Ska/engarchive/tests/test_comps.py index c1752e66..cdaa2367 100644 --- a/Ska/engarchive/tests/test_comps.py +++ b/Ska/engarchive/tests/test_comps.py @@ -102,27 +102,27 @@ def test_mups_valve(): colnames = ['times', 'vals', 'bads', 'vals_raw', 'vals_nan', 'vals_corr', 'vals_model', 'source'] - dat = fetch.MSID('comp_PM2THV1T', '2020:001', '2020:010') + dat = fetch.MSID('PM2THV1T_clean', '2020:001', '2020:010') assert len(dat.vals) == 36661 assert np.count_nonzero(dat.source != 0) == 34499 assert dat.colnames == colnames for attr in colnames: assert len(dat.vals) == len(getattr(dat, attr)) - dat = fetch.Msid('comp_PM2THV1T', '2020:001', '2020:010') + dat = fetch.Msid('PM2THV1T_clean', '2020:001', '2020:010') assert len(dat.vals) == 34499 # Some bad values assert dat.colnames == colnames for attr in colnames: if attr != 'bads': assert len(dat.vals) == len(getattr(dat, attr)) - dat = fetch.MSID('comp_PM1THV2T', '2020:001', '2020:010') + dat = fetch.MSID('PM1THV2T_clean', '2020:001', '2020:010') assert len(dat.vals) == 36661 # Same as PM2THV1T assert dat.colnames == colnames for attr in colnames: assert len(dat.vals) == len(getattr(dat, attr)) - dat = fetch.Msid('comp_PM1THV2T', '2020:001', '2020:010') + dat = fetch.Msid('pm1thv2t_clean', '2020:001', '2020:010') assert len(dat.vals) == 36240 # Some bad values assert len(dat.source) == 36240 # Filtering applies to sources assert dat.colnames == colnames From 3e48f2db1362eff44ed4c8c19588a97725f9f2ff Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Sun, 15 Mar 2020 11:01:25 -0400 Subject: [PATCH 04/19] More tidy --- Ska/engarchive/derived/comps.py | 7 ------- 1 file changed, 7 deletions(-) diff --git a/Ska/engarchive/derived/comps.py b/Ska/engarchive/derived/comps.py index d19887af..fc0ab055 100644 --- a/Ska/engarchive/derived/comps.py +++ b/Ska/engarchive/derived/comps.py @@ -14,18 +14,11 @@ class ComputedMsid: def __init_subclass__(cls, **kwargs): super().__init_subclass__(**kwargs) - # Validate class name, msid_attrs, msid_match - if not cls.__name__.startswith('Comp_'): - raise ValueError(f'comp class name {cls.__name__} must start with Comp_') - cls.msid_attrs = ComputedMsid.msid_attrs + cls.extra_msid_attrs if not hasattr(cls, 'msid_match'): raise ValueError(f'comp {cls.__name__} must define msid_match') - # Force match to include entire line - cls.msid_match = cls.msid_match + '$' - cls.msid_classes.append(cls) @classmethod From e9a662b1047f7e13348ef9fdbc361aa4290b0e44 Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Sun, 15 Mar 2020 11:04:41 -0400 Subject: [PATCH 05/19] Flake8 fixes --- Ska/engarchive/derived/mups_valve.py | 11 ++++++++--- 1 file changed, 8 insertions(+), 3 deletions(-) diff --git a/Ska/engarchive/derived/mups_valve.py b/Ska/engarchive/derived/mups_valve.py index 51602295..e71350bf 100644 --- a/Ska/engarchive/derived/mups_valve.py +++ b/Ska/engarchive/derived/mups_valve.py @@ -91,13 +91,18 @@ def c2f(degc): count_to_degf = interp1d(pp_counts, pp_temps) degf_to_counts = interp1d(pp_temps, pp_counts) + # Define MUPS valve thermistor voltage point pair calibration table. -count_to_volts = lambda counts: counts / 256 * 5.12 -volts_to_counts = lambda volts: volts / 5.12 * 256 +def count_to_volts(counts): + return counts / 256 * 5.12 + + +def volts_to_counts(volts): + return volts / 5.12 * 256 + # Voltage and Temperature, with and without resistor (see flight note 447 for # the source of these numbers). - volt_with_resistor = [4.153325779, 3.676396578, 3.175100371, 2.587948965, 2.435, 2.025223702, 1.538506813, 1.148359251, 0.63128179, 0.354868907, 0.208375569] volt_without_resistor = [28.223, 15, 9.1231, 5.5228, 4.87, 3.467, 2.249, From 0c16991f3624b0d81c380c86924a8833383565ed Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Mon, 16 Mar 2020 10:02:28 -0400 Subject: [PATCH 06/19] Add Comp_KadiCommandState and docs --- Ska/engarchive/derived/comps.py | 83 +++++++++++++++++++++++++++++- Ska/engarchive/tests/test_comps.py | 17 ++++++ 2 files changed, 98 insertions(+), 2 deletions(-) diff --git a/Ska/engarchive/derived/comps.py b/Ska/engarchive/derived/comps.py index fc0ab055..b8d0a3eb 100644 --- a/Ska/engarchive/derived/comps.py +++ b/Ska/engarchive/derived/comps.py @@ -1,5 +1,16 @@ +# Licensed under a 3-clause BSD style license - see LICENSE.rst +""" +Support computed MSIDs in the cheta archive. + +- Base class ComputedMsid for user-generated comps. +- Cleaned MUPS valve temperatures MSIDs: '(pm2thv1t|pm1thv2t)_clean'. +- Commanded states 'cmd_state__
' for any kadi commanded state value. +""" + import re +import numpy as np + class ComputedMsid: # Global dict of registered computed MSIDs @@ -66,17 +77,85 @@ def get_msid_attrs(self, start, stop, msid, msid_args): class Comp_MUPS_Valve_Temp_Clean(ComputedMsid): + """Computed MSID for cleaned MUPS valve temps PM2THV1T, PM1THV2T + + This uses the cleaning method demonstrated in the following notebook + to return a version of the MUPS valve temperature comprised of + telemetry values that are consistent with a thermal model. + + https://nbviewer.jupyter.org/urls/cxc.cfa.harvard.edu/mta/ASPECT/ipynb/misc/mups-valve-xija-filtering.ipynb + + Allowed MSIDs are 'pm2thv1t_clean' and 'pm1thv2t_clean' (as always case is + not important). + """ msid_match = r'(pm2thv1t|pm1thv2t)_clean' extra_msid_attrs = ('vals_raw', 'vals_nan', 'vals_corr', 'vals_model', 'source') - def get_msid_attrs(self, start, stop, msid, msid_args): + def get_msid_attrs(self, tstart, tstop, msid, msid_args): + """Get attributes for computed MSID: ``vals``, ``bads``, ``times`` + + :param tstart: start time (CXC secs) + :param tstop: stop time (CXC secs) + :param msid: full MSID name e.g. pm2thv1t_clean + :param msid_args: tuple of regex match groups (msid_name,) + + :returns: dict of MSID attributes + """ from .mups_valve import fetch_clean_msid # Get cleaned MUPS valve temperature data as an MSID object - dat = fetch_clean_msid(msid_args[0], start, stop, + dat = fetch_clean_msid(msid_args[0], tstart, tstop, dt_thresh=5.0, median=7, model_spec=None) # Convert to dict as required by the get_msids_attrs API msid_attrs = {attr: getattr(dat, attr) for attr in self.msid_attrs} return msid_attrs + + +class Comp_KadiCommandState(ComputedMsid): + """Computed MSID for kadi dynamic commanded states. + + The MSID here takes the form ``cmd_state__
`` where: + - ``state_key`` is a valid commanded state key such as ``pitch`` or + ``pcad_mode`` or ``acisfp_temp``. + - ``dt`` is the sampling time expressed as a multiple of 1.025 sec + frames. + + Example MSID names:: + + 'cmd_state_pcad_mode_1': sample ``pcad_mode`` every 1.025 secs + 'cmd_state_acisfp_temp_32': sample ``acisfp_temp`` every 32.8 secs + """ + msid_match = r'cmd_state_(\w+)_(\d+)' + + def get_msid_attrs(self, start, stop, msid, msid_args): + """Get attributes for computed MSID: ``vals``, ``bads``, ``times`` + + :param tstart: start time (CXC secs) + :param tstop: stop time (CXC secs) + :param msid: full MSID name e.g. cmd_state_pitch_clean + :param msid_args: tuple of regex match groups: (state_key, dt) + + :returns: dict of MSID attributes + """ + from kadi.commands.states import get_states + from Chandra.Time import date2secs + + state_key = msid_args[0] + dt = 1.025 * int(msid_args[1]) + states = get_states(start, stop, state_keys=[state_key]) + + tstart = date2secs(states['datestart'][0]) + tstops = date2secs(states['datestop']) + + times = np.arange(tstart, tstops[-1], dt) + vals = states[state_key].view(np.ndarray) + + indexes = np.searchsorted(tstops, times) + + out = {'vals': vals[indexes], + 'times': times, + 'bads': np.zeros(len(times), dtype=bool)} + + return out diff --git a/Ska/engarchive/tests/test_comps.py b/Ska/engarchive/tests/test_comps.py index cdaa2367..de35b0b4 100644 --- a/Ska/engarchive/tests/test_comps.py +++ b/Ska/engarchive/tests/test_comps.py @@ -129,3 +129,20 @@ def test_mups_valve(): for attr in colnames: if attr != 'bads': assert len(dat.vals) == len(getattr(dat, attr)) + + +def test_cmd_states(): + start, stop = '2020:002:08:00:00', '2020:002:10:00:00' + dat = fetch.Msid('cmd_state_pitch_1000', start, stop) + exp_vals = np.array([55.99128956, 55.8747053, 55.8747053, 90.66266599, + 159.06945155, 173.11528258, 173.11528258, 173.11528258]) + assert np.allclose(dat.vals, exp_vals) + assert type(dat.vals) is np.ndarray + assert np.allclose(np.diff(dat.times), 1025.0) + assert not np.any(dat.bads) + + dat = fetch.Msid('cmd_state_pcad_mode_1000', start, stop) + exp_vals = np.array(['NPNT', 'NPNT', 'NPNT', 'NMAN', 'NMAN', 'NPNT', 'NPNT', 'NPNT']) + assert np.all(dat.vals == exp_vals) + assert type(dat.vals) is np.ndarray + assert np.allclose(np.diff(dat.times), 1025.0) From e9aa79bf5908992af2430bae9802bee401f85e4a Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Mon, 16 Mar 2020 10:05:13 -0400 Subject: [PATCH 07/19] Use tstart, tstop instead of start, stop --- Ska/engarchive/derived/comps.py | 4 ++-- Ska/engarchive/tests/test_comps.py | 8 ++++---- 2 files changed, 6 insertions(+), 6 deletions(-) diff --git a/Ska/engarchive/derived/comps.py b/Ska/engarchive/derived/comps.py index b8d0a3eb..60b184dc 100644 --- a/Ska/engarchive/derived/comps.py +++ b/Ska/engarchive/derived/comps.py @@ -129,7 +129,7 @@ class Comp_KadiCommandState(ComputedMsid): """ msid_match = r'cmd_state_(\w+)_(\d+)' - def get_msid_attrs(self, start, stop, msid, msid_args): + def get_msid_attrs(self, tstart, tstop, msid, msid_args): """Get attributes for computed MSID: ``vals``, ``bads``, ``times`` :param tstart: start time (CXC secs) @@ -144,7 +144,7 @@ def get_msid_attrs(self, start, stop, msid, msid_args): state_key = msid_args[0] dt = 1.025 * int(msid_args[1]) - states = get_states(start, stop, state_keys=[state_key]) + states = get_states(tstart, tstop, state_keys=[state_key]) tstart = date2secs(states['datestart'][0]) tstops = date2secs(states['datestop']) diff --git a/Ska/engarchive/tests/test_comps.py b/Ska/engarchive/tests/test_comps.py index de35b0b4..158e1d59 100644 --- a/Ska/engarchive/tests/test_comps.py +++ b/Ska/engarchive/tests/test_comps.py @@ -24,8 +24,8 @@ class Comp_Val_Plus_Five(ComputedMsid): """Silly base comp to add 5 to the value""" msid_match = r'comp_(\w+)_plus_five' - def get_msid_attrs(self, start, stop, msid, msid_args): - dat = self.fetch_eng.MSID(msid_args[0], start, stop) + def get_msid_attrs(self, tstart, tstop, msid, msid_args): + dat = self.fetch_eng.MSID(msid_args[0], tstart, tstop) out = {'vals': dat.vals + 5, 'bads': dat.bads, @@ -46,9 +46,9 @@ class Comp_CSS1_NPM_SUN(ComputedMsid, DerivedParameter): max_gap = 10.0 msid_match = 'comp_css1_npm_sun' - def get_msid_attrs(self, start, stop, msid, msid_args): + def get_msid_attrs(self, tstart, tstop, msid, msid_args): # Get an interpolated MSIDset for rootparams - msids = self.fetch(start, stop) + msids = self.fetch(tstart, tstop) # Do the computation and set bad values npm_sun = ((msids['aopcadmd'].vals == 'NPNT') & From 09732cdb596c89c00b0729c4731dcc9bd4d4f4a7 Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Tue, 17 Mar 2020 07:41:01 -0400 Subject: [PATCH 08/19] Support 5min, daily stats for computed MSIDs --- Ska/engarchive/derived/comps.py | 153 ++++++++++++++++++++++++++++- Ska/engarchive/fetch.py | 4 +- Ska/engarchive/tests/test_comps.py | 34 +++++++ 3 files changed, 183 insertions(+), 8 deletions(-) diff --git a/Ska/engarchive/derived/comps.py b/Ska/engarchive/derived/comps.py index 60b184dc..7736df59 100644 --- a/Ska/engarchive/derived/comps.py +++ b/Ska/engarchive/derived/comps.py @@ -11,6 +11,107 @@ import numpy as np +from Chandra.Time import DateTime + + +def calc_stats_vals(msid, rows, indexes, interval): + """ + Compute statistics values for ``msid`` over specified intervals. + This is a very slightly modified version of the same function from + update_archive.py. However, cannot directly import that because + it has side effects that break everything, probably related to + enabling caching. + + The mods here are basically to take out handling of state codes + and turn a warning about negative dts into an exception. + + :param msid: Msid object (filter_bad=True) + :param rows: Msid row indices corresponding to stat boundaries + :param indexes: Universal index values for stat (row times // dt) + :param interval: interval name (5min or daily) + """ + import scipy.stats + + quantiles = (1, 5, 16, 50, 84, 95, 99) + n_out = len(rows) - 1 + + # Check if data type is "numeric". Boolean values count as numeric, + # partly for historical reasons, in that they support funcs like + # mean (with implicit conversion to float). + msid_dtype = msid.vals.dtype + msid_is_numeric = issubclass(msid_dtype.type, (np.number, np.bool_)) + + # If MSID data is unicode, then for stats purposes cast back to bytes + # by creating the output array as a like-sized S-type array. + if msid_dtype.kind == 'U': + msid_dtype = re.sub(r'U', 'S', msid.vals.dtype.str) + + # Predeclare numpy arrays of correct type and sufficient size for accumulating results. + out = {} + out['index'] = np.ndarray((n_out,), dtype=np.int32) + out['n'] = np.ndarray((n_out,), dtype=np.int32) + out['val'] = np.ndarray((n_out,), dtype=msid_dtype) + + if msid_is_numeric: + out['min'] = np.ndarray((n_out,), dtype=msid_dtype) + out['max'] = np.ndarray((n_out,), dtype=msid_dtype) + out['mean'] = np.ndarray((n_out,), dtype=np.float32) + + if interval == 'daily': + out['std'] = np.ndarray((n_out,), dtype=msid_dtype) + for quantile in quantiles: + out['p{:02d}'.format(quantile)] = np.ndarray((n_out,), dtype=msid_dtype) + + i = 0 + for row0, row1, index in zip(rows[:-1], rows[1:], indexes[:-1]): + vals = msid.vals[row0:row1] + times = msid.times[row0:row1] + + n_vals = len(vals) + if n_vals > 0: + out['index'][i] = index + out['n'][i] = n_vals + out['val'][i] = vals[n_vals // 2] + if msid_is_numeric: + if n_vals <= 2: + dts = np.ones(n_vals, dtype=np.float64) + else: + dts = np.empty(n_vals, dtype=np.float64) + dts[0] = times[1] - times[0] + dts[-1] = times[-1] - times[-2] + dts[1:-1] = ((times[1:-1] - times[:-2]) + + (times[2:] - times[1:-1])) / 2.0 + negs = dts < 0.0 + if np.any(negs): + times_dts = [(DateTime(t).date, dt) + for t, dt in zip(times[negs], dts[negs])] + raise ValueError('WARNING - negative dts in {} at {}' + .format(msid.MSID, times_dts)) + + # Clip to range 0.001 to 300.0. The low bound is just there + # for data with identical time stamps. This shouldn't happen + # but in practice might. The 300.0 represents 5 minutes and + # is the largest normal time interval. Data near large gaps + # will get a weight of 5 mins. + dts.clip(0.001, 300.0, out=dts) + sum_dts = np.sum(dts) + + out['min'][i] = np.min(vals) + out['max'][i] = np.max(vals) + out['mean'][i] = np.sum(dts * vals) / sum_dts + if interval == 'daily': + # biased weighted estimator of variance (N should be big enough) + # http://en.wikipedia.org/wiki/Mean_square_weighted_deviation + sigma_sq = np.sum(dts * (vals - out['mean'][i]) ** 2) / sum_dts + out['std'][i] = np.sqrt(sigma_sq) + quant_vals = scipy.stats.mstats.mquantiles(vals, np.array(quantiles) / 100.0) + for quant_val, quantile in zip(quant_vals, quantiles): + out['p%02d' % quantile][i] = quant_val + + i += 1 + + return np.rec.fromarrays([x[:i] for x in out.values()], names=list(out.keys())) + class ComputedMsid: # Global dict of registered computed MSIDs @@ -56,25 +157,67 @@ def fetch_cxc(self): from .. import fetch_cxc return fetch_cxc - def __call__(self, start, stop, msid): + def __call__(self, tstart, tstop, msid, interval=None): match = re.match(self.msid_match, msid, re.IGNORECASE) if not match: raise RuntimeError(f'unexpected mismatch of {msid} with {self.msid_match}') match_args = [arg.lower() for arg in match.groups()] - msid_attrs = self.get_msid_attrs(start, stop, msid.lower(), match_args) - if set(msid_attrs) != set(self.msid_attrs): - raise ValueError(f'computed class did not return expected attributes') + if interval is None: + msid_attrs = self.get_msid_attrs(tstart, tstop, msid.lower(), match_args) + else: + msid_attrs = self.get_stats_attrs(tstart, tstop, msid.lower(), match_args, interval) return msid_attrs - def get_msid_attrs(self, start, stop, msid, msid_args): + def get_msid_attrs(self, tstart, tstop, msid, msid_args): """Get the attributes required for this MSID. TODO: detailed docs here since this is the main user-defined method """ raise NotImplementedError() + def get_stats_attrs(self, tstart, tstop, msid, match_args, interval): + from ..fetch import _plural + + class MsidStub: + def __init__(self, kwargs): + self.state_codes = None + for key, val in kwargs.items(): + setattr(self, key, val) + + # Replicate a stripped-down version of processing in update_archive. + # This produces a recarray with columns that correspond to the raw + # stats HDF5 files. + dt = {'5min': 328, + 'daily': 86400}[interval] + index0 = int(np.floor(tstart / dt)) + index1 = int(np.ceil(tstop / dt)) + tstart = (index0 - 1) * dt + tstop = (index1 + 1) * dt + msid_obj = self.fetch_eng.Msid(msid, tstart, tstop) + indexes = np.arange(index0, index1 + 1) + times = indexes * dt # This is the *start* time of each bin + + if len(times) > 0: + rows = np.searchsorted(msid_obj.times, times) + vals_stats = calc_stats_vals(msid_obj, rows, indexes, interval) + else: + raise ValueError() + + # Replicate the name munging that fetch does going from the HDF5 columns + # to what is seen in a stats fetch query. + out = {} + for key in vals_stats.dtype.names: + out_key = _plural(key) if key != 'n' else 'samples' + out[out_key] = vals_stats[key] + out['times'] = (vals_stats['index'] + 0.5) * dt + out['bads'] = np.zeros(len(vals_stats), dtype=bool) + out['midvals'] = out['vals'] + out['vals'] = out['means'] + + return out + class Comp_MUPS_Valve_Temp_Clean(ComputedMsid): """Computed MSID for cleaned MUPS valve temps PM2THV1T, PM1THV2T diff --git a/Ska/engarchive/fetch.py b/Ska/engarchive/fetch.py index eb21b589..faddd66f 100644 --- a/Ska/engarchive/fetch.py +++ b/Ska/engarchive/fetch.py @@ -585,8 +585,6 @@ def _get_data(self): comp_cls = ComputedMsid.get_matching_comp_cls(self.msid) if comp_cls: - if self.stat: - raise ValueError('stats are not supported for computed MSIDs') self._get_comp_data(comp_cls) return @@ -636,7 +634,7 @@ def _get_comp_data(self, comp_cls): logger.info(f'Getting comp data for {self.msid}') # Do computation. This returns a dict of MSID attribute values. - dat = comp_cls()(self.tstart, self.tstop, self.msid) + dat = comp_cls()(self.tstart, self.tstop, self.msid, self.stat) # Allow upstream class to be a bit sloppy on times and include samples # outside the time range. This can happen with classes that inherit diff --git a/Ska/engarchive/tests/test_comps.py b/Ska/engarchive/tests/test_comps.py index 158e1d59..31b00cfc 100644 --- a/Ska/engarchive/tests/test_comps.py +++ b/Ska/engarchive/tests/test_comps.py @@ -20,6 +20,19 @@ HAS_MAUDE = True +class Comp_Passthru(ComputedMsid): + """Pass MSID through unchanged (for checking that stats work)""" + msid_match = r'passthru_(\w+)' + + def get_msid_attrs(self, tstart, tstop, msid, msid_args): + dat = self.fetch_eng.MSID(msid_args[0], tstart, tstop) + + out = {'vals': dat.vals, + 'bads': dat.bads, + 'times': dat.times} + return out + + class Comp_Val_Plus_Five(ComputedMsid): """Silly base comp to add 5 to the value""" msid_match = r'comp_(\w+)_plus_five' @@ -146,3 +159,24 @@ def test_cmd_states(): assert np.all(dat.vals == exp_vals) assert type(dat.vals) is np.ndarray assert np.allclose(np.diff(dat.times), 1025.0) + + +@pytest.mark.parametrize('stat', ['5min', 'daily']) +def test_stats(stat): + start, stop = '2020:001', '2020:010' + + dat = fetch.Msid('pitch', start, stop, stat=stat) + datc = fetch.Msid('passthru_pitch', start, stop, stat=stat) + + for attr in datc.colnames: + val = getattr(dat, attr) + valc = getattr(datc, attr) + if attr == 'bads': + assert val == valc + continue + assert val.dtype == valc.dtype + if val.dtype.kind == 'f': + assert np.allclose(val, valc) + else: + assert np.all(val == valc) + From 988b51395dd29cad89a080191d9f599567ad0d18 Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Mon, 6 Apr 2020 18:12:05 -0400 Subject: [PATCH 09/19] Add temperature unit to MUPS valve comp msid --- Ska/engarchive/derived/comps.py | 2 +- Ska/engarchive/derived/mups_valve.py | 30 +++++++++++++++++++++------- Ska/engarchive/tests/test_comps.py | 7 ------- 3 files changed, 24 insertions(+), 15 deletions(-) diff --git a/Ska/engarchive/derived/comps.py b/Ska/engarchive/derived/comps.py index 7736df59..16521df9 100644 --- a/Ska/engarchive/derived/comps.py +++ b/Ska/engarchive/derived/comps.py @@ -248,7 +248,7 @@ def get_msid_attrs(self, tstart, tstop, msid, msid_args): # Get cleaned MUPS valve temperature data as an MSID object dat = fetch_clean_msid(msid_args[0], tstart, tstop, - dt_thresh=5.0, median=7, model_spec=None) + dt_thresh=5.0, median=7, model_spec=None, unit='degc') # Convert to dict as required by the get_msids_attrs API msid_attrs = {attr: getattr(dat, attr) for attr in self.msid_attrs} diff --git a/Ska/engarchive/derived/mups_valve.py b/Ska/engarchive/derived/mups_valve.py index e71350bf..1e8abf2d 100644 --- a/Ska/engarchive/derived/mups_valve.py +++ b/Ska/engarchive/derived/mups_valve.py @@ -82,6 +82,11 @@ def c2f(degc): return degf +def f2c(degf): + degc = (degf - 32) / 1.8 + return degc + + # Define MUPS valve thermistor point pair calibration table (from TDB). pp_counts = np.array([0, 27, 36, 44, 55, 70, 90, 118, 175, 195, 210, 219, 226, 231, 235, 255]) pp_temps = np.array([369.5305, 263.32577, 239.03652, 222.30608, 203.6944, 183.2642, 161.0796, @@ -187,7 +192,8 @@ def select_using_model(data1, data2, model, dt_thresh, out, source): source[ii] = 0 -def fetch_clean_msid(msid, start, stop=None, dt_thresh=5.0, median=7, model_spec=None): +def fetch_clean_msid(msid, start, stop=None, dt_thresh=5.0, median=7, + model_spec=None, unit='degf'): """Fetch a cleaned version of telemetry for ``msid``. If not supplied the model spec will be downloaded from github using this URL: @@ -212,6 +218,7 @@ def fetch_clean_msid(msid, start, stop=None, dt_thresh=5.0, median=7, model_spec :param dt_thresh: tolerance for matching model to data in degF (default=5 degF) :param median: length of median filter (default=7, use 0 to disable) :param model_spec: file name or URL containing relevant xija model spec + :param unit: temperature unit: 'degf' (default) or 'degc' :returns: fetch.Msid object """ @@ -220,6 +227,8 @@ def fetch_clean_msid(msid, start, stop=None, dt_thresh=5.0, median=7, model_spec f'mups_valve/{msid}_spec.json') model_spec = download_file(url, cache=True) + # NOTE: all temperatures are in degF unless otherwise noted + dat = fetch_eng.Msid(msid, start, stop) t_obs = dat.vals @@ -227,29 +236,36 @@ def fetch_clean_msid(msid, start, stop=None, dt_thresh=5.0, median=7, model_spec t_obs = median_filter(t_obs, size=median) t_corr = get_corr_mups_temp(t_obs) + # Model is computed in degC mdl = XijaModel(model_spec=model_spec, start=dat.times[0] - 400000, stop=dat.times[-1]) - mdl.comp['mups0'].set_data(75) + mdl.comp['mups0'].set_data(75) # degC mdl.make() mdl.calc() + # Interpolate model temperatures (converted to degF with c2f) onto telemetry times t_model = np.interp(xp=mdl.times, fp=c2f(mdl.comp[msid].mvals), x=dat.times) out = np.zeros_like(dat.vals) source = np.zeros(len(dat.vals), dtype=int) select_using_model(t_obs, t_corr, t_model, dt_thresh=dt_thresh, out=out, source=source) - dat.vals_raw = dat.vals + if unit not in ('degf', 'degc'): + raise ValueError('unit must be either degf or degc') + + convert = (lambda x: x) if unit == 'degf' else f2c + + dat.vals_raw = convert(dat.vals) - dat.vals = out - dat.vals_nan = out.copy() + dat.vals = convert(out) + dat.vals_nan = convert(out).copy() dat.bads = (source == 0) dat.vals_nan[dat.bads] = np.nan dat.source = source - dat.vals_corr = t_corr - dat.vals_model = t_model + dat.vals_corr = convert(t_corr) + dat.vals_model = convert(t_model) return dat diff --git a/Ska/engarchive/tests/test_comps.py b/Ska/engarchive/tests/test_comps.py index 31b00cfc..f3b58ea8 100644 --- a/Ska/engarchive/tests/test_comps.py +++ b/Ska/engarchive/tests/test_comps.py @@ -104,13 +104,6 @@ def test_simple_comp_with_maude(): assert np.all(dat1.bads == dat2.bads) -def test_comp_with_stat(): - with pytest.raises(ValueError, - match='stats are not supported for computed MSIDs'): - fetch_eng.Msid('comp_tephin_plus_five', - '2020:001', '2020:010', stat='5min') - - def test_mups_valve(): colnames = ['times', 'vals', 'bads', 'vals_raw', 'vals_nan', 'vals_corr', 'vals_model', 'source'] From e3455544e5b53e5f77fb23fc60b0193b6425b25a Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Mon, 6 Apr 2020 18:20:04 -0400 Subject: [PATCH 10/19] Remove orphaned MsidStub class (what was it?) --- Ska/engarchive/derived/comps.py | 6 ------ 1 file changed, 6 deletions(-) diff --git a/Ska/engarchive/derived/comps.py b/Ska/engarchive/derived/comps.py index 16521df9..20fe9e44 100644 --- a/Ska/engarchive/derived/comps.py +++ b/Ska/engarchive/derived/comps.py @@ -180,12 +180,6 @@ def get_msid_attrs(self, tstart, tstop, msid, msid_args): def get_stats_attrs(self, tstart, tstop, msid, match_args, interval): from ..fetch import _plural - class MsidStub: - def __init__(self, kwargs): - self.state_codes = None - for key, val in kwargs.items(): - setattr(self, key, val) - # Replicate a stripped-down version of processing in update_archive. # This produces a recarray with columns that correspond to the raw # stats HDF5 files. From 773e66a113884a6f5b0290f9ff443e4721767279 Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Tue, 7 Apr 2020 18:06:43 -0400 Subject: [PATCH 11/19] Remove unneccesary msid_attrs and extra_msid_attrs class attributes --- Ska/engarchive/derived/comps.py | 20 +++++++++----------- 1 file changed, 9 insertions(+), 11 deletions(-) diff --git a/Ska/engarchive/derived/comps.py b/Ska/engarchive/derived/comps.py index 20fe9e44..161406ac 100644 --- a/Ska/engarchive/derived/comps.py +++ b/Ska/engarchive/derived/comps.py @@ -117,17 +117,9 @@ class ComputedMsid: # Global dict of registered computed MSIDs msid_classes = [] - # Standard base MSID attributes that must be provided - msid_attrs = ('times', 'vals', 'bads') - - # Extra MSID attributes that are provided beyond times, vals, bads - extra_msid_attrs = () - def __init_subclass__(cls, **kwargs): super().__init_subclass__(**kwargs) - cls.msid_attrs = ComputedMsid.msid_attrs + cls.extra_msid_attrs - if not hasattr(cls, 'msid_match'): raise ValueError(f'comp {cls.__name__} must define msid_match') @@ -165,6 +157,10 @@ def __call__(self, tstart, tstop, msid, interval=None): if interval is None: msid_attrs = self.get_msid_attrs(tstart, tstop, msid.lower(), match_args) + for attr in ('vals', 'bads', 'times'): + if attr not in msid_attrs: + raise ValueError(f'computed MSID {self.__class__.__name__} failed ' + f'to set required attribute {attr}') else: msid_attrs = self.get_stats_attrs(tstart, tstop, msid.lower(), match_args, interval) @@ -226,7 +222,6 @@ class Comp_MUPS_Valve_Temp_Clean(ComputedMsid): not important). """ msid_match = r'(pm2thv1t|pm1thv2t)_clean' - extra_msid_attrs = ('vals_raw', 'vals_nan', 'vals_corr', 'vals_model', 'source') def get_msid_attrs(self, tstart, tstop, msid, msid_args): """Get attributes for computed MSID: ``vals``, ``bads``, ``times`` @@ -244,8 +239,11 @@ def get_msid_attrs(self, tstart, tstop, msid, msid_args): dat = fetch_clean_msid(msid_args[0], tstart, tstop, dt_thresh=5.0, median=7, model_spec=None, unit='degc') - # Convert to dict as required by the get_msids_attrs API - msid_attrs = {attr: getattr(dat, attr) for attr in self.msid_attrs} + # Convert to dict as required by the get_msids_attrs API. `fetch_clean_msid` + # returns an MSID object with the following attrs. + attrs = ('vals', 'times', 'bads', 'vals_raw', 'vals_nan', + 'vals_corr', 'vals_model', 'source') + msid_attrs = {attr: getattr(dat, attr) for attr in attrs} return msid_attrs From 3fe7e8fd0a65b480c87727c7c5a7fbee8c67cf55 Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Tue, 7 Apr 2020 18:07:20 -0400 Subject: [PATCH 12/19] Add some docs --- Ska/engarchive/derived/comps.py | 39 ++++++++++++++++++++++++++++++ Ska/engarchive/tests/test_comps.py | 3 +-- doc/pseudo_msids.rst | 36 ++++++++++++++++++++++++++- 3 files changed, 75 insertions(+), 3 deletions(-) diff --git a/Ska/engarchive/derived/comps.py b/Ska/engarchive/derived/comps.py index 161406ac..5e0c864a 100644 --- a/Ska/engarchive/derived/comps.py +++ b/Ska/engarchive/derived/comps.py @@ -13,6 +13,8 @@ from Chandra.Time import DateTime +__all__ = ['ComputedMsid', 'Comp_MUPS_Valve_Temp_Clean', 'Comp_KadiCommandState'] + def calc_stats_vals(msid, rows, indexes, interval): """ @@ -114,10 +116,45 @@ def calc_stats_vals(msid, rows, indexes, interval): class ComputedMsid: + """Base class for cheta computed MSID. + + Sub-classes must define at least the following: + + * ``msid_match`` class attribute as a regex to match for the MSID. + * ``get_msid_attrs()`` method to perform the computation and return + a dict with the result. + + Example:: + + class Comp_Val_Plus_Offset(ComputedMsid): + ''' + Computed MSID to add an integer offset to MSID value. + + MSID format is "comp__plus_". + ''' + # Regex to match including arguments and + msid_match = r'comp_(\w+)_plus_(\d+)' + + def get_msid_attrs(self, tstart, tstop, msid, msid_args): + offset = int(msid_args[1]) + dat = self.fetch_eng.MSID(msid_args[0], tstart, tstop) + + # Must return a dict with at least `vals`, `times` and `bads`. + # Additional attributes are allowed and will be set on the + # final MSID object. + out = {'vals': dat.vals + int(msid_args[1]), + 'bads': dat.bads, + 'times': dat.times, + 'vals_raw': dat.vals, + 'offset': offset} + return out + """ # Global dict of registered computed MSIDs msid_classes = [] def __init_subclass__(cls, **kwargs): + """Validate and register ComputedMSID subclass. + """ super().__init_subclass__(**kwargs) if not hasattr(cls, 'msid_match'): @@ -134,6 +171,8 @@ def get_matching_comp_cls(cls, msid): return None + # These three properties are provided as a convenience because the module + # itself cannot import fetch because this is circular. @property def fetch_eng(self): from .. import fetch_eng diff --git a/Ska/engarchive/tests/test_comps.py b/Ska/engarchive/tests/test_comps.py index f3b58ea8..3752d1a8 100644 --- a/Ska/engarchive/tests/test_comps.py +++ b/Ska/engarchive/tests/test_comps.py @@ -105,7 +105,7 @@ def test_simple_comp_with_maude(): def test_mups_valve(): - colnames = ['times', 'vals', 'bads', 'vals_raw', + colnames = ['vals', 'times', 'bads', 'vals_raw', 'vals_nan', 'vals_corr', 'vals_model', 'source'] dat = fetch.MSID('PM2THV1T_clean', '2020:001', '2020:010') @@ -172,4 +172,3 @@ def test_stats(stat): assert np.allclose(val, valc) else: assert np.all(val == valc) - diff --git a/doc/pseudo_msids.rst b/doc/pseudo_msids.rst index 8483420c..84f31b68 100644 --- a/doc/pseudo_msids.rst +++ b/doc/pseudo_msids.rst @@ -8,7 +8,8 @@ stream are also available in the archive. These are: * ACIS DEA housekeeping: status from the DEA (including detector focal plane temperature) * Ephemeris: predictive and definitive orbital (Chandra), solar, and lunar ephemeris values * SIM telemetry: SIM position and moving status -* Derived parameters: values computed from other MSIDs in the archive +* Derived parameters: values computed from other MSIDs in the archive and stored in archive files +* Computed parameters: values computed on the fly from other MSIDs ACIS DEA housekeeping -------------------------------------------------- @@ -656,3 +657,36 @@ Thermal :members: :undoc-members: +Computed MSIDs +-------------- + +Cheta provides support for on-the-fly computed MSIDs with the following features: + +* Designed for simplicity and ease of use by non-expert Ska3 users. +* Following a simple recipe in user code, new MSIDs are automatically registered to fetch. +* Computed MSIDs can included embedded parameters to allow further customization or application of a function to any other MSID. +* Support for 5-minute and daily stats also included. + +What's the advantage? +^^^^^^^^^^^^^^^^^^^^^ + +If you have a function of MSID values, what is the advantage of going through +this formalism to create a computed MSID instead of just using the function +output directly? + +* It gives you the entire :func:`~Ska.engarchive.fetch.MSID` API! You get for free all the + `fetch bling Date: Tue, 7 Apr 2020 18:08:44 -0400 Subject: [PATCH 13/19] Fix flake8 oddity --- Ska/engarchive/derived/comps.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Ska/engarchive/derived/comps.py b/Ska/engarchive/derived/comps.py index 5e0c864a..a5d26af7 100644 --- a/Ska/engarchive/derived/comps.py +++ b/Ska/engarchive/derived/comps.py @@ -133,7 +133,7 @@ class Comp_Val_Plus_Offset(ComputedMsid): MSID format is "comp__plus_". ''' # Regex to match including arguments and - msid_match = r'comp_(\w+)_plus_(\d+)' + msid_match = r'comp_(\w+)_plus_(\d+)' # noqa def get_msid_attrs(self, tstart, tstop, msid, msid_args): offset = int(msid_args[1]) From 18629750bbc36615a22bd33e2ce63d2116a6ebe3 Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Wed, 8 Apr 2020 13:49:54 -0400 Subject: [PATCH 14/19] Add more docs --- doc/fetch_tutorial.rst | 10 ++++- doc/pseudo_msids.rst | 87 +++++++++++++++++++++++++++++++++++++----- 2 files changed, 86 insertions(+), 11 deletions(-) diff --git a/doc/fetch_tutorial.rst b/doc/fetch_tutorial.rst index b6615322..2084ce02 100644 --- a/doc/fetch_tutorial.rst +++ b/doc/fetch_tutorial.rst @@ -1547,7 +1547,7 @@ local or remote. In order to use this option, the user must have a special key file ``ska_remote_access.json``placed at the root of the local Python installation folder. This is placed in the directory shown with ``import sys; print(sys.prefix)``. To get a copy of this file - contact Mark Baski or Tom Aldcroft. +contact Mark Baski or Tom Aldcroft. Remote access is controlled as follows: @@ -1559,3 +1559,11 @@ Remote access is controlled as follows: is *enabled* unless the system finds a local engineering data archive. It looks for data in either ``$SKA/data/eng_archive`` or ``$ENG_ARCHIVE``, where those refer to user-defined environment variables. + +Local cheta archive +=================== + +Instructions for creating, using, and maintaining a local cheta archive (on your +laptop for instance) are found at the `Tutorial for installing and maintaining a +cheta telemetry archive +`_. diff --git a/doc/pseudo_msids.rst b/doc/pseudo_msids.rst index 84f31b68..c018ccb0 100644 --- a/doc/pseudo_msids.rst +++ b/doc/pseudo_msids.rst @@ -663,10 +663,16 @@ Computed MSIDs Cheta provides support for on-the-fly computed MSIDs with the following features: * Designed for simplicity and ease of use by non-expert Ska3 users. -* Following a simple recipe in user code, new MSIDs are automatically registered to fetch. -* Computed MSIDs can included embedded parameters to allow further customization or application of a function to any other MSID. +* Following a simple recipe in user code, new MSIDs are automatically registered + to fetch. +* Computed MSIDs can included embedded parameters to allow further customization + or application of a function to any other MSID. * Support for 5-minute and daily stats also included. +See the `DAWG computed MSIDs notebook +`_ +for a longer introduction to the topic with a number of examples. + What's the advantage? ^^^^^^^^^^^^^^^^^^^^^ @@ -674,19 +680,80 @@ If you have a function of MSID values, what is the advantage of going through this formalism to create a computed MSID instead of just using the function output directly? -* It gives you the entire :func:`~Ska.engarchive.fetch.MSID` API! You get for free all the - `fetch bling `_ like plotting, selecting intervals, kadi event integration, interpolation. -* If your computed MSID is useful enough it will be trivial to add to the released ``cheta`` - for other users. +* If your computed MSID is useful to the community it is simple to add to the + released ``cheta`` package for other users. * It allows creation of arbitrary MSIDs that can be used in xija models. Examples: - ``pm2tv1t_clean`` as a ``Node`` component - ``cmd_state_acisfp_temp_32`` as a telemetry input (``TelemData`` component). - Built-in computed MSIDs and API - ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ - - .. automodule:: Ska.engarchive.derived.thermal +Example +^^^^^^^ + +This example is not terribly useful but illustrates the key concepts and +requirements for defining a computed MSID:: + + from cheta.derived.comps import ComputedMsid # Inherit from this class + + # Class name is arbitrary, but by convention start with `Comp_` + class Comp_Val_Plus_Offset(ComputedMsid): + """ + Computed MSID to add an integer offset to MSID value. + + MSID format is "_plus_", where is an existing MSID + in the archive and is an integer offset. + + """ + msid_match = r'(\w+)_plus_(\d+)' + + # `msid_match` is a class attribute that defines a regular expresion to + # match for this computed MSID. This must be defined and it must be + # unambiguous (not matching an existing MSID or other computed MSID). + # + # The two groups in parentheses specify the arguments and . + # These are passed to `get_msid_attrs` as msid_args[0] and msid_args[1]. + # The \w symbol means to match a-z, A-Z, 0-9 and underscore (_). + # The \d symbol means to match digits 0-9. + + def get_msid_attrs(self, tstart, tstop, msid, msid_args): + """ + Get attributes for computed MSID: ``vals``, ``bads``, ``times``, + ``raw_vals``, and ``offset``. The first three must always be + provided. + + :param tstart: start time (CXC secs) + :param tstop: stop time (CXC secs) + :param msid: full MSID name e.g. tephin_plus_5 + :param msid_args: tuple of regex match groups (msid_name,) + :returns: dict of MSID attributes + """ + # Process the arguments parsed from the MSID + msid msid_args[0] + offset = int(msid_args[1]) + + # Get the raw telemetry value in engineering units + dat = self.fetch_eng.MSID(msid, tstart, tstop) + + # Do the computation + vals = dat.vals + offset + + # Return a dict with at least `vals`, `times` and `bads`. + # Additional attributes are allowed and will be set on the + # final MSID object. + out = {'vals': vals, + 'bads': dat.bads, + 'times': dat.times, + 'vals_raw': dat.vals, # Provide original values without offset + 'offset': offset # Provide the offset for reference + } + return out + +Built-in computed MSIDs and API +^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ + + .. automodule:: Ska.engarchive.derived.comps :members: :undoc-members: From 570a79cabb8251a0b3c9ed0dbe57084aad3ceeb1 Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Wed, 8 Apr 2020 14:02:21 -0400 Subject: [PATCH 15/19] Fix indentation --- doc/pseudo_msids.rst | 112 +++++++++++++++++++++---------------------- 1 file changed, 56 insertions(+), 56 deletions(-) diff --git a/doc/pseudo_msids.rst b/doc/pseudo_msids.rst index c018ccb0..b07c21a3 100644 --- a/doc/pseudo_msids.rst +++ b/doc/pseudo_msids.rst @@ -686,8 +686,8 @@ output directly? * If your computed MSID is useful to the community it is simple to add to the released ``cheta`` package for other users. * It allows creation of arbitrary MSIDs that can be used in xija models. Examples: - - ``pm2tv1t_clean`` as a ``Node`` component - - ``cmd_state_acisfp_temp_32`` as a telemetry input (``TelemData`` component). + * ``pm2tv1t_clean`` as a ``Node`` component + * ``cmd_state_acisfp_temp_32`` as a telemetry input (``TelemData`` component). Example ^^^^^^^ @@ -695,60 +695,60 @@ Example This example is not terribly useful but illustrates the key concepts and requirements for defining a computed MSID:: - from cheta.derived.comps import ComputedMsid # Inherit from this class - - # Class name is arbitrary, but by convention start with `Comp_` - class Comp_Val_Plus_Offset(ComputedMsid): - """ - Computed MSID to add an integer offset to MSID value. - - MSID format is "_plus_", where is an existing MSID - in the archive and is an integer offset. - - """ - msid_match = r'(\w+)_plus_(\d+)' - - # `msid_match` is a class attribute that defines a regular expresion to - # match for this computed MSID. This must be defined and it must be - # unambiguous (not matching an existing MSID or other computed MSID). - # - # The two groups in parentheses specify the arguments and . - # These are passed to `get_msid_attrs` as msid_args[0] and msid_args[1]. - # The \w symbol means to match a-z, A-Z, 0-9 and underscore (_). - # The \d symbol means to match digits 0-9. - - def get_msid_attrs(self, tstart, tstop, msid, msid_args): - """ - Get attributes for computed MSID: ``vals``, ``bads``, ``times``, - ``raw_vals``, and ``offset``. The first three must always be - provided. - - :param tstart: start time (CXC secs) - :param tstop: stop time (CXC secs) - :param msid: full MSID name e.g. tephin_plus_5 - :param msid_args: tuple of regex match groups (msid_name,) - :returns: dict of MSID attributes - """ - # Process the arguments parsed from the MSID - msid msid_args[0] - offset = int(msid_args[1]) - - # Get the raw telemetry value in engineering units - dat = self.fetch_eng.MSID(msid, tstart, tstop) - - # Do the computation - vals = dat.vals + offset - - # Return a dict with at least `vals`, `times` and `bads`. - # Additional attributes are allowed and will be set on the - # final MSID object. - out = {'vals': vals, - 'bads': dat.bads, - 'times': dat.times, - 'vals_raw': dat.vals, # Provide original values without offset - 'offset': offset # Provide the offset for reference - } - return out + from cheta.derived.comps import ComputedMsid # Inherit from this class + + # Class name is arbitrary, but by convention start with `Comp_` + class Comp_Val_Plus_Offset(ComputedMsid): + """ + Computed MSID to add an integer offset to MSID value. + + MSID format is "_plus_", where is an existing MSID + in the archive and is an integer offset. + + """ + msid_match = r'(\w+)_plus_(\d+)' + + # `msid_match` is a class attribute that defines a regular expresion to + # match for this computed MSID. This must be defined and it must be + # unambiguous (not matching an existing MSID or other computed MSID). + # + # The two groups in parentheses specify the arguments and . + # These are passed to `get_msid_attrs` as msid_args[0] and msid_args[1]. + # The \w symbol means to match a-z, A-Z, 0-9 and underscore (_). + # The \d symbol means to match digits 0-9. + + def get_msid_attrs(self, tstart, tstop, msid, msid_args): + """ + Get attributes for computed MSID: ``vals``, ``bads``, ``times``, + ``raw_vals``, and ``offset``. The first three must always be + provided. + + :param tstart: start time (CXC secs) + :param tstop: stop time (CXC secs) + :param msid: full MSID name e.g. tephin_plus_5 + :param msid_args: tuple of regex match groups (msid_name,) + :returns: dict of MSID attributes + """ + # Process the arguments parsed from the MSID + msid msid_args[0] + offset = int(msid_args[1]) + + # Get the raw telemetry value in engineering units + dat = self.fetch_eng.MSID(msid, tstart, tstop) + + # Do the computation + vals = dat.vals + offset + + # Return a dict with at least `vals`, `times` and `bads`. + # Additional attributes are allowed and will be set on the + # final MSID object. + out = {'vals': vals, + 'bads': dat.bads, + 'times': dat.times, + 'vals_raw': dat.vals, # Provide original values without offset + 'offset': offset # Provide the offset for reference + } + return out Built-in computed MSIDs and API ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ From b5dec2330b333e80a36f69ba0dad3888b1b47033 Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Wed, 8 Apr 2020 14:08:44 -0400 Subject: [PATCH 16/19] Even more docs --- Ska/engarchive/derived/comps.py | 14 ++++++++++++-- 1 file changed, 12 insertions(+), 2 deletions(-) diff --git a/Ska/engarchive/derived/comps.py b/Ska/engarchive/derived/comps.py index a5d26af7..dcd9b196 100644 --- a/Ska/engarchive/derived/comps.py +++ b/Ska/engarchive/derived/comps.py @@ -5,6 +5,8 @@ - Base class ComputedMsid for user-generated comps. - Cleaned MUPS valve temperatures MSIDs: '(pm2thv1t|pm1thv2t)_clean'. - Commanded states 'cmd_state__
' for any kadi commanded state value. + +See: https://nbviewer.jupyter.org/urls/cxc.harvard.edu/mta/ASPECT/ipynb/misc/DAWG-mups-valve-xija-filtering.ipynb """ import re @@ -208,9 +210,16 @@ def __call__(self, tstart, tstop, msid, interval=None): def get_msid_attrs(self, tstart, tstop, msid, msid_args): """Get the attributes required for this MSID. - TODO: detailed docs here since this is the main user-defined method + Get attributes for computed MSID, which must include at least + ``vals``, ``bads``, ``times``, and may include additional attributes. + + :param tstart: start time (CXC secs) + :param tstop: stop time (CXC secs) + :param msid: full MSID name e.g. tephin_plus_5 + :param msid_args: tuple of regex match groups (msid_name,) + :returns: dict of MSID attributes """ - raise NotImplementedError() + raise NotImplementedError('sub-class must implement get_msid_attrs()') def get_stats_attrs(self, tstart, tstop, msid, match_args, interval): from ..fetch import _plural @@ -256,6 +265,7 @@ class Comp_MUPS_Valve_Temp_Clean(ComputedMsid): telemetry values that are consistent with a thermal model. https://nbviewer.jupyter.org/urls/cxc.cfa.harvard.edu/mta/ASPECT/ipynb/misc/mups-valve-xija-filtering.ipynb + https://nbviewer.jupyter.org/urls/cxc.harvard.edu/mta/ASPECT/ipynb/misc/DAWG-mups-valve-xija-filtering.ipynb Allowed MSIDs are 'pm2thv1t_clean' and 'pm1thv2t_clean' (as always case is not important). From 2114701389e0d6d2130ca2fdab3b8bf3680ac159 Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Wed, 8 Apr 2020 17:37:09 -0400 Subject: [PATCH 17/19] Remove unused functions --- Ska/engarchive/derived/mups_valve.py | 10 ++-------- 1 file changed, 2 insertions(+), 8 deletions(-) diff --git a/Ska/engarchive/derived/mups_valve.py b/Ska/engarchive/derived/mups_valve.py index 1e8abf2d..94460f7c 100644 --- a/Ska/engarchive/derived/mups_valve.py +++ b/Ska/engarchive/derived/mups_valve.py @@ -112,14 +112,8 @@ def volts_to_counts(volts): 1.538506813, 1.148359251, 0.63128179, 0.354868907, 0.208375569] volt_without_resistor = [28.223, 15, 9.1231, 5.5228, 4.87, 3.467, 2.249, 1.5027, 0.7253, 0.38276, 0.21769] -temp_without_resistor = [50, 77, 100, 125, 130, 150, 175, 200, 250, 300, 350] - -volt_without_resistor_to_temp_without_resistor = interp1d(volt_without_resistor, - temp_without_resistor) -temp_without_resistor_to_volt_without_resistor = interp1d(temp_without_resistor, - volt_without_resistor) -volt_with_resistor_to_volt_without_resistor = interp1d(volt_with_resistor, volt_without_resistor) -volt_without_resistor_to_volt_with_resistor = interp1d(volt_without_resistor, volt_with_resistor) +volt_without_resistor_to_volt_with_resistor = interp1d(volt_without_resistor, + volt_with_resistor) def get_corr_mups_temp(temp): From 81a9bd3f114558c37d6a7b806fc52eb494c46939 Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Thu, 9 Apr 2020 06:54:08 -0400 Subject: [PATCH 18/19] Add support for units --- Ska/engarchive/derived/comps.py | 155 ++++++++++++++++++++------- Ska/engarchive/derived/mups_valve.py | 24 ++--- Ska/engarchive/fetch.py | 17 +-- Ska/engarchive/tests/test_comps.py | 67 ++++++++---- Ska/engarchive/units.py | 8 ++ doc/pseudo_msids.rst | 43 +++++++- 6 files changed, 226 insertions(+), 88 deletions(-) diff --git a/Ska/engarchive/derived/comps.py b/Ska/engarchive/derived/comps.py index dcd9b196..80c4ddf1 100644 --- a/Ska/engarchive/derived/comps.py +++ b/Ska/engarchive/derived/comps.py @@ -12,9 +12,10 @@ import re import numpy as np - from Chandra.Time import DateTime +from ..units import converters as unit_converter_funcs + __all__ = ['ComputedMsid', 'Comp_MUPS_Valve_Temp_Clean', 'Comp_KadiCommandState'] @@ -33,6 +34,8 @@ def calc_stats_vals(msid, rows, indexes, interval): :param rows: Msid row indices corresponding to stat boundaries :param indexes: Universal index values for stat (row times // dt) :param interval: interval name (5min or daily) + + :returns: np.recarray of stats values """ import scipy.stats @@ -126,37 +129,23 @@ class ComputedMsid: * ``get_msid_attrs()`` method to perform the computation and return a dict with the result. - Example:: - - class Comp_Val_Plus_Offset(ComputedMsid): - ''' - Computed MSID to add an integer offset to MSID value. - - MSID format is "comp__plus_". - ''' - # Regex to match including arguments and - msid_match = r'comp_(\w+)_plus_(\d+)' # noqa - - def get_msid_attrs(self, tstart, tstop, msid, msid_args): - offset = int(msid_args[1]) - dat = self.fetch_eng.MSID(msid_args[0], tstart, tstop) - - # Must return a dict with at least `vals`, `times` and `bads`. - # Additional attributes are allowed and will be set on the - # final MSID object. - out = {'vals': dat.vals + int(msid_args[1]), - 'bads': dat.bads, - 'times': dat.times, - 'vals_raw': dat.vals, - 'offset': offset} - return out + Optionally: + + * ``units`` attribute to specify unit handling. + + See the fetch tutorial Computed MSIDs section for details. """ # Global dict of registered computed MSIDs msid_classes = [] + # Base units specification (None implies no unit handling) + units = None + + def __init__(self, unit_system='eng'): + self.unit_system = unit_system + def __init_subclass__(cls, **kwargs): - """Validate and register ComputedMSID subclass. - """ + """Validate and register ComputedMSID subclass.""" super().__init_subclass__(**kwargs) if not hasattr(cls, 'msid_match'): @@ -166,6 +155,11 @@ def __init_subclass__(cls, **kwargs): @classmethod def get_matching_comp_cls(cls, msid): + """Get computed classes that match ``msid`` + + :param msid: str, input msid + :returns: first ComputedMsid subclass that matches ``msid`` or None + """ for comp_cls in ComputedMsid.msid_classes: match = re.match(comp_cls.msid_match + '$', msid, re.IGNORECASE) if match: @@ -173,46 +167,102 @@ def get_matching_comp_cls(cls, msid): return None - # These three properties are provided as a convenience because the module + # These four properties are provided as a convenience because the module # itself cannot import fetch because this is circular. @property def fetch_eng(self): + """Fetch in TDB engineering units like DEGF""" from .. import fetch_eng return fetch_eng @property def fetch_sci(self): + """Fetch in scientific units like DEGC""" from .. import fetch_sci return fetch_sci @property def fetch_cxc(self): - from .. import fetch_cxc - return fetch_cxc + """Fetch in CXC (FITS standard) units like K""" + from .. import fetch + return fetch + + @property + def fetch_sys(self): + """Fetch in the unit system specified for the class""" + fetch = getattr(self, f'fetch_{self.unit_system}') + return fetch def __call__(self, tstart, tstop, msid, interval=None): + """Emulate the fetch.MSID() API, but return a dict of MSID attributes. + + The returned dict turned into a proper MSID object by the upstream caller + `fetch.MSID._get_comp_data()`. + + :param tstart: float, start time (CXC seconds) + :param tstop: float, stop time (CXC seconds) + :param msid: str, MSID name + :param interval: str or None, stats interval (None, '5min', 'daily') + + :returns: dict of MSID attributes including 'times', 'vals', 'bads' + """ + # Parse any arguments from the input `msid` match = re.match(self.msid_match, msid, re.IGNORECASE) if not match: raise RuntimeError(f'unexpected mismatch of {msid} with {self.msid_match}') match_args = [arg.lower() for arg in match.groups()] if interval is None: + # Call the actual user-supplied work method to compute the MSID values msid_attrs = self.get_msid_attrs(tstart, tstop, msid.lower(), match_args) - for attr in ('vals', 'bads', 'times'): + + for attr in ('vals', 'bads', 'times', 'unit'): if attr not in msid_attrs: raise ValueError(f'computed MSID {self.__class__.__name__} failed ' f'to set required attribute {attr}') + + # Presence of a non-None `units` class attribute means that the MSID has + # units that should be converted to `self.unit_system` if required, where + # unit_system is 'cxc', 'sci', or 'eng'. + if self.units is not None: + msid_attrs = self.convert_units(msid_attrs) else: - msid_attrs = self.get_stats_attrs(tstart, tstop, msid.lower(), match_args, interval) + msid_attrs = self.get_stats_attrs(tstart, tstop, msid.lower(), match_args, + interval) return msid_attrs + def convert_units(self, msid_attrs): + """ + Convert required elements of ``msid_attrs`` to ``self.unit_system``. + + Unit_system can be one of 'cxc', 'sci', 'eng'. + + :param msid_attrs: dict, input MSID attributes + :param unit_system: str, unit system + + :returns: dict, converted MSID attributes + """ + unit_current = self.units[self.units['internal_system']] + unit_new = self.units[self.unit_system] + + out = msid_attrs.copy() + out['unit'] = unit_new + + if unit_current != unit_new: + for attr in self.units['convert_attrs']: + out[attr] = unit_converter_funcs[unit_current, unit_new](msid_attrs[attr]) + + return out + def get_msid_attrs(self, tstart, tstop, msid, msid_args): """Get the attributes required for this MSID. Get attributes for computed MSID, which must include at least ``vals``, ``bads``, ``times``, and may include additional attributes. + This MUST be supplied by sub-classes. + :param tstart: start time (CXC secs) :param tstop: stop time (CXC secs) :param msid: full MSID name e.g. tephin_plus_5 @@ -222,6 +272,16 @@ def get_msid_attrs(self, tstart, tstop, msid, msid_args): raise NotImplementedError('sub-class must implement get_msid_attrs()') def get_stats_attrs(self, tstart, tstop, msid, match_args, interval): + """Get 5-min or daily stats attributes. + + This is normally not overridden by sub-classes. + + :param tstart: start time (CXC secs) + :param tstop: stop time (CXC secs) + :param msid: full MSID name e.g. tephin_plus_5 + :param msid_args: tuple of regex match groups (msid_name,) + :returns: dict of MSID attributes + """ from ..fetch import _plural # Replicate a stripped-down version of processing in update_archive. @@ -233,7 +293,9 @@ def get_stats_attrs(self, tstart, tstop, msid, match_args, interval): index1 = int(np.ceil(tstop / dt)) tstart = (index0 - 1) * dt tstop = (index1 + 1) * dt - msid_obj = self.fetch_eng.Msid(msid, tstart, tstop) + + msid_obj = self.fetch_sys.Msid(msid, tstart, tstop) + indexes = np.arange(index0, index1 + 1) times = indexes * dt # This is the *start* time of each bin @@ -253,9 +315,14 @@ def get_stats_attrs(self, tstart, tstop, msid, match_args, interval): out['bads'] = np.zeros(len(vals_stats), dtype=bool) out['midvals'] = out['vals'] out['vals'] = out['means'] + out['unit'] = msid_obj.unit return out +############################################################################ +# Built-in computed MSIDs +############################################################################ + class Comp_MUPS_Valve_Temp_Clean(ComputedMsid): """Computed MSID for cleaned MUPS valve temps PM2THV1T, PM1THV2T @@ -272,8 +339,17 @@ class Comp_MUPS_Valve_Temp_Clean(ComputedMsid): """ msid_match = r'(pm2thv1t|pm1thv2t)_clean' + units = { + 'internal_system': 'eng', # Unit system for attrs from get_msid_attrs() + 'eng': 'DEGF', # Units for eng, sci, cxc systems + 'sci': 'DEGC', + 'cxc': 'K', + # Attributes that need conversion + 'convert_attrs': ['vals', 'vals_raw', 'vals_nan', 'vals_corr', 'vals_model'] + } + def get_msid_attrs(self, tstart, tstop, msid, msid_args): - """Get attributes for computed MSID: ``vals``, ``bads``, ``times`` + """Get attributes for computed MSID: ``vals``, ``bads``, ``times``, ``unit`` :param tstart: start time (CXC secs) :param tstop: stop time (CXC secs) @@ -286,13 +362,14 @@ def get_msid_attrs(self, tstart, tstop, msid, msid_args): # Get cleaned MUPS valve temperature data as an MSID object dat = fetch_clean_msid(msid_args[0], tstart, tstop, - dt_thresh=5.0, median=7, model_spec=None, unit='degc') + dt_thresh=5.0, median=7, model_spec=None) # Convert to dict as required by the get_msids_attrs API. `fetch_clean_msid` # returns an MSID object with the following attrs. attrs = ('vals', 'times', 'bads', 'vals_raw', 'vals_nan', 'vals_corr', 'vals_model', 'source') msid_attrs = {attr: getattr(dat, attr) for attr in attrs} + msid_attrs['unit'] = 'DEGF' return msid_attrs @@ -301,9 +378,10 @@ class Comp_KadiCommandState(ComputedMsid): """Computed MSID for kadi dynamic commanded states. The MSID here takes the form ``cmd_state__
`` where: - - ``state_key`` is a valid commanded state key such as ``pitch`` or + + * ``state_key`` is a valid commanded state key such as ``pitch`` or ``pcad_mode`` or ``acisfp_temp``. - - ``dt`` is the sampling time expressed as a multiple of 1.025 sec + * ``dt`` is the sampling time expressed as a multiple of 1.025 sec frames. Example MSID names:: @@ -340,6 +418,7 @@ def get_msid_attrs(self, tstart, tstop, msid, msid_args): out = {'vals': vals[indexes], 'times': times, - 'bads': np.zeros(len(times), dtype=bool)} + 'bads': np.zeros(len(times), dtype=bool), + 'unit': None} return out diff --git a/Ska/engarchive/derived/mups_valve.py b/Ska/engarchive/derived/mups_valve.py index 94460f7c..d6f49008 100644 --- a/Ska/engarchive/derived/mups_valve.py +++ b/Ska/engarchive/derived/mups_valve.py @@ -82,11 +82,6 @@ def c2f(degc): return degf -def f2c(degf): - degc = (degf - 32) / 1.8 - return degc - - # Define MUPS valve thermistor point pair calibration table (from TDB). pp_counts = np.array([0, 27, 36, 44, 55, 70, 90, 118, 175, 195, 210, 219, 226, 231, 235, 255]) pp_temps = np.array([369.5305, 263.32577, 239.03652, 222.30608, 203.6944, 183.2642, 161.0796, @@ -186,8 +181,7 @@ def select_using_model(data1, data2, model, dt_thresh, out, source): source[ii] = 0 -def fetch_clean_msid(msid, start, stop=None, dt_thresh=5.0, median=7, - model_spec=None, unit='degf'): +def fetch_clean_msid(msid, start, stop=None, dt_thresh=5.0, median=7, model_spec=None): """Fetch a cleaned version of telemetry for ``msid``. If not supplied the model spec will be downloaded from github using this URL: @@ -212,7 +206,6 @@ def fetch_clean_msid(msid, start, stop=None, dt_thresh=5.0, median=7, :param dt_thresh: tolerance for matching model to data in degF (default=5 degF) :param median: length of median filter (default=7, use 0 to disable) :param model_spec: file name or URL containing relevant xija model spec - :param unit: temperature unit: 'degf' (default) or 'degc' :returns: fetch.Msid object """ @@ -245,21 +238,16 @@ def fetch_clean_msid(msid, start, stop=None, dt_thresh=5.0, median=7, source = np.zeros(len(dat.vals), dtype=int) select_using_model(t_obs, t_corr, t_model, dt_thresh=dt_thresh, out=out, source=source) - if unit not in ('degf', 'degc'): - raise ValueError('unit must be either degf or degc') - - convert = (lambda x: x) if unit == 'degf' else f2c - - dat.vals_raw = convert(dat.vals) + dat.vals_raw = dat.vals - dat.vals = convert(out) - dat.vals_nan = convert(out).copy() + dat.vals = out + dat.vals_nan = out.copy() dat.bads = (source == 0) dat.vals_nan[dat.bads] = np.nan dat.source = source - dat.vals_corr = convert(t_corr) - dat.vals_model = convert(t_model) + dat.vals_corr = t_corr + dat.vals_model = t_model return dat diff --git a/Ska/engarchive/fetch.py b/Ska/engarchive/fetch.py index faddd66f..16ac817e 100644 --- a/Ska/engarchive/fetch.py +++ b/Ska/engarchive/fetch.py @@ -631,23 +631,28 @@ def _get_data(self): self._get_msid_data_from_maude(*args) def _get_comp_data(self, comp_cls): - logger.info(f'Getting comp data for {self.msid}') + logger.info(f'Getting computed values for {self.msid}') # Do computation. This returns a dict of MSID attribute values. - dat = comp_cls()(self.tstart, self.tstop, self.msid, self.stat) + attrs = comp_cls(self.units['system'])(self.tstart, self.tstop, self.msid, self.stat) # Allow upstream class to be a bit sloppy on times and include samples # outside the time range. This can happen with classes that inherit # from DerivedParameter. - ok = (dat['times'] >= self.tstart) & (dat['times'] <= self.tstop) + ok = (attrs['times'] >= self.tstart) & (attrs['times'] <= self.tstop) all_ok = np.all(ok) + # List of "colnames", which is the ndarray attributes. There can be + # non-ndarray attributes that get returned, including typically 'unit'. + self.colnames = [attr for attr, val in attrs.items() + if (isinstance(val, np.ndarray) + and len(val) == len(attrs['times']))] + # Apply attributes to self - self.colnames = list(dat) - for attr, val in dat.items(): + for attr, val in attrs.items(): if (not all_ok and isinstance(val, np.ndarray) - and len(val) == len(dat['times'])): + and len(val) == len(attrs['times'])): val = val[ok] setattr(self, attr, val) diff --git a/Ska/engarchive/tests/test_comps.py b/Ska/engarchive/tests/test_comps.py index 3752d1a8..a6fda013 100644 --- a/Ska/engarchive/tests/test_comps.py +++ b/Ska/engarchive/tests/test_comps.py @@ -5,7 +5,7 @@ import numpy as np import pytest -from .. import fetch_eng, fetch +from .. import fetch_eng, fetch_sci, fetch as fetch_cxc from ..derived.base import DerivedParameter from ..derived.comps import ComputedMsid @@ -29,7 +29,8 @@ def get_msid_attrs(self, tstart, tstop, msid, msid_args): out = {'vals': dat.vals, 'bads': dat.bads, - 'times': dat.times} + 'times': dat.times, + 'unit': dat.unit} return out @@ -38,11 +39,12 @@ class Comp_Val_Plus_Five(ComputedMsid): msid_match = r'comp_(\w+)_plus_five' def get_msid_attrs(self, tstart, tstop, msid, msid_args): - dat = self.fetch_eng.MSID(msid_args[0], tstart, tstop) + dat = self.fetch_sys.MSID(msid_args[0], tstart, tstop) out = {'vals': dat.vals + 5, 'bads': dat.bads, - 'times': dat.times} + 'times': dat.times, + 'unit': dat.unit} return out @@ -71,7 +73,8 @@ def get_msid_attrs(self, tstart, tstop, msid, msid_args): out = {'vals': css1_npm_sun, 'times': msids.times, - 'bads': msids.bads} + 'bads': msids.bads, + 'unit': None} return out @@ -82,53 +85,74 @@ def test_comp_from_derived_parameter(): dat1 = fetch_eng.Msid('comp_css1_npm_sun', '2020:001', '2020:010') dat2 = fetch_eng.Msid('dp_css1_npm_sun', '2020:001', '2020:010') - for attr in ('vals', 'times', 'bads'): + for attr in ('vals', 'times', 'bads', 'unit'): assert np.all(getattr(dat1, attr) == getattr(dat2, attr)) -def test_simple_comp(): - dat1 = fetch_eng.Msid('tephin', '2020:001', '2020:010') - dat2 = fetch_eng.Msid('comp_tephin_plus_five', '2020:001', '2020:010') +@pytest.mark.parametrize('fetch,unit', [(fetch_eng, 'DEGF'), + (fetch_sci, 'DEGC'), + (fetch_cxc, 'K')]) +def test_simple_comp(fetch, unit): + dat1 = fetch.Msid('tephin', '2020:001', '2020:010') + dat2 = fetch.Msid('comp_tephin_plus_five', '2020:001', '2020:010') assert np.all(dat1.vals + 5 == dat2.vals) assert np.all(dat1.times == dat2.times) assert np.all(dat1.bads == dat2.bads) + assert dat1.unit == dat2.unit + assert dat2.unit == unit @pytest.mark.skipif("not HAS_MAUDE") -def test_simple_comp_with_maude(): +@pytest.mark.parametrize('fetch,unit', [(fetch_eng, 'DEGF'), + (fetch_sci, 'DEGC'), + (fetch_cxc, 'K')]) +def test_simple_comp_with_maude(fetch, unit): with fetch.data_source('maude'): - dat1 = fetch_eng.Msid('tephin', '2020:001', '2020:003') - dat2 = fetch_eng.Msid('comp_tephin_plus_five', '2020:001', '2020:003') + dat1 = fetch.Msid('tephin', '2020:001', '2020:003') + dat2 = fetch.Msid('comp_tephin_plus_five', '2020:001', '2020:003') assert np.all(dat1.vals + 5 == dat2.vals) assert np.all(dat1.times == dat2.times) assert np.all(dat1.bads == dat2.bads) + assert dat1.unit == dat2.unit + assert dat2.unit == unit def test_mups_valve(): colnames = ['vals', 'times', 'bads', 'vals_raw', 'vals_nan', 'vals_corr', 'vals_model', 'source'] - dat = fetch.MSID('PM2THV1T_clean', '2020:001', '2020:010') + dat = fetch_eng.MSID('PM2THV1T_clean', '2020:001', '2020:010') + assert dat.unit == 'DEGF' assert len(dat.vals) == 36661 - assert np.count_nonzero(dat.source != 0) == 34499 + ok = dat.source != 0 + # Temps are reasonable for degF + assert np.all((dat.vals[ok] > 55) & (dat.vals[ok] < 220)) + assert np.count_nonzero(ok) == 34499 assert dat.colnames == colnames for attr in colnames: assert len(dat.vals) == len(getattr(dat, attr)) - dat = fetch.Msid('PM2THV1T_clean', '2020:001', '2020:010') + dat = fetch_sci.Msid('PM2THV1T_clean', '2020:001', '2020:010') + assert dat.unit == 'DEGC' + ok = dat.source != 0 + # Temps are reasonable for degC + assert np.all((dat.vals[ok] > 10) & (dat.vals[ok] < 110)) assert len(dat.vals) == 34499 # Some bad values assert dat.colnames == colnames for attr in colnames: if attr != 'bads': assert len(dat.vals) == len(getattr(dat, attr)) - dat = fetch.MSID('PM1THV2T_clean', '2020:001', '2020:010') + dat = fetch_cxc.MSID('PM1THV2T_clean', '2020:001', '2020:010') + ok = dat.source != 0 + # Temps are reasonable for K + assert np.all((dat.vals[ok] > 280) & (dat.vals[ok] < 380)) assert len(dat.vals) == 36661 # Same as PM2THV1T assert dat.colnames == colnames for attr in colnames: assert len(dat.vals) == len(getattr(dat, attr)) - dat = fetch.Msid('pm1thv2t_clean', '2020:001', '2020:010') + dat = fetch_eng.Msid('pm1thv2t_clean', '2020:001', '2020:010') assert len(dat.vals) == 36240 # Some bad values assert len(dat.source) == 36240 # Filtering applies to sources assert dat.colnames == colnames @@ -139,15 +163,16 @@ def test_mups_valve(): def test_cmd_states(): start, stop = '2020:002:08:00:00', '2020:002:10:00:00' - dat = fetch.Msid('cmd_state_pitch_1000', start, stop) + dat = fetch_eng.Msid('cmd_state_pitch_1000', start, stop) exp_vals = np.array([55.99128956, 55.8747053, 55.8747053, 90.66266599, 159.06945155, 173.11528258, 173.11528258, 173.11528258]) assert np.allclose(dat.vals, exp_vals) assert type(dat.vals) is np.ndarray assert np.allclose(np.diff(dat.times), 1025.0) assert not np.any(dat.bads) + assert dat.unit is None - dat = fetch.Msid('cmd_state_pcad_mode_1000', start, stop) + dat = fetch_eng.Msid('cmd_state_pcad_mode_1000', start, stop) exp_vals = np.array(['NPNT', 'NPNT', 'NPNT', 'NMAN', 'NMAN', 'NPNT', 'NPNT', 'NPNT']) assert np.all(dat.vals == exp_vals) assert type(dat.vals) is np.ndarray @@ -158,8 +183,8 @@ def test_cmd_states(): def test_stats(stat): start, stop = '2020:001', '2020:010' - dat = fetch.Msid('pitch', start, stop, stat=stat) - datc = fetch.Msid('passthru_pitch', start, stop, stat=stat) + dat = fetch_eng.Msid('pitch', start, stop, stat=stat) + datc = fetch_eng.Msid('passthru_pitch', start, stop, stat=stat) for attr in datc.colnames: val = getattr(dat, attr) diff --git a/Ska/engarchive/units.py b/Ska/engarchive/units.py index 221b5ad6..9ec8d157 100644 --- a/Ska/engarchive/units.py +++ b/Ska/engarchive/units.py @@ -93,6 +93,13 @@ def F_to_C(vals, delta_val=False): return (vals - 32.0) / 1.8 +def C_to_F(vals, delta_val=False): + if delta_val: + return vals * 1.8 + else: + return vals * 1.8 + 32 + + def C_to_K(vals, delta_val=False): if delta_val: return vals @@ -171,6 +178,7 @@ def divide(scale_factor, decimals=None): ('mm', 'TSCSTEP'): mult(1.0 / 0.00251431530156, decimals=3), ('mm', 'FASTEP'): mm_to_FASTEP, ('PWM', 'PWMSTEP'): mult(16), + ('DEGC', 'DEGF'): C_to_F, # Eng units to CXC or Sci ('DEGC', 'K'): C_to_K, diff --git a/doc/pseudo_msids.rst b/doc/pseudo_msids.rst index b07c21a3..503770f6 100644 --- a/doc/pseudo_msids.rst +++ b/doc/pseudo_msids.rst @@ -720,8 +720,8 @@ requirements for defining a computed MSID:: def get_msid_attrs(self, tstart, tstop, msid, msid_args): """ Get attributes for computed MSID: ``vals``, ``bads``, ``times``, - ``raw_vals``, and ``offset``. The first three must always be - provided. + ``unit``, ``raw_vals``, and ``offset``. The first four must always + be provided. :param tstart: start time (CXC secs) :param tstop: stop time (CXC secs) @@ -733,23 +733,56 @@ requirements for defining a computed MSID:: msid msid_args[0] offset = int(msid_args[1]) - # Get the raw telemetry value in engineering units - dat = self.fetch_eng.MSID(msid, tstart, tstop) + # Get the raw telemetry value in user-requested unit system + dat = self.fetch_sys.MSID(msid, tstart, tstop) # Do the computation vals = dat.vals + offset - # Return a dict with at least `vals`, `times` and `bads`. + # Return a dict with at least `vals`, `times`, `bads`, and `unit`. # Additional attributes are allowed and will be set on the # final MSID object. out = {'vals': vals, 'bads': dat.bads, 'times': dat.times, + 'unit': dat.unit, 'vals_raw': dat.vals, # Provide original values without offset 'offset': offset # Provide the offset for reference } return out +Units +^^^^^ + +Computed MSIDs should support units where applicable. In the example above +this was done by using ``self.fetch_sys`` in order to get the original data +in the user-requested system. In other words, if the user did a call +``dat = fetch_sci.Msid('msid_plus_8', '2010:001', '2010:002')``, then +``self.fetch_sys`` would translate to ``self.fetch_sci``. + +In some cases the unit handling may require additional specification. This +can happen if the computation needs to be done in a particular unit, as is +the case for the built-in +:class:`~Ska.engarchive.derived.comps.Comp_MUPS_Valve_Temp_Clean` class. +Here the class must define an additional ``units`` attribute with the +following structure:: + + units = { + # Unit system for attrs from get_msid_attrs(), one of 'eng', 'sci', 'cxc' + 'internal_system': 'eng', + + # Units for eng, sci, cxc systems + 'eng': 'DEGF', + 'sci': 'DEGC', + 'cxc': 'K', + + # Attributes that need conversion. At least `vals` but maybe others. + 'convert_attrs': ['vals'] + } + +The specified units must all be convertable using functions defined in the +``converters`` dict in the ``Ska.engarchive.units`` module. + Built-in computed MSIDs and API ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ From 77cff94b5a76fac49ba43c533700e170c9aa2968 Mon Sep 17 00:00:00 2001 From: Tom Aldcroft Date: Thu, 9 Apr 2020 06:56:26 -0400 Subject: [PATCH 19/19] Fix flake8 --- Ska/engarchive/derived/comps.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Ska/engarchive/derived/comps.py b/Ska/engarchive/derived/comps.py index 80c4ddf1..0d2a6f2c 100644 --- a/Ska/engarchive/derived/comps.py +++ b/Ska/engarchive/derived/comps.py @@ -7,7 +7,7 @@ - Commanded states 'cmd_state__
' for any kadi commanded state value. See: https://nbviewer.jupyter.org/urls/cxc.harvard.edu/mta/ASPECT/ipynb/misc/DAWG-mups-valve-xija-filtering.ipynb -""" +""" # noqa import re @@ -378,7 +378,7 @@ class Comp_KadiCommandState(ComputedMsid): """Computed MSID for kadi dynamic commanded states. The MSID here takes the form ``cmd_state__
`` where: - + * ``state_key`` is a valid commanded state key such as ``pitch`` or ``pcad_mode`` or ``acisfp_temp``. * ``dt`` is the sampling time expressed as a multiple of 1.025 sec