`get_calalign_offsets` fails on subsets, and assigns `after_caldb` in the wrong row order

[javierggt]

A report from Claudio. I have not verified it yet.

get_calalign_offsets fails on subsets, and assigns after_caldb in the wrong row order

Two independent bugs in astromon.utils.get_calalign_offsets. Both are invisible when the
function is called the way celmon calls it (the full cross-match table, which happens to
come out sorted by (obsid, x_id)), and both appear as soon as it is called with anything
else.

Environment: astromon 1.3.2, astropy 7.2.0, numpy 2.3.5, python 3.13.11 (ska3-aca).

---

Bug 1 — IndexError when every caldb_version has the same number of components

This bug is reportedly fixed in #50.

Reproduction

Self-contained, no database or CALDB access beyond the default CALALIGN directory:

from astromon import utils
from astropy.table import Table
from cxotime import CxoTime


def matches(caldb_versions):
    return Table(
        {
            "obsid": [1, 2],
            "x_id": [1, 1],
            "detector": ["ACIS-S", "ACIS-S"],
            "time": CxoTime(["2023:001", "2023:002"]),
            "caldb_version": caldb_versions,
        }
    )


for label, versions in [
    ("mixed  ['4.10.0', '4.9.0.1']", ["4.10.0", "4.9.0.1"]),
    ("uniform ['4.10.0', '4.10.0']", ["4.10.0", "4.10.0"]),
]:
    try:
        utils.get_calalign_offsets(matches(versions))
        print(f"{label}  -> ok")
    except Exception as exc:
        print(f"{label}  -> {type(exc).__name__}: {exc}")

Output:

mixed  ['4.10.0', '4.9.0.1']  -> ok
uniform ['4.10.0', '4.10.0']  -> IndexError: too many indices for array: array is 1-dimensional, but 2 were indexed

Traceback:

  File "astromon/utils.py", line 505, in get_calalign_offsets
    actual = match_w_cal[[r["cdbv"] >= r["cav"].tolist() for r in match_w_cal]]
  File "astropy/table/table.py", line 2132, in __getitem__
    return self._new_from_slice(item)
IndexError: too many indices for array: array is 1-dimensional, but 2 were indexed

Cause

utils.py:481-486 stores the parsed versions as table columns:

match_w_cal["cdbv"] = [
    [int(f) for f in s.split(".")] for s in match_w_cal["caldb_version"]
]

The comment at utils.py:503-504 states the assumption:

# NOTE: cdbv's dtype is "object" because versions have different lengthd
# but cav's dtype is array, so need to convert to list

That only holds when the version strings have differing numbers of dot-separated
components. When they are all the same length, astropy stores a 2-D integer column instead
of an object column of lists. r["cdbv"] is then an ndarray, so r["cdbv"] >= r["cav"].tolist()
broadcasts element-wise and returns an array rather than a bool. The list comprehension at
utils.py:505 becomes a list of arrays, i.e. a 2-D mask, and indexing the table with it
raises.

The dtype therefore depends on the contents of the input table, which is why this is
invisible on the full cross-match table (it always contains both 3-component versions like
4.12.0 and 4-component ones like 4.11.0.1) and fails on short time slices:

input	n	caldb_version component counts	result
full mission (mta selection)	1377	[3, 4]	works
last 5 years	278	[3, 4]	works
last 1 year	13	[3]	IndexError

The 5-year window works only by accident; once the last 4-component CalDB version falls out
of it, it will start failing too.

Suggested fix

Don't round-trip the parsed versions through table columns — compare plain Python tuples,
which also gives the right ordering semantics for versions of different lengths
((4, 10, 0) < (4, 10, 0, 1)):

def _version_tuple(version):
    return tuple(int(field) for field in version.split("."))

cdbv = [_version_tuple(v) for v in match_w_cal["caldb_version"]]
cav = [_version_tuple(v) for v in match_w_cal["calalign_version"]]

actual = match_w_cal[[c >= v for c, v in zip(cdbv, cav, strict=True)]]
reference = match_w_cal[[tuple(ref_calalign) >= v for v in cav]]

(ref_calalign is built as a list of ints at utils.py:493-501, so it needs the explicit
tuple() — comparing a list with a tuple raises TypeError.)

