413 self, exposure, psfCandidateList, metadata=None, flagKey=None
415 """Determine a Piff PSF model for an exposure given a list of PSF
420 exposure : `lsst.afw.image.Exposure`
421 Exposure containing the PSF candidates.
422 psfCandidateList : `list` of `lsst.meas.algorithms.PsfCandidate`
423 A sequence of PSF candidates typically obtained by detecting sources
424 and then running them through a star selector.
425 metadata : `lsst.daf.base import PropertyList` or `None`, optional
426 A home for interesting tidbits of information.
427 flagKey : `str` or `None`, optional
428 Schema key used to mark sources actually used in PSF determination.
432 psf : `lsst.meas.extensions.piff.PiffPsf`
433 The measured PSF model.
435 Unused by this PsfDeterminer.
439 if self.config.stampSize:
440 stampSize = self.config.stampSize
441 if stampSize > psfCandidateList[0].getWidth():
442 self.log.warning(
"stampSize is larger than the PSF candidate size. Using candidate size.")
443 stampSize = psfCandidateList[0].getWidth()
445 self.log.debug(
"stampSize not set. Using candidate size.")
446 stampSize = psfCandidateList[0].getWidth()
448 scale = exposure.getWcs().getPixelScale(exposure.getBBox().getCenter()).asArcseconds()
450 match self.config.useCoordinates:
452 detector = exposure.getDetector()
453 pix_to_field = detector.getTransform(PIXELS, FIELD_ANGLE)
459 skyOrigin = exposure.getWcs().getSkyOrigin()
460 ra = skyOrigin.getLongitude().asDegrees()
461 dec = skyOrigin.getLatitude().asDegrees()
462 pointing = galsim.CelestialCoord(
468 gswcs = galsim.PixelScale(scale)
472 for candidate
in psfCandidateList:
473 cmi = candidate.getMaskedImage(stampSize, stampSize)
475 fracGood = np.sum(good)/good.size
476 if fracGood < self.config.minimumUnmaskedFraction:
481 bds = galsim.BoundsI(
482 galsim.PositionI(*bbox.getMin()),
483 galsim.PositionI(*bbox.getMax())
485 gsImage = galsim.Image(bds, wcs=gswcs, dtype=float)
486 gsImage.array[:] = cmi.image.array
487 gsWeight = galsim.Image(bds, wcs=gswcs, dtype=float)
488 gsWeight.array[:] = weight
490 source = candidate.getSource()
491 starId = source.getId()
492 image_pos = galsim.PositionD(source.getX(), source.getY())
494 data = piff.StarData(
500 star = piff.Star(data,
None)
501 star.data.properties[
'starId'] = starId
504 if self.config.piffPsfConfigYaml
is None:
507 if hasattr(exposure.filter,
"bandLabel"):
508 band = exposure.filter.bandLabel
511 spatialOrder = self.config.spatialOrderPerBand.get(band, self.config.spatialOrder)
516 'scale': scale * self.config.samplingSize,
517 'size': self.config.modelSize,
518 'interp': self.config.interpolant,
519 'centered': self.config.piffPixelGridFitCenter,
522 'type':
'BasisPolynomial',
523 'order': spatialOrder,
524 'solver': self.config.piffBasisPolynomialSolver,
528 'nsigma': self.config.outlierNSigma,
529 'max_remove': self.config.outlierMaxRemove,
531 'max_iter': self.config.piffMaxIter
534 piffConfig = yaml.safe_load(self.config.piffPsfConfigYaml)
536 def _get_threshold(nth_order):
538 freeParameters = ((nth_order + 1) * (nth_order + 2)) // 2
539 return freeParameters
541 if piffConfig[
'interp'][
'type']
in [
'BasisPolynomial',
'Polynomial']:
542 threshold = _get_threshold(piffConfig[
'interp'][
'order'])
543 if len(stars) < threshold:
544 if self.config.zerothOrderInterpNotEnoughStars:
546 "Only %d stars found, "
547 "but %d required. Using zeroth order interpolation."%((len(stars), threshold))
549 piffConfig[
'interp'][
'order'] = 0
552 piffConfig[
'max_iter'] = 1
555 num_good_stars=len(stars),
556 minimum_dof=threshold,
557 poly_ndim=piffConfig[
'interp'][
'order'],
561 piffResult = piff.PSF.process(piffConfig)
564 piffResult.fit(stars, wcs, pointing, logger=self.
piffLogger)
566 nUsedStars = len([s
for s
in piffResult.stars
if not s.is_flagged
and not s.is_reserve])
568 if piffConfig[
'interp'][
'type']
in [
'BasisPolynomial',
'Polynomial']:
569 threshold = _get_threshold(piffConfig[
'interp'][
'order'])
570 if nUsedStars < threshold:
571 if self.config.zerothOrderInterpNotEnoughStars:
573 "Only %d after outlier rejection, "
574 "but %d required. Using zeroth order interpolation."%((nUsedStars, threshold))
576 piffConfig[
'interp'][
'order'] = 0
579 piffConfig[
'max_iter'] = 1
580 piffResult.fit(stars, wcs, pointing, logger=self.
piffLogger)
581 nUsedStars = len(stars)
584 num_good_stars=nUsedStars,
585 minimum_dof=threshold,
586 poly_ndim=piffConfig[
'interp'][
'order'],
589 drawSize = 2*np.floor(0.5*stampSize/self.config.samplingSize) + 1
591 used_image_starId = {s.data.properties[
'starId']
for s
in piffResult.stars
592 if not s.is_flagged
and not s.is_reserve}
594 assert len(used_image_starId) == nUsedStars,
"Star IDs are not unique"
597 for candidate
in psfCandidateList:
598 source = candidate.getSource()
599 starId = source.getId()
600 if starId
in used_image_starId:
601 source.set(flagKey,
True)
603 if metadata
is not None:
604 metadata[
"spatialFitChi2"] = piffResult.chisq
605 metadata[
"numAvailStars"] = len(stars)
606 metadata[
"numGoodStars"] = nUsedStars
607 metadata[
"avgX"] = np.mean([s.x
for s
in piffResult.stars
608 if not s.is_flagged
and not s.is_reserve])
609 metadata[
"avgY"] = np.mean([s.y
for s
in piffResult.stars
610 if not s.is_flagged
and not s.is_reserve])
612 if not self.config.debugStarData:
613 for star
in piffResult.stars:
616 del star.fit.params_var
621 del star.data.orig_weight
623 return PiffPsf(drawSize, drawSize, piffResult),
None