N_EPOCH counts only PSF-validated epochs

astra.yaml

@@ -733,7 +733,9 @@ analyses:
           admits its upper endpoint, column 2080. [HARDCODED] The comparisons
           are strict in make_post_process; CCD_SIZE sets only the bounds. A
           WCS inversion failure skips that CCD, lowering N_EPOCH, which caps
-          how many exposures enter each galaxy's multi-epoch fit.
+          how many exposures enter each galaxy's multi-epoch fit. This geometric
+          count remains internal; the final catalogue's N_EPOCH counts only
+          epochs with a validated PSF.
           Values:
           SEXTRACTOR_RUNNER.CCD_SIZE = 33,2080,1,4612;
           SEXTRACTOR_RUNNER.MAKE_POST_PROCESS = True.
@@ -1173,6 +1175,12 @@ analyses:
           22-star floor to that accepted count. In the science path the CCD's
           epoch is dropped for every object on it, and an object left with no
           epoch has no shape.
+          The final catalogue's N_EPOCH counts the epochs with an
+          interpolated PSF (for psfex, those passing this gate; for mccd,
+          those with a model). Only this count is written as N_EPOCH: one column
+          avoids downstream cuts on the geometric count by mistake.
+          PSF-gate rejections are studied per CCD in validation output and
+          the exposure maps (#887), rather than per object.
           The pipeline has no minimum-epoch floor: NGMIX_N_EPOCH records what
           survived and epoch-count cuts happen downstream. Guinot+22 discards
           the CCD from PSF estimation rather than gating at interpolation.

example/cfis_image_sims/final_cat.param

@@ -31,7 +31,8 @@ NGMIX_MCAL_TYPES_FAIL
 NGMIX_G1_PSF_ORIG_NOSHEAR
 NGMIX_G2_PSF_ORIG_NOSHEAR
 
-# Number of epochs (exposures)
+# Number of epochs: N_EPOCH with a validated PSF model,
+# NGMIX_N_EPOCH fitted by ngmix
 N_EPOCH
 NGMIX_N_EPOCH
 

src/shapepipe/modules/make_cat_package/__init__.py

@@ -21,7 +21,9 @@
 previous module runs. This gathers all relevant information on the measured
 galaxies for weak-lensing post-processing. This includes galaxy detection and
 basic measurement parameters, the PSF model at galaxy positions, and the
-shape measurement. Every detected object is kept: the catalogue carries no
+shape measurement. ``N_EPOCH`` counts the object's epochs with a validated
+PSF model, the epochs ngmix can fit.
+The geometric count of covering CCDs is used internally, not written out. Every detected object is kept: the catalogue carries no
 star/galaxy classification, which is done downstream. Each object is
 tagged with its source tile via the ``TILE_ID`` column and carries the
 survey-wide object ID ``TILE_UNIQUE_ID = tile_id * 10**6 + NUMBER``
@@ -43,7 +45,8 @@
     (``HSM_*_PSF_n``, ``EXP_ID_n``, ``CCD_n``) when ``SAVE_PSF_DATA`` is
     ``True``; unfilled slots hold the family's sentinel, and an object with
     more epochs than slots raises an error. A fixed count gives every tile
-    the same schema. Default is the tile's maximum ``N_EPOCH`` plus one
+    the same schema. Default is the sexcat's maximum geometric ``N_EPOCH`` plus
+    one
 
 """
 

src/shapepipe/modules/make_cat_package/make_cat.py

@@ -81,6 +81,12 @@ def read_sextractor_data(sexcat_path, remove_vignet=True):
     ``NUMBER``). The tile is read from the catalogue's file name, e.g.
     ``sexcat-301-279.fits``.
 
+    The SExtractor catalogue's ``N_EPOCH`` counts the exposure CCDs whose
+    footprint holds the object, before any PSF model exists; it is used
+    internally for slot sizing and overwritten before output. The final catalogue's ``N_EPOCH`` counts only the
+    epochs with a validated PSF and is written by
+    :meth:`SaveCatalogue.process` in ``psf`` mode.
+
     Parameters
     ----------
     sexcat_path : str
@@ -214,6 +220,7 @@ def process(
         cat_path=None,
         moments=False,
         n_epoch_slots=None,
+        epoch_slots=True,
     ):
         """Process Catalogue.
 
@@ -227,7 +234,10 @@ def process(
             Option to run ``ngmix`` mode with moments
         n_epoch_slots : int, optional
             Number of per-epoch slots in ``psf`` mode; if ``None``, the
-            tile's ``max(N_EPOCH) + 1``
+            sexcat's geometric ``max(N_EPOCH) + 1``
+        epoch_slots : bool, optional
+            In ``psf`` mode, write the per-epoch column families as well as
+            ``N_EPOCH``; default ``True``
 
         Returns
         --------
@@ -242,7 +252,7 @@ def process(
         if mode == "ngmix":
             err_msg = self._save_ngmix_data(cat_path, moments)
         elif mode == "psf":
-            self._save_psf_data(cat_path, n_epoch_slots)
+            self._save_psf_data(cat_path, n_epoch_slots, epoch_slots)
         else:
             err_msg = (
                 f"Invalid process mode ({mode}) for "
@@ -549,14 +559,24 @@ def _save_ngmix_data(self, ngmix_cat_path, moments=False):
 
         return None
 
-    def _save_psf_data(self, galaxy_psf_path, n_epoch_slots=None):
+    def _save_psf_data(
+        self, galaxy_psf_path, n_epoch_slots=None, epoch_slots=True
+    ):
         """Save PSF data.
 
-        Save the PSF catalogue into the final one, as per-epoch column
+        Save the PSF catalogue into the final one: the per-object epoch
+        count ``N_EPOCH`` and, if ``epoch_slots``, the per-epoch column
         families ``HSM_*_PSF_n`` (``psf_shape_cols``), ``EXP_ID_n`` and
-        ``CCD_n``. Slot ``n`` of every
-        family refers to the same epoch; slots no epoch fills keep the
-        family's sentinel.
+        ``CCD_n``. Slot ``n`` of every family refers to the same epoch;
+        slots no epoch fills keep the family's sentinel.
+
+        ``N_EPOCH`` counts the object's entries in the PSF catalogue, i.e.
+        the epochs whose CCD passed PSF-model validation and so have an
+        interpolated PSF; ngmix fits no other epoch. It is at most
+        the sexcat's geometric ``N_EPOCH``, the count of CCDs whose footprint holds
+        the object, and 0 for an object the PSF catalogue marks
+        ``"empty"``. An epoch whose PSF shape fit failed
+        (``HSM_FLAG_PSF != 0``) still counts: its PSF exists.
 
         @sc [label:schema] psf-epoch-slot-columns
         The per-epoch families are ``psf_shape_cols`` plus ``EXP_ID``/``CCD``,
@@ -570,7 +590,10 @@ def _save_psf_data(self, galaxy_psf_path, n_epoch_slots=None):
             Path to the PSF catalogue to save
         n_epoch_slots : int, optional
             Number of slots written per family; if ``None``, the tile's
-            ``max(N_EPOCH) + 1``
+            geometric ``max(N_EPOCH) + 1`` from the sexcat, which bounds every object's PSF
+            epoch count
+        epoch_slots : bool, optional
+            Write the per-epoch column families; default ``True``
 
         Raises
         ------
@@ -580,12 +603,15 @@ def _save_psf_data(self, galaxy_psf_path, n_epoch_slots=None):
         """
         galaxy_psf_cat = SqliteDict(galaxy_psf_path)
 
-        n_epoch = self._final_cat["N_EPOCH"]
-        if n_epoch_slots is None:
-            n_slots = np.max(n_epoch) + 1
+        n_obj = len(self._obj_id)
+        n_epoch = np.zeros(n_obj, dtype="int32")
+
+        if not epoch_slots:
+            n_slots = 0
+        elif n_epoch_slots is None:
+            n_slots = np.max(self._final_cat["N_EPOCH"]) + 1
         else:
             n_slots = n_epoch_slots
-        n_obj = len(self._obj_id)
 
         # Per-epoch PSF shape columns copied from the producer's SHAPES dict:
         # (column, empty-slot fill, dtype). Fills are out of physical range
@@ -613,15 +639,27 @@ def _save_psf_data(self, galaxy_psf_path, n_epoch_slots=None):
 
         for idx, id_tmp in enumerate(self._obj_id):
 
-            obj_epochs = galaxy_psf_cat[str(id_tmp)]
+            try:
+                obj_epochs = galaxy_psf_cat[str(id_tmp)]
+            except KeyError:
+                galaxy_psf_cat.close()
+                raise KeyError(
+                    f"Object {id_tmp} is missing from {galaxy_psf_path}; "
+                    "the galaxy PSF catalogue and the SExtractor catalogue "
+                    "do not match"
+                ) from None
             if obj_epochs == "empty":
                 continue
 
+            n_epoch[idx] = len(obj_epochs)
+            if not epoch_slots:
+                continue
+
             if len(obj_epochs) > n_slots:
                 galaxy_psf_cat.close()
                 raise ValueError(
-                    f"Object {id_tmp} has {len(obj_epochs)} PSF epochs"
-                    + f" (N_EPOCH={n_epoch[idx]}), more than the {n_slots}"
+                    f"Object {id_tmp} has N_EPOCH={n_epoch[idx]} PSF epochs,"
+                    + f" more than the {n_slots}"
                     + f" per-epoch slots (N_EPOCH_SLOTS={n_epoch_slots})"
                 )
 
@@ -649,3 +687,5 @@ def _save_psf_data(self, galaxy_psf_path, n_epoch_slots=None):
                     )
 
         galaxy_psf_cat.close()
+
+        self._output_dict = {"N_EPOCH": n_epoch, **self._output_dict}

src/shapepipe/modules/make_cat_runner.py

@@ -76,8 +76,14 @@ def make_cat_runner(
             w_log.info(err_msg)
             return None, None
 
-    if save_psf:
-        sc_inst.process("psf", galaxy_psf_path, n_epoch_slots=n_epoch_slots)
+    # N_EPOCH (epochs with a validated PSF) always; the per-epoch PSF
+    # columns only when SAVE_PSF_DATA is set.
+    sc_inst.process(
+        "psf",
+        galaxy_psf_path,
+        n_epoch_slots=n_epoch_slots,
+        epoch_slots=save_psf,
+    )
 
     # Optional per-band external healsparse mask lookup (UNIONS-WL/spherex#38):
     # add one MASK_<BAND> column per band, queried at each object's world