DM-39592: Add support for Auxtel/LATISS and add associated tests.

.gitignore

@@ -17,6 +17,7 @@ doc/*.tag
 doc/*.inc
 doc/doxygen.conf
 tests/.tests
+tests/TestFgcm*
 version.py
 bin/
 pytest_session.txt

python/lsst/fgcmcal/fgcmFitCycle.py

@@ -442,6 +442,13 @@ class FgcmFitCycleConfig(pipeBase.PipelineTaskConfig,
         dtype=float,
         default=5.0,
     )
+    superStarPlotCcdResiduals = pexConfig.Field(
+        doc="If plotting is enabled, should per-detector residuals be plotted? "
+            "This may produce a lot of output, and should be used only for "
+            "debugging purposes.",
+        dtype=bool,
+        default=False,
+    )
     focalPlaneSigmaClip = pexConfig.Field(
         doc="Number of sigma to clip outliers per focal-plane.",
         dtype=float,
@@ -542,6 +549,13 @@ class FgcmFitCycleConfig(pipeBase.PipelineTaskConfig,
         dtype=float,
         default=None,
     )
+    mirrorArea = pexConfig.Field(
+        doc="Mirror area (square meters) of telescope.  If not set, will "
+            "try to estimate from camera.telescopeDiameter.",
+        dtype=float,
+        default=None,
+        optional=True,
+    )
     defaultCameraOrientation = pexConfig.Field(
         doc="Default camera orientation for QA plots.",
         dtype=float,

python/lsst/fgcmcal/fgcmMakeLut.py

@@ -112,6 +112,12 @@ class FgcmMakeLutConnections(pipeBase.PipelineTaskConnections,
         multiple=True,
     )
 
+    def __init__(self, *, config=None):
+        if not config.doOpticsTransmission:
+            del self.transmission_optics
+        if not config.doSensorTransmission:
+            del self.transmission_sensor
+
 
 class FgcmMakeLutParametersConfig(pexConfig.Config):
     """Config for parameters if atmosphereTableName not available"""
@@ -261,6 +267,16 @@ class FgcmMakeLutConfig(pipeBase.PipelineTaskConfig,
         default=None,
         optional=True,
     )
+    doOpticsTransmission = pexConfig.Field(
+        doc="Include optics transmission?",
+        dtype=bool,
+        default=True,
+    )
+    doSensorTransmission = pexConfig.Field(
+        doc="Include sensor transmission?",
+        dtype=bool,
+        default=True,
+    )
     parameters = pexConfig.ConfigField(
         doc="Atmosphere parameters (required if no atmosphereTableName)",
         dtype=FgcmMakeLutParametersConfig,
@@ -308,11 +324,17 @@ def __init__(self, initInputs=None, **kwargs):
     def runQuantum(self, butlerQC, inputRefs, outputRefs):
         camera = butlerQC.get(inputRefs.camera)
 
-        opticsHandle = butlerQC.get(inputRefs.transmission_optics)
+        if self.config.doOpticsTransmission:
+            opticsHandle = butlerQC.get(inputRefs.transmission_optics)
+        else:
+            opticsHandle = None
 
-        sensorHandles = butlerQC.get(inputRefs.transmission_sensor)
-        sensorHandleDict = {sensorHandle.dataId.byName()['detector']: sensorHandle for
-                            sensorHandle in sensorHandles}
+        if self.config.doSensorTransmission:
+            sensorHandles = butlerQC.get(inputRefs.transmission_sensor)
+            sensorHandleDict = {sensorHandle.dataId.byName()['detector']: sensorHandle for
+                                sensorHandle in sensorHandles}
+        else:
+            sensorHandleDict = {}
 
         filterHandles = butlerQC.get(inputRefs.transmission_filter)
         filterHandleDict = {filterHandle.dataId['physical_filter']: filterHandle for
@@ -528,11 +550,35 @@ def _loadThroughputs(self, camera, opticsHandle, sensorHandleDict, filterHandleD
         ValueError : Raised if configured filter name does not match any of the
             available filter transmission curves.
         """
-        self._opticsTransmission = opticsHandle.get()
+        if self.config.doOpticsTransmission:
+            self._opticsTransmission = opticsHandle.get()
+        else:
+            self._opticsTransmission = TransmissionCurve.makeSpatiallyConstant(
+                throughput=np.ones(100),
+                wavelengths=np.linspace(
+                    self.config.parameters.lambdaRange[0],
+                    self.config.parameters.lambdaRange[1],
+                    100,
+                ),
+                throughputAtMin=1.0,
+                throughputAtMax=1.0,
+            )
 
         self._sensorsTransmission = {}
         for detector in camera:
-            self._sensorsTransmission[detector.getId()] = sensorHandleDict[detector.getId()].get()
+            if self.config.doSensorTransmission:
+                self._sensorsTransmission[detector.getId()] = sensorHandleDict[detector.getId()].get()
+            else:
+                self._sensorsTransmission[detector.getId()] = TransmissionCurve.makeSpatiallyConstant(
+                    throughput=np.ones(100),
+                    wavelengths=np.linspace(
+                        self.config.parameters.lambdaRange[0],
+                        self.config.parameters.lambdaRange[1],
+                        100,
+                    ),
+                    throughputAtMin=1.0,
+                    throughputAtMax=1.0,
+                )
 
         self._filtersTransmission = {}
         for physicalFilter in self.config.physicalFilters:

python/lsst/fgcmcal/utilities.py

@@ -99,7 +99,11 @@ def makeConfigDict(config, log, camera, maxIter,
     # TODO: Having direct access to the mirror area from the camera would be
     #  useful.  See DM-16489.
     # Mirror area in cm**2
-    mirrorArea = np.pi*(camera.telescopeDiameter*100./2.)**2.
+    if config.mirrorArea is None:
+        mirrorArea = np.pi*(camera.telescopeDiameter*100./2.)**2.
+    else:
+        # Convert to square cm.
+        mirrorArea = config.mirrorArea * 100.**2.
 
     # Get approximate average camera gain:
     gains = [amp.getGain() for detector in camera for amp in detector.getAmplifiers()]
@@ -146,6 +150,7 @@ def makeConfigDict(config, log, camera, maxIter,
                   'superStarSubCCDChebyshevOrder': config.superStarSubCcdChebyshevOrder,
                   'superStarSubCCDTriangular': config.superStarSubCcdTriangular,
                   'superStarSigmaClip': config.superStarSigmaClip,
+                  'superStarPlotCCDResiduals': config.superStarPlotCcdResiduals,
                   'focalPlaneSigmaClip': config.focalPlaneSigmaClip,
                   'ccdGraySubCCDDict': dict(config.ccdGraySubCcdDict),
                   'ccdGraySubCCDChebyshevOrder': config.ccdGraySubCcdChebyshevOrder,
@@ -573,7 +578,7 @@ def makeZptSchema(superStarChebyshevSize, zptChebyshevSize):
                             '(value set by deltaMagBkgOffsetPercentile) calibration '
                             'stars.'))
     zptSchema.addField('exptime', type=np.float32, doc='Exposure time')
-    zptSchema.addField('filtername', type=str, size=10, doc='Filter name')
+    zptSchema.addField('filtername', type=str, size=30, doc='Filter name')
 
     return zptSchema
 
@@ -644,7 +649,7 @@ def makeAtmSchema():
 
     atmSchema = afwTable.Schema()
 
-    atmSchema.addField('visit', type=np.int32, doc='Visit number')
+    atmSchema.addField('visit', type=np.int64, doc='Visit number')
     atmSchema.addField('pmb', type=np.float64, doc='Barometric pressure (mb)')
     atmSchema.addField('pwv', type=np.float64, doc='Water vapor (mm)')
     atmSchema.addField('tau', type=np.float64, doc='Aerosol optical depth')

setup.cfg

@@ -5,6 +5,7 @@ exclude =
 	doc/conf.py,
 	**/*/__init__.py,
 	tests/config/*,
-	cookbook
+	cookbook,
+	tests/TestFgcm*
 
 [tool:pytest]

tests/config/colortermsLatiss.py

@@ -0,0 +1,37 @@
+from lsst.pipe.tasks.colorterms import Colorterm, ColortermDict
+
+
+config.data = {
+    "atlas_refcat2*": ColortermDict(data={
+        "SDSSg_65mm~empty": Colorterm(
+            primary="g",
+            secondary="r",
+            c0=-0.09034144345111599,
+            c1=0.1710923238086337,
+            c2=-0.038260355621929296,
+        ),
+        "SDSSr_65mm~empty": Colorterm(
+            primary="r",
+            secondary="i",
+            c0=0.0073632488906825045,
+            c1=-0.026620900037027242,
+            c2=-0.03203533692013322,
+        ),
+        "SDSSi_65mm~empty": Colorterm(
+            primary="i",
+            secondary="r",
+            c0=0.016940180565664747,
+            c1=0.0610018330811135,
+            c2=-0.0722575356707918,
+        ),
+        # The following two are blank until we have data to measure them.
+        "SDSSz_65mm~empty": Colorterm(
+            primary="z",
+            secondary="z",
+        ),
+        "SDSSy_65mm~empty": Colorterm(
+            primary="y",
+            secondary="y",
+        ),
+    }),
+}

tests/config/fgcmBuildFromIsolatedStarsHsc.py

@@ -25,9 +25,9 @@
 config.physicalFilterMap = HSC_FILTER_DEFINITIONS.physical_to_band
 config.doSubtractLocalBackground = True
 config.sourceSelector["science"].flags.bad.append("localBackground_flag")
-config.fgcmLoadReferenceCatalog.load(os.path.join(configDir, 'filterMap.py'))
+config.fgcmLoadReferenceCatalog.load(os.path.join(configDir, 'filterMapHsc.py'))
 config.fgcmLoadReferenceCatalog.applyColorTerms = True
-config.fgcmLoadReferenceCatalog.colorterms.load(os.path.join(configDir, 'colorterms.py'))
+config.fgcmLoadReferenceCatalog.colorterms.load(os.path.join(configDir, 'colortermsHsc.py'))
 config.fgcmLoadReferenceCatalog.referenceSelector.doSignalToNoise = True
 config.fgcmLoadReferenceCatalog.referenceSelector.signalToNoise.fluxField = 'i_flux'
 config.fgcmLoadReferenceCatalog.referenceSelector.signalToNoise.errField = 'i_fluxErr'

tests/config/fgcmBuildFromIsolatedStarsLatiss.py

@@ -0,0 +1,35 @@
+import os
+from lsst.obs.lsst.filters import LATISS_FILTER_DEFINITIONS
+
+# The filterMap and bands are for the small subset of bands used in the tests
+config.physicalFilterMap = {
+    "SDSSg_65mm~empty": "g",
+    "SDSSr_65mm~empty": "r",
+    "SDSSi_65mm~empty": "i",
+}
+config.requiredBands = ["g", "r"]
+config.primaryBands = ["r"]
+# The coarseNside is set appropriate to the area of the test data
+config.coarseNside = 64
+# The tests are done with only the brightest reference stars
+config.fgcmLoadReferenceCatalog.referenceSelector.signalToNoise.minimum = 50.0
+
+config.instFluxField = "apFlux_35_0_instFlux"
+config.apertureInnerInstFluxField = "apFlux_35_0_instFlux"
+config.apertureOuterInstFluxField = "apFlux_50_0_instFlux"
+
+config.minPerBand = 2
+config.connections.ref_cat = "atlas_refcat2_20220201"
+
+configDir = os.path.join(os.path.dirname(__file__))
+config.physicalFilterMap = LATISS_FILTER_DEFINITIONS.physical_to_band
+config.doSubtractLocalBackground = True
+config.sourceSelector["science"].flags.bad.append("localBackground_flag")
+config.sourceSelector["science"].signalToNoise.fluxField = "apFlux_35_0_instFlux"
+config.sourceSelector["science"].signalToNoise.errField = "apFlux_35_0_instFluxErr"
+config.fgcmLoadReferenceCatalog.load(os.path.join(configDir, "filterMapLatiss.py"))
+config.fgcmLoadReferenceCatalog.applyColorTerms = True
+config.fgcmLoadReferenceCatalog.colorterms.load(os.path.join(configDir, "colortermsLatiss.py"))
+config.fgcmLoadReferenceCatalog.referenceSelector.doSignalToNoise = True
+config.fgcmLoadReferenceCatalog.referenceSelector.signalToNoise.fluxField = "i_flux"
+config.fgcmLoadReferenceCatalog.referenceSelector.signalToNoise.errField = "i_fluxErr"

tests/config/fgcmBuildStarsTableHsc.py

@@ -29,9 +29,9 @@
 config.sourceSelector["science"].flags.bad.append("localBackground_flag")
 # Our test files do not have detect_isPrimary in the columns.
 config.sourceSelector["science"].doRequirePrimary = False
-config.fgcmLoadReferenceCatalog.load(os.path.join(configDir, 'filterMap.py'))
+config.fgcmLoadReferenceCatalog.load(os.path.join(configDir, 'filterMapHsc.py'))
 config.fgcmLoadReferenceCatalog.applyColorTerms = True
-config.fgcmLoadReferenceCatalog.colorterms.load(os.path.join(configDir, 'colorterms.py'))
+config.fgcmLoadReferenceCatalog.colorterms.load(os.path.join(configDir, 'colortermsHsc.py'))
 config.fgcmLoadReferenceCatalog.referenceSelector.doSignalToNoise = True
 config.fgcmLoadReferenceCatalog.referenceSelector.signalToNoise.fluxField = 'i_flux'
 config.fgcmLoadReferenceCatalog.referenceSelector.signalToNoise.errField = 'i_fluxErr'

tests/config/fgcmFitCycleLatiss.py

@@ -0,0 +1,108 @@
+import lsst.fgcmcal as fgcmcal
+
+config.outfileBase = "TestFgcm"
+# Use these bands and fit them.
+# i band does not have any observations in tests, this ensures
+# that the code performs properly when there are missing bands.
+config.bands = ["g", "r", "i"]
+config.fitBands = ["g", "r", "i"]
+from lsst.obs.lsst.filters import LATISS_FILTER_DEFINITIONS
+config.physicalFilterMap = LATISS_FILTER_DEFINITIONS.physical_to_band
+# Only require g, r observations for a star to be a calibration star.
+config.requiredBands = ["g", "r"]
+# Do 5 iterations in multi-cycle run mode.
+config.maxIterBeforeFinalCycle = 5
+config.nCore = 1
+config.cycleNumber = 0
+config.utBoundary = 0.0
+config.washMjds = (0.0, )
+# For tests, define 1 observing epoch that encompasses everything.
+config.epochMjds = (0.0, 100000.0)
+config.coatingMjds = []
+config.latitude = -30.2333
+# This is pi*(1.2/2.)**2.
+config.mirrorArea = 1.13097
+config.defaultCameraOrientation = 0.0
+config.brightObsGrayMax = 0.5
+config.expGrayInitialCut = -0.5
+config.expGrayPhotometricCutDict = {"g": -0.5, "r": -0.5, "i": -0.5}
+config.expGrayHighCutDict = {"g": 0.2, "r": 0.2, "i": 0.2}
+config.expVarGrayPhotometricCutDict = {"g": 0.1**2.,
+                                       "r": 0.1**2.,
+                                       "i": 0.1**2.}
+# For tests, make a broad cut for outliers.
+config.autoPhotometricCutNSig = 5.0
+config.autoHighCutNSig = 5.0
+# Fit aperture corrections with only 2 bins to exercise the code.
+config.aperCorrFitNBins = 2
+config.aperCorrInputSlopeDict = {"g": -1.0,
+                                 "r": -1.0,
+                                 "i": -1.0}
+# Define the band to SED constants approximately so they work
+# for data that only has r, r observations.
+config.sedboundaryterms = fgcmcal.SedboundarytermDict()
+config.sedboundaryterms.data["gr"] = fgcmcal.Sedboundaryterm(primary="g",
+                                                             secondary="r")
+config.sedterms = fgcmcal.SedtermDict()
+config.sedterms.data["g"] = fgcmcal.Sedterm(primaryTerm="gr", secondaryTerm=None,
+                                            extrapolated=False, constant=0.0)
+config.sedterms.data["r"] = fgcmcal.Sedterm(primaryTerm="gr", secondaryTerm=None,
+                                            extrapolated=False, constant=1.0)
+config.sedterms.data["i"] = fgcmcal.Sedterm(primaryTerm="gr", secondaryTerm=None,
+                                            extrapolated=False, constant=0.75)
+# Define good stars with an g-r color cut.
+config.starColorCuts = ("g, r, -0.50, 2.25",)
+config.refStarColorCuts = ("g, r, -0.50, 2.25",)
+config.useExposureReferenceOffset = True
+config.precomputeSuperStarInitialCycle = False
+config.superStarSubCcdDict = {"g": True,
+                              "r": True,
+                              "i": True}
+config.superStarPlotCcdResiduals = True
+# Allow calibration to work with just 1 exposure on a night.
+config.minExpPerNight = 1
+# Allow calibration to work with very few stars per exposure.
+config.minStarPerExp = 5
+# Allow calibration to work with small number of stars in processing batches.
+config.nStarPerRun = 50
+config.nExpPerRun = 2
+# Define g-r color as the primary way to split by color.
+config.colorSplitBands = ["g", "r"]
+config.freezeStdAtmosphere = True
+# For tests, do low-order per-ccd polynomial.
+config.superStarSubCcdChebyshevOrder = 1
+config.ccdGraySubCcdDict = {"g": False,
+                            "r": False,
+                            "i": False}
+config.ccdGrayFocalPlaneDict = {"g": False,
+                                "r": False,
+                                "i": False}
+config.ccdGrayFocalPlaneFitMinCcd = 1
+config.ccdGrayFocalPlaneChebyshevOrder = 1
+# Do not model the magnitude errors (use errors as reported).
+config.modelMagErrors = False
+# Fix the sigma_cal calibration noise to 0.003 mag.
+config.sigmaCalRange = (0.003, 0.003)
+# Do not fit instrumental parameters (mirror decay) per band.
+config.instrumentParsPerBand = False
+# Set the random seed for repeatability in fits.
+config.randomSeed = 12345
+# Do not use star repeatability metrics for selecting exposures.
+# (Instead, use exposure repeatability metrics).
+config.useRepeatabilityForExpGrayCutsDict = {"g": False,
+                                             "r": False,
+                                             "i": False}
+config.sigFgcmMaxEGrayDict = {"g": 0.1,
+                              "r": 0.1,
+                              "i": 0.1}
+config.approxThroughputDict = {"g": 1.0,
+                              "r": 1.0,
+                              "i": 1.0}
+
+config.deltaAperFitPerCcdNx = 2
+config.deltaAperFitPerCcdNy = 2
+config.deltaAperInnerRadiusArcsec = 2.04
+config.deltaAperOuterRadiusArcsec = 2.89
+config.doComputeDeltaAperPerVisit = True
+config.doComputeDeltaAperMap = True
+config.doComputeDeltaAperPerCcd = True

tests/config/fgcmMakeLutLatiss.py

@@ -0,0 +1,2 @@
+config.physicalFilters = ['SDSSg_65mm~empty', 'SDSSr_65mm~empty', 'SDSSi_65mm~empty']
+config.atmosphereTableName = 'fgcm_atm_lsst2_test'

tests/config/fgcmOutputProductsHsc.py

@@ -15,5 +15,5 @@
 config.photoCal.photoCatName = 'ps1_pv3_3pi_20170110'
 
 configDir = os.path.join(os.path.dirname(__file__))
-config.photoCal.colorterms.load(os.path.join(configDir, 'colorterms.py'))
+config.photoCal.colorterms.load(os.path.join(configDir, 'colortermsHsc.py'))
 config.connections.refCat = 'ps1_pv3_3pi_20170110'