Skip to content
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
75 changes: 63 additions & 12 deletions python/lsst/ip/diffim/modelPsfMatch.py
Original file line number Diff line number Diff line change
Expand Up @@ -45,6 +45,55 @@
sigma2fwhm = 2.*np.sqrt(2.*np.log(2.))


def _safeEstimateRoughShape(psf, position, jitter=2.0, maxTries=5, log=None):

"""
Try to compute the shape of `psf` at `position`.

If it fails, perturb the position until success or maxTries.

Parameters
----------
psf : `lsst.afw.detection.Psf`
The PSF object.
position : `lsst.geom.Point2D`
The nominal position at which to evaluate the PSF shape.
jitter : `float`
Size of random offset in pixels to try when retrying.
maxTries : `int`
Maximum number of attempts (including the original).
log : `lsst.utils.logging.LsstLogAdapter` optional
logger to use for logging retries.

Returns
-------
shape : lsst.afw.geom.ellipses.Quadrupole
The shape of the PSF at the (possibly perturbed) position.
"""
try:
return psf.computeShape(position)
except pexExceptions.Exception as err:
if log:
log.info(f"safeEstimateRoughShape: initial evaluation of PSF failed with message: {str(err)}."
" Retrying nearby")
firstErr = err

# random offsets must be deterministic
offsets = np.random.default_rng(seed=40).uniform(-jitter, jitter, size=(maxTries - 1, 2))
candidates = [position + geom.Extent2D(dx, dy) for dx, dy in offsets]
for i, candidate in enumerate(candidates):
try:
shape = psf.computeShape(candidate)
if log:
log.info(f"safeEstimateRoughShape succeeded on try {i + 2} at position {candidate}")
return shape
except pexExceptions.Exception:
continue

# If nothing worked, raise AlgorithmError from original position
raise PsfComputeShapeError(position) from firstErr


def nextOddInteger(x):
nextInt = int(np.ceil(x))
return nextInt + 1 if nextInt%2 == 0 else nextInt
Expand All @@ -61,14 +110,16 @@ def metadata(self) -> dict:

class PsfComputeShapeError(pipeBase.AlgorithmError):
def __init__(self, position):
message = f"Unable to compute the FWHM of the science Psf at {position}"
message = f"Unable to evaluate Psf at ({float(position.getX()):.4f}, {float(position.getY()):.4f})"
super().__init__(message)
self.position = position

@property
def metadata(self) -> dict:
return {
"position": self.position,
# must have values `str`, `int`, `float`, `bool`or nested-dict
"position": {"x": float(self.position.getX()),
"y": float(self.position.getY())}
}


Expand Down Expand Up @@ -173,15 +224,12 @@ def run(self, exposure, referencePsfModel, kernelSum=1.0):
# exposure's bounding box in DM-32756.
sciAvgPos = exposure.getPsf().getAveragePosition()
modelAvgPos = referencePsfModel.getAveragePosition()
try:
fwhmScience = exposure.getPsf().computeShape(sciAvgPos).getDeterminantRadius()*sigma2fwhm
except pexExceptions.RangeError:
raise WarpedPsfTransformTooBigError(
f"Unable to compute the FWHM of the science Psf at {sciAvgPos}"
"due to an unexpectedly large transform."
)
Comment on lines -178 to -182
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You don't seem to be catching this exception anymore. Is there a reason that we don't need to worry about this? Is pexExcepetions.RangeError a subset of pexExceptions.Exception?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

pexExceptions.Exception IS a catch all:

>>> import lsst.pex.exceptions as pexExceptions
>>> 
>>> try:
...     raise pexExceptions.RangeError("blah")
... except pexExceptions.Exception as err:
...     print("caught it")
... 
caught it

One of the reasons I sent this to you to review is I was wondering if you were attached to WarpedPsfTransformTooBigError specifically.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It's used in MakePsfMatchedWarpTask to catch errors and make them annotated errors (part of DM-50702). That's a ticket that I completed but without a lot of background. Your PR catches PsfComputeShapeError too so the same failure modes should still return an annotated error, the question is do we want to persist a more informative error about why it failed or not. I have no opinion on this.

except pexExceptions.Exception:
raise PsfComputeShapeError(sciAvgPos)

# TODO If DM-52537 eliminates the need for _safeEstimateRoughShape,
# the remove it, and replace next line with a direct call to
# computeShape in a try/except block.
fwhmScience = _safeEstimateRoughShape(exposure.getPsf(), sciAvgPos,
log=self.log).getDeterminantRadius()*sigma2fwhm
fwhmModel = referencePsfModel.computeShape(modelAvgPos).getDeterminantRadius()*sigma2fwhm

basisList = makeKernelBasisList(self.kConfig, fwhmScience, fwhmModel, metadata=self.metadata)
Expand Down Expand Up @@ -281,7 +329,10 @@ def _buildCellSet(self, exposure, referencePsfModel):
posY = sizeCellY*row + sizeCellY//2 + scienceY0
for col in range(nCellX):
posX = sizeCellX*col + sizeCellX//2 + scienceX0
widthS, heightS = sciencePsfModel.computeBBox(geom.Point2D(posX, posY)).getDimensions()
try:
widthS, heightS = sciencePsfModel.computeBBox(geom.Point2D(posX, posY)).getDimensions()
except pexExceptions.InvalidParameterError as err:
raise PsfComputeShapeError(geom.Point2D(posX, posY)) from err
widthList.append(widthS)
heightList.append(heightS)

Expand Down