diff --git a/acis_thermal_check/acis_obs.py b/acis_thermal_check/acis_obs.py index 5396c0c6..e030e69a 100755 --- a/acis_thermal_check/acis_obs.py +++ b/acis_thermal_check/acis_obs.py @@ -1,5 +1,6 @@ import numpy as np from acis_thermal_check.utils import mylog +from kadi.commands.states import decode_power def who_in_fp(simpos=80655): @@ -86,21 +87,11 @@ def fetch_ocat_data(obsid_list): if got_table: tab = ascii.read(resp.text, header_start=0, data_start=2) tab.sort("OBSID") - # We figure out the CCD count from the table by finding out - # which ccds were on, optional, or dropped, and then - # subtracting off the dropped chip count entry in the table - ccd_count = np.zeros(tab["S3"].size, dtype='int') - for a, r in zip(["I", "S"], [range(4), range(6)]): - for i in r: - ccd = np.ma.filled(tab[f"{a}{i}"].data) - ccd_count += (ccd == "Y").astype('int') - ccd_count += (ccd == "D").astype('int') - ccd_count += np.char.startswith(ccd, "O").astype('int') - ccd_count -= tab["DROPPED_CHIP_CNT"].data.astype('int') # Now we have to find all of the obsids in each sequence and then # compute the complete exposure for each sequence - seq_nums = list(tab["SEQ_NUM"].data.astype("str")) - seq_num_list = ",".join([seq_num for seq_num in seq_nums if seq_num != " "]) + seq_nums = np.unique([str(sn) for sn in tab["SEQ_NUM"].data.astype("str") + if sn != " "]) + seq_num_list = ",".join(seq_nums) obsids = tab["OBSID"].data.astype("int") cnt_rate = tab["EST_CNT_RATE"].data.astype("float64") params = {"seqNum": seq_num_list} @@ -121,15 +112,18 @@ def fetch_ocat_data(obsid_list): return None tab_seq = ascii.read(resp.text, header_start=0, data_start=2) app_exp = np.zeros_like(cnt_rate) - for row in tab_seq: - i = seq_nums.index(str(row["SEQ_NUM"])) - app_exp[i] += np.float64(row["APP_EXP"]) + seq_nums = tab_seq["SEQ_NUM"].data.astype("str") + for i, row in enumerate(tab): + sn = str(row["SEQ_NUM"]) + if sn == " ": + continue + j = np.where(sn == seq_nums)[0] + app_exp[i] += np.float64(tab_seq["APP_EXP"][j]).sum() app_exp *= 1000.0 table_dict = {"obsid": np.array(obsids), "grating": tab["GRAT"].data, - "ccd_count": ccd_count, - "S3": np.ma.filled(tab["S3"].data), - "num_counts": cnt_rate*app_exp} + "cnt_rate": cnt_rate, + "app_exp": app_exp} else: # We weren't able to get a valid table for some reason, so # we cannot check for -109 data, but we proceed with the @@ -210,12 +204,18 @@ def find_obsid_intervals(cmd_states, load_start): datestart = eachstate['datestart'] tstart = eachstate['tstart'] + # Process the power command which turns things on + if pow_cmd.startswith("WSPOW") and pow_cmd != 'WSPOW00000' \ + and firstpow: + ccds = decode_power(pow_cmd)['ccds'].replace(" ", ",")[:-1] + # Process the first XTZ0000005 line you see if pow_cmd in ['XTZ0000005', 'XCZ0000005'] and \ (xtztime is None and firstpow): xtztime = eachstate['tstart'] # MUST fix the instrument now instrument = who_in_fp(eachstate['simpos']) + ccd_count = eachstate["ccd_count"] # Process the first AA00000000 line you see if pow_cmd == 'AA00000000' and firstpow: @@ -234,9 +234,11 @@ def find_obsid_intervals(cmd_states, load_start): "datestop": datestop, "tstart": tstart, "tstop": tstop, + "ccds": ccds, "start_science": xtztime, "obsid": eachstate['obsid'], - "instrument": instrument} + "instrument": instrument, + "ccd_count": ccd_count} obsid_interval_list.append(obsid_dict) # now clear out the data values @@ -295,15 +297,30 @@ def acis_filter(obsid_interval_list): acis_i = [] cold_ecs = [] + mylog.debug("OBSID\tCNT_RATE\tAPP_EXP\tNUM_CTS\tGRATING\tCCDS") for eachobs in obsid_interval_list: + # First we check that we got ocat data using "grating" + hot_acis = False if "grating" in eachobs: - hetg = eachobs["grating"] == "HETG" - s3_only = eachobs["S3"] == "Y" and eachobs["ccd_count"] == 1 - hot_acis = hetg or (eachobs["num_counts"] < 300.0 and s3_only) - else: - hot_acis = False + eachobs["num_counts"] = int(eachobs["cnt_rate"] * eachobs["app_exp"]) + # First check to see if this is an S3 observation + mylog.debug(f"{eachobs['obsid']}\t{eachobs['cnt_rate']}\t" + f"{eachobs['app_exp']*1.0e-3}\t{eachobs['num_counts']}\t" + f"{eachobs['grating']}\t{eachobs['ccds']}") + if eachobs["ccd_count"] <= 2: + # S3 with low counts + low_ct_s3 = eachobs["num_counts"] < 300 and "S3" in eachobs["ccds"] + # Is there another chip on? Make sure it's not S1 + if eachobs["ccd_count"] == 2: + low_ct_s3 &= "S1" not in eachobs["ccds"] + else: + # All higher ccd counts are invalid + low_ct_s3 = False + # Also check grating status + hot_acis = (eachobs["grating"] == "HETG") or low_ct_s3 + eachobs['hot_acis'] = hot_acis if hot_acis: - acis_hot.append(eachobs) + acis_hot.append(eachobs) else: if eachobs["instrument"] == "ACIS-S": acis_s.append(eachobs) diff --git a/acis_thermal_check/apps/acisfp_check.py b/acis_thermal_check/apps/acisfp_check.py index 425670b1..ed9fa9b1 100755 --- a/acis_thermal_check/apps/acisfp_check.py +++ b/acis_thermal_check/apps/acisfp_check.py @@ -137,13 +137,6 @@ def make_prediction_plots(self, outdir, states, temps, load_start): those from the commanded states data structure called "states" """ - # extract the OBSID's from the commanded states. NOTE: this contains all - # observations including ECS runs and HRC observations - observation_intervals = find_obsid_intervals(states, load_start) - - # Filter out any HRC science observations BUT keep ACIS ECS observations - self.acis_and_ecs_obs = hrc_science_obs_filter(observation_intervals) - # create an empty dictionary called plots to contain the returned # figures, axes 1 and axes 2 of the PredictPlot class plots = {} @@ -192,6 +185,9 @@ def make_prediction_plots(self, outdir, states, temps, load_start): textypos[i], fontsize[i], plot_start) # These next lines are dummies so we can get the obsids in the legend + plots[name].ax.errorbar([0.0, 0.0], [1.0, 1.0], xerr=1.0, + lw=2, xlolims=True, color='blue', + capsize=4, capthick=2, label='Cold ECS') plots[name].ax.errorbar([0.0, 0.0], [1.0, 1.0], xerr=1.0, lw=2, xlolims=True, color='red', capsize=4, capthick=2, label='ACIS-I') @@ -199,8 +195,8 @@ def make_prediction_plots(self, outdir, states, temps, load_start): lw=2, xlolims=True, color='green', capsize=4, capthick=2, label='ACIS-S') plots[name].ax.errorbar([0.0, 0.0], [1.0, 1.0], xerr=1.0, - lw=2, xlolims=True, color='blue', - capsize=4, capthick=2, label='ECS') + lw=2, xlolims=True, color='saddlebrown', + capsize=4, capthick=2, label='Hot ACIS-S') # Make the legend on the temperature plot plots[name].ax.legend(bbox_to_anchor=(0.15, 0.99), @@ -230,7 +226,7 @@ def make_prediction_plots(self, outdir, states, temps, load_start): return plots - def make_prediction_viols(self, temps, load_start): + def make_prediction_viols(self, temps, states, load_start): """ Find limit violations where predicted temperature is above the red minus margin. @@ -252,9 +248,16 @@ def make_prediction_viols(self, temps, load_start): - science_viols """ + # extract the OBSID's from the commanded states. NOTE: this contains all + # observations including ECS runs and HRC observations + observation_intervals = find_obsid_intervals(states, load_start) + + # Filter out any HRC science observations BUT keep ACIS ECS observations + self.acis_and_ecs_obs = hrc_science_obs_filter(observation_intervals) + times = self.predict_model.times - mylog.info(f"\nMAKE VIOLS Checking for limit violations in " + mylog.info(f"MAKE VIOLS Checking for limit violations in " f"{len(self.acis_and_ecs_obs)} total science observations") viols = {} @@ -278,7 +281,7 @@ def make_prediction_viols(self, temps, load_start): acis_hot_limit = self.limits["acis_hot"].value cold_ecs_limit = self.limits["cold_ecs"].value - mylog.info(f'\n\nACIS-I Science ({acis_i_limit} C) violations') + mylog.info(f'ACIS-I Science ({acis_i_limit} C) violations') # Create the violation data structure. acis_i_viols = self.search_obsids_for_viols("ACIS-I", @@ -293,7 +296,7 @@ def make_prediction_viols(self, temps, load_start): # science run. These are load killers # ------------------------------------------------------------ # - mylog.info(f'\n\nACIS-S Science ({acis_s_limit} C) violations') + mylog.info(f'ACIS-S Science ({acis_s_limit} C) violations') acis_s_viols = self.search_obsids_for_viols("ACIS-S", acis_s_limit, ACIS_S_obs, temp, times, load_start) @@ -306,7 +309,7 @@ def make_prediction_viols(self, temps, load_start): # science run which can run hot. These are load killers # ------------------------------------------------------------ # - mylog.info(f'\n\nACIS-S Science ({acis_hot_limit} C) violations') + mylog.info(f'ACIS-S Science ({acis_hot_limit} C) violations') acis_hot_viols = self.search_obsids_for_viols("Hot ACIS-S", acis_hot_limit, ACIS_hot_obs, temp, times, load_start) @@ -317,7 +320,7 @@ def make_prediction_viols(self, temps, load_start): # ------------------------------------------------------------ # Science Orbit ECS -119.5 violations; -119.5 violation check # ------------------------------------------------------------ - mylog.info(f'\n\nScience Orbit ECS ({cold_ecs_limit} C) violations') + mylog.info(f'Science Orbit ECS ({cold_ecs_limit} C) violations') ecs_viols = self.search_obsids_for_viols("Science Orbit ECS", cold_ecs_limit, sci_ecs_obs, temp, times, load_start) @@ -328,7 +331,7 @@ def make_prediction_viols(self, temps, load_start): # Store all obsids which can go to -109 C for obs in ACIS_hot_obs: - self.acis_hot_obs.append(obs["obsid"]) + self.acis_hot_obs.append(obs) return viols @@ -416,6 +419,7 @@ def draw_obsids(obs_list, plots, msid, ypos, endcapstart, endcapstop, obsid = eachobservation['obsid'] in_fp = eachobservation['instrument'] + hot_acis = eachobservation['hot_acis'] if obsid > 60000: # ECS observations during the science orbit are colored blue @@ -426,7 +430,7 @@ def draw_obsids(obs_list, plots, msid, ypos, endcapstart, endcapstop, if in_fp == "ACIS-I": color = 'red' else: - color = 'green' + color = 'saddlebrown' if hot_acis else 'green' obsid_txt = str(obsid) # If this is an ECS measurement in the science orbit mark diff --git a/acis_thermal_check/main.py b/acis_thermal_check/main.py index 2be64a94..5961869d 100644 --- a/acis_thermal_check/main.py +++ b/acis_thermal_check/main.py @@ -339,12 +339,13 @@ def make_week_predict(self, tstart, tstop, tlm, T_init, model_spec, plt.rc("grid", linewidth=1.5) temps = {self.name: model.comp[self.msid].mvals} - # make_prediction_plots runs the validation of the model - # against previous telemetry - plots = self.make_prediction_plots(outdir, states, temps, tstart) # make_prediction_viols determines the violations and prints them out - viols = self.make_prediction_viols(temps, tstart) + viols = self.make_prediction_viols(temps, states, tstart) + + # make_prediction_plots runs the validation of the model + # against previous telemetry + plots = self.make_prediction_plots(outdir, states, temps, tstart) # write_states writes the commanded states to states.dat self.write_states(outdir, states) @@ -504,7 +505,7 @@ def _make_prediction_viols(self, times, temp, load_start, limit, lim_name, return viols - def make_prediction_viols(self, temps, load_start): + def make_prediction_viols(self, temps, states, load_start): """ Find limit violations where predicted temperature is above the specified limits. @@ -513,6 +514,8 @@ def make_prediction_viols(self, temps, load_start): ---------- temps : dict of NumPy arrays NumPy arrays corresponding to the modeled temperatures + states : NumPy record array + Commanded states load_start : float The start time of the load, used so that we only report violations for times later than this time for the model diff --git a/acis_thermal_check/templates/index_template.rst b/acis_thermal_check/templates/index_template.rst index 0817a70e..9293fecd 100644 --- a/acis_thermal_check/templates/index_template.rst +++ b/acis_thermal_check/templates/index_template.rst @@ -34,9 +34,15 @@ States ``_ {% if proc.msid == "FPTEMP" %} "Hot" ACIS Observations (-109 C limit) -------------------------------------- -{% for obsid in acis_hot_obs %} -| {{ obsid }} +{% if acis_hot_obs|length > 0 %} +===== ================= ================== ======= +Obsid CCDs # of counts in seq Grating +===== ================= ================== ======= +{% for eachobs in acis_hot_obs %} +{{ eachobs.obsid }} {{"{0: <17}".format(eachobs.ccds)}} {{"{0: <18}".format(eachobs.num_counts)}} {{eachobs.grating}} {% endfor %} +===== ================= ================== ======= +{% endif %} {% endif %} {% for key in viols.keys() %} diff --git a/acis_thermal_check/tests/acisfp/answers/FEB2122A_hot_acis.json b/acis_thermal_check/tests/acisfp/answers/FEB2122A_hot_acis.json new file mode 100644 index 00000000..f1c51583 --- /dev/null +++ b/acis_thermal_check/tests/acisfp/answers/FEB2122A_hot_acis.json @@ -0,0 +1,26 @@ +[ + [ + "24244", + "S0,S1,S2,S3,S4,S5", + "765000", + "HETG" + ], + [ + "25556", + "S3", + "3", + "NONE" + ], + [ + "26335", + "S1,S2,S3,S4", + "765000", + "HETG" + ], + [ + "26334", + "S3", + "29", + "NONE" + ] +] \ No newline at end of file diff --git a/acis_thermal_check/tests/acisfp/answers/FEB2822A_hot_acis.json b/acis_thermal_check/tests/acisfp/answers/FEB2822A_hot_acis.json new file mode 100644 index 00000000..fb0155c2 --- /dev/null +++ b/acis_thermal_check/tests/acisfp/answers/FEB2822A_hot_acis.json @@ -0,0 +1,38 @@ +[ + [ + "25333", + "S3", + "5", + "NONE" + ], + [ + "23653", + "S2,S3", + "260", + "NONE" + ], + [ + "24705", + "S2,S3", + "500000", + "HETG" + ], + [ + "26344", + "S2,S3", + "500000", + "HETG" + ], + [ + "26345", + "S2,S3", + "500000", + "HETG" + ], + [ + "25291", + "S3", + "61", + "NONE" + ] +] \ No newline at end of file diff --git a/acis_thermal_check/tests/acisfp/answers/MAR0722A_hot_acis.json b/acis_thermal_check/tests/acisfp/answers/MAR0722A_hot_acis.json new file mode 100644 index 00000000..b0df4501 --- /dev/null +++ b/acis_thermal_check/tests/acisfp/answers/MAR0722A_hot_acis.json @@ -0,0 +1,32 @@ +[ + [ + "24283", + "S3", + "50", + "NONE" + ], + [ + "25404", + "S1,S2,S3,S4", + "50000", + "HETG" + ], + [ + "26356", + "S0,S1,S2,S3,S4,S5", + "50000", + "HETG" + ], + [ + "26352", + "S3", + "50", + "NONE" + ], + [ + "25266", + "S3", + "16", + "NONE" + ] +] \ No newline at end of file diff --git a/acis_thermal_check/tests/acisfp/test_acisfp_hot_acis.py b/acis_thermal_check/tests/acisfp/test_acisfp_hot_acis.py new file mode 100644 index 00000000..cb92f3b8 --- /dev/null +++ b/acis_thermal_check/tests/acisfp/test_acisfp_hot_acis.py @@ -0,0 +1,24 @@ +from acis_thermal_check.apps.acisfp_check import ACISFPCheck +from acis_thermal_check.tests.regression_testing import \ + RegressionTester, tests_path + + +def test_FEB2122_hot_acis(answer_store, test_root): + answer_data = tests_path / "acisfp/answers/FEB2122A_hot_acis.json" + fp_rt = RegressionTester(ACISFPCheck, test_root=test_root, sub_dir='hot_acis') + fp_rt.check_hot_acis_reporting("FEB2122A", answer_data, + answer_store=answer_store) + + +def test_FEB2822_hot_acis(answer_store, test_root): + answer_data = tests_path / "acisfp/answers/FEB2822A_hot_acis.json" + fp_rt = RegressionTester(ACISFPCheck, test_root=test_root, sub_dir='hot_acis') + fp_rt.check_hot_acis_reporting("FEB2822A", answer_data, + answer_store=answer_store) + + +def test_MAR0722_hot_acis(answer_store, test_root): + answer_data = tests_path / "acisfp/answers/MAR0722A_hot_acis.json" + fp_rt = RegressionTester(ACISFPCheck, test_root=test_root, sub_dir='hot_acis') + fp_rt.check_hot_acis_reporting("MAR0722A", answer_data, + answer_store=answer_store) diff --git a/acis_thermal_check/tests/regression_testing.py b/acis_thermal_check/tests/regression_testing.py index 37b4ad87..7d0844d0 100644 --- a/acis_thermal_check/tests/regression_testing.py +++ b/acis_thermal_check/tests/regression_testing.py @@ -415,3 +415,28 @@ def check_violation_reporting(self, load_week, viol_json, with open(viol_json, "w") as f: json.dump(viol_data, f, indent=4) + def check_hot_acis_reporting(self, load_week, hot_json, + answer_store=False): + import json + self.run_model(load_week) + out_dir = self.outdir / load_week / self.name + index_rst = out_dir / "index.rst" + hot_data = [] + with open(index_rst, 'r') as myfile: + read_hot_data = False + for i, line in enumerate(myfile.readlines()): + words = line.strip().split() + if line.startswith("Obsid"): + read_hot_data = True + elif read_hot_data: + if len(words) == 4 and not line.startswith("=="): + hot_data.append(words) + elif len(words) == 0: + read_hot_data = False + if answer_store: + with open(hot_json, "w") as f: + json.dump(hot_data, f, indent=4) + else: + with open(hot_json, "r") as f: + hot_data_stored = json.load(f) + assert_array_equal(hot_data, hot_data_stored)