348 minMatchedPairs, matchTolerance, sourceFluxField, verbose):
349 """Implementation of matching sources to position reference objects
350
351 Unlike matchObjectsToSources, this method does not check if the sources
352 are suitable.
353
354 Parameters
355 ----------
356 refCat : `lsst.afw.table.SimpleCatalog`
357 catalog of position reference objects that overlap an exposure
358 sourceCat : `lsst.afw.table.SourceCatalog`
359 catalog of sources found on the exposure
360 wcs : `lsst.afw.geom.SkyWcs`
361 estimated WCS of exposure
362 refFluxField : `str`
363 field of refCat to use for flux
364 numUsableSources : `int`
365 number of usable sources (sources with known centroid that are not
366 near the edge, but may be saturated)
367 minMatchedPairs : `int`
368 minimum number of matches
369 matchTolerance : `lsst.meas.astrom.MatchTolerancePessimistic`
370 a MatchTolerance object containing variables specifying matcher
371 tolerances and state from possible previous runs.
372 sourceFluxField : `str`
373 Name of the flux field in the source catalog.
374 verbose : `bool`
375 Set true to print diagnostic information to std::cout
376
377 Returns
378 -------
379 result :
380 Results struct with components:
381
382 - ``matches`` : a list the matches found
383 (`list` of `lsst.afw.table.ReferenceMatch`).
384 - ``matchTolerance`` : MatchTolerance containing updated values from
385 this fit iteration (`lsst.meas.astrom.MatchTolerancePessimistic`)
386 """
387
388
389
390
391
392
393
394 src_array = np.empty((len(sourceCat), 4), dtype=np.float64)
395 for src_idx, srcObj in enumerate(sourceCat):
396 coord = wcs.pixelToSky(srcObj.getCentroid())
397 theta = np.pi / 2 - coord.getLatitude().asRadians()
398 phi = coord.getLongitude().asRadians()
399 flux = srcObj[sourceFluxField]
400 src_array[src_idx, :] = \
401 self._latlong_flux_to_xyz_mag(theta, phi, flux)
402
403 if matchTolerance.PPMbObj is None or \
404 matchTolerance.autoMaxMatchDist is None:
405
406
407
408 ref_array = np.empty((len(refCat), 4), dtype=np.float64)
409 for ref_idx, refObj in enumerate(refCat):
410 theta = np.pi / 2 - refObj.getDec().asRadians()
411 phi = refObj.getRa().asRadians()
412 flux = refObj[refFluxField]
413 ref_array[ref_idx, :] = \
414 self._latlong_flux_to_xyz_mag(theta, phi, flux)
415
416 matchTolerance.PPMbObj = PessimisticPatternMatcherB(
417 ref_array[:, :3], self.log)
418 self.log.debug("Computing source statistics...")
419 maxMatchDistArcSecSrc = self._get_pair_pattern_statistics(
420 src_array)
421 self.log.debug("Computing reference statistics...")
422 maxMatchDistArcSecRef = self._get_pair_pattern_statistics(
423 ref_array)
424 maxMatchDistArcSec = np.max((
425 self.config.minMatchDistPixels
426 * wcs.getPixelScale().asArcseconds(),
427 np.min((maxMatchDistArcSecSrc,
428 maxMatchDistArcSecRef))))
430 maxMatchDistArcSec, geom.arcseconds)
431
432
433
434 if matchTolerance.maxShift is None:
435 maxShiftArcseconds = (self.config.maxOffsetPix
436 * wcs.getPixelScale().asArcseconds())
437 else:
438
439
440 maxShiftArcseconds = np.max(
441 (matchTolerance.maxShift.asArcseconds(),
442 self.config.minMatchDistPixels
443 * wcs.getPixelScale().asArcseconds()))
444
445
446
447
448
449 if matchTolerance.maxMatchDist is None:
450 matchTolerance.maxMatchDist = matchTolerance.autoMaxMatchDist
451 else:
452 maxMatchDistArcSec = np.max(
453 (self.config.minMatchDistPixels
454 * wcs.getPixelScale().asArcseconds(),
455 np.min((matchTolerance.maxMatchDist.asArcseconds(),
456 matchTolerance.autoMaxMatchDist.asArcseconds()))))
457
458
459
460
461
462 numConsensus = self.config.numPatternConsensus
463 if len(refCat) < self.config.numRefRequireConsensus:
464 minObjectsForConsensus = \
465 self.config.numBrightStars + \
466 self.config.numPointsForShapeAttempt
467 if len(refCat) < minObjectsForConsensus or \
468 len(sourceCat) < minObjectsForConsensus:
469 numConsensus = 1
470
471 self.log.debug("Current tol maxDist: %.4f arcsec",
472 maxMatchDistArcSec)
473 self.log.debug("Current shift: %.4f arcsec",
474 maxShiftArcseconds)
475
476 match_found = False
477
478 for soften_dist in range(self.config.matcherIterations):
479 if soften_dist == 0 and \
480 matchTolerance.lastMatchedPattern is not None:
481
482
483
484
485 run_n_consent = 1
486 else:
487
488
489 run_n_consent = numConsensus
490
491
492 matcher_struct = matchTolerance.PPMbObj.match(
493 source_array=src_array,
494 n_check=self.config.numPointsForShapeAttempt,
495 n_match=self.config.numPointsForShape,
496 n_agree=run_n_consent,
497 max_n_patterns=self.config.numBrightStars,
498 max_shift=maxShiftArcseconds,
499 max_rotation=self.config.maxRotationDeg,
500 max_dist=maxMatchDistArcSec * 2. ** soften_dist,
501 min_matches=minMatchedPairs,
502 pattern_skip_array=np.array(
503 matchTolerance.failedPatternList)
504 )
505
506 if soften_dist == 0 and \
507 len(matcher_struct.match_ids) == 0 and \
508 matchTolerance.lastMatchedPattern is not None:
509
510
511
512
513
514
515 matchTolerance.failedPatternList.append(
516 matchTolerance.lastMatchedPattern)
517 matchTolerance.lastMatchedPattern = None
518 maxShiftArcseconds = \
519 self.config.maxOffsetPix * wcs.getPixelScale().asArcseconds()
520 elif len(matcher_struct.match_ids) > 0:
521
522
523 matchTolerance.maxShift = \
524 matcher_struct.shift * geom.arcseconds
525 matchTolerance.lastMatchedPattern = \
526 matcher_struct.pattern_idx
527 match_found = True
528 break
529
530
531 if not match_found:
532 return pipeBase.Struct(
533 matches=[],
534 matchTolerance=matchTolerance,
535 )
536
537
538
539
540
541
542
543
544
545
546 distances_arcsec = np.degrees(matcher_struct.distances_rad) * 3600
547 dist_cut_arcsec = np.max(
548 (np.degrees(matcher_struct.max_dist_rad) * 3600,
549 self.config.minMatchDistPixels * wcs.getPixelScale().asArcseconds()))
550
551
552
553 matches = []
554 for match_id_pair, dist_arcsec in zip(matcher_struct.match_ids,
555 distances_arcsec):
556 if dist_arcsec < dist_cut_arcsec:
558 match.first = refCat[int(match_id_pair[1])]
559 match.second = sourceCat[int(match_id_pair[0])]
560
561
562
563 match.distance = match.first.getCoord().separation(
564 match.second.getCoord()).asArcseconds()
565 matches.append(match)
566
567 return pipeBase.Struct(
568 matches=matches,
569 matchTolerance=matchTolerance,
570 )
571
Tag types used to declare specialized field types.
A class representing an angle.