23 """Make a look-up-table (LUT) for FGCM calibration.
25 This task computes a look-up-table for the range in expected atmosphere
26 variation and variation in instrumental throughput (as tracked by the
27 transmission_filter products). By pre-computing linearized integrals,
28 the FGCM fit is orders of magnitude faster for stars with a broad range
29 of colors and observing bands, yielding precision at the 1-2 mmag level.
31 Computing a LUT requires running MODTRAN or with a pre-generated
32 atmosphere table packaged with fgcm.
45 from lsst.utils.timer
import timeMethod
46 from .utilities
import lookupStaticCalibrations
50 __all__ = [
'FgcmMakeLutParametersConfig',
'FgcmMakeLutConfig',
'FgcmMakeLutTask',
55 dimensions=(
'instrument',),
57 camera = connectionTypes.PrerequisiteInput(
58 doc=
"Camera instrument",
60 storageClass=
"Camera",
61 dimensions=(
"instrument",),
62 lookupFunction=lookupStaticCalibrations,
66 transmission_optics = connectionTypes.PrerequisiteInput(
67 doc=
"Optics transmission curve information",
68 name=
"transmission_optics",
69 storageClass=
"TransmissionCurve",
70 dimensions=(
"instrument",),
71 lookupFunction=lookupStaticCalibrations,
76 transmission_sensor = connectionTypes.PrerequisiteInput(
77 doc=
"Sensor transmission curve information",
78 name=
"transmission_sensor",
79 storageClass=
"TransmissionCurve",
80 dimensions=(
"instrument",
"detector",),
81 lookupFunction=lookupStaticCalibrations,
87 transmission_filter = connectionTypes.PrerequisiteInput(
88 doc=
"Filter transmission curve information",
89 name=
"transmission_filter",
90 storageClass=
"TransmissionCurve",
91 dimensions=(
"band",
"instrument",
"physical_filter",),
92 lookupFunction=lookupStaticCalibrations,
98 fgcmLookUpTable = connectionTypes.Output(
99 doc=(
"Atmosphere + instrument look-up-table for FGCM throughput and "
100 "chromatic corrections."),
101 name=
"fgcmLookUpTable",
102 storageClass=
"Catalog",
103 dimensions=(
"instrument",),
108 """Config for parameters if atmosphereTableName not available"""
111 elevation = pexConfig.Field(
112 doc=
"Telescope elevation (m)",
116 pmbRange = pexConfig.ListField(
117 doc=(
"Barometric Pressure range (millibar) "
118 "Recommended range depends on the site."),
122 pmbSteps = pexConfig.Field(
123 doc=
"Barometric Pressure number of steps",
127 pwvRange = pexConfig.ListField(
128 doc=(
"Precipitable Water Vapor range (mm) "
129 "Recommended range depends on the site."),
133 pwvSteps = pexConfig.Field(
134 doc=
"Precipitable Water Vapor number of steps",
138 o3Range = pexConfig.ListField(
139 doc=
"Ozone range (dob)",
141 default=[220.0, 310.0],
143 o3Steps = pexConfig.Field(
144 doc=
"Ozone number of steps",
148 tauRange = pexConfig.ListField(
149 doc=
"Aerosol Optical Depth range (unitless)",
151 default=[0.002, 0.35],
153 tauSteps = pexConfig.Field(
154 doc=
"Aerosol Optical Depth number of steps",
158 alphaRange = pexConfig.ListField(
159 doc=
"Aerosol alpha range (unitless)",
163 alphaSteps = pexConfig.Field(
164 doc=
"Aerosol alpha number of steps",
168 zenithRange = pexConfig.ListField(
169 doc=
"Zenith angle range (degree)",
173 zenithSteps = pexConfig.Field(
174 doc=
"Zenith angle number of steps",
180 pmbStd = pexConfig.Field(
181 doc=(
"Standard Atmosphere pressure (millibar); "
182 "Recommended default depends on the site."),
186 pwvStd = pexConfig.Field(
187 doc=(
"Standard Atmosphere PWV (mm); "
188 "Recommended default depends on the site."),
192 o3Std = pexConfig.Field(
193 doc=
"Standard Atmosphere O3 (dob)",
197 tauStd = pexConfig.Field(
198 doc=
"Standard Atmosphere aerosol optical depth",
202 alphaStd = pexConfig.Field(
203 doc=
"Standard Atmosphere aerosol alpha",
207 airmassStd = pexConfig.Field(
208 doc=(
"Standard Atmosphere airmass; "
209 "Recommended default depends on the survey strategy."),
213 lambdaNorm = pexConfig.Field(
214 doc=
"Aerosol Optical Depth normalization wavelength (Angstrom)",
218 lambdaStep = pexConfig.Field(
219 doc=
"Wavelength step for generating atmospheres (nm)",
223 lambdaRange = pexConfig.ListField(
224 doc=
"Wavelength range for LUT (Angstrom)",
226 default=[3000.0, 11000.0],
231 pipelineConnections=FgcmMakeLutConnections):
232 """Config for FgcmMakeLutTask"""
233 physicalFilters = pexConfig.ListField(
234 doc=
"List of physicalFilter labels to generate look-up table.",
238 stdPhysicalFilterOverrideMap = pexConfig.DictField(
239 doc=(
"Override mapping from physical filter labels to 'standard' physical "
240 "filter labels. The 'standard' physical filter defines the transmission "
241 "curve that the FGCM standard bandpass will be based on. "
242 "Any filter not listed here will be mapped to "
243 "itself (e.g. g->g or HSC-G->HSC-G). Use this override for cross-"
244 "filter calibration such as HSC-R->HSC-R2 and HSC-I->HSC-I2."),
249 atmosphereTableName = pexConfig.Field(
250 doc=
"FGCM name or filename of precomputed atmospheres",
255 parameters = pexConfig.ConfigField(
256 doc=
"Atmosphere parameters (required if no atmosphereTableName)",
257 dtype=FgcmMakeLutParametersConfig,
263 Validate the config parameters.
265 This method behaves differently from the parent validate in the case
266 that atmosphereTableName is set. In this case, the config values
267 for standard values, step sizes, and ranges are loaded
268 directly from the specified atmosphereTableName.
271 self._fields[
'physicalFilters'].
validate(self)
272 self._fields[
'stdPhysicalFilterOverrideMap'].
validate(self)
276 self._fields[
'parameters'].
validate(self)
280 """Subclass of TaskRunner for fgcmMakeLutTask
282 fgcmMakeLutTask.run() takes one argument, the butler, and
283 does not run on any data in the repository.
284 This runner does not use any parallelization.
290 Return a list with one element, the butler.
292 return [parsedCmd.butler]
298 butler: `lsst.daf.persistence.Butler`
302 exitStatus: `list` with `pipeBase.Struct`
303 exitStatus (0: success; 1: failure)
305 task = self.TaskClass(config=self.config, log=self.log)
309 task.runDataRef(butler)
312 task.runDataRef(butler)
313 except Exception
as e:
315 task.log.fatal(
"Failed: %s" % e)
316 if not isinstance(e, pipeBase.TaskError):
317 traceback.print_exc(file=sys.stderr)
319 task.writeMetadata(butler)
322 return [pipeBase.Struct(exitStatus=exitStatus)]
326 Run the task, with no multiprocessing
330 parsedCmd: ArgumentParser parsed command line
335 if self.precall(parsedCmd):
338 resultList = self(targetList[0])
345 Make Look-Up Table for FGCM.
347 This task computes a look-up-table for the range in expected atmosphere
348 variation and variation in instrumental throughput (as tracked by the
349 transmission_filter products). By pre-computing linearized integrals,
350 the FGCM fit is orders of magnitude faster for stars with a broad range
351 of colors and observing bands, yielding precision at the 1-2 mmag level.
353 Computing a LUT requires running MODTRAN or with a pre-generated
354 atmosphere table packaged with fgcm.
357 ConfigClass = FgcmMakeLutConfig
358 RunnerClass = FgcmMakeLutRunner
359 _DefaultName =
"fgcmMakeLut"
361 def __init__(self, butler=None, initInputs=None, **kwargs):
365 def _getMetadataName(self):
371 Make a Look-Up Table for FGCM
375 butler: `lsst.daf.persistence.Butler`
379 ValueError : Raised if configured filter name does not match any of the
380 available filter transmission curves.
382 camera = butler.get(
'camera')
383 opticsDataRef = butler.dataRef(
'transmission_optics')
385 sensorDataRefDict = {}
386 for detector
in camera:
387 sensorDataRefDict[detector.getId()] = butler.dataRef(
'transmission_sensor',
388 dataId={
'ccd': detector.getId()})
390 filterDataRefDict = {}
391 for physicalFilter
in self.config.physicalFilters:
395 dataRef = butler.dataRef(
'transmission_filter', filter=physicalFilter)
396 if not dataRef.datasetExists():
397 raise ValueError(f
"Could not find transmission for filter {physicalFilter}.")
398 filterDataRefDict[physicalFilter] = dataRef
404 butler.put(lutCat,
'fgcmLookUpTable')
407 camera = butlerQC.get(inputRefs.camera)
409 opticsDataRef = butlerQC.get(inputRefs.transmission_optics)
411 sensorRefs = butlerQC.get(inputRefs.transmission_sensor)
412 sensorDataRefDict = {sensorRef.dataId.byName()[
'detector']: sensorRef
for
413 sensorRef
in sensorRefs}
415 filterRefs = butlerQC.get(inputRefs.transmission_filter)
416 filterDataRefDict = {filterRef.dataId[
'physical_filter']: filterRef
for
417 filterRef
in filterRefs}
423 butlerQC.put(lutCat, outputRefs.fgcmLookUpTable)
425 def _fgcmMakeLut(self, camera, opticsDataRef, sensorDataRefDict,
428 Make a FGCM Look-up Table
432 camera : `lsst.afw.cameraGeom.Camera`
433 Camera from the butler.
434 opticsDataRef : `lsst.daf.persistence.ButlerDataRef` or
435 `lsst.daf.butler.DeferredDatasetHandle`
436 Reference to optics transmission curve.
437 sensorDataRefDict : `dict` of [`int`, `lsst.daf.persistence.ButlerDataRef` or
438 `lsst.daf.butler.DeferredDatasetHandle`]
439 Dictionary of references to sensor transmission curves. Key will
441 filterDataRefDict : `dict` of [`str`, `lsst.daf.persistence.ButlerDataRef` or
442 `lsst.daf.butler.DeferredDatasetHandle`]
443 Dictionary of references to filter transmission curves. Key will
444 be physical filter label.
448 fgcmLookUpTable : `BaseCatalog`
449 The FGCM look-up table.
453 self.log.
info(
"Found %d ccds for look-up table" % (nCcd))
464 self.log.
info(
"Making the LUT maker object")
472 throughputLambda = np.arange(self.
fgcmLutMakerfgcmLutMaker.lambdaRange[0],
476 self.log.
info(
"Built throughput lambda, %.1f-%.1f, step %.2f" %
477 (throughputLambda[0], throughputLambda[-1],
478 throughputLambda[1] - throughputLambda[0]))
481 for i, physicalFilter
in enumerate(self.config.physicalFilters):
483 tDict[
'LAMBDA'] = throughputLambda
484 for ccdIndex, detector
in enumerate(camera):
485 tDict[ccdIndex] = self.
_getThroughputDetector_getThroughputDetector(detector, physicalFilter, throughputLambda)
486 throughputDict[physicalFilter] = tDict
489 self.
fgcmLutMakerfgcmLutMaker.setThroughputs(throughputDict)
492 self.log.
info(
"Making LUT")
499 physicalFilterString = comma.join(self.config.physicalFilters)
502 atmosphereTableName =
'NoTableWasUsed'
503 if self.config.atmosphereTableName
is not None:
504 atmosphereTableName = self.config.atmosphereTableName
506 lutSchema = self.
_makeLutSchema_makeLutSchema(physicalFilterString, stdPhysicalFilterString,
509 lutCat = self.
_makeLutCat_makeLutCat(lutSchema, physicalFilterString,
510 stdPhysicalFilterString, atmosphereTableName)
513 def _getStdPhysicalFilterList(self):
514 """Get the standard physical filter lists from config.physicalFilters
515 and config.stdPhysicalFilterOverrideMap
519 stdPhysicalFilters : `list`
521 override = self.config.stdPhysicalFilterOverrideMap
522 return [override.get(physicalFilter, physicalFilter)
for
523 physicalFilter
in self.config.physicalFilters]
525 def _createLutConfig(self, nCcd):
527 Create the fgcmLut config dictionary
532 Number of CCDs in the camera
537 lutConfig[
'logger'] = self.log
538 lutConfig[
'filterNames'] = self.config.physicalFilters
540 lutConfig[
'nCCD'] = nCcd
543 if self.config.atmosphereTableName
is not None:
544 lutConfig[
'atmosphereTableName'] = self.config.atmosphereTableName
547 lutConfig[
'elevation'] = self.config.parameters.elevation
548 lutConfig[
'pmbRange'] = self.config.parameters.pmbRange
549 lutConfig[
'pmbSteps'] = self.config.parameters.pmbSteps
550 lutConfig[
'pwvRange'] = self.config.parameters.pwvRange
551 lutConfig[
'pwvSteps'] = self.config.parameters.pwvSteps
552 lutConfig[
'o3Range'] = self.config.parameters.o3Range
553 lutConfig[
'o3Steps'] = self.config.parameters.o3Steps
554 lutConfig[
'tauRange'] = self.config.parameters.tauRange
555 lutConfig[
'tauSteps'] = self.config.parameters.tauSteps
556 lutConfig[
'alphaRange'] = self.config.parameters.alphaRange
557 lutConfig[
'alphaSteps'] = self.config.parameters.alphaSteps
558 lutConfig[
'zenithRange'] = self.config.parameters.zenithRange
559 lutConfig[
'zenithSteps'] = self.config.parameters.zenithSteps
560 lutConfig[
'pmbStd'] = self.config.parameters.pmbStd
561 lutConfig[
'pwvStd'] = self.config.parameters.pwvStd
562 lutConfig[
'o3Std'] = self.config.parameters.o3Std
563 lutConfig[
'tauStd'] = self.config.parameters.tauStd
564 lutConfig[
'alphaStd'] = self.config.parameters.alphaStd
565 lutConfig[
'airmassStd'] = self.config.parameters.airmassStd
566 lutConfig[
'lambdaRange'] = self.config.parameters.lambdaRange
567 lutConfig[
'lambdaStep'] = self.config.parameters.lambdaStep
568 lutConfig[
'lambdaNorm'] = self.config.parameters.lambdaNorm
572 def _loadThroughputs(self, camera, opticsDataRef, sensorDataRefDict, filterDataRefDict):
573 """Internal method to load throughput data for filters
577 camera: `lsst.afw.cameraGeom.Camera`
578 Camera from the butler
579 opticsDataRef : `lsst.daf.persistence.ButlerDataRef` or
580 `lsst.daf.butler.DeferredDatasetHandle`
581 Reference to optics transmission curve.
582 sensorDataRefDict : `dict` of [`int`, `lsst.daf.persistence.ButlerDataRef` or
583 `lsst.daf.butler.DeferredDatasetHandle`]
584 Dictionary of references to sensor transmission curves. Key will
586 filterDataRefDict : `dict` of [`str`, `lsst.daf.persistence.ButlerDataRef` or
587 `lsst.daf.butler.DeferredDatasetHandle`]
588 Dictionary of references to filter transmission curves. Key will
589 be physical filter label.
593 ValueError : Raised if configured filter name does not match any of the
594 available filter transmission curves.
599 for detector
in camera:
600 self.
_sensorsTransmission_sensorsTransmission[detector.getId()] = sensorDataRefDict[detector.getId()].get()
603 for physicalFilter
in self.config.physicalFilters:
604 self.
_filtersTransmission_filtersTransmission[physicalFilter] = filterDataRefDict[physicalFilter].get()
606 def _getThroughputDetector(self, detector, physicalFilter, throughputLambda):
607 """Internal method to get throughput for a detector.
609 Returns the throughput at the center of the detector for a given filter.
613 detector: `lsst.afw.cameraGeom._detector.Detector`
615 physicalFilter: `str`
616 Physical filter label
617 throughputLambda: `np.array(dtype=np.float64)`
618 Wavelength steps (Angstrom)
622 throughput: `np.array(dtype=np.float64)`
623 Throughput (max 1.0) at throughputLambda
626 c = detector.getCenter(afwCameraGeom.FOCAL_PLANE)
627 c.scale(1.0/detector.getPixelSize()[0])
630 wavelengths=throughputLambda)
632 throughput *= self.
_sensorsTransmission_sensorsTransmission[detector.getId()].sampleAt(position=c,
633 wavelengths=throughputLambda)
635 throughput *= self.
_filtersTransmission_filtersTransmission[physicalFilter].sampleAt(position=c,
636 wavelengths=throughputLambda)
639 throughput = np.clip(throughput, 0.0, 1.0)
643 def _makeLutSchema(self, physicalFilterString, stdPhysicalFilterString,
644 atmosphereTableName):
650 physicalFilterString: `str`
651 Combined string of all the physicalFilters
652 stdPhysicalFilterString: `str`
653 Combined string of all the standard physicalFilters
654 atmosphereTableName: `str`
655 Name of the atmosphere table used to generate LUT
659 lutSchema: `afwTable.schema`
664 lutSchema.addField(
'tablename', type=str, doc=
'Atmosphere table name',
665 size=len(atmosphereTableName))
666 lutSchema.addField(
'elevation', type=float, doc=
"Telescope elevation used for LUT")
667 lutSchema.addField(
'physicalFilters', type=str, doc=
'physicalFilters in LUT',
668 size=len(physicalFilterString))
669 lutSchema.addField(
'stdPhysicalFilters', type=str, doc=
'Standard physicalFilters in LUT',
670 size=len(stdPhysicalFilterString))
671 lutSchema.addField(
'pmb', type=
'ArrayD', doc=
'Barometric Pressure',
673 lutSchema.addField(
'pmbFactor', type=
'ArrayD', doc=
'PMB scaling factor',
675 lutSchema.addField(
'pmbElevation', type=np.float64, doc=
'PMB Scaling at elevation')
676 lutSchema.addField(
'pwv', type=
'ArrayD', doc=
'Preciptable Water Vapor',
678 lutSchema.addField(
'o3', type=
'ArrayD', doc=
'Ozone',
680 lutSchema.addField(
'tau', type=
'ArrayD', doc=
'Aerosol optical depth',
682 lutSchema.addField(
'lambdaNorm', type=np.float64, doc=
'AOD wavelength')
683 lutSchema.addField(
'alpha', type=
'ArrayD', doc=
'Aerosol alpha',
685 lutSchema.addField(
'zenith', type=
'ArrayD', doc=
'Zenith angle',
687 lutSchema.addField(
'nCcd', type=np.int32, doc=
'Number of CCDs')
690 lutSchema.addField(
'pmbStd', type=np.float64, doc=
'PMB Standard')
691 lutSchema.addField(
'pwvStd', type=np.float64, doc=
'PWV Standard')
692 lutSchema.addField(
'o3Std', type=np.float64, doc=
'O3 Standard')
693 lutSchema.addField(
'tauStd', type=np.float64, doc=
'Tau Standard')
694 lutSchema.addField(
'alphaStd', type=np.float64, doc=
'Alpha Standard')
695 lutSchema.addField(
'zenithStd', type=np.float64, doc=
'Zenith angle Standard')
696 lutSchema.addField(
'lambdaRange', type=
'ArrayD', doc=
'Wavelength range',
698 lutSchema.addField(
'lambdaStep', type=np.float64, doc=
'Wavelength step')
699 lutSchema.addField(
'lambdaStd', type=
'ArrayD', doc=
'Standard Wavelength',
701 lutSchema.addField(
'lambdaStdFilter', type=
'ArrayD', doc=
'Standard Wavelength (raw)',
703 lutSchema.addField(
'i0Std', type=
'ArrayD', doc=
'I0 Standard',
705 lutSchema.addField(
'i1Std', type=
'ArrayD', doc=
'I1 Standard',
707 lutSchema.addField(
'i10Std', type=
'ArrayD', doc=
'I10 Standard',
709 lutSchema.addField(
'i2Std', type=
'ArrayD', doc=
'I2 Standard',
711 lutSchema.addField(
'lambdaB', type=
'ArrayD', doc=
'Wavelength for passband (no atm)',
713 lutSchema.addField(
'atmLambda', type=
'ArrayD', doc=
'Atmosphere wavelengths (Angstrom)',
715 lutSchema.addField(
'atmStdTrans', type=
'ArrayD', doc=
'Standard Atmosphere Throughput',
719 lutSchema.addField(
'luttype', type=str, size=20, doc=
'Look-up table type')
720 lutSchema.addField(
'lut', type=
'ArrayF', doc=
'Look-up table for luttype',
725 def _makeLutCat(self, lutSchema, physicalFilterString, stdPhysicalFilterString,
726 atmosphereTableName):
732 lutSchema: `afwTable.schema`
734 physicalFilterString: `str`
735 Combined string of all the physicalFilters
736 stdPhysicalFilterString: `str`
737 Combined string of all the standard physicalFilters
738 atmosphereTableName: `str`
739 Name of the atmosphere table used to generate LUT
743 lutCat: `afwTable.BaseCatalog`
744 Lut catalog for persistence
752 lutCat.table.preallocate(14)
755 rec = lutCat.addNew()
757 rec[
'tablename'] = atmosphereTableName
758 rec[
'elevation'] = self.
fgcmLutMakerfgcmLutMaker.atmosphereTable.elevation
759 rec[
'physicalFilters'] = physicalFilterString
760 rec[
'stdPhysicalFilters'] = stdPhysicalFilterString
762 rec[
'pmbFactor'][:] = self.
fgcmLutMakerfgcmLutMaker.pmbFactor
763 rec[
'pmbElevation'] = self.
fgcmLutMakerfgcmLutMaker.pmbElevation
767 rec[
'lambdaNorm'] = self.
fgcmLutMakerfgcmLutMaker.lambdaNorm
776 rec[
'alphaStd'] = self.
fgcmLutMakerfgcmLutMaker.alphaStd
777 rec[
'zenithStd'] = self.
fgcmLutMakerfgcmLutMaker.zenithStd
778 rec[
'lambdaRange'][:] = self.
fgcmLutMakerfgcmLutMaker.lambdaRange
779 rec[
'lambdaStep'] = self.
fgcmLutMakerfgcmLutMaker.lambdaStep
780 rec[
'lambdaStd'][:] = self.
fgcmLutMakerfgcmLutMaker.lambdaStd
781 rec[
'lambdaStdFilter'][:] = self.
fgcmLutMakerfgcmLutMaker.lambdaStdFilter
786 rec[
'lambdaB'][:] = self.
fgcmLutMakerfgcmLutMaker.lambdaB
787 rec[
'atmLambda'][:] = self.
fgcmLutMakerfgcmLutMaker.atmLambda
788 rec[
'atmStdTrans'][:] = self.
fgcmLutMakerfgcmLutMaker.atmStdTrans
790 rec[
'luttype'] =
'I0'
791 rec[
'lut'][:] = self.
fgcmLutMakerfgcmLutMaker.lut[
'I0'].flatten()
794 rec = lutCat.addNew()
795 rec[
'luttype'] =
'I1'
796 rec[
'lut'][:] = self.
fgcmLutMakerfgcmLutMaker.lut[
'I1'].flatten()
798 derivTypes = [
'D_PMB',
'D_LNPWV',
'D_O3',
'D_LNTAU',
'D_ALPHA',
'D_SECZENITH',
799 'D_PMB_I1',
'D_LNPWV_I1',
'D_O3_I1',
'D_LNTAU_I1',
'D_ALPHA_I1',
801 for derivType
in derivTypes:
802 rec = lutCat.addNew()
803 rec[
'luttype'] = derivType
804 rec[
'lut'][:] = self.
fgcmLutMakerfgcmLutMaker.lutDeriv[derivType].flatten()
Defines the fields and offsets for a table.
def getTargetList(parsedCmd)
def __call__(self, butler)
def _makeLutCat(self, lutSchema, physicalFilterString, stdPhysicalFilterString, atmosphereTableName)
def _getStdPhysicalFilterList(self)
def _createLutConfig(self, nCcd)
def runDataRef(self, butler)
def _fgcmMakeLut(self, camera, opticsDataRef, sensorDataRefDict, filterDataRefDict)
def _makeLutSchema(self, physicalFilterString, stdPhysicalFilterString, atmosphereTableName)
def _getThroughputDetector(self, detector, physicalFilter, throughputLambda)
def runQuantum(self, butlerQC, inputRefs, outputRefs)
def _loadThroughputs(self, camera, opticsDataRef, sensorDataRefDict, filterDataRefDict)
def __init__(self, butler=None, initInputs=None, **kwargs)