lsst.jointcal  7.5-hsc+1
jointcal.py
Go to the documentation of this file.
1 # This file is part of jointcal.
2 #
3 # Developed for the LSST Data Management System.
4 # This product includes software developed by the LSST Project
5 # (https://www.lsst.org).
6 # See the COPYRIGHT file at the top-level directory of this distribution
7 # for details of code ownership.
8 #
9 # This program is free software: you can redistribute it and/or modify
10 # it under the terms of the GNU General Public License as published by
11 # the Free Software Foundation, either version 3 of the License, or
12 # (at your option) any later version.
13 #
14 # This program is distributed in the hope that it will be useful,
15 # but WITHOUT ANY WARRANTY; without even the implied warranty of
16 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 # GNU General Public License for more details.
18 #
19 # You should have received a copy of the GNU General Public License
20 # along with this program. If not, see <https://www.gnu.org/licenses/>.
21 
22 import collections
23 import numpy as np
24 import astropy.units as u
25 
26 import lsst.geom
27 import lsst.utils
28 import lsst.pex.config as pexConfig
29 import lsst.pipe.base as pipeBase
30 from lsst.afw.image import fluxErrFromABMagErr
31 import lsst.pex.exceptions as pexExceptions
32 import lsst.afw.table
34 from lsst.pipe.tasks.colorterms import ColortermLibrary
35 from lsst.verify import Job, Measurement
36 
37 from lsst.meas.algorithms import LoadIndexedReferenceObjectsTask, ReferenceSourceSelectorTask
38 from lsst.meas.algorithms.sourceSelector import sourceSelectorRegistry
39 
40 from .dataIds import PerTractCcdDataIdContainer
41 
42 import lsst.jointcal
43 from lsst.jointcal import MinimizeResult
44 
45 __all__ = ["JointcalConfig", "JointcalRunner", "JointcalTask"]
46 
47 Photometry = collections.namedtuple('Photometry', ('fit', 'model'))
48 Astrometry = collections.namedtuple('Astrometry', ('fit', 'model', 'sky_to_tan_projection'))
49 
50 
51 # TODO: move this to MeasurementSet in lsst.verify per DM-12655.
52 def add_measurement(job, name, value):
53  meas = Measurement(job.metrics[name], value)
54  job.measurements.insert(meas)
55 
56 
57 class JointcalRunner(pipeBase.ButlerInitializedTaskRunner):
58  """Subclass of TaskRunner for jointcalTask
59 
60  jointcalTask.runDataRef() takes a number of arguments, one of which is a list of dataRefs
61  extracted from the command line (whereas most CmdLineTasks' runDataRef methods take
62  single dataRef, are are called repeatedly). This class transforms the processed
63  arguments generated by the ArgumentParser into the arguments expected by
64  Jointcal.runDataRef().
65 
66  See pipeBase.TaskRunner for more information.
67  """
68 
69  @staticmethod
70  def getTargetList(parsedCmd, **kwargs):
71  """
72  Return a list of tuples per tract, each containing (dataRefs, kwargs).
73 
74  Jointcal operates on lists of dataRefs simultaneously.
75  """
76  kwargs['profile_jointcal'] = parsedCmd.profile_jointcal
77  kwargs['butler'] = parsedCmd.butler
78 
79  # organize data IDs by tract
80  refListDict = {}
81  for ref in parsedCmd.id.refList:
82  refListDict.setdefault(ref.dataId["tract"], []).append(ref)
83  # we call runDataRef() once with each tract
84  result = [(refListDict[tract], kwargs) for tract in sorted(refListDict.keys())]
85  return result
86 
87  def __call__(self, args):
88  """
89  Parameters
90  ----------
91  args
92  Arguments for Task.runDataRef()
93 
94  Returns
95  -------
96  pipe.base.Struct
97  if self.doReturnResults is False:
98 
99  - ``exitStatus``: 0 if the task completed successfully, 1 otherwise.
100 
101  if self.doReturnResults is True:
102 
103  - ``result``: the result of calling jointcal.runDataRef()
104  - ``exitStatus``: 0 if the task completed successfully, 1 otherwise.
105  """
106  exitStatus = 0 # exit status for shell
107 
108  # NOTE: cannot call self.makeTask because that assumes args[0] is a single dataRef.
109  dataRefList, kwargs = args
110  butler = kwargs.pop('butler')
111  task = self.TaskClass(config=self.config, log=self.log, butler=butler)
112  result = None
113  try:
114  result = task.runDataRef(dataRefList, **kwargs)
115  exitStatus = result.exitStatus
116  job_path = butler.get('verify_job_filename')
117  result.job.write(job_path[0])
118  except Exception as e: # catch everything, sort it out later.
119  if self.doRaise:
120  raise e
121  else:
122  exitStatus = 1
123  eName = type(e).__name__
124  tract = dataRefList[0].dataId['tract']
125  task.log.fatal("Failed processing tract %s, %s: %s", tract, eName, e)
126 
127  # Put the butler back into kwargs for the other Tasks.
128  kwargs['butler'] = butler
129  if self.doReturnResults:
130  return pipeBase.Struct(result=result, exitStatus=exitStatus)
131  else:
132  return pipeBase.Struct(exitStatus=exitStatus)
133 
134 
135 class JointcalConfig(pexConfig.Config):
136  """Configuration for JointcalTask"""
137 
138  doAstrometry = pexConfig.Field(
139  doc="Fit astrometry and write the fitted result.",
140  dtype=bool,
141  default=True
142  )
143  doPhotometry = pexConfig.Field(
144  doc="Fit photometry and write the fitted result.",
145  dtype=bool,
146  default=True
147  )
148  coaddName = pexConfig.Field(
149  doc="Type of coadd, typically deep or goodSeeing",
150  dtype=str,
151  default="deep"
152  )
153  positionErrorPedestal = pexConfig.Field(
154  doc="Systematic term to apply to the measured position error (pixels)",
155  dtype=float,
156  default=0.02,
157  )
158  photometryErrorPedestal = pexConfig.Field(
159  doc="Systematic term to apply to the measured error on flux or magnitude as a "
160  "fraction of source flux or magnitude delta (e.g. 0.05 is 5% of flux or +50 millimag).",
161  dtype=float,
162  default=0.0,
163  )
164  # TODO: DM-6885 matchCut should be an afw.geom.Angle
165  matchCut = pexConfig.Field(
166  doc="Matching radius between fitted and reference stars (arcseconds)",
167  dtype=float,
168  default=3.0,
169  )
170  minMeasurements = pexConfig.Field(
171  doc="Minimum number of associated measured stars for a fitted star to be included in the fit",
172  dtype=int,
173  default=2,
174  )
175  minMeasuredStarsPerCcd = pexConfig.Field(
176  doc="Minimum number of measuredStars per ccdImage before printing warnings",
177  dtype=int,
178  default=100,
179  )
180  minRefStarsPerCcd = pexConfig.Field(
181  doc="Minimum number of measuredStars per ccdImage before printing warnings",
182  dtype=int,
183  default=30,
184  )
185  allowLineSearch = pexConfig.Field(
186  doc="Allow a line search during minimization, if it is reasonable for the model"
187  " (models with a significant non-linear component, e.g. constrainedPhotometry).",
188  dtype=bool,
189  default=False
190  )
191  astrometrySimpleOrder = pexConfig.Field(
192  doc="Polynomial order for fitting the simple astrometry model.",
193  dtype=int,
194  default=3,
195  )
196  astrometryChipOrder = pexConfig.Field(
197  doc="Order of the per-chip transform for the constrained astrometry model.",
198  dtype=int,
199  default=1,
200  )
201  astrometryVisitOrder = pexConfig.Field(
202  doc="Order of the per-visit transform for the constrained astrometry model.",
203  dtype=int,
204  default=5,
205  )
206  useInputWcs = pexConfig.Field(
207  doc="Use the input calexp WCSs to initialize a SimpleAstrometryModel.",
208  dtype=bool,
209  default=True,
210  )
211  astrometryModel = pexConfig.ChoiceField(
212  doc="Type of model to fit to astrometry",
213  dtype=str,
214  default="constrained",
215  allowed={"simple": "One polynomial per ccd",
216  "constrained": "One polynomial per ccd, and one polynomial per visit"}
217  )
218  photometryModel = pexConfig.ChoiceField(
219  doc="Type of model to fit to photometry",
220  dtype=str,
221  default="constrainedMagnitude",
222  allowed={"simpleFlux": "One constant zeropoint per ccd and visit, fitting in flux space.",
223  "constrainedFlux": "Constrained zeropoint per ccd, and one polynomial per visit,"
224  " fitting in flux space.",
225  "simpleMagnitude": "One constant zeropoint per ccd and visit,"
226  " fitting in magnitude space.",
227  "constrainedMagnitude": "Constrained zeropoint per ccd, and one polynomial per visit,"
228  " fitting in magnitude space.",
229  }
230  )
231  applyColorTerms = pexConfig.Field(
232  doc="Apply photometric color terms to reference stars?"
233  "Requires that colorterms be set to a ColortermLibrary",
234  dtype=bool,
235  default=False
236  )
237  colorterms = pexConfig.ConfigField(
238  doc="Library of photometric reference catalog name to color term dict.",
239  dtype=ColortermLibrary,
240  )
241  photometryVisitOrder = pexConfig.Field(
242  doc="Order of the per-visit polynomial transform for the constrained photometry model.",
243  dtype=int,
244  default=7,
245  )
246  photometryDoRankUpdate = pexConfig.Field(
247  doc="Do the rank update step during minimization. "
248  "Skipping this can help deal with models that are too non-linear.",
249  dtype=bool,
250  default=True,
251  )
252  astrometryDoRankUpdate = pexConfig.Field(
253  doc="Do the rank update step during minimization (should not change the astrometry fit). "
254  "Skipping this can help deal with models that are too non-linear.",
255  dtype=bool,
256  default=True,
257  )
258  outlierRejectSigma = pexConfig.Field(
259  doc="How many sigma to reject outliers at during minimization.",
260  dtype=float,
261  default=5.0,
262  )
263  maxPhotometrySteps = pexConfig.Field(
264  doc="Maximum number of minimize iterations to take when fitting photometry.",
265  dtype=int,
266  default=20,
267  )
268  maxAstrometrySteps = pexConfig.Field(
269  doc="Maximum number of minimize iterations to take when fitting photometry.",
270  dtype=int,
271  default=20,
272  )
273  astrometryRefObjLoader = pexConfig.ConfigurableField(
274  target=LoadIndexedReferenceObjectsTask,
275  doc="Reference object loader for astrometric fit",
276  )
277  photometryRefObjLoader = pexConfig.ConfigurableField(
278  target=LoadIndexedReferenceObjectsTask,
279  doc="Reference object loader for photometric fit",
280  )
281  sourceSelector = sourceSelectorRegistry.makeField(
282  doc="How to select sources for cross-matching",
283  default="astrometry"
284  )
285  astrometryReferenceSelector = pexConfig.ConfigurableField(
286  target=ReferenceSourceSelectorTask,
287  doc="How to down-select the loaded astrometry reference catalog.",
288  )
289  photometryReferenceSelector = pexConfig.ConfigurableField(
290  target=ReferenceSourceSelectorTask,
291  doc="How to down-select the loaded photometry reference catalog.",
292  )
293  astrometryReferenceErr = pexConfig.Field(
294  doc="Uncertainty on reference catalog coordinates [mas] to use in place of the `coord_*_err` fields."
295  " If None, then raise an exception if the reference catalog is missing coordinate errors."
296  " If specified, overrides any existing `coord_*_err` values.",
297  dtype=float,
298  default=None,
299  optional=True
300  )
301  writeInitMatrix = pexConfig.Field(
302  dtype=bool,
303  doc="Write the pre/post-initialization Hessian and gradient to text files, for debugging."
304  "The output files will be of the form 'astrometry_preinit-mat.txt', in the current directory."
305  "Note that these files are the dense versions of the matrix, and so may be very large.",
306  default=False
307  )
308  writeChi2FilesInitialFinal = pexConfig.Field(
309  dtype=bool,
310  doc="Write .csv files containing the contributions to chi2 for the initialization and final fit.",
311  default=False
312  )
313  writeChi2FilesOuterLoop = pexConfig.Field(
314  dtype=bool,
315  doc="Write .csv files containing the contributions to chi2 for the outer fit loop.",
316  default=False
317  )
318  sourceFluxType = pexConfig.Field(
319  dtype=str,
320  doc="Source flux field to use in source selection and to get fluxes from the catalog.",
321  default='Calib'
322  )
323 
324  def validate(self):
325  super().validate()
326  if self.applyColorTerms and len(self.colorterms.data) == 0:
327  msg = "applyColorTerms=True requires the `colorterms` field be set to a ColortermLibrary."
328  raise pexConfig.FieldValidationError(JointcalConfig.colorterms, self, msg)
329 
330  def setDefaults(self):
331  # Use science source selector which can filter on extendedness, SNR, and whether blended
332  self.sourceSelector.name = 'science'
333  # Use only stars because aperture fluxes of galaxies are biased and depend on seeing
334  self.sourceSelector['science'].doUnresolved = True
335  # with dependable signal to noise ratio.
336  self.sourceSelector['science'].doSignalToNoise = True
337  # Min SNR must be > 0 because jointcal cannot handle negative fluxes,
338  # and S/N > 10 to use sources that are not too faint, and thus better measured.
339  self.sourceSelector['science'].signalToNoise.minimum = 10.
340  # Base SNR on CalibFlux because that is the flux jointcal that fits and must be positive
341  fluxField = f"slot_{self.sourceFluxType}Flux_instFlux"
342  self.sourceSelector['science'].signalToNoise.fluxField = fluxField
343  self.sourceSelector['science'].signalToNoise.errField = fluxField + "Err"
344  # Do not trust blended sources' aperture fluxes which also depend on seeing
345  self.sourceSelector['science'].doIsolated = True
346  # Do not trust either flux or centroid measurements with flags,
347  # chosen from the usual QA flags for stars)
348  self.sourceSelector['science'].doFlags = True
349  badFlags = ['base_PixelFlags_flag_edge', 'base_PixelFlags_flag_saturated',
350  'base_PixelFlags_flag_interpolatedCenter', 'base_SdssCentroid_flag',
351  'base_PsfFlux_flag', 'base_PixelFlags_flag_suspectCenter']
352  self.sourceSelector['science'].flags.bad = badFlags
353 
354 
355 class JointcalTask(pipeBase.CmdLineTask):
356  """Jointly astrometrically and photometrically calibrate a group of images."""
357 
358  ConfigClass = JointcalConfig
359  RunnerClass = JointcalRunner
360  _DefaultName = "jointcal"
361 
362  def __init__(self, butler=None, profile_jointcal=False, **kwargs):
363  """
364  Instantiate a JointcalTask.
365 
366  Parameters
367  ----------
368  butler : `lsst.daf.persistence.Butler`
369  The butler is passed to the refObjLoader constructor in case it is
370  needed. Ignored if the refObjLoader argument provides a loader directly.
371  Used to initialize the astrometry and photometry refObjLoaders.
372  profile_jointcal : `bool`
373  Set to True to profile different stages of this jointcal run.
374  """
375  pipeBase.CmdLineTask.__init__(self, **kwargs)
376  self.profile_jointcal = profile_jointcal
377  self.makeSubtask("sourceSelector")
378  if self.config.doAstrometry:
379  self.makeSubtask('astrometryRefObjLoader', butler=butler)
380  self.makeSubtask("astrometryReferenceSelector")
381  else:
383  if self.config.doPhotometry:
384  self.makeSubtask('photometryRefObjLoader', butler=butler)
385  self.makeSubtask("photometryReferenceSelector")
386  else:
388 
389  # To hold various computed metrics for use by tests
390  self.job = Job.load_metrics_package(subset='jointcal')
391 
392  # We don't currently need to persist the metadata.
393  # If we do in the future, we will have to add appropriate dataset templates
394  # to each obs package (the metadata template should look like `jointcal_wcs`).
395  def _getMetadataName(self):
396  return None
397 
398  @classmethod
399  def _makeArgumentParser(cls):
400  """Create an argument parser"""
401  parser = pipeBase.ArgumentParser(name=cls._DefaultName)
402  parser.add_argument("--profile_jointcal", default=False, action="store_true",
403  help="Profile steps of jointcal separately.")
404  parser.add_id_argument("--id", "calexp", help="data ID, e.g. --id visit=6789 ccd=0..9",
405  ContainerClass=PerTractCcdDataIdContainer)
406  return parser
407 
408  def _build_ccdImage(self, dataRef, associations, jointcalControl):
409  """
410  Extract the necessary things from this dataRef to add a new ccdImage.
411 
412  Parameters
413  ----------
414  dataRef : `lsst.daf.persistence.ButlerDataRef`
415  DataRef to extract info from.
416  associations : `lsst.jointcal.Associations`
417  Object to add the info to, to construct a new CcdImage
418  jointcalControl : `jointcal.JointcalControl`
419  Control object for associations management
420 
421  Returns
422  ------
423  namedtuple
424  ``wcs``
425  The TAN WCS of this image, read from the calexp
426  (`lsst.afw.geom.SkyWcs`).
427  ``key``
428  A key to identify this dataRef by its visit and ccd ids
429  (`namedtuple`).
430  ``filter``
431  This calexp's filter (`str`).
432  """
433  if "visit" in dataRef.dataId.keys():
434  visit = dataRef.dataId["visit"]
435  else:
436  visit = dataRef.getButler().queryMetadata("calexp", ("visit"), dataRef.dataId)[0]
437 
438  src = dataRef.get("src", flags=lsst.afw.table.SOURCE_IO_NO_FOOTPRINTS, immediate=True)
439 
440  visitInfo = dataRef.get('calexp_visitInfo')
441  detector = dataRef.get('calexp_detector')
442  ccdId = detector.getId()
443  photoCalib = dataRef.get('calexp_photoCalib')
444  tanWcs = dataRef.get('calexp_wcs')
445  bbox = dataRef.get('calexp_bbox')
446  filt = dataRef.get('calexp_filter')
447  filterName = filt.getName()
448 
449  goodSrc = self.sourceSelector.run(src)
450 
451  if len(goodSrc.sourceCat) == 0:
452  self.log.warn("No sources selected in visit %s ccd %s", visit, ccdId)
453  else:
454  self.log.info("%d sources selected in visit %d ccd %d", len(goodSrc.sourceCat), visit, ccdId)
455  associations.createCcdImage(goodSrc.sourceCat,
456  tanWcs,
457  visitInfo,
458  bbox,
459  filterName,
460  photoCalib,
461  detector,
462  visit,
463  ccdId,
464  jointcalControl)
465 
466  Result = collections.namedtuple('Result_from_build_CcdImage', ('wcs', 'key', 'filter'))
467  Key = collections.namedtuple('Key', ('visit', 'ccd'))
468  return Result(tanWcs, Key(visit, ccdId), filterName)
469 
470  @pipeBase.timeMethod
471  def runDataRef(self, dataRefs, profile_jointcal=False):
472  """
473  Jointly calibrate the astrometry and photometry across a set of images.
474 
475  Parameters
476  ----------
477  dataRefs : `list` of `lsst.daf.persistence.ButlerDataRef`
478  List of data references to the exposures to be fit.
479  profile_jointcal : `bool`
480  Profile the individual steps of jointcal.
481 
482  Returns
483  -------
484  result : `lsst.pipe.base.Struct`
485  Struct of metadata from the fit, containing:
486 
487  ``dataRefs``
488  The provided data references that were fit (with updated WCSs)
489  ``oldWcsList``
490  The original WCS from each dataRef
491  ``metrics``
492  Dictionary of internally-computed metrics for testing/validation.
493  """
494  if len(dataRefs) == 0:
495  raise ValueError('Need a non-empty list of data references!')
496 
497  exitStatus = 0 # exit status for shell
498 
499  sourceFluxField = "slot_%sFlux" % (self.config.sourceFluxType,)
500  jointcalControl = lsst.jointcal.JointcalControl(sourceFluxField)
501  associations = lsst.jointcal.Associations()
502 
503  visit_ccd_to_dataRef = {}
504  oldWcsList = []
505  filters = []
506  load_cat_prof_file = 'jointcal_build_ccdImage.prof' if profile_jointcal else ''
507  with pipeBase.cmdLineTask.profile(load_cat_prof_file):
508  # We need the bounding-box of the focal plane for photometry visit models.
509  # NOTE: we only need to read it once, because its the same for all exposures of a camera.
510  camera = dataRefs[0].get('camera', immediate=True)
511  self.focalPlaneBBox = camera.getFpBBox()
512  for ref in dataRefs:
513  result = self._build_ccdImage(ref, associations, jointcalControl)
514  oldWcsList.append(result.wcs)
515  visit_ccd_to_dataRef[result.key] = ref
516  filters.append(result.filter)
517  filters = collections.Counter(filters)
518 
519  associations.computeCommonTangentPoint()
520 
521  boundingCircle = associations.computeBoundingCircle()
522  center = lsst.geom.SpherePoint(boundingCircle.getCenter())
523  radius = lsst.geom.Angle(boundingCircle.getOpeningAngle().asRadians(), lsst.geom.radians)
524 
525  # Determine a default filter associated with the catalog. See DM-9093
526  defaultFilter = filters.most_common(1)[0][0]
527  self.log.debug("Using %s band for reference flux", defaultFilter)
528 
529  # TODO: need a better way to get the tract.
530  tract = dataRefs[0].dataId['tract']
531 
532  if self.config.doAstrometry:
533  astrometry = self._do_load_refcat_and_fit(associations, defaultFilter, center, radius,
534  name="astrometry",
535  refObjLoader=self.astrometryRefObjLoader,
536  referenceSelector=self.astrometryReferenceSelector,
537  fit_function=self._fit_astrometry,
538  profile_jointcal=profile_jointcal,
539  tract=tract)
540  self._write_astrometry_results(associations, astrometry.model, visit_ccd_to_dataRef)
541  else:
542  astrometry = Astrometry(None, None, None)
543 
544  if self.config.doPhotometry:
545  photometry = self._do_load_refcat_and_fit(associations, defaultFilter, center, radius,
546  name="photometry",
547  refObjLoader=self.photometryRefObjLoader,
548  referenceSelector=self.photometryReferenceSelector,
549  fit_function=self._fit_photometry,
550  profile_jointcal=profile_jointcal,
551  tract=tract,
552  filters=filters,
553  reject_bad_fluxes=True)
554  self._write_photometry_results(associations, photometry.model, visit_ccd_to_dataRef)
555  else:
556  photometry = Photometry(None, None)
557 
558  return pipeBase.Struct(dataRefs=dataRefs,
559  oldWcsList=oldWcsList,
560  job=self.job,
561  astrometryRefObjLoader=self.astrometryRefObjLoader,
562  photometryRefObjLoader=self.photometryRefObjLoader,
563  defaultFilter=defaultFilter,
564  exitStatus=exitStatus)
565 
566  def _do_load_refcat_and_fit(self, associations, defaultFilter, center, radius,
567  name="", refObjLoader=None, referenceSelector=None,
568  filters=[], fit_function=None,
569  tract=None, profile_jointcal=False, match_cut=3.0,
570  reject_bad_fluxes=False):
571  """Load reference catalog, perform the fit, and return the result.
572 
573  Parameters
574  ----------
575  associations : `lsst.jointcal.Associations`
576  The star/reference star associations to fit.
577  defaultFilter : `str`
578  filter to load from reference catalog.
579  center : `lsst.afw.geom.SpherePoint`
580  ICRS center of field to load from reference catalog.
581  radius : `lsst.afw.geom.Angle`
582  On-sky radius to load from reference catalog.
583  name : `str`
584  Name of thing being fit: "Astrometry" or "Photometry".
585  refObjLoader : `lsst.meas.algorithms.LoadReferenceObjectsTask`
586  Reference object loader to load from for fit.
587  filters : `list` of `str`, optional
588  List of filters to load from the reference catalog.
589  fit_function : callable
590  Function to call to perform fit (takes associations object).
591  tract : `str`
592  Name of tract currently being fit.
593  profile_jointcal : `bool`, optional
594  Separately profile the fitting step.
595  match_cut : `float`, optional
596  Radius in arcseconds to find cross-catalog matches to during
597  associations.associateCatalogs.
598  reject_bad_fluxes : `bool`, optional
599  Reject refCat sources with NaN/inf flux or NaN/0 fluxErr.
600 
601  Returns
602  -------
603  result : `Photometry` or `Astrometry`
604  Result of `fit_function()`
605  """
606  self.log.info("====== Now processing %s...", name)
607  # TODO: this should not print "trying to invert a singular transformation:"
608  # if it does that, something's not right about the WCS...
609  associations.associateCatalogs(match_cut)
610  add_measurement(self.job, 'jointcal.associated_%s_fittedStars' % name,
611  associations.fittedStarListSize())
612 
613  applyColorterms = False if name == "Astrometry" else self.config.applyColorTerms
614  if name == "Astrometry":
615  referenceSelector = self.config.astrometryReferenceSelector
616  elif name == "Photometry":
617  referenceSelector = self.config.photometryReferenceSelector
618  refCat, fluxField = self._load_reference_catalog(refObjLoader, referenceSelector,
619  center, radius, defaultFilter,
620  applyColorterms=applyColorterms)
621 
622  if self.config.astrometryReferenceErr is None:
623  refCoordErr = float('nan')
624  else:
625  refCoordErr = self.config.astrometryReferenceErr
626 
627  associations.collectRefStars(refCat,
628  self.config.matchCut*lsst.geom.arcseconds,
629  fluxField,
630  refCoordinateErr=refCoordErr,
631  rejectBadFluxes=reject_bad_fluxes)
632  add_measurement(self.job, 'jointcal.collected_%s_refStars' % name,
633  associations.refStarListSize())
634 
635  associations.prepareFittedStars(self.config.minMeasurements)
636 
637  self._check_star_lists(associations, name)
638  add_measurement(self.job, 'jointcal.selected_%s_refStars' % name,
639  associations.nFittedStarsWithAssociatedRefStar())
640  add_measurement(self.job, 'jointcal.selected_%s_fittedStars' % name,
641  associations.fittedStarListSize())
642  add_measurement(self.job, 'jointcal.selected_%s_ccdImages' % name,
643  associations.nCcdImagesValidForFit())
644 
645  load_cat_prof_file = 'jointcal_fit_%s.prof'%name if profile_jointcal else ''
646  dataName = "{}_{}".format(tract, defaultFilter)
647  with pipeBase.cmdLineTask.profile(load_cat_prof_file):
648  result = fit_function(associations, dataName)
649  # TODO DM-12446: turn this into a "butler save" somehow.
650  # Save reference and measurement chi2 contributions for this data
651  if self.config.writeChi2FilesInitialFinal:
652  baseName = f"{name}_final_chi2-{dataName}"
653  result.fit.saveChi2Contributions(baseName+"{type}")
654 
655  return result
656 
657  def _load_reference_catalog(self, refObjLoader, referenceSelector, center, radius, filterName,
658  applyColorterms=False):
659  """Load the necessary reference catalog sources, convert fluxes to
660  correct units, and apply color term corrections if requested.
661 
662  Parameters
663  ----------
664  refObjLoader : `lsst.meas.algorithms.LoadReferenceObjectsTask`
665  The reference catalog loader to use to get the data.
666  referenceSelector : `lsst.meas.algorithms.ReferenceSourceSelectorTask`
667  Source selector to apply to loaded reference catalog.
668  center : `lsst.geom.SpherePoint`
669  The center around which to load sources.
670  radius : `lsst.geom.Angle`
671  The radius around ``center`` to load sources in.
672  filterName : `str`
673  The name of the camera filter to load fluxes for.
674  applyColorterms : `bool`
675  Apply colorterm corrections to the refcat for ``filterName``?
676 
677  Returns
678  -------
679  refCat : `lsst.afw.table.SimpleCatalog`
680  The loaded reference catalog.
681  fluxField : `str`
682  The name of the reference catalog flux field appropriate for ``filterName``.
683  """
684  skyCircle = refObjLoader.loadSkyCircle(center,
685  radius,
686  filterName)
687 
688  selected = referenceSelector.run(skyCircle.refCat)
689  # Need memory contiguity to get reference filters as a vector.
690  if not selected.sourceCat.isContiguous():
691  refCat = selected.sourceCat.copy(deep=True)
692  else:
693  refCat = selected.sourceCat
694 
695  if self.config.astrometryReferenceErr is None and 'coord_ra_err' not in refCat.schema:
696  msg = ("Reference catalog does not contain coordinate errors, "
697  "and config.astrometryReferenceErr not supplied.")
698  raise pexConfig.FieldValidationError(JointcalConfig.astrometryReferenceErr,
699  self.config,
700  msg)
701 
702  if self.config.astrometryReferenceErr is not None and 'coord_ra_err' in refCat.schema:
703  self.log.warn("Overriding reference catalog coordinate errors with %f/coordinate [mas]",
704  self.config.astrometryReferenceErr)
705 
706  if applyColorterms:
707  try:
708  refCatName = refObjLoader.ref_dataset_name
709  except AttributeError:
710  # NOTE: we need this try:except: block in place until we've completely removed a.net support.
711  raise RuntimeError("Cannot perform colorterm corrections with a.net refcats.")
712  self.log.info("Applying color terms for filterName=%r reference catalog=%s",
713  filterName, refCatName)
714  colorterm = self.config.colorterms.getColorterm(
715  filterName=filterName, photoCatName=refCatName, doRaise=True)
716 
717  refMag, refMagErr = colorterm.getCorrectedMagnitudes(refCat, filterName)
718  refCat[skyCircle.fluxField] = u.Magnitude(refMag, u.ABmag).to_value(u.nJy)
719  # TODO: I didn't want to use this, but I'll deal with it in DM-16903
720  refCat[skyCircle.fluxField+'Err'] = fluxErrFromABMagErr(refMagErr, refMag) * 1e9
721 
722  return refCat, skyCircle.fluxField
723 
724  def _check_star_lists(self, associations, name):
725  # TODO: these should be len(blah), but we need this properly wrapped first.
726  if associations.nCcdImagesValidForFit() == 0:
727  raise RuntimeError('No images in the ccdImageList!')
728  if associations.fittedStarListSize() == 0:
729  raise RuntimeError('No stars in the {} fittedStarList!'.format(name))
730  if associations.refStarListSize() == 0:
731  raise RuntimeError('No stars in the {} reference star list!'.format(name))
732 
733  def _logChi2AndValidate(self, associations, fit, model, chi2Label="Model",
734  writeChi2Name=None):
735  """Compute chi2, log it, validate the model, and return chi2.
736 
737  Parameters
738  ----------
739  associations : `lsst.jointcal.Associations`
740  The star/reference star associations to fit.
741  fit : `lsst.jointcal.FitterBase`
742  The fitter to use for minimization.
743  model : `lsst.jointcal.Model`
744  The model being fit.
745  chi2Label : str, optional
746  Label to describe the chi2 (e.g. "Initialized", "Final").
747  writeChi2Name : `str`, optional
748  Filename prefix to write the chi2 contributions to.
749  Do not supply an extension: an appropriate one will be added.
750 
751  Returns
752  -------
753  chi2: `lsst.jointcal.Chi2Accumulator`
754  The chi2 object for the current fitter and model.
755 
756  Raises
757  ------
758  FloatingPointError
759  Raised if chi2 is infinite or NaN.
760  ValueError
761  Raised if the model is not valid.
762  """
763  if writeChi2Name is not None:
764  fit.saveChi2Contributions(writeChi2Name+"{type}")
765  self.log.info("Wrote chi2 contributions files: %s", writeChi2Name)
766 
767  chi2 = fit.computeChi2()
768  self.log.info("%s %s", chi2Label, chi2)
769  self._check_stars(associations)
770  if not np.isfinite(chi2.chi2):
771  raise FloatingPointError(f'{chi2Label} chi2 is invalid: {chi2}')
772  if not model.validate(associations.getCcdImageList(), chi2.ndof):
773  raise ValueError("Model is not valid: check log messages for warnings.")
774  return chi2
775 
776  def _fit_photometry(self, associations, dataName=None):
777  """
778  Fit the photometric data.
779 
780  Parameters
781  ----------
782  associations : `lsst.jointcal.Associations`
783  The star/reference star associations to fit.
784  dataName : `str`
785  Name of the data being processed (e.g. "1234_HSC-Y"), for
786  identifying debugging files.
787 
788  Returns
789  -------
790  fit_result : `namedtuple`
791  fit : `lsst.jointcal.PhotometryFit`
792  The photometric fitter used to perform the fit.
793  model : `lsst.jointcal.PhotometryModel`
794  The photometric model that was fit.
795  """
796  self.log.info("=== Starting photometric fitting...")
797 
798  # TODO: should use pex.config.RegistryField here (see DM-9195)
799  if self.config.photometryModel == "constrainedFlux":
800  model = lsst.jointcal.ConstrainedFluxModel(associations.getCcdImageList(),
801  self.focalPlaneBBox,
802  visitOrder=self.config.photometryVisitOrder,
803  errorPedestal=self.config.photometryErrorPedestal)
804  # potentially nonlinear problem, so we may need a line search to converge.
805  doLineSearch = self.config.allowLineSearch
806  elif self.config.photometryModel == "constrainedMagnitude":
807  model = lsst.jointcal.ConstrainedMagnitudeModel(associations.getCcdImageList(),
808  self.focalPlaneBBox,
809  visitOrder=self.config.photometryVisitOrder,
810  errorPedestal=self.config.photometryErrorPedestal)
811  # potentially nonlinear problem, so we may need a line search to converge.
812  doLineSearch = self.config.allowLineSearch
813  elif self.config.photometryModel == "simpleFlux":
814  model = lsst.jointcal.SimpleFluxModel(associations.getCcdImageList(),
815  errorPedestal=self.config.photometryErrorPedestal)
816  doLineSearch = False # purely linear in model parameters, so no line search needed
817  elif self.config.photometryModel == "simpleMagnitude":
818  model = lsst.jointcal.SimpleMagnitudeModel(associations.getCcdImageList(),
819  errorPedestal=self.config.photometryErrorPedestal)
820  doLineSearch = False # purely linear in model parameters, so no line search needed
821 
822  fit = lsst.jointcal.PhotometryFit(associations, model)
823  # TODO DM-12446: turn this into a "butler save" somehow.
824  # Save reference and measurement chi2 contributions for this data
825  if self.config.writeChi2FilesInitialFinal:
826  baseName = f"photometry_initial_chi2-{dataName}"
827  else:
828  baseName = None
829  self._logChi2AndValidate(associations, fit, model, "Initialized", writeChi2Name=baseName)
830 
831  def getChi2Name(whatToFit):
832  if self.config.writeChi2FilesOuterLoop:
833  return f"photometry_init-%s_chi2-{dataName}" % whatToFit
834  else:
835  return None
836 
837  # The constrained model needs the visit transform fit first; the chip
838  # transform is initialized from the singleFrame PhotoCalib, so it's close.
839  dumpMatrixFile = "photometry_preinit" if self.config.writeInitMatrix else ""
840  if self.config.photometryModel.startswith("constrained"):
841  # no line search: should be purely (or nearly) linear,
842  # and we want a large step size to initialize with.
843  fit.minimize("ModelVisit", dumpMatrixFile=dumpMatrixFile)
844  self._logChi2AndValidate(associations, fit, model, writeChi2Name=getChi2Name("ModelVisit"))
845  dumpMatrixFile = "" # so we don't redo the output on the next step
846 
847  fit.minimize("Model", doLineSearch=doLineSearch, dumpMatrixFile=dumpMatrixFile)
848  self._logChi2AndValidate(associations, fit, model, writeChi2Name=getChi2Name("Model"))
849 
850  fit.minimize("Fluxes") # no line search: always purely linear.
851  self._logChi2AndValidate(associations, fit, model, writeChi2Name=getChi2Name("Fluxes"))
852 
853  fit.minimize("Model Fluxes", doLineSearch=doLineSearch)
854  self._logChi2AndValidate(associations, fit, model, "Fit prepared",
855  writeChi2Name=getChi2Name("ModelFluxes"))
856 
857  model.freezeErrorTransform()
858  self.log.debug("Photometry error scales are frozen.")
859 
860  chi2 = self._iterate_fit(associations,
861  fit,
862  self.config.maxPhotometrySteps,
863  "photometry",
864  "Model Fluxes",
865  doRankUpdate=self.config.photometryDoRankUpdate,
866  doLineSearch=doLineSearch,
867  dataName=dataName)
868 
869  add_measurement(self.job, 'jointcal.photometry_final_chi2', chi2.chi2)
870  add_measurement(self.job, 'jointcal.photometry_final_ndof', chi2.ndof)
871  return Photometry(fit, model)
872 
873  def _fit_astrometry(self, associations, dataName=None):
874  """
875  Fit the astrometric data.
876 
877  Parameters
878  ----------
879  associations : `lsst.jointcal.Associations`
880  The star/reference star associations to fit.
881  dataName : `str`
882  Name of the data being processed (e.g. "1234_HSC-Y"), for
883  identifying debugging files.
884 
885  Returns
886  -------
887  fit_result : `namedtuple`
888  fit : `lsst.jointcal.AstrometryFit`
889  The astrometric fitter used to perform the fit.
890  model : `lsst.jointcal.AstrometryModel`
891  The astrometric model that was fit.
892  sky_to_tan_projection : `lsst.jointcal.ProjectionHandler`
893  The model for the sky to tangent plane projection that was used in the fit.
894  """
895 
896  self.log.info("=== Starting astrometric fitting...")
897 
898  associations.deprojectFittedStars()
899 
900  # NOTE: need to return sky_to_tan_projection so that it doesn't get garbage collected.
901  # TODO: could we package sky_to_tan_projection and model together so we don't have to manage
902  # them so carefully?
903  sky_to_tan_projection = lsst.jointcal.OneTPPerVisitHandler(associations.getCcdImageList())
904 
905  if self.config.astrometryModel == "constrained":
906  model = lsst.jointcal.ConstrainedAstrometryModel(associations.getCcdImageList(),
907  sky_to_tan_projection,
908  chipOrder=self.config.astrometryChipOrder,
909  visitOrder=self.config.astrometryVisitOrder)
910  elif self.config.astrometryModel == "simple":
911  model = lsst.jointcal.SimpleAstrometryModel(associations.getCcdImageList(),
912  sky_to_tan_projection,
913  self.config.useInputWcs,
914  nNotFit=0,
915  order=self.config.astrometrySimpleOrder)
916 
917  fit = lsst.jointcal.AstrometryFit(associations, model, self.config.positionErrorPedestal)
918  # TODO DM-12446: turn this into a "butler save" somehow.
919  # Save reference and measurement chi2 contributions for this data
920  if self.config.writeChi2FilesInitialFinal:
921  baseName = f"astrometry_initial_chi2-{dataName}"
922  else:
923  baseName = None
924  self._logChi2AndValidate(associations, fit, model, "Initial", writeChi2Name=baseName)
925 
926  def getChi2Name(whatToFit):
927  if self.config.writeChi2FilesOuterLoop:
928  return f"astrometry_init-%s_chi2-{dataName}" % whatToFit
929  else:
930  return None
931 
932  dumpMatrixFile = "astrometry_preinit" if self.config.writeInitMatrix else ""
933  # The constrained model needs the visit transform fit first; the chip
934  # transform is initialized from the detector's cameraGeom, so it's close.
935  if self.config.astrometryModel == "constrained":
936  fit.minimize("DistortionsVisit", dumpMatrixFile=dumpMatrixFile)
937  self._logChi2AndValidate(associations, fit, model, writeChi2Name=getChi2Name("DistortionsVisit"))
938  dumpMatrixFile = "" # so we don't redo the output on the next step
939 
940  fit.minimize("Distortions", dumpMatrixFile=dumpMatrixFile)
941  self._logChi2AndValidate(associations, fit, model, writeChi2Name=getChi2Name("Distortions"))
942 
943  fit.minimize("Positions")
944  self._logChi2AndValidate(associations, fit, model, writeChi2Name=getChi2Name("Positions"))
945 
946  fit.minimize("Distortions Positions")
947  self._logChi2AndValidate(associations, fit, model, "Fit prepared",
948  writeChi2Name=getChi2Name("DistortionsPositions"))
949 
950  chi2 = self._iterate_fit(associations,
951  fit,
952  self.config.maxAstrometrySteps,
953  "astrometry",
954  "Distortions Positions",
955  doRankUpdate=self.config.astrometryDoRankUpdate,
956  dataName=dataName)
957 
958  add_measurement(self.job, 'jointcal.astrometry_final_chi2', chi2.chi2)
959  add_measurement(self.job, 'jointcal.astrometry_final_ndof', chi2.ndof)
960 
961  return Astrometry(fit, model, sky_to_tan_projection)
962 
963  def _check_stars(self, associations):
964  """Count measured and reference stars per ccd and warn/log them."""
965  for ccdImage in associations.getCcdImageList():
966  nMeasuredStars, nRefStars = ccdImage.countStars()
967  self.log.debug("ccdImage %s has %s measured and %s reference stars",
968  ccdImage.getName(), nMeasuredStars, nRefStars)
969  if nMeasuredStars < self.config.minMeasuredStarsPerCcd:
970  self.log.warn("ccdImage %s has only %s measuredStars (desired %s)",
971  ccdImage.getName(), nMeasuredStars, self.config.minMeasuredStarsPerCcd)
972  if nRefStars < self.config.minRefStarsPerCcd:
973  self.log.warn("ccdImage %s has only %s RefStars (desired %s)",
974  ccdImage.getName(), nRefStars, self.config.minRefStarsPerCcd)
975 
976  def _iterate_fit(self, associations, fitter, max_steps, name, whatToFit,
977  dataName="",
978  doRankUpdate=True,
979  doLineSearch=False):
980  """Run fitter.minimize up to max_steps times, returning the final chi2.
981 
982  Parameters
983  ----------
984  associations : `lsst.jointcal.Associations`
985  The star/reference star associations to fit.
986  fitter : `lsst.jointcal.FitterBase`
987  The fitter to use for minimization.
988  max_steps : `int`
989  Maximum number of steps to run outlier rejection before declaring
990  convergence failure.
991  name : {'photometry' or 'astrometry'}
992  What type of data are we fitting (for logs and debugging files).
993  whatToFit : `str`
994  Passed to ``fitter.minimize()`` to define the parameters to fit.
995  dataName : `str`, optional
996  Descriptive name for this dataset (e.g. tract and filter),
997  for debugging.
998  doRankUpdate : `bool`, optional
999  Do an Eigen rank update during minimization, or recompute the full
1000  matrix and gradient?
1001  doLineSearch : `bool`, optional
1002  Do a line search for the optimum step during minimization?
1003 
1004  Returns
1005  -------
1006  chi2: `lsst.jointcal.Chi2Statistic`
1007  The final chi2 after the fit converges, or is forced to end.
1008 
1009  Raises
1010  ------
1011  FloatingPointError
1012  Raised if the fitter fails with a non-finite value.
1013  RuntimeError
1014  Raised if the fitter fails for some other reason;
1015  log messages will provide further details.
1016  """
1017  dumpMatrixFile = "%s_postinit" % name if self.config.writeInitMatrix else ""
1018  for i in range(max_steps):
1019  if self.config.writeChi2FilesOuterLoop:
1020  writeChi2Name = f"{name}_iterate_{i}_chi2-{dataName}"
1021  else:
1022  writeChi2Name = None
1023  result = fitter.minimize(whatToFit,
1024  self.config.outlierRejectSigma,
1025  doRankUpdate=doRankUpdate,
1026  doLineSearch=doLineSearch,
1027  dumpMatrixFile=dumpMatrixFile)
1028  dumpMatrixFile = "" # clear it so we don't write the matrix again.
1029  chi2 = self._logChi2AndValidate(associations, fitter, fitter.getModel(),
1030  writeChi2Name=writeChi2Name)
1031 
1032  if result == MinimizeResult.Converged:
1033  if doRankUpdate:
1034  self.log.debug("fit has converged - no more outliers - redo minimization "
1035  "one more time in case we have lost accuracy in rank update.")
1036  # Redo minimization one more time in case we have lost accuracy in rank update
1037  result = fitter.minimize(whatToFit, self.config.outlierRejectSigma)
1038  chi2 = self._logChi2AndValidate(associations, fitter, fitter.getModel(), "Fit completed")
1039 
1040  # log a message for a large final chi2, TODO: DM-15247 for something better
1041  if chi2.chi2/chi2.ndof >= 4.0:
1042  self.log.error("Potentially bad fit: High chi-squared/ndof.")
1043 
1044  break
1045  elif result == MinimizeResult.Chi2Increased:
1046  self.log.warn("still some outliers but chi2 increases - retry")
1047  elif result == MinimizeResult.NonFinite:
1048  filename = "{}_failure-nonfinite_chi2-{}.csv".format(name, dataName)
1049  # TODO DM-12446: turn this into a "butler save" somehow.
1050  fitter.saveChi2Contributions(filename)
1051  msg = "Nonfinite value in chi2 minimization, cannot complete fit. Dumped star tables to: {}"
1052  raise FloatingPointError(msg.format(filename))
1053  elif result == MinimizeResult.Failed:
1054  raise RuntimeError("Chi2 minimization failure, cannot complete fit.")
1055  else:
1056  raise RuntimeError("Unxepected return code from minimize().")
1057  else:
1058  self.log.error("%s failed to converge after %d steps"%(name, max_steps))
1059 
1060  return chi2
1061 
1062  def _write_astrometry_results(self, associations, model, visit_ccd_to_dataRef):
1063  """
1064  Write the fitted astrometric results to a new 'jointcal_wcs' dataRef.
1065 
1066  Parameters
1067  ----------
1068  associations : `lsst.jointcal.Associations`
1069  The star/reference star associations to fit.
1070  model : `lsst.jointcal.AstrometryModel`
1071  The astrometric model that was fit.
1072  visit_ccd_to_dataRef : `dict` of Key: `lsst.daf.persistence.ButlerDataRef`
1073  Dict of ccdImage identifiers to dataRefs that were fit.
1074  """
1075 
1076  ccdImageList = associations.getCcdImageList()
1077  for ccdImage in ccdImageList:
1078  # TODO: there must be a better way to identify this ccdImage than a visit,ccd pair?
1079  ccd = ccdImage.ccdId
1080  visit = ccdImage.visit
1081  dataRef = visit_ccd_to_dataRef[(visit, ccd)]
1082  self.log.info("Updating WCS for visit: %d, ccd: %d", visit, ccd)
1083  skyWcs = model.makeSkyWcs(ccdImage)
1084  try:
1085  dataRef.put(skyWcs, 'jointcal_wcs')
1086  except pexExceptions.Exception as e:
1087  self.log.fatal('Failed to write updated Wcs: %s', str(e))
1088  raise e
1089 
1090  def _write_photometry_results(self, associations, model, visit_ccd_to_dataRef):
1091  """
1092  Write the fitted photometric results to a new 'jointcal_photoCalib' dataRef.
1093 
1094  Parameters
1095  ----------
1096  associations : `lsst.jointcal.Associations`
1097  The star/reference star associations to fit.
1098  model : `lsst.jointcal.PhotometryModel`
1099  The photoometric model that was fit.
1100  visit_ccd_to_dataRef : `dict` of Key: `lsst.daf.persistence.ButlerDataRef`
1101  Dict of ccdImage identifiers to dataRefs that were fit.
1102  """
1103 
1104  ccdImageList = associations.getCcdImageList()
1105  for ccdImage in ccdImageList:
1106  # TODO: there must be a better way to identify this ccdImage than a visit,ccd pair?
1107  ccd = ccdImage.ccdId
1108  visit = ccdImage.visit
1109  dataRef = visit_ccd_to_dataRef[(visit, ccd)]
1110  self.log.info("Updating PhotoCalib for visit: %d, ccd: %d", visit, ccd)
1111  photoCalib = model.toPhotoCalib(ccdImage)
1112  try:
1113  dataRef.put(photoCalib, 'jointcal_photoCalib')
1114  except pexExceptions.Exception as e:
1115  self.log.fatal('Failed to write updated PhotoCalib: %s', str(e))
1116  raise e
def runDataRef(self, dataRefs, profile_jointcal=False)
Definition: jointcal.py:471
def _load_reference_catalog(self, refObjLoader, referenceSelector, center, radius, filterName, applyColorterms=False)
Definition: jointcal.py:658
def _build_ccdImage(self, dataRef, associations, jointcalControl)
Definition: jointcal.py:408
def _fit_photometry(self, associations, dataName=None)
Definition: jointcal.py:776
def getTargetList(parsedCmd, kwargs)
Definition: jointcal.py:70
def _fit_astrometry(self, associations, dataName=None)
Definition: jointcal.py:873
def _check_star_lists(self, associations, name)
Definition: jointcal.py:724
The class that implements the relations between MeasuredStar and FittedStar.
Definition: Associations.h:54
A projection handler in which all CCDs from the same visit have the same tangent point.
def _logChi2AndValidate(self, associations, fit, model, chi2Label="Model", writeChi2Name=None)
Definition: jointcal.py:734
A model where there is one independent transform per CcdImage.
def _iterate_fit(self, associations, fitter, max_steps, name, whatToFit, dataName="", doRankUpdate=True, doLineSearch=False)
Definition: jointcal.py:979
def _write_photometry_results(self, associations, model, visit_ccd_to_dataRef)
Definition: jointcal.py:1090
def _check_stars(self, associations)
Definition: jointcal.py:963
Class that handles the photometric least squares problem.
Definition: PhotometryFit.h:45
Class that handles the astrometric least squares problem.
Definition: AstrometryFit.h:78
def add_measurement(job, name, value)
Definition: jointcal.py:52
def _do_load_refcat_and_fit(self, associations, defaultFilter, center, radius, name="", refObjLoader=None, referenceSelector=None, filters=[], fit_function=None, tract=None, profile_jointcal=False, match_cut=3.0, reject_bad_fluxes=False)
Definition: jointcal.py:570
A multi-component model, fitting mappings for sensors and visits simultaneously.
def _write_astrometry_results(self, associations, model, visit_ccd_to_dataRef)
Definition: jointcal.py:1062
def __init__(self, butler=None, profile_jointcal=False, kwargs)
Definition: jointcal.py:362