Merge pull request #133 from lsst/tickets/DM-47919

DM-47919: Add configuration for cutting visits from the fit by PSF FWHM

python/lsst/fgcmcal/fgcmBuildStarsBase.py

@@ -33,6 +33,7 @@
 from lsst.meas.algorithms.sourceSelector import sourceSelectorRegistry
 
 from .fgcmLoadReferenceCatalog import FgcmLoadReferenceCatalogTask
+from .utilities import computeReferencePixelScale
 
 import fgcm
 
@@ -259,13 +260,15 @@ def fgcmMakeVisitCatalog(self, camera, groupedHandles):
         visitCat['used'] = 0
         visitCat['sources_read'] = False
 
+        defaultPixelScale = computeReferencePixelScale(camera)
+
         # No matter what, fill the catalog. This will check if it was
         # already read.
-        self._fillVisitCatalog(visitCat, groupedHandles)
+        self._fillVisitCatalog(visitCat, groupedHandles, defaultPixelScale)
 
         return visitCat
 
-    def _fillVisitCatalog(self, visitCat, groupedHandles):
+    def _fillVisitCatalog(self, visitCat, groupedHandles, defaultPixelScale):
         """
         Fill the visit catalog with visit metadata
 
@@ -275,6 +278,8 @@ def _fillVisitCatalog(self, visitCat, groupedHandles):
             Visit catalog.  See _makeFgcmVisitSchema() for schema definition.
         groupedHandles : `dict` [`list` [`lsst.daf.butler.DeferredDatasetHandle`]]
             Dataset handles, grouped by visit.
+        defaultPixelScale : `float`
+            Default pixel scale to use if not in visit summary (arcsecond/pixel).
         """
 
         # Guarantee that these are sorted.
@@ -292,15 +297,24 @@ def _fillVisitCatalog(self, visitCat, groupedHandles):
 
             visitInfo = summaryRow.getVisitInfo()
             physicalFilter = summaryRow['physical_filter']
-            # Compute the median psf sigma if possible
-            goodSigma, = np.where(summary['psfSigma'] > 0)
+            # Compute the median psf sigma and fwhm if possible.
+            if 'pixelScale' in summary.schema:
+                # This is not available in the older test summaries
+                pixelScales = summary['pixelScale']
+            else:
+                pixelScales = np.full(len(summary['psfSigma']), defaultPixelScale)
+            psfSigmas = summary['psfSigma']
+            goodSigma, = np.where((np.nan_to_num(psfSigmas) > 0) & (np.nan_to_num(pixelScales) > 0))
             if goodSigma.size > 2:
-                psfSigma = np.median(summary['psfSigma'][goodSigma])
+                psfSigma = np.median(psfSigmas[goodSigma])
+                psfFwhm = np.median(psfSigmas[goodSigma] * pixelScales[goodSigma]) * np.sqrt(8.*np.log(2.))
             elif goodSigma.size > 0:
-                psfSigma = summary['psfSigma'][goodSigma[0]]
+                psfSigma = psfSigmas[goodSigma[0]]
+                psfFwhm = psfSigmas[goodSigma[0]] * pixelScales[goodSigma[0]] * np.sqrt(8.)*np.log(2.)
             else:
                 self.log.warning("Could not find any good summary psfSigma for visit %d", visit)
                 psfSigma = 0.0
+                psfFwhm = 0.0
             # Compute median background if possible
             goodBackground, = np.where(np.nan_to_num(summary['skyBg']) > 0.0)
             if goodBackground.size > 2:
@@ -332,6 +346,7 @@ def _fillVisitCatalog(self, visitCat, groupedHandles):
             # Median delta aperture, to be measured from stars
             rec['deltaAper'] = 0.0
             rec['psfSigma'] = psfSigma
+            rec['psfFwhm'] = psfFwhm
             rec['skyBackground'] = skyBackground
             rec['used'] = 1
 
