Skip to content

Commit

Permalink
Implement full compliance with maude_decom API
Browse files Browse the repository at this point in the history
  • Loading branch information
@taldcroft
taldcroft committed Oct 4, 2024
1 parent 1f31698 commit afd867d
Show file tree
Hide file tree
Showing 2 changed files with 180 additions and 103 deletions.
195 changes: 128 additions & 67 deletions mica/archive/aca_l0.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -17,7 +17,7 @@
import Ska.File
import ska_dbi
import tables
from astropy.table import Table, vstack
from astropy.table import MaskedColumn, Table, vstack
from Chandra.Time import DateTime
from chandra_aca.aca_image import ACAImage
from cxotime import CxoTimeLike
Expand Down Expand Up @@ -62,6 +62,8 @@
("MNF", ">i4"),
("END_INTEG_TIME", ">f8"),
("INTEG", ">f4"),
("PIXTLM", "<U4"),
("BGDTYP", "<U4"),
("GLBSTAT", "|u1"),
("COMMCNT", "|u1"),
("COMMPROG", "|u1"),
Expand Down Expand Up @@ -119,6 +121,27 @@
cda_table="cda_aca0.h5",
)

# Names of bits in GLBSTAT and IMGSTAT fields, in order from most to least significant
GLBSTAT_NAMES = [
"HIGH_BGD",
"RAM_FAIL",
"ROM_FAIL",
"POWER_FAIL",
"CAL_FAIL",
"COMM_CHECKSUM_FAIL",
"RESET",
"SYNTAX_ERROR",
]

IMGSTAT_NAMES = [
"SAT_PIXEL",
"DEF_PIXEL",
"QUAD_BOUND",
"COMMON_COL",
"MULTI_STAR",
"ION_RAD",
]


