lsst.meas.astrom  16.0-21-gdae7b8c+3
matchPessimisticB.py
Go to the documentation of this file.
1 
2 import numpy as np
3 from scipy.spatial import cKDTree
4 from scipy.stats import sigmaclip
5 
6 import lsst.pex.config as pexConfig
7 import lsst.pipe.base as pipeBase
8 import lsst.afw.geom as afwgeom
9 import lsst.afw.table as afwTable
10 from lsst.meas.algorithms.sourceSelector import sourceSelectorRegistry
11 
12 from .matchOptimisticB import MatchTolerance
13 
14 from .pessimistic_pattern_matcher_b_3D import PessimisticPatternMatcherB
15 
16 __all__ = ["MatchPessimisticBTask", "MatchPessimisticBConfig",
17  "MatchTolerancePessimistic"]
18 
19 
21  """Stores match tolerances for use in AstrometryTask and later
22  iterations of the matcher.
23 
24  MatchPessimisticBTask relies on several state variables to be
25  preserved over different iterations in the
26  AstrometryTask.matchAndFitWcs loop of AstrometryTask.
27 
28  Parameters
29  ----------
30  maxMatchDist : `lsst.geom.Angle`
31  Maximum distance to consider a match from the previous match/fit
32  iteration.
33  autoMaxMatchDist : `lsst.geom.Angle`
34  Automated estimation of the maxMatchDist from the sky statistics of the
35  source and reference catalogs.
36  maxShift : `lsst.geom.Angle`
37  Maximum shift found in the previous match/fit cycle.
38  lastMatchedPattern : `int`
39  Index of the last source pattern that was matched into the reference
40  data.
41  failedPatternList : `list` of `int`
42  Previous matches were found to be false positives.
43  PPMbObj : `lsst.meas.astrom.PessimisticPatternMatcherB`
44  Initialized Pessimistic pattern matcher object. Storing this prevents
45  the need for recalculation of the searchable distances in the PPMB.
46  """
47 
48  def __init__(self, maxMatchDist=None, autoMaxMatchDist=None,
49  maxShift=None, lastMatchedPattern=None,
50  failedPatternList=None, PPMbObj=None):
51  self.maxMatchDist = maxMatchDist
52  self.autoMaxMatchDist = autoMaxMatchDist
53  self.maxShift = maxShift
54  self.lastMatchedPattern = lastMatchedPattern
55  self.PPMbObj = PPMbObj
56  if failedPatternList is None:
58  else:
59  self.failedPatternList = failedPatternList
60 
61 
62 class MatchPessimisticBConfig(pexConfig.Config):
63  """Configuration for MatchPessimisticBTask
64  """
65  numBrightStars = pexConfig.RangeField(
66  doc="Number of bright stars to use. Sets the max number of patterns "
67  "that can be tested.",
68  dtype=int,
69  default=200,
70  min=2,
71  )
72  minMatchedPairs = pexConfig.RangeField(
73  doc="Minimum number of matched pairs; see also minFracMatchedPairs.",
74  dtype=int,
75  default=30,
76  min=2,
77  )
78  minFracMatchedPairs = pexConfig.RangeField(
79  doc="Minimum number of matched pairs as a fraction of the smaller of "
80  "the number of reference stars or the number of good sources; "
81  "the actual minimum is the smaller of this value or "
82  "minMatchedPairs.",
83  dtype=float,
84  default=0.3,
85  min=0,
86  max=1,
87  )
88  matcherIterations = pexConfig.RangeField(
89  doc="Number of softening iterations in matcher.",
90  dtype=int,
91  default=3,
92  min=1,
93  )
94  maxOffsetPix = pexConfig.RangeField(
95  doc="Maximum allowed shift of WCS, due to matching (pixel). "
96  "When changing this value, the "
97  "LoadReferenceObjectsConfig.pixelMargin should also be updated.",
98  dtype=int,
99  default=300,
100  max=4000,
101  )
102  maxRotationDeg = pexConfig.RangeField(
103  doc="Rotation angle allowed between sources and position reference "
104  "objects (degrees).",
105  dtype=float,
106  default=1.0,
107  max=6.0,
108  )
109  numPointsForShape = pexConfig.Field(
110  doc="Number of points to define a shape for matching.",
111  dtype=int,
112  default=6,
113  )
114  numPointsForShapeAttempt = pexConfig.Field(
115  doc="Number of points to try for creating a shape. This value should "
116  "be greater than or equal to numPointsForShape. Besides "
117  "loosening the signal to noise cut in the matcherSourceSelector, "
118  "increasing this number will solve CCDs where no match was found.",
119  dtype=int,
120  default=6,
121  )
122  minMatchDistPixels = pexConfig.RangeField(
123  doc="Distance in units of pixels to always consider a source-"
124  "reference pair a match. This prevents the astrometric fitter "
125  "from over-fitting and removing stars that should be matched and "
126  "allows for inclusion of new matches as the wcs improves.",
127  dtype=float,
128  default=1.0,
129  min=0.0,
130  max=6.0,
131  )
132  numPatternConsensus = pexConfig.Field(
133  doc="Number of implied shift/rotations from patterns that must agree "
134  "before it a given shift/rotation is accepted. This is only used "
135  "after the first softening iteration fails and if both the "
136  "number of reference and source objects is greater than "
137  "numBrightStars.",
138  dtype=int,
139  default=3,
140  )
141  sourceSelector = sourceSelectorRegistry.makeField(
142  doc="How to select sources for cross-matching. The default "
143  "matcherSourceSelector removes objects with low S/N, bad "
144  "saturated objects, edge objects, and interpolated objects.",
145  default="matcherPessimistic"
146  )
147 
148  def setDefaults(self):
149  sourceSelector = self.sourceSelector["matcherPessimistic"]
150  sourceSelector.setDefaults()
151 
152  def validate(self):
153  pexConfig.Config.validate(self)
155  raise ValueError("numPointsForShapeAttempt must be greater than "
156  "or equal to numPointsForShape.")
157  if self.numPointsForShape > self.numBrightStars:
158  raise ValueError("numBrightStars must be greater than "
159  "numPointsForShape.")
160 
161 
162 # The following block adds links to this task from the Task Documentation page.
163 # \addtogroup LSST_task_documentation
164 # \{
165 # \page measAstrom_MatchPessimisticBTask
166 # \ref MatchPessimisticBTask "MatchPessimisticBTask"
167 # Match sources to reference objects
168 # \}
169 
170 
171 class MatchPessimisticBTask(pipeBase.Task):
172  """Match sources to reference objects.
173  """
174 
175  ConfigClass = MatchPessimisticBConfig
176  _DefaultName = "matchObjectsToSources"
177 
178  def __init__(self, **kwargs):
179  pipeBase.Task.__init__(self, **kwargs)
180  self.makeSubtask("sourceSelector")
181 
182  @pipeBase.timeMethod
183  def matchObjectsToSources(self, refCat, sourceCat, wcs, refFluxField,
184  match_tolerance=None):
185  """Match sources to position reference stars
186 
187  refCat : `lsst.afw.table.SimpleCatalog`
188  catalog of reference objects that overlap the exposure; reads
189  fields for:
190 
191  - coord
192  - the specified flux field
193 
194  sourceCat : `lsst.afw.table.SourceCatalog`
195  catalog of sources found on an exposure; Please check the required
196  fields of your specified source selector that the correct flags are present.
197  wcs : `lsst.afw.geom.SkyWcs`
198  estimated WCS
199  refFluxField : `str`
200  field of refCat to use for flux
201  match_tolerance : `lsst.meas.astrom.MatchTolerancePessimistic`
202  is a MatchTolerance class object or `None`. This this class is used
203  to communicate state between AstrometryTask and MatcherTask.
204  AstrometryTask will also set the MatchTolerance class variable
205  maxMatchDist based on the scatter AstrometryTask has found after
206  fitting for the wcs.
207 
208  Returns
209  -------
210  result : `lsst.pipe.base.Struct`
211  Result struct with components:
212 
213  - ``matches`` : source to reference matches found (`list` of
214  `lsst.afw.table.ReferenceMatch`)
215  - ``usableSourcCat`` : a catalog of sources potentially usable for
216  matching and WCS fitting (`lsst.afw.table.SourceCatalog`).
217  - ``match_tolerance`` : a MatchTolerance object containing the
218  resulting state variables from the match
219  (`lsst.meas.astrom.MatchTolerancePessimistic`).
220  """
221  import lsstDebug
222  debug = lsstDebug.Info(__name__)
223 
224  # If we get an empty tolerance struct create the variables we need for
225  # this matcher.
226  if match_tolerance is None:
227  match_tolerance = MatchTolerancePessimistic()
228 
229  # usableSourceCat: sources that are good but may be saturated
230  numSources = len(sourceCat)
231  selectedSources = self.sourceSelector.run(sourceCat)
232  goodSourceCat = selectedSources.sourceCat
233  numUsableSources = len(goodSourceCat)
234  self.log.info("Purged %d sources, leaving %d good sources" %
235  (numSources - numUsableSources, numUsableSources))
236 
237  if len(goodSourceCat) == 0:
238  raise pipeBase.TaskError("No sources are good")
239 
240  # avoid accidentally using sourceCat; use goodSourceCat from now on
241  del sourceCat
242 
243  minMatchedPairs = min(self.config.minMatchedPairs,
244  int(self.config.minFracMatchedPairs *
245  min([len(refCat), len(goodSourceCat)])))
246 
247  doMatchReturn = self._doMatch(
248  refCat=refCat,
249  sourceCat=goodSourceCat,
250  wcs=wcs,
251  refFluxField=refFluxField,
252  numUsableSources=numUsableSources,
253  minMatchedPairs=minMatchedPairs,
254  match_tolerance=match_tolerance,
255  sourceFluxField=self.sourceSelector.fluxField,
256  verbose=debug.verbose,
257  )
258  matches = doMatchReturn.matches
259  match_tolerance = doMatchReturn.match_tolerance
260 
261  if len(matches) == 0:
262  raise RuntimeError("Unable to match sources")
263 
264  self.log.info("Matched %d sources" % len(matches))
265  if len(matches) < minMatchedPairs:
266  self.log.warn("Number of matches is smaller than request")
267 
268  return pipeBase.Struct(
269  matches=matches,
270  usableSourceCat=goodSourceCat,
271  match_tolerance=match_tolerance,
272  )
273 
274  @pipeBase.timeMethod
275  def _doMatch(self, refCat, sourceCat, wcs, refFluxField, numUsableSources,
276  minMatchedPairs, match_tolerance, sourceFluxField, verbose):
277  """Implementation of matching sources to position reference stars
278 
279  Unlike matchObjectsToSources, this method does not check if the sources
280  are suitable.
281 
282  Parameters
283  ----------
284  refCat : `lsst.afw.table.SimpleCatalog`
285  catalog of position reference stars that overlap an exposure
286  sourceCat : `lsst.afw.table.SourceCatalog`
287  catalog of sources found on the exposure
288  wcs : `lsst.afw.geom.SkyWcs`
289  estimated WCS of exposure
290  refFluxField : `str`
291  field of refCat to use for flux
292  numUsableSources : `int`
293  number of usable sources (sources with known centroid that are not
294  near the edge, but may be saturated)
295  minMatchedPairs : `int`
296  minimum number of matches
297  match_tolerance : `lsst.meas.astrom.MatchTolerancePessimistic`
298  a MatchTolerance object containing variables specifying matcher
299  tolerances and state from possible previous runs.
300  sourceFluxField : `str`
301  Name of the flux field in the source catalog.
302  verbose : `bool`
303  Set true to print diagnostic information to std::cout
304 
305  Returns
306  -------
307  result :
308  Results struct with components:
309 
310  - ``matches`` : a list the matches found
311  (`list` of `lsst.afw.table.ReferenceMatch`).
312  - ``match_tolerance`` : MatchTolerance containing updated values from
313  this fit iteration (`lsst.meas.astrom.MatchTolerancePessimistic`)
314  """
315 
316  # Load the source and reference catalog as spherical points
317  # in numpy array. We do this rather than relying on internal
318  # lsst C objects for simplicity and because we require
319  # objects contiguous in memory. We need to do these slightly
320  # differently for the reference and source cats as they are
321  # different catalog objects with different fields.
322  src_array = np.empty((len(sourceCat), 4), dtype=np.float64)
323  for src_idx, srcObj in enumerate(sourceCat):
324  coord = wcs.pixelToSky(srcObj.getCentroid())
325  theta = np.pi / 2 - coord.getLatitude().asRadians()
326  phi = coord.getLongitude().asRadians()
327  flux = srcObj.getPsfInstFlux()
328  src_array[src_idx, :] = \
329  self._latlong_flux_to_xyz_mag(theta, phi, flux)
330 
331  if match_tolerance.PPMbObj is None or \
332  match_tolerance.autoMaxMatchDist is None:
333  # The reference catalog is fixed per AstrometryTask so we only
334  # create the data needed if this is the first step in the match
335  # fit cycle.
336  ref_array = np.empty((len(refCat), 4), dtype=np.float64)
337  for ref_idx, refObj in enumerate(refCat):
338  theta = np.pi / 2 - refObj.getDec().asRadians()
339  phi = refObj.getRa().asRadians()
340  flux = refObj[refFluxField]
341  ref_array[ref_idx, :] = \
342  self._latlong_flux_to_xyz_mag(theta, phi, flux)
343  # Create our matcher object.
344  match_tolerance.PPMbObj = PessimisticPatternMatcherB(
345  ref_array[:, :3], self.log)
346  self.log.debug("Computing source statistics...")
347  maxMatchDistArcSecSrc = self._get_pair_pattern_statistics(
348  src_array)
349  self.log.debug("Computing reference statistics...")
350  maxMatchDistArcSecRef = self._get_pair_pattern_statistics(
351  ref_array)
352  maxMatchDistArcSec = np.min((maxMatchDistArcSecSrc,
353  maxMatchDistArcSecRef))
354  match_tolerance.autoMaxMatchDist = afwgeom.Angle(
355  maxMatchDistArcSec, afwgeom.arcseconds)
356 
357  # Set configurable defaults when we encounter None type or set
358  # state based on previous run of AstrometryTask._matchAndFitWcs.
359  if match_tolerance.maxShift is None:
360  maxShiftArcseconds = (self.config.maxOffsetPix *
361  wcs.getPixelScale().asArcseconds())
362  else:
363  # We don't want to clamp down too hard on the allowed shift so
364  # we test that the smallest we ever allow is the pixel scale.
365  maxShiftArcseconds = np.max(
366  (match_tolerance.maxShift.asArcseconds(),
367  self.config.minMatchDistPixels *
368  wcs.getPixelScale().asArcseconds()))
369 
370  # If our tolerances are not set from a previous run, estimate a
371  # starting tolerance guess from the statistics of patterns we can
372  # create on both the source and reference catalog. We use the smaller
373  # of the two.
374  if match_tolerance.maxMatchDist is None:
375  match_tolerance.maxMatchDist = match_tolerance.autoMaxMatchDist
376  else:
377  maxMatchDistArcSec = np.max(
378  (self.config.minMatchDistPixels *
379  wcs.getPixelScale().asArcseconds(),
380  np.min((match_tolerance.maxMatchDist.asArcseconds(),
381  match_tolerance.autoMaxMatchDist.asArcseconds()))))
382 
383  # Make sure the data we are considering is dense enough to require
384  # the consensus mode of the matcher. If not default to Optimistic
385  # pattern matcher behavior.
386  numConsensus = self.config.numPatternConsensus
387  minObjectsForConsensus = \
388  self.config.numBrightStars + self.config.numPointsForShapeAttempt
389  if len(refCat) < minObjectsForConsensus or \
390  len(sourceCat) < minObjectsForConsensus:
391  numConsensus = 1
392 
393  self.log.debug("Current tol maxDist: %.4f arcsec" %
394  maxMatchDistArcSec)
395  self.log.debug("Current shift: %.4f arcsec" %
396  maxShiftArcseconds)
397 
398  match_found = False
399  # Start the iteration over our tolerances.
400  for soften_dist in range(self.config.matcherIterations):
401  for soften_pattern in range(self.config.matcherIterations):
402  if soften_pattern == 0 and soften_dist == 0 and \
403  match_tolerance.lastMatchedPattern is not None:
404  # If we are on the first, most stringent tolerance,
405  # and have already found a match, the matcher should behave
406  # like an optimistic pattern matcher. Exiting at the first
407  # match.
408  run_n_consent = 1
409  else:
410  # If we fail or this is the first match attempt, set the
411  # pattern consensus to the specified config value.
412  run_n_consent = numConsensus
413  # We double the match dist tolerance each round and add 1 to the
414  # to the number of candidate spokes to check.
415  matcher_struct = match_tolerance.PPMbObj.match(
416  source_array=src_array,
417  n_check=self.config.numPointsForShapeAttempt + soften_pattern,
418  n_match=self.config.numPointsForShape,
419  n_agree=run_n_consent,
420  max_n_patterns=self.config.numBrightStars,
421  max_shift=maxShiftArcseconds,
422  max_rotation=self.config.maxRotationDeg,
423  max_dist=maxMatchDistArcSec * 2. ** soften_dist,
424  min_matches=minMatchedPairs,
425  pattern_skip_array=np.array(
426  match_tolerance.failedPatternList)
427  )
428 
429  if soften_pattern == 0 and soften_dist == 0 and \
430  len(matcher_struct.match_ids) == 0 and \
431  match_tolerance.lastMatchedPattern is not None:
432  # If we found a pattern on a previous match-fit iteration and
433  # can't find an optimistic match on our first try with the
434  # tolerances as found in the previous match-fit,
435  # the match we found in the last iteration was likely bad. We
436  # append the bad match's index to the a list of
437  # patterns/matches to skip on subsequent iterations.
438  match_tolerance.failedPatternList.append(
439  match_tolerance.lastMatchedPattern)
440  match_tolerance.lastMatchedPattern = None
441  maxShiftArcseconds = \
442  self.config.maxOffsetPix * wcs.getPixelScale().asArcseconds()
443  elif len(matcher_struct.match_ids) > 0:
444  # Match found, save a bit a state regarding this pattern
445  # in the match tolerance class object and exit.
446  match_tolerance.maxShift = \
447  matcher_struct.shift * afwgeom.arcseconds
448  match_tolerance.lastMatchedPattern = \
449  matcher_struct.pattern_idx
450  match_found = True
451  break
452  if match_found:
453  break
454 
455  # If we didn't find a match, exit early.
456  if not match_found:
457  return pipeBase.Struct(
458  matches=[],
459  match_tolerance=match_tolerance,
460  )
461 
462  # The matcher returns all the nearest neighbors that agree between
463  # the reference and source catalog. For the current astrometric solver
464  # we need to remove as many false positives as possible before sending
465  # the matches off to the solver. The low value of 100 and high value of
466  # 2 are the low number of sigma and high respectively. The exact values
467  # were found after testing on data of various reference/source
468  # densities and astrometric distortion quality, specifically the
469  # visits: HSC (3358), DECam (406285, 410827),
470  # CFHT (793169, 896070, 980526).
471  distances_arcsec = np.degrees(matcher_struct.distances_rad) * 3600
472  clip_max_dist = np.max(
473  (sigmaclip(distances_arcsec, low=100, high=2)[-1],
474  self.config.minMatchDistPixels * wcs.getPixelScale().asArcseconds())
475  )
476  # Assuming the number of matched objects surviving the clip_max_dist
477  # cut if greater the requested min number of pairs, we select the
478  # smaller of the current maxMatchDist or the sigma clipped distance.
479  if not np.isfinite(clip_max_dist):
480  clip_max_dist = maxMatchDistArcSec
481 
482  if clip_max_dist < maxMatchDistArcSec and \
483  len(distances_arcsec[distances_arcsec < clip_max_dist]) < \
484  minMatchedPairs:
485  dist_cut_arcsec = maxMatchDistArcSec
486  else:
487  dist_cut_arcsec = np.min((clip_max_dist, maxMatchDistArcSec))
488 
489  # A match has been found, return our list of matches and
490  # return.
491  matches = []
492  for match_id_pair, dist_arcsec in zip(matcher_struct.match_ids,
493  distances_arcsec):
494  if dist_arcsec < dist_cut_arcsec:
495  match = afwTable.ReferenceMatch()
496  match.first = refCat[match_id_pair[1]]
497  match.second = sourceCat[match_id_pair[0]]
498  # We compute the true distance along and sphere instead
499  # and store it in units of arcseconds. The previous
500  # distances we used were approximate.
501  match.distance = match.first.getCoord().separation(
502  match.second.getCoord()).asArcseconds()
503  matches.append(match)
504 
505  return pipeBase.Struct(
506  matches=matches,
507  match_tolerance=match_tolerance,
508  )
509 
510  def _latlong_flux_to_xyz_mag(self, theta, phi, flux):
511  """Convert angles theta and phi and a flux into unit sphere
512  x, y, z, and a relative magnitude.
513 
514  Takes in a afw catalog object and converts the catalog object RA, DECs
515  to points on the unit sphere. Also converts the flux into a simple,
516  non-zero-pointed magnitude for relative sorting.
517 
518  Parameters
519  ----------
520  theta : `float`
521  Angle from the north pole (z axis) of the sphere
522  phi : `float`
523  Rotation around the sphere
524 
525  Return
526  ------
527  output_array : `numpy.ndarray`, (N, 4)
528  Spherical unit vector x, y, z with flux.
529  """
530  output_array = np.empty(4, dtype=np.float64)
531  output_array[0] = np.sin(theta)*np.cos(phi)
532  output_array[1] = np.sin(theta)*np.sin(phi)
533  output_array[2] = np.cos(theta)
534  if flux > 0:
535  output_array[3] = -2.5 * np.log10(flux)
536  else:
537  # Set flux to a very faint mag if its for some reason it
538  # does not exist
539  output_array[3] = 99.
540 
541  return output_array
542 
543  def _get_pair_pattern_statistics(self, cat_array):
544  """ Compute the tolerances for the matcher automatically by comparing
545  pinwheel patterns as we would in the matcher.
546 
547  We test how similar the patterns we can create from a given set of
548  objects by computing the spoke lengths for each pattern and sorting
549  them from smallest to largest. The match tolerance is the average
550  distance per spoke between the closest two patterns in the sorted
551  spoke length space.
552 
553  Parameters
554  ----------
555  cat_array : `numpy.ndarray`, (N, 3)
556  array of 3 vectors representing the x, y, z position of catalog
557  objects on the unit sphere.
558 
559  Returns
560  -------
561  dist_tol : `float`
562  Suggested max match tolerance distance calculated from comparisons
563  between pinwheel patterns used in optimistic/pessimistic pattern
564  matcher.
565  """
566 
567  self.log.debug("Starting automated tolerance calculation...")
568 
569  # Create an empty array of all the patterns we possibly make
570  # sorting from brightest to faintest.
571  pattern_array = np.empty(
572  (cat_array.shape[0] - self.config.numPointsForShape,
573  self.config.numPointsForShape - 1))
574  flux_args_array = np.argsort(cat_array[:, -1])
575 
576  # Sort our input array.
577  tmp_sort_array = cat_array[flux_args_array]
578 
579  # Start making patterns.
580  for start_idx in range(cat_array.shape[0] -
581  self.config.numPointsForShape):
582  pattern_points = tmp_sort_array[start_idx:start_idx +
583  self.config.numPointsForShape, :-1]
584  pattern_delta = pattern_points[1:, :] - pattern_points[0, :]
585  pattern_array[start_idx, :] = np.sqrt(
586  pattern_delta[:, 0] ** 2 +
587  pattern_delta[:, 1] ** 2 +
588  pattern_delta[:, 2] ** 2)
589 
590  # When we store the length of each spoke in our pattern we
591  # sort from shortest to longest so we have a defined space
592  # to compare them in.
593  pattern_array[start_idx, :] = pattern_array[
594  start_idx, np.argsort(pattern_array[start_idx, :])]
595 
596  # Create a searchable tree object of the patterns and find
597  # for any given pattern the closest pattern in the sorted
598  # spoke length space.
599  dist_tree = cKDTree(
600  pattern_array[:, :(self.config.numPointsForShape - 1)])
601  dist_nearest_array, ids = dist_tree.query(
602  pattern_array[:, :(self.config.numPointsForShape - 1)], k=2)
603  dist_nearest_array = dist_nearest_array[:, 1]
604  dist_nearest_array.sort()
605 
606  # We use the two closest patterns to set our tolerance.
607  dist_idx = 0
608  dist_tol = (np.degrees(dist_nearest_array[dist_idx]) * 3600. /
609  (self.config.numPointsForShape - 1.))
610 
611  self.log.debug("Automated tolerance")
612  self.log.debug("\tdistance/match tol: %.4f [arcsec]" % dist_tol)
613 
614  return dist_tol
def __init__(self, maxMatchDist=None, autoMaxMatchDist=None, maxShift=None, lastMatchedPattern=None, failedPatternList=None, PPMbObj=None)
def _doMatch(self, refCat, sourceCat, wcs, refFluxField, numUsableSources, minMatchedPairs, match_tolerance, sourceFluxField, verbose)
def matchObjectsToSources(self, refCat, sourceCat, wcs, refFluxField, match_tolerance=None)