diff --git a/scripts/calibration/params.py b/scripts/calibration/params.py index ea871119..7f6d7c41 100644 --- a/scripts/calibration/params.py +++ b/scripts/calibration/params.py @@ -143,7 +143,6 @@ "NUMBER", "IMAFLAGS_ISO", "FLAGS", - "NGMIX_MCAL_FLAGS", "NGMIX_MCAL_TYPES_FAIL", "N_EPOCH", "NGMIX_N_EPOCH", @@ -152,6 +151,8 @@ add_cols_pre_cal_format["TILE_ID"] = "A7" add_cols_pre_cal_format["NUMBER"] = "J" +# Metacal bitmasks hold ngmix flag bits 0-15; bit 15 overflows signed 16-bit I. +add_cols_pre_cal_format["NGMIX_MCAL_FLAGS"] = "J" # Create key names for metacal information prefix = "NGMIX" @@ -162,7 +163,7 @@ add_cols_pre_cal.append(f"{prefix}_{center}_{suffix}") for suffix in suffixes: - add_cols_pre_cal_format[f"FLAGS_{suffix}"] = "I" + add_cols_pre_cal_format[f"{prefix}_FLAGS_{suffix}"] = "J" # Catalog parameters diff --git a/src/sp_validation/catalog.py b/src/sp_validation/catalog.py index e39d2076..76022243 100644 --- a/src/sp_validation/catalog.py +++ b/src/sp_validation/catalog.py @@ -439,6 +439,13 @@ def write_shape_catalog( Write catalogue with galaxy shapes = shear estimates. + @sc [label:convention] metacal-flag-width + The data and image-simulation parameter files supply FITS ``J`` (int32) + for metacal bitmasks, so all ngmix flag bits (0-15; ``ZERO_DOF`` = 2**15 + overflows a signed 16-bit ``I``) survive FITS/HDF5 output; ``JointCat`` + never narrows integers. + The number of failed metacal types is a count in [0, 5], not a bitmask. + Parameters ---------- output_path : str diff --git a/src/sp_validation/catalog_builders.py b/src/sp_validation/catalog_builders.py index 78dc7838..b6abafc0 100644 --- a/src/sp_validation/catalog_builders.py +++ b/src/sp_validation/catalog_builders.py @@ -437,26 +437,19 @@ def dtype_out(self, name, dtype_in): output dtype """ - # Specify columns for which original (high-precision) format - # needs to be kept and not reduced to lower precision - cols_keep_dtype = [ - "RA", - "Dec", - "FLAGS", - "IMAFLAGS_ISO", - "NUMBER", - ] if dtype_in.kind == "U": # Transform unicode to string of equal length return np.dtype(f"S{dtype_in.itemsize // 4}") - if self._params["reduce_mem"] == False: - return dtype_in - elif name not in cols_keep_dtype: - if dtype_in.kind == "f" and dtype_in.itemsize == 8: - return np.float32 - if dtype_in.kind == "i" and dtype_in.itemsize == 4: - return np.int8 + # reduce_mem narrows float64 to float32, except coordinates. Integer + # columns (bitmasks, IDs, counts) are never narrowed. + if ( + self._params["reduce_mem"] + and dtype_in.kind == "f" + and dtype_in.itemsize == 8 + and name not in ("RA", "Dec") + ): + return np.dtype(np.float32) return dtype_in diff --git a/src/sp_validation/tests/test_catalog_flag_roundtrip.py b/src/sp_validation/tests/test_catalog_flag_roundtrip.py new file mode 100644 index 00000000..ea3cac97 --- /dev/null +++ b/src/sp_validation/tests/test_catalog_flag_roundtrip.py @@ -0,0 +1,80 @@ +"""Comprehensive catalogues preserve metacal failure bits exactly.""" + +import runpy +from pathlib import Path + +import h5py +import numpy as np +import pytest +from astropy.io import fits + +from sp_validation.catalog import write_shape_catalog +from sp_validation.catalog_builders import JointCat + +ROOT = Path(__file__).resolve().parents[3] + + +@pytest.mark.parametrize( + "params_path", + ["scripts/calibration/params.py", "workflow/image_sims/params_im_sim.py"], + ids=["data", "image-sims"], +) +@pytest.mark.parametrize("extension", [".fits", ".hdf5"]) +@pytest.mark.parametrize("reduce_mem", [False, True]) +def test_metacal_flags_roundtrip(tmp_path, params_path, extension, reduce_mem): + """Contract metacal-flag-width: neither output nor merging loses flag bits. + + Values are ngmix flag bits: none, LM_FUNC_NOTFINITE (2**12), ZERO_DOF + (2**15), and every bit 0-15 at once. Float64 inputs match ShapePipe's final + catalogue. Bit 15 overflows signed int16, and int8 wraps everything above + 127, so the columns must stay int32. + """ + with np.printoptions(): + params = runpy.run_path(str(ROOT / params_path)) + flags = np.array([0, 2**12, 2**15, 2**15 | 2**12 | 8, 2**16 - 1], dtype=np.float64) + bit_columns = ["NGMIX_MCAL_FLAGS"] + [ + f"NGMIX_FLAGS_{suffix}" for suffix in ("NOSHEAR", "1P", "1M", "2P", "2M") + ] + columns = {name: flags for name in bit_columns} + columns["NGMIX_MCAL_TYPES_FAIL"] = np.array([0, 1, 5, 2, 5]) + assert set(columns) <= set(params["add_cols_pre_cal"]) + path = tmp_path / f"comprehensive{extension}" + write_shape_catalog( + str(path), + np.zeros(5), + np.zeros(5), + np.ones(5), + add_cols=columns, + add_cols_format=params["add_cols_pre_cal_format"], + ) + if extension == ".fits": + written = fits.getdata(path, 1) + else: + with h5py.File(path, "r") as catalog: + written = catalog["data"][:] + + builder = JointCat() + builder._params["reduce_mem"] = reduce_mem + for name, expected in columns.items(): + values = written[name] + np.testing.assert_array_equal(values, expected, err_msg=name) + if name in bit_columns: + assert values.dtype.kind == "i" and values.dtype.itemsize == 4, name + reduced = values.astype(builder.dtype_out(name, values.dtype)) + np.testing.assert_array_equal(reduced, expected, err_msg=name) + + +@pytest.mark.parametrize("reduce_mem", [False, True]) +def test_reduce_mem_never_narrows_integers(reduce_mem): + """reduce_mem narrows float64 (except RA/Dec) and leaves integers intact.""" + builder = JointCat() + builder._params["reduce_mem"] = reduce_mem + for dtype in (">i4", "i8", "f8")) == np.dtype(">f8") + expected = np.float32 if reduce_mem else np.dtype(">f8") + assert builder.dtype_out("NGMIX_G1_NOSHEAR", np.dtype(">f8")) == expected diff --git a/workflow/image_sims/params_im_sim.py b/workflow/image_sims/params_im_sim.py index 1cd395c6..b7db399e 100644 --- a/workflow/image_sims/params_im_sim.py +++ b/workflow/image_sims/params_im_sim.py @@ -148,7 +148,6 @@ for key in ( "NUMBER", "FLAGS", - "NGMIX_MCAL_FLAGS", "NGMIX_MCAL_TYPES_FAIL", "N_EPOCH", "NGMIX_N_EPOCH", @@ -157,6 +156,8 @@ add_cols_pre_cal_format["TILE_ID"] = "A7" add_cols_pre_cal_format["NUMBER"] = "J" +# Metacal bitmasks hold ngmix flag bits 0-15; bit 15 overflows signed 16-bit I. +add_cols_pre_cal_format["NGMIX_MCAL_FLAGS"] = "J" # Create key names for metacal information prefix = "NGMIX" @@ -167,7 +168,7 @@ add_cols_pre_cal.append(f"{prefix}_{center}_{suffix}") for suffix in suffixes: - add_cols_pre_cal_format[f"FLAGS_{suffix}"] = "I" + add_cols_pre_cal_format[f"{prefix}_FLAGS_{suffix}"] = "J" # Catalog parameters