23"""Perform a single fit cycle of FGCM.
25This task runs a single "fit cycle" of fgcm. Prior to running this task
26one must run both fgcmMakeLut (to construct the atmosphere and instrumental
27look-up-table) and fgcmBuildStars (to extract visits and star observations
30The fgcmFitCycle is meant to be run multiple times, and is tracked by the
31'cycleNumber'. After each run of the fit cycle, diagnostic plots should
32be inspected to set parameters for outlier rejection on the following
33cycle. Please see the fgcmcal Cookbook for details.
45from .utilities
import makeConfigDict, translateFgcmLut, translateVisitCatalog
46from .utilities
import extractReferenceMags
47from .utilities
import makeZptSchema, makeZptCat
48from .utilities
import makeAtmSchema, makeAtmCat, makeStdSchema, makeStdCat
49from .sedterms
import SedboundarytermDict, SedtermDict
50from .focalPlaneProjector
import FocalPlaneProjector
54__all__ = [
'FgcmFitCycleConfig',
'FgcmFitCycleTask']
56MULTIPLE_CYCLES_MAX = 10
60 dimensions=(
"instrument",),
61 defaultTemplates={
"previousCycleNumber":
"-1",
63 camera = connectionTypes.PrerequisiteInput(
64 doc=
"Camera instrument",
66 storageClass=
"Camera",
67 dimensions=(
"instrument",),
71 fgcmLookUpTable = connectionTypes.PrerequisiteInput(
72 doc=(
"Atmosphere + instrument look-up-table for FGCM throughput and "
73 "chromatic corrections."),
74 name=
"fgcmLookUpTable",
75 storageClass=
"Catalog",
76 dimensions=(
"instrument",),
80 fgcmVisitCatalog = connectionTypes.Input(
81 doc=
"Catalog of visit information for fgcm",
82 name=
"fgcmVisitCatalog",
83 storageClass=
"Catalog",
84 dimensions=(
"instrument",),
88 fgcmStarObservationsParquet = connectionTypes.Input(
89 doc=(
"Catalog of star observations for fgcm, in parquet format. "
90 "Used if useParquetCatalogFormat is True."),
91 name=
"fgcm_star_observations",
92 storageClass=
"ArrowAstropy",
93 dimensions=(
"instrument",),
97 fgcmStarIdsParquet = connectionTypes.Input(
98 doc=(
"Catalog of fgcm calibration star IDs, in parquet format. "
99 "Used if useParquetCatalogFormat is True."),
100 name=
"fgcm_star_ids",
101 storageClass=
"ArrowAstropy",
102 dimensions=(
"instrument",),
106 fgcmReferenceStarsParquet = connectionTypes.Input(
107 doc=(
"Catalog of fgcm-matched reference stars, in parquet format. "
108 "Used if useParquetCatalogFormat is True."),
109 name=
"fgcm_reference_stars",
110 storageClass=
"ArrowAstropy",
111 dimensions=(
"instrument",),
115 fgcmStarObservations = connectionTypes.Input(
116 doc=(
"Catalog of star observations for fgcm; old format. "
117 "Used if useParquetCatalogFormat is False."),
118 name=
"fgcmStarObservations",
119 storageClass=
"Catalog",
120 dimensions=(
"instrument",),
124 fgcmStarIds = connectionTypes.Input(
125 doc=(
"Catalog of fgcm calibration star IDs. "
126 "Used if useParquetCatalogFormat is False."),
128 storageClass=
"Catalog",
129 dimensions=(
"instrument",),
133 fgcmStarIndices = connectionTypes.Input(
134 doc=(
"Catalog of fgcm calibration star indices; old format."
135 "Used if useParquetCatalogFormat is False."),
136 name=
"fgcmStarIndices",
137 storageClass=
"Catalog",
138 dimensions=(
"instrument",),
142 fgcmReferenceStars = connectionTypes.Input(
143 doc=(
"Catalog of fgcm-matched reference stars; old format."
144 "Used if useParquetCatalogFormat is False."),
145 name=
"fgcmReferenceStars",
146 storageClass=
"Catalog",
147 dimensions=(
"instrument",),
151 def __init__(self, *, config=None):
152 super().__init__(config=config)
157 if str(int(config.connections.cycleNumber)) != config.connections.cycleNumber:
158 raise ValueError(
"cycleNumber must be of integer format")
159 if str(int(config.connections.previousCycleNumber)) != config.connections.previousCycleNumber:
160 raise ValueError(
"previousCycleNumber must be of integer format")
164 if int(config.connections.previousCycleNumber) != (int(config.connections.cycleNumber) - 1):
165 raise ValueError(
"previousCycleNumber must be 1 less than cycleNumber")
167 inputAndOutputConnections = [
168 (
"FitParameters",
"Catalog",
"Catalog of fgcm fit parameters."),
169 (
"FlaggedStars",
"Catalog",
"Catalog of flagged stars for fgcm calibration."),
171 multicycleOutputConnections = [
172 (
"OutputConfig",
"Config",
"Configuration for next fgcm fit cycle."),
174 optionalZpOutputConnections = [
175 (
"Zeropoints",
"Catalog",
"Catalog of fgcm zeropoint data."),
176 (
"AtmosphereParameters",
"Catalog",
"Catalog of atmospheric fit parameters."),
178 optionalStarOutputConnections = [
179 (
"StandardStars",
"SimpleCatalog",
"Catalog of standard star magnitudes."),
184 physical_filters = []
186 for key, val
in config.physicalFilterMap.items():
189 physical_filters.append(key.replace(
"-",
"_").replace(
" ",
"_").replace(
"~",
"_"))
191 epochs = [f
"epoch{i}" for i
in range(len(config.epochMjds))]
195 (
"Zeropoints_Plot",
"Plot",
"Plot of fgcm zeropoints."),
196 (
"ExpgrayDeep_Plot",
"Plot",
"Plot of gray term per exposure per time for deep fields."),
197 (
"NightlyAlpha_Plot",
"Plot",
"Plot of nightly AOD alpha term."),
198 (
"NightlyTau_Plot",
"Plot",
"Plot of nightly aerosol optical depth (tau)."),
199 (
"NightlyPwv_Plot",
"Plot",
"Plot of nightly water vapor."),
200 (
"NightlyO3_Plot",
"Plot",
"Plot of nightly ozone."),
201 (
"FilterOffsets_Plot",
"Plot",
"Plot of in-band filter offsets."),
202 (
"AbsThroughputs_Plot",
"Plot",
"Plot of absolute throughput fractions."),
203 (
"QESysWashesInitial_Plot",
"Plot",
"Plot of initial system QE with mirror washes."),
204 (
"QESysWashesFinal_Plot",
"Plot",
"Plot of final system QE with mirror washes."),
205 (
"RpwvVsRpwvInput_Plot",
207 "Plot of change in per-visit ``retrieved`` PWV from previous fit cycle."),
208 (
"RpwvVsRpwvSmooth_Plot",
210 "Plot of per-visit ``retrieved`` PWV vs. smoothed PWV."),
211 (
"ModelPwvVsRpwv_Plot",
213 "Plot of model PWV vs. per-visit ``retrieved`` PWV."),
216 "Plot of chisq as a function of iteration."),
219 plotConnections.extend(
221 (f
"Apercorr_{band}_Plot",
"Plot",
"Plot of fgcm aperture corrections."),
222 (f
"EpsilonGlobal_{band}_Plot",
224 "Plot of global background over/undersubtraction."),
225 (f
"EpsilonMap_{band}_Plot",
227 "Map of spatially varying background over/undersubtraction."),
228 (f
"ExpgrayInitial_{band}_Plot",
230 "Histogram of initial gray term per exposure."),
231 (f
"CompareRedblueExpgray_{band}_Plot",
233 "Plot of red/blue split gray term per exposure"),
234 (f
"Expgray_{band}_Plot",
"Plot",
"Histogram of gray term per exposure."),
235 (f
"ExpgrayAirmass_{band}_Plot",
237 "Plot of exposure gray term as a function of airmass."),
238 (f
"ExpgrayCompareMjdRedblue_{band}_Plot",
240 "Plot of red/blue split gray term per exposure as a function of time."),
241 (f
"ExpgrayUT_{band}_Plot",
243 "Plot of grey term per exposure as a function of time of night."),
244 (f
"ExpgrayCompareBands_{band}_Plot",
246 "Plot of gray term per exposure between bands nearby in time."),
247 (f
"ExpgrayReference_{band}_Plot",
249 "Histogram of gray term per exposure compared to reference mags."),
250 (f
"QESysRefstarsStdInitial_{band}_Plot",
252 "Plot of reference mag - calibrated (standard) mag vs. time (before fit)."),
253 (f
"QESysRefstarsStdFinal_{band}_Plot",
255 "Plot of reference mag - calibrated (standard) mag vs. time (after fit)."),
256 (f
"QESysRefstarsObsInitial_{band}_Plot",
258 "Plot of reference mag - observed (instrumental) mag vs. time (before fit)."),
259 (f
"QESysRefstarsObsFinal_{band}_Plot",
261 "Plot of reference mag - observed (instrumental) mag vs. time (after fit)."),
262 (f
"ModelMagerrInitial_{band}_Plot",
264 "Plots for magnitude error model, initial estimate."),
265 (f
"ModelMagerrPostfit_{band}_Plot",
267 "Plots for magnitude error model, after fitting."),
268 (f
"SigmaFgcmAllStars_{band}_Plot",
270 "Histograms for intrinsic scatter for all bright stars."),
271 (f
"SigmaFgcmReservedStars_{band}_Plot",
273 "Histograms for intrinsic scatter for reserved bright stars."),
274 (f
"SigmaFgcmReservedStarsCrunched_{band}_Plot",
276 "Histograms for intrinsic scatter for reserved bright stars (after gray correction)."),
277 (f
"SigmaCal_{band}_Plot",
279 "Plot showing scatter as a function of systematic error floor."),
280 (f
"SigmaRef_{band}_Plot",
282 "Histograms of scatter with respect to reference stars."),
283 (f
"RefResidVsColorAll_{band}_Plot",
285 "Plot of reference star residuals vs. color (all stars)."),
286 (f
"RefResidVsColorCut_{band}_Plot",
288 "Plot of reference star residuals vs. color (reference star color cuts)."),
291 for physical_filter
in physical_filters:
292 plotConnections.extend(
294 (f
"I1R1_{physical_filter}_Plot",
"Plot",
"Plot of fgcm R1 vs. I1."),
295 (f
"I1_{physical_filter}_Plot",
"Plot",
"Focal plane map of fgcm I1."),
296 (f
"R1_{physical_filter}_Plot",
"Plot",
"Focal plane map of fgcm R1."),
297 (f
"R1mI1_{physical_filter}_Plot",
"Plot",
"Focal plane map of fgcm R1 - I1."),
298 (f
"R1mI1_vs_mjd_{physical_filter}_Plot",
"Plot",
"R1 - I1 residuals vs. mjd."),
299 (f
"CompareRedblueMirrorchrom_{physical_filter}_Plot",
301 "Comparison of mirror chromaticity residuals for red/blue stars."),
302 (f
"CcdChromaticity_{physical_filter}_Plot",
304 "Focal plane map of fgcm ccd chromaticity."),
305 (f
"EpsilonDetector_{physical_filter}_Plot",
307 "Focal plane map of background over/undersubtraction."),
308 (f
"EpsilonDetectorMatchscale_{physical_filter}_Plot",
310 "Focal plane map of background over/undersubtraction."),
314 plotConnections.extend(
317 f
"Superstar_{physical_filter}_{epoch}_Plot",
319 "Plot of illumination Correction.",
324 if config.doMultipleCycles:
328 for cycle
in range(config.multipleCyclesFinalCycleNumber):
329 outputConnections = copy.copy(inputAndOutputConnections)
330 outputConnections.extend(multicycleOutputConnections)
331 if config.outputZeropointsBeforeFinalCycle:
332 outputConnections.extend(optionalZpOutputConnections)
333 if config.outputStandardsBeforeFinalCycle:
334 outputConnections.extend(optionalStarOutputConnections)
342 if cycle == (config.multipleCyclesFinalCycleNumber - 1) \
343 or config.doPlotsBeforeFinalCycles:
344 outputConnections.extend(plotConnections)
346 for (name, storageClass, doc)
in outputConnections:
347 connectionName = f
"fgcm_Cycle{cycle}_{name}"
348 storageName = connectionName
349 outConnection = connectionTypes.Output(
351 storageClass=storageClass,
353 dimensions=(
"instrument",),
355 setattr(self, connectionName, outConnection)
358 outputConnections = copy.copy(inputAndOutputConnections)
359 outputConnections.extend(multicycleOutputConnections)
360 outputConnections.extend(optionalZpOutputConnections)
361 outputConnections.extend(optionalStarOutputConnections)
363 outputConnections.extend(plotConnections)
364 for (name, storageClass, doc)
in outputConnections:
365 connectionName = f
"fgcm_Cycle{config.multipleCyclesFinalCycleNumber}_{name}"
366 storageName = connectionName
367 outConnection = connectionTypes.Output(
369 storageClass=storageClass,
371 dimensions=(
"instrument",),
373 setattr(self, connectionName, outConnection)
376 if config.cycleNumber > 0:
377 inputConnections = copy.copy(inputAndOutputConnections)
379 inputConnections = []
380 outputConnections = copy.copy(inputAndOutputConnections)
385 if config.isFinalCycle
or config.outputZeropointsBeforeFinalCycle:
386 outputConnections.extend(optionalZpOutputConnections)
387 if config.isFinalCycle
or config.outputStandardsBeforeFinalCycle:
388 outputConnections.extend(optionalStarOutputConnections)
391 outputConnections.extend(plotConnections)
393 for (name, storageClass, doc)
in inputConnections:
394 connectionName = f
"fgcm{name}Input"
395 storageName = f
"fgcm_Cycle{config.cycleNumber - 1}_{name}"
396 inConnection = connectionTypes.PrerequisiteInput(
398 storageClass=storageClass,
400 dimensions=(
"instrument",),
402 setattr(self, connectionName, inConnection)
404 for (name, storageClass, doc)
in outputConnections:
405 connectionName = f
"fgcm{name}"
406 storageName = f
"fgcm_Cycle{config.cycleNumber}_{name}"
408 if storageClass ==
"Plot":
409 connectionName = storageName
410 outConnection = connectionTypes.Output(
412 storageClass=storageClass,
414 dimensions=(
"instrument",),
416 setattr(self, connectionName, outConnection)
418 if not config.doReferenceCalibration:
419 self.inputs.remove(
"fgcmReferenceStars")
420 self.inputs.remove(
"fgcmReferenceStarsParquet")
422 if config.useParquetCatalogFormat:
423 self.inputs.remove(
"fgcmStarObservations")
424 self.inputs.remove(
"fgcmStarIds")
425 self.inputs.remove(
"fgcmStarIndices")
426 if config.doReferenceCalibration:
427 self.inputs.remove(
"fgcmReferenceStars")
429 self.inputs.remove(
"fgcmStarObservationsParquet")
430 self.inputs.remove(
"fgcmStarIdsParquet")
431 if config.doReferenceCalibration:
432 self.inputs.remove(
"fgcmReferenceStarsParquet")
435class FgcmFitCycleConfig(pipeBase.PipelineTaskConfig,
436 pipelineConnections=FgcmFitCycleConnections):
437 """Config for FgcmFitCycle"""
439 doMultipleCycles = pexConfig.Field(
440 doc=
"Run multiple fit cycles in one task",
444 useParquetCatalogFormat = pexConfig.Field(
445 doc=
"Use parquet catalog format?",
449 multipleCyclesFinalCycleNumber = pexConfig.RangeField(
450 doc=(
"Final cycle number in multiple cycle mode. The initial cycle "
451 "is 0, with limited parameters fit. The next cycle is 1 with "
452 "full parameter fit. The final cycle is a clean-up with no "
453 "parameters fit. There will be a total of "
454 "(multipleCycleFinalCycleNumber + 1) cycles run, and the final "
455 "cycle number cannot be less than 2."),
459 max=MULTIPLE_CYCLES_MAX,
462 bands = pexConfig.ListField(
463 doc=
"Bands to run calibration",
467 fitBands = pexConfig.ListField(
468 doc=(
"Bands to use in atmospheric fit. The bands not listed here will have "
469 "the atmosphere constrained from the 'fitBands' on the same night. "
470 "Must be a subset of `config.bands`"),
474 requiredBands = pexConfig.ListField(
475 doc=(
"Bands that are required for a star to be considered a calibration star. "
476 "Must be a subset of `config.bands`"),
480 physicalFilterMap = pexConfig.DictField(
481 doc=
"Mapping from 'physicalFilter' to band.",
486 doReferenceCalibration = pexConfig.Field(
487 doc=
"Use reference catalog as additional constraint on calibration",
491 refStarSnMin = pexConfig.Field(
492 doc=
"Reference star signal-to-noise minimum to use in calibration. Set to <=0 for no cut.",
496 refStarOutlierNSig = pexConfig.Field(
497 doc=(
"Number of sigma compared to average mag for reference star to be considered an outlier. "
498 "Computed per-band, and if it is an outlier in any band it is rejected from fits."),
502 applyRefStarColorCuts = pexConfig.Field(
503 doc=(
"Apply color cuts defined in ``starColorCuts`` to reference stars? "
504 "These cuts are in addition to any cuts defined in ``refStarColorCuts``"),
508 refStarMaxFracUse = pexConfig.Field(
509 doc=(
"Maximum fraction of reference stars to use in the fit. Remainder will "
510 "be used only for validation."),
514 useExposureReferenceOffset = pexConfig.Field(
515 doc=(
"Use per-exposure (visit) offsets between calibrated stars and reference stars "
516 "for final zeropoints? This may help uniformity for disjoint surveys."),
520 nCore = pexConfig.Field(
521 doc=
"Number of cores to use",
524 deprecated=
"Number of cores is deprecated as a config, and will be removed after v27. "
525 "Please use ``pipetask run --cores-per-quantum`` instead.",
527 nStarPerRun = pexConfig.Field(
528 doc=
"Number of stars to run in each chunk",
532 nExpPerRun = pexConfig.Field(
533 doc=
"Number of exposures to run in each chunk",
537 reserveFraction = pexConfig.Field(
538 doc=
"Fraction of stars to reserve for testing",
542 freezeStdAtmosphere = pexConfig.Field(
543 doc=
"Freeze atmosphere parameters to standard (for testing)",
547 precomputeSuperStarInitialCycle = pexConfig.Field(
548 doc=
"Precompute superstar flat for initial cycle",
552 superStarSubCcdDict = pexConfig.DictField(
553 doc=(
"Per-band specification on whether to compute superstar flat on sub-ccd scale. "
554 "Must have one entry per band."),
559 superStarSubCcdChebyshevOrder = pexConfig.Field(
560 doc=(
"Order of the 2D chebyshev polynomials for sub-ccd superstar fit. "
561 "Global default is first-order polynomials, and should be overridden "
562 "on a camera-by-camera basis depending on the ISR."),
566 superStarSubCcdTriangular = pexConfig.Field(
567 doc=(
"Should the sub-ccd superstar chebyshev matrix be triangular to "
568 "suppress high-order cross terms?"),
572 superStarSigmaClip = pexConfig.Field(
573 doc=
"Number of sigma to clip outliers when selecting for superstar flats",
577 superStarPlotCcdResiduals = pexConfig.Field(
578 doc=
"If plotting is enabled, should per-detector residuals be plotted? "
579 "This may produce a lot of output, and should be used only for "
580 "debugging purposes.",
584 focalPlaneSigmaClip = pexConfig.Field(
585 doc=
"Number of sigma to clip outliers per focal-plane.",
589 ccdGraySubCcdDict = pexConfig.DictField(
590 doc=(
"Per-band specification on whether to compute achromatic per-ccd residual "
591 "('ccd gray') on a sub-ccd scale."),
596 ccdGraySubCcdChebyshevOrder = pexConfig.Field(
597 doc=
"Order of the 2D chebyshev polynomials for sub-ccd gray fit.",
601 ccdGraySubCcdTriangular = pexConfig.Field(
602 doc=(
"Should the sub-ccd gray chebyshev matrix be triangular to "
603 "suppress high-order cross terms?"),
607 ccdGrayFocalPlaneDict = pexConfig.DictField(
608 doc=(
"Per-band specification on whether to compute focal-plane residual "
609 "('ccd gray') corrections."),
614 ccdGrayFocalPlaneFitMinCcd = pexConfig.Field(
615 doc=(
"Minimum number of 'good' CCDs required to perform focal-plane "
616 "gray corrections. If there are fewer good CCDs then the gray "
617 "correction is computed per-ccd."),
621 ccdGrayFocalPlaneChebyshevOrder = pexConfig.Field(
622 doc=
"Order of the 2D chebyshev polynomials for focal plane fit.",
626 cycleNumber = pexConfig.Field(
627 doc=(
"FGCM fit cycle number. This is automatically incremented after each run "
628 "and stage of outlier rejection. See cookbook for details."),
632 isFinalCycle = pexConfig.Field(
633 doc=(
"Is this the final cycle of the fitting? Will automatically compute final "
634 "selection of stars and photometric exposures, and will output zeropoints "
635 "and standard stars for use in fgcmOutputProducts"),
639 maxIterBeforeFinalCycle = pexConfig.Field(
640 doc=(
"Maximum fit iterations, prior to final cycle. The number of iterations "
641 "will always be 0 in the final cycle for cleanup and final selection."),
645 deltaMagBkgOffsetPercentile = pexConfig.Field(
646 doc=(
"Percentile brightest stars on a visit/ccd to use to compute net "
647 "offset from local background subtraction."),
651 deltaMagBkgPerCcd = pexConfig.Field(
652 doc=(
"Compute net offset from local background subtraction per-ccd? "
653 "Otherwise, use computation per visit."),
657 utBoundary = pexConfig.Field(
658 doc=
"Boundary (in UTC) from day-to-day",
662 washMjds = pexConfig.ListField(
663 doc=
"Mirror wash MJDs",
667 epochMjds = pexConfig.ListField(
668 doc=
"Epoch boundaries in MJD",
672 minObsPerBand = pexConfig.Field(
673 doc=
"Minimum good observations per band",
679 latitude = pexConfig.Field(
680 doc=
"Observatory latitude",
684 mirrorArea = pexConfig.Field(
685 doc=
"Mirror area (square meters) of telescope. If not set, will "
686 "try to estimate from camera.telescopeDiameter.",
691 cameraGain = pexConfig.Field(
692 doc=
"Average camera gain. If not set, will use the median of the "
693 "camera model/detector/amplifier gains.",
698 defaultCameraOrientation = pexConfig.Field(
699 doc=
"Default camera orientation for QA plots.",
703 brightObsGrayMax = pexConfig.Field(
704 doc=
"Maximum gray extinction to be considered bright observation",
708 minStarPerCcd = pexConfig.Field(
709 doc=(
"Minimum number of good stars per CCD to be used in calibration fit. "
710 "CCDs with fewer stars will have their calibration estimated from other "
711 "CCDs in the same visit, with zeropoint error increased accordingly."),
715 minCcdPerExp = pexConfig.Field(
716 doc=(
"Minimum number of good CCDs per exposure/visit to be used in calibration fit. "
717 "Visits with fewer good CCDs will have CCD zeropoints estimated where possible."),
721 maxCcdGrayErr = pexConfig.Field(
722 doc=
"Maximum error on CCD gray offset to be considered photometric",
726 minStarPerExp = pexConfig.Field(
727 doc=(
"Minimum number of good stars per exposure/visit to be used in calibration fit. "
728 "Visits with fewer good stars will have CCD zeropoints estimated where possible."),
732 minExpPerNight = pexConfig.Field(
733 doc=
"Minimum number of good exposures/visits to consider a partly photometric night",
737 expGrayInitialCut = pexConfig.Field(
738 doc=(
"Maximum exposure/visit gray value for initial selection of possible photometric "
743 expGrayPhotometricCutDict = pexConfig.DictField(
744 doc=(
"Per-band specification on maximum (negative) achromatic exposure residual "
745 "('gray term') for a visit to be considered photometric. Must have one "
746 "entry per band. Broad-band filters should be -0.05."),
751 expGrayHighCutDict = pexConfig.DictField(
752 doc=(
"Per-band specification on maximum (positive) achromatic exposure residual "
753 "('gray term') for a visit to be considered photometric. Must have one "
754 "entry per band. Broad-band filters should be 0.2."),
759 expGrayRecoverCut = pexConfig.Field(
760 doc=(
"Maximum (negative) exposure gray to be able to recover bad ccds via interpolation. "
761 "Visits with more gray extinction will only get CCD zeropoints if there are "
762 "sufficient star observations (minStarPerCcd) on that CCD."),
766 expVarGrayPhotometricCutDict = pexConfig.DictField(
767 doc=(
"Per-band specification on maximum exposure variance to be considered possibly "
768 "photometric. Must have one entry per band. Broad-band filters should be "
774 expGrayErrRecoverCut = pexConfig.Field(
775 doc=(
"Maximum exposure gray error to be able to recover bad ccds via interpolation. "
776 "Visits with more gray variance will only get CCD zeropoints if there are "
777 "sufficient star observations (minStarPerCcd) on that CCD."),
781 aperCorrFitNBins = pexConfig.Field(
782 doc=(
"Number of aperture bins used in aperture correction fit. When set to 0"
783 "no fit will be performed, and the config.aperCorrInputSlopes will be "
784 "used if available."),
788 aperCorrInputSlopeDict = pexConfig.DictField(
789 doc=(
"Per-band specification of aperture correction input slope parameters. These "
790 "are used on the first fit iteration, and aperture correction parameters will "
791 "be updated from the data if config.aperCorrFitNBins > 0. It is recommended "
792 "to set this when there is insufficient data to fit the parameters (e.g. "
798 sedboundaryterms = pexConfig.ConfigField(
799 doc=
"Mapping from bands to SED boundary term names used is sedterms.",
800 dtype=SedboundarytermDict,
802 sedterms = pexConfig.ConfigField(
803 doc=
"Mapping from terms to bands for fgcm linear SED approximations.",
806 sigFgcmMaxErr = pexConfig.Field(
807 doc=
"Maximum mag error for fitting sigma_FGCM",
811 sigFgcmMaxEGrayDict = pexConfig.DictField(
812 doc=(
"Per-band specification for maximum (absolute) achromatic residual (gray value) "
813 "for observations in sigma_fgcm (raw repeatability). Broad-band filters "
819 ccdGrayMaxStarErr = pexConfig.Field(
820 doc=(
"Maximum error on a star observation to use in ccd gray (achromatic residual) "
825 approxThroughputDict = pexConfig.DictField(
826 doc=(
"Per-band specification of the approximate overall throughput at the start of "
827 "calibration observations. Must have one entry per band. Typically should "
833 sigmaCalRange = pexConfig.ListField(
834 doc=
"Allowed range for systematic error floor estimation",
836 default=(0.001, 0.003),
838 sigmaCalFitPercentile = pexConfig.ListField(
839 doc=
"Magnitude percentile range to fit systematic error floor",
841 default=(0.05, 0.15),
843 sigmaCalPlotPercentile = pexConfig.ListField(
844 doc=
"Magnitude percentile range to plot systematic error floor",
846 default=(0.05, 0.95),
848 sigma0Phot = pexConfig.Field(
849 doc=
"Systematic error floor for all zeropoints",
853 mapLongitudeRef = pexConfig.Field(
854 doc=
"Reference longitude for plotting maps",
858 mapNSide = pexConfig.Field(
859 doc=
"Healpix nside for plotting maps",
863 outfileBase = pexConfig.Field(
864 doc=
"Filename start for plot output files",
868 starColorCuts = pexConfig.ListField(
869 doc=(
"Encoded star-color cuts (using calibration star colors). "
870 "This is a list with each entry a string of the format "
871 "``band1,band2,low,high`` such that only stars of color "
872 "low < band1 - band2 < high will be used for calibration."),
874 default=(
"NO_DATA",),
876 refStarColorCuts = pexConfig.ListField(
877 doc=(
"Encoded star color cuts specifically to apply to reference stars. "
878 "This is a list with each entry a string of the format "
879 "``band1,band2,low,high`` such that only stars of color "
880 "low < band1 - band2 < high will be used as reference stars."),
882 default=(
"NO_DATA",),
884 colorSplitBands = pexConfig.ListField(
885 doc=
"Band names to use to split stars by color. Must have 2 entries.",
890 modelMagErrors = pexConfig.Field(
891 doc=
"Should FGCM model the magnitude errors from sky/fwhm? (False means trust inputs)",
895 useQuadraticPwv = pexConfig.Field(
896 doc=
"Model PWV with a quadratic term for variation through the night?",
900 instrumentParsPerBand = pexConfig.Field(
901 doc=(
"Model instrumental parameters per band? "
902 "Otherwise, instrumental parameters (QE changes with time) are "
903 "shared among all bands."),
907 instrumentSlopeMinDeltaT = pexConfig.Field(
908 doc=(
"Minimum time change (in days) between observations to use in constraining "
909 "instrument slope."),
913 fitMirrorChromaticity = pexConfig.Field(
914 doc=
"Fit (intraband) mirror chromatic term?",
918 fitCcdChromaticityDict = pexConfig.DictField(
919 doc=
"Specification on whether to compute first-order quantum efficiency (QE) "
920 "adjustments. Key is band, and value will be True or False. Any band "
921 "not explicitly specified will default to False.",
926 coatingMjds = pexConfig.ListField(
927 doc=
"Mirror coating dates in MJD",
931 outputStandardsBeforeFinalCycle = pexConfig.Field(
932 doc=
"Output standard stars prior to final cycle? Used in debugging.",
936 outputZeropointsBeforeFinalCycle = pexConfig.Field(
937 doc=
"Output standard stars prior to final cycle? Used in debugging.",
941 useRepeatabilityForExpGrayCutsDict = pexConfig.DictField(
942 doc=(
"Per-band specification on whether to use star repeatability (instead of exposures) "
943 "for computing photometric cuts. Recommended for tract mode or bands with few visits."),
948 autoPhotometricCutNSig = pexConfig.Field(
949 doc=(
"Number of sigma for automatic computation of (low) photometric cut. "
950 "Cut is based on exposure gray width (per band), unless "
951 "useRepeatabilityForExpGrayCuts is set, in which case the star "
952 "repeatability is used (also per band)."),
956 autoHighCutNSig = pexConfig.Field(
957 doc=(
"Number of sigma for automatic computation of (high) outlier cut. "
958 "Cut is based on exposure gray width (per band), unless "
959 "useRepeatabilityForExpGrayCuts is set, in which case the star "
960 "repeatability is used (also per band)."),
964 quietMode = pexConfig.Field(
965 doc=
"Be less verbose with logging.",
969 doPlots = pexConfig.Field(
970 doc=
"Make fgcm QA plots.",
974 doPlotsBeforeFinalCycles = pexConfig.Field(
975 doc=
"Make fgcm QA plots before the final two fit cycles? This only applies in"
976 "multi-cycle mode, and if doPlots is True. These are typically the most"
977 "important QA plots to inspect.",
981 randomSeed = pexConfig.Field(
982 doc=
"Random seed for fgcm for consistency in tests.",
987 deltaAperFitMinNgoodObs = pexConfig.Field(
988 doc=
"Minimum number of good observations to use mean delta-aper values in fits.",
992 deltaAperFitPerCcdNx = pexConfig.Field(
993 doc=(
"Number of x bins per ccd when computing delta-aper background offsets. "
994 "Only used when ``doComputeDeltaAperPerCcd`` is True."),
998 deltaAperFitPerCcdNy = pexConfig.Field(
999 doc=(
"Number of y bins per ccd when computing delta-aper background offsets. "
1000 "Only used when ``doComputeDeltaAperPerCcd`` is True."),
1004 deltaAperFitSpatialNside = pexConfig.Field(
1005 doc=
"Healpix nside to compute spatial delta-aper background offset maps.",
1009 deltaAperInnerRadiusArcsec = pexConfig.Field(
1010 doc=(
"Inner radius used to compute deltaMagAper (arcseconds). "
1011 "Must be positive and less than ``deltaAperOuterRadiusArcsec`` if "
1012 "any of ``doComputeDeltaAperPerVisit``, ``doComputeDeltaAperPerStar``, "
1013 "``doComputeDeltaAperMap``, ``doComputeDeltaAperPerCcd`` are set."),
1017 deltaAperOuterRadiusArcsec = pexConfig.Field(
1018 doc=(
"Outer radius used to compute deltaMagAper (arcseconds). "
1019 "Must be positive and greater than ``deltaAperInnerRadiusArcsec`` if "
1020 "any of ``doComputeDeltaAperPerVisit``, ``doComputeDeltaAperPerStar``, "
1021 "``doComputeDeltaAperMap``, ``doComputeDeltaAperPerCcd`` are set."),
1025 doComputeDeltaAperPerVisit = pexConfig.Field(
1026 doc=(
"Do the computation of delta-aper background offsets per visit? "
1027 "Note: this option can be very slow when there are many visits."),
1031 doComputeDeltaAperPerStar = pexConfig.Field(
1032 doc=
"Do the computation of delta-aper mean values per star?",
1036 doComputeDeltaAperMap = pexConfig.Field(
1037 doc=(
"Do the computation of delta-aper spatial maps? "
1038 "This is only used if ``doComputeDeltaAperPerStar`` is True,"),
1042 doComputeDeltaAperPerCcd = pexConfig.Field(
1043 doc=
"Do the computation of per-ccd delta-aper background offsets?",
1051 if self.connections.previousCycleNumber != str(self.cycleNumber - 1):
1052 msg =
"cycleNumber in template must be connections.previousCycleNumber + 1"
1053 raise RuntimeError(msg)
1054 if self.connections.cycleNumber != str(self.cycleNumber):
1055 msg =
"cycleNumber in template must be equal to connections.cycleNumber"
1056 raise RuntimeError(msg)
1058 for band
in self.fitBands:
1059 if band
not in self.bands:
1060 msg =
'fitBand %s not in bands' % (band)
1061 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.fitBands, self, msg)
1062 for band
in self.requiredBands:
1063 if band
not in self.bands:
1064 msg =
'requiredBand %s not in bands' % (band)
1065 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.requiredBands, self, msg)
1066 for band
in self.colorSplitBands:
1067 if band
not in self.bands:
1068 msg =
'colorSplitBand %s not in bands' % (band)
1069 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.colorSplitBands, self, msg)
1070 for band
in self.bands:
1071 if band
not in self.superStarSubCcdDict:
1072 msg =
'band %s not in superStarSubCcdDict' % (band)
1073 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.superStarSubCcdDict,
1075 if band
not in self.ccdGraySubCcdDict:
1076 msg =
'band %s not in ccdGraySubCcdDict' % (band)
1077 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.ccdGraySubCcdDict,
1079 if band
not in self.expGrayPhotometricCutDict:
1080 msg =
'band %s not in expGrayPhotometricCutDict' % (band)
1081 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.expGrayPhotometricCutDict,
1083 if band
not in self.expGrayHighCutDict:
1084 msg =
'band %s not in expGrayHighCutDict' % (band)
1085 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.expGrayHighCutDict,
1087 if band
not in self.expVarGrayPhotometricCutDict:
1088 msg =
'band %s not in expVarGrayPhotometricCutDict' % (band)
1089 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.expVarGrayPhotometricCutDict,
1091 if band
not in self.sigFgcmMaxEGrayDict:
1092 msg =
'band %s not in sigFgcmMaxEGrayDict' % (band)
1093 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.sigFgcmMaxEGrayDict,
1095 if band
not in self.approxThroughputDict:
1096 msg =
'band %s not in approxThroughputDict' % (band)
1097 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.approxThroughputDict,
1099 if band
not in self.useRepeatabilityForExpGrayCutsDict:
1100 msg =
'band %s not in useRepeatabilityForExpGrayCutsDict' % (band)
1101 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.useRepeatabilityForExpGrayCutsDict,
1104 if self.doComputeDeltaAperPerVisit
or self.doComputeDeltaAperMap \
1105 or self.doComputeDeltaAperPerCcd:
1106 if self.deltaAperInnerRadiusArcsec <= 0.0:
1107 msg =
'deltaAperInnerRadiusArcsec must be positive if deltaAper computations are turned on.'
1108 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.deltaAperInnerRadiusArcsec,
1110 if self.deltaAperOuterRadiusArcsec <= 0.0:
1111 msg =
'deltaAperOuterRadiusArcsec must be positive if deltaAper computations are turned on.'
1112 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.deltaAperOuterRadiusArcsec,
1114 if self.deltaAperOuterRadiusArcsec <= self.deltaAperInnerRadiusArcsec:
1115 msg = (
'deltaAperOuterRadiusArcsec must be greater than deltaAperInnerRadiusArcsec if '
1116 'deltaAper computations are turned on.')
1117 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.deltaAperOuterRadiusArcsec,
1121class FgcmFitCycleTask(pipeBase.PipelineTask):
1123 Run Single fit cycle for FGCM global calibration
1126 ConfigClass = FgcmFitCycleConfig
1127 _DefaultName =
"fgcmFitCycle"
1129 def __init__(self, initInputs=None, **kwargs):
1130 super().__init__(**kwargs)
1132 def runQuantum(self, butlerQC, inputRefs, outputRefs):
1133 camera = butlerQC.get(inputRefs.camera)
1135 nCore = butlerQC.resources.num_cores
1139 handleDict[
'fgcmLookUpTable'] = butlerQC.get(inputRefs.fgcmLookUpTable)
1140 handleDict[
'fgcmVisitCatalog'] = butlerQC.get(inputRefs.fgcmVisitCatalog)
1142 if self.config.useParquetCatalogFormat:
1143 handleDict[
'fgcmStarObservations'] = butlerQC.get(inputRefs.fgcmStarObservationsParquet)
1144 handleDict[
'fgcmStarIds'] = butlerQC.get(inputRefs.fgcmStarIdsParquet)
1145 if self.config.doReferenceCalibration:
1146 handleDict[
'fgcmReferenceStars'] = butlerQC.get(inputRefs.fgcmReferenceStarsParquet)
1148 handleDict[
'fgcmStarObservations'] = butlerQC.get(inputRefs.fgcmStarObservations)
1149 handleDict[
'fgcmStarIds'] = butlerQC.get(inputRefs.fgcmStarIds)
1150 handleDict[
'fgcmStarIndices'] = butlerQC.get(inputRefs.fgcmStarIndices)
1151 if self.config.doReferenceCalibration:
1152 handleDict[
'fgcmReferenceStars'] = butlerQC.get(inputRefs.fgcmReferenceStars)
1153 if self.config.cycleNumber > 0:
1154 handleDict[
'fgcmFlaggedStars'] = butlerQC.get(inputRefs.fgcmFlaggedStarsInput)
1155 handleDict[
'fgcmFitParameters'] = butlerQC.get(inputRefs.fgcmFitParametersInput)
1157 fgcmDatasetDict =
None
1158 if self.config.doMultipleCycles:
1160 config = copy.copy(self.config)
1161 config.update(cycleNumber=0)
1162 for cycle
in range(self.config.multipleCyclesFinalCycleNumber + 1):
1163 if cycle == self.config.multipleCyclesFinalCycleNumber:
1164 config.update(isFinalCycle=
True)
1167 handleDict[
'fgcmFlaggedStars'] = fgcmDatasetDict[
'fgcmFlaggedStars']
1168 handleDict[
'fgcmFitParameters'] = fgcmDatasetDict[
'fgcmFitParameters']
1173 for outputRefName
in outputRefs.keys():
1174 if outputRefName.endswith(
"Plot")
and f
"Cycle{cycle}" in outputRefName:
1175 ref = getattr(outputRefs, outputRefName)
1176 plotHandleDict[outputRefName] = ref
1178 fgcmDatasetDict, config = self._fgcmFitCycle(
1182 plotHandleDict=plotHandleDict,
1186 butlerQC.put(fgcmDatasetDict[
'fgcmFitParameters'],
1187 getattr(outputRefs, f
'fgcm_Cycle{cycle}_FitParameters'))
1188 butlerQC.put(fgcmDatasetDict[
'fgcmFlaggedStars'],
1189 getattr(outputRefs, f
'fgcm_Cycle{cycle}_FlaggedStars'))
1190 butlerQC.put(config,
1191 getattr(outputRefs, f
'fgcm_Cycle{cycle}_OutputConfig'))
1192 if self.outputZeropoints:
1193 butlerQC.put(fgcmDatasetDict[
'fgcmZeropoints'],
1194 getattr(outputRefs, f
'fgcm_Cycle{cycle}_Zeropoints'))
1195 butlerQC.put(fgcmDatasetDict[
'fgcmAtmosphereParameters'],
1196 getattr(outputRefs, f
'fgcm_Cycle{cycle}_AtmosphereParameters'))
1197 if self.outputStandards:
1198 butlerQC.put(fgcmDatasetDict[
'fgcmStandardStars'],
1199 getattr(outputRefs, f
'fgcm_Cycle{cycle}_StandardStars'))
1206 for outputRefName
in outputRefs.keys():
1207 if outputRefName.endswith(
"Plot")
and f
"Cycle{self.config.cycleNumber}" in outputRefName:
1208 ref = getattr(outputRefs, outputRefName)
1209 plotHandleDict[outputRefName] = ref
1211 fgcmDatasetDict, _ = self._fgcmFitCycle(
1216 plotHandleDict=plotHandleDict,
1220 butlerQC.put(fgcmDatasetDict[
'fgcmFitParameters'], outputRefs.fgcmFitParameters)
1221 butlerQC.put(fgcmDatasetDict[
'fgcmFlaggedStars'], outputRefs.fgcmFlaggedStars)
1222 if self.outputZeropoints:
1223 butlerQC.put(fgcmDatasetDict[
'fgcmZeropoints'], outputRefs.fgcmZeropoints)
1224 butlerQC.put(fgcmDatasetDict[
'fgcmAtmosphereParameters'], outputRefs.fgcmAtmosphereParameters)
1225 if self.outputStandards:
1226 butlerQC.put(fgcmDatasetDict[
'fgcmStandardStars'], outputRefs.fgcmStandardStars)
1233 plotHandleDict=None,
1243 camera : `lsst.afw.cameraGeom.Camera`
1245 All handles are `lsst.daf.butler.DeferredDatasetHandle`
1246 handle dictionary with keys:
1248 ``"fgcmLookUpTable"``
1249 handle for the FGCM look-up table.
1250 ``"fgcmVisitCatalog"``
1251 handle for visit summary catalog.
1252 ``"fgcmStarObservations"``
1253 handle for star observation catalog.
1255 handle for star id catalog.
1256 ``"fgcmStarIndices"``
1257 handle for star index catalog.
1258 ``"fgcmReferenceStars"``
1259 handle for matched reference star catalog.
1260 ``"fgcmFlaggedStars"``
1261 handle for flagged star catalog.
1262 ``"fgcmFitParameters"``
1263 handle for fit parameter catalog.
1264 butlerQC : `lsst.pipe.base.QuantumContext`, optional
1265 Quantum context used for serializing plots.
1266 plotHandleDict : `dict` [`str`, `lsst.daf.butler.DatasetRef`], optional
1267 Dictionary of plot dataset refs, keyed by plot name.
1268 config : `lsst.pex.config.Config`, optional
1269 Configuration to use to override self.config.
1270 nCore : `int`, optional
1271 Number of cores to use during fitting.
1272 multiCycle : `bool`, optional
1273 Is this part of a multicycle run?
1277 fgcmDatasetDict : `dict`
1278 Dictionary of datasets to persist.
1280 if config
is not None:
1283 _config = self.config
1286 self.maxIter = _config.maxIterBeforeFinalCycle
1287 self.outputStandards = _config.outputStandardsBeforeFinalCycle
1288 self.outputZeropoints = _config.outputZeropointsBeforeFinalCycle
1289 self.resetFitParameters =
True
1291 if _config.isFinalCycle:
1296 self.outputStandards =
True
1297 self.outputZeropoints =
True
1298 self.resetFitParameters =
False
1300 lutCat = handleDict[
'fgcmLookUpTable'].get()
1301 fgcmLut, lutIndexVals, lutStd = translateFgcmLut(lutCat,
1302 dict(_config.physicalFilterMap))
1306 doPlots = _config.doPlots
1307 if doPlots
and multiCycle:
1308 if _config.cycleNumber < (_config.multipleCyclesFinalCycleNumber - 1) \
1309 and not _config.doPlotsBeforeFinalCycles:
1312 configDict = makeConfigDict(_config, self.log, camera,
1313 self.maxIter, self.resetFitParameters,
1314 self.outputZeropoints,
1315 lutIndexVals[0][
'FILTERNAMES'],
1320 visitCat = handleDict[
'fgcmVisitCatalog'].get()
1321 fgcmExpInfo = translateVisitCatalog(visitCat)
1325 self.config.defaultCameraOrientation)
1327 noFitsDict = {
'lutIndex': lutIndexVals,
1329 'expInfo': fgcmExpInfo,
1330 'focalPlaneProjector': focalPlaneProjector}
1333 fgcmFitCycle = fgcm.FgcmFitCycle(
1336 noFitsDict=noFitsDict,
1339 plotHandleDict=plotHandleDict,
1343 if (fgcmFitCycle.initialCycle):
1345 fgcmPars = fgcm.FgcmParameters.newParsWithArrays(fgcmFitCycle.fgcmConfig,
1349 plotHandleDict=plotHandleDict)
1352 parCat = handleDict[
'fgcmFitParameters']
1354 parCat = handleDict[
'fgcmFitParameters'].get()
1355 inParInfo, inParams, inSuperStar = self._loadParameters(parCat)
1357 fgcmPars = fgcm.FgcmParameters.loadParsWithArrays(fgcmFitCycle.fgcmConfig,
1363 plotHandleDict=plotHandleDict)
1366 fgcmStars = fgcm.FgcmStars(fgcmFitCycle.fgcmConfig, butlerQC=butlerQC, plotHandleDict=plotHandleDict)
1368 starObs = handleDict[
'fgcmStarObservations'].get()
1369 starIds = handleDict[
'fgcmStarIds'].get()
1370 if not self.config.useParquetCatalogFormat:
1371 starIndices = handleDict[
'fgcmStarIndices'].get()
1376 if 'fgcmFlaggedStars' in handleDict:
1378 flaggedStars = handleDict[
'fgcmFlaggedStars']
1380 flaggedStars = handleDict[
'fgcmFlaggedStars'].get()
1381 flagId = flaggedStars[
'objId'][:]
1382 flagFlag = flaggedStars[
'objFlag'][:]
1385 elif self.config.useParquetCatalogFormat:
1391 (flagged,) = (starIds[
'obj_flag'] > 0).nonzero()
1392 flagId = starIds[
'fgcm_id'][flagged]
1393 flagFlag = starIds[
'obj_flag'][flagged]
1398 if _config.doReferenceCalibration:
1399 refStars = handleDict[
'fgcmReferenceStars'].get()
1401 refMag, refMagErr = extractReferenceMags(refStars,
1403 _config.physicalFilterMap)
1405 refId = refStars[
'fgcm_id'][:]
1415 if self.config.useParquetCatalogFormat:
1416 visitIndex = np.searchsorted(fgcmExpInfo[
'VISIT'], starObs[
'visit'])
1418 visitIndex = np.searchsorted(fgcmExpInfo[
'VISIT'], starObs[
'visit'][starIndices[
'obsIndex']])
1427 if self.config.useParquetCatalogFormat:
1430 fgcmStars.loadStars(fgcmPars,
1431 starObs[
'visit'][:],
1432 starObs[
'detector'][:],
1435 starObs[
'inst_mag'][:],
1436 starObs[
'inst_mag_err'][:],
1437 fgcmExpInfo[
'FILTERNAME'][visitIndex],
1438 starIds[
'fgcm_id'][:],
1441 starIds[
'obs_arr_index'][:],
1442 starIds[
'n_obs'][:],
1443 obsX=starObs[
'x'][:],
1444 obsY=starObs[
'y'][:],
1445 obsDeltaMagBkg=starObs[
'delta_mag_bkg'][:],
1446 obsDeltaAper=starObs[
'delta_mag_aper'][:],
1449 refMagErr=refMagErr,
1457 conv = starObs[0][
'ra'].asDegrees() / float(starObs[0][
'ra'])
1459 fgcmStars.loadStars(fgcmPars,
1460 starObs[
'visit'][starIndices[
'obsIndex']],
1461 starObs[
'ccd'][starIndices[
'obsIndex']],
1462 starObs[
'ra'][starIndices[
'obsIndex']] * conv,
1463 starObs[
'dec'][starIndices[
'obsIndex']] * conv,
1464 starObs[
'instMag'][starIndices[
'obsIndex']],
1465 starObs[
'instMagErr'][starIndices[
'obsIndex']],
1466 fgcmExpInfo[
'FILTERNAME'][visitIndex],
1467 starIds[
'fgcm_id'][:],
1470 starIds[
'obsArrIndex'][:],
1472 obsX=starObs[
'x'][starIndices[
'obsIndex']],
1473 obsY=starObs[
'y'][starIndices[
'obsIndex']],
1474 obsDeltaMagBkg=starObs[
'deltaMagBkg'][starIndices[
'obsIndex']],
1475 obsDeltaAper=starObs[
'deltaMagAper'][starIndices[
'obsIndex']],
1476 psfCandidate=starObs[
'psf_candidate'][starIndices[
'obsIndex']],
1479 refMagErr=refMagErr,
1496 fgcmFitCycle.setLUT(fgcmLut)
1497 fgcmFitCycle.setStars(fgcmStars, fgcmPars)
1498 fgcmFitCycle.setPars(fgcmPars)
1501 fgcmFitCycle.finishSetup()
1510 fgcmDatasetDict = self._makeFgcmOutputDatasets(fgcmFitCycle)
1515 updatedPhotometricCutDict = {b: float(fgcmFitCycle.updatedPhotometricCut[i])
for
1516 i, b
in enumerate(_config.bands)}
1517 updatedHighCutDict = {band: float(fgcmFitCycle.updatedHighCut[i])
for
1518 i, band
in enumerate(_config.bands)}
1520 outConfig = copy.copy(_config)
1521 outConfig.update(cycleNumber=(_config.cycleNumber + 1),
1522 precomputeSuperStarInitialCycle=
False,
1523 freezeStdAtmosphere=
False,
1524 expGrayPhotometricCutDict=updatedPhotometricCutDict,
1525 expGrayHighCutDict=updatedHighCutDict)
1527 outConfig.connections.update(previousCycleNumber=str(_config.cycleNumber),
1528 cycleNumber=str(_config.cycleNumber + 1))
1531 configFileName =
'%s_cycle%02d_config.py' % (outConfig.outfileBase,
1532 outConfig.cycleNumber)
1533 outConfig.save(configFileName)
1535 if _config.isFinalCycle == 1:
1537 self.log.info(
"Everything is in place to run fgcmOutputProducts.py")
1539 self.log.info(
"Saved config for next cycle to %s" % (configFileName))
1540 self.log.info(
"Be sure to look at:")
1541 self.log.info(
" config.expGrayPhotometricCut")
1542 self.log.info(
" config.expGrayHighCut")
1543 self.log.info(
"If you are satisfied with the fit, please set:")
1544 self.log.info(
" config.isFinalCycle = True")
1546 fgcmFitCycle.freeSharedMemory()
1548 return fgcmDatasetDict, outConfig
1550 def _loadParameters(self, parCat):
1552 Load FGCM parameters from a previous fit cycle
1556 parCat : `lsst.afw.table.BaseCatalog`
1557 Parameter catalog in afw table form.
1561 inParInfo: `numpy.ndarray`
1562 Numpy array parameter information formatted for input to fgcm
1563 inParameters: `numpy.ndarray`
1564 Numpy array parameter values formatted for input to fgcm
1565 inSuperStar: `numpy.array`
1566 Superstar flat formatted for input to fgcm
1568 parLutFilterNames = np.array(parCat[0][
'lutFilterNames'].split(
','))
1569 parFitBands = np.array(parCat[0][
'fitBands'].split(
','))
1571 inParInfo = np.zeros(1, dtype=[(
'NCCD',
'i4'),
1572 (
'LUTFILTERNAMES', parLutFilterNames.dtype.str,
1573 (parLutFilterNames.size, )),
1574 (
'FITBANDS', parFitBands.dtype.str, (parFitBands.size, )),
1575 (
'LNTAUUNIT',
'f8'),
1576 (
'LNTAUSLOPEUNIT',
'f8'),
1577 (
'ALPHAUNIT',
'f8'),
1578 (
'LNPWVUNIT',
'f8'),
1579 (
'LNPWVSLOPEUNIT',
'f8'),
1580 (
'LNPWVQUADRATICUNIT',
'f8'),
1581 (
'LNPWVGLOBALUNIT',
'f8'),
1583 (
'QESYSUNIT',
'f8'),
1584 (
'FILTEROFFSETUNIT',
'f8'),
1585 (
'HASEXTERNALPWV',
'i2'),
1586 (
'HASEXTERNALTAU',
'i2')])
1587 inParInfo[
'NCCD'] = parCat[
'nCcd']
1588 inParInfo[
'LUTFILTERNAMES'][:] = parLutFilterNames
1589 inParInfo[
'FITBANDS'][:] = parFitBands
1590 inParInfo[
'HASEXTERNALPWV'] = parCat[
'hasExternalPwv']
1591 inParInfo[
'HASEXTERNALTAU'] = parCat[
'hasExternalTau']
1593 inParams = np.zeros(1, dtype=[(
'PARALPHA',
'f8', (parCat[
'parAlpha'].size, )),
1594 (
'PARO3',
'f8', (parCat[
'parO3'].size, )),
1595 (
'PARLNTAUINTERCEPT',
'f8',
1596 (parCat[
'parLnTauIntercept'].size, )),
1597 (
'PARLNTAUSLOPE',
'f8',
1598 (parCat[
'parLnTauSlope'].size, )),
1599 (
'PARLNPWVINTERCEPT',
'f8',
1600 (parCat[
'parLnPwvIntercept'].size, )),
1601 (
'PARLNPWVSLOPE',
'f8',
1602 (parCat[
'parLnPwvSlope'].size, )),
1603 (
'PARLNPWVQUADRATIC',
'f8',
1604 (parCat[
'parLnPwvQuadratic'].size, )),
1605 (
'PARQESYSINTERCEPT',
'f8',
1606 (parCat[
'parQeSysIntercept'].size, )),
1607 (
'COMPQESYSSLOPE',
'f8',
1608 (parCat[
'compQeSysSlope'].size, )),
1609 (
'PARFILTEROFFSET',
'f8',
1610 (parCat[
'parFilterOffset'].size, )),
1611 (
'PARFILTEROFFSETFITFLAG',
'i2',
1612 (parCat[
'parFilterOffsetFitFlag'].size, )),
1613 (
'PARRETRIEVEDLNPWVSCALE',
'f8'),
1614 (
'PARRETRIEVEDLNPWVOFFSET',
'f8'),
1615 (
'PARRETRIEVEDLNPWVNIGHTLYOFFSET',
'f8',
1616 (parCat[
'parRetrievedLnPwvNightlyOffset'].size, )),
1617 (
'COMPABSTHROUGHPUT',
'f8',
1618 (parCat[
'compAbsThroughput'].size, )),
1619 (
'COMPREFOFFSET',
'f8',
1620 (parCat[
'compRefOffset'].size, )),
1621 (
'COMPREFSIGMA',
'f8',
1622 (parCat[
'compRefSigma'].size, )),
1623 (
'COMPMIRRORCHROMATICITY',
'f8',
1624 (parCat[
'compMirrorChromaticity'].size, )),
1625 (
'MIRRORCHROMATICITYPIVOT',
'f8',
1626 (parCat[
'mirrorChromaticityPivot'].size, )),
1627 (
'COMPCCDCHROMATICITY',
'f8',
1628 (parCat[
'compCcdChromaticity'].size, )),
1629 (
'COMPMEDIANSEDSLOPE',
'f8',
1630 (parCat[
'compMedianSedSlope'].size, )),
1631 (
'COMPAPERCORRPIVOT',
'f8',
1632 (parCat[
'compAperCorrPivot'].size, )),
1633 (
'COMPAPERCORRSLOPE',
'f8',
1634 (parCat[
'compAperCorrSlope'].size, )),
1635 (
'COMPAPERCORRSLOPEERR',
'f8',
1636 (parCat[
'compAperCorrSlopeErr'].size, )),
1637 (
'COMPAPERCORRRANGE',
'f8',
1638 (parCat[
'compAperCorrRange'].size, )),
1639 (
'COMPMODELERREXPTIMEPIVOT',
'f8',
1640 (parCat[
'compModelErrExptimePivot'].size, )),
1641 (
'COMPMODELERRFWHMPIVOT',
'f8',
1642 (parCat[
'compModelErrFwhmPivot'].size, )),
1643 (
'COMPMODELERRSKYPIVOT',
'f8',
1644 (parCat[
'compModelErrSkyPivot'].size, )),
1645 (
'COMPMODELERRPARS',
'f8',
1646 (parCat[
'compModelErrPars'].size, )),
1647 (
'COMPEXPGRAY',
'f8',
1648 (parCat[
'compExpGray'].size, )),
1649 (
'COMPVARGRAY',
'f8',
1650 (parCat[
'compVarGray'].size, )),
1651 (
'COMPEXPDELTAMAGBKG',
'f8',
1652 (parCat[
'compExpDeltaMagBkg'].size, )),
1653 (
'COMPNGOODSTARPEREXP',
'i4',
1654 (parCat[
'compNGoodStarPerExp'].size, )),
1655 (
'COMPEXPREFOFFSET',
'f8',
1656 (parCat[
'compExpRefOffset'].size, )),
1657 (
'COMPSIGFGCM',
'f8',
1658 (parCat[
'compSigFgcm'].size, )),
1659 (
'COMPSIGMACAL',
'f8',
1660 (parCat[
'compSigmaCal'].size, )),
1661 (
'COMPRETRIEVEDLNPWV',
'f8',
1662 (parCat[
'compRetrievedLnPwv'].size, )),
1663 (
'COMPRETRIEVEDLNPWVRAW',
'f8',
1664 (parCat[
'compRetrievedLnPwvRaw'].size, )),
1665 (
'COMPRETRIEVEDLNPWVFLAG',
'i2',
1666 (parCat[
'compRetrievedLnPwvFlag'].size, )),
1667 (
'COMPRETRIEVEDTAUNIGHT',
'f8',
1668 (parCat[
'compRetrievedTauNight'].size, )),
1669 (
'COMPEPSILON',
'f8',
1670 (parCat[
'compEpsilon'].size, )),
1671 (
'COMPMEDDELTAAPER',
'f8',
1672 (parCat[
'compMedDeltaAper'].size, )),
1673 (
'COMPGLOBALEPSILON',
'f4',
1674 (parCat[
'compGlobalEpsilon'].size, )),
1675 (
'COMPEPSILONMAP',
'f4',
1676 (parCat[
'compEpsilonMap'].size, )),
1677 (
'COMPEPSILONNSTARMAP',
'i4',
1678 (parCat[
'compEpsilonNStarMap'].size, )),
1679 (
'COMPEPSILONCCDMAP',
'f4',
1680 (parCat[
'compEpsilonCcdMap'].size, )),
1681 (
'COMPEPSILONCCDNSTARMAP',
'i4',
1682 (parCat[
'compEpsilonCcdNStarMap'].size, ))])
1684 inParams[
'PARALPHA'][:] = parCat[
'parAlpha'][0, :]
1685 inParams[
'PARO3'][:] = parCat[
'parO3'][0, :]
1686 inParams[
'PARLNTAUINTERCEPT'][:] = parCat[
'parLnTauIntercept'][0, :]
1687 inParams[
'PARLNTAUSLOPE'][:] = parCat[
'parLnTauSlope'][0, :]
1688 inParams[
'PARLNPWVINTERCEPT'][:] = parCat[
'parLnPwvIntercept'][0, :]
1689 inParams[
'PARLNPWVSLOPE'][:] = parCat[
'parLnPwvSlope'][0, :]
1690 inParams[
'PARLNPWVQUADRATIC'][:] = parCat[
'parLnPwvQuadratic'][0, :]
1691 inParams[
'PARQESYSINTERCEPT'][:] = parCat[
'parQeSysIntercept'][0, :]
1692 inParams[
'COMPQESYSSLOPE'][:] = parCat[
'compQeSysSlope'][0, :]
1693 inParams[
'PARFILTEROFFSET'][:] = parCat[
'parFilterOffset'][0, :]
1694 inParams[
'PARFILTEROFFSETFITFLAG'][:] = parCat[
'parFilterOffsetFitFlag'][0, :]
1695 inParams[
'PARRETRIEVEDLNPWVSCALE'] = parCat[
'parRetrievedLnPwvScale']
1696 inParams[
'PARRETRIEVEDLNPWVOFFSET'] = parCat[
'parRetrievedLnPwvOffset']
1697 inParams[
'PARRETRIEVEDLNPWVNIGHTLYOFFSET'][:] = parCat[
'parRetrievedLnPwvNightlyOffset'][0, :]
1698 inParams[
'COMPABSTHROUGHPUT'][:] = parCat[
'compAbsThroughput'][0, :]
1699 inParams[
'COMPREFOFFSET'][:] = parCat[
'compRefOffset'][0, :]
1700 inParams[
'COMPREFSIGMA'][:] = parCat[
'compRefSigma'][0, :]
1701 inParams[
'COMPMIRRORCHROMATICITY'][:] = parCat[
'compMirrorChromaticity'][0, :]
1702 inParams[
'MIRRORCHROMATICITYPIVOT'][:] = parCat[
'mirrorChromaticityPivot'][0, :]
1703 inParams[
'COMPCCDCHROMATICITY'][:] = parCat[
'compCcdChromaticity'][0, :]
1704 inParams[
'COMPMEDIANSEDSLOPE'][:] = parCat[
'compMedianSedSlope'][0, :]
1705 inParams[
'COMPAPERCORRPIVOT'][:] = parCat[
'compAperCorrPivot'][0, :]
1706 inParams[
'COMPAPERCORRSLOPE'][:] = parCat[
'compAperCorrSlope'][0, :]
1707 inParams[
'COMPAPERCORRSLOPEERR'][:] = parCat[
'compAperCorrSlopeErr'][0, :]
1708 inParams[
'COMPAPERCORRRANGE'][:] = parCat[
'compAperCorrRange'][0, :]
1709 inParams[
'COMPMODELERREXPTIMEPIVOT'][:] = parCat[
'compModelErrExptimePivot'][0, :]
1710 inParams[
'COMPMODELERRFWHMPIVOT'][:] = parCat[
'compModelErrFwhmPivot'][0, :]
1711 inParams[
'COMPMODELERRSKYPIVOT'][:] = parCat[
'compModelErrSkyPivot'][0, :]
1712 inParams[
'COMPMODELERRPARS'][:] = parCat[
'compModelErrPars'][0, :]
1713 inParams[
'COMPEXPGRAY'][:] = parCat[
'compExpGray'][0, :]
1714 inParams[
'COMPVARGRAY'][:] = parCat[
'compVarGray'][0, :]
1715 inParams[
'COMPEXPDELTAMAGBKG'][:] = parCat[
'compExpDeltaMagBkg'][0, :]
1716 inParams[
'COMPNGOODSTARPEREXP'][:] = parCat[
'compNGoodStarPerExp'][0, :]
1717 inParams[
'COMPEXPREFOFFSET'][:] = parCat[
'compExpRefOffset'][0, :]
1718 inParams[
'COMPSIGFGCM'][:] = parCat[
'compSigFgcm'][0, :]
1719 inParams[
'COMPSIGMACAL'][:] = parCat[
'compSigmaCal'][0, :]
1720 inParams[
'COMPRETRIEVEDLNPWV'][:] = parCat[
'compRetrievedLnPwv'][0, :]
1721 inParams[
'COMPRETRIEVEDLNPWVRAW'][:] = parCat[
'compRetrievedLnPwvRaw'][0, :]
1722 inParams[
'COMPRETRIEVEDLNPWVFLAG'][:] = parCat[
'compRetrievedLnPwvFlag'][0, :]
1723 inParams[
'COMPRETRIEVEDTAUNIGHT'][:] = parCat[
'compRetrievedTauNight'][0, :]
1724 inParams[
'COMPEPSILON'][:] = parCat[
'compEpsilon'][0, :]
1725 inParams[
'COMPMEDDELTAAPER'][:] = parCat[
'compMedDeltaAper'][0, :]
1726 inParams[
'COMPGLOBALEPSILON'][:] = parCat[
'compGlobalEpsilon'][0, :]
1727 inParams[
'COMPEPSILONMAP'][:] = parCat[
'compEpsilonMap'][0, :]
1728 inParams[
'COMPEPSILONNSTARMAP'][:] = parCat[
'compEpsilonNStarMap'][0, :]
1729 inParams[
'COMPEPSILONCCDMAP'][:] = parCat[
'compEpsilonCcdMap'][0, :]
1730 inParams[
'COMPEPSILONCCDNSTARMAP'][:] = parCat[
'compEpsilonCcdNStarMap'][0, :]
1732 inSuperStar = np.zeros(parCat[
'superstarSize'][0, :], dtype=
'f8')
1733 inSuperStar[:, :, :, :] = parCat[
'superstar'][0, :].reshape(inSuperStar.shape)
1735 return (inParInfo, inParams, inSuperStar)
1737 def _makeFgcmOutputDatasets(self, fgcmFitCycle):
1739 Persist FGCM datasets through the butler.
1743 fgcmFitCycle: `lsst.fgcm.FgcmFitCycle`
1744 Fgcm Fit cycle object
1746 fgcmDatasetDict = {}
1749 parInfo, pars = fgcmFitCycle.fgcmPars.parsToArrays()
1754 lutFilterNameString = comma.join([n.decode(
'utf-8')
1755 for n
in parInfo[
'LUTFILTERNAMES'][0]])
1756 fitBandString = comma.join([n.decode(
'utf-8')
1757 for n
in parInfo[
'FITBANDS'][0]])
1759 parSchema = self._makeParSchema(parInfo, pars, fgcmFitCycle.fgcmPars.parSuperStarFlat,
1760 lutFilterNameString, fitBandString)
1761 parCat = self._makeParCatalog(parSchema, parInfo, pars,
1762 fgcmFitCycle.fgcmPars.parSuperStarFlat,
1763 lutFilterNameString, fitBandString)
1765 fgcmDatasetDict[
'fgcmFitParameters'] = parCat
1770 flagStarSchema = self._makeFlagStarSchema()
1771 flagStarStruct = fgcmFitCycle.fgcmStars.getFlagStarIndices()
1772 flagStarCat = self._makeFlagStarCat(flagStarSchema, flagStarStruct)
1774 fgcmDatasetDict[
'fgcmFlaggedStars'] = flagStarCat
1777 if self.outputZeropoints:
1778 superStarChebSize = fgcmFitCycle.fgcmZpts.zpStruct[
'FGCM_FZPT_SSTAR_CHEB'].shape[1]
1779 zptChebSize = fgcmFitCycle.fgcmZpts.zpStruct[
'FGCM_FZPT_CHEB'].shape[1]
1781 zptSchema = makeZptSchema(superStarChebSize, zptChebSize)
1782 zptCat = makeZptCat(zptSchema, fgcmFitCycle.fgcmZpts.zpStruct)
1784 fgcmDatasetDict[
'fgcmZeropoints'] = zptCat
1788 atmSchema = makeAtmSchema()
1789 atmCat = makeAtmCat(atmSchema, fgcmFitCycle.fgcmZpts.atmStruct)
1791 fgcmDatasetDict[
'fgcmAtmosphereParameters'] = atmCat
1794 if self.outputStandards:
1795 stdStruct, goodBands = fgcmFitCycle.fgcmStars.retrieveStdStarCatalog(fgcmFitCycle.fgcmPars)
1796 stdSchema = makeStdSchema(len(goodBands))
1797 stdCat = makeStdCat(stdSchema, stdStruct, goodBands)
1799 fgcmDatasetDict[
'fgcmStandardStars'] = stdCat
1801 return fgcmDatasetDict
1803 def _makeParSchema(self, parInfo, pars, parSuperStarFlat,
1804 lutFilterNameString, fitBandString):
1806 Make the parameter persistence schema
1810 parInfo: `numpy.ndarray`
1811 Parameter information returned by fgcm
1812 pars: `numpy.ndarray`
1813 Parameter values returned by fgcm
1814 parSuperStarFlat: `numpy.array`
1815 Superstar flat values returned by fgcm
1816 lutFilterNameString: `str`
1817 Combined string of all the lutFilterNames
1818 fitBandString: `str`
1819 Combined string of all the fitBands
1823 parSchema: `afwTable.schema`
1829 parSchema.addField(
'nCcd', type=np.int32, doc=
'Number of CCDs')
1830 parSchema.addField(
'lutFilterNames', type=str, doc=
'LUT Filter names in parameter file',
1831 size=len(lutFilterNameString))
1832 parSchema.addField(
'fitBands', type=str, doc=
'Bands that were fit',
1833 size=len(fitBandString))
1834 parSchema.addField(
'lnTauUnit', type=np.float64, doc=
'Step units for ln(AOD)')
1835 parSchema.addField(
'lnTauSlopeUnit', type=np.float64,
1836 doc=
'Step units for ln(AOD) slope')
1837 parSchema.addField(
'alphaUnit', type=np.float64, doc=
'Step units for alpha')
1838 parSchema.addField(
'lnPwvUnit', type=np.float64, doc=
'Step units for ln(pwv)')
1839 parSchema.addField(
'lnPwvSlopeUnit', type=np.float64,
1840 doc=
'Step units for ln(pwv) slope')
1841 parSchema.addField(
'lnPwvQuadraticUnit', type=np.float64,
1842 doc=
'Step units for ln(pwv) quadratic term')
1843 parSchema.addField(
'lnPwvGlobalUnit', type=np.float64,
1844 doc=
'Step units for global ln(pwv) parameters')
1845 parSchema.addField(
'o3Unit', type=np.float64, doc=
'Step units for O3')
1846 parSchema.addField(
'qeSysUnit', type=np.float64, doc=
'Step units for mirror gray')
1847 parSchema.addField(
'filterOffsetUnit', type=np.float64, doc=
'Step units for filter offset')
1848 parSchema.addField(
'hasExternalPwv', type=np.int32, doc=
'Parameters fit using external pwv')
1849 parSchema.addField(
'hasExternalTau', type=np.int32, doc=
'Parameters fit using external tau')
1852 parSchema.addField(
'parAlpha', type=
'ArrayD', doc=
'Alpha parameter vector',
1853 size=pars[
'PARALPHA'].size)
1854 parSchema.addField(
'parO3', type=
'ArrayD', doc=
'O3 parameter vector',
1855 size=pars[
'PARO3'].size)
1856 parSchema.addField(
'parLnTauIntercept', type=
'ArrayD',
1857 doc=
'ln(Tau) intercept parameter vector',
1858 size=pars[
'PARLNTAUINTERCEPT'].size)
1859 parSchema.addField(
'parLnTauSlope', type=
'ArrayD',
1860 doc=
'ln(Tau) slope parameter vector',
1861 size=pars[
'PARLNTAUSLOPE'].size)
1862 parSchema.addField(
'parLnPwvIntercept', type=
'ArrayD', doc=
'ln(pwv) intercept parameter vector',
1863 size=pars[
'PARLNPWVINTERCEPT'].size)
1864 parSchema.addField(
'parLnPwvSlope', type=
'ArrayD', doc=
'ln(pwv) slope parameter vector',
1865 size=pars[
'PARLNPWVSLOPE'].size)
1866 parSchema.addField(
'parLnPwvQuadratic', type=
'ArrayD', doc=
'ln(pwv) quadratic parameter vector',
1867 size=pars[
'PARLNPWVQUADRATIC'].size)
1868 parSchema.addField(
'parQeSysIntercept', type=
'ArrayD', doc=
'Mirror gray intercept parameter vector',
1869 size=pars[
'PARQESYSINTERCEPT'].size)
1870 parSchema.addField(
'compQeSysSlope', type=
'ArrayD', doc=
'Mirror gray slope parameter vector',
1871 size=pars[0][
'COMPQESYSSLOPE'].size)
1872 parSchema.addField(
'parFilterOffset', type=
'ArrayD', doc=
'Filter offset parameter vector',
1873 size=pars[
'PARFILTEROFFSET'].size)
1874 parSchema.addField(
'parFilterOffsetFitFlag', type=
'ArrayI', doc=
'Filter offset parameter fit flag',
1875 size=pars[
'PARFILTEROFFSETFITFLAG'].size)
1876 parSchema.addField(
'parRetrievedLnPwvScale', type=np.float64,
1877 doc=
'Global scale for retrieved ln(pwv)')
1878 parSchema.addField(
'parRetrievedLnPwvOffset', type=np.float64,
1879 doc=
'Global offset for retrieved ln(pwv)')
1880 parSchema.addField(
'parRetrievedLnPwvNightlyOffset', type=
'ArrayD',
1881 doc=
'Nightly offset for retrieved ln(pwv)',
1882 size=pars[
'PARRETRIEVEDLNPWVNIGHTLYOFFSET'].size)
1883 parSchema.addField(
'compAbsThroughput', type=
'ArrayD',
1884 doc=
'Absolute throughput (relative to transmission curves)',
1885 size=pars[
'COMPABSTHROUGHPUT'].size)
1886 parSchema.addField(
'compRefOffset', type=
'ArrayD',
1887 doc=
'Offset between reference stars and calibrated stars',
1888 size=pars[
'COMPREFOFFSET'].size)
1889 parSchema.addField(
'compRefSigma', type=
'ArrayD',
1890 doc=
'Width of reference star/calibrated star distribution',
1891 size=pars[
'COMPREFSIGMA'].size)
1892 parSchema.addField(
'compMirrorChromaticity', type=
'ArrayD',
1893 doc=
'Computed mirror chromaticity terms',
1894 size=pars[
'COMPMIRRORCHROMATICITY'].size)
1895 parSchema.addField(
'mirrorChromaticityPivot', type=
'ArrayD',
1896 doc=
'Mirror chromaticity pivot mjd',
1897 size=pars[
'MIRRORCHROMATICITYPIVOT'].size)
1898 parSchema.addField(
'compCcdChromaticity', type=
'ArrayD',
1899 doc=
'Computed CCD chromaticity terms',
1900 size=pars[
'COMPCCDCHROMATICITY'].size)
1901 parSchema.addField(
'compMedianSedSlope', type=
'ArrayD',
1902 doc=
'Computed median SED slope (per band)',
1903 size=pars[
'COMPMEDIANSEDSLOPE'].size)
1904 parSchema.addField(
'compAperCorrPivot', type=
'ArrayD', doc=
'Aperture correction pivot',
1905 size=pars[
'COMPAPERCORRPIVOT'].size)
1906 parSchema.addField(
'compAperCorrSlope', type=
'ArrayD', doc=
'Aperture correction slope',
1907 size=pars[
'COMPAPERCORRSLOPE'].size)
1908 parSchema.addField(
'compAperCorrSlopeErr', type=
'ArrayD', doc=
'Aperture correction slope error',
1909 size=pars[
'COMPAPERCORRSLOPEERR'].size)
1910 parSchema.addField(
'compAperCorrRange', type=
'ArrayD', doc=
'Aperture correction range',
1911 size=pars[
'COMPAPERCORRRANGE'].size)
1912 parSchema.addField(
'compModelErrExptimePivot', type=
'ArrayD', doc=
'Model error exptime pivot',
1913 size=pars[
'COMPMODELERREXPTIMEPIVOT'].size)
1914 parSchema.addField(
'compModelErrFwhmPivot', type=
'ArrayD', doc=
'Model error fwhm pivot',
1915 size=pars[
'COMPMODELERRFWHMPIVOT'].size)
1916 parSchema.addField(
'compModelErrSkyPivot', type=
'ArrayD', doc=
'Model error sky pivot',
1917 size=pars[
'COMPMODELERRSKYPIVOT'].size)
1918 parSchema.addField(
'compModelErrPars', type=
'ArrayD', doc=
'Model error parameters',
1919 size=pars[
'COMPMODELERRPARS'].size)
1920 parSchema.addField(
'compExpGray', type=
'ArrayD', doc=
'Computed exposure gray',
1921 size=pars[
'COMPEXPGRAY'].size)
1922 parSchema.addField(
'compVarGray', type=
'ArrayD', doc=
'Computed exposure variance',
1923 size=pars[
'COMPVARGRAY'].size)
1924 parSchema.addField(
'compExpDeltaMagBkg', type=
'ArrayD',
1925 doc=
'Computed exposure offset due to background',
1926 size=pars[
'COMPEXPDELTAMAGBKG'].size)
1927 parSchema.addField(
'compNGoodStarPerExp', type=
'ArrayI',
1928 doc=
'Computed number of good stars per exposure',
1929 size=pars[
'COMPNGOODSTARPEREXP'].size)
1930 parSchema.addField(
'compExpRefOffset', type=
'ArrayD',
1931 doc=
'Computed per-visit median offset between standard stars and ref stars.',
1932 size=pars[
'COMPEXPREFOFFSET'].size)
1933 parSchema.addField(
'compSigFgcm', type=
'ArrayD', doc=
'Computed sigma_fgcm (intrinsic repeatability)',
1934 size=pars[
'COMPSIGFGCM'].size)
1935 parSchema.addField(
'compSigmaCal', type=
'ArrayD', doc=
'Computed sigma_cal (systematic error floor)',
1936 size=pars[
'COMPSIGMACAL'].size)
1937 parSchema.addField(
'compRetrievedLnPwv', type=
'ArrayD', doc=
'Retrieved ln(pwv) (smoothed)',
1938 size=pars[
'COMPRETRIEVEDLNPWV'].size)
1939 parSchema.addField(
'compRetrievedLnPwvRaw', type=
'ArrayD', doc=
'Retrieved ln(pwv) (raw)',
1940 size=pars[
'COMPRETRIEVEDLNPWVRAW'].size)
1941 parSchema.addField(
'compRetrievedLnPwvFlag', type=
'ArrayI', doc=
'Retrieved ln(pwv) Flag',
1942 size=pars[
'COMPRETRIEVEDLNPWVFLAG'].size)
1943 parSchema.addField(
'compRetrievedTauNight', type=
'ArrayD', doc=
'Retrieved tau (per night)',
1944 size=pars[
'COMPRETRIEVEDTAUNIGHT'].size)
1945 parSchema.addField(
'compEpsilon', type=
'ArrayD',
1946 doc=
'Computed epsilon background offset per visit (nJy/arcsec2)',
1947 size=pars[
'COMPEPSILON'].size)
1948 parSchema.addField(
'compMedDeltaAper', type=
'ArrayD',
1949 doc=
'Median delta mag aper per visit',
1950 size=pars[
'COMPMEDDELTAAPER'].size)
1951 parSchema.addField(
'compGlobalEpsilon', type=
'ArrayD',
1952 doc=
'Computed epsilon bkg offset (global) (nJy/arcsec2)',
1953 size=pars[
'COMPGLOBALEPSILON'].size)
1954 parSchema.addField(
'compEpsilonMap', type=
'ArrayD',
1955 doc=
'Computed epsilon maps (nJy/arcsec2)',
1956 size=pars[
'COMPEPSILONMAP'].size)
1957 parSchema.addField(
'compEpsilonNStarMap', type=
'ArrayI',
1958 doc=
'Number of stars per pixel in computed epsilon maps',
1959 size=pars[
'COMPEPSILONNSTARMAP'].size)
1960 parSchema.addField(
'compEpsilonCcdMap', type=
'ArrayD',
1961 doc=
'Computed epsilon ccd maps (nJy/arcsec2)',
1962 size=pars[
'COMPEPSILONCCDMAP'].size)
1963 parSchema.addField(
'compEpsilonCcdNStarMap', type=
'ArrayI',
1964 doc=
'Number of stars per ccd bin in epsilon ccd maps',
1965 size=pars[
'COMPEPSILONCCDNSTARMAP'].size)
1967 parSchema.addField(
'superstarSize', type=
'ArrayI', doc=
'Superstar matrix size',
1969 parSchema.addField(
'superstar', type=
'ArrayD', doc=
'Superstar matrix (flattened)',
1970 size=parSuperStarFlat.size)
1974 def _makeParCatalog(self, parSchema, parInfo, pars, parSuperStarFlat,
1975 lutFilterNameString, fitBandString):
1977 Make the FGCM parameter catalog for persistence
1981 parSchema: `lsst.afw.table.Schema`
1982 Parameter catalog schema
1983 pars: `numpy.ndarray`
1984 FGCM parameters to put into parCat
1985 parSuperStarFlat: `numpy.array`
1986 FGCM superstar flat array to put into parCat
1987 lutFilterNameString: `str`
1988 Combined string of all the lutFilterNames
1989 fitBandString: `str`
1990 Combined string of all the fitBands
1994 parCat: `afwTable.BasicCatalog`
1995 Atmosphere and instrumental model parameter catalog for persistence
2003 rec = parCat.addNew()
2006 rec[
'nCcd'] = parInfo[
'NCCD'][0]
2007 rec[
'lutFilterNames'] = lutFilterNameString
2008 rec[
'fitBands'] = fitBandString
2010 rec[
'hasExternalPwv'] = 0
2011 rec[
'hasExternalTau'] = 0
2015 scalarNames = [
'parRetrievedLnPwvScale',
'parRetrievedLnPwvOffset']
2017 arrNames = [
'parAlpha',
'parO3',
'parLnTauIntercept',
'parLnTauSlope',
2018 'parLnPwvIntercept',
'parLnPwvSlope',
'parLnPwvQuadratic',
2019 'parQeSysIntercept',
'compQeSysSlope',
2020 'parRetrievedLnPwvNightlyOffset',
'compAperCorrPivot',
2021 'parFilterOffset',
'parFilterOffsetFitFlag',
2022 'compAbsThroughput',
'compRefOffset',
'compRefSigma',
2023 'compMirrorChromaticity',
'mirrorChromaticityPivot',
'compCcdChromaticity',
2024 'compAperCorrSlope',
'compAperCorrSlopeErr',
'compAperCorrRange',
2025 'compModelErrExptimePivot',
'compModelErrFwhmPivot',
2026 'compModelErrSkyPivot',
'compModelErrPars',
2027 'compExpGray',
'compVarGray',
'compNGoodStarPerExp',
'compSigFgcm',
2028 'compSigmaCal',
'compExpDeltaMagBkg',
'compMedianSedSlope',
2029 'compRetrievedLnPwv',
'compRetrievedLnPwvRaw',
'compRetrievedLnPwvFlag',
2030 'compRetrievedTauNight',
'compEpsilon',
'compMedDeltaAper',
2031 'compGlobalEpsilon',
'compEpsilonMap',
'compEpsilonNStarMap',
2032 'compEpsilonCcdMap',
'compEpsilonCcdNStarMap',
'compExpRefOffset']
2034 for scalarName
in scalarNames:
2035 rec[scalarName] = pars[scalarName.upper()][0]
2037 for arrName
in arrNames:
2038 rec[arrName][:] = np.atleast_1d(pars[0][arrName.upper()])[:]
2041 rec[
'superstarSize'][:] = parSuperStarFlat.shape
2042 rec[
'superstar'][:] = parSuperStarFlat.ravel()
2046 def _makeFlagStarSchema(self):
2048 Make the flagged-stars schema
2052 flagStarSchema: `lsst.afw.table.Schema`
2057 flagStarSchema.addField(
'objId', type=np.int32, doc=
'FGCM object id')
2058 flagStarSchema.addField(
'objFlag', type=np.int32, doc=
'FGCM object flag')
2060 return flagStarSchema
2062 def _makeFlagStarCat(self, flagStarSchema, flagStarStruct):
2064 Make the flagged star catalog for persistence
2068 flagStarSchema: `lsst.afw.table.Schema`
2070 flagStarStruct: `numpy.ndarray`
2071 Flagged star structure from fgcm
2075 flagStarCat: `lsst.afw.table.BaseCatalog`
2076 Flagged star catalog for persistence
2080 flagStarCat.resize(flagStarStruct.size)
2082 flagStarCat[
'objId'][:] = flagStarStruct[
'OBJID']
2083 flagStarCat[
'objFlag'][:] = flagStarStruct[
'OBJFLAG']
Defines the fields and offsets for a table.