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 defaultCameraOrientation = pexConfig.Field(
692 doc=
"Default camera orientation for QA plots.",
696 brightObsGrayMax = pexConfig.Field(
697 doc=
"Maximum gray extinction to be considered bright observation",
701 minStarPerCcd = pexConfig.Field(
702 doc=(
"Minimum number of good stars per CCD to be used in calibration fit. "
703 "CCDs with fewer stars will have their calibration estimated from other "
704 "CCDs in the same visit, with zeropoint error increased accordingly."),
708 minCcdPerExp = pexConfig.Field(
709 doc=(
"Minimum number of good CCDs per exposure/visit to be used in calibration fit. "
710 "Visits with fewer good CCDs will have CCD zeropoints estimated where possible."),
714 maxCcdGrayErr = pexConfig.Field(
715 doc=
"Maximum error on CCD gray offset to be considered photometric",
719 minStarPerExp = pexConfig.Field(
720 doc=(
"Minimum number of good stars per exposure/visit to be used in calibration fit. "
721 "Visits with fewer good stars will have CCD zeropoints estimated where possible."),
725 minExpPerNight = pexConfig.Field(
726 doc=
"Minimum number of good exposures/visits to consider a partly photometric night",
730 expGrayInitialCut = pexConfig.Field(
731 doc=(
"Maximum exposure/visit gray value for initial selection of possible photometric "
736 expGrayPhotometricCutDict = pexConfig.DictField(
737 doc=(
"Per-band specification on maximum (negative) achromatic exposure residual "
738 "('gray term') for a visit to be considered photometric. Must have one "
739 "entry per band. Broad-band filters should be -0.05."),
744 expGrayHighCutDict = pexConfig.DictField(
745 doc=(
"Per-band specification on maximum (positive) achromatic exposure residual "
746 "('gray term') for a visit to be considered photometric. Must have one "
747 "entry per band. Broad-band filters should be 0.2."),
752 expGrayRecoverCut = pexConfig.Field(
753 doc=(
"Maximum (negative) exposure gray to be able to recover bad ccds via interpolation. "
754 "Visits with more gray extinction will only get CCD zeropoints if there are "
755 "sufficient star observations (minStarPerCcd) on that CCD."),
759 expVarGrayPhotometricCutDict = pexConfig.DictField(
760 doc=(
"Per-band specification on maximum exposure variance to be considered possibly "
761 "photometric. Must have one entry per band. Broad-band filters should be "
767 expGrayErrRecoverCut = pexConfig.Field(
768 doc=(
"Maximum exposure gray error to be able to recover bad ccds via interpolation. "
769 "Visits with more gray variance will only get CCD zeropoints if there are "
770 "sufficient star observations (minStarPerCcd) on that CCD."),
774 aperCorrFitNBins = pexConfig.Field(
775 doc=(
"Number of aperture bins used in aperture correction fit. When set to 0"
776 "no fit will be performed, and the config.aperCorrInputSlopes will be "
777 "used if available."),
781 aperCorrInputSlopeDict = pexConfig.DictField(
782 doc=(
"Per-band specification of aperture correction input slope parameters. These "
783 "are used on the first fit iteration, and aperture correction parameters will "
784 "be updated from the data if config.aperCorrFitNBins > 0. It is recommended "
785 "to set this when there is insufficient data to fit the parameters (e.g. "
791 sedboundaryterms = pexConfig.ConfigField(
792 doc=
"Mapping from bands to SED boundary term names used is sedterms.",
793 dtype=SedboundarytermDict,
795 sedterms = pexConfig.ConfigField(
796 doc=
"Mapping from terms to bands for fgcm linear SED approximations.",
799 sigFgcmMaxErr = pexConfig.Field(
800 doc=
"Maximum mag error for fitting sigma_FGCM",
804 sigFgcmMaxEGrayDict = pexConfig.DictField(
805 doc=(
"Per-band specification for maximum (absolute) achromatic residual (gray value) "
806 "for observations in sigma_fgcm (raw repeatability). Broad-band filters "
812 ccdGrayMaxStarErr = pexConfig.Field(
813 doc=(
"Maximum error on a star observation to use in ccd gray (achromatic residual) "
818 approxThroughputDict = pexConfig.DictField(
819 doc=(
"Per-band specification of the approximate overall throughput at the start of "
820 "calibration observations. Must have one entry per band. Typically should "
826 sigmaCalRange = pexConfig.ListField(
827 doc=
"Allowed range for systematic error floor estimation",
829 default=(0.001, 0.003),
831 sigmaCalFitPercentile = pexConfig.ListField(
832 doc=
"Magnitude percentile range to fit systematic error floor",
834 default=(0.05, 0.15),
836 sigmaCalPlotPercentile = pexConfig.ListField(
837 doc=
"Magnitude percentile range to plot systematic error floor",
839 default=(0.05, 0.95),
841 sigma0Phot = pexConfig.Field(
842 doc=
"Systematic error floor for all zeropoints",
846 mapLongitudeRef = pexConfig.Field(
847 doc=
"Reference longitude for plotting maps",
851 mapNSide = pexConfig.Field(
852 doc=
"Healpix nside for plotting maps",
856 outfileBase = pexConfig.Field(
857 doc=
"Filename start for plot output files",
861 starColorCuts = pexConfig.ListField(
862 doc=(
"Encoded star-color cuts (using calibration star colors). "
863 "This is a list with each entry a string of the format "
864 "``band1,band2,low,high`` such that only stars of color "
865 "low < band1 - band2 < high will be used for calibration."),
867 default=(
"NO_DATA",),
869 refStarColorCuts = pexConfig.ListField(
870 doc=(
"Encoded star color cuts specifically to apply to reference stars. "
871 "This is a list with each entry a string of the format "
872 "``band1,band2,low,high`` such that only stars of color "
873 "low < band1 - band2 < high will be used as reference stars."),
875 default=(
"NO_DATA",),
877 colorSplitBands = pexConfig.ListField(
878 doc=
"Band names to use to split stars by color. Must have 2 entries.",
883 modelMagErrors = pexConfig.Field(
884 doc=
"Should FGCM model the magnitude errors from sky/fwhm? (False means trust inputs)",
888 useQuadraticPwv = pexConfig.Field(
889 doc=
"Model PWV with a quadratic term for variation through the night?",
893 instrumentParsPerBand = pexConfig.Field(
894 doc=(
"Model instrumental parameters per band? "
895 "Otherwise, instrumental parameters (QE changes with time) are "
896 "shared among all bands."),
900 instrumentSlopeMinDeltaT = pexConfig.Field(
901 doc=(
"Minimum time change (in days) between observations to use in constraining "
902 "instrument slope."),
906 fitMirrorChromaticity = pexConfig.Field(
907 doc=
"Fit (intraband) mirror chromatic term?",
911 fitCcdChromaticityDict = pexConfig.DictField(
912 doc=
"Specification on whether to compute first-order quantum efficiency (QE) "
913 "adjustments. Key is band, and value will be True or False. Any band "
914 "not explicitly specified will default to False.",
919 coatingMjds = pexConfig.ListField(
920 doc=
"Mirror coating dates in MJD",
924 outputStandardsBeforeFinalCycle = pexConfig.Field(
925 doc=
"Output standard stars prior to final cycle? Used in debugging.",
929 outputZeropointsBeforeFinalCycle = pexConfig.Field(
930 doc=
"Output standard stars prior to final cycle? Used in debugging.",
934 useRepeatabilityForExpGrayCutsDict = pexConfig.DictField(
935 doc=(
"Per-band specification on whether to use star repeatability (instead of exposures) "
936 "for computing photometric cuts. Recommended for tract mode or bands with few visits."),
941 autoPhotometricCutNSig = pexConfig.Field(
942 doc=(
"Number of sigma for automatic computation of (low) photometric cut. "
943 "Cut is based on exposure gray width (per band), unless "
944 "useRepeatabilityForExpGrayCuts is set, in which case the star "
945 "repeatability is used (also per band)."),
949 autoHighCutNSig = pexConfig.Field(
950 doc=(
"Number of sigma for automatic computation of (high) outlier cut. "
951 "Cut is based on exposure gray width (per band), unless "
952 "useRepeatabilityForExpGrayCuts is set, in which case the star "
953 "repeatability is used (also per band)."),
957 quietMode = pexConfig.Field(
958 doc=
"Be less verbose with logging.",
962 doPlots = pexConfig.Field(
963 doc=
"Make fgcm QA plots.",
967 doPlotsBeforeFinalCycles = pexConfig.Field(
968 doc=
"Make fgcm QA plots before the final two fit cycles? This only applies in"
969 "multi-cycle mode, and if doPlots is True. These are typically the most"
970 "important QA plots to inspect.",
974 randomSeed = pexConfig.Field(
975 doc=
"Random seed for fgcm for consistency in tests.",
980 deltaAperFitMinNgoodObs = pexConfig.Field(
981 doc=
"Minimum number of good observations to use mean delta-aper values in fits.",
985 deltaAperFitPerCcdNx = pexConfig.Field(
986 doc=(
"Number of x bins per ccd when computing delta-aper background offsets. "
987 "Only used when ``doComputeDeltaAperPerCcd`` is True."),
991 deltaAperFitPerCcdNy = pexConfig.Field(
992 doc=(
"Number of y bins per ccd when computing delta-aper background offsets. "
993 "Only used when ``doComputeDeltaAperPerCcd`` is True."),
997 deltaAperFitSpatialNside = pexConfig.Field(
998 doc=
"Healpix nside to compute spatial delta-aper background offset maps.",
1002 deltaAperInnerRadiusArcsec = pexConfig.Field(
1003 doc=(
"Inner radius used to compute deltaMagAper (arcseconds). "
1004 "Must be positive and less than ``deltaAperOuterRadiusArcsec`` if "
1005 "any of ``doComputeDeltaAperPerVisit``, ``doComputeDeltaAperPerStar``, "
1006 "``doComputeDeltaAperMap``, ``doComputeDeltaAperPerCcd`` are set."),
1010 deltaAperOuterRadiusArcsec = pexConfig.Field(
1011 doc=(
"Outer radius used to compute deltaMagAper (arcseconds). "
1012 "Must be positive and greater than ``deltaAperInnerRadiusArcsec`` if "
1013 "any of ``doComputeDeltaAperPerVisit``, ``doComputeDeltaAperPerStar``, "
1014 "``doComputeDeltaAperMap``, ``doComputeDeltaAperPerCcd`` are set."),
1018 doComputeDeltaAperPerVisit = pexConfig.Field(
1019 doc=(
"Do the computation of delta-aper background offsets per visit? "
1020 "Note: this option can be very slow when there are many visits."),
1024 doComputeDeltaAperPerStar = pexConfig.Field(
1025 doc=
"Do the computation of delta-aper mean values per star?",
1029 doComputeDeltaAperMap = pexConfig.Field(
1030 doc=(
"Do the computation of delta-aper spatial maps? "
1031 "This is only used if ``doComputeDeltaAperPerStar`` is True,"),
1035 doComputeDeltaAperPerCcd = pexConfig.Field(
1036 doc=
"Do the computation of per-ccd delta-aper background offsets?",
1044 if self.connections.previousCycleNumber != str(self.cycleNumber - 1):
1045 msg =
"cycleNumber in template must be connections.previousCycleNumber + 1"
1046 raise RuntimeError(msg)
1047 if self.connections.cycleNumber != str(self.cycleNumber):
1048 msg =
"cycleNumber in template must be equal to connections.cycleNumber"
1049 raise RuntimeError(msg)
1051 for band
in self.fitBands:
1052 if band
not in self.bands:
1053 msg =
'fitBand %s not in bands' % (band)
1054 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.fitBands, self, msg)
1055 for band
in self.requiredBands:
1056 if band
not in self.bands:
1057 msg =
'requiredBand %s not in bands' % (band)
1058 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.requiredBands, self, msg)
1059 for band
in self.colorSplitBands:
1060 if band
not in self.bands:
1061 msg =
'colorSplitBand %s not in bands' % (band)
1062 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.colorSplitBands, self, msg)
1063 for band
in self.bands:
1064 if band
not in self.superStarSubCcdDict:
1065 msg =
'band %s not in superStarSubCcdDict' % (band)
1066 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.superStarSubCcdDict,
1068 if band
not in self.ccdGraySubCcdDict:
1069 msg =
'band %s not in ccdGraySubCcdDict' % (band)
1070 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.ccdGraySubCcdDict,
1072 if band
not in self.expGrayPhotometricCutDict:
1073 msg =
'band %s not in expGrayPhotometricCutDict' % (band)
1074 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.expGrayPhotometricCutDict,
1076 if band
not in self.expGrayHighCutDict:
1077 msg =
'band %s not in expGrayHighCutDict' % (band)
1078 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.expGrayHighCutDict,
1080 if band
not in self.expVarGrayPhotometricCutDict:
1081 msg =
'band %s not in expVarGrayPhotometricCutDict' % (band)
1082 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.expVarGrayPhotometricCutDict,
1084 if band
not in self.sigFgcmMaxEGrayDict:
1085 msg =
'band %s not in sigFgcmMaxEGrayDict' % (band)
1086 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.sigFgcmMaxEGrayDict,
1088 if band
not in self.approxThroughputDict:
1089 msg =
'band %s not in approxThroughputDict' % (band)
1090 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.approxThroughputDict,
1092 if band
not in self.useRepeatabilityForExpGrayCutsDict:
1093 msg =
'band %s not in useRepeatabilityForExpGrayCutsDict' % (band)
1094 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.useRepeatabilityForExpGrayCutsDict,
1097 if self.doComputeDeltaAperPerVisit
or self.doComputeDeltaAperMap \
1098 or self.doComputeDeltaAperPerCcd:
1099 if self.deltaAperInnerRadiusArcsec <= 0.0:
1100 msg =
'deltaAperInnerRadiusArcsec must be positive if deltaAper computations are turned on.'
1101 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.deltaAperInnerRadiusArcsec,
1103 if self.deltaAperOuterRadiusArcsec <= 0.0:
1104 msg =
'deltaAperOuterRadiusArcsec must be positive if deltaAper computations are turned on.'
1105 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.deltaAperOuterRadiusArcsec,
1107 if self.deltaAperOuterRadiusArcsec <= self.deltaAperInnerRadiusArcsec:
1108 msg = (
'deltaAperOuterRadiusArcsec must be greater than deltaAperInnerRadiusArcsec if '
1109 'deltaAper computations are turned on.')
1110 raise pexConfig.FieldValidationError(FgcmFitCycleConfig.deltaAperOuterRadiusArcsec,
1114class FgcmFitCycleTask(pipeBase.PipelineTask):
1116 Run Single fit cycle for FGCM global calibration
1119 ConfigClass = FgcmFitCycleConfig
1120 _DefaultName =
"fgcmFitCycle"
1122 def __init__(self, initInputs=None, **kwargs):
1123 super().__init__(**kwargs)
1125 def runQuantum(self, butlerQC, inputRefs, outputRefs):
1126 camera = butlerQC.get(inputRefs.camera)
1128 nCore = butlerQC.resources.num_cores
1132 handleDict[
'fgcmLookUpTable'] = butlerQC.get(inputRefs.fgcmLookUpTable)
1133 handleDict[
'fgcmVisitCatalog'] = butlerQC.get(inputRefs.fgcmVisitCatalog)
1135 if self.config.useParquetCatalogFormat:
1136 handleDict[
'fgcmStarObservations'] = butlerQC.get(inputRefs.fgcmStarObservationsParquet)
1137 handleDict[
'fgcmStarIds'] = butlerQC.get(inputRefs.fgcmStarIdsParquet)
1138 if self.config.doReferenceCalibration:
1139 handleDict[
'fgcmReferenceStars'] = butlerQC.get(inputRefs.fgcmReferenceStarsParquet)
1141 handleDict[
'fgcmStarObservations'] = butlerQC.get(inputRefs.fgcmStarObservations)
1142 handleDict[
'fgcmStarIds'] = butlerQC.get(inputRefs.fgcmStarIds)
1143 handleDict[
'fgcmStarIndices'] = butlerQC.get(inputRefs.fgcmStarIndices)
1144 if self.config.doReferenceCalibration:
1145 handleDict[
'fgcmReferenceStars'] = butlerQC.get(inputRefs.fgcmReferenceStars)
1146 if self.config.cycleNumber > 0:
1147 handleDict[
'fgcmFlaggedStars'] = butlerQC.get(inputRefs.fgcmFlaggedStarsInput)
1148 handleDict[
'fgcmFitParameters'] = butlerQC.get(inputRefs.fgcmFitParametersInput)
1150 fgcmDatasetDict =
None
1151 if self.config.doMultipleCycles:
1153 config = copy.copy(self.config)
1154 config.update(cycleNumber=0)
1155 for cycle
in range(self.config.multipleCyclesFinalCycleNumber + 1):
1156 if cycle == self.config.multipleCyclesFinalCycleNumber:
1157 config.update(isFinalCycle=
True)
1160 handleDict[
'fgcmFlaggedStars'] = fgcmDatasetDict[
'fgcmFlaggedStars']
1161 handleDict[
'fgcmFitParameters'] = fgcmDatasetDict[
'fgcmFitParameters']
1166 for outputRefName
in outputRefs.keys():
1167 if outputRefName.endswith(
"Plot")
and f
"Cycle{cycle}" in outputRefName:
1168 ref = getattr(outputRefs, outputRefName)
1169 plotHandleDict[outputRefName] = ref
1171 fgcmDatasetDict, config = self._fgcmFitCycle(
1175 plotHandleDict=plotHandleDict,
1179 butlerQC.put(fgcmDatasetDict[
'fgcmFitParameters'],
1180 getattr(outputRefs, f
'fgcm_Cycle{cycle}_FitParameters'))
1181 butlerQC.put(fgcmDatasetDict[
'fgcmFlaggedStars'],
1182 getattr(outputRefs, f
'fgcm_Cycle{cycle}_FlaggedStars'))
1183 butlerQC.put(config,
1184 getattr(outputRefs, f
'fgcm_Cycle{cycle}_OutputConfig'))
1185 if self.outputZeropoints:
1186 butlerQC.put(fgcmDatasetDict[
'fgcmZeropoints'],
1187 getattr(outputRefs, f
'fgcm_Cycle{cycle}_Zeropoints'))
1188 butlerQC.put(fgcmDatasetDict[
'fgcmAtmosphereParameters'],
1189 getattr(outputRefs, f
'fgcm_Cycle{cycle}_AtmosphereParameters'))
1190 if self.outputStandards:
1191 butlerQC.put(fgcmDatasetDict[
'fgcmStandardStars'],
1192 getattr(outputRefs, f
'fgcm_Cycle{cycle}_StandardStars'))
1199 for outputRefName
in outputRefs.keys():
1200 if outputRefName.endswith(
"Plot")
and f
"Cycle{self.config.cycleNumber}" in outputRefName:
1201 ref = getattr(outputRefs, outputRefName)
1202 plotHandleDict[outputRefName] = ref
1204 fgcmDatasetDict, _ = self._fgcmFitCycle(
1209 plotHandleDict=plotHandleDict,
1213 butlerQC.put(fgcmDatasetDict[
'fgcmFitParameters'], outputRefs.fgcmFitParameters)
1214 butlerQC.put(fgcmDatasetDict[
'fgcmFlaggedStars'], outputRefs.fgcmFlaggedStars)
1215 if self.outputZeropoints:
1216 butlerQC.put(fgcmDatasetDict[
'fgcmZeropoints'], outputRefs.fgcmZeropoints)
1217 butlerQC.put(fgcmDatasetDict[
'fgcmAtmosphereParameters'], outputRefs.fgcmAtmosphereParameters)
1218 if self.outputStandards:
1219 butlerQC.put(fgcmDatasetDict[
'fgcmStandardStars'], outputRefs.fgcmStandardStars)
1226 plotHandleDict=None,
1236 camera : `lsst.afw.cameraGeom.Camera`
1238 All handles are `lsst.daf.butler.DeferredDatasetHandle`
1239 handle dictionary with keys:
1241 ``"fgcmLookUpTable"``
1242 handle for the FGCM look-up table.
1243 ``"fgcmVisitCatalog"``
1244 handle for visit summary catalog.
1245 ``"fgcmStarObservations"``
1246 handle for star observation catalog.
1248 handle for star id catalog.
1249 ``"fgcmStarIndices"``
1250 handle for star index catalog.
1251 ``"fgcmReferenceStars"``
1252 handle for matched reference star catalog.
1253 ``"fgcmFlaggedStars"``
1254 handle for flagged star catalog.
1255 ``"fgcmFitParameters"``
1256 handle for fit parameter catalog.
1257 butlerQC : `lsst.pipe.base.QuantumContext`, optional
1258 Quantum context used for serializing plots.
1259 plotHandleDict : `dict` [`str`, `lsst.daf.butler.DatasetRef`], optional
1260 Dictionary of plot dataset refs, keyed by plot name.
1261 config : `lsst.pex.config.Config`, optional
1262 Configuration to use to override self.config.
1263 nCore : `int`, optional
1264 Number of cores to use during fitting.
1265 multiCycle : `bool`, optional
1266 Is this part of a multicycle run?
1270 fgcmDatasetDict : `dict`
1271 Dictionary of datasets to persist.
1273 if config
is not None:
1276 _config = self.config
1279 self.maxIter = _config.maxIterBeforeFinalCycle
1280 self.outputStandards = _config.outputStandardsBeforeFinalCycle
1281 self.outputZeropoints = _config.outputZeropointsBeforeFinalCycle
1282 self.resetFitParameters =
True
1284 if _config.isFinalCycle:
1289 self.outputStandards =
True
1290 self.outputZeropoints =
True
1291 self.resetFitParameters =
False
1293 lutCat = handleDict[
'fgcmLookUpTable'].get()
1294 fgcmLut, lutIndexVals, lutStd = translateFgcmLut(lutCat,
1295 dict(_config.physicalFilterMap))
1299 doPlots = _config.doPlots
1300 if doPlots
and multiCycle:
1301 if _config.cycleNumber < (_config.multipleCyclesFinalCycleNumber - 1) \
1302 and not _config.doPlotsBeforeFinalCycles:
1305 configDict = makeConfigDict(_config, self.log, camera,
1306 self.maxIter, self.resetFitParameters,
1307 self.outputZeropoints,
1308 lutIndexVals[0][
'FILTERNAMES'],
1313 visitCat = handleDict[
'fgcmVisitCatalog'].get()
1314 fgcmExpInfo = translateVisitCatalog(visitCat)
1318 self.config.defaultCameraOrientation)
1320 noFitsDict = {
'lutIndex': lutIndexVals,
1322 'expInfo': fgcmExpInfo,
1323 'focalPlaneProjector': focalPlaneProjector}
1326 fgcmFitCycle = fgcm.FgcmFitCycle(
1329 noFitsDict=noFitsDict,
1332 plotHandleDict=plotHandleDict,
1336 if (fgcmFitCycle.initialCycle):
1338 fgcmPars = fgcm.FgcmParameters.newParsWithArrays(fgcmFitCycle.fgcmConfig,
1342 plotHandleDict=plotHandleDict)
1345 parCat = handleDict[
'fgcmFitParameters']
1347 parCat = handleDict[
'fgcmFitParameters'].get()
1348 inParInfo, inParams, inSuperStar = self._loadParameters(parCat)
1350 fgcmPars = fgcm.FgcmParameters.loadParsWithArrays(fgcmFitCycle.fgcmConfig,
1356 plotHandleDict=plotHandleDict)
1359 fgcmStars = fgcm.FgcmStars(fgcmFitCycle.fgcmConfig, butlerQC=butlerQC, plotHandleDict=plotHandleDict)
1361 starObs = handleDict[
'fgcmStarObservations'].get()
1362 starIds = handleDict[
'fgcmStarIds'].get()
1363 if not self.config.useParquetCatalogFormat:
1364 starIndices = handleDict[
'fgcmStarIndices'].get()
1369 if 'fgcmFlaggedStars' in handleDict:
1371 flaggedStars = handleDict[
'fgcmFlaggedStars']
1373 flaggedStars = handleDict[
'fgcmFlaggedStars'].get()
1374 flagId = flaggedStars[
'objId'][:]
1375 flagFlag = flaggedStars[
'objFlag'][:]
1378 elif self.config.useParquetCatalogFormat:
1384 (flagged,) = (starIds[
'obj_flag'] > 0).nonzero()
1385 flagId = starIds[
'fgcm_id'][flagged]
1386 flagFlag = starIds[
'obj_flag'][flagged]
1391 if _config.doReferenceCalibration:
1392 refStars = handleDict[
'fgcmReferenceStars'].get()
1394 refMag, refMagErr = extractReferenceMags(refStars,
1396 _config.physicalFilterMap)
1398 refId = refStars[
'fgcm_id'][:]
1408 if self.config.useParquetCatalogFormat:
1409 visitIndex = np.searchsorted(fgcmExpInfo[
'VISIT'], starObs[
'visit'])
1411 visitIndex = np.searchsorted(fgcmExpInfo[
'VISIT'], starObs[
'visit'][starIndices[
'obsIndex']])
1420 if self.config.useParquetCatalogFormat:
1423 fgcmStars.loadStars(fgcmPars,
1424 starObs[
'visit'][:],
1425 starObs[
'detector'][:],
1428 starObs[
'inst_mag'][:],
1429 starObs[
'inst_mag_err'][:],
1430 fgcmExpInfo[
'FILTERNAME'][visitIndex],
1431 starIds[
'fgcm_id'][:],
1434 starIds[
'obs_arr_index'][:],
1435 starIds[
'n_obs'][:],
1436 obsX=starObs[
'x'][:],
1437 obsY=starObs[
'y'][:],
1438 obsDeltaMagBkg=starObs[
'delta_mag_bkg'][:],
1439 obsDeltaAper=starObs[
'delta_mag_aper'][:],
1442 refMagErr=refMagErr,
1450 conv = starObs[0][
'ra'].asDegrees() / float(starObs[0][
'ra'])
1452 fgcmStars.loadStars(fgcmPars,
1453 starObs[
'visit'][starIndices[
'obsIndex']],
1454 starObs[
'ccd'][starIndices[
'obsIndex']],
1455 starObs[
'ra'][starIndices[
'obsIndex']] * conv,
1456 starObs[
'dec'][starIndices[
'obsIndex']] * conv,
1457 starObs[
'instMag'][starIndices[
'obsIndex']],
1458 starObs[
'instMagErr'][starIndices[
'obsIndex']],
1459 fgcmExpInfo[
'FILTERNAME'][visitIndex],
1460 starIds[
'fgcm_id'][:],
1463 starIds[
'obsArrIndex'][:],
1465 obsX=starObs[
'x'][starIndices[
'obsIndex']],
1466 obsY=starObs[
'y'][starIndices[
'obsIndex']],
1467 obsDeltaMagBkg=starObs[
'deltaMagBkg'][starIndices[
'obsIndex']],
1468 obsDeltaAper=starObs[
'deltaMagAper'][starIndices[
'obsIndex']],
1469 psfCandidate=starObs[
'psf_candidate'][starIndices[
'obsIndex']],
1472 refMagErr=refMagErr,
1489 fgcmFitCycle.setLUT(fgcmLut)
1490 fgcmFitCycle.setStars(fgcmStars, fgcmPars)
1491 fgcmFitCycle.setPars(fgcmPars)
1494 fgcmFitCycle.finishSetup()
1503 fgcmDatasetDict = self._makeFgcmOutputDatasets(fgcmFitCycle)
1508 updatedPhotometricCutDict = {b: float(fgcmFitCycle.updatedPhotometricCut[i])
for
1509 i, b
in enumerate(_config.bands)}
1510 updatedHighCutDict = {band: float(fgcmFitCycle.updatedHighCut[i])
for
1511 i, band
in enumerate(_config.bands)}
1513 outConfig = copy.copy(_config)
1514 outConfig.update(cycleNumber=(_config.cycleNumber + 1),
1515 precomputeSuperStarInitialCycle=
False,
1516 freezeStdAtmosphere=
False,
1517 expGrayPhotometricCutDict=updatedPhotometricCutDict,
1518 expGrayHighCutDict=updatedHighCutDict)
1520 outConfig.connections.update(previousCycleNumber=str(_config.cycleNumber),
1521 cycleNumber=str(_config.cycleNumber + 1))
1524 configFileName =
'%s_cycle%02d_config.py' % (outConfig.outfileBase,
1525 outConfig.cycleNumber)
1526 outConfig.save(configFileName)
1528 if _config.isFinalCycle == 1:
1530 self.log.info(
"Everything is in place to run fgcmOutputProducts.py")
1532 self.log.info(
"Saved config for next cycle to %s" % (configFileName))
1533 self.log.info(
"Be sure to look at:")
1534 self.log.info(
" config.expGrayPhotometricCut")
1535 self.log.info(
" config.expGrayHighCut")
1536 self.log.info(
"If you are satisfied with the fit, please set:")
1537 self.log.info(
" config.isFinalCycle = True")
1539 fgcmFitCycle.freeSharedMemory()
1541 return fgcmDatasetDict, outConfig
1543 def _loadParameters(self, parCat):
1545 Load FGCM parameters from a previous fit cycle
1549 parCat : `lsst.afw.table.BaseCatalog`
1550 Parameter catalog in afw table form.
1554 inParInfo: `numpy.ndarray`
1555 Numpy array parameter information formatted for input to fgcm
1556 inParameters: `numpy.ndarray`
1557 Numpy array parameter values formatted for input to fgcm
1558 inSuperStar: `numpy.array`
1559 Superstar flat formatted for input to fgcm
1561 parLutFilterNames = np.array(parCat[0][
'lutFilterNames'].split(
','))
1562 parFitBands = np.array(parCat[0][
'fitBands'].split(
','))
1564 inParInfo = np.zeros(1, dtype=[(
'NCCD',
'i4'),
1565 (
'LUTFILTERNAMES', parLutFilterNames.dtype.str,
1566 (parLutFilterNames.size, )),
1567 (
'FITBANDS', parFitBands.dtype.str, (parFitBands.size, )),
1568 (
'LNTAUUNIT',
'f8'),
1569 (
'LNTAUSLOPEUNIT',
'f8'),
1570 (
'ALPHAUNIT',
'f8'),
1571 (
'LNPWVUNIT',
'f8'),
1572 (
'LNPWVSLOPEUNIT',
'f8'),
1573 (
'LNPWVQUADRATICUNIT',
'f8'),
1574 (
'LNPWVGLOBALUNIT',
'f8'),
1576 (
'QESYSUNIT',
'f8'),
1577 (
'FILTEROFFSETUNIT',
'f8'),
1578 (
'HASEXTERNALPWV',
'i2'),
1579 (
'HASEXTERNALTAU',
'i2')])
1580 inParInfo[
'NCCD'] = parCat[
'nCcd']
1581 inParInfo[
'LUTFILTERNAMES'][:] = parLutFilterNames
1582 inParInfo[
'FITBANDS'][:] = parFitBands
1583 inParInfo[
'HASEXTERNALPWV'] = parCat[
'hasExternalPwv']
1584 inParInfo[
'HASEXTERNALTAU'] = parCat[
'hasExternalTau']
1586 inParams = np.zeros(1, dtype=[(
'PARALPHA',
'f8', (parCat[
'parAlpha'].size, )),
1587 (
'PARO3',
'f8', (parCat[
'parO3'].size, )),
1588 (
'PARLNTAUINTERCEPT',
'f8',
1589 (parCat[
'parLnTauIntercept'].size, )),
1590 (
'PARLNTAUSLOPE',
'f8',
1591 (parCat[
'parLnTauSlope'].size, )),
1592 (
'PARLNPWVINTERCEPT',
'f8',
1593 (parCat[
'parLnPwvIntercept'].size, )),
1594 (
'PARLNPWVSLOPE',
'f8',
1595 (parCat[
'parLnPwvSlope'].size, )),
1596 (
'PARLNPWVQUADRATIC',
'f8',
1597 (parCat[
'parLnPwvQuadratic'].size, )),
1598 (
'PARQESYSINTERCEPT',
'f8',
1599 (parCat[
'parQeSysIntercept'].size, )),
1600 (
'COMPQESYSSLOPE',
'f8',
1601 (parCat[
'compQeSysSlope'].size, )),
1602 (
'PARFILTEROFFSET',
'f8',
1603 (parCat[
'parFilterOffset'].size, )),
1604 (
'PARFILTEROFFSETFITFLAG',
'i2',
1605 (parCat[
'parFilterOffsetFitFlag'].size, )),
1606 (
'PARRETRIEVEDLNPWVSCALE',
'f8'),
1607 (
'PARRETRIEVEDLNPWVOFFSET',
'f8'),
1608 (
'PARRETRIEVEDLNPWVNIGHTLYOFFSET',
'f8',
1609 (parCat[
'parRetrievedLnPwvNightlyOffset'].size, )),
1610 (
'COMPABSTHROUGHPUT',
'f8',
1611 (parCat[
'compAbsThroughput'].size, )),
1612 (
'COMPREFOFFSET',
'f8',
1613 (parCat[
'compRefOffset'].size, )),
1614 (
'COMPREFSIGMA',
'f8',
1615 (parCat[
'compRefSigma'].size, )),
1616 (
'COMPMIRRORCHROMATICITY',
'f8',
1617 (parCat[
'compMirrorChromaticity'].size, )),
1618 (
'MIRRORCHROMATICITYPIVOT',
'f8',
1619 (parCat[
'mirrorChromaticityPivot'].size, )),
1620 (
'COMPCCDCHROMATICITY',
'f8',
1621 (parCat[
'compCcdChromaticity'].size, )),
1622 (
'COMPMEDIANSEDSLOPE',
'f8',
1623 (parCat[
'compMedianSedSlope'].size, )),
1624 (
'COMPAPERCORRPIVOT',
'f8',
1625 (parCat[
'compAperCorrPivot'].size, )),
1626 (
'COMPAPERCORRSLOPE',
'f8',
1627 (parCat[
'compAperCorrSlope'].size, )),
1628 (
'COMPAPERCORRSLOPEERR',
'f8',
1629 (parCat[
'compAperCorrSlopeErr'].size, )),
1630 (
'COMPAPERCORRRANGE',
'f8',
1631 (parCat[
'compAperCorrRange'].size, )),
1632 (
'COMPMODELERREXPTIMEPIVOT',
'f8',
1633 (parCat[
'compModelErrExptimePivot'].size, )),
1634 (
'COMPMODELERRFWHMPIVOT',
'f8',
1635 (parCat[
'compModelErrFwhmPivot'].size, )),
1636 (
'COMPMODELERRSKYPIVOT',
'f8',
1637 (parCat[
'compModelErrSkyPivot'].size, )),
1638 (
'COMPMODELERRPARS',
'f8',
1639 (parCat[
'compModelErrPars'].size, )),
1640 (
'COMPEXPGRAY',
'f8',
1641 (parCat[
'compExpGray'].size, )),
1642 (
'COMPVARGRAY',
'f8',
1643 (parCat[
'compVarGray'].size, )),
1644 (
'COMPEXPDELTAMAGBKG',
'f8',
1645 (parCat[
'compExpDeltaMagBkg'].size, )),
1646 (
'COMPNGOODSTARPEREXP',
'i4',
1647 (parCat[
'compNGoodStarPerExp'].size, )),
1648 (
'COMPEXPREFOFFSET',
'f8',
1649 (parCat[
'compExpRefOffset'].size, )),
1650 (
'COMPSIGFGCM',
'f8',
1651 (parCat[
'compSigFgcm'].size, )),
1652 (
'COMPSIGMACAL',
'f8',
1653 (parCat[
'compSigmaCal'].size, )),
1654 (
'COMPRETRIEVEDLNPWV',
'f8',
1655 (parCat[
'compRetrievedLnPwv'].size, )),
1656 (
'COMPRETRIEVEDLNPWVRAW',
'f8',
1657 (parCat[
'compRetrievedLnPwvRaw'].size, )),
1658 (
'COMPRETRIEVEDLNPWVFLAG',
'i2',
1659 (parCat[
'compRetrievedLnPwvFlag'].size, )),
1660 (
'COMPRETRIEVEDTAUNIGHT',
'f8',
1661 (parCat[
'compRetrievedTauNight'].size, )),
1662 (
'COMPEPSILON',
'f8',
1663 (parCat[
'compEpsilon'].size, )),
1664 (
'COMPMEDDELTAAPER',
'f8',
1665 (parCat[
'compMedDeltaAper'].size, )),
1666 (
'COMPGLOBALEPSILON',
'f4',
1667 (parCat[
'compGlobalEpsilon'].size, )),
1668 (
'COMPEPSILONMAP',
'f4',
1669 (parCat[
'compEpsilonMap'].size, )),
1670 (
'COMPEPSILONNSTARMAP',
'i4',
1671 (parCat[
'compEpsilonNStarMap'].size, )),
1672 (
'COMPEPSILONCCDMAP',
'f4',
1673 (parCat[
'compEpsilonCcdMap'].size, )),
1674 (
'COMPEPSILONCCDNSTARMAP',
'i4',
1675 (parCat[
'compEpsilonCcdNStarMap'].size, ))])
1677 inParams[
'PARALPHA'][:] = parCat[
'parAlpha'][0, :]
1678 inParams[
'PARO3'][:] = parCat[
'parO3'][0, :]
1679 inParams[
'PARLNTAUINTERCEPT'][:] = parCat[
'parLnTauIntercept'][0, :]
1680 inParams[
'PARLNTAUSLOPE'][:] = parCat[
'parLnTauSlope'][0, :]
1681 inParams[
'PARLNPWVINTERCEPT'][:] = parCat[
'parLnPwvIntercept'][0, :]
1682 inParams[
'PARLNPWVSLOPE'][:] = parCat[
'parLnPwvSlope'][0, :]
1683 inParams[
'PARLNPWVQUADRATIC'][:] = parCat[
'parLnPwvQuadratic'][0, :]
1684 inParams[
'PARQESYSINTERCEPT'][:] = parCat[
'parQeSysIntercept'][0, :]
1685 inParams[
'COMPQESYSSLOPE'][:] = parCat[
'compQeSysSlope'][0, :]
1686 inParams[
'PARFILTEROFFSET'][:] = parCat[
'parFilterOffset'][0, :]
1687 inParams[
'PARFILTEROFFSETFITFLAG'][:] = parCat[
'parFilterOffsetFitFlag'][0, :]
1688 inParams[
'PARRETRIEVEDLNPWVSCALE'] = parCat[
'parRetrievedLnPwvScale']
1689 inParams[
'PARRETRIEVEDLNPWVOFFSET'] = parCat[
'parRetrievedLnPwvOffset']
1690 inParams[
'PARRETRIEVEDLNPWVNIGHTLYOFFSET'][:] = parCat[
'parRetrievedLnPwvNightlyOffset'][0, :]
1691 inParams[
'COMPABSTHROUGHPUT'][:] = parCat[
'compAbsThroughput'][0, :]
1692 inParams[
'COMPREFOFFSET'][:] = parCat[
'compRefOffset'][0, :]
1693 inParams[
'COMPREFSIGMA'][:] = parCat[
'compRefSigma'][0, :]
1694 inParams[
'COMPMIRRORCHROMATICITY'][:] = parCat[
'compMirrorChromaticity'][0, :]
1695 inParams[
'MIRRORCHROMATICITYPIVOT'][:] = parCat[
'mirrorChromaticityPivot'][0, :]
1696 inParams[
'COMPCCDCHROMATICITY'][:] = parCat[
'compCcdChromaticity'][0, :]
1697 inParams[
'COMPMEDIANSEDSLOPE'][:] = parCat[
'compMedianSedSlope'][0, :]
1698 inParams[
'COMPAPERCORRPIVOT'][:] = parCat[
'compAperCorrPivot'][0, :]
1699 inParams[
'COMPAPERCORRSLOPE'][:] = parCat[
'compAperCorrSlope'][0, :]
1700 inParams[
'COMPAPERCORRSLOPEERR'][:] = parCat[
'compAperCorrSlopeErr'][0, :]
1701 inParams[
'COMPAPERCORRRANGE'][:] = parCat[
'compAperCorrRange'][0, :]
1702 inParams[
'COMPMODELERREXPTIMEPIVOT'][:] = parCat[
'compModelErrExptimePivot'][0, :]
1703 inParams[
'COMPMODELERRFWHMPIVOT'][:] = parCat[
'compModelErrFwhmPivot'][0, :]
1704 inParams[
'COMPMODELERRSKYPIVOT'][:] = parCat[
'compModelErrSkyPivot'][0, :]
1705 inParams[
'COMPMODELERRPARS'][:] = parCat[
'compModelErrPars'][0, :]
1706 inParams[
'COMPEXPGRAY'][:] = parCat[
'compExpGray'][0, :]
1707 inParams[
'COMPVARGRAY'][:] = parCat[
'compVarGray'][0, :]
1708 inParams[
'COMPEXPDELTAMAGBKG'][:] = parCat[
'compExpDeltaMagBkg'][0, :]
1709 inParams[
'COMPNGOODSTARPEREXP'][:] = parCat[
'compNGoodStarPerExp'][0, :]
1710 inParams[
'COMPEXPREFOFFSET'][:] = parCat[
'compExpRefOffset'][0, :]
1711 inParams[
'COMPSIGFGCM'][:] = parCat[
'compSigFgcm'][0, :]
1712 inParams[
'COMPSIGMACAL'][:] = parCat[
'compSigmaCal'][0, :]
1713 inParams[
'COMPRETRIEVEDLNPWV'][:] = parCat[
'compRetrievedLnPwv'][0, :]
1714 inParams[
'COMPRETRIEVEDLNPWVRAW'][:] = parCat[
'compRetrievedLnPwvRaw'][0, :]
1715 inParams[
'COMPRETRIEVEDLNPWVFLAG'][:] = parCat[
'compRetrievedLnPwvFlag'][0, :]
1716 inParams[
'COMPRETRIEVEDTAUNIGHT'][:] = parCat[
'compRetrievedTauNight'][0, :]
1717 inParams[
'COMPEPSILON'][:] = parCat[
'compEpsilon'][0, :]
1718 inParams[
'COMPMEDDELTAAPER'][:] = parCat[
'compMedDeltaAper'][0, :]
1719 inParams[
'COMPGLOBALEPSILON'][:] = parCat[
'compGlobalEpsilon'][0, :]
1720 inParams[
'COMPEPSILONMAP'][:] = parCat[
'compEpsilonMap'][0, :]
1721 inParams[
'COMPEPSILONNSTARMAP'][:] = parCat[
'compEpsilonNStarMap'][0, :]
1722 inParams[
'COMPEPSILONCCDMAP'][:] = parCat[
'compEpsilonCcdMap'][0, :]
1723 inParams[
'COMPEPSILONCCDNSTARMAP'][:] = parCat[
'compEpsilonCcdNStarMap'][0, :]
1725 inSuperStar = np.zeros(parCat[
'superstarSize'][0, :], dtype=
'f8')
1726 inSuperStar[:, :, :, :] = parCat[
'superstar'][0, :].reshape(inSuperStar.shape)
1728 return (inParInfo, inParams, inSuperStar)
1730 def _makeFgcmOutputDatasets(self, fgcmFitCycle):
1732 Persist FGCM datasets through the butler.
1736 fgcmFitCycle: `lsst.fgcm.FgcmFitCycle`
1737 Fgcm Fit cycle object
1739 fgcmDatasetDict = {}
1742 parInfo, pars = fgcmFitCycle.fgcmPars.parsToArrays()
1747 lutFilterNameString = comma.join([n.decode(
'utf-8')
1748 for n
in parInfo[
'LUTFILTERNAMES'][0]])
1749 fitBandString = comma.join([n.decode(
'utf-8')
1750 for n
in parInfo[
'FITBANDS'][0]])
1752 parSchema = self._makeParSchema(parInfo, pars, fgcmFitCycle.fgcmPars.parSuperStarFlat,
1753 lutFilterNameString, fitBandString)
1754 parCat = self._makeParCatalog(parSchema, parInfo, pars,
1755 fgcmFitCycle.fgcmPars.parSuperStarFlat,
1756 lutFilterNameString, fitBandString)
1758 fgcmDatasetDict[
'fgcmFitParameters'] = parCat
1763 flagStarSchema = self._makeFlagStarSchema()
1764 flagStarStruct = fgcmFitCycle.fgcmStars.getFlagStarIndices()
1765 flagStarCat = self._makeFlagStarCat(flagStarSchema, flagStarStruct)
1767 fgcmDatasetDict[
'fgcmFlaggedStars'] = flagStarCat
1770 if self.outputZeropoints:
1771 superStarChebSize = fgcmFitCycle.fgcmZpts.zpStruct[
'FGCM_FZPT_SSTAR_CHEB'].shape[1]
1772 zptChebSize = fgcmFitCycle.fgcmZpts.zpStruct[
'FGCM_FZPT_CHEB'].shape[1]
1774 zptSchema = makeZptSchema(superStarChebSize, zptChebSize)
1775 zptCat = makeZptCat(zptSchema, fgcmFitCycle.fgcmZpts.zpStruct)
1777 fgcmDatasetDict[
'fgcmZeropoints'] = zptCat
1781 atmSchema = makeAtmSchema()
1782 atmCat = makeAtmCat(atmSchema, fgcmFitCycle.fgcmZpts.atmStruct)
1784 fgcmDatasetDict[
'fgcmAtmosphereParameters'] = atmCat
1787 if self.outputStandards:
1788 stdStruct, goodBands = fgcmFitCycle.fgcmStars.retrieveStdStarCatalog(fgcmFitCycle.fgcmPars)
1789 stdSchema = makeStdSchema(len(goodBands))
1790 stdCat = makeStdCat(stdSchema, stdStruct, goodBands)
1792 fgcmDatasetDict[
'fgcmStandardStars'] = stdCat
1794 return fgcmDatasetDict
1796 def _makeParSchema(self, parInfo, pars, parSuperStarFlat,
1797 lutFilterNameString, fitBandString):
1799 Make the parameter persistence schema
1803 parInfo: `numpy.ndarray`
1804 Parameter information returned by fgcm
1805 pars: `numpy.ndarray`
1806 Parameter values returned by fgcm
1807 parSuperStarFlat: `numpy.array`
1808 Superstar flat values returned by fgcm
1809 lutFilterNameString: `str`
1810 Combined string of all the lutFilterNames
1811 fitBandString: `str`
1812 Combined string of all the fitBands
1816 parSchema: `afwTable.schema`
1822 parSchema.addField(
'nCcd', type=np.int32, doc=
'Number of CCDs')
1823 parSchema.addField(
'lutFilterNames', type=str, doc=
'LUT Filter names in parameter file',
1824 size=len(lutFilterNameString))
1825 parSchema.addField(
'fitBands', type=str, doc=
'Bands that were fit',
1826 size=len(fitBandString))
1827 parSchema.addField(
'lnTauUnit', type=np.float64, doc=
'Step units for ln(AOD)')
1828 parSchema.addField(
'lnTauSlopeUnit', type=np.float64,
1829 doc=
'Step units for ln(AOD) slope')
1830 parSchema.addField(
'alphaUnit', type=np.float64, doc=
'Step units for alpha')
1831 parSchema.addField(
'lnPwvUnit', type=np.float64, doc=
'Step units for ln(pwv)')
1832 parSchema.addField(
'lnPwvSlopeUnit', type=np.float64,
1833 doc=
'Step units for ln(pwv) slope')
1834 parSchema.addField(
'lnPwvQuadraticUnit', type=np.float64,
1835 doc=
'Step units for ln(pwv) quadratic term')
1836 parSchema.addField(
'lnPwvGlobalUnit', type=np.float64,
1837 doc=
'Step units for global ln(pwv) parameters')
1838 parSchema.addField(
'o3Unit', type=np.float64, doc=
'Step units for O3')
1839 parSchema.addField(
'qeSysUnit', type=np.float64, doc=
'Step units for mirror gray')
1840 parSchema.addField(
'filterOffsetUnit', type=np.float64, doc=
'Step units for filter offset')
1841 parSchema.addField(
'hasExternalPwv', type=np.int32, doc=
'Parameters fit using external pwv')
1842 parSchema.addField(
'hasExternalTau', type=np.int32, doc=
'Parameters fit using external tau')
1845 parSchema.addField(
'parAlpha', type=
'ArrayD', doc=
'Alpha parameter vector',
1846 size=pars[
'PARALPHA'].size)
1847 parSchema.addField(
'parO3', type=
'ArrayD', doc=
'O3 parameter vector',
1848 size=pars[
'PARO3'].size)
1849 parSchema.addField(
'parLnTauIntercept', type=
'ArrayD',
1850 doc=
'ln(Tau) intercept parameter vector',
1851 size=pars[
'PARLNTAUINTERCEPT'].size)
1852 parSchema.addField(
'parLnTauSlope', type=
'ArrayD',
1853 doc=
'ln(Tau) slope parameter vector',
1854 size=pars[
'PARLNTAUSLOPE'].size)
1855 parSchema.addField(
'parLnPwvIntercept', type=
'ArrayD', doc=
'ln(pwv) intercept parameter vector',
1856 size=pars[
'PARLNPWVINTERCEPT'].size)
1857 parSchema.addField(
'parLnPwvSlope', type=
'ArrayD', doc=
'ln(pwv) slope parameter vector',
1858 size=pars[
'PARLNPWVSLOPE'].size)
1859 parSchema.addField(
'parLnPwvQuadratic', type=
'ArrayD', doc=
'ln(pwv) quadratic parameter vector',
1860 size=pars[
'PARLNPWVQUADRATIC'].size)
1861 parSchema.addField(
'parQeSysIntercept', type=
'ArrayD', doc=
'Mirror gray intercept parameter vector',
1862 size=pars[
'PARQESYSINTERCEPT'].size)
1863 parSchema.addField(
'compQeSysSlope', type=
'ArrayD', doc=
'Mirror gray slope parameter vector',
1864 size=pars[0][
'COMPQESYSSLOPE'].size)
1865 parSchema.addField(
'parFilterOffset', type=
'ArrayD', doc=
'Filter offset parameter vector',
1866 size=pars[
'PARFILTEROFFSET'].size)
1867 parSchema.addField(
'parFilterOffsetFitFlag', type=
'ArrayI', doc=
'Filter offset parameter fit flag',
1868 size=pars[
'PARFILTEROFFSETFITFLAG'].size)
1869 parSchema.addField(
'parRetrievedLnPwvScale', type=np.float64,
1870 doc=
'Global scale for retrieved ln(pwv)')
1871 parSchema.addField(
'parRetrievedLnPwvOffset', type=np.float64,
1872 doc=
'Global offset for retrieved ln(pwv)')
1873 parSchema.addField(
'parRetrievedLnPwvNightlyOffset', type=
'ArrayD',
1874 doc=
'Nightly offset for retrieved ln(pwv)',
1875 size=pars[
'PARRETRIEVEDLNPWVNIGHTLYOFFSET'].size)
1876 parSchema.addField(
'compAbsThroughput', type=
'ArrayD',
1877 doc=
'Absolute throughput (relative to transmission curves)',
1878 size=pars[
'COMPABSTHROUGHPUT'].size)
1879 parSchema.addField(
'compRefOffset', type=
'ArrayD',
1880 doc=
'Offset between reference stars and calibrated stars',
1881 size=pars[
'COMPREFOFFSET'].size)
1882 parSchema.addField(
'compRefSigma', type=
'ArrayD',
1883 doc=
'Width of reference star/calibrated star distribution',
1884 size=pars[
'COMPREFSIGMA'].size)
1885 parSchema.addField(
'compMirrorChromaticity', type=
'ArrayD',
1886 doc=
'Computed mirror chromaticity terms',
1887 size=pars[
'COMPMIRRORCHROMATICITY'].size)
1888 parSchema.addField(
'mirrorChromaticityPivot', type=
'ArrayD',
1889 doc=
'Mirror chromaticity pivot mjd',
1890 size=pars[
'MIRRORCHROMATICITYPIVOT'].size)
1891 parSchema.addField(
'compCcdChromaticity', type=
'ArrayD',
1892 doc=
'Computed CCD chromaticity terms',
1893 size=pars[
'COMPCCDCHROMATICITY'].size)
1894 parSchema.addField(
'compMedianSedSlope', type=
'ArrayD',
1895 doc=
'Computed median SED slope (per band)',
1896 size=pars[
'COMPMEDIANSEDSLOPE'].size)
1897 parSchema.addField(
'compAperCorrPivot', type=
'ArrayD', doc=
'Aperture correction pivot',
1898 size=pars[
'COMPAPERCORRPIVOT'].size)
1899 parSchema.addField(
'compAperCorrSlope', type=
'ArrayD', doc=
'Aperture correction slope',
1900 size=pars[
'COMPAPERCORRSLOPE'].size)
1901 parSchema.addField(
'compAperCorrSlopeErr', type=
'ArrayD', doc=
'Aperture correction slope error',
1902 size=pars[
'COMPAPERCORRSLOPEERR'].size)
1903 parSchema.addField(
'compAperCorrRange', type=
'ArrayD', doc=
'Aperture correction range',
1904 size=pars[
'COMPAPERCORRRANGE'].size)
1905 parSchema.addField(
'compModelErrExptimePivot', type=
'ArrayD', doc=
'Model error exptime pivot',
1906 size=pars[
'COMPMODELERREXPTIMEPIVOT'].size)
1907 parSchema.addField(
'compModelErrFwhmPivot', type=
'ArrayD', doc=
'Model error fwhm pivot',
1908 size=pars[
'COMPMODELERRFWHMPIVOT'].size)
1909 parSchema.addField(
'compModelErrSkyPivot', type=
'ArrayD', doc=
'Model error sky pivot',
1910 size=pars[
'COMPMODELERRSKYPIVOT'].size)
1911 parSchema.addField(
'compModelErrPars', type=
'ArrayD', doc=
'Model error parameters',
1912 size=pars[
'COMPMODELERRPARS'].size)
1913 parSchema.addField(
'compExpGray', type=
'ArrayD', doc=
'Computed exposure gray',
1914 size=pars[
'COMPEXPGRAY'].size)
1915 parSchema.addField(
'compVarGray', type=
'ArrayD', doc=
'Computed exposure variance',
1916 size=pars[
'COMPVARGRAY'].size)
1917 parSchema.addField(
'compExpDeltaMagBkg', type=
'ArrayD',
1918 doc=
'Computed exposure offset due to background',
1919 size=pars[
'COMPEXPDELTAMAGBKG'].size)
1920 parSchema.addField(
'compNGoodStarPerExp', type=
'ArrayI',
1921 doc=
'Computed number of good stars per exposure',
1922 size=pars[
'COMPNGOODSTARPEREXP'].size)
1923 parSchema.addField(
'compExpRefOffset', type=
'ArrayD',
1924 doc=
'Computed per-visit median offset between standard stars and ref stars.',
1925 size=pars[
'COMPEXPREFOFFSET'].size)
1926 parSchema.addField(
'compSigFgcm', type=
'ArrayD', doc=
'Computed sigma_fgcm (intrinsic repeatability)',
1927 size=pars[
'COMPSIGFGCM'].size)
1928 parSchema.addField(
'compSigmaCal', type=
'ArrayD', doc=
'Computed sigma_cal (systematic error floor)',
1929 size=pars[
'COMPSIGMACAL'].size)
1930 parSchema.addField(
'compRetrievedLnPwv', type=
'ArrayD', doc=
'Retrieved ln(pwv) (smoothed)',
1931 size=pars[
'COMPRETRIEVEDLNPWV'].size)
1932 parSchema.addField(
'compRetrievedLnPwvRaw', type=
'ArrayD', doc=
'Retrieved ln(pwv) (raw)',
1933 size=pars[
'COMPRETRIEVEDLNPWVRAW'].size)
1934 parSchema.addField(
'compRetrievedLnPwvFlag', type=
'ArrayI', doc=
'Retrieved ln(pwv) Flag',
1935 size=pars[
'COMPRETRIEVEDLNPWVFLAG'].size)
1936 parSchema.addField(
'compRetrievedTauNight', type=
'ArrayD', doc=
'Retrieved tau (per night)',
1937 size=pars[
'COMPRETRIEVEDTAUNIGHT'].size)
1938 parSchema.addField(
'compEpsilon', type=
'ArrayD',
1939 doc=
'Computed epsilon background offset per visit (nJy/arcsec2)',
1940 size=pars[
'COMPEPSILON'].size)
1941 parSchema.addField(
'compMedDeltaAper', type=
'ArrayD',
1942 doc=
'Median delta mag aper per visit',
1943 size=pars[
'COMPMEDDELTAAPER'].size)
1944 parSchema.addField(
'compGlobalEpsilon', type=
'ArrayD',
1945 doc=
'Computed epsilon bkg offset (global) (nJy/arcsec2)',
1946 size=pars[
'COMPGLOBALEPSILON'].size)
1947 parSchema.addField(
'compEpsilonMap', type=
'ArrayD',
1948 doc=
'Computed epsilon maps (nJy/arcsec2)',
1949 size=pars[
'COMPEPSILONMAP'].size)
1950 parSchema.addField(
'compEpsilonNStarMap', type=
'ArrayI',
1951 doc=
'Number of stars per pixel in computed epsilon maps',
1952 size=pars[
'COMPEPSILONNSTARMAP'].size)
1953 parSchema.addField(
'compEpsilonCcdMap', type=
'ArrayD',
1954 doc=
'Computed epsilon ccd maps (nJy/arcsec2)',
1955 size=pars[
'COMPEPSILONCCDMAP'].size)
1956 parSchema.addField(
'compEpsilonCcdNStarMap', type=
'ArrayI',
1957 doc=
'Number of stars per ccd bin in epsilon ccd maps',
1958 size=pars[
'COMPEPSILONCCDNSTARMAP'].size)
1960 parSchema.addField(
'superstarSize', type=
'ArrayI', doc=
'Superstar matrix size',
1962 parSchema.addField(
'superstar', type=
'ArrayD', doc=
'Superstar matrix (flattened)',
1963 size=parSuperStarFlat.size)
1967 def _makeParCatalog(self, parSchema, parInfo, pars, parSuperStarFlat,
1968 lutFilterNameString, fitBandString):
1970 Make the FGCM parameter catalog for persistence
1974 parSchema: `lsst.afw.table.Schema`
1975 Parameter catalog schema
1976 pars: `numpy.ndarray`
1977 FGCM parameters to put into parCat
1978 parSuperStarFlat: `numpy.array`
1979 FGCM superstar flat array to put into parCat
1980 lutFilterNameString: `str`
1981 Combined string of all the lutFilterNames
1982 fitBandString: `str`
1983 Combined string of all the fitBands
1987 parCat: `afwTable.BasicCatalog`
1988 Atmosphere and instrumental model parameter catalog for persistence
1996 rec = parCat.addNew()
1999 rec[
'nCcd'] = parInfo[
'NCCD'][0]
2000 rec[
'lutFilterNames'] = lutFilterNameString
2001 rec[
'fitBands'] = fitBandString
2003 rec[
'hasExternalPwv'] = 0
2004 rec[
'hasExternalTau'] = 0
2008 scalarNames = [
'parRetrievedLnPwvScale',
'parRetrievedLnPwvOffset']
2010 arrNames = [
'parAlpha',
'parO3',
'parLnTauIntercept',
'parLnTauSlope',
2011 'parLnPwvIntercept',
'parLnPwvSlope',
'parLnPwvQuadratic',
2012 'parQeSysIntercept',
'compQeSysSlope',
2013 'parRetrievedLnPwvNightlyOffset',
'compAperCorrPivot',
2014 'parFilterOffset',
'parFilterOffsetFitFlag',
2015 'compAbsThroughput',
'compRefOffset',
'compRefSigma',
2016 'compMirrorChromaticity',
'mirrorChromaticityPivot',
'compCcdChromaticity',
2017 'compAperCorrSlope',
'compAperCorrSlopeErr',
'compAperCorrRange',
2018 'compModelErrExptimePivot',
'compModelErrFwhmPivot',
2019 'compModelErrSkyPivot',
'compModelErrPars',
2020 'compExpGray',
'compVarGray',
'compNGoodStarPerExp',
'compSigFgcm',
2021 'compSigmaCal',
'compExpDeltaMagBkg',
'compMedianSedSlope',
2022 'compRetrievedLnPwv',
'compRetrievedLnPwvRaw',
'compRetrievedLnPwvFlag',
2023 'compRetrievedTauNight',
'compEpsilon',
'compMedDeltaAper',
2024 'compGlobalEpsilon',
'compEpsilonMap',
'compEpsilonNStarMap',
2025 'compEpsilonCcdMap',
'compEpsilonCcdNStarMap',
'compExpRefOffset']
2027 for scalarName
in scalarNames:
2028 rec[scalarName] = pars[scalarName.upper()][0]
2030 for arrName
in arrNames:
2031 rec[arrName][:] = np.atleast_1d(pars[0][arrName.upper()])[:]
2034 rec[
'superstarSize'][:] = parSuperStarFlat.shape
2035 rec[
'superstar'][:] = parSuperStarFlat.ravel()
2039 def _makeFlagStarSchema(self):
2041 Make the flagged-stars schema
2045 flagStarSchema: `lsst.afw.table.Schema`
2050 flagStarSchema.addField(
'objId', type=np.int32, doc=
'FGCM object id')
2051 flagStarSchema.addField(
'objFlag', type=np.int32, doc=
'FGCM object flag')
2053 return flagStarSchema
2055 def _makeFlagStarCat(self, flagStarSchema, flagStarStruct):
2057 Make the flagged star catalog for persistence
2061 flagStarSchema: `lsst.afw.table.Schema`
2063 flagStarStruct: `numpy.ndarray`
2064 Flagged star structure from fgcm
2068 flagStarCat: `lsst.afw.table.BaseCatalog`
2069 Flagged star catalog for persistence
2073 flagStarCat.resize(flagStarStruct.size)
2075 flagStarCat[
'objId'][:] = flagStarStruct[
'OBJID']
2076 flagStarCat[
'objFlag'][:] = flagStarStruct[
'OBJFLAG']
Defines the fields and offsets for a table.