Note match_w_cal.sort(keys=["obsid", "x_id", "cav", "start"], reverse=True) at
utils.py:490 sorts on cav, so it would need a sortable column (or a sort key derived
from the tuples) if the columns are dropped.

---

Bug 2 — after_caldb is assigned in the input row order, not the output row order

Reproduction

from astromon import utils
from astropy.table import Table
from cxotime import CxoTime

# Deliberately NOT sorted by (obsid, x_id). Row 3 is left in place so that
# `np.all(all_matches["obsid"] != result["obsid"])` is False and the internal guard
# at utils.py:540 does not fire.
matches = Table(
    {
        "obsid": [2, 1, 3],
        "x_id": [1, 1, 1],
        "detector": ["ACIS-S", "ACIS-S", "ACIS-S"],
        "time": CxoTime(["2023:001", "2015:001", "2023:001"]),
        "caldb_version": ["4.12.0", "4.11.0.1", "4.12.0"],
    }
)

result = utils.get_calalign_offsets(matches)
time_of = dict(zip(matches["obsid"], matches["time"].date, strict=True))
result["obs_date"] = [time_of[obsid] for obsid in result["obsid"]]
result["expected"] = CxoTime(result["since"]) < CxoTime(result["obs_date"])
result["ok"] = result["after_caldb"] == result["expected"]
print(result[["obsid", "since", "obs_date", "after_caldb", "expected", "ok"]])

Output — the first two rows are both wrong:

obsid        since               obs_date       after_caldb expected   ok
----- ------------------- --------------------- ----------- -------- -----
    1 2022-06-28T14:00:00 2015:001:00:00:00.000        True    False False
    2 2022-06-28T14:00:00 2023:001:00:00:00.000       False     True False
    3 2022-06-28T14:00:00 2023:001:00:00:00.000        True     True  True

OBSID 1 was observed in 2015, well before the reference CALALIGN was added in 2022, yet it
is reported as after_caldb=True; OBSID 2, observed in 2023, is reported as False.

Cause

utils.py:534-546:

result = join(
    all_matches[["obsid", "x_id"]],
    hstack([actual[["obsid", "x_id"] + cols], reference[ref_cols + ["since"]]]),
    keys=["obsid", "x_id"],
)
...
result["after_caldb"] = result["since"] < all_matches["time"]

result comes out of join sorted by (obsid, x_id), but all_matches["time"] is still in
the caller's row order. The two are zipped positionally, so each row gets some other
observation's time whenever the input is not already key-sorted.

The returned obsid/x_id columns are correct, and so are calalign_dy, ref_calalign_dy
and friends (they travel through the join). Only after_caldb is scrambled — which means a
caller that carefully re-joins the result on (obsid, x_id) still gets a corrupted
after_caldb, because the damage is done before the function returns.

Impact on real data

Full mta cross-match table (the celmon selection: 1377 matches after dropping
caldb_version == "0.0"), calling the function once with the table sorted by
(obsid, x_id) and once sorted by time:

matches: 1377
after_caldb differs between the two input orders: 52
  key_sorted  input -> 0 rows wrong
  time_sorted input -> 52 rows wrong

so roughly 4% of the rows get the wrong after_caldb as soon as the caller sorts by time.

celmon is not affected today, because db.get_cross_matches returns rows grouped by
(obsid, x_id) and celmon does not re-sort before calling this. But that is a coincidence of
the two functions, not something either one guarantees.

Suggested fix

Carry time through the join instead of reading it from the input table:

result = join(
    all_matches[["obsid", "x_id", "time"]],
    hstack([actual[["obsid", "x_id"] + cols], reference[ref_cols + ["since"]]]),
    keys=["obsid", "x_id"],
)
...
result["after_caldb"] = CxoTime(result["since"]) < result["time"]

---

Related: the internal consistency guards never fire

The guards meant to catch exactly this kind of misalignment use np.all where they want
np.any, so they only trigger when every single row is misplaced:

• utils.py:529 if np.all(reference["obsid"] != actual["obsid"]):
• utils.py:531 if np.all(reference["x_id"] != actual["x_id"]):
• utils.py:540 if np.all(all_matches["obsid"] != result["obsid"]):
• utils.py:542 if np.all(all_matches["x_id"] != result["x_id"]):

With np.any the 3-row reproduction above would have raised instead of returning wrong
values. (The 2-row version of it does raise, since there every row is out of place — which
is a good illustration of how narrow the current check is.)