Skip to content
New issue

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

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

Already on GitHub? Sign in to your account

Let 2-CCD (S3+FI) observations with less than 300 counts in the sequence go to -109 C #49

Merged
merged 16 commits into from
Mar 22, 2022
Merged
69 changes: 43 additions & 26 deletions acis_thermal_check/acis_obs.py
Original file line number Diff line number Diff line change
@@ -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):
Expand Down Expand Up @@ -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}
Expand All @@ -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
Expand Down Expand Up @@ -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:
Expand All @@ -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
Expand Down Expand Up @@ -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)
Expand Down
38 changes: 21 additions & 17 deletions acis_thermal_check/apps/acisfp_check.py
Original file line number Diff line number Diff line change
Expand Up @@ -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 = {}
Expand Down Expand Up @@ -192,15 +185,18 @@ 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')
plots[name].ax.errorbar([0.0, 0.0], [1.0, 1.0], xerr=1.0,
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),
Expand Down Expand Up @@ -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.
Expand All @@ -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 = {}
Expand All @@ -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",
Expand All @@ -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)
Expand All @@ -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)
Expand All @@ -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)
Expand All @@ -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

Expand Down Expand Up @@ -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
Expand All @@ -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
Expand Down
13 changes: 8 additions & 5 deletions acis_thermal_check/main.py
Original file line number Diff line number Diff line change
Expand Up @@ -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)
Expand Down Expand Up @@ -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.
Expand All @@ -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
Expand Down
10 changes: 8 additions & 2 deletions acis_thermal_check/templates/index_template.rst
Original file line number Diff line number Diff line change
Expand Up @@ -34,9 +34,15 @@ States `<states.dat>`_
{% 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() %}
Expand Down
26 changes: 26 additions & 0 deletions acis_thermal_check/tests/acisfp/answers/FEB2122A_hot_acis.json
Original file line number Diff line number Diff line change
@@ -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"
]
]
38 changes: 38 additions & 0 deletions acis_thermal_check/tests/acisfp/answers/FEB2822A_hot_acis.json
Original file line number Diff line number Diff line change
@@ -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"
]
]
Loading