@@ -559,7 +574,8 @@ def _makeFgcmVisitSchema(self, nCcd):
         schema.addField('mjd', type=np.float64, doc="MJD of visit")
         schema.addField('exptime', type=np.float32, doc="Exposure time")
         schema.addField('pmb', type=np.float32, doc="Pressure (millibar)")
-        schema.addField('psfSigma', type=np.float32, doc="PSF sigma (reference CCD)")
+        schema.addField('psfSigma', type=np.float32, doc="PSF sigma (median); pixels")
+        schema.addField('psfFwhm', type=np.float32, doc="PSF FWHM (median); arcseconds")
         schema.addField('deltaAper', type=np.float32, doc="Delta-aperture")
         schema.addField('skyBackground', type=np.float32, doc="Sky background (ADU) (reference CCD)")
         # the following field is not used yet

python/lsst/fgcmcal/fgcmFitCycle.py

@@ -772,6 +772,15 @@ class FgcmFitCycleConfig(pipeBase.PipelineTaskConfig,
         dtype=float,
         default=-0.25,
     )
+    expFwhmCutDict = pexConfig.DictField(
+        doc=("Per-band specification on maximum exposure FWHM (arcseconds) that will "
+             "be considered for the model fit. Exposures with median FWHM larger "
+             "than this threshold will get zeropoints based on matching to good "
+             "stars."),
+        keytype=str,
+        itemtype=float,
+        default={},
+    )
     expGrayPhotometricCutDict = pexConfig.DictField(
         doc=("Per-band specification on maximum (negative) achromatic exposure residual "
              "('gray term') for a visit to be considered photometric.  Must have one "

python/lsst/fgcmcal/utilities.py

@@ -138,7 +138,7 @@ def makeConfigDict(config, log, camera, maxIter,
                   'expField': FGCM_EXP_FIELD,
                   'ccdField': FGCM_CCD_FIELD,
                   'seeingField': 'DELTA_APER',
-                  'fwhmField': 'PSFSIGMA',
+                  'fwhmField': 'PSFFWHM',
                   'skyBrightnessField': 'SKYBACKGROUND',
                   'deepFlag': 'DEEPFLAG',  # unused
                   'bands': list(config.bands),
@@ -183,6 +183,7 @@ def makeConfigDict(config, log, camera, maxIter,
                   'minStarPerExp': config.minStarPerExp,
                   'minExpPerNight': config.minExpPerNight,
                   'expGrayInitialCut': config.expGrayInitialCut,
+                  'expFwhmCutDict': dict(config.expFwhmCutDict),
                   'expGrayPhotometricCutDict': dict(config.expGrayPhotometricCutDict),
                   'expGrayHighCutDict': dict(config.expGrayHighCutDict),
                   'expGrayRecoverCut': config.expGrayRecoverCut,
@@ -404,7 +405,7 @@ def translateVisitCatalog(visitCat):
     fgcmExpInfo = np.zeros(len(visitCat), dtype=[('VISIT', 'i8'),
                                                  ('MJD', 'f8'),
                                                  ('EXPTIME', 'f8'),
-                                                 ('PSFSIGMA', 'f8'),
+                                                 ('PSFFWHM', 'f8'),
                                                  ('DELTA_APER', 'f8'),
                                                  ('SKYBACKGROUND', 'f8'),
                                                  ('DEEPFLAG', 'i2'),
@@ -423,7 +424,7 @@ def translateVisitCatalog(visitCat):
     fgcmExpInfo['TELDEC'][:] = visitCat['teldec']
     fgcmExpInfo['TELROT'][:] = visitCat['telrot']
     fgcmExpInfo['PMB'][:] = visitCat['pmb']
-    fgcmExpInfo['PSFSIGMA'][:] = visitCat['psfSigma']
+    fgcmExpInfo['PSFFWHM'][:] = visitCat['psfFwhm']
     fgcmExpInfo['DELTA_APER'][:] = visitCat['deltaAper']
     fgcmExpInfo['SKYBACKGROUND'][:] = visitCat['skyBackground']
     # Note that we have to go through asAstropy() to get a string