def get_options():
parser = argparse.ArgumentParser(
Expand All @@ -145,95 +168,129 @@ def get_options():
return opt


def get_aca_images(start: CxoTimeLike, stop: CxoTimeLike, imgsize=None, columns=None):
def _extract_bits_as_uint8(arr: np.ndarray, bit_start: int, bit_stop:int) -> np.ndarray:
"""Extract UINT8's from a 2-D array of uint8 bit values into a new array.
Input array is Nx8 where N is the number of rows and 8 is the number of bits.
It must be in little-endian order with the most significant bit in the rightmost
position of the array.
"""
Get ACA images in a table from mica l0 archive for a given time range.
# This returns a shape of (N, 1) so we need to squeeze it.
out = np.packbits(arr[:, bit_start:bit_stop], axis=1, bitorder="little")
return out.squeeze()


def get_aca_images(
start: CxoTimeLike,
stop: CxoTimeLike,
):
"""
Get ACA images from mica L0 archive in same format as maude_decom.get_aca_images.
This gets a table of images from the mica L0 file archive that includes the same
output columns as ``chandra_aca.maude_decom.get_aca_images()``. The images are 8x8
and are centered in the 8x8 image, with masks indicating empty pixels (for 4x4 or
6x6 data).
This is a convenience method to get a table of images that should be very similar to
that from maude_decom.get_aca_images. These mica images are centered in the 8x8 images
via the centered_8x8=True option to aca_l0.get_slot_data to match the maude images.
The units of the image data are DN, so you need to multiply by 5.0 / INTEG to get
e-/sec.
This function returns additional columns that are not in the maude images: -
'HD3TLM*': uint8, raw 8-bit header-3 telemetry values - 'IMGSIZE': int32, 4, 6, or 8
Parameters
----------
start : CxoTimeLike
Start time of the time range.
stop : CxoTimeLike
Stop time of the time range.
imgsize : list of int, optional
List of integers of desired image sizes (passed to get_slot_data which defaults
to [4, 6, 8] if None).
columns : list of str, optional
List of desired columns in the ACA0 telemetry. If None defaults to a reasonable
minimal list ('TIME', 'QUALITY', 'MJF', 'MNF', 'IMGFUNC1', 'IMGSIZE',
'IMGRAW', 'IMGROW0', 'IMGCOL0', 'BGDAVG', 'INTEG').
Returns
-------
astropy.table.Table of ACA images and accompanying data.
Note that IMGRAW is renamed to IMG.
IMGNUM is included (same as slot number)
IMGROW0_8X8 and IMGCOL0_8X8 are added to the table to match maude_decom.
astropy.table.Table
ACA images and accompanying data.
"""
if imgsize is None:
imgsize = [4, 6, 8]

if columns is None:
columns = [
"TIME",
"QUALITY",
"MJF",
"MNF",
"IMGFUNC1",
"IMGRAW",
"IMGSIZE",
"IMGROW0",
"IMGCOL0",
"BGDAVG",
"INTEG",
]
else:
# Require that these columns are present
for col in ["TIME", "IMGRAW", "IMGROW0", "IMGCOL0"]:
if col not in columns:
columns.append(col)

img_data = []
slot_data_list = []
for slot in range(8):
slot_data_raw = get_slot_data(
start,
stop,
imgsize=imgsize,
slot=slot,
columns=columns,
centered_8x8=True,
)
# Convert from numpy array to astropy table
slot_data = Table(slot_data_raw)

# Rename and rejigger colums to match maude
IMGROW0_8X8 = slot_data["IMGROW0"].copy()
IMGCOL0_8X8 = slot_data["IMGCOL0"].copy()
IMGROW0_8X8[slot_data["IMGSIZE"] == 6] -= 1
IMGCOL0_8X8[slot_data["IMGSIZE"] == 6] -= 1
IMGROW0_8X8[slot_data["IMGSIZE"] == 4] -= 2
IMGCOL0_8X8[slot_data["IMGSIZE"] == 4] -= 2
slot_data["IMGROW0_8X8"] = IMGROW0_8X8
slot_data["IMGCOL0_8X8"] = IMGCOL0_8X8
slot_data.rename_column("IMGRAW", "IMG")
# Convert from numpy structured array to astropy table, keeping only good data
ok = slot_data_raw["QUALITY"] == 0
slot_data = Table(slot_data_raw[ok])


# Add slot number if there are any rows (if statement not needed after
# https://github.com/astropy/astropy/pull/17102).
# if len(slot_data) > 0:
# slot_data["IMGNUM"] = slot

for name in ["IMGFID", "IMGFUNC", "IMGNUM"]:
slot_data[name] = slot_data[f"{name}1"]
for ii in (1, 2, 3, 4):
name_ii = f"{name}{ii}"
if name_ii in slot_data.colnames:
del slot_data[name_ii]

slot_data_list.append(slot_data)

# Combine all slot data into one output table and sort by time and slot
out = vstack(slot_data_list)
out.sort(keys=["TIME", "IMGNUM"])

# Rename and rejigger colums to match chandra_aca.maude_decum.get_aca_images
is_6x6 = out["IMGSIZE"] == 6
is_4x4 = out["IMGSIZE"] == 4
for rc in ["ROW", "COL"]:
# IMGROW/COL0 are row/col of lower/left of 4x4, 6x6, or 8x8 image
# Make these new columns:
# IMGROW/COL_A1: row/col of pixel A1 of 4x4, 6x6, or 8x8 image
# IMGROW/COL0_8x8: row/col of lower/left of 8x8 image with smaller img centered
out[f"IMG{rc}_A1"] = out[f"IMG{rc}0"]
out[f"IMG{rc}_A1"][is_6x6] += 1
out[f"IMG{rc}0_8X8"] = out[f"IMG{rc}0"]
out[f"IMG{rc}0_8X8"][is_6x6] -= 1
out[f"IMG{rc}0_8X8"][is_4x4] -= 2

out["IMGTYPE"] = np.full(shape=len(out), fill_value=4) # 8x8
out["IMGTYPE"][is_6x6] = 1 # 6x6
out["IMGTYPE"][is_4x4] = 0 # 4x4

out.rename_column("IMGRAW", "IMG")
out.rename_column("BGDTYP", "AABGDTYP")
out.rename_column("PIXTLM", "AAPIXTLM")
del out["FILENAME"]
del out["QUALITY"]
out["VCDUCTR"] = out["MJF"] * 128 + out["MNF"]

# Split uint8 GLBSTAT and IMGSTAT into individual bits. The stat_bit_names are
# in MSB order so we need to reverse them below.
for stat_name, stat_bit_names in (
("GLBSTAT", GLBSTAT_NAMES),
("IMGSTAT", IMGSTAT_NAMES),
):
for bit, name in zip(range(8), reversed(stat_bit_names)):
out[name] = (out[stat_name] & (1 << bit)) > 0

# Add slot number if there are any rows
if len(slot_data) > 0:
slot_data["IMGNUM"] = slot
img_data.append(slot_data)
# Extract fields from the mica L0 8-bit COMMCNT value which is really three MSIDs.
# We need to use .data attribute of MaskedColumn to get a numpy masked array. The
# unpackbits function fails on MaskedColumn because it looks for a _mask attribute.
bits = np.unpackbits(out["COMMCNT"].data).reshape(-1, 8)
bits_rev = bits[:, ::-1] # MSB is on the right (index=7)
out["COMMCNT"] = _extract_bits_as_uint8(bits_rev, 2, 8)
out["COMMCNT_CHECKSUM_FAIL"] = bits_rev[:, 0] > 0
out["COMMCNT_SYNTAX_ERROR"] = bits_rev[:, 1] > 0

combined_data = vstack(img_data)
# Extract fields from the mica L0 8-bit COMMPROG value which is really two MSIDs
bits = np.unpackbits(out["COMMPROG"].data).reshape(-1, 8)
bits_rev = bits[:, ::-1]
out["COMMPROG"] = _extract_bits_as_uint8(bits_rev, 2, 8)
out["COMMPROG_REPEAT"] = _extract_bits_as_uint8(bits_rev, 0, 2)

# And sort if there are any rows
if len(combined_data) > 0:
combined_data.sort(keys=["TIME", "IMGNUM"])
return combined_data
return out


def get_slot_data(
Expand Down Expand Up @@ -310,6 +367,10 @@ def get_slot_data(
continue
if fname in chunk.dtype.names:
all_rows[fname][idx0:idx1] = chunk[fname]
elif fname == "PIXTLM":
all_rows[fname][idx0:idx1] = "ORIG"
elif fname == "BGDTYP":
all_rows[fname][idx0:idx1] = "FLAT"
if "IMGSIZE" in columns:
all_rows["IMGSIZE"][idx0:idx1] = f_imgsize
if "FILENAME" in columns:
Expand Down
88 changes: 52 additions & 36 deletions mica/archive/tests/test_aca_l0.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3,12 +3,13 @@

import os

import astropy.units as u
import numpy as np
import pytest
from astropy.table import Table
from chandra_aca import maude_decom
from cxotime import CxoTime
from Ska.Numpy import interpolate
from ska_numpy import interpolate

from mica.archive import aca_l0, asp_l1

Expand Down Expand Up @@ -39,47 +40,62 @@ def test_get_aca_images():
"""
Confirm mica and maude images sources agree for a small time range.
"""
tstart = CxoTime("2012:180").secs
imgs_table_maude = maude_decom.get_aca_images(tstart, tstart + 30)
imgs_table_maude.sort(keys=["TIME", "IMGNUM"])
imgs_table_mica = aca_l0.get_aca_images(tstart, tstart + 30)
assert len(imgs_table_maude) == len(imgs_table_mica)
for col in ["IMG", "IMGROW0_8X8", "IMGCOL0_8X8", "TIME"]:
assert np.allclose(
imgs_table_maude[col], imgs_table_mica[col], atol=0.1, rtol=0
)
for col in ["IMGNUM"]:
assert np.all(imgs_table_maude[col] == imgs_table_maude[col])
# Includes 4x4, 6x6, 8x8 images, commanding, NPNT, NMAN, and fids
start = CxoTime("2012:270:02:25:00")
stop = CxoTime("2012:270:02:46:00")
imgs_maude = maude_decom.get_aca_images(start, stop)
imgs_maude.sort(keys=["TIME", "IMGNUM"])
imgs_mica = aca_l0.get_aca_images(start, stop)

assert len(imgs_maude) == len(imgs_mica)

# Test data set includes commanding, searches, and all data types
assert set(imgs_mica["IMGSIZE"]) == {4, 6, 8}
assert set(imgs_mica["COMMPROG"]) == {0, 1, 3, 5, 7, 9, 11, 13, 15, 19, 21, 22}
assert set(imgs_mica["COMMCNT"]) == {0, 24, 13, 25}
assert set(imgs_mica["COMMPROG_REPEAT"]) == {0, 1}

for colname in imgs_maude.colnames:
if colname == "END_INTEG_TIME":
# Currently this is off by INTEG/2
continue
if imgs_maude[colname].dtype.kind == "f":
assert np.allclose(
imgs_mica[colname], imgs_maude[colname], rtol=0, atol=1e-3
)
else:
assert np.all(imgs_mica[colname] == imgs_maude[colname])

assert set(imgs_mica.colnames) - set(imgs_maude.colnames) == {
"HD3TLM63",
"HD3TLM77",
"HD3TLM66",
"HD3TLM74",
"HD3TLM64",
"HD3TLM76",
"IMGSIZE",
"HD3TLM65",
"HD3TLM72",
"HD3TLM73",
"HD3TLM75",
"HD3TLM62",
"HD3TLM67",
}


@pytest.mark.skipif(not has_l0_2012_archive, reason="Test requires 2012 L0 archive")
def test_get_aca_images_empty():
"""
Confirm that get_aca_images returns an empty table when no images are found
Confirm that get_aca_images returns a zero-length table when no images are found
"""
tstart = CxoTime("2012:180").secs
# select on a non-sensical imgsize to ensure no images are found
imgs_table_mica = aca_l0.get_aca_images(tstart, tstart + 30, imgsize=[7])

# Assert it has no rows
assert len(imgs_table_mica) == 0
# assert that it is an astropy table
assert type(imgs_table_mica) is Table
assert imgs_table_mica.colnames == [
"TIME",
"QUALITY",
"MJF",
"MNF",
"INTEG",
"IMGFUNC1",
"IMGROW0",
"IMGCOL0",
"BGDAVG",
"IMG",
"IMGSIZE",
"IMGROW0_8X8",
"IMGCOL0_8X8",
]
start = CxoTime("2012:270")
imgs_mica = aca_l0.get_aca_images(start, start)
imgs_maude = maude_decom.get_aca_images(start, start + 10 * u.s)

# zero-length table with the required columns to match maude_decom
assert len(imgs_mica) == 0
assert type(imgs_mica) is Table
assert set(imgs_maude.colnames).issubset(imgs_mica.colnames)


has_l0_2007_archive = os.path.exists(os.path.join(aca_l0.CONFIG["data_root"], "2007"))
Expand Down

0 comments on commit afd867d

Please sign in to comment.