lsst.jointcal  17.0.1-15-g1718489+9
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.utils
27 import lsst.pex.config as pexConfig
28 import lsst.pipe.base as pipeBase
29 from lsst.afw.image import fluxErrFromABMagErr
30 import lsst.afw.geom as afwGeom
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  # Use external reference catalogs handled by LSST stack mechanism
522  # Get the bounding box overlapping all associated images
523  # ==> This is probably a bad idea to do it this way <== To be improved
524  bbox = associations.getRaDecBBox()
525  bboxCenter = bbox.getCenter()
526  center = afwGeom.SpherePoint(bboxCenter[0], bboxCenter[1], afwGeom.degrees)
527  bboxMax = bbox.getMax()
528  corner = afwGeom.SpherePoint(bboxMax[0], bboxMax[1], afwGeom.degrees)
529  radius = center.separation(corner).asRadians()
530 
531  # Get astrometry_net_data path
532  anDir = lsst.utils.getPackageDir('astrometry_net_data')
533  if anDir is None:
534  raise RuntimeError("astrometry_net_data is not setup")
535 
536  # Determine a default filter associated with the catalog. See DM-9093
537  defaultFilter = filters.most_common(1)[0][0]
538  self.log.debug("Using %s band for reference flux", defaultFilter)
539 
540  # TODO: need a better way to get the tract.
541  tract = dataRefs[0].dataId['tract']
542 
543  if self.config.doAstrometry:
544  astrometry = self._do_load_refcat_and_fit(associations, defaultFilter, center, radius,
545  name="astrometry",
546  refObjLoader=self.astrometryRefObjLoader,
547  referenceSelector=self.astrometryReferenceSelector,
548  fit_function=self._fit_astrometry,
549  profile_jointcal=profile_jointcal,
550  tract=tract)
551  self._write_astrometry_results(associations, astrometry.model, visit_ccd_to_dataRef)
552  else:
553  astrometry = Astrometry(None, None, None)
554 
555  if self.config.doPhotometry:
556  photometry = self._do_load_refcat_and_fit(associations, defaultFilter, center, radius,
557  name="photometry",
558  refObjLoader=self.photometryRefObjLoader,
559  referenceSelector=self.photometryReferenceSelector,
560  fit_function=self._fit_photometry,
561  profile_jointcal=profile_jointcal,
562  tract=tract,
563  filters=filters,
564  reject_bad_fluxes=True)
565  self._write_photometry_results(associations, photometry.model, visit_ccd_to_dataRef)
566  else:
567  photometry = Photometry(None, None)
568 
569  return pipeBase.Struct(dataRefs=dataRefs,
570  oldWcsList=oldWcsList,
571  job=self.job,
572  astrometryRefObjLoader=self.astrometryRefObjLoader,
573  photometryRefObjLoader=self.photometryRefObjLoader,
574  defaultFilter=defaultFilter,
575  exitStatus=exitStatus)
576 
577  def _do_load_refcat_and_fit(self, associations, defaultFilter, center, radius,
578  name="", refObjLoader=None, referenceSelector=None,
579  filters=[], fit_function=None,
580  tract=None, profile_jointcal=False, match_cut=3.0,
581  reject_bad_fluxes=False):
582  """Load reference catalog, perform the fit, and return the result.
583 
584  Parameters
585  ----------
586  associations : `lsst.jointcal.Associations`
587  The star/reference star associations to fit.
588  defaultFilter : `str`
589  filter to load from reference catalog.
590  center : `lsst.afw.geom.SpherePoint`
591  ICRS center of field to load from reference catalog.
592  radius : `lsst.afw.geom.Angle`
593  On-sky radius to load from reference catalog.
594  name : `str`
595  Name of thing being fit: "Astrometry" or "Photometry".
596  refObjLoader : `lsst.meas.algorithms.LoadReferenceObjectsTask`
597  Reference object loader to load from for fit.
598  filters : `list` of `str`, optional
599  List of filters to load from the reference catalog.
600  fit_function : callable
601  Function to call to perform fit (takes associations object).
602  tract : `str`
603  Name of tract currently being fit.
604  profile_jointcal : `bool`, optional
605  Separately profile the fitting step.
606  match_cut : `float`, optional
607  Radius in arcseconds to find cross-catalog matches to during
608  associations.associateCatalogs.
609  reject_bad_fluxes : `bool`, optional
610  Reject refCat sources with NaN/inf flux or NaN/0 fluxErr.
611 
612  Returns
613  -------
614  result : `Photometry` or `Astrometry`
615  Result of `fit_function()`
616  """
617  self.log.info("====== Now processing %s...", name)
618  # TODO: this should not print "trying to invert a singular transformation:"
619  # if it does that, something's not right about the WCS...
620  associations.associateCatalogs(match_cut)
621  add_measurement(self.job, 'jointcal.associated_%s_fittedStars' % name,
622  associations.fittedStarListSize())
623 
624  applyColorterms = False if name == "Astrometry" else self.config.applyColorTerms
625  if name == "Astrometry":
626  referenceSelector = self.config.astrometryReferenceSelector
627  elif name == "Photometry":
628  referenceSelector = self.config.photometryReferenceSelector
629  refCat, fluxField = self._load_reference_catalog(refObjLoader, referenceSelector,
630  center, radius, defaultFilter,
631  applyColorterms=applyColorterms)
632 
633  if self.config.astrometryReferenceErr is None:
634  refCoordErr = float('nan')
635  else:
636  refCoordErr = self.config.astrometryReferenceErr
637 
638  associations.collectRefStars(refCat,
639  self.config.matchCut*afwGeom.arcseconds,
640  fluxField,
641  refCoordinateErr=refCoordErr,
642  rejectBadFluxes=reject_bad_fluxes)
643  add_measurement(self.job, 'jointcal.collected_%s_refStars' % name,
644  associations.refStarListSize())
645 
646  associations.prepareFittedStars(self.config.minMeasurements)
647 
648  self._check_star_lists(associations, name)
649  add_measurement(self.job, 'jointcal.selected_%s_refStars' % name,
650  associations.nFittedStarsWithAssociatedRefStar())
651  add_measurement(self.job, 'jointcal.selected_%s_fittedStars' % name,
652  associations.fittedStarListSize())
653  add_measurement(self.job, 'jointcal.selected_%s_ccdImages' % name,
654  associations.nCcdImagesValidForFit())
655 
656  load_cat_prof_file = 'jointcal_fit_%s.prof'%name if profile_jointcal else ''
657  dataName = "{}_{}".format(tract, defaultFilter)
658  with pipeBase.cmdLineTask.profile(load_cat_prof_file):
659  result = fit_function(associations, dataName)
660  # TODO DM-12446: turn this into a "butler save" somehow.
661  # Save reference and measurement chi2 contributions for this data
662  if self.config.writeChi2FilesInitialFinal:
663  baseName = f"{name}_final_chi2-{dataName}"
664  result.fit.saveChi2Contributions(baseName+"{type}")
665 
666  return result
667 
668  def _load_reference_catalog(self, refObjLoader, referenceSelector, center, radius, filterName,
669  applyColorterms=False):
670  """Load the necessary reference catalog sources, convert fluxes to
671  correct units, and apply color term corrections if requested.
672 
673  Parameters
674  ----------
675  refObjLoader : `lsst.meas.algorithms.LoadReferenceObjectsTask`
676  The reference catalog loader to use to get the data.
677  referenceSelector : `lsst.meas.algorithms.ReferenceSourceSelectorTask`
678  Source selector to apply to loaded reference catalog.
679  center : `lsst.geom.SpherePoint`
680  The center around which to load sources.
681  radius : `lsst.geom.Angle`
682  The radius around ``center`` to load sources in.
683  filterName : `str`
684  The name of the camera filter to load fluxes for.
685  applyColorterms : `bool`
686  Apply colorterm corrections to the refcat for ``filterName``?
687 
688  Returns
689  -------
690  refCat : `lsst.afw.table.SimpleCatalog`
691  The loaded reference catalog.
692  fluxField : `str`
693  The name of the reference catalog flux field appropriate for ``filterName``.
694  """
695  skyCircle = refObjLoader.loadSkyCircle(center,
696  afwGeom.Angle(radius, afwGeom.radians),
697  filterName)
698 
699  selected = referenceSelector.run(skyCircle.refCat)
700  # Need memory contiguity to get reference filters as a vector.
701  if not selected.sourceCat.isContiguous():
702  refCat = selected.sourceCat.copy(deep=True)
703  else:
704  refCat = selected.sourceCat
705 
706  if self.config.astrometryReferenceErr is None and 'coord_ra_err' not in refCat.schema:
707  msg = ("Reference catalog does not contain coordinate errors, "
708  "and config.astrometryReferenceErr not supplied.")
709  raise pexConfig.FieldValidationError(JointcalConfig.astrometryReferenceErr,
710  self.config,
711  msg)
712 
713  if self.config.astrometryReferenceErr is not None and 'coord_ra_err' in refCat.schema:
714  self.log.warn("Overriding reference catalog coordinate errors with %f/coordinate [mas]",
715  self.config.astrometryReferenceErr)
716 
717  if applyColorterms:
718  try:
719  refCatName = refObjLoader.ref_dataset_name
720  except AttributeError:
721  # NOTE: we need this try:except: block in place until we've completely removed a.net support.
722  raise RuntimeError("Cannot perform colorterm corrections with a.net refcats.")
723  self.log.info("Applying color terms for filterName=%r reference catalog=%s",
724  filterName, refCatName)
725  colorterm = self.config.colorterms.getColorterm(
726  filterName=filterName, photoCatName=refCatName, doRaise=True)
727 
728  refMag, refMagErr = colorterm.getCorrectedMagnitudes(refCat, filterName)
729  refCat[skyCircle.fluxField] = u.Magnitude(refMag, u.ABmag).to_value(u.nJy)
730  # TODO: I didn't want to use this, but I'll deal with it in DM-16903
731  refCat[skyCircle.fluxField+'Err'] = fluxErrFromABMagErr(refMagErr, refMag) * 1e9
732 
733  return refCat, skyCircle.fluxField
734 
735  def _check_star_lists(self, associations, name):
736  # TODO: these should be len(blah), but we need this properly wrapped first.
737  if associations.nCcdImagesValidForFit() == 0:
738  raise RuntimeError('No images in the ccdImageList!')
739  if associations.fittedStarListSize() == 0:
740  raise RuntimeError('No stars in the {} fittedStarList!'.format(name))
741  if associations.refStarListSize() == 0:
742  raise RuntimeError('No stars in the {} reference star list!'.format(name))
743 
744  def _logChi2AndValidate(self, associations, fit, model, chi2Label="Model",
745  writeChi2Name=None):
746  """Compute chi2, log it, validate the model, and return chi2.
747 
748  Parameters
749  ----------
750  associations : `lsst.jointcal.Associations`
751  The star/reference star associations to fit.
752  fit : `lsst.jointcal.FitterBase`
753  The fitter to use for minimization.
754  model : `lsst.jointcal.Model`
755  The model being fit.
756  chi2Label : str, optional
757  Label to describe the chi2 (e.g. "Initialized", "Final").
758  writeChi2Name : `str`, optional
759  Filename prefix to write the chi2 contributions to.
760  Do not supply an extension: an appropriate one will be added.
761 
762  Returns
763  -------
764  chi2: `lsst.jointcal.Chi2Accumulator`
765  The chi2 object for the current fitter and model.
766 
767  Raises
768  ------
769  FloatingPointError
770  Raised if chi2 is infinite or NaN.
771  ValueError
772  Raised if the model is not valid.
773  """
774  if writeChi2Name is not None:
775  fit.saveChi2Contributions(writeChi2Name+"{type}")
776  self.log.info("Wrote chi2 contributions files: %s", writeChi2Name)
777 
778  chi2 = fit.computeChi2()
779  self.log.info("%s %s", chi2Label, chi2)
780  self._check_stars(associations)
781  if not np.isfinite(chi2.chi2):
782  raise FloatingPointError(f'{chi2Label} chi2 is invalid: {chi2}')
783  if not model.validate(associations.getCcdImageList(), chi2.ndof):
784  raise ValueError("Model is not valid: check log messages for warnings.")
785  return chi2
786 
787  def _fit_photometry(self, associations, dataName=None):
788  """
789  Fit the photometric data.
790 
791  Parameters
792  ----------
793  associations : `lsst.jointcal.Associations`
794  The star/reference star associations to fit.
795  dataName : `str`
796  Name of the data being processed (e.g. "1234_HSC-Y"), for
797  identifying debugging files.
798 
799  Returns
800  -------
801  fit_result : `namedtuple`
802  fit : `lsst.jointcal.PhotometryFit`
803  The photometric fitter used to perform the fit.
804  model : `lsst.jointcal.PhotometryModel`
805  The photometric model that was fit.
806  """
807  self.log.info("=== Starting photometric fitting...")
808 
809  # TODO: should use pex.config.RegistryField here (see DM-9195)
810  if self.config.photometryModel == "constrainedFlux":
811  model = lsst.jointcal.ConstrainedFluxModel(associations.getCcdImageList(),
812  self.focalPlaneBBox,
813  visitOrder=self.config.photometryVisitOrder,
814  errorPedestal=self.config.photometryErrorPedestal)
815  # potentially nonlinear problem, so we may need a line search to converge.
816  doLineSearch = self.config.allowLineSearch
817  elif self.config.photometryModel == "constrainedMagnitude":
818  model = lsst.jointcal.ConstrainedMagnitudeModel(associations.getCcdImageList(),
819  self.focalPlaneBBox,
820  visitOrder=self.config.photometryVisitOrder,
821  errorPedestal=self.config.photometryErrorPedestal)
822  # potentially nonlinear problem, so we may need a line search to converge.
823  doLineSearch = self.config.allowLineSearch
824  elif self.config.photometryModel == "simpleFlux":
825  model = lsst.jointcal.SimpleFluxModel(associations.getCcdImageList(),
826  errorPedestal=self.config.photometryErrorPedestal)
827  doLineSearch = False # purely linear in model parameters, so no line search needed
828  elif self.config.photometryModel == "simpleMagnitude":
829  model = lsst.jointcal.SimpleMagnitudeModel(associations.getCcdImageList(),
830  errorPedestal=self.config.photometryErrorPedestal)
831  doLineSearch = False # purely linear in model parameters, so no line search needed
832 
833  fit = lsst.jointcal.PhotometryFit(associations, model)
834  # TODO DM-12446: turn this into a "butler save" somehow.
835  # Save reference and measurement chi2 contributions for this data
836  if self.config.writeChi2FilesInitialFinal:
837  baseName = f"photometry_initial_chi2-{dataName}"
838  else:
839  baseName = None
840  self._logChi2AndValidate(associations, fit, model, "Initialized", writeChi2Name=baseName)
841 
842  def getChi2Name(whatToFit):
843  if self.config.writeChi2FilesOuterLoop:
844  return f"photometry_init-%s_chi2-{dataName}" % whatToFit
845  else:
846  return None
847 
848  # The constrained model needs the visit transform fit first; the chip
849  # transform is initialized from the singleFrame PhotoCalib, so it's close.
850  dumpMatrixFile = "photometry_preinit" if self.config.writeInitMatrix else ""
851  if self.config.photometryModel.startswith("constrained"):
852  # no line search: should be purely (or nearly) linear,
853  # and we want a large step size to initialize with.
854  fit.minimize("ModelVisit", dumpMatrixFile=dumpMatrixFile)
855  self._logChi2AndValidate(associations, fit, model, writeChi2Name=getChi2Name("ModelVisit"))
856  dumpMatrixFile = "" # so we don't redo the output on the next step
857 
858  fit.minimize("Model", doLineSearch=doLineSearch, dumpMatrixFile=dumpMatrixFile)
859  self._logChi2AndValidate(associations, fit, model, writeChi2Name=getChi2Name("Model"))
860 
861  fit.minimize("Fluxes") # no line search: always purely linear.
862  self._logChi2AndValidate(associations, fit, model, writeChi2Name=getChi2Name("Fluxes"))
863 
864  fit.minimize("Model Fluxes", doLineSearch=doLineSearch)
865  self._logChi2AndValidate(associations, fit, model, "Fit prepared",
866  writeChi2Name=getChi2Name("ModelFluxes"))
867 
868  model.freezeErrorTransform()
869  self.log.debug("Photometry error scales are frozen.")
870 
871  chi2 = self._iterate_fit(associations,
872  fit,
873  self.config.maxPhotometrySteps,
874  "photometry",
875  "Model Fluxes",
876  doRankUpdate=self.config.photometryDoRankUpdate,
877  doLineSearch=doLineSearch,
878  dataName=dataName)
879 
880  add_measurement(self.job, 'jointcal.photometry_final_chi2', chi2.chi2)
881  add_measurement(self.job, 'jointcal.photometry_final_ndof', chi2.ndof)
882  return Photometry(fit, model)
883 
884  def _fit_astrometry(self, associations, dataName=None):
885  """
886  Fit the astrometric data.
887 
888  Parameters
889  ----------
890  associations : `lsst.jointcal.Associations`
891  The star/reference star associations to fit.
892  dataName : `str`
893  Name of the data being processed (e.g. "1234_HSC-Y"), for
894  identifying debugging files.
895 
896  Returns
897  -------
898  fit_result : `namedtuple`
899  fit : `lsst.jointcal.AstrometryFit`
900  The astrometric fitter used to perform the fit.
901  model : `lsst.jointcal.AstrometryModel`
902  The astrometric model that was fit.
903  sky_to_tan_projection : `lsst.jointcal.ProjectionHandler`
904  The model for the sky to tangent plane projection that was used in the fit.
905  """
906 
907  self.log.info("=== Starting astrometric fitting...")
908 
909  associations.deprojectFittedStars()
910 
911  # NOTE: need to return sky_to_tan_projection so that it doesn't get garbage collected.
912  # TODO: could we package sky_to_tan_projection and model together so we don't have to manage
913  # them so carefully?
914  sky_to_tan_projection = lsst.jointcal.OneTPPerVisitHandler(associations.getCcdImageList())
915 
916  if self.config.astrometryModel == "constrained":
917  model = lsst.jointcal.ConstrainedAstrometryModel(associations.getCcdImageList(),
918  sky_to_tan_projection,
919  chipOrder=self.config.astrometryChipOrder,
920  visitOrder=self.config.astrometryVisitOrder)
921  elif self.config.astrometryModel == "simple":
922  model = lsst.jointcal.SimpleAstrometryModel(associations.getCcdImageList(),
923  sky_to_tan_projection,
924  self.config.useInputWcs,
925  nNotFit=0,
926  order=self.config.astrometrySimpleOrder)
927 
928  fit = lsst.jointcal.AstrometryFit(associations, model, self.config.positionErrorPedestal)
929  # TODO DM-12446: turn this into a "butler save" somehow.
930  # Save reference and measurement chi2 contributions for this data
931  if self.config.writeChi2FilesInitialFinal:
932  baseName = f"astrometry_initial_chi2-{dataName}"
933  else:
934  baseName = None
935  self._logChi2AndValidate(associations, fit, model, "Initial", writeChi2Name=baseName)
936 
937  def getChi2Name(whatToFit):
938  if self.config.writeChi2FilesOuterLoop:
939  return f"astrometry_init-%s_chi2-{dataName}" % whatToFit
940  else:
941  return None
942 
943  dumpMatrixFile = "astrometry_preinit" if self.config.writeInitMatrix else ""
944  # The constrained model needs the visit transform fit first; the chip
945  # transform is initialized from the detector's cameraGeom, so it's close.
946  if self.config.astrometryModel == "constrained":
947  fit.minimize("DistortionsVisit", dumpMatrixFile=dumpMatrixFile)
948  self._logChi2AndValidate(associations, fit, model, writeChi2Name=getChi2Name("DistortionsVisit"))
949  dumpMatrixFile = "" # so we don't redo the output on the next step
950 
951  fit.minimize("Distortions", dumpMatrixFile=dumpMatrixFile)
952  self._logChi2AndValidate(associations, fit, model, writeChi2Name=getChi2Name("Distortions"))
953 
954  fit.minimize("Positions")
955  self._logChi2AndValidate(associations, fit, model, writeChi2Name=getChi2Name("Positions"))
956 
957  fit.minimize("Distortions Positions")
958  self._logChi2AndValidate(associations, fit, model, "Fit prepared",
959  writeChi2Name=getChi2Name("DistortionsPositions"))
960 
961  chi2 = self._iterate_fit(associations,
962  fit,
963  self.config.maxAstrometrySteps,
964  "astrometry",
965  "Distortions Positions",
966  doRankUpdate=self.config.astrometryDoRankUpdate,
967  dataName=dataName)
968 
969  add_measurement(self.job, 'jointcal.astrometry_final_chi2', chi2.chi2)
970  add_measurement(self.job, 'jointcal.astrometry_final_ndof', chi2.ndof)
971 
972  return Astrometry(fit, model, sky_to_tan_projection)
973 
974  def _check_stars(self, associations):
975  """Count measured and reference stars per ccd and warn/log them."""
976  for ccdImage in associations.getCcdImageList():
977  nMeasuredStars, nRefStars = ccdImage.countStars()
978  self.log.debug("ccdImage %s has %s measured and %s reference stars",
979  ccdImage.getName(), nMeasuredStars, nRefStars)
980  if nMeasuredStars < self.config.minMeasuredStarsPerCcd:
981  self.log.warn("ccdImage %s has only %s measuredStars (desired %s)",
982  ccdImage.getName(), nMeasuredStars, self.config.minMeasuredStarsPerCcd)
983  if nRefStars < self.config.minRefStarsPerCcd:
984  self.log.warn("ccdImage %s has only %s RefStars (desired %s)",
985  ccdImage.getName(), nRefStars, self.config.minRefStarsPerCcd)
986 
987  def _iterate_fit(self, associations, fitter, max_steps, name, whatToFit,
988  dataName="",
989  doRankUpdate=True,
990  doLineSearch=False):
991  """Run fitter.minimize up to max_steps times, returning the final chi2.
992 
993  Parameters
994  ----------
995  associations : `lsst.jointcal.Associations`
996  The star/reference star associations to fit.
997  fitter : `lsst.jointcal.FitterBase`
998  The fitter to use for minimization.
999  max_steps : `int`
1000  Maximum number of steps to run outlier rejection before declaring
1001  convergence failure.
1002  name : {'photometry' or 'astrometry'}
1003  What type of data are we fitting (for logs and debugging files).
1004  whatToFit : `str`
1005  Passed to ``fitter.minimize()`` to define the parameters to fit.
1006  dataName : `str`, optional
1007  Descriptive name for this dataset (e.g. tract and filter),
1008  for debugging.
1009  doRankUpdate : `bool`, optional
1010  Do an Eigen rank update during minimization, or recompute the full
1011  matrix and gradient?
1012  doLineSearch : `bool`, optional
1013  Do a line search for the optimum step during minimization?
1014 
1015  Returns
1016  -------
1017  chi2: `lsst.jointcal.Chi2Statistic`
1018  The final chi2 after the fit converges, or is forced to end.
1019 
1020  Raises
1021  ------
1022  FloatingPointError
1023  Raised if the fitter fails with a non-finite value.
1024  RuntimeError
1025  Raised if the fitter fails for some other reason;
1026  log messages will provide further details.
1027  """
1028  dumpMatrixFile = "%s_postinit" % name if self.config.writeInitMatrix else ""
1029  for i in range(max_steps):
1030  if self.config.writeChi2FilesOuterLoop:
1031  writeChi2Name = f"{name}_iterate_{i}_chi2-{dataName}"
1032  else:
1033  writeChi2Name = None
1034  result = fitter.minimize(whatToFit,
1035  self.config.outlierRejectSigma,
1036  doRankUpdate=doRankUpdate,
1037  doLineSearch=doLineSearch,
1038  dumpMatrixFile=dumpMatrixFile)
1039  dumpMatrixFile = "" # clear it so we don't write the matrix again.
1040  chi2 = self._logChi2AndValidate(associations, fitter, fitter.getModel(),
1041  writeChi2Name=writeChi2Name)
1042 
1043  if result == MinimizeResult.Converged:
1044  if doRankUpdate:
1045  self.log.debug("fit has converged - no more outliers - redo minimization "
1046  "one more time in case we have lost accuracy in rank update.")
1047  # Redo minimization one more time in case we have lost accuracy in rank update
1048  result = fitter.minimize(whatToFit, self.config.outlierRejectSigma)
1049  chi2 = self._logChi2AndValidate(associations, fitter, fitter.getModel(), "Fit completed")
1050 
1051  # log a message for a large final chi2, TODO: DM-15247 for something better
1052  if chi2.chi2/chi2.ndof >= 4.0:
1053  self.log.error("Potentially bad fit: High chi-squared/ndof.")
1054 
1055  break
1056  elif result == MinimizeResult.Chi2Increased:
1057  self.log.warn("still some outliers but chi2 increases - retry")
1058  elif result == MinimizeResult.NonFinite:
1059  filename = "{}_failure-nonfinite_chi2-{}.csv".format(name, dataName)
1060  # TODO DM-12446: turn this into a "butler save" somehow.
1061  fitter.saveChi2Contributions(filename)
1062  msg = "Nonfinite value in chi2 minimization, cannot complete fit. Dumped star tables to: {}"
1063  raise FloatingPointError(msg.format(filename))
1064  elif result == MinimizeResult.Failed:
1065  raise RuntimeError("Chi2 minimization failure, cannot complete fit.")
1066  else:
1067  raise RuntimeError("Unxepected return code from minimize().")
1068  else:
1069  self.log.error("%s failed to converge after %d steps"%(name, max_steps))
1070 
1071  return chi2
1072 
1073  def _write_astrometry_results(self, associations, model, visit_ccd_to_dataRef):
1074  """
1075  Write the fitted astrometric results to a new 'jointcal_wcs' dataRef.
1076 
1077  Parameters
1078  ----------
1079  associations : `lsst.jointcal.Associations`
1080  The star/reference star associations to fit.
1081  model : `lsst.jointcal.AstrometryModel`
1082  The astrometric model that was fit.
1083  visit_ccd_to_dataRef : `dict` of Key: `lsst.daf.persistence.ButlerDataRef`
1084  Dict of ccdImage identifiers to dataRefs that were fit.
1085  """
1086 
1087  ccdImageList = associations.getCcdImageList()
1088  for ccdImage in ccdImageList:
1089  # TODO: there must be a better way to identify this ccdImage than a visit,ccd pair?
1090  ccd = ccdImage.ccdId
1091  visit = ccdImage.visit
1092  dataRef = visit_ccd_to_dataRef[(visit, ccd)]
1093  self.log.info("Updating WCS for visit: %d, ccd: %d", visit, ccd)
1094  skyWcs = model.makeSkyWcs(ccdImage)
1095  try:
1096  dataRef.put(skyWcs, 'jointcal_wcs')
1097  except pexExceptions.Exception as e:
1098  self.log.fatal('Failed to write updated Wcs: %s', str(e))
1099  raise e
1100 
1101  def _write_photometry_results(self, associations, model, visit_ccd_to_dataRef):
1102  """
1103  Write the fitted photometric results to a new 'jointcal_photoCalib' dataRef.
1104 
1105  Parameters
1106  ----------
1107  associations : `lsst.jointcal.Associations`
1108  The star/reference star associations to fit.
1109  model : `lsst.jointcal.PhotometryModel`
1110  The photoometric model that was fit.
1111  visit_ccd_to_dataRef : `dict` of Key: `lsst.daf.persistence.ButlerDataRef`
1112  Dict of ccdImage identifiers to dataRefs that were fit.
1113  """
1114 
1115  ccdImageList = associations.getCcdImageList()
1116  for ccdImage in ccdImageList:
1117  # TODO: there must be a better way to identify this ccdImage than a visit,ccd pair?
1118  ccd = ccdImage.ccdId
1119  visit = ccdImage.visit
1120  dataRef = visit_ccd_to_dataRef[(visit, ccd)]
1121  self.log.info("Updating PhotoCalib for visit: %d, ccd: %d", visit, ccd)
1122  photoCalib = model.toPhotoCalib(ccdImage)
1123  try:
1124  dataRef.put(photoCalib, 'jointcal_photoCalib')
1125  except pexExceptions.Exception as e:
1126  self.log.fatal('Failed to write updated PhotoCalib: %s', str(e))
1127  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:669
def _build_ccdImage(self, dataRef, associations, jointcalControl)
Definition: jointcal.py:408
def _fit_photometry(self, associations, dataName=None)
Definition: jointcal.py:787
def getTargetList(parsedCmd, kwargs)
Definition: jointcal.py:70
def _fit_astrometry(self, associations, dataName=None)
Definition: jointcal.py:884
def _check_star_lists(self, associations, name)
Definition: jointcal.py:735
The class that implements the relations between MeasuredStar and FittedStar.
Definition: Associations.h:53
A projection handler in which all CCDs from the same visit have the same tangent point.
std::string getPackageDir(std::string const &packageName)
def _logChi2AndValidate(self, associations, fit, model, chi2Label="Model", writeChi2Name=None)
Definition: jointcal.py:745
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:990
def _write_photometry_results(self, associations, model, visit_ccd_to_dataRef)
Definition: jointcal.py:1101
def _check_stars(self, associations)
Definition: jointcal.py:974
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:581
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:1073
def __init__(self, butler=None, profile_jointcal=False, kwargs)
Definition: jointcal.py:362