23"""Make a look-up-table (LUT) for FGCM calibration.
25This task computes a look-up-table for the range in expected atmosphere
26variation and variation in instrumental throughput (as tracked by the
27transmission_filter products). By pre-computing linearized integrals,
28the FGCM fit is orders of magnitude faster
for stars
with a broad range
29of colors
and observing bands, yielding precision at the 1-2 mmag level.
31Computing a LUT requires running MODTRAN
or with a pre-generated
32atmosphere table packaged
with fgcm.
37import lsst.pex.config as pexConfig
38import lsst.pipe.base as pipeBase
39from lsst.pipe.base import connectionTypes
40import lsst.afw.table as afwTable
41import lsst.afw.cameraGeom as afwCameraGeom
42from .utilities import lookupStaticCalibrations
46__all__ = ['FgcmMakeLutParametersConfig', 'FgcmMakeLutConfig', 'FgcmMakeLutTask']
49class FgcmMakeLutConnections(pipeBase.PipelineTaskConnections,
50 dimensions=('instrument',),
52 camera = connectionTypes.PrerequisiteInput(
53 doc=
"Camera instrument",
55 storageClass=
"Camera",
56 dimensions=(
"instrument",),
57 lookupFunction=lookupStaticCalibrations,
61 transmission_optics = connectionTypes.PrerequisiteInput(
62 doc=
"Optics transmission curve information",
63 name=
"transmission_optics",
64 storageClass=
"TransmissionCurve",
65 dimensions=(
"instrument",),
66 lookupFunction=lookupStaticCalibrations,
71 transmission_sensor = connectionTypes.PrerequisiteInput(
72 doc=
"Sensor transmission curve information",
73 name=
"transmission_sensor",
74 storageClass=
"TransmissionCurve",
75 dimensions=(
"instrument",
"detector",),
76 lookupFunction=lookupStaticCalibrations,
82 transmission_filter = connectionTypes.PrerequisiteInput(
83 doc=
"Filter transmission curve information",
84 name=
"transmission_filter",
85 storageClass=
"TransmissionCurve",
86 dimensions=(
"band",
"instrument",
"physical_filter",),
87 lookupFunction=lookupStaticCalibrations,
93 fgcmLookUpTable = connectionTypes.Output(
94 doc=(
"Atmosphere + instrument look-up-table for FGCM throughput and "
95 "chromatic corrections."),
96 name=
"fgcmLookUpTable",
97 storageClass=
"Catalog",
98 dimensions=(
"instrument",),
103 """Config for parameters if atmosphereTableName not available"""
106 elevation = pexConfig.Field(
107 doc=
"Telescope elevation (m)",
111 pmbRange = pexConfig.ListField(
112 doc=(
"Barometric Pressure range (millibar) "
113 "Recommended range depends on the site."),
117 pmbSteps = pexConfig.Field(
118 doc=
"Barometric Pressure number of steps",
122 pwvRange = pexConfig.ListField(
123 doc=(
"Precipitable Water Vapor range (mm) "
124 "Recommended range depends on the site."),
128 pwvSteps = pexConfig.Field(
129 doc=
"Precipitable Water Vapor number of steps",
133 o3Range = pexConfig.ListField(
134 doc=
"Ozone range (dob)",
136 default=[220.0, 310.0],
138 o3Steps = pexConfig.Field(
139 doc=
"Ozone number of steps",
143 tauRange = pexConfig.ListField(
144 doc=
"Aerosol Optical Depth range (unitless)",
146 default=[0.002, 0.35],
148 tauSteps = pexConfig.Field(
149 doc=
"Aerosol Optical Depth number of steps",
153 alphaRange = pexConfig.ListField(
154 doc=
"Aerosol alpha range (unitless)",
158 alphaSteps = pexConfig.Field(
159 doc=
"Aerosol alpha number of steps",
163 zenithRange = pexConfig.ListField(
164 doc=
"Zenith angle range (degree)",
168 zenithSteps = pexConfig.Field(
169 doc=
"Zenith angle number of steps",
175 pmbStd = pexConfig.Field(
176 doc=(
"Standard Atmosphere pressure (millibar); "
177 "Recommended default depends on the site."),
181 pwvStd = pexConfig.Field(
182 doc=(
"Standard Atmosphere PWV (mm); "
183 "Recommended default depends on the site."),
187 o3Std = pexConfig.Field(
188 doc=
"Standard Atmosphere O3 (dob)",
192 tauStd = pexConfig.Field(
193 doc=
"Standard Atmosphere aerosol optical depth",
197 alphaStd = pexConfig.Field(
198 doc=
"Standard Atmosphere aerosol alpha",
202 airmassStd = pexConfig.Field(
203 doc=(
"Standard Atmosphere airmass; "
204 "Recommended default depends on the survey strategy."),
208 lambdaNorm = pexConfig.Field(
209 doc=
"Aerosol Optical Depth normalization wavelength (Angstrom)",
213 lambdaStep = pexConfig.Field(
214 doc=
"Wavelength step for generating atmospheres (nm)",
218 lambdaRange = pexConfig.ListField(
219 doc=
"Wavelength range for LUT (Angstrom)",
221 default=[3000.0, 11000.0],
226 pipelineConnections=FgcmMakeLutConnections):
227 """Config for FgcmMakeLutTask"""
228 physicalFilters = pexConfig.ListField(
229 doc=
"List of physicalFilter labels to generate look-up table.",
233 stdPhysicalFilterOverrideMap = pexConfig.DictField(
234 doc=(
"Override mapping from physical filter labels to 'standard' physical "
235 "filter labels. The 'standard' physical filter defines the transmission "
236 "curve that the FGCM standard bandpass will be based on. "
237 "Any filter not listed here will be mapped to "
238 "itself (e.g. g->g or HSC-G->HSC-G). Use this override for cross-"
239 "filter calibration such as HSC-R->HSC-R2 and HSC-I->HSC-I2."),
244 atmosphereTableName = pexConfig.Field(
245 doc=
"FGCM name or filename of precomputed atmospheres",
250 parameters = pexConfig.ConfigField(
251 doc=
"Atmosphere parameters (required if no atmosphereTableName)",
252 dtype=FgcmMakeLutParametersConfig,
258 Validate the config parameters.
260 This method behaves differently from the parent validate
in the case
261 that atmosphereTableName
is set. In this case, the config values
262 for standard values, step sizes,
and ranges are loaded
263 directly
from the specified atmosphereTableName.
266 self._fields[
'physicalFilters'].
validate(self)
267 self._fields[
'stdPhysicalFilterOverrideMap'].
validate(self)
271 self._fields[
'parameters'].
validate(self)
276 Make Look-Up Table for FGCM.
278 This task computes a look-up-table
for the range
in expected atmosphere
279 variation
and variation
in instrumental throughput (
as tracked by the
280 transmission_filter products). By pre-computing linearized integrals,
281 the FGCM fit
is orders of magnitude faster
for stars
with a broad range
282 of colors
and observing bands, yielding precision at the 1-2 mmag level.
284 Computing a LUT requires running MODTRAN
or with a pre-generated
285 atmosphere table packaged
with fgcm.
288 ConfigClass = FgcmMakeLutConfig
289 _DefaultName = "fgcmMakeLut"
291 def __init__(self, butler=None, initInputs=None, **kwargs):
295 camera = butlerQC.get(inputRefs.camera)
297 opticsHandle = butlerQC.get(inputRefs.transmission_optics)
299 sensorHandles = butlerQC.get(inputRefs.transmission_sensor)
300 sensorHandleDict = {sensorHandle.dataId.byName()[
'detector']: sensorHandle
for
301 sensorHandle
in sensorHandles}
303 filterHandles = butlerQC.get(inputRefs.transmission_filter)
304 filterHandleDict = {filterHandle.dataId[
'physical_filter']: filterHandle
for
305 filterHandle
in filterHandles}
311 butlerQC.put(lutCat, outputRefs.fgcmLookUpTable)
313 def _fgcmMakeLut(self, camera, opticsHandle, sensorHandleDict,
316 Make a FGCM Look-up Table
321 Camera from the butler.
322 opticsHandle : `lsst.daf.butler.DeferredDatasetHandle`
323 Reference to optics transmission curve.
324 sensorHandleDict : `dict` of [`int`, `lsst.daf.butler.DeferredDatasetHandle`]
325 Dictionary of references to sensor transmission curves. Key will
327 filterHandleDict : `dict` of [`str`, `lsst.daf.butler.DeferredDatasetHandle`]
328 Dictionary of references to filter transmission curves. Key will
329 be physical filter label.
333 fgcmLookUpTable : `BaseCatalog`
334 The FGCM look-up table.
338 self.log.
info(
"Found %d ccds for look-up table" % (nCcd))
349 self.log.
info(
"Making the LUT maker object")
357 throughputLambda = np.arange(self.
fgcmLutMakerfgcmLutMaker.lambdaRange[0],
361 self.log.
info(
"Built throughput lambda, %.1f-%.1f, step %.2f" %
362 (throughputLambda[0], throughputLambda[-1],
363 throughputLambda[1] - throughputLambda[0]))
366 for i, physicalFilter
in enumerate(self.config.physicalFilters):
368 tDict[
'LAMBDA'] = throughputLambda
369 for ccdIndex, detector
in enumerate(camera):
370 tDict[ccdIndex] = self.
_getThroughputDetector_getThroughputDetector(detector, physicalFilter, throughputLambda)
371 throughputDict[physicalFilter] = tDict
374 self.
fgcmLutMakerfgcmLutMaker.setThroughputs(throughputDict)
377 self.log.
info(
"Making LUT")
384 physicalFilterString = comma.join(self.config.physicalFilters)
387 atmosphereTableName =
'NoTableWasUsed'
388 if self.config.atmosphereTableName
is not None:
389 atmosphereTableName = self.config.atmosphereTableName
391 lutSchema = self.
_makeLutSchema_makeLutSchema(physicalFilterString, stdPhysicalFilterString,
394 lutCat = self.
_makeLutCat_makeLutCat(lutSchema, physicalFilterString,
395 stdPhysicalFilterString, atmosphereTableName)
398 def _getStdPhysicalFilterList(self):
399 """Get the standard physical filter lists from config.physicalFilters
400 and config.stdPhysicalFilterOverrideMap
404 stdPhysicalFilters : `list`
406 override = self.config.stdPhysicalFilterOverrideMap
407 return [override.get(physicalFilter, physicalFilter)
for
408 physicalFilter
in self.config.physicalFilters]
410 def _createLutConfig(self, nCcd):
412 Create the fgcmLut config dictionary
417 Number of CCDs in the camera
422 lutConfig[
'logger'] = self.log
423 lutConfig[
'filterNames'] = self.config.physicalFilters
425 lutConfig[
'nCCD'] = nCcd
428 if self.config.atmosphereTableName
is not None:
429 lutConfig[
'atmosphereTableName'] = self.config.atmosphereTableName
432 lutConfig[
'elevation'] = self.config.parameters.elevation
433 lutConfig[
'pmbRange'] = self.config.parameters.pmbRange
434 lutConfig[
'pmbSteps'] = self.config.parameters.pmbSteps
435 lutConfig[
'pwvRange'] = self.config.parameters.pwvRange
436 lutConfig[
'pwvSteps'] = self.config.parameters.pwvSteps
437 lutConfig[
'o3Range'] = self.config.parameters.o3Range
438 lutConfig[
'o3Steps'] = self.config.parameters.o3Steps
439 lutConfig[
'tauRange'] = self.config.parameters.tauRange
440 lutConfig[
'tauSteps'] = self.config.parameters.tauSteps
441 lutConfig[
'alphaRange'] = self.config.parameters.alphaRange
442 lutConfig[
'alphaSteps'] = self.config.parameters.alphaSteps
443 lutConfig[
'zenithRange'] = self.config.parameters.zenithRange
444 lutConfig[
'zenithSteps'] = self.config.parameters.zenithSteps
445 lutConfig[
'pmbStd'] = self.config.parameters.pmbStd
446 lutConfig[
'pwvStd'] = self.config.parameters.pwvStd
447 lutConfig[
'o3Std'] = self.config.parameters.o3Std
448 lutConfig[
'tauStd'] = self.config.parameters.tauStd
449 lutConfig[
'alphaStd'] = self.config.parameters.alphaStd
450 lutConfig[
'airmassStd'] = self.config.parameters.airmassStd
451 lutConfig[
'lambdaRange'] = self.config.parameters.lambdaRange
452 lutConfig[
'lambdaStep'] = self.config.parameters.lambdaStep
453 lutConfig[
'lambdaNorm'] = self.config.parameters.lambdaNorm
457 def _loadThroughputs(self, camera, opticsHandle, sensorHandleDict, filterHandleDict):
458 """Internal method to load throughput data for filters
463 Camera from the butler
464 opticsHandle : `lsst.daf.butler.DeferredDatasetHandle`
465 Reference to optics transmission curve.
466 sensorHandleDict : `dict` of [`int`, `lsst.daf.butler.DeferredDatasetHandle`]
467 Dictionary of references to sensor transmission curves. Key will
469 filterHandleDict : `dict` of [`str`, `lsst.daf.butler.DeferredDatasetHandle`]
470 Dictionary of references to filter transmission curves. Key will
471 be physical filter label.
475 ValueError : Raised
if configured filter name does
not match any of the
476 available filter transmission curves.
481 for detector
in camera:
482 self.
_sensorsTransmission_sensorsTransmission[detector.getId()] = sensorHandleDict[detector.getId()].get()
485 for physicalFilter
in self.config.physicalFilters:
486 self.
_filtersTransmission_filtersTransmission[physicalFilter] = filterHandleDict[physicalFilter].get()
488 def _getThroughputDetector(self, detector, physicalFilter, throughputLambda):
489 """Internal method to get throughput for a detector.
491 Returns the throughput at the center of the detector for a given filter.
495 detector: `lsst.afw.cameraGeom._detector.Detector`
497 physicalFilter: `str`
498 Physical filter label
499 throughputLambda: `np.array(dtype=np.float64)`
500 Wavelength steps (Angstrom)
504 throughput: `np.array(dtype=np.float64)`
505 Throughput (max 1.0) at throughputLambda
508 c = detector.getCenter(afwCameraGeom.FOCAL_PLANE)
509 c.scale(1.0/detector.getPixelSize()[0])
512 wavelengths=throughputLambda)
514 throughput *= self.
_sensorsTransmission_sensorsTransmission[detector.getId()].sampleAt(position=c,
515 wavelengths=throughputLambda)
517 throughput *= self.
_filtersTransmission_filtersTransmission[physicalFilter].sampleAt(position=c,
518 wavelengths=throughputLambda)
521 throughput = np.clip(throughput, 0.0, 1.0)
525 def _makeLutSchema(self, physicalFilterString, stdPhysicalFilterString,
526 atmosphereTableName):
532 physicalFilterString: `str`
533 Combined string of all the physicalFilters
534 stdPhysicalFilterString: `str`
535 Combined string of all the standard physicalFilters
536 atmosphereTableName: `str`
537 Name of the atmosphere table used to generate LUT
541 lutSchema: `afwTable.schema`
546 lutSchema.addField('tablename', type=str, doc=
'Atmosphere table name',
547 size=len(atmosphereTableName))
548 lutSchema.addField(
'elevation', type=float, doc=
"Telescope elevation used for LUT")
549 lutSchema.addField(
'physicalFilters', type=str, doc=
'physicalFilters in LUT',
550 size=len(physicalFilterString))
551 lutSchema.addField(
'stdPhysicalFilters', type=str, doc=
'Standard physicalFilters in LUT',
552 size=len(stdPhysicalFilterString))
553 lutSchema.addField(
'pmb', type=
'ArrayD', doc=
'Barometric Pressure',
555 lutSchema.addField(
'pmbFactor', type=
'ArrayD', doc=
'PMB scaling factor',
557 lutSchema.addField(
'pmbElevation', type=np.float64, doc=
'PMB Scaling at elevation')
558 lutSchema.addField(
'pwv', type=
'ArrayD', doc=
'Preciptable Water Vapor',
560 lutSchema.addField(
'o3', type=
'ArrayD', doc=
'Ozone',
562 lutSchema.addField(
'tau', type=
'ArrayD', doc=
'Aerosol optical depth',
564 lutSchema.addField(
'lambdaNorm', type=np.float64, doc=
'AOD wavelength')
565 lutSchema.addField(
'alpha', type=
'ArrayD', doc=
'Aerosol alpha',
567 lutSchema.addField(
'zenith', type=
'ArrayD', doc=
'Zenith angle',
569 lutSchema.addField(
'nCcd', type=np.int32, doc=
'Number of CCDs')
572 lutSchema.addField(
'pmbStd', type=np.float64, doc=
'PMB Standard')
573 lutSchema.addField(
'pwvStd', type=np.float64, doc=
'PWV Standard')
574 lutSchema.addField(
'o3Std', type=np.float64, doc=
'O3 Standard')
575 lutSchema.addField(
'tauStd', type=np.float64, doc=
'Tau Standard')
576 lutSchema.addField(
'alphaStd', type=np.float64, doc=
'Alpha Standard')
577 lutSchema.addField(
'zenithStd', type=np.float64, doc=
'Zenith angle Standard')
578 lutSchema.addField(
'lambdaRange', type=
'ArrayD', doc=
'Wavelength range',
580 lutSchema.addField(
'lambdaStep', type=np.float64, doc=
'Wavelength step')
581 lutSchema.addField(
'lambdaStd', type=
'ArrayD', doc=
'Standard Wavelength',
583 lutSchema.addField(
'lambdaStdFilter', type=
'ArrayD', doc=
'Standard Wavelength (raw)',
585 lutSchema.addField(
'i0Std', type=
'ArrayD', doc=
'I0 Standard',
587 lutSchema.addField(
'i1Std', type=
'ArrayD', doc=
'I1 Standard',
589 lutSchema.addField(
'i10Std', type=
'ArrayD', doc=
'I10 Standard',
591 lutSchema.addField(
'i2Std', type=
'ArrayD', doc=
'I2 Standard',
593 lutSchema.addField(
'lambdaB', type=
'ArrayD', doc=
'Wavelength for passband (no atm)',
595 lutSchema.addField(
'atmLambda', type=
'ArrayD', doc=
'Atmosphere wavelengths (Angstrom)',
597 lutSchema.addField(
'atmStdTrans', type=
'ArrayD', doc=
'Standard Atmosphere Throughput',
601 lutSchema.addField(
'luttype', type=str, size=20, doc=
'Look-up table type')
602 lutSchema.addField(
'lut', type=
'ArrayF', doc=
'Look-up table for luttype',
607 def _makeLutCat(self, lutSchema, physicalFilterString, stdPhysicalFilterString,
608 atmosphereTableName):
614 lutSchema: `afwTable.schema`
616 physicalFilterString: `str`
617 Combined string of all the physicalFilters
618 stdPhysicalFilterString: `str`
619 Combined string of all the standard physicalFilters
620 atmosphereTableName: `str`
621 Name of the atmosphere table used to generate LUT
626 Lut catalog for persistence
634 lutCat.table.preallocate(14)
637 rec = lutCat.addNew()
639 rec[
'tablename'] = atmosphereTableName
640 rec[
'elevation'] = self.
fgcmLutMakerfgcmLutMaker.atmosphereTable.elevation
641 rec[
'physicalFilters'] = physicalFilterString
642 rec[
'stdPhysicalFilters'] = stdPhysicalFilterString
644 rec[
'pmbFactor'][:] = self.
fgcmLutMakerfgcmLutMaker.pmbFactor
645 rec[
'pmbElevation'] = self.
fgcmLutMakerfgcmLutMaker.pmbElevation
649 rec[
'lambdaNorm'] = self.
fgcmLutMakerfgcmLutMaker.lambdaNorm
658 rec[
'alphaStd'] = self.
fgcmLutMakerfgcmLutMaker.alphaStd
659 rec[
'zenithStd'] = self.
fgcmLutMakerfgcmLutMaker.zenithStd
660 rec[
'lambdaRange'][:] = self.
fgcmLutMakerfgcmLutMaker.lambdaRange
661 rec[
'lambdaStep'] = self.
fgcmLutMakerfgcmLutMaker.lambdaStep
662 rec[
'lambdaStd'][:] = self.
fgcmLutMakerfgcmLutMaker.lambdaStd
663 rec[
'lambdaStdFilter'][:] = self.
fgcmLutMakerfgcmLutMaker.lambdaStdFilter
668 rec[
'lambdaB'][:] = self.
fgcmLutMakerfgcmLutMaker.lambdaB
669 rec[
'atmLambda'][:] = self.
fgcmLutMakerfgcmLutMaker.atmLambda
670 rec[
'atmStdTrans'][:] = self.
fgcmLutMakerfgcmLutMaker.atmStdTrans
672 rec[
'luttype'] =
'I0'
673 rec[
'lut'][:] = self.
fgcmLutMakerfgcmLutMaker.lut[
'I0'].flatten()
676 rec = lutCat.addNew()
677 rec[
'luttype'] =
'I1'
678 rec[
'lut'][:] = self.
fgcmLutMakerfgcmLutMaker.lut[
'I1'].flatten()
680 derivTypes = [
'D_PMB',
'D_LNPWV',
'D_O3',
'D_LNTAU',
'D_ALPHA',
'D_SECZENITH',
681 'D_PMB_I1',
'D_LNPWV_I1',
'D_O3_I1',
'D_LNTAU_I1',
'D_ALPHA_I1',
683 for derivType
in derivTypes:
684 rec = lutCat.addNew()
685 rec[
'luttype'] = derivType
686 rec[
'lut'][:] = self.
fgcmLutMakerfgcmLutMaker.lutDeriv[derivType].flatten()
An immutable representation of a camera.
Defines the fields and offsets for a table.
def _makeLutCat(self, lutSchema, physicalFilterString, stdPhysicalFilterString, atmosphereTableName)
def _fgcmMakeLut(self, camera, opticsHandle, sensorHandleDict, filterHandleDict)
def _getStdPhysicalFilterList(self)
def _createLutConfig(self, nCcd)
def _loadThroughputs(self, camera, opticsHandle, sensorHandleDict, filterHandleDict)
def _makeLutSchema(self, physicalFilterString, stdPhysicalFilterString, atmosphereTableName)
def _getThroughputDetector(self, detector, physicalFilter, throughputLambda)
def runQuantum(self, butlerQC, inputRefs, outputRefs)
def __init__(self, butler=None, initInputs=None, **kwargs)