diff --git a/python/lsst/pipe/tasks/calibrateImage.py b/python/lsst/pipe/tasks/calibrateImage.py index cc21d3e82..3d7fadedc 100644 --- a/python/lsst/pipe/tasks/calibrateImage.py +++ b/python/lsst/pipe/tasks/calibrateImage.py @@ -225,6 +225,51 @@ class CalibrateImageConfig(pipeBase.PipelineTaskConfig, pipelineConnections=Cali doc="Task to normalize the calibration flux (e.g. compensated tophats) " "for the bright stars used for psf estimation.", ) + max_unnorm_psf_ellipticity_per_band = pexConfig.DictField( + keytype=str, + itemtype=float, + default={ + "u": 3.8, + "g": 3.8, + "r": 3.8, + "i": 3.8, + "z": 3.8, + "y": 3.8, + }, + doc="Maximum unnormalized ellipticity (defined as hypot(Ixx - Iyy, 2Ixy)) of the PSF model " + "deemed good enough for further consideration. Values above this threshold raise " + "UnprocessableDataError.", + ) + max_unnorm_psf_ellipticity_fallback = pexConfig.Field( + dtype=float, + default=3.8, + doc="Fallback maximum unnormalized ellipticity (defined as hypot(Ixx - Iyy, 2Ixy)) of the " + "PSF model deemed good enough for further consideration if the current band is not in " + "the config.max_unnorm-psf_ellipticity_per_band dict. Values above this threshold " + "raise UnprocessableDataError.", + ) + max_psf_ellipticity_per_band = pexConfig.DictField( + keytype=str, + itemtype=float, + default={ + "u": 0.33, + "g": 0.32, + "r": 0.35, + "i": 0.35, + "z": 0.37, + "y": 0.32, + }, + doc="Value of the PSF model ellipticity deemed good enough for further consideration, " + "regardless of the value of the unnormalized PSF model ellipticity. Values above " + "this threshold raise UnprocessableDataError.", + ) + max_psf_ellipticity_fallback = pexConfig.Field( + dtype=float, + default=0.35, + doc="Fallback maximum ellipticity of the PSF model deemed good enough for further " + "consideration if the current band is not in the config.max_psf_ellipticity_per_band " + "dict. Values above this threshold raise UnprocessableDataError.", + ) # TODO DM-39203: we can remove aperture correction from this task once we are # using the shape-based star/galaxy code. @@ -630,6 +675,13 @@ def run(self, *, exposures, id_generator=None, result=None): ``photometry_matches`` Reference catalog stars matches used in the photometric fit. (`list` [`lsst.afw.table.ReferenceMatch`] or `lsst.afw.table.BaseCatalog`) + + Raises + ------ + UnprocessableDataError + Raised if the unnormalized model PSF ellipticity is greater than + max_unnorm_psf_ellipticity or the model PSF ellipticity is greater + than max_psf_ellipticity. """ if result is None: result = pipeBase.Struct() @@ -732,34 +784,84 @@ def _compute_psf(self, exposure, id_generator): cell_set : `lsst.afw.math.SpatialCellSet` PSF candidates returned by the psf determiner. """ - def log_psf(msg, addToMetadata=False): + def log_and_validate_psf(msg, addToMetadata=False, doValidation=False): """Log the parameters of the psf and background, with a prepended message. There is also the option to add the PSF sigma to the task metadata. + Also validate the PSF model based on its ellipticity paramerters + and the thresholds set in the config (see docs for descriptions). + If any model value exceeds its associated threshold, an + UnprocessableDataError is raised as this exposure is now deemed + too low quality for further consideration. + Parameters ---------- - msg : `str` - Message to prepend the log info with. - addToMetadata : `bool`, optional - Whether to add the final psf sigma value to the task metadata - (the default is False). + msg : `str` + Message to prepend the log info with. + addToMetadata : `bool`, optional + Whether to add the final psf sigma value to the task metadata + (the default is False). + doValidation : `bool` + Whether to perform the PSF model ellipticity based validation. + + Raises + ------ + UnprocessableDataError + Raised if the unnormalized model PSF ellipticity is greater than + max_unnorm_psf_ellipticity or the model PSF ellipticity is greater + than max_psf_ellipticity. """ - position = exposure.psf.getAveragePosition() - sigma = exposure.psf.computeShape(position).getDeterminantRadius() - dimensions = exposure.psf.computeImage(position).getDimensions() + psf = exposure.psf + average_position = psf.getAveragePosition() + psf_shape = psf.computeShape(average_position) + psf_sigma = psf_shape.getDeterminantRadius() + psf_e1 = (psf_shape.getIxx() - psf_shape.getIyy())/(psf_shape.getIxx() + psf_shape.getIyy()) + psf_e2 = 2.0*psf_shape.getIxy()/(psf_shape.getIxx() + psf_shape.getIyy()) + psf_e = np.sqrt(psf_e1**2.0 + psf_e2**2.0) + psf_unnorm_e = np.hypot(psf_shape.getIxx() - psf_shape.getIyy(), 2.0*psf_shape.getIxy()) + dimensions = psf.computeImage(average_position).getDimensions() median_background = np.median(background.getImage().array) - self.log.info("%s sigma=%0.4f, dimensions=%s; median background=%0.2f", - msg, sigma, dimensions, median_background) + self.log.info("%s sigma=%0.4f, psfE=%.3f, unNormPsfE=%.2f, dimensions=%s, median " + "background=%0.2f", + msg, psf_sigma, psf_e, psf_unnorm_e, dimensions, median_background) + + # Validate the PSF model based on its ellipticity properties + band = exposure.filter.bandLabel + if band in self.config.max_unnorm_psf_ellipticity_per_band: + max_unnorm_psf_ellipticity = self.config.max_unnorm_psf_ellipticity_per_band[band] + else: + max_unnorm_psf_ellipticity = self.config.max_unnorm_psf_ellipticity_fallback + self.log.warning( + f"Band {band} was not included in self.config.max_unnorm_psf_ellipticity_per_band. " + f"Setting max_unnorm_psf_ellipticity to fallback value of {max_unnorm_psf_ellipticity}." + ) + if band in self.config.max_psf_ellipticity_per_band: + max_psf_ellipticity = self.config.max_psf_ellipticity_per_band[band] + else: + max_psf_ellipticity = self.config.max_psf_ellipticity_fallback + self.log.warning( + f"Band {band} was not included in self.config.max_psf_ellipticity_per_band. " + f"Setting max_unnorm_psf_ellipticity to fallback value of {max_psf_ellipticity:.2f}." + ) + + if psf_unnorm_e > max_unnorm_psf_ellipticity or psf_e > max_psf_ellipticity: + raise pipeBase.UnprocessableDataError( + "Either the unnormalized model PSF ellipticity is greater than the maximum allowed " + f"for band {band} ({max_unnorm_psf_ellipticity:.2f}) or the model PSF ellipticity " + f"is greater than the maximum allowed({max_psf_ellipticity:.2f} )" + f"(unNormPsfE={psf_unnorm_e:.2f}, psfE={psf_e:.2f})" + ) + if addToMetadata: - self.metadata["final_psf_sigma"] = sigma + self.metadata["final_psf_sigma"] = psf_sigma self.log.info("First pass detection with Guassian PSF FWHM=%s pixels", self.config.install_simple_psf.fwhm) self.install_simple_psf.run(exposure=exposure) background = self.psf_subtract_background.run(exposure=exposure).background - log_psf("Initial PSF:") + log_and_validate_psf("Initial PSF:") self.psf_repair.run(exposure=exposure, keepCRs=True) table = afwTable.SourceTable.make(self.psf_schema, id_generator.make_table_id_factory()) @@ -776,7 +878,7 @@ def log_psf(msg, addToMetadata=False): # repair/detect/measure/measure_psf step: this can help it converge. self.install_simple_psf.run(exposure=exposure) - log_psf("Rerunning with simple PSF:") + log_and_validate_psf("Rerunning with simple PSF:") # TODO investigation: Should we only re-run repair here, to use the # new PSF? Maybe we *do* need to re-run measurement with PsfFlux, to # use the fitted PSF? @@ -795,7 +897,7 @@ def log_psf(msg, addToMetadata=False): self.psf_source_measurement.run(detections.sources, exposure) psf_result = self.psf_measure_psf.run(exposure=exposure, sources=detections.sources) - log_psf("Final PSF:", addToMetadata=True) + log_and_validate_psf("Final PSF:", addToMetadata=True, doValidation=True) # Final repair with final PSF, removing cosmic rays this time. self.psf_repair.run(exposure=exposure) diff --git a/python/lsst/pipe/tasks/characterizeImage.py b/python/lsst/pipe/tasks/characterizeImage.py index 6dbc428e6..f206717f8 100644 --- a/python/lsst/pipe/tasks/characterizeImage.py +++ b/python/lsst/pipe/tasks/characterizeImage.py @@ -131,6 +131,52 @@ class CharacterizeImageConfig(pipeBase.PipelineTaskConfig, "estimate PSF. If useSimplePsf is True then 2 should be plenty; " "otherwise more may be wanted.", ) + maxUnNormPsfEllipticityPerBand = pexConfig.DictField( + keytype=str, + itemtype=float, + default={ + "u": 3.8, + "g": 3.8, + "r": 3.8, + "i": 3.8, + "z": 3.8, + "y": 3.8, + }, + doc="Maximum unnormalized ellipticity (defined as hypot(Ixx - Iyy, 2Ixy)) of the PSF model " + "deemed good enough for further consideration. Values above this threshold raise " + "UnprocessableDataError.", + ) + maxUnNormPsfEllipticityFallback = pexConfig.Field( + dtype=float, + default=3.8, + doc="Fallback maximum unnormalized ellipticity (defined as hypot(Ixx - Iyy, 2Ixy)) of the " + "PSF model deemed good enough for further consideration if the current band is not in " + "the config.maxUnNormPsfEllipticityPerBand dict. Values above this threshold " + "raise UnprocessableDataError.", + ) + maxPsfEllipticityPerBand = pexConfig.DictField( + keytype=str, + itemtype=float, + default={ + "u": 0.33, + "g": 0.32, + "r": 0.35, + "i": 0.35, + "z": 0.37, + "y": 0.32, + }, + doc="Value of the PSF model ellipticity deemed good enough for further consideration, " + "regardless of the value of the unnormalized PSF model ellipticity. Values above " + "this threshold raise UnprocessableDataError.", + ) + maxPsfEllipticityFallback = pexConfig.Field( + dtype=float, + default=0.35, + doc="Fallback maximum ellipticity of the PSF model deemed good enough for further " + "consideration if the current band is not in the config.maxPsfEllipticityPerBand " + "dict. Values above this threshold raise UnprocessableDataError.", + ) + background = pexConfig.ConfigurableField( target=SubtractBackgroundTask, doc="Configuration for initial background estimation", @@ -174,8 +220,9 @@ class CharacterizeImageConfig(pipeBase.PipelineTaskConfig, target=ApplyApCorrTask, doc="Subtask to apply aperture corrections" ) - # If doApCorr is False, and the exposure does not have apcorrections already applied, the - # active plugins in catalogCalculation almost certainly should not contain the characterization plugin + # If doApCorr is False, and the exposure does not have apcorrections + # already applied, the active plugins in catalogCalculation almost + # certainly should not contain the characterization plugin. catalogCalculation = pexConfig.ConfigurableField( target=CatalogCalculationTask, doc="Subtask to run catalogCalculation plugins on catalog" @@ -241,8 +288,8 @@ class CharacterizeImageConfig(pipeBase.PipelineTaskConfig, def setDefaults(self): super().setDefaults() # Just detect bright stars. - # The thresholdValue sets the minimum flux in a pixel to be included in the - # footprint, while peaks are only detected when they are above + # The thresholdValue sets the minimum flux in a pixel to be included + # in the footprint, while peaks are only detected when they are above # thresholdValue * includeThresholdMultiplier. The low thresholdValue # ensures that the footprints are large enough for the noise replacer # to mask out faint undetected neighbors that are not to be measured. @@ -250,18 +297,19 @@ def setDefaults(self): self.detection.includeThresholdMultiplier = 10.0 # do not deblend, as it makes a mess self.doDeblend = False - # measure and apply aperture correction; note: measuring and applying aperture - # correction are disabled until the final measurement, after PSF is measured + # Measure and apply aperture correction; note: measuring and applying + # aperture correction are disabled until the final measurement, after + # PSF is measured. self.doApCorr = True - # During characterization, we don't have full source measurement information, - # so must do the aperture correction with only psf stars, combined with the - # default signal-to-noise cuts in MeasureApCorrTask. + # During characterization, we don't have full source measurement + # information, so must do the aperture correction with only psf stars, + # combined with the default signal-to-noise cuts in MeasureApCorrTask. selector = self.measureApCorr.sourceSelector["science"] selector.doUnresolved = False selector.flags.good = ["calib_psf_used"] selector.flags.bad = [] - # minimal set of measurements needed to determine PSF + # Minimal set of measurements needed to determine PSF. self.measurement.plugins.names = [ "base_PixelFlags", "base_SdssCentroid", @@ -278,7 +326,7 @@ def validate(self): if self.doApCorr and not self.measurePsf: raise RuntimeError("Must measure PSF to measure aperture correction, " "because flags determined by PSF measurement are used to identify " - "sources used to measure aperture correction") + "sources used to measure aperture correction.") class CharacterizeImageTask(pipeBase.PipelineTask): @@ -305,21 +353,25 @@ class CharacterizeImageTask(pipeBase.PipelineTask): CharacterizeImageTask has a debug dictionary with the following keys: frame - int: if specified, the frame of first debug image displayed (defaults to 1) + int: if specified, the frame of first debug image displayed (defaults + to 1). repair_iter - bool; if True display image after each repair in the measure PSF loop + bool; if True display image after each repair in the measure PSF loop. background_iter - bool; if True display image after each background subtraction in the measure PSF loop + bool; if True display image after each background subtraction in the + measure PSF loop. measure_iter - bool; if True display image and sources at the end of each iteration of the measure PSF loop - See `~lsst.meas.astrom.displayAstrometry` for the meaning of the various symbols. + bool; if True display image and sources at the end of each iteration + of the measure PSF loop. See `~lsst.meas.astrom.displayAstrometry` for + the meaning of the various symbols. psf bool; if True display image and sources after PSF is measured; - this will be identical to the final image displayed by measure_iter if measure_iter is true + this will be identical to the final image displayed by measure_iter if + measure_iter is true. repair - bool; if True display image and sources after final repair + bool; if True display image and sources after final repair. measure - bool; if True display image and sources after final measurement + bool; if True display image and sources after final measurement. """ ConfigClass = CharacterizeImageConfig @@ -367,8 +419,10 @@ def run(self, exposure, background=None, idGenerator=None): """Characterize a science image. Peforms the following operations: - - Iterate the following config.psfIterations times, or once if config.doMeasurePsf false: - - detect and measure sources and estimate PSF (see detectMeasureAndEstimatePsf for details) + - Iterate the following config.psfIterations times, or once if + config.doMeasurePsf false: + - detect and measure sources and estimate PSF (see + detectMeasureAndEstimatePsf for details) - interpolate over cosmic rays - perform final measurement @@ -393,7 +447,8 @@ def run(self, exposure, background=None, idGenerator=None): ``background`` Model of subtracted background (`lsst.afw.math.BackgroundList`). ``psfCellSet`` - Spatial cells of PSF candidates (`lsst.afw.math.SpatialCellSet`). + Spatial cells of PSF candidates + (`lsst.afw.math.SpatialCellSet`). ``characterized`` Another reference to ``exposure`` for compatibility. ``backgroundModel`` @@ -403,6 +458,10 @@ def run(self, exposure, background=None, idGenerator=None): ------ RuntimeError Raised if PSF sigma is NaN. + UnprocessableDataError + Raised if the unnormalized model PSF ellipticity is greater than + maxUnNormPsfEllipticity or the model PSF ellipticity is greater + than maxPsfEllipticity. """ self._frame = self._initialFrame # reset debug display frame @@ -423,31 +482,63 @@ def run(self, exposure, background=None, idGenerator=None): idGenerator=idGenerator, background=background, ) - psf = dmeRes.exposure.getPsf() # Just need a rough estimate; average positions are fine psfAvgPos = psf.getAveragePosition() - psfSigma = psf.computeShape(psfAvgPos).getDeterminantRadius() + psfShape = psf.computeShape(psfAvgPos) + psfSigma = psfShape.getDeterminantRadius() + psfE1 = (psfShape.getIxx() - psfShape.getIyy())/(psfShape.getIxx() + psfShape.getIyy()) + psfE2 = 2.0*psfShape.getIxy()/(psfShape.getIxx() + psfShape.getIyy()) + psfE = np.sqrt(psfE1**2.0 + psfE2**2.0) + unNormPsfE = np.hypot(psfShape.getIxx() - psfShape.getIyy(), 2.0*psfShape.getIxy()) + psfDimensions = psf.computeImage(psfAvgPos).getDimensions() medBackground = np.median(dmeRes.background.getImage().getArray()) - self.log.info("iter %s; PSF sigma=%0.4f, dimensions=%s; median background=%0.2f", - i + 1, psfSigma, psfDimensions, medBackground) + self.log.info( + "iter %s: PSF sigma=%0.4f, psfE=%.3f, unNormPsfE=%.2f, dimensions=%s, " + "median background=%0.2f", + i + 1, psfSigma, psfE, unNormPsfE, psfDimensions, medBackground) if np.isnan(psfSigma): raise RuntimeError("PSF sigma is NaN, cannot continue PSF determination.") + band = exposure.filter.bandLabel + if band in self.config.maxUnNormPsfEllipticityPerBand: + maxUnNormPsfEllipticity = self.config.maxUnNormPsfEllipticityPerBand[band] + else: + maxUnNormPsfEllipticity = self.config.maxUnNormPsfEllipticityFallback + self.log.warning( + f"Band {band} was not included in self.config.maxUnNormPsfEllipticityPerBand. " + f"Setting maxUnNormPsfEllipticity to fallback value of {maxUnNormPsfEllipticity}." + ) + if band in self.config.maxPsfEllipticityPerBand: + maxPsfEllipticity = self.config.maxPsfEllipticityPerBand[band] + else: + maxPsfEllipticity = self.config.maxPsfEllipticityFallback + self.log.warning( + f"Band {band} was not included in self.config.maxPsfEllipticityPerBand. " + f"Setting maxUnNormPsfEllipticity to fallback value of {maxPsfEllipticity:.2f}." + ) + + if unNormPsfE > maxUnNormPsfEllipticity or psfE > maxPsfEllipticity: + raise pipeBase.UnprocessableDataError( + "Either the unnormalized model PSF ellipticity is greater than the maximum allowed " + f"for band {band} ({maxUnNormPsfEllipticity:.2f}) or the model PSF ellipticity " + f"is greater than the maximum allowed ({maxPsfEllipticity:.2f}) " + f"(unNormPsfE={unNormPsfE:.2f}, psfE={psfE:.2f})" + ) self.display("psf", exposure=dmeRes.exposure, sourceCat=dmeRes.sourceCat) - # perform final repair with final PSF + # Perform final repair with final PSF. self.repair.run(exposure=dmeRes.exposure) self.display("repair", exposure=dmeRes.exposure, sourceCat=dmeRes.sourceCat) - # mask streaks + # Mask streaks. # TODO: Remove in DM-44658, streak masking to happen only in ip_diffim if self.config.doMaskStreaks: _ = self.maskStreaks.run(dmeRes.exposure) - # perform final measurement with final PSF, including measuring and applying aperture correction, - # if wanted + # Perform final measurement with final PSF, including measuring and + # applying aperture correction, if wanted. self.measurement.run(measCat=dmeRes.sourceCat, exposure=dmeRes.exposure, exposureId=idGenerator.catalog_id) @@ -475,7 +566,8 @@ def run(self, exposure, background=None, idGenerator=None): # downstream. dmeRes.exposure.info.setApCorrMap(None) else: - # Need to merge the aperture correction map from the normalization. + # Need to merge the aperture correction map from the + # normalization. if normApCorrMap: for key in normApCorrMap: apCorrMap[key] = normApCorrMap[key] @@ -503,12 +595,13 @@ def detectMeasureAndEstimatePsf(self, exposure, idGenerator, background): Performs the following operations: - if config.doMeasurePsf or not exposure.hasPsf(): - - - install a simple PSF model (replacing the existing one, if need be) + - install a simple PSF model (replacing the existing one, if + need be) - interpolate over cosmic rays with keepCRs=True - estimate background and subtract it from the exposure - - detect, deblend and measure sources, and subtract a refined background model; + - detect, deblend and measure sources, and subtract a refined + background model; - if config.doMeasurePsf: - measure PSF @@ -533,19 +626,21 @@ def detectMeasureAndEstimatePsf(self, exposure, idGenerator, background): ``background`` Model of subtracted background (`lsst.afw.math.BackgroundList`). ``psfCellSet`` - Spatial cells of PSF candidates (`lsst.afw.math.SpatialCellSet`). + Spatial cells of PSF candidates + (`lsst.afw.math.SpatialCellSet`). Raises ------ LengthError Raised if there are too many CR pixels. """ - # install a simple PSF model, if needed or wanted + # Install a simple PSF model, if needed or wanted. if not exposure.hasPsf() or (self.config.doMeasurePsf and self.config.useSimplePsf): self.log.info("PSF estimation initialized with 'simple' PSF") self.installSimplePsf.run(exposure=exposure) - # run repair, but do not interpolate over cosmic rays (do that elsewhere, with the final PSF model) + # Run repair, but do not interpolate over cosmic rays (do that + # elsewhere, with the final PSF model). if self.config.requireCrForPsf: self.repair.run(exposure=exposure, keepCRs=True) else: @@ -572,7 +667,7 @@ def detectMeasureAndEstimatePsf(self, exposure, idGenerator, background): if self.config.doDeblend: self.deblend.run(exposure=exposure, sources=sourceCat) - # We need the output catalog to be contiguous for further processing. + # The output catalog needs to be contiguous for further processing. if not sourceCat.isContiguous(): sourceCat = sourceCat.copy(deep=True) diff --git a/tests/test_calibrate.py b/tests/test_calibrate.py index b5bfef2e4..a5b2ca716 100755 --- a/tests/test_calibrate.py +++ b/tests/test_calibrate.py @@ -202,8 +202,10 @@ class CalibrateTaskTestCaseNoButler(lsst.utils.tests.TestCase): def testNoAperCorrMap(self): expPath = os.path.join(TESTDIR, "data", "v695833-e0-c000-a00.sci.fits") exposure = lsst.afw.image.ExposureF(expPath) + band = exposure.filter.bandLabel charImConfig = CharacterizeImageConfig() + charImConfig.maxUnNormPsfEllipticityPerBand[band] = 3.0 charImConfig.measurePsf.psfDeterminer = 'piff' charImConfig.measurePsf.psfDeterminer['piff'].spatialOrder = 0 charImConfig.measureApCorr.sourceSelector["science"].doSignalToNoise = False diff --git a/tests/test_psfCandidateSelection.py b/tests/test_psfCandidateSelection.py index fd5bc9328..f1e7fa5f3 100755 --- a/tests/test_psfCandidateSelection.py +++ b/tests/test_psfCandidateSelection.py @@ -36,6 +36,7 @@ def setUp(self): # Load sample input from disk expPath = os.path.join(TESTDIR, "data", "v695833-e0-c000-a00.sci.fits") self.exposure = afwImage.ExposureF(expPath) + self.band = self.exposure.filter.bandLabel def tearDown(self): del self.exposure @@ -44,6 +45,7 @@ def testFlags(self): # Test that all of the flags are defined and there is no reservation by default # also test that used sources are a subset of candidate sources config = CharacterizeImageConfig() + config.maxUnNormPsfEllipticityPerBand[self.band] = 3.0 config.measurePsf.psfDeterminer = 'piff' config.measurePsf.psfDeterminer['piff'].spatialOrder = 0 config.measureApCorr.sourceSelector["science"].doSignalToNoise = False diff --git a/tests/test_unNormPsfEllipticity.py b/tests/test_unNormPsfEllipticity.py new file mode 100755 index 000000000..0c9bb05a6 --- /dev/null +++ b/tests/test_unNormPsfEllipticity.py @@ -0,0 +1,104 @@ +# This file is part of pipe_tasks. +# +# Developed for the LSST Data Management System. +# This product includes software developed by the LSST Project +# (https://www.lsst.org). +# See the COPYRIGHT file at the top-level directory of this distribution +# for details of code ownership. +# +# This program is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# This program is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with this program. If not, see . + +import os +import unittest + +import lsst.afw.image as afwImage +import lsst.meas.extensions.piff.piffPsfDeterminer +import lsst.pipe.base as pipeBase +import lsst.utils.tests +from lsst.pipe.tasks.characterizeImage import CharacterizeImageTask, CharacterizeImageConfig + +TESTDIR = os.path.abspath(os.path.dirname(__file__)) + + +class UnNormPsfEllipticityTestCase(lsst.utils.tests.TestCase): + + def setUp(self): + # Load sample input from disk + expPath = os.path.join(TESTDIR, "data", "v695833-e0-c000-a00.sci.fits") + self.exposure = afwImage.ExposureF(expPath) + self.band = self.exposure.filter.bandLabel + + def tearDown(self): + del self.exposure + + def testUnNormPsfEllipticity(self): + """Tests that UnprocessableDataError is raised when the unnormalized + PSF model ellipticity excedes config.maxUnNormPsfEllipticity. + """ + self._checkUnNormPsfEllipticity(allowPass=True) + self._checkUnNormPsfEllipticity(allowPass=False) + + def testUnNormPsfEllipticityFallback(self): + """Tests that the fallback value is set when band is missing from + config.maxUnNormPsfEllipticityPerBand. + """ + charImConfig = CharacterizeImageConfig() + charImConfig.measurePsf.psfDeterminer = 'piff' + charImConfig.measurePsf.psfDeterminer['piff'].spatialOrder = 0 + charImConfig.measureApCorr.sourceSelector["science"].doSignalToNoise = False + # Pop the band entry from the config dict to impose need for fallback. + tempDict = charImConfig.maxUnNormPsfEllipticityPerBand + tempDict.pop(self.band) + charImConfig.maxUnNormPsfEllipticityPerBand = tempDict + charImTask = CharacterizeImageTask(config=charImConfig) + charImTask.run(self.exposure) + with self.assertRaises(pipeBase.UnprocessableDataError): + charImConfig.maxUnNormPsfEllipticityFallback = 0.0 + charImTask = CharacterizeImageTask(config=charImConfig) + charImTask.run(self.exposure) + + def _checkUnNormPsfEllipticity(self, allowPass): + """Check unnormalized model PSF ellipticity threshold functionality. + + Parameters + ---------- + allowPass : `bool` + Whether to update from the default config to allow this exporsure + to pass the threshold check. + """ + charImConfig = CharacterizeImageConfig() + charImConfig.measurePsf.psfDeterminer = 'piff' + charImConfig.measurePsf.psfDeterminer['piff'].spatialOrder = 0 + charImConfig.measureApCorr.sourceSelector["science"].doSignalToNoise = False + if not allowPass: + charImConfig.maxUnNormPsfEllipticityPerBand[self.band] = 0.0 + charImTask = CharacterizeImageTask(config=charImConfig) + if allowPass: + charImTask.run(self.exposure) + else: + with self.assertRaises(pipeBase.UnprocessableDataError): + charImTask.run(self.exposure) + + +class MemoryTester(lsst.utils.tests.MemoryTestCase): + pass + + +def setup_module(module): + lsst.utils.tests.init() + + +if __name__ == "__main__": + lsst.utils.tests.init() + unittest.main()