29import astropy.units
as u
34import lsst.pipe.base
as pipeBase
38from lsst.pipe.base
import Instrument
39from lsst.pipe.tasks.colorterms
import ColortermLibrary
40from lsst.verify
import Job, Measurement
43 LoadReferenceObjectsConfig)
49__all__ = [
"JointcalConfig",
"JointcalTask"]
51Photometry = collections.namedtuple(
'Photometry', (
'fit',
'model'))
52Astrometry = collections.namedtuple(
'Astrometry', (
'fit',
'model',
'sky_to_tan_projection'))
57 meas = Measurement(job.metrics[name], value)
58 job.measurements.insert(meas)
62 """Lookup function that asserts/hopes that a static calibration dataset
63 exists in a particular collection, since this task can
't provide a single
64 date/time to use to search for one properly.
66 This
is mostly useful
for the ``camera`` dataset,
in cases where the task
's
67 quantum dimensions do *not* include something temporal, like ``exposure``
72 datasetType : `lsst.daf.butler.DatasetType`
73 Type of dataset being searched
for.
74 registry : `lsst.daf.butler.Registry`
75 Data repository registry to search.
76 quantumDataId : `lsst.daf.butler.DataCoordinate`
77 Data ID of the quantum this camera should match.
78 collections : `Iterable` [ `str` ]
79 Collections that should be searched - but this lookup function works
80 by ignoring this
in favor of a more-
or-less hard-coded value.
84 refs : `Iterator` [ `lsst.daf.butler.DatasetRef` ]
85 Iterator over dataset references; should have only one element.
89 This implementation duplicates one
in fgcmcal,
and is at least quite
90 similar to another
in cp_pipe. This duplicate has the most documentation.
91 Fixing this
is DM-29661.
93 instrument = Instrument.fromName(quantumDataId["instrument"], registry)
94 unboundedCollection = instrument.makeUnboundedCalibrationRunName()
95 return registry.queryDatasets(datasetType,
97 collections=[unboundedCollection],
102 """Lookup function that finds all refcats for all visits that overlap a
103 tract, rather than just the refcats that directly overlap the tract.
107 datasetType : `lsst.daf.butler.DatasetType`
108 Type of dataset being searched for.
109 registry : `lsst.daf.butler.Registry`
110 Data repository registry to search.
111 quantumDataId : `lsst.daf.butler.DataCoordinate`
112 Data ID of the quantum; expected to be something we can use
as a
113 constraint to query
for overlapping visits.
114 collections : `Iterable` [ `str` ]
115 Collections to search.
119 refs : `Iterator` [ `lsst.daf.butler.DatasetRef` ]
120 Iterator over refcat references.
129 for visit_data_id
in set(registry.queryDataIds(
"visit", dataId=quantumDataId).expanded()):
131 registry.queryDatasets(
133 collections=collections,
134 dataId=visit_data_id,
142 dimensions=(
"skymap",
"tract",
"instrument",
"physical_filter")):
143 """Middleware input/output connections for jointcal data."""
144 inputCamera = pipeBase.connectionTypes.PrerequisiteInput(
145 doc=
"The camera instrument that took these observations.",
147 storageClass=
"Camera",
148 dimensions=(
"instrument",),
150 lookupFunction=lookupStaticCalibrations,
152 inputSourceTableVisit = pipeBase.connectionTypes.Input(
153 doc=
"Source table in parquet format, per visit",
154 name=
"sourceTable_visit",
155 storageClass=
"DataFrame",
156 dimensions=(
"instrument",
"visit"),
160 inputVisitSummary = pipeBase.connectionTypes.Input(
161 doc=(
"Per-visit consolidated exposure metadata built from calexps. "
162 "These catalogs use detector id for the id and must be sorted for "
163 "fast lookups of a detector."),
165 storageClass=
"ExposureCatalog",
166 dimensions=(
"instrument",
"visit"),
170 astrometryRefCat = pipeBase.connectionTypes.PrerequisiteInput(
171 doc=
"The astrometry reference catalog to match to loaded input catalog sources.",
172 name=
"gaia_dr2_20200414",
173 storageClass=
"SimpleCatalog",
174 dimensions=(
"skypix",),
177 lookupFunction=lookupVisitRefCats,
179 photometryRefCat = pipeBase.connectionTypes.PrerequisiteInput(
180 doc=
"The photometry reference catalog to match to loaded input catalog sources.",
181 name=
"ps1_pv3_3pi_20170110",
182 storageClass=
"SimpleCatalog",
183 dimensions=(
"skypix",),
186 lookupFunction=lookupVisitRefCats,
189 outputWcs = pipeBase.connectionTypes.Output(
190 doc=(
"Per-tract, per-visit world coordinate systems derived from the fitted model."
191 " These catalogs only contain entries for detectors with an output, and use"
192 " the detector id for the catalog id, sorted on id for fast lookups of a detector."),
193 name=
"jointcalSkyWcsCatalog",
194 storageClass=
"ExposureCatalog",
195 dimensions=(
"instrument",
"visit",
"skymap",
"tract"),
198 outputPhotoCalib = pipeBase.connectionTypes.Output(
199 doc=(
"Per-tract, per-visit photometric calibrations derived from the fitted model."
200 " These catalogs only contain entries for detectors with an output, and use"
201 " the detector id for the catalog id, sorted on id for fast lookups of a detector."),
202 name=
"jointcalPhotoCalibCatalog",
203 storageClass=
"ExposureCatalog",
204 dimensions=(
"instrument",
"visit",
"skymap",
"tract"),
212 for name
in (
"astrometry",
"photometry"):
213 vars()[f
"{name}_matched_fittedStars"] = pipeBase.connectionTypes.Output(
214 doc=f
"The number of cross-matched fittedStars for {name}",
215 name=f
"metricvalue_jointcal_{name}_matched_fittedStars",
216 storageClass=
"MetricValue",
217 dimensions=(
"skymap",
"tract",
"instrument",
"physical_filter"),
219 vars()[f
"{name}_collected_refStars"] = pipeBase.connectionTypes.Output(
220 doc=f
"The number of {name} reference stars drawn from the reference catalog, before matching.",
221 name=f
"metricvalue_jointcal_{name}_collected_refStars",
222 storageClass=
"MetricValue",
223 dimensions=(
"skymap",
"tract",
"instrument",
"physical_filter"),
225 vars()[f
"{name}_prepared_refStars"] = pipeBase.connectionTypes.Output(
226 doc=f
"The number of {name} reference stars matched to fittedStars.",
227 name=f
"metricvalue_jointcal_{name}_prepared_refStars",
228 storageClass=
"MetricValue",
229 dimensions=(
"skymap",
"tract",
"instrument",
"physical_filter"),
231 vars()[f
"{name}_prepared_fittedStars"] = pipeBase.connectionTypes.Output(
232 doc=f
"The number of cross-matched fittedStars after cleanup, for {name}.",
233 name=f
"metricvalue_jointcal_{name}_prepared_fittedStars",
234 storageClass=
"MetricValue",
235 dimensions=(
"skymap",
"tract",
"instrument",
"physical_filter"),
237 vars()[f
"{name}_prepared_ccdImages"] = pipeBase.connectionTypes.Output(
238 doc=f
"The number of ccdImages that will be fit for {name}, after cleanup.",
239 name=f
"metricvalue_jointcal_{name}_prepared_ccdImages",
240 storageClass=
"MetricValue",
241 dimensions=(
"skymap",
"tract",
"instrument",
"physical_filter"),
243 vars()[f
"{name}_final_chi2"] = pipeBase.connectionTypes.Output(
244 doc=f
"The final chi2 of the {name} fit.",
245 name=f
"metricvalue_jointcal_{name}_final_chi2",
246 storageClass=
"MetricValue",
247 dimensions=(
"skymap",
"tract",
"instrument",
"physical_filter"),
249 vars()[f
"{name}_final_ndof"] = pipeBase.connectionTypes.Output(
250 doc=f
"The number of degrees of freedom of the fitted {name}.",
251 name=f
"metricvalue_jointcal_{name}_final_ndof",
252 storageClass=
"MetricValue",
253 dimensions=(
"skymap",
"tract",
"instrument",
"physical_filter"),
264 if not config.doAstrometry:
265 self.prerequisiteInputs.
remove(
"astrometryRefCat")
266 self.outputs.
remove(
"outputWcs")
267 for key
in list(self.outputs):
268 if "metricvalue_jointcal_astrometry" in key:
270 if not config.doPhotometry:
271 self.prerequisiteInputs.
remove(
"photometryRefCat")
272 self.outputs.
remove(
"outputPhotoCalib")
273 for key
in list(self.outputs):
274 if "metricvalue_jointcal_photometry" in key:
279 pipelineConnections=JointcalTaskConnections):
280 """Configuration for JointcalTask"""
282 doAstrometry = pexConfig.Field(
283 doc=
"Fit astrometry and write the fitted result.",
287 doPhotometry = pexConfig.Field(
288 doc=
"Fit photometry and write the fitted result.",
292 sourceFluxType = pexConfig.Field(
294 doc=
"Source flux field to use in source selection and to get fluxes from the catalog.",
295 default=
'apFlux_12_0'
297 positionErrorPedestal = pexConfig.Field(
298 doc=
"Systematic term to apply to the measured position error (pixels)",
302 photometryErrorPedestal = pexConfig.Field(
303 doc=
"Systematic term to apply to the measured error on flux or magnitude as a "
304 "fraction of source flux or magnitude delta (e.g. 0.05 is 5% of flux or +50 millimag).",
309 matchCut = pexConfig.Field(
310 doc=
"Matching radius between fitted and reference stars (arcseconds)",
314 minMeasurements = pexConfig.Field(
315 doc=
"Minimum number of associated measured stars for a fitted star to be included in the fit",
319 minMeasuredStarsPerCcd = pexConfig.Field(
320 doc=
"Minimum number of measuredStars per ccdImage before printing warnings",
324 minRefStarsPerCcd = pexConfig.Field(
325 doc=
"Minimum number of measuredStars per ccdImage before printing warnings",
329 allowLineSearch = pexConfig.Field(
330 doc=
"Allow a line search during minimization, if it is reasonable for the model"
331 " (models with a significant non-linear component, e.g. constrainedPhotometry).",
335 astrometrySimpleOrder = pexConfig.Field(
336 doc=
"Polynomial order for fitting the simple astrometry model.",
340 astrometryChipOrder = pexConfig.Field(
341 doc=
"Order of the per-chip transform for the constrained astrometry model.",
345 astrometryVisitOrder = pexConfig.Field(
346 doc=
"Order of the per-visit transform for the constrained astrometry model.",
350 useInputWcs = pexConfig.Field(
351 doc=
"Use the input calexp WCSs to initialize a SimpleAstrometryModel.",
355 astrometryModel = pexConfig.ChoiceField(
356 doc=
"Type of model to fit to astrometry",
358 default=
"constrained",
359 allowed={
"simple":
"One polynomial per ccd",
360 "constrained":
"One polynomial per ccd, and one polynomial per visit"}
362 photometryModel = pexConfig.ChoiceField(
363 doc=
"Type of model to fit to photometry",
365 default=
"constrainedMagnitude",
366 allowed={
"simpleFlux":
"One constant zeropoint per ccd and visit, fitting in flux space.",
367 "constrainedFlux":
"Constrained zeropoint per ccd, and one polynomial per visit,"
368 " fitting in flux space.",
369 "simpleMagnitude":
"One constant zeropoint per ccd and visit,"
370 " fitting in magnitude space.",
371 "constrainedMagnitude":
"Constrained zeropoint per ccd, and one polynomial per visit,"
372 " fitting in magnitude space.",
375 applyColorTerms = pexConfig.Field(
376 doc=
"Apply photometric color terms to reference stars?"
377 "Requires that colorterms be set to a ColortermLibrary",
381 colorterms = pexConfig.ConfigField(
382 doc=
"Library of photometric reference catalog name to color term dict.",
383 dtype=ColortermLibrary,
385 photometryVisitOrder = pexConfig.Field(
386 doc=
"Order of the per-visit polynomial transform for the constrained photometry model.",
390 photometryDoRankUpdate = pexConfig.Field(
391 doc=(
"Do the rank update step during minimization. "
392 "Skipping this can help deal with models that are too non-linear."),
396 astrometryDoRankUpdate = pexConfig.Field(
397 doc=(
"Do the rank update step during minimization (should not change the astrometry fit). "
398 "Skipping this can help deal with models that are too non-linear."),
402 outlierRejectSigma = pexConfig.Field(
403 doc=
"How many sigma to reject outliers at during minimization.",
407 astrometryOutlierRelativeTolerance = pexConfig.Field(
408 doc=(
"Convergence tolerance for outlier rejection threshold when fitting astrometry. Iterations will "
409 "stop when the fractional change in the chi2 cut level is below this value. If tolerance is set "
410 "to zero, iterations will continue until there are no more outliers. We suggest a value of 0.002"
411 "as a balance between a shorter minimization runtime and achieving a final fitted model that is"
412 "close to the solution found when removing all outliers."),
416 maxPhotometrySteps = pexConfig.Field(
417 doc=
"Maximum number of minimize iterations to take when fitting photometry.",
421 maxAstrometrySteps = pexConfig.Field(
422 doc=
"Maximum number of minimize iterations to take when fitting astrometry.",
426 astrometryRefObjLoader = pexConfig.ConfigField(
427 dtype=LoadReferenceObjectsConfig,
428 doc=
"Reference object loader for astrometric fit",
430 photometryRefObjLoader = pexConfig.ConfigField(
431 dtype=LoadReferenceObjectsConfig,
432 doc=
"Reference object loader for photometric fit",
434 sourceSelector = sourceSelectorRegistry.makeField(
435 doc=
"How to select sources for cross-matching",
438 astrometryReferenceSelector = pexConfig.ConfigurableField(
439 target=ReferenceSourceSelectorTask,
440 doc=
"How to down-select the loaded astrometry reference catalog.",
442 photometryReferenceSelector = pexConfig.ConfigurableField(
443 target=ReferenceSourceSelectorTask,
444 doc=
"How to down-select the loaded photometry reference catalog.",
446 astrometryReferenceErr = pexConfig.Field(
447 doc=(
"Uncertainty on reference catalog coordinates [mas] to use in place of the `coord_*Err` fields. "
448 "If None, then raise an exception if the reference catalog is missing coordinate errors. "
449 "If specified, overrides any existing `coord_*Err` values."),
456 writeInitMatrix = pexConfig.Field(
458 doc=(
"Write the pre/post-initialization Hessian and gradient to text files, for debugging. "
459 "Output files will be written to `config.debugOutputPath` and will "
460 "be of the form 'astrometry_[pre|post]init-TRACT-FILTER-mat.txt'. "
461 "Note that these files are the dense versions of the matrix, and so may be very large."),
464 writeChi2FilesInitialFinal = pexConfig.Field(
466 doc=(
"Write .csv files containing the contributions to chi2 for the initialization and final fit. "
467 "Output files will be written to `config.debugOutputPath` and will "
468 "be of the form `astrometry_[initial|final]_chi2-TRACT-FILTER."),
471 writeChi2FilesOuterLoop = pexConfig.Field(
473 doc=(
"Write .csv files containing the contributions to chi2 for the outer fit loop. "
474 "Output files will be written to `config.debugOutputPath` and will "
475 "be of the form `astrometry_init-NN_chi2-TRACT-FILTER`."),
478 writeInitialModel = pexConfig.Field(
480 doc=(
"Write the pre-initialization model to text files, for debugging. "
481 "Output files will be written to `config.debugOutputPath` and will be "
482 "of the form `initial_astrometry_model-TRACT_FILTER.txt`."),
485 debugOutputPath = pexConfig.Field(
488 doc=(
"Path to write debug output files to. Used by "
489 "`writeInitialModel`, `writeChi2Files*`, `writeInitMatrix`.")
491 detailedProfile = pexConfig.Field(
494 doc=
"Output separate profiling information for different parts of jointcal, e.g. data read, fitting"
500 msg =
"applyColorTerms=True requires the `colorterms` field be set to a ColortermLibrary."
501 raise pexConfig.FieldValidationError(JointcalConfig.colorterms, self, msg)
503 msg = (
"Only doing astrometry, but Colorterms are not applied for astrometry;"
504 "applyColorTerms=True will be ignored.")
505 logging.getLogger(
"lsst.jointcal").warning(msg)
519 self.
sourceSelector[
"science"].signalToNoise.fluxField = f
"{self.sourceFluxType}_instFlux"
520 self.
sourceSelector[
"science"].signalToNoise.errField = f
"{self.sourceFluxType}_instFluxErr"
523 self.
sourceSelector[
"science"].isolated.parentName =
"parentSourceId"
524 self.
sourceSelector[
"science"].isolated.nChildName =
"deblend_nChild"
528 badFlags = [
"pixelFlags_edge",
529 "pixelFlags_saturated",
530 "pixelFlags_interpolatedCenter",
531 "pixelFlags_interpolated",
532 "pixelFlags_crCenter",
534 "hsmPsfMoments_flag",
535 f
"{self.sourceFluxType}_flag",
539 self.
sourceSelector[
"science"].requireFiniteRaDec.raColName =
"ra"
540 self.
sourceSelector[
"science"].requireFiniteRaDec.decColName =
"decl"
549 """Write model to outfile."""
550 with open(filename,
"w")
as file:
551 file.write(repr(model))
552 log.info(
"Wrote %s to file: %s", model, filename)
555@dataclasses.dataclass
557 """The input data jointcal needs for each detector/visit."""
559 """The visit identifier of this exposure."""
561 """The catalog derived from this exposure."""
563 """The VisitInfo of this exposure."""
565 """The detector of this exposure."""
567 """The photometric calibration of this exposure."""
569 """The WCS of this exposure."""
571 """The bounding box of this exposure."""
573 """The filter of this exposure."""
577 """Astrometricly and photometricly calibrate across multiple visits of the
581 ConfigClass = JointcalConfig
582 _DefaultName = "jointcal"
586 self.makeSubtask(
"sourceSelector")
587 if self.config.doAstrometry:
588 self.makeSubtask(
"astrometryReferenceSelector")
591 if self.config.doPhotometry:
592 self.makeSubtask(
"photometryReferenceSelector")
597 self.
job = Job.load_metrics_package(subset=
'jointcal')
602 inputs = butlerQC.get(inputRefs)
604 tract = butlerQC.quantum.dataId[
'tract']
605 if self.config.doAstrometry:
607 dataIds=[ref.datasetRef.dataId
608 for ref
in inputRefs.astrometryRefCat],
609 refCats=inputs.pop(
'astrometryRefCat'),
610 config=self.config.astrometryRefObjLoader,
611 name=self.config.connections.astrometryRefCat,
613 if self.config.doPhotometry:
615 dataIds=[ref.datasetRef.dataId
616 for ref
in inputRefs.photometryRefCat],
617 refCats=inputs.pop(
'photometryRefCat'),
618 config=self.config.photometryRefObjLoader,
619 name=self.config.connections.photometryRefCat,
621 outputs = self.
run(**inputs, tract=tract)
623 if self.config.doAstrometry:
624 self.
_put_output(butlerQC, outputs.outputWcs, outputRefs.outputWcs,
625 inputs[
'inputCamera'],
"setWcs")
626 if self.config.doPhotometry:
627 self.
_put_output(butlerQC, outputs.outputPhotoCalib, outputRefs.outputPhotoCalib,
628 inputs[
'inputCamera'],
"setPhotoCalib")
630 def _put_metrics(self, butlerQC, job, outputRefs):
631 """Persist all measured metrics stored in a job.
635 butlerQC : `lsst.pipe.base.ButlerQuantumContext`
636 A butler which is specialized to operate
in the context of a
637 `lsst.daf.butler.Quantum`; This
is the input to `runQuantum`.
638 job : `lsst.verify.job.Job`
639 Measurements of metrics to persist.
640 outputRefs : `list` [`lsst.pipe.base.connectionTypes.OutputQuantizedConnection`]
641 The DatasetRefs to persist the data to.
643 for key
in job.measurements.keys():
644 butlerQC.put(job.measurements[key], getattr(outputRefs, key.fqn.replace(
'jointcal.',
'')))
646 def _put_output(self, butlerQC, outputs, outputRefs, camera, setter):
647 """Persist the output datasets to their appropriate datarefs.
651 butlerQC : `lsst.pipe.base.ButlerQuantumContext`
652 A butler which is specialized to operate
in the context of a
653 `lsst.daf.butler.Quantum`; This
is the input to `runQuantum`.
656 The fitted objects to persist.
657 outputRefs : `list` [`lsst.pipe.base.connectionTypes.OutputQuantizedConnection`]
658 The DatasetRefs to persist the data to.
660 The camera
for this instrument, to get detector ids
from.
662 The method to call on the ExposureCatalog to set each output.
665 schema.addField('visit', type=
'L', doc=
'Visit number')
667 def new_catalog(visit, size):
668 """Return an catalog ready to be filled with appropriate output."""
671 catalog[
'visit'] = visit
673 metadata.add(
"COMMENT",
"Catalog id is detector id, sorted.")
674 metadata.add(
"COMMENT",
"Only detectors with data have entries.")
678 detectors_per_visit = collections.defaultdict(int)
681 detectors_per_visit[key[0]] += 1
683 for ref
in outputRefs:
684 visit = ref.dataId[
'visit']
685 catalog = new_catalog(visit, detectors_per_visit[visit])
688 for detector
in camera:
689 detectorId = detector.getId()
690 key = (ref.dataId[
'visit'], detectorId)
691 if key
not in outputs:
693 self.log.debug(
"No %s output for detector %s in visit %s",
694 setter[3:], detectorId, visit)
697 catalog[i].setId(detectorId)
698 getattr(catalog[i], setter)(outputs[key])
702 butlerQC.put(catalog, ref)
703 self.log.info(
"Wrote %s detectors to %s", i, ref)
705 def run(self, inputSourceTableVisit, inputVisitSummary, inputCamera, tract=None):
711 sourceFluxField =
"flux"
715 oldWcsList, bands = self.
_load_data(inputSourceTableVisit,
721 boundingCircle, center, radius, defaultFilter, epoch = self.
_prep_sky(associations, bands)
723 if self.config.doAstrometry:
727 referenceSelector=self.astrometryReferenceSelector,
731 astrometry_output = self.
_make_output(associations.getCcdImageList(),
735 astrometry_output =
None
737 if self.config.doPhotometry:
741 referenceSelector=self.photometryReferenceSelector,
745 reject_bad_fluxes=
True)
746 photometry_output = self.
_make_output(associations.getCcdImageList(),
750 photometry_output =
None
752 return pipeBase.Struct(outputWcs=astrometry_output,
753 outputPhotoCalib=photometry_output,
758 def _load_data(self, inputSourceTableVisit, inputVisitSummary, associations,
759 jointcalControl, camera):
760 """Read the data that jointcal needs to run.
762 Modifies ``associations`` in-place
with the loaded data.
766 inputSourceTableVisit : `list` [`lsst.daf.butler.DeferredDatasetHandle`]
767 References to visit-level DataFrames to load the catalog data
from.
768 inputVisitSummary : `list` [`lsst.daf.butler.DeferredDatasetHandle`]
769 Visit-level exposure summary catalog
with metadata.
771 Object to add the loaded data to by constructing new CcdImages.
773 Control object
for C++ associations management.
775 Camera object
for detector geometry.
780 The original WCS of the input data, to aid
in writing tests.
781 bands : `list` [`str`]
782 The filter bands of each input dataset.
786 load_cat_profile_file = 'jointcal_load_data.prof' if self.config.detailedProfile
else ''
787 with lsst.utils.timer.profile(load_cat_profile_file):
790 catalogMap = {ref.dataId[
'visit']: i
for i, ref
in enumerate(inputSourceTableVisit)}
791 detectorDict = {detector.getId(): detector
for detector
in camera}
795 for visitSummaryRef
in inputVisitSummary:
796 visitSummary = visitSummaryRef.get()
798 dataRef = inputSourceTableVisit[catalogMap[visitSummaryRef.dataId[
'visit']]]
800 inColumns = dataRef.get(component=
'columns')
804 visitCatalog = dataRef.get(parameters={
'columns': columns})
806 selected = self.sourceSelector.
run(visitCatalog)
807 if len(selected.sourceCat) == 0:
808 self.log.warning(
"No sources selected in visit %s. Skipping...",
809 visitSummary[
"visit"][0])
813 detectors = {id: index
for index, id
in enumerate(visitSummary[
'id'])}
814 for id, index
in detectors.items():
819 self.config.sourceFluxType,
825 if result
is not None:
826 oldWcsList.append(result.wcs)
828 filters.append(data.filter)
829 if len(filters) == 0:
830 raise RuntimeError(
"No data to process: did source selector remove all sources?")
831 filters = collections.Counter(filters)
833 return oldWcsList, filters
835 def _make_one_input_data(self, visitRecord, catalog, detectorDict):
836 """Return a data structure for this detector+visit."""
839 visitInfo=visitRecord.getVisitInfo(),
840 detector=detectorDict[visitRecord.getId()],
841 photoCalib=visitRecord.getPhotoCalib(),
842 wcs=visitRecord.getWcs(),
843 bbox=visitRecord.getBBox(),
847 physical=visitRecord[
'physical_filter']))
849 def _build_ccdImage(self, data, associations, jointcalControl):
851 Extract the necessary things from this catalog+metadata to add a new
856 data : `JointcalInputData`
857 The loaded input data.
859 Object to add the info to, to construct a new CcdImage
861 Control object
for associations management
867 The TAN WCS of this image, read
from the calexp
870 A key to identify this dataRef by its visit
and ccd ids
873 if there are no sources
in the loaded catalog.
875 if len(data.catalog) == 0:
876 self.log.warning(
"No sources selected in visit %s ccd %s", data.visit, data.detector.getId())
879 associations.createCcdImage(data.catalog,
883 data.filter.physicalLabel,
887 data.detector.getId(),
890 Result = collections.namedtuple(
'Result_from_build_CcdImage', (
'wcs',
'key'))
891 Key = collections.namedtuple(
'Key', (
'visit',
'ccd'))
892 return Result(data.wcs, Key(data.visit, data.detector.getId()))
894 def _getDebugPath(self, filename):
895 """Constructs a path to filename using the configured debug path.
897 return os.path.join(self.config.debugOutputPath, filename)
899 def _prep_sky(self, associations, filters):
900 """Prepare on-sky and other data that must be computed after data has
903 associations.computeCommonTangentPoint()
905 boundingCircle = associations.computeBoundingCircle()
907 radius = lsst.geom.Angle(boundingCircle.getOpeningAngle().asRadians(), lsst.geom.radians)
909 self.log.info(f"Data has center={center} with radius={radius.asDegrees()} degrees.")
912 defaultFilter = filters.most_common(1)[0][0]
913 self.log.debug(
"Using '%s' filter for reference flux", defaultFilter.physicalLabel)
917 associations.setEpoch(epoch.jyear)
919 return boundingCircle, center, radius, defaultFilter, epoch
921 def _get_refcat_coordinate_error_override(self, refCat, name):
922 """Check whether we should override the refcat coordinate errors, and
923 return the overridden error
if necessary.
928 The reference catalog to check
for a ``coord_raErr`` field.
930 Whether we are doing
"astrometry" or "photometry".
934 refCoordErr : `float`
935 The refcat coordinate error to use,
or NaN
if we are
not overriding
940 lsst.pex.config.FieldValidationError
941 Raised
if the refcat does
not contain coordinate errors
and
942 ``config.astrometryReferenceErr``
is not set.
946 if name.lower() ==
"photometry":
947 if 'coord_raErr' not in refCat.schema:
952 if self.config.astrometryReferenceErr
is None and 'coord_raErr' not in refCat.schema:
953 msg = (
"Reference catalog does not contain coordinate errors, "
954 "and config.astrometryReferenceErr not supplied.")
955 raise pexConfig.FieldValidationError(JointcalConfig.astrometryReferenceErr,
959 if self.config.astrometryReferenceErr
is not None and 'coord_raErr' in refCat.schema:
960 self.log.warning(
"Overriding reference catalog coordinate errors with %f/coordinate [mas]",
961 self.config.astrometryReferenceErr)
963 if self.config.astrometryReferenceErr
is None:
966 return self.config.astrometryReferenceErr
968 def _compute_proper_motion_epoch(self, ccdImageList):
969 """Return the proper motion correction epoch of the provided images.
974 The images to compute the appropriate epoch for.
978 epoch : `astropy.time.Time`
979 The date to use
for proper motion corrections.
981 return astropy.time.Time(np.mean([ccdImage.getEpoch()
for ccdImage
in ccdImageList]),
985 def _do_load_refcat_and_fit(self, associations, defaultFilter, center, radius,
986 tract="", match_cut=3.0,
987 reject_bad_fluxes=False, *,
988 name="", refObjLoader=None, referenceSelector=None,
989 fit_function=None, epoch=None):
990 """Load reference catalog, perform the fit, and return the result.
995 The star/reference star associations to fit.
997 filter to load from reference catalog.
999 ICRS center of field to load
from reference catalog.
1001 On-sky radius to load
from reference catalog.
1003 Name of thing being fit:
"astrometry" or "photometry".
1004 refObjLoader : `lsst.meas.algorithms.ReferenceObjectLoader`
1005 Reference object loader to use to load a reference catalog.
1006 referenceSelector : `lsst.meas.algorithms.ReferenceSourceSelectorTask`
1007 Selector to use to pick objects
from the loaded reference catalog.
1008 fit_function : callable
1009 Function to call to perform fit (takes Associations object).
1010 tract : `str`, optional
1011 Name of tract currently being fit.
1012 match_cut : `float`, optional
1013 Radius
in arcseconds to find cross-catalog matches to during
1014 associations.associateCatalogs.
1015 reject_bad_fluxes : `bool`, optional
1016 Reject refCat sources
with NaN/inf flux
or NaN/0 fluxErr.
1017 epoch : `astropy.time.Time`, optional
1018 Epoch to which to correct refcat proper motion
and parallax,
1019 or `
None` to
not apply such corrections.
1023 result : `Photometry`
or `Astrometry`
1024 Result of `fit_function()`
1026 self.log.info("====== Now processing %s...", name)
1029 associations.associateCatalogs(match_cut)
1031 associations.fittedStarListSize())
1033 applyColorterms =
False if name.lower() ==
"astrometry" else self.config.applyColorTerms
1035 center, radius, defaultFilter,
1036 applyColorterms=applyColorterms,
1040 associations.collectRefStars(refCat,
1041 self.config.matchCut*lsst.geom.arcseconds,
1043 refCoordinateErr=refCoordErr,
1044 rejectBadFluxes=reject_bad_fluxes)
1046 associations.refStarListSize())
1048 associations.prepareFittedStars(self.config.minMeasurements)
1052 associations.nFittedStarsWithAssociatedRefStar())
1054 associations.fittedStarListSize())
1056 associations.nCcdImagesValidForFit())
1058 fit_profile_file =
'jointcal_fit_%s.prof'%name
if self.config.detailedProfile
else ''
1059 dataName =
"{}_{}".format(tract, defaultFilter.physicalLabel)
1060 with lsst.utils.timer.profile(fit_profile_file):
1061 result = fit_function(associations, dataName)
1064 if self.config.writeChi2FilesInitialFinal:
1065 baseName = self.
_getDebugPath(f
"{name}_final_chi2-{dataName}")
1066 result.fit.saveChi2Contributions(baseName+
"{type}")
1067 self.log.info(
"Wrote chi2 contributions files: %s", baseName)
1071 def _load_reference_catalog(self, refObjLoader, referenceSelector, center, radius, filterLabel,
1072 applyColorterms=False, epoch=None):
1073 """Load the necessary reference catalog sources, convert fluxes to
1074 correct units, and apply color term corrections
if requested.
1078 refObjLoader : `lsst.meas.algorithms.ReferenceObjectLoader`
1079 The reference catalog loader to use to get the data.
1080 referenceSelector : `lsst.meas.algorithms.ReferenceSourceSelectorTask`
1081 Source selector to apply to loaded reference catalog.
1083 The center around which to load sources.
1085 The radius around ``center`` to load sources
in.
1087 The camera filter to load fluxes
for.
1088 applyColorterms : `bool`
1089 Apply colorterm corrections to the refcat
for ``filterName``?
1090 epoch : `astropy.time.Time`, optional
1091 Epoch to which to correct refcat proper motion
and parallax,
1092 or `
None` to
not apply such corrections.
1097 The loaded reference catalog.
1099 The name of the reference catalog flux field appropriate
for ``filterName``.
1101 skyCircle = refObjLoader.loadSkyCircle(center,
1103 filterLabel.bandLabel,
1106 selected = referenceSelector.run(skyCircle.refCat)
1108 if not selected.sourceCat.isContiguous():
1109 refCat = selected.sourceCat.copy(deep=
True)
1111 refCat = selected.sourceCat
1114 refCatName = refObjLoader.name
1115 self.log.info(
"Applying color terms for physical filter=%r reference catalog=%s",
1116 filterLabel.physicalLabel, refCatName)
1117 colorterm = self.config.colorterms.getColorterm(filterLabel.physicalLabel,
1121 refMag, refMagErr = colorterm.getCorrectedMagnitudes(refCat)
1122 refCat[skyCircle.fluxField] = u.Magnitude(refMag, u.ABmag).to_value(u.nJy)
1124 refCat[skyCircle.fluxField+
'Err'] = fluxErrFromABMagErr(refMagErr, refMag) * 1e9
1126 return refCat, skyCircle.fluxField
1128 def _check_star_lists(self, associations, name):
1130 if associations.nCcdImagesValidForFit() == 0:
1131 raise RuntimeError(
'No images in the ccdImageList!')
1132 if associations.fittedStarListSize() == 0:
1133 raise RuntimeError(
'No stars in the {} fittedStarList!'.format(name))
1134 if associations.refStarListSize() == 0:
1135 raise RuntimeError(
'No stars in the {} reference star list!'.format(name))
1137 def _logChi2AndValidate(self, associations, fit, model, chi2Label, writeChi2Name=None):
1138 """Compute chi2, log it, validate the model, and return chi2.
1143 The star/reference star associations to fit.
1145 The fitter to use for minimization.
1146 model : `lsst.jointcal.Model`
1147 The model being fit.
1149 Label to describe the chi2 (e.g.
"Initialized",
"Final").
1150 writeChi2Name : `str`, optional
1151 Filename prefix to write the chi2 contributions to.
1152 Do
not supply an extension: an appropriate one will be added.
1157 The chi2 object
for the current fitter
and model.
1162 Raised
if chi2
is infinite
or NaN.
1164 Raised
if the model
is not valid.
1166 if writeChi2Name
is not None:
1168 fit.saveChi2Contributions(fullpath+
"{type}")
1169 self.log.info(
"Wrote chi2 contributions files: %s", fullpath)
1171 chi2 = fit.computeChi2()
1172 self.log.info(
"%s %s", chi2Label, chi2)
1174 if not np.isfinite(chi2.chi2):
1175 raise FloatingPointError(f
'{chi2Label} chi2 is invalid: {chi2}')
1176 if not model.validate(associations.getCcdImageList(), chi2.ndof):
1177 raise ValueError(
"Model is not valid: check log messages for warnings.")
1180 def _fit_photometry(self, associations, dataName=None):
1182 Fit the photometric data.
1187 The star/reference star associations to fit.
1189 Name of the data being processed (e.g. "1234_HSC-Y"),
for
1190 identifying debugging files.
1194 fit_result : `namedtuple`
1196 The photometric fitter used to perform the fit.
1198 The photometric model that was fit.
1200 self.log.info("=== Starting photometric fitting...")
1203 if self.config.photometryModel ==
"constrainedFlux":
1206 visitOrder=self.config.photometryVisitOrder,
1207 errorPedestal=self.config.photometryErrorPedestal)
1209 doLineSearch = self.config.allowLineSearch
1210 elif self.config.photometryModel ==
"constrainedMagnitude":
1213 visitOrder=self.config.photometryVisitOrder,
1214 errorPedestal=self.config.photometryErrorPedestal)
1216 doLineSearch = self.config.allowLineSearch
1217 elif self.config.photometryModel ==
"simpleFlux":
1219 errorPedestal=self.config.photometryErrorPedestal)
1220 doLineSearch =
False
1221 elif self.config.photometryModel ==
"simpleMagnitude":
1223 errorPedestal=self.config.photometryErrorPedestal)
1224 doLineSearch =
False
1229 if self.config.writeChi2FilesInitialFinal:
1230 baseName = f
"photometry_initial_chi2-{dataName}"
1233 if self.config.writeInitialModel:
1234 fullpath = self.
_getDebugPath(f
"initial_photometry_model-{dataName}.txt")
1238 def getChi2Name(whatToFit):
1239 if self.config.writeChi2FilesOuterLoop:
1240 return f
"photometry_init-%s_chi2-{dataName}" % whatToFit
1246 if self.config.writeInitMatrix:
1247 dumpMatrixFile = self.
_getDebugPath(f
"photometry_preinit-{dataName}")
1250 if self.config.photometryModel.startswith(
"constrained"):
1253 fit.minimize(
"ModelVisit", dumpMatrixFile=dumpMatrixFile)
1255 writeChi2Name=getChi2Name(
"ModelVisit"))
1258 fit.minimize(
"Model", doLineSearch=doLineSearch, dumpMatrixFile=dumpMatrixFile)
1260 writeChi2Name=getChi2Name(
"Model"))
1262 fit.minimize(
"Fluxes")
1264 writeChi2Name=getChi2Name(
"Fluxes"))
1266 fit.minimize(
"Model Fluxes", doLineSearch=doLineSearch)
1268 writeChi2Name=getChi2Name(
"ModelFluxes"))
1270 model.freezeErrorTransform()
1271 self.log.debug(
"Photometry error scales are frozen.")
1275 self.config.maxPhotometrySteps,
1278 doRankUpdate=self.config.photometryDoRankUpdate,
1279 doLineSearch=doLineSearch,
1286 def _fit_astrometry(self, associations, dataName=None):
1288 Fit the astrometric data.
1293 The star/reference star associations to fit.
1295 Name of the data being processed (e.g. "1234_HSC-Y"),
for
1296 identifying debugging files.
1300 fit_result : `namedtuple`
1302 The astrometric fitter used to perform the fit.
1304 The astrometric model that was fit.
1306 The model
for the sky to tangent plane projection that was used
in the fit.
1309 self.log.info("=== Starting astrometric fitting...")
1311 associations.deprojectFittedStars()
1318 if self.config.astrometryModel ==
"constrained":
1320 sky_to_tan_projection,
1321 chipOrder=self.config.astrometryChipOrder,
1322 visitOrder=self.config.astrometryVisitOrder)
1323 elif self.config.astrometryModel ==
"simple":
1325 sky_to_tan_projection,
1326 self.config.useInputWcs,
1328 order=self.config.astrometrySimpleOrder)
1333 if self.config.writeChi2FilesInitialFinal:
1334 baseName = f
"astrometry_initial_chi2-{dataName}"
1337 if self.config.writeInitialModel:
1338 fullpath = self.
_getDebugPath(f
"initial_astrometry_model-{dataName}.txt")
1342 def getChi2Name(whatToFit):
1343 if self.config.writeChi2FilesOuterLoop:
1344 return f
"astrometry_init-%s_chi2-{dataName}" % whatToFit
1348 if self.config.writeInitMatrix:
1349 dumpMatrixFile = self.
_getDebugPath(f
"astrometry_preinit-{dataName}")
1354 if self.config.astrometryModel ==
"constrained":
1355 fit.minimize(
"DistortionsVisit", dumpMatrixFile=dumpMatrixFile)
1357 writeChi2Name=getChi2Name(
"DistortionsVisit"))
1360 fit.minimize(
"Distortions", dumpMatrixFile=dumpMatrixFile)
1362 writeChi2Name=getChi2Name(
"Distortions"))
1364 fit.minimize(
"Positions")
1366 writeChi2Name=getChi2Name(
"Positions"))
1368 fit.minimize(
"Distortions Positions")
1370 writeChi2Name=getChi2Name(
"DistortionsPositions"))
1374 self.config.maxAstrometrySteps,
1376 "Distortions Positions",
1377 sigmaRelativeTolerance=self.config.astrometryOutlierRelativeTolerance,
1378 doRankUpdate=self.config.astrometryDoRankUpdate,
1384 return Astrometry(fit, model, sky_to_tan_projection)
1386 def _check_stars(self, associations):
1387 """Count measured and reference stars per ccd and warn/log them."""
1388 for ccdImage
in associations.getCcdImageList():
1389 nMeasuredStars, nRefStars = ccdImage.countStars()
1390 self.log.debug(
"ccdImage %s has %s measured and %s reference stars",
1391 ccdImage.getName(), nMeasuredStars, nRefStars)
1392 if nMeasuredStars < self.config.minMeasuredStarsPerCcd:
1393 self.log.warning(
"ccdImage %s has only %s measuredStars (desired %s)",
1394 ccdImage.getName(), nMeasuredStars, self.config.minMeasuredStarsPerCcd)
1395 if nRefStars < self.config.minRefStarsPerCcd:
1396 self.log.warning(
"ccdImage %s has only %s RefStars (desired %s)",
1397 ccdImage.getName(), nRefStars, self.config.minRefStarsPerCcd)
1399 def _iterate_fit(self, associations, fitter, max_steps, name, whatToFit,
1401 sigmaRelativeTolerance=0,
1403 doLineSearch=False):
1404 """Run fitter.minimize up to max_steps times, returning the final chi2.
1409 The star/reference star associations to fit.
1411 The fitter to use for minimization.
1413 Maximum number of steps to run outlier rejection before declaring
1414 convergence failure.
1415 name : {
'photometry' or 'astrometry'}
1416 What type of data are we fitting (
for logs
and debugging files).
1418 Passed to ``fitter.minimize()`` to define the parameters to fit.
1419 dataName : `str`, optional
1420 Descriptive name
for this dataset (e.g. tract
and filter),
1422 sigmaRelativeTolerance : `float`, optional
1423 Convergence tolerance
for the fractional change
in the chi2 cut
1424 level
for determining outliers. If set to zero, iterations will
1425 continue until there are no outliers.
1426 doRankUpdate : `bool`, optional
1427 Do an Eigen rank update during minimization,
or recompute the full
1428 matrix
and gradient?
1429 doLineSearch : `bool`, optional
1430 Do a line search
for the optimum step during minimization?
1435 The final chi2 after the fit converges,
or is forced to end.
1440 Raised
if the fitter fails
with a non-finite value.
1442 Raised
if the fitter fails
for some other reason;
1443 log messages will provide further details.
1445 if self.config.writeInitMatrix:
1446 dumpMatrixFile = self.
_getDebugPath(f
"{name}_postinit-{dataName}")
1450 oldChi2.chi2 = float(
"inf")
1451 for i
in range(max_steps):
1452 if self.config.writeChi2FilesOuterLoop:
1453 writeChi2Name = f
"{name}_iterate_{i}_chi2-{dataName}"
1455 writeChi2Name =
None
1456 result = fitter.minimize(whatToFit,
1457 self.config.outlierRejectSigma,
1458 sigmaRelativeTolerance=sigmaRelativeTolerance,
1459 doRankUpdate=doRankUpdate,
1460 doLineSearch=doLineSearch,
1461 dumpMatrixFile=dumpMatrixFile)
1464 f
"Fit iteration {i}", writeChi2Name=writeChi2Name)
1466 if result == MinimizeResult.Converged:
1468 self.log.debug(
"fit has converged - no more outliers - redo minimization "
1469 "one more time in case we have lost accuracy in rank update.")
1471 result = fitter.minimize(whatToFit, self.config.outlierRejectSigma,
1472 sigmaRelativeTolerance=sigmaRelativeTolerance)
1476 if chi2.chi2/chi2.ndof >= 4.0:
1477 self.log.error(
"Potentially bad fit: High chi-squared/ndof.")
1480 elif result == MinimizeResult.Chi2Increased:
1481 self.log.warning(
"Still some outliers remaining but chi2 increased - retry")
1483 chi2Ratio = chi2.chi2 / oldChi2.chi2
1485 self.log.warning(
'Significant chi2 increase by a factor of %.4g / %.4g = %.4g',
1486 chi2.chi2, oldChi2.chi2, chi2Ratio)
1493 msg = (
"Large chi2 increase between steps: fit likely cannot converge."
1494 " Try setting one or more of the `writeChi2*` config fields and looking"
1495 " at how individual star chi2-values evolve during the fit.")
1496 raise RuntimeError(msg)
1498 elif result == MinimizeResult.NonFinite:
1499 filename = self.
_getDebugPath(
"{}_failure-nonfinite_chi2-{}.csv".format(name, dataName))
1501 fitter.saveChi2Contributions(filename+
"{type}")
1502 msg =
"Nonfinite value in chi2 minimization, cannot complete fit. Dumped star tables to: {}"
1503 raise FloatingPointError(msg.format(filename))
1504 elif result == MinimizeResult.Failed:
1505 raise RuntimeError(
"Chi2 minimization failure, cannot complete fit.")
1507 raise RuntimeError(
"Unxepected return code from minimize().")
1509 self.log.error(
"%s failed to converge after %d steps"%(name, max_steps))
1513 def _make_output(self, ccdImageList, model, func):
1514 """Return the internal jointcal models converted to the afw
1515 structures that will be saved to disk.
1520 The list of CcdImages to get the output for.
1524 The name of the function to call on ``model`` to get the converted
1531 The data to be saved, keyed on (visit, detector).
1534 for ccdImage
in ccdImageList:
1535 ccd = ccdImage.ccdId
1536 visit = ccdImage.visit
1537 self.log.debug(
"%s for visit: %d, ccd: %d", func, visit, ccd)
1538 output[(visit, ccd)] = getattr(model, func)(ccdImage)
1543 """Return an afw SourceTable to use as a base for creating the
1544 SourceCatalog to insert values from the dataFrame into.
1549 Table
with schema
and slots to use to make SourceCatalogs.
1552 schema.addField("centroid_x",
"D")
1553 schema.addField(
"centroid_y",
"D")
1554 schema.addField(
"centroid_xErr",
"F")
1555 schema.addField(
"centroid_yErr",
"F")
1556 schema.addField(
"shape_xx",
"D")
1557 schema.addField(
"shape_yy",
"D")
1558 schema.addField(
"shape_xy",
"D")
1559 schema.addField(
"flux_instFlux",
"D")
1560 schema.addField(
"flux_instFluxErr",
"D")
1562 table.defineCentroid(
"centroid")
1563 table.defineShape(
"shape")
1569 Get the sourceTable_visit columns to load from the catalogs.
1574 List of columns known to be available
in the sourceTable_visit.
1575 config : `JointcalConfig`
1576 A filled-
in config to to help define column names.
1577 sourceSelector : `lsst.meas.algorithms.BaseSourceSelectorTask`
1578 A configured source selector to define column names to load.
1583 List of columns to read
from sourceTable_visit.
1585 Name of the ixx/iyy/ixy columns.
1587 columns = ['visit',
'detector',
1588 'sourceId',
'x',
'xErr',
'y',
'yErr',
1589 config.sourceFluxType +
'_instFlux', config.sourceFluxType +
'_instFluxErr']
1591 if 'ixx' in inColumns:
1593 ixxColumns = [
'ixx',
'iyy',
'ixy']
1596 ixxColumns = [
'Ixx',
'Iyy',
'Ixy']
1597 columns.extend(ixxColumns)
1599 if sourceSelector.config.doFlags:
1600 columns.extend(sourceSelector.config.flags.bad)
1601 if sourceSelector.config.doUnresolved:
1602 columns.append(sourceSelector.config.unresolved.name)
1603 if sourceSelector.config.doIsolated:
1604 columns.append(sourceSelector.config.isolated.parentName)
1605 columns.append(sourceSelector.config.isolated.nChildName)
1606 if sourceSelector.config.doRequireFiniteRaDec:
1607 columns.append(sourceSelector.config.requireFiniteRaDec.raColName)
1608 columns.append(sourceSelector.config.requireFiniteRaDec.decColName)
1609 if sourceSelector.config.doRequirePrimary:
1610 columns.append(sourceSelector.config.requirePrimary.primaryColName)
1612 return columns, ixxColumns
1616 ixxColumns, sourceFluxType, log):
1617 """Return an afw SourceCatalog extracted from a visit-level dataframe,
1618 limited to just one detector.
1623 Table factory to use to make the SourceCatalog that will be
1624 populated with data
from ``visitCatalog``.
1625 visitCatalog : `pandas.DataFrame`
1626 DataFrame to extract a detector catalog
from.
1628 Numeric id of the detector to extract
from ``visitCatalog``.
1629 ixxColumns : `list` [`str`]
1630 Names of the ixx/iyy/ixy columns
in the catalog.
1631 sourceFluxType : `str`
1632 Name of the catalog field to load instFluxes
from.
1633 log : `logging.Logger`
1634 Logging instance to log to.
1639 Detector-level catalog extracted
from ``visitCatalog``,
or `
None`
1640 if there was no data to load.
1643 mapping = {
'x':
'centroid_x',
1645 'xErr':
'centroid_xErr',
1646 'yErr':
'centroid_yErr',
1647 ixxColumns[0]:
'shape_xx',
1648 ixxColumns[1]:
'shape_yy',
1649 ixxColumns[2]:
'shape_xy',
1650 f
'{sourceFluxType}_instFlux':
'flux_instFlux',
1651 f
'{sourceFluxType}_instFluxErr':
'flux_instFluxErr',
1655 matched = visitCatalog[
'detector'] == detectorId
1659 catalog.resize(sum(matched))
1660 view = visitCatalog.loc[matched]
1661 catalog[
'id'] = view.index
1662 for dfCol, afwCol
in mapping.items():
1663 catalog[afwCol] = view[dfCol]
1665 log.debug(
"%d sources selected in visit %d detector %d",
1667 view[
'visit'].iloc[0],
static Schema makeMinimalSchema()
static std::shared_ptr< SourceTable > make(Schema const &schema, std::shared_ptr< IdFactory > const &idFactory)
static Schema makeMinimalSchema()
The class that implements the relations between MeasuredStar and FittedStar.
Class that handles the astrometric least squares problem.
Interface between AstrometryFit and the combinations of Mappings from pixels to some tangent plane (a...
Handler of an actual image from a single CCD.
Base class for Chi2Statistic and Chi2List, to allow addEntry inside Fitter for either class.
Simple structure to accumulate chi2 and ndof.
A multi-component model, fitting mappings for sensors and visits simultaneously.
A projection handler in which all CCDs from the same visit have the same tangent point.
Class that handles the photometric least squares problem.
A model where there is one independent transform per CcdImage.
def __init__(self, *config=None)
def _load_data(self, inputSourceTableVisit, inputVisitSummary, associations, jointcalControl, camera)
def _compute_proper_motion_epoch(self, ccdImageList)
def _get_refcat_coordinate_error_override(self, refCat, name)
def _getDebugPath(self, filename)
def _check_star_lists(self, associations, name)
def _make_one_input_data(self, visitRecord, catalog, detectorDict)
def _iterate_fit(self, associations, fitter, max_steps, name, whatToFit, dataName="", sigmaRelativeTolerance=0, doRankUpdate=True, doLineSearch=False)
def _check_stars(self, associations)
def _prep_sky(self, associations, filters)
def runQuantum(self, butlerQC, inputRefs, outputRefs)
def _do_load_refcat_and_fit(self, associations, defaultFilter, center, radius, tract="", match_cut=3.0, reject_bad_fluxes=False, *name="", refObjLoader=None, referenceSelector=None, fit_function=None, epoch=None)
def _put_metrics(self, butlerQC, job, outputRefs)
def _build_ccdImage(self, data, associations, jointcalControl)
def __init__(self, **kwargs)
def run(self, inputSourceTableVisit, inputVisitSummary, inputCamera, tract=None)
def _fit_photometry(self, associations, dataName=None)
def _logChi2AndValidate(self, associations, fit, model, chi2Label, writeChi2Name=None)
def _make_output(self, ccdImageList, model, func)
def _put_output(self, butlerQC, outputs, outputRefs, camera, setter)
def _load_reference_catalog(self, refObjLoader, referenceSelector, center, radius, filterLabel, applyColorterms=False, epoch=None)
def _fit_astrometry(self, associations, dataName=None)
def extract_detector_catalog_from_visit_catalog(table, visitCatalog, detectorId, ixxColumns, sourceFluxType, log)
def lookupStaticCalibrations(datasetType, registry, quantumDataId, collections)
def add_measurement(job, name, value)
def get_sourceTable_visit_columns(inColumns, config, sourceSelector)
def lookupVisitRefCats(datasetType, registry, quantumDataId, collections)
def writeModel(model, filename, log)
This is a virtual class that allows a lot of freedom in the choice of the projection from "Sky" (wher...