Change name of "used" field and add ability to reserve candidates

Flag marks sources using calib_photometryUsed flag: changed from photometryStandard.

Add feature to reserve a fraction of candidates selected from the initial match.
Mark those calib_photometryReserved.

Mark the remainder of the candidates calib_photometryCandidate.

Unit test added to testPhotoCal.py

python/lsst/pipe/tasks/photoCal.py

@@ -27,6 +27,7 @@
 from builtins import range
 
 import math
+import random
 import sys
 
 import numpy as np
@@ -67,10 +68,16 @@ class PhotoCalConfig(RefMatchConfig):
         default=22.0,
         doc="Don't use objects fainter than this magnitude",
     )
-    doWriteOutput = pexConf.Field(
-        dtype=bool,
-        default=True,
-        doc="Write a field name astrom_usedByPhotoCal to the schema",
+    reserveFraction = pexConf.Field(
+        dtype=float,
+        doc="Fraction of candidates to reserve from fitting; none if <= 0",
+        default=-1.0,
+    )
+    reserveSeed = pexConf.Field(
+        dtype = int,
+        doc = "This number will be multiplied by the exposure ID "
+        "to set the random seed for reserving candidates",
+        default = 1,
     )
     fluxField = pexConf.Field(
         dtype=str,
@@ -290,14 +297,20 @@ def DebugInfo(name):
     def __init__(self, refObjLoader, schema=None, **kwds):
         """!Create the photometric calibration task.  See PhotoCalTask.init for documentation
         """
-        RefMatchTask.__init__(self, refObjLoader, schema=schema, **kwds)
+        RefMatchTask.__init__(self, refObjLoader, schema=None, **kwds)
         self.scatterPlot = None
         self.fig = None
-        if self.config.doWriteOutput:
-            self.outputField = schema.addField("photocal_photometricStandard", type="Flag",
-                                               doc="set if source was used in photometric calibration")
+        if schema is not None:
+            self.usedKey = schema.addField("calib_photometryUsed", type="Flag",
+                                           doc="set if source was used in photometric calibration")
+            self.candidateKey = schema.addField("calib_photometryCandidate", type="Flag",
+                                                doc="set if source was a candidate for use in calibration")
+            self.reservedKey = schema.addField("calib_photometryReserved", type="Flag",
+                                               doc="set if source was reserved, so not used in calibration")
         else:
-            self.outputField = None
+            self.usedKey = None
+            self.candidateKey = None
+            self.reservedKey = None
 
     def getSourceKeys(self, schema):
         """!Return a struct containing the source catalog keys for fields used by PhotoCalTask.
@@ -365,6 +378,10 @@ def selectMatches(self, matches, sourceKeys, filterName, frame=None):
         if len(matches) == 0:
             raise ValueError("No input matches")
 
+        for m in matches:
+            if self.candidateKey is not None:
+                m.second.set(self.candidateKey, True)
+
         # Only use stars for which the flags indicate the photometry is good.
         afterFlagCutInd = [i for i, m in enumerate(matches) if checkSourceFlags(m.second, sourceKeys)]
         afterFlagCut = [matches[i] for i in afterFlagCutInd]
@@ -456,8 +473,8 @@ def selectMatches(self, matches, sourceKeys, filterName, frame=None):
 
         result = []
         for m in matches:
-            if self.outputField is not None:
-                m.second.set(self.outputField, True)
+            if self.usedKey is not None:
+                m.second.set(self.usedKey, True)
             result.append(m)
         return result
 
@@ -577,7 +594,7 @@ def extractMagArrays(self, matches, filterName, sourceKeys):
         )
 
     @pipeBase.timeMethod
-    def run(self, exposure, sourceCat):
+    def run(self, exposure, sourceCat, expId=0):
         """!Do photometric calibration - select matches to use and (possibly iteratively) compute
         the zero point.
 
@@ -648,21 +665,41 @@ def run(self, exposure, sourceCat):
             frame = None
 
         res = self.loadAndMatch(exposure, sourceCat)
+
+        #from res.matches, reserve a fraction of the population and mark the sources as reserved
+
+        if self.config.reserveFraction > 0:
+            random.seed(self.config.reserveSeed*expId)
+            reserveList = random.sample(res.matches,
+                                        int((self.config.reserveFraction)*len(res.matches)))
+
+            for candidate in reserveList:
+                res.matches.remove(candidate)
+
+            if reserveList and self.reservedKey is not None:
+                for candidate in reserveList:
+                    candidate.second.set(self.reservedKey, True)
+
         matches = res.matches
+        for m in matches:
+            if self.candidateKey is not None:
+                m.second.set(self.candidateKey, True)
+
         filterName = exposure.getFilter().getName()
         sourceKeys = self.getSourceKeys(matches[0].second.schema)
+
         matches = self.selectMatches(matches=matches, sourceKeys=sourceKeys, filterName=filterName,
                                      frame=frame)
         arrays = self.extractMagArrays(matches=matches, filterName=filterName, sourceKeys=sourceKeys)
 
-        if matches and self.outputField:
+        if matches and self.usedKey:
             try:
                 # matches[].second is a measured source, wherein we wish to set outputField.
                 # Check that the field is present in the Sources schema.
-                matches[0].second.getSchema().find(self.outputField)
+                matches[0].second.getSchema().find(self.usedKey)
             except:
-                raise RuntimeError("sources' schema does not contain the used-in-calibration flag \"%s\"" %
-                                   self.outputField)
+                raise RuntimeError("sources' schema does not contain the calib_photometryUsed flag \"%s\"" %
+                                   self.usedKey)
 
         # Fit for zeropoint.  We can run the code more than once, so as to
         # give good stars that got clipped by a bad first guess a second

tests/testPhotoCal.py

@@ -92,7 +92,6 @@ def setUp(self):
         # The test and associated data have been prepared on the basis that we
         # use the PsfFlux to perform photometry.
         self.config.fluxField = "base_PsfFlux_flux"
-        self.config.doWriteOutput = False    # schema is fixed because we already loaded the data
 
     def tearDown(self):
         del self.srcCat
@@ -102,7 +101,7 @@ def tearDown(self):
 
     def _runTask(self):
         """All the common setup to actually test the results"""
-        task = PhotoCalTask(self.refObjLoader, config=self.config, schema=self.srcCat.schema)
+        task = PhotoCalTask(self.refObjLoader, config=self.config, schema=None)
         pCal = task.run(exposure=self.exposure, sourceCat=self.srcCat)
         matches = pCal.matches
         print("Ref flux fields list =", pCal.arrays.refFluxFieldList)
@@ -125,6 +124,29 @@ def _runTask(self):
         self.zp = pCal.calib.getMagnitude(1.)
         self.fitdiff = pCal.arrays.srcMag + self.zp - pCal.arrays.refMag
 
+    def testFlags(self):
+        """test that all the calib_photometry flags are set to reasonable values"""
+        schema = self.srcCat.schema
+        task = PhotoCalTask(self.refObjLoader, config=self.config, schema=schema)
+        mapper = afwTable.SchemaMapper(self.srcCat.schema, schema)
+        cat = afwTable.SourceCatalog(schema)
+        for name in self.srcCat.schema.getNames():
+            mapper.addMapping(self.srcCat.schema.find(name).key)
+        cat.extend(self.srcCat, mapper=mapper)
+
+        #   test that by default, no stars are reserved and used < candidates
+        task.run(exposure=self.exposure, sourceCat=cat)
+        candidates = cat.get("calib_photometryCandidate").sum()
+        self.assertEqual(cat.get("calib_photometryReserved").sum(), 0)
+        self.assertLessEqual(cat.get("calib_photometryUsed").sum(), candidates)
+
+        #   set the reserve fraction, and see if the right proportion are reserved.
+        self.config.reserveFraction = .3
+        task.run(exposure=self.exposure, sourceCat=cat)
+        reserved = cat.get("calib_photometryReserved").sum()
+        self.assertEqual(reserved, int(.3 * candidates))
+        self.assertLessEqual(cat.get("calib_photometryUsed").sum(), (candidates - reserved))
+
     def testZeroPoint(self):
         """ Test to see if we can compute a photometric zeropoint given a reference task"""
         self._runTask()