diff --git a/astra.yaml b/astra.yaml index d93c2067c..bcd074da0 100644 --- a/astra.yaml +++ b/astra.yaml @@ -469,22 +469,19 @@ analyses: release: on the DR5 image of a tile only ~95% of detections pair. Unpaired detections still stop the run beyond both 2% and 20 objects, preserving the check for other pixels. External-only - losses stop beyond both 20 objects and 2% of the external count - plus an additive 50-object allowance (only when pairs exist). - Four production edge tiles have 49, 59, 55 and 101 external-only - objects among 1680-2861 objects, concentrated in bright-star - halos/spikes at normal positive weight; only 0-3 detections lack - partners and paired positions agree to 7e-5 px. Local correlated - deblending losses do not scale with total tile population. The - allowance is empirical, rounded above the largest observed excess - over 2% (44 objects), not a binomial confidence bound. It weakens - sensitivity to modest external-side detection-config losses on - sparse tiles; there are too few uncensored sparse tiles to - calibrate a false-rejection rate. A zero-weight exemption would - not explain these losses, and spatial exemptions would add mask - and clustering choices. Only the acceptance gate changes: no - external-only object gains a shape and every previously accepted - catalogue retains identical membership and measurements. + losses never reject a tile: warn and record, then cut downstream. + The warning threshold is max(20, 0.02 N + 50), without the 50 if + no pairs exist; it no longer shares the SExtractor guard settings. + Sparse edge tiles 058.220 and 207.238 exceed the former budget, + despite almost every detection pairing. Local correlated + deblending losses need not scale with tile population. + Every joined tile records the five match counts, external paired + fraction and warning budget in its detection catalogue headers, + then in tile_detect and the published tile_make_cat.json manifest + beside final_cat. Join its unit (RRR.DDD) as a float to TILE_ID in + the final catalogue to select the population downstream. + No external-only object gains a shape; previously accepted + catalogue rows retain identical membership and measurements. PSF stars come from the exposure-level SExtractor run in both options. Image simulations, which have no UNIONS catalogue, keep SExtractor's own NUMBER. @@ -492,8 +489,7 @@ analyses: input_types.data.tile_detection = unions_catalogue; config_tile_Sx.ini#SEXTRACTOR_RUNNER.MATCH_RADIUS = 1.0; config_tile_Sx.ini#SEXTRACTOR_RUNNER.MATCH_MIN_FRACTION = 0.98; - config_tile_Sx.ini#SEXTRACTOR_RUNNER.MATCH_TOLERATED_UNPAIRED = 20; - config_tile_Sx.ini#SEXTRACTOR_RUNNER.MATCH_EXTERNAL_UNPAIRED_ALLOWANCE = 50. + config_tile_Sx.ini#SEXTRACTOR_RUNNER.MATCH_TOLERATED_UNPAIRED = 20. default: unions_catalogue options: unions_catalogue: diff --git a/src/shapepipe/modules/make_cat_package/make_cat.py b/src/shapepipe/modules/make_cat_package/make_cat.py index a04d4c52c..1220ffa28 100644 --- a/src/shapepipe/modules/make_cat_package/make_cat.py +++ b/src/shapepipe/modules/make_cat_package/make_cat.py @@ -364,7 +364,7 @@ def _save_ngmix_data(self, ngmix_cat_path, moments=False): # A low match fraction can be valid (e.g. missing PSFs). Warn but # continue; unmatched detections retain their pre-filled sentinels. # Completeness gates reject premature or incomplete merges. - if len(ngmix_n_epoch) / self._cat_size_target < 0.1: + if self._cat_size_target and len(ngmix_n_epoch) / self._cat_size_target < 0.1: warning = ( f"Merged shape catalogue {ngmix_cat_path} contains only " f"{len(ngmix_n_epoch)} of {self._cat_size_target} target " @@ -638,7 +638,7 @@ def _save_psf_data( if not epoch_slots: n_slots = 0 elif n_epoch_slots is None: - n_slots = np.max(self._final_cat["N_EPOCH"]) + 1 + n_slots = int(np.max(self._final_cat["N_EPOCH"], initial=0)) + 1 else: n_slots = n_epoch_slots diff --git a/src/shapepipe/modules/sextractor_package/match_catalogue.py b/src/shapepipe/modules/sextractor_package/match_catalogue.py index 57ab7b59c..1424b9bdd 100644 --- a/src/shapepipe/modules/sextractor_package/match_catalogue.py +++ b/src/shapepipe/modules/sextractor_package/match_catalogue.py @@ -16,6 +16,8 @@ """ +import logging + import numpy as np from astropy.io import ascii as asc from astropy.io import fits @@ -30,6 +32,21 @@ # The largest NUMBER the int32 NUMBER and SEG_VIGNET columns hold. MAX_NUMBER = np.iinfo(np.int32).max +# Carried through tile_detect into the published tile_make_cat.json manifest. +MATCH_HEADER = dict( + n_sextractor="MCHNSEX", n_external="MCHNEXT", n_paired="MCHNPAIR", + n_dropped="MCHNDROP", n_external_only="MCHNONLY", + external_paired_fraction="MCHFRAC", external_warning_budget="MCHBUDG", +) + + +def read_match_record(cat_path): + """Read join diagnostics, or None for a catalogue with no external join.""" + header = fits.getheader(cat_path, "LDAC_OBJECTS") + if MATCH_HEADER["n_sextractor"] not in header: + return None + return {name: header[key] for name, key in MATCH_HEADER.items()} + def mutual_nearest(x_a, y_a, x_b, y_b, radius): """One-to-one pairs of mutual nearest neighbours closer than ``radius``. @@ -97,8 +114,7 @@ def relabel(seg_vignets, old_number, new_number): def match_catalogue(cat_path, ext_cat_path, radius=1.0, min_fraction=0.98, - tolerated_unpaired=20, w_log=None, - external_unpaired_allowance=50): + tolerated_unpaired=20, w_log=None): """Take membership and ``NUMBER`` from an external catalogue. Each SExtractor row is paired with its mutual nearest neighbour in the @@ -119,15 +135,11 @@ def match_catalogue(cat_path, ext_cat_path, radius=1.0, min_fraction=0.98, the DR6 catalogue) only ~95% pair. The guard is there to catch that: the SExtractor-side guard stops when unpaired detections exceed both ``(1 - min_fraction)`` of that side and ``tolerated_unpaired``. - External-only losses instead get an additive local deblending allowance - on top of the fractional budget (and retain the absolute tolerance). - Sparse DR6 edge tiles can lose dozens of children in one bright-star - halo while almost every detection pairs; those losses do not scale with - the tile's total population. This is an empirical allowance, not an - independent-object confidence bound. It weakens sensitivity to modest - external-side detection-config losses on sparse tiles, but does not - relax the SExtractor-side check for other pixels. With no pairs there - is no evidence for a local loss, so no additive allowance is granted. + External-only losses never reject a tile: warn above max(20, 0.02 N + + 50), the former rejection budget (without the 50 when no pairs exist). + Record every tile's counts, external paired fraction and warning budget + in LDAC_OBJECTS headers for publication in the make-cat manifest, so + analysts can cut downstream. Sparse edge losses need not scale with N. Parameters ---------- @@ -139,31 +151,27 @@ def match_catalogue(cat_path, ext_cat_path, radius=1.0, min_fraction=0.98, radius : float, optional Largest separation of a pair, in pixels; default 1 min_fraction : float, optional - Baseline paired fraction on each side; default 0.98. External-only - losses additionally receive ``external_unpaired_allowance``. + Minimum SExtractor paired fraction; default 0.98. tolerated_unpaired : int, optional - Number of unpaired rows on either side accepted whatever the + Number of unpaired SExtractor rows accepted whatever the fraction; default 20 w_log : logging.Logger, optional Pipeline logger - external_unpaired_allowance : int, optional - Additive external-only loss budget, in objects; default 50. Only - applied when at least one pair exists; zero restores the old guard. - Returns ------- dict ``n_sextractor``, ``n_external``, ``n_paired``, ``n_dropped`` (SExtractor rows without a partner) and ``n_external_only`` - (external objects without a SExtractor row) + (external objects without a SExtractor row), plus + ``external_paired_fraction`` (zero for an empty external catalogue) + and ``external_warning_budget`` Raises ------ ValueError If the external ``NUMBER`` repeats or is not an integer in ``[1, MAX_NUMBER]``, or if unpaired rows exceed the absolute tolerance - and the side's loss budget: fractional only for SExtractor, - fractional plus the local allowance for the external catalogue + and the SExtractor fractional loss budget @sc [decision:detection.tile_detection] """ @@ -194,10 +202,10 @@ def match_catalogue(cat_path, ext_cat_path, radius=1.0, min_fraction=0.98, n_external_only=len(ext) - len(i_sex), ) summary = ", ".join(f"{k}={v}" for k, v in counts.items()) - def too_many_unpaired(n_total, allowance=0): + def too_many_unpaired(n_total): unpaired = n_total - len(i_sex) return (unpaired > tolerated_unpaired - and unpaired > (1 - min_fraction) * n_total + allowance) + and unpaired > (1 - min_fraction) * n_total) fraction = len(i_sex) / max(len(data), 1) if too_many_unpaired(len(data)): @@ -210,16 +218,15 @@ def too_many_unpaired(n_total, allowance=0): + " release of the catalogue (DR6)." ) fraction_ext = len(i_sex) / max(len(ext), 1) - allowance = external_unpaired_allowance if len(i_sex) else 0 - if too_many_unpaired(len(ext), allowance): - raise ValueError( + budget = max(20, 0.02 * len(ext) + (50 if len(i_sex) else 0)) + counts.update(external_paired_fraction=fraction_ext, + external_warning_budget=budget) + if counts["n_external_only"] > budget: + (w_log or logging.getLogger(__name__)).warning( f"Only {fraction_ext:.4f} of the objects of {ext_cat_path} pair" - + f" with a row of {cat_path} within {radius} px (baseline" - + f" {min_fraction} beyond {tolerated_unpaired} unpaired;" - + f" additive external allowance {allowance};" - + f" {summary}). SExtractor finds fewer objects" - + " than the catalogue: check that the detection configuration" - + " is the catalogue's." + + f" with a row of {cat_path} within {radius} px (external-only" + + f" warning budget {budget}; {summary}). Recording for a" + + " downstream catalogue cut; the tile is not rejected." ) old_number = np.asarray(data["NUMBER"]) @@ -241,6 +248,8 @@ def too_many_unpaired(n_total, allowance=0): new = fits.BinTableHDU.from_columns( columns, header=objects.header, name="LDAC_OBJECTS", ) + for name, key in MATCH_HEADER.items(): + new.header[key] = counts[name] fits.HDUList( [new if h is objects else h for h in hdus] ).writeto(cat_path, overwrite=True) diff --git a/src/shapepipe/modules/sextractor_package/sextractor_script.py b/src/shapepipe/modules/sextractor_package/sextractor_script.py index e41ddad57..f44ce68e4 100644 --- a/src/shapepipe/modules/sextractor_package/sextractor_script.py +++ b/src/shapepipe/modules/sextractor_package/sextractor_script.py @@ -76,7 +76,7 @@ def seg_vignet_column(seg_vignets): return fits.Column( name="SEG_VIGNET", format=f"{ny * nx}J", - array=seg_vignets.astype(np.int32, copy=False).reshape(n_obj, -1), + array=seg_vignets.astype(np.int32, copy=False).reshape(n_obj, ny * nx), dim=f"({nx},{ny})", ) @@ -193,7 +193,8 @@ def add_seg_vignet(cat_path, seg_path, w_log=None): ).writeto(cat_path, overwrite=True) if w_log: centre = size // 2 - own = np.mean(seg_vignets[:, centre, centre] == data["NUMBER"]) + own = (np.mean(seg_vignets[:, centre, centre] == data["NUMBER"]) + if len(col) else np.nan) w_log.info( f"SEG_VIGNET cut from {seg_path} for {len(col)} objects;" + f" centre label = NUMBER for {own:.4f}" @@ -430,7 +431,7 @@ def make_post_process(cat_path, f_wcs_path, pos_params, ccd_size, w_log=None): # (laid out as FITSCatalogue.add_col lays it out), then the EPOCH HDUs. new = fits.BinTableHDU.from_columns( objects.data.columns + fits.ColDefs([column("N_EPOCH", n_epoch)]), - name="LDAC_OBJECTS", + header=objects.header, name="LDAC_OBJECTS", ) fits.HDUList( [new if h is objects else h for h in hdus] + epochs diff --git a/src/shapepipe/modules/sextractor_runner.py b/src/shapepipe/modules/sextractor_runner.py index 1b4754a8b..2a0df4713 100644 --- a/src/shapepipe/modules/sextractor_runner.py +++ b/src/shapepipe/modules/sextractor_runner.py @@ -177,10 +177,6 @@ def sextractor_runner( tolerated_unpaired=config.getint( module_config_sec, "MATCH_TOLERATED_UNPAIRED" ), - external_unpaired_allowance=config.getint( - module_config_sec, "MATCH_EXTERNAL_UNPAIRED_ALLOWANCE", - fallback=50, - ), w_log=w_log, ) diff --git a/src/shapepipe/modules/vignetmaker_package/vignetmaker.py b/src/shapepipe/modules/vignetmaker_package/vignetmaker.py index f0d7d9dd9..cf4173971 100644 --- a/src/shapepipe/modules/vignetmaker_package/vignetmaker.py +++ b/src/shapepipe/modules/vignetmaker_package/vignetmaker.py @@ -335,6 +335,8 @@ def _get_stamp_me(self, image_dirs, image_pattern): + f" HDU #{hdu_index}" ) if len(exp_names) == 0: + if len(all_id) == 0: + continue # A valid zero-object tile has empty EPOCH HDUs. raise IndexError( f"EXP_NAME list is empty in file {self._galcat_path}" + f" HDU #{hdu_index}" diff --git a/tests/module/test_match_catalogue.py b/tests/module/test_match_catalogue.py index d56767d3a..b17a21e67 100644 --- a/tests/module/test_match_catalogue.py +++ b/tests/module/test_match_catalogue.py @@ -73,7 +73,10 @@ def test_rows_take_the_external_number_one_to_one(tmp_path): expected[order] = numbers assert list(data["NUMBER"]) == list(expected) assert counts == dict(n_sextractor=5, n_external=6, n_paired=5, - n_dropped=0, n_external_only=1) + n_dropped=0, n_external_only=1, + external_paired_fraction=5 / 6, + external_warning_budget=50.12) + assert mc.read_match_record(cat) == counts # Every measurement is SExtractor's, row for row. np.testing.assert_array_equal(data["XWIN_IMAGE"], X + 0.01) np.testing.assert_array_equal(data["VIGNET"][2].ravel(), @@ -116,17 +119,14 @@ def test_pairs_lie_within_the_radius(tmp_path): assert counts["n_dropped"] == 1 and counts["n_external_only"] == 1 -def test_too_few_catalogue_objects_paired_stop_the_run(tmp_path): - """Fewer detections than catalogue objects (a drifted detection - configuration) fails, though every detection pairs.""" +def test_external_record_is_independent_of_sextractor_guard(tmp_path): cat = _sexcat(tmp_path / "sexcat.fits", X[:4], Y[:4]) ext = _external(tmp_path / "ext.cat", [1, 2, 3, 4, 5], X, Y) - with pytest.raises(ValueError, match="fewer objects"): - mc.match_catalogue(cat, ext, min_fraction=0.99, tolerated_unpaired=0, - external_unpaired_allowance=0) - assert list(_objects(cat)["NUMBER"]) == [1, 2, 3, 4] - counts = mc.match_catalogue(cat, ext, min_fraction=0.8) + counts = mc.match_catalogue(cat, ext, min_fraction=0.99, + tolerated_unpaired=0) assert counts["n_paired"] == 4 and counts["n_external_only"] == 1 + assert counts["external_warning_budget"] == 50.1 + assert mc.read_match_record(cat) == counts def test_too_few_pairs_stop_the_run(tmp_path): @@ -162,6 +162,8 @@ def test_small_tile_with_few_unpaired_children_passes(tmp_path): (2490, 2431, 2431), # 577.246 (1680, 1625, 1625), # 603.227 (2861, 2760, 2760), # 618.220 + (845, 774, 774), # 058.220, real b619 counts; formerly rejected + (2400, 2296, 2296), # 207.238-like: 0.9567; synthetic counts ]) @pytest.mark.decision("detection.tile_detection") def test_production_sparse_tiles_keep_only_the_pairs( @@ -180,7 +182,10 @@ def test_production_sparse_tiles_keep_only_the_pairs( counts = mc.match_catalogue(cat, ext) assert counts == dict(n_external=n_ext, n_sextractor=n_sex, n_paired=n_paired, n_dropped=n_sex - n_paired, - n_external_only=n_ext - n_paired) + n_external_only=n_ext - n_paired, + external_paired_fraction=n_paired / n_ext, + external_warning_budget=0.02 * n_ext + 50) + assert mc.read_match_record(cat) == counts data = _objects(cat) np.testing.assert_array_equal(data["NUMBER"], np.arange(10001, 10001 + n_paired)) @@ -211,49 +216,49 @@ def test_large_tile_with_unpaired_children_passes(tmp_path): assert counts["n_external_only"] == 30 -@pytest.mark.parametrize("side", ["sextractor", "external"]) -def test_other_pixels_stop_the_run(tmp_path, side): - """5% unpaired on either side (a DR5 image: ~95% pair) stops the run with - the default guard, far beyond the 20 tolerated.""" +def test_other_pixels_stop_the_run(tmp_path): + """5% of our detections unpaired still stops the run.""" x, y = _grid(2000) - n_sex, n_ext = (2000, 1900) if side == "sextractor" else (1900, 2000) + n_sex, n_ext = 2000, 1900 cat = _sexcat(tmp_path / "sexcat.fits", x[:n_sex], y[:n_sex]) ext = _external(tmp_path / "ext.cat", np.arange(1, n_ext + 1), x[:n_ext], y[:n_ext]) - with pytest.raises(ValueError, match="DR6" if side == "sextractor" - else "fewer objects"): + with pytest.raises(ValueError, match="DR6"): mc.match_catalogue(cat, ext) -@pytest.mark.parametrize("n_ext,missing,passes", [ - (2000, 89, True), (2000, 91, False), - (40000, 849, True), (40000, 851, False), +@pytest.mark.parametrize("n_ext,missing,warns", [ + (2000, 89, False), (2000, 91, True), + (40000, 849, False), (40000, 851, True), + (845, 71, True), (2400, 104, True), ]) @pytest.mark.decision("detection.tile_detection") -def test_external_loss_budget_boundary(tmp_path, n_ext, missing, passes): - """The local allowance is bounded; a detection-config regression fails. - - Missing rows occupy a single grid region: clustering grants no exemption. - """ +def test_external_loss_budget_boundary(tmp_path, caplog, n_ext, missing, warns): + """The former rejection budget only controls a warning, never rows.""" x, y = _grid(n_ext) cat = _sexcat(tmp_path / "sexcat.fits", x[:-missing], y[:-missing]) ext = _external(tmp_path / "ext.cat", np.arange(1, n_ext + 1), x, y) - if passes: - assert mc.match_catalogue(cat, ext)["n_external_only"] == missing - else: - before = (tmp_path / "sexcat.fits").read_bytes() - with pytest.raises(ValueError, match="fewer objects"): - mc.match_catalogue(cat, ext) - assert (tmp_path / "sexcat.fits").read_bytes() == before + counts = mc.match_catalogue(cat, ext) + assert counts["n_external_only"] == missing + assert mc.read_match_record(cat) == counts + assert ("downstream catalogue cut" in caplog.text) == warns + if warns: + assert f"n_external_only={missing}" in caplog.text + assert len(_objects(cat)) == n_ext - missing -def test_empty_detection_does_not_get_a_deblending_allowance(tmp_path): - """No paired population is not evidence for a local deblending loss.""" +def test_empty_detection_warns_and_records(tmp_path, caplog): + """Even no pairs is a downstream decision, not an external-side veto.""" x, y = _grid(30) cat = _sexcat(tmp_path / "sexcat.fits", [], []) ext = _external(tmp_path / "ext.cat", np.arange(1, 31), x, y) - with pytest.raises(ValueError, match="fewer objects"): - mc.match_catalogue(cat, ext) + counts = mc.match_catalogue(cat, ext) + assert counts["n_external_only"] == 30 + assert counts["external_paired_fraction"] == 0 + assert counts["external_warning_budget"] == 20 + assert mc.read_match_record(cat) == counts + assert "downstream catalogue cut" in caplog.text + assert len(_objects(cat)) == 0 @pytest.mark.decision("detection.tile_detection") @@ -265,19 +270,11 @@ def test_previously_accepted_catalogue_is_byte_identical(tmp_path): ext = _external(tmp_path / "ext.cat", np.arange(10001, 12001), x, y) old = _sexcat(tmp_path / "old.fits", x[:1997], y[:1997], seg) new = _sexcat(tmp_path / "new.fits", x[:1997], y[:1997], seg) - old_counts = mc.match_catalogue(old, ext, external_unpaired_allowance=0) + old_counts = mc.match_catalogue(old, ext) assert mc.match_catalogue(new, ext) == old_counts assert (tmp_path / "old.fits").read_bytes() == (tmp_path / "new.fits").read_bytes() -def test_external_allowance_can_be_disabled(tmp_path): - x, y = _grid(2000) - cat = _sexcat(tmp_path / "sexcat.fits", x[:1911], y[:1911]) - ext = _external(tmp_path / "ext.cat", np.arange(1, 2001), x, y) - with pytest.raises(ValueError, match="fewer objects"): - mc.match_catalogue(cat, ext, external_unpaired_allowance=0) - - def test_seg_vignet_is_relabelled(tmp_path): """SEG_VIGNET labels follow NUMBER; a dropped row's footprint is -1.""" x, y = X[:3], Y[:3] diff --git a/tests/module/test_match_record_publication.py b/tests/module/test_match_record_publication.py new file mode 100644 index 000000000..03573ace4 --- /dev/null +++ b/tests/module/test_match_record_publication.py @@ -0,0 +1,60 @@ +"""Join diagnostics reach the manifest published beside each tile catalogue.""" + +import importlib +import json +from pathlib import Path +import subprocess + +import pytest + +from shapepipe.modules.sextractor_package import match_catalogue as mc +from tests.module.test_match_catalogue import _external, _grid, _sexcat + + +@pytest.mark.parametrize("n_ext,n_pair", [(845, 774), (2400, 2296), (100, 100), (30, 0)]) +def test_match_record_reaches_published_manifest(tmp_path, monkeypatch, n_ext, n_pair): + scripts = Path(__file__).resolve().parents[2] / "workflow" / "scripts" + monkeypatch.syspath_prepend(str(scripts)) + completeness = importlib.import_module("completeness") + unit = "058.220" + run = tmp_path / unit + detect_out = run / "output/run_sp_tile_Sx/sextractor_runner/output" + detect_out.mkdir(parents=True) + x, y = _grid(n_ext) + cat = _sexcat(detect_out / "sexcat-058-220.fits", x[:n_pair], y[:n_pair]) + ext = _external(tmp_path / "external.cat", range(1, n_ext + 1), x, y) + expected = mc.match_catalogue(cat, ext) + (detect_out / "background.fits").touch() + detected, ok = completeness.build_manifest("tile_detect", run, unit) + assert ok + assert detected["catalogue_match"] == expected + mdir = run / "manifests" + mdir.mkdir() + (mdir / "tile_detect.json").write_text(json.dumps(detected)) + # Chunk aggregation is orthogonal; exercise the real make-cat manifest + # builder with its already-tested epoch-cuts reader replaced here. + monkeypatch.setattr(completeness, "gather_manifests", lambda *args: {}) + monkeypatch.setenv("NGMIX_N_CHUNKS", "1") + output = run / "output/run_sp_tile_Mc/make_cat_runner/output" + output.mkdir(parents=True) + (output / "final_cat-058-220.hdf5").touch() + made, ok = completeness.build_manifest("tile_make_cat", run, unit) + assert ok + manifest = mdir / "tile_make_cat.json" + manifest.write_text(json.dumps(made)) + published = tmp_path / "products/tiles/05/058.220/tile_make_cat.json" + published.parent.mkdir(parents=True) + # Execute the publication commands from the actual rule, not a test copy. + rules = (scripts.parent / "rules/tile.smk").read_text() + start = rules.index("' cp -f {output.manifest}") + stop = rules.index("'", rules.index("' mv -f", start) + 1) + import ast + commands = "".join(ast.literal_eval(line.strip()) + for line in rules[start:stop + 1].splitlines()) + commands = commands.replace("{output.manifest}", str(manifest)) + commands = commands.replace("{output.stats}", str(published)) + subprocess.run(["bash", "-ec", commands], check=True) + record = json.loads(published.read_text()) + assert record["unit"] == unit + assert record["catalogue_match"] == expected + assert float(record["unit"]) == 58.220 # final catalogue TILE_ID diff --git a/tests/module/test_sextractor_post_process.py b/tests/module/test_sextractor_post_process.py index 868559e36..2952d66f0 100644 --- a/tests/module/test_sextractor_post_process.py +++ b/tests/module/test_sextractor_post_process.py @@ -118,6 +118,15 @@ def test_post_process_scans_all_ccds_despite_tile_id_key(tmp_path): ) cat_path = tmp_path / "sexcat.fits" _write_sex_ldac(cat_path, exp_names, positions) + from shapepipe.modules.sextractor_package.match_catalogue import ( + MATCH_HEADER, read_match_record, + ) + record = dict(n_sextractor=2, n_external=3, n_paired=2, n_dropped=0, + n_external_only=1, external_paired_fraction=2 / 3, + external_warning_budget=50.06) + with fits.open(cat_path, mode="update") as hdus: + for name, key in MATCH_HEADER.items(): + hdus["LDAC_OBJECTS"].header[key] = record[name] sextractor_script.make_post_process( str(cat_path), @@ -126,6 +135,7 @@ def test_post_process_scans_all_ccds_despite_tile_id_key(tmp_path): ["0", str(CCD_NPIX), "0", str(CCD_NPIX)], ) + assert read_match_record(cat_path) == record with fits.open(cat_path) as hdu_list: n_epoch = hdu_list["LDAC_OBJECTS"].data["N_EPOCH"] ccd_n = { diff --git a/tests/module/test_zero_object_tile.py b/tests/module/test_zero_object_tile.py new file mode 100644 index 000000000..e65d6e078 --- /dev/null +++ b/tests/module/test_zero_object_tile.py @@ -0,0 +1,197 @@ +"""A zero-detection tile traverses the real tile modules without fake rows.""" + +import importlib +import json +from pathlib import Path +import shutil +import sqlite3 +import subprocess +import sys + +from astropy.io import fits +import h5py +import numpy as np +import pytest +from sqlitedict import SqliteDict + +from shapepipe.modules.sextractor_package import sextractor_script as sx +from shapepipe.modules.sextractor_package import match_catalogue as mc +from shapepipe.modules.merge_headers_package.merge_headers import merge_headers +from shapepipe.modules.vignetmaker_package.vignetmaker import VignetMaker +from shapepipe.modules.psfex_interp_package.psfex_interp import PSFExInterpolator +from shapepipe.modules.ngmix_runner import ngmix_runner +from shapepipe.modules.merge_sep_cats_package.merge_sep_cats import MergeSep +from shapepipe.modules.make_cat_runner import make_cat_runner +from shapepipe.pipeline.config import CustomParser +from tests.module.test_sextractor_post_process import ( + _write_exposure_headers, _write_sex_ldac, +) +from tests.module.test_match_catalogue import _external, _grid +from tests.module.test_ngmix_empty_tile import _RecordingLogger, _runner_config + + +@pytest.mark.parametrize("epoch_slots", [None, 12]) +def test_zero_object_tile_modules_and_publication( + tmp_path, monkeypatch, epoch_slots, *, tile="177.257", n_chunks=2, + products_dir=None, exposure="2360341"): + scripts = Path(__file__).resolve().parents[2] / "workflow/scripts" + monkeypatch.syspath_prepend(str(scripts)) + complete = importlib.import_module("completeness") + number = "-" + tile.replace(".", "-") + products_dir = products_dir or tmp_path / "products" + run = tmp_path / tile + log = _RecordingLogger() + + def output(stage, runner): + path = run / "output" / stage / runner / "output" + path.mkdir(parents=True, exist_ok=True) + return path + + mdir = run / "manifests" + mdir.mkdir(parents=True, exist_ok=True) + + def manifest(stage, subdir=None, suffix=""): + record, ok = complete.build_manifest(stage, run, tile, subdir) + assert ok, record["failures"] + (mdir / f"{stage}{suffix}.json").write_text(json.dumps(record)) + return record + + # Fir 177.257's one exposure, but synthetic zero-row SExtractor output. + exp = exposure + npy = tmp_path / f"headers-{exp}.npy" + _write_exposure_headers(npy) + merge_headers([[str(npy)]], str(tmp_path), tile_number=number) + headers = str(tmp_path / f"log_exp_headers{number}.sqlite") + detected = output("run_sp_tile_Sx", "sextractor_runner") + cat = detected / f"sexcat{number}.fits" + _write_sex_ldac(cat, [exp], np.empty((0, 2))) + with fits.open(cat, mode="update") as hdus: + cols = hdus["LDAC_OBJECTS"].columns + fits.ColDefs([ + fits.Column(name=name, format="D", array=np.empty(0)) + for name in ["X_IMAGE", "Y_IMAGE", "X_IMAGE_DBL", "Y_IMAGE_DBL"] + ] + [fits.Column(name="VIGNET", format="9E", dim="(3,3)", + array=np.empty((0, 3, 3), np.float32))]) + hdus[2] = fits.BinTableHDU.from_columns(cols, name="LDAC_OBJECTS") + seg = tmp_path / "seg.fits" + fits.writeto(seg, np.zeros((10, 10), np.int32)) + sx.add_seg_vignet(str(cat), str(seg), w_log=log) + with fits.open(cat) as hdus: + assert hdus[2].data["SEG_VIGNET"].shape == (0, 3, 3) + assert hdus[2].data["SEG_VIGNET"].dtype.kind == "i" + x, y = _grid(30) + ext = _external(tmp_path / "external.cat", range(1, 31), x, y) + counts = mc.match_catalogue(str(cat), ext, w_log=log) + sx.make_post_process(str(cat), headers, ["XWIN_WORLD", "YWIN_WORLD"], + ["0", "100", "0", "100"]) + assert mc.read_match_record(cat) == counts + (detected / "background.fits").touch() + assert manifest("tile_detect")["catalogue_match"] == counts + with fits.open(cat) as hdus: + assert len(hdus["EPOCH_0"].data) == 0 + assert len(hdus["LDAC_OBJECTS"].data) == 0 + + # Real PSF interpolation and all four multi-epoch vignette stores. + psfout = output("run_sp_tile_PiViVi", "psfex_interp_runner") + psf = PSFExInterpolator(None, str(cat), str(psfout), number, log, + pos_params=["XWIN_WORLD", "YWIN_WORLD"]) + psf.process_me([str(tmp_path)], "psf", headers) + viout = output("run_sp_tile_PiViVi", "vignetmaker_runner") + vi = VignetMaker(str(cat), "SPHE", ["XWIN_WORLD", "YWIN_WORLD"], + str(viout), number, log) + vi.process_me([[str(tmp_path)]], ["image", "weight", "flag"], headers, 1) + bgout = viout + vi = VignetMaker(str(cat), "SPHE", ["XWIN_WORLD", "YWIN_WORLD"], + str(bgout), number, log) + vi.process_me([[str(tmp_path)]], ["exp_background", "exp_background_rms"], headers, 1) + stores = [viout / f"{name}_vignet{number}.sqlite" + for name in ["image", "weight", "flag"]] + psfpath = psfout / f"galaxy_psf{number}.sqlite" + bgpath = bgout / f"exp_background_vignet{number}.sqlite" + for path in [*stores, psfpath, bgpath]: + with SqliteDict(str(path)) as db: + assert len(db) == 0 + monkeypatch.setenv("SP_PSF", "psfex") + manifest("tile_vignets") + + # Run both chunks, merge their five empty shear HDUs, assemble HDF5. + for k in range(1, n_chunks + 1): + subdir = f"run_sp_tile_ngmix_Ng{k}u" + ngout = output(subdir, "ngmix_runner") + log.messages.clear() + config = _runner_config() + config.set("NGMIX_RUNNER", "BLEND_HANDLING", "uberseg") + config.set("NGMIX_RUNNER", "ID_OBJ_MIN", str(k)) + config.set("NGMIX_RUNNER", "ID_OBJ_MAX", str(k)) + assert ngmix_runner( + [str(cat), str(stores[0]), str(bgpath), str(psfpath), + str(stores[1]), str(stores[2]), headers], + {"output": str(ngout)}, number, config, "NGMIX_RUNNER", log, + ) == (None, None) + logs = ngout.parent / "logs" + logs.mkdir() + (logs / f"process{number}.log").write_text("\n".join(log.messages)) + manifest("tile_ngmix", subdir, f"_{k}") + with fits.open(ngout / f"ngmix{number}.fits") as hdus: + assert len(hdus) == 6 + assert all(len(h.data) == 0 for h in hdus[1:]) + merged = output("run_sp_tile_Ms", "merge_sep_cats_runner") + first = run / "output/run_sp_tile_ngmix_Ng1u/ngmix_runner/output" / f"ngmix{number}.fits" + MergeSep([str(first)], number, ["ngmix"], [".fits"], str(merged), + n_chunks, "error", log).process() + manifest("tile_merge_cats") + finalout = output("run_sp_tile_Mc", "make_cat_runner") + config = CustomParser() + options = {"SHAPE_MEASUREMENT_TYPE": "ngmix", "SAVE_PSF_DATA": "True"} + if epoch_slots is not None: + options["N_EPOCH_SLOTS"] = str(epoch_slots) + config.read_dict({"MAKE_CAT_RUNNER": options}) + assert make_cat_runner( + [str(cat), str(psfpath), str(merged / f"ngmix{number}.fits")], + {"output": str(finalout)}, number, config, "MAKE_CAT_RUNNER", log, + ) == (None, None) + final = finalout / f"final_cat{number}.hdf5" + with h5py.File(final) as h5: + assert "TILE_ID" in h5 and "NGMIX_MCAL_FLAGS" in h5 + assert len(h5) > 50 + assert all(d.shape[0] == 0 for d in h5.values()) + monkeypatch.setenv("NGMIX_N_CHUNKS", str(n_chunks)) + made = manifest("tile_make_cat") + assert made["catalogue_match"] == counts + assert made["epoch_cuts"]["totals"] == dict( + considered=0, masked_fraction=0, central_veto=0, failed=0, objects_emptied=0) + # Same persistent products as the publication rule (whose commands are + # exercised separately by test_match_record_publication). + published = products_dir / "tiles" / tile[:2] / tile + published.mkdir(parents=True) + shutil.copyfile(final, published / f"final_cat-{tile}.hdf5") + shutil.copyfile(mdir / "tile_make_cat.json", published / "tile_make_cat.json") + assert json.loads((published / "tile_make_cat.json").read_text())["catalogue_match"]["n_sextractor"] == 0 + + # Campaign merge uses the shipped reader and reconciler, also used by + # sp_latest; the tile remains present as a zero-length structured dataset. + tiles = tmp_path / "tiles.txt" + tiles.write_text(tile + "\n") + index = tmp_path / "index.sqlite" + with sqlite3.connect(index) as db: + db.execute("CREATE TABLE tiles(tile_id TEXT PRIMARY KEY, ra_dir TEXT, n_exp INTEGER)") + db.execute("INSERT INTO tiles VALUES (?, ?, 1)", (tile, tile.split(".")[0])) + with h5py.File(final) as h5: + columns = list(h5) + params = tmp_path / "final_cat.param" + params.write_text("\n".join(columns) + "\n") + campaign_cat = products_dir / "final_cat_empty.hdf5" + argv = [sys.executable, str(scripts / "merge_final_cat.py"), + "--products-dir", str(products_dir), "--tile-list", str(tiles), + "--index-db", str(index), "--output", str(campaign_cat), + "--campaign", "empty", "--param-file", str(params)] + merged = subprocess.run(argv, capture_output=True, text=True) + assert merged.returncode == 0, merged.stdout + merged.stderr + with h5py.File(campaign_cat) as h5: + data = h5[f"patches/empty/{tile}"] + assert data.shape == (0,) + assert set(data.dtype.names) == set(columns) + assert h5.attrs["n_tiles"] == 1 + # An unchanged campaign must also pass the incremental-merge path. + rerun = subprocess.run(argv, capture_output=True, text=True) + assert rerun.returncode == 0, rerun.stderr + assert "unchanged" in rerun.stdout diff --git a/tests/unit/test_ngmix_range.py b/tests/unit/test_ngmix_range.py index fe2bea719..53a4be32a 100644 --- a/tests/unit/test_ngmix_range.py +++ b/tests/unit/test_ngmix_range.py @@ -44,8 +44,11 @@ def _load(): def assert_tiles(ranges, n_obj, n_chunks): """Assert the partition invariant, in the form each failure mode takes.""" assert len(ranges) == n_chunks - assert ranges[0][0] == 1 - assert ranges[-1][1] == n_obj + if n_obj: + assert ranges[0][0] == 1 + assert ranges[-1][1] == n_obj + else: + assert ranges == [(2, 1)] * n_chunks for (lo, hi), (prev_lo, prev_hi) in zip(ranges[1:], ranges): assert lo == prev_hi + 1 assert hi >= lo - 1, "an empty range is (hi + 1, hi), never narrower" @@ -64,7 +67,7 @@ def assert_tiles(ranges, n_obj, n_chunks): @settings(deadline=None) @given( - n_obj=st.integers(min_value=1, max_value=400), + n_obj=st.integers(min_value=0, max_value=400), n_chunks=st.integers(min_value=1, max_value=16), data=st.data(), ) @@ -82,7 +85,7 @@ def test_ranges_tile_the_catalogue(n_obj, n_chunks, data): @settings(deadline=None) @given( - n_obj=st.integers(min_value=1, max_value=60), + n_obj=st.integers(min_value=0, max_value=60), n_chunks=st.integers(min_value=1, max_value=8), data=st.data(), ) @@ -116,7 +119,7 @@ def _roundtrip(epochs, n_chunks): @settings(deadline=None) @given( - n_obj=st.integers(min_value=1, max_value=200), + n_obj=st.integers(min_value=0, max_value=200), n_chunks=st.integers(min_value=1, max_value=12), data=st.data(), ) @@ -190,10 +193,13 @@ def test_fewer_objects_than_chunks_pads_with_empty_ranges(): assert_tiles(ngmix_range.row_ranges([2, 5, 1], 6), 3, 6) -def test_zero_objects_is_fatal(): - """No split is meaningful, and (1, 0) is ngmix's unbounded sentinel.""" - with pytest.raises(ValueError, match="zero objects"): - ngmix_range.row_ranges([], 8) +@pytest.mark.parametrize("n_chunks", [1, 8, 16]) +def test_zero_objects_preserves_configured_empty_chunks(n_chunks): + """All scheduled readers get a bounded empty range and can write products.""" + doc = _roundtrip([], n_chunks) + assert doc["n_obj"] == 0 and doc["n_chunks"] == n_chunks + ranges = [ngmix_range.chunk_range(doc, k) for k in range(1, n_chunks + 1)] + assert_tiles(ranges, 0, n_chunks) def test_zero_chunks_is_fatal(): diff --git a/tests/workflow/test_zero_object_tile.py b/tests/workflow/test_zero_object_tile.py new file mode 100644 index 000000000..ba1578cbb --- /dev/null +++ b/tests/workflow/test_zero_object_tile.py @@ -0,0 +1,90 @@ +"""A zero-row catalogue uses the resolved DAG's real partition commands. + +Planning alone missed Fir 177.257: tile_vignets' prelude runs ngmix_range +before the modules, and it used to reject the otherwise valid empty tile. +Execute that writer and every scheduled reader from the resolved job params, +then feed their exported bounds to the real module-to-publication regression. +""" + +import json +import os +from pathlib import Path +import subprocess + +import h5py + +from tests.module import test_zero_object_tile as modules +from tests.workflow.harness import Campaign, resolve + + +def test_zero_object_tile_dag_partition_through_publication(tmp_path, monkeypatch): + campaign = Campaign(tmp_path / "campaign", "data", "psfex") + campaign.config["ngmix_chunks"] = 16 # Fir's configured static fan-out + campaign.write_config() + tile = "124.456" # one exposure in the campaign fixture + run = campaign.run_dir / "tiles" / tile[:2] / tile + local = tmp_path / "node-local" + local.mkdir() + selected = lambda jobs: [j for j in jobs if j.wildcards.tile == tile] + seen = [] + real_runner = modules.ngmix_runner + + with resolve(campaign, monkeypatch) as dag: + vignets, = selected(dag.jobs_for("tile_vignets")) + chunks = selected(dag.jobs_for("tile_ngmix")) + merged, = selected(dag.jobs_for("tile_merge_cats")) + made, = selected(dag.jobs_for("tile_make_cat")) + assert len(chunks) == 16 + assert set(map(str, merged.input.manifests)) == { + str(j.output.manifest) for j in chunks} + assert str(merged.output.manifest) in map(str, made.input) + assert all(str(vignets.output.manifest) in map(str, j.input) + for j in chunks) + assert str(campaign.final_cat(tile)) in map(str, made.output) + writer, = [line for line in vignets.params.pre.splitlines() + if "ngmix_range.py" in line and "--write" in line] + readers = {} + for job in chunks: + lookup, = [line for line in job.params.pre.splitlines() + if "ngmix_range.py" in line and "--read" in line] + assert 'eval "$ngmix_range_out"' in job.params.pre + readers[int(job.wildcards.chunk)] = lookup + + def scheduled_runner(inputs, dirs, number, config, section, log): + # The fixture has now produced its real zero-row LDAC/EPOCH tables. + # Run the exact shell command Snakemake resolved for tile_vignets, + # not partition() directly or a test copy of its CLI arguments. + env = dict(os.environ, SP_RUN=str(run), SP_LOCAL=str(local)) + if not seen: + written = subprocess.run(["bash", "-ec", writer], env=env, + capture_output=True, text=True) + assert written.returncode == 0, written.stdout + written.stderr + doc = json.loads((local / "ngmix_ranges.json").read_text()) + assert doc["n_obj"] == 0 and doc["n_chunks"] == 16 + assert len(doc["chunks"]) == len(chunks) + k = len(seen) + 1 + command = (readers[k] + '\neval "$ngmix_range_out"\n' + + 'printf "%s %s" "$NGMIX_ROW_MIN" "$NGMIX_ROW_MAX"') + read = subprocess.run(["bash", "-ec", command], env=env, + capture_output=True, text=True) + assert read.returncode == 0, read.stdout + read.stderr + lo, hi = map(int, read.stdout.split()) + assert (lo, hi) == (2, 1) # bounded empty interval, not unbounded 0 + config.set(section, "ID_OBJ_MIN", str(lo)) + config.set(section, "ID_OBJ_MAX", str(hi)) + seen.append(k) + return real_runner(inputs, dirs, number, config, section, log) + + monkeypatch.setattr(modules, "ngmix_runner", scheduled_runner) + modules.test_zero_object_tile_modules_and_publication( + run.parent, monkeypatch, 12, tile=tile, n_chunks=len(chunks), + products_dir=campaign.products_dir, exposure=campaign.ready[tile][0]) + assert seen == list(range(1, 17)) + for job in [vignets, *chunks, merged, made]: + manifest = Path(job.output.manifest) + assert json.loads(manifest.read_text())["status"] == "complete" + with h5py.File(campaign.final_cat(tile)) as cat: + assert all(column.shape[0] == 0 for column in cat.values()) + stats = json.loads(Path(made.output.stats).read_text()) + assert stats["catalogue_match"]["n_sextractor"] == 0 + assert stats["epoch_cuts"]["n_chunks"] == 16 diff --git a/workflow/README.md b/workflow/README.md index 0d9b0eb1c..897e2a204 100644 --- a/workflow/README.md +++ b/workflow/README.md @@ -376,6 +376,21 @@ profiles/nibi/config.yaml SLURM executor; apptainer SDM; per-user jobs cap; kee This requires an ngmix build with PR #922's end-of-loop `epoch cuts:` summary, including on empty chunks. Older code without that summary fails the ngmix completeness check and leaves scratch intact. + Every tile joined to the external catalogue also publishes `catalogue_match` in this same `tile_make_cat.json` (under `/tiles///`). + It contains integer `n_sextractor`, `n_external`, `n_paired`, `n_external_only`, `n_dropped`, plus `external_paired_fraction` (zero for an empty external catalogue) and `external_warning_budget`. + The budget is `max(20, 0.02 * n_external + 50)` when any pairs exist, otherwise `max(20, 0.02 * n_external)`; exceeding it only warns. + The SExtractor-side mismatch still fails hard. + The record travels from detection FITS headers through `tile_detect.json` before scratch cleanup, without changing catalogue rows. + Join the manifest's `unit` converted to float to the final catalogue's `TILE_ID` to cut downstream; records exist for passing and warning tiles alike, including zero-object tiles. + Runs without an external join omit `catalogue_match`. + A zero-detection tile is valid: detection writes zero-row `LDAC_OBJECTS`/`EPOCH_*` tables and an explicitly shaped `SEG_VIGNET`; pairing records zero pairs and warns when above budget. + PSF interpolation and vignette extraction write empty SQLite stores, each ngmix chunk and their merge write five zero-row shear HDUs with zero epoch-cut counts, and make-cat publishes an HDF5 catalogue with the normal columns and zero rows. + Completeness counts these files, not their rows; the campaign merge retains the tile as a zero-length structured dataset rather than losing its identity. + The workflow keeps its configured static ngmix fan-out even for zero objects: `ngmix_range.py` writes a normal ranges document with every chunk set to the positive empty interval `(2, 1)`. + This avoids the unbounded `ROW_MAX <= 0` sentinel and lets every scheduled chunk publish its normal empty product, so merge and completeness need no separate zero-chunk graph. + Nonzero partitions are unchanged. + For example, real b619 tile 058.220 has `(n_sextractor, n_external, n_paired, n_external_only, n_dropped) = (774, 845, 774, 71, 0)`, fraction `774/845 = 0.9159763313609467`, budget `66.9`, and now warns rather than failing. + For an **idle existing campaign** with kept logs, use `workflow/scripts/backfill_epoch_cuts.py --run-dir --products-dir --n-chunks --dry-run` first. Alternatively pass `--config ` to resolve products and chunk count from the run's frozen `workflow/config.yaml` plus its run config. The script skips tiles without a nonempty persisted catalogue and refuses missing or ambiguous chunk evidence. diff --git a/workflow/config/cfis/config_tile_Sx.ini b/workflow/config/cfis/config_tile_Sx.ini index 8d11674ef..d30ebff2a 100644 --- a/workflow/config/cfis/config_tile_Sx.ini +++ b/workflow/config/cfis/config_tile_Sx.ini @@ -148,20 +148,13 @@ MATCH_RADIUS = 1.0 # SExtractor-only losses must exceed both the fractional budget and the # absolute tolerance to stop. This catches other pixels (DR5 pairs ~95%). -# External-only losses also get the additive allowance below: sparse DR6 -# tiles lose localized children in star halos, even at positive weight. -# This is an empirical loss budget, not an independent-object probability. +# External-only losses warn and are recorded for downstream cuts. # @sc [decision:detection.tile_detection] MATCH_MIN_FRACTION = 0.98 # @sc [decision:detection.tile_detection] MATCH_TOLERATED_UNPAIRED = 20 -# Extra external-only objects allowed on top of (1 - MATCH_MIN_FRACTION)*N. -# Requires at least one pair; zero restores the original two-sided guard. -# @sc [decision:detection.tile_detection] -MATCH_EXTERNAL_UNPAIRED_ALLOWANCE = 50 - ## Post-processing # Necessary for tiles, to enable multi-exposure processing diff --git a/workflow/scripts/completeness.py b/workflow/scripts/completeness.py index 7596b7b78..84955a37a 100644 --- a/workflow/scripts/completeness.py +++ b/workflow/scripts/completeness.py @@ -343,6 +343,24 @@ def build_manifest(stage, run_dir, unit, stage_subdir=None): "warn": warn, "status": status, "reasons": scrape_reasons(stage_dir, runner), }) + # Capture the join before detection scratch is reclaimed, and carry it + # into the existing published make-cat manifest beside final_cat. + if ok and stage == "tile_detect": + from shapepipe.modules.sextractor_package.match_catalogue import ( + read_match_record, + ) + cats = sorted(stage_dir.glob("sextractor_runner/output/sexcat*.fits")) + if len(cats) == 1: + record = read_match_record(cats[0]) + if record is not None: + manifest["catalogue_match"] = record + if ok and stage == "tile_make_cat": + detect = run_dir / "manifests" / "tile_detect.json" + if detect.exists(): + record = json.loads(detect.read_text()).get("catalogue_match") + if record is not None: + manifest["catalogue_match"] = record + # Capture ngmix's diagnostics while its temp() directory still exists. # Make-cat consumes only structured records, never the retained log copies. if ok and stage in ("tile_ngmix", "tile_make_cat"): diff --git a/workflow/scripts/ngmix_range.py b/workflow/scripts/ngmix_range.py index 173cb51f7..6d3fac443 100644 --- a/workflow/scripts/ngmix_range.py +++ b/workflow/scripts/ngmix_range.py @@ -73,19 +73,20 @@ def row_ranges(epochs, n_chunks: int) -> list[tuple[int, int]]: which is arbitrary but fixed, and fixed is the property that matters. With ``n_obj < n_chunks`` the first ``n_obj`` chunks take one object each - and the remainder are empty, written ``(n_obj + 1, n_obj)``. This preserves - ``ranges[k][0] == ranges[k-1][1] + 1`` and keeps the upper bound positive - as required by the module's row-partition contract. + and the remainder are empty, written ``(max(1, n_obj) + 1, max(1, n_obj))``. + For nonempty catalogues this preserves adjacency to the last object. + With no objects every configured chunk receives ``(2, 1)``: a positive, + inverted interval that selects nothing, never ngmix's unbounded sentinel. + Keeping the configured chunk count lets the static DAG, merge and + completeness use their normal empty-product paths without a second graph. """ if n_chunks < 1: raise ValueError(f"n_chunks must be >= 1, got {n_chunks}") n_obj = len(epochs) - if n_obj < 1: - raise ValueError("cannot split a catalogue of zero objects") - if n_obj < n_chunks: + empty_hi = max(1, n_obj) return ([(k, k) for k in range(1, n_obj + 1)] - + [(n_obj + 1, n_obj)] * (n_chunks - n_obj)) + + [(empty_hi + 1, empty_hi)] * (n_chunks - n_obj)) weights = [MILLI_EPOCH * int(e) + ALPHA_MILLI_EPOCHS for e in epochs] cap = _min_feasible_capacity(weights, n_chunks)