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