diff --git a/python/lsst/pipe/tasks/calibrateImage.py b/python/lsst/pipe/tasks/calibrateImage.py index 4275453642..8bd0c588d6 100644 --- a/python/lsst/pipe/tasks/calibrateImage.py +++ b/python/lsst/pipe/tasks/calibrateImage.py @@ -26,7 +26,7 @@ import lsst.meas.algorithms import lsst.meas.algorithms.installGaussianPsf import lsst.meas.algorithms.measureApCorr -from lsst.meas.algorithms import sourceSelector +import lsst.meas.base import lsst.meas.astrom import lsst.meas.deblender import lsst.meas.extensions.shapeHSM @@ -144,6 +144,9 @@ class CalibrateImageConfig(pipeBase.PipelineTaskConfig, pipelineConnections=Cali optional=True ) + # To generate catalog ids consistently across subtasks. + id_generator = lsst.meas.base.DetectorVisitIdGeneratorConfig.make_field() + # subtasks used during psf characterization install_simple_psf = pexConfig.ConfigurableField( target=lsst.meas.algorithms.installGaussianPsf.InstallGaussianPsfTask, @@ -183,6 +186,10 @@ class CalibrateImageConfig(pipeBase.PipelineTaskConfig, pipelineConnections=Cali target=lsst.meas.algorithms.SourceDetectionTask, doc="Task to detect stars to return in the output catalog." ) + star_sky_sources = pexConfig.ConfigurableField( + target=lsst.meas.algorithms.SkySourcesTask, + doc="Task to generate sky sources ('empty' regions where there are no detections).", + ) star_mask_streaks = pexConfig.ConfigurableField( target=maskStreaks.MaskStreaksTask, doc="Task for masking streaks. Adds a STREAK mask plane to an exposure.", @@ -308,15 +315,20 @@ def setDefaults(self): self.star_selector["science"].doSignalToNoise = True self.star_selector["science"].doIsolated = True self.star_selector["science"].signalToNoise.minimum = 10.0 + # Keep sky sources in the output catalog, even though they aren't + # wanted for calibration. + self.star_selector["science"].doSkySources = True # Use the affine WCS fitter (assumes we have a good camera geometry). self.astrometry.wcsFitter.retarget(lsst.meas.astrom.FitAffineWcsTask) # phot_g_mean is the primary Gaia band for all input bands. self.astrometry_ref_loader.anyFilterMapsToThis = "phot_g_mean" - # Do not subselect during fitting; we already selected good stars. - self.astrometry.sourceSelector = "null" - self.photometry.match.sourceSelection.retarget(sourceSelector.NullSourceSelectorTask) + # Only reject sky sources; we already selected good stars. + self.astrometry.sourceSelector["science"].doFlags = True + self.astrometry.sourceSelector["science"].flags.bad = ["sky_source"] + self.photometry.match.sourceSelection.doFlags = True + self.photometry.match.sourceSelection.flags.bad = ["sky_source"] # All sources should be good for PSF summary statistics. # TODO: These should both be changed to calib_psf_used with DM-41640. @@ -356,6 +368,7 @@ def __init__(self, initial_stars_schema=None, **kwargs): initial_stars_schema = afwTable.SourceTable.makeMinimalSchema() afwTable.CoordKey.addErrorFields(initial_stars_schema) self.makeSubtask("star_detection", schema=initial_stars_schema) + self.makeSubtask("star_sky_sources", schema=initial_stars_schema) self.makeSubtask("star_mask_streaks") self.makeSubtask("star_deblend", schema=initial_stars_schema) self.makeSubtask("star_measurement", schema=initial_stars_schema) @@ -374,6 +387,7 @@ def __init__(self, initial_stars_schema=None, **kwargs): def runQuantum(self, butlerQC, inputRefs, outputRefs): inputs = butlerQC.get(inputRefs) + id_generator = self.config.id_generator.apply(butlerQC.quantum.dataId) astrometry_loader = lsst.meas.algorithms.ReferenceObjectLoader( dataIds=[ref.datasetRef.dataId for ref in inputRefs.astrometry_ref_cat], @@ -389,12 +403,12 @@ def runQuantum(self, butlerQC, inputRefs, outputRefs): config=self.config.photometry_ref_loader, log=self.log) self.photometry.match.setRefObjLoader(photometry_loader) - outputs = self.run(**inputs) + outputs = self.run(id_generator=id_generator, **inputs) butlerQC.put(outputs, outputRefs) @timeMethod - def run(self, *, exposure): + def run(self, *, exposure, id_generator=None): """Find stars and perform psf measurement, then do a deeper detection and measurement and calibrate astrometry and photometry from that. @@ -403,6 +417,8 @@ def run(self, *, exposure): exposure : `lsst.afw.image.Exposure` Post-ISR exposure, with an initial WCS, VisitInfo, and Filter. Modified in-place during processing. + id_generator : `lsst.meas.base.IdGenerator`, optional + Object that generates source IDs and provides random seeds. Returns ------- @@ -432,11 +448,14 @@ def run(self, *, exposure): Reference catalog stars matches used in the photometric fit. (`list` [`lsst.afw.table.ReferenceMatch`] or `lsst.afw.table.BaseCatalog`) """ + if id_generator is None: + id_generator = lsst.meas.base.IdGenerator() + psf_stars, background, candidates = self._compute_psf(exposure) self._measure_aperture_correction(exposure, psf_stars) - stars = self._find_stars(exposure, background) + stars = self._find_stars(exposure, background, id_generator) astrometry_matches, astrometry_meta = self._fit_astrometry(exposure, stars) stars, photometry_matches, photometry_meta, photo_calib = self._fit_photometry(exposure, stars) @@ -546,7 +565,7 @@ def _measure_aperture_correction(self, exposure, bright_sources): result = self.measure_aperture_correction.run(exposure, bright_sources) exposure.setApCorrMap(result.apCorrMap) - def _find_stars(self, exposure, background): + def _find_stars(self, exposure, background, id_generator): """Detect stars on an exposure that has a PSF model, and measure their PSF, circular aperture, compensated gaussian fluxes. @@ -557,6 +576,8 @@ def _find_stars(self, exposure, background): background : `lsst.afw.math.BackgroundList` Background that was fit to the exposure during detection; modified in-place during subsequent detection. + id_generator : `lsst.meas.base.IdGenerator` + Object that generates source IDs and provides random seeds. Returns ------- @@ -569,6 +590,7 @@ def _find_stars(self, exposure, background): # measurement uses the most accurate background-subtraction. detections = self.star_detection.run(table=table, exposure=exposure, background=background) sources = detections.sources + self.star_sky_sources.run(exposure.mask, sources, seed=id_generator.catalog_id) # Mask streaks self.star_mask_streaks.run(exposure) diff --git a/tests/test_calibrateImage.py b/tests/test_calibrateImage.py index 9fca3c9c3b..35a54e3819 100644 --- a/tests/test_calibrateImage.py +++ b/tests/test_calibrateImage.py @@ -100,6 +100,8 @@ def setUp(self): # We don't have many sources, so have to fit simpler models. self.config.psf_detection.background.approxOrderX = 1 self.config.star_detection.background.approxOrderX = 1 + # Only insert 2 sky sources, for simplicity. + self.config.star_sky_sources.nSources = 2 # Use PCA psf fitter, as psfex fails if there are only 4 stars. self.config.psf_measure_psf.psfDeterminer = 'pca' # We don't have many test points, so can't match on complicated shapes. @@ -111,6 +113,9 @@ def setUp(self): # will use those fluxes here, and hopefully can remove this. self.config.astrometry.magnitudeOutlierRejectionNSigma = 9.0 + # find_stars needs an id generator. + self.id_generator = lsst.meas.base.IdGenerator() + # Something about this test dataset prefers the older fluxRatio here. self.config.star_catalog_calculation.plugins['base_ClassificationExtendedness'].fluxRatio = 0.925 @@ -120,7 +125,7 @@ def test_run(self): calibrate = CalibrateImageTask(config=self.config) calibrate.astrometry.setRefObjLoader(self.ref_loader) calibrate.photometry.match.setRefObjLoader(self.ref_loader) - result = calibrate.run(exposure=self.exposure) + result = calibrate.run(exposure=self.exposure, id_generator=self.id_generator) # Background should have 4 elements: 3 from compute_psf and one from # re-estimation during source detection. @@ -196,14 +201,15 @@ def test_find_stars(self): sources, background, candidates = calibrate._compute_psf(self.exposure) calibrate._measure_aperture_correction(self.exposure, sources) - stars = calibrate._find_stars(self.exposure, background) + stars = calibrate._find_stars(self.exposure, background, self.id_generator) # Background should have 4 elements: 3 from compute_psf and one from # re-estimation during source detection. self.assertEqual(len(background), 4) - # Only psf-like sources with S/N>10 should be in the output catalog. - self.assertEqual(len(stars), 4) + # Only psf-like sources with S/N>10 should be in the output catalog, + # plus two sky sources + self.assertEqual(len(stars), 6) self.assertTrue(sources.isContiguous()) # Sort in order of brightness, to easily compare with expected positions. sources.sort(sources.getPsfFluxSlot().getMeasKey()) @@ -219,12 +225,13 @@ def test_astrometry(self): calibrate.astrometry.setRefObjLoader(self.ref_loader) sources, background, candidates = calibrate._compute_psf(self.exposure) calibrate._measure_aperture_correction(self.exposure, sources) - stars = calibrate._find_stars(self.exposure, background) + stars = calibrate._find_stars(self.exposure, background, self.id_generator) calibrate._fit_astrometry(self.exposure, stars) # Check that we got reliable matches with the truth coordinates. - fitted = SkyCoord(stars['coord_ra'], stars['coord_dec'], unit="radian") + sky = stars["sky_source"] + fitted = SkyCoord(stars[~sky]['coord_ra'], stars[~sky]['coord_dec'], unit="radian") truth = SkyCoord(self.truth_cat['coord_ra'], self.truth_cat['coord_dec'], unit="radian") idx, d2d, _ = fitted.match_to_catalog_sky(truth) np.testing.assert_array_less(d2d.to_value(u.milliarcsecond), 30.0) @@ -238,7 +245,7 @@ def test_photometry(self): calibrate.photometry.match.setRefObjLoader(self.ref_loader) sources, background, candidates = calibrate._compute_psf(self.exposure) calibrate._measure_aperture_correction(self.exposure, sources) - stars = calibrate._find_stars(self.exposure, background) + stars = calibrate._find_stars(self.exposure, background, self.id_generator) calibrate._fit_astrometry(self.exposure, stars) stars, matches, meta, photoCalib = calibrate._fit_photometry(self.exposure, stars) @@ -252,15 +259,21 @@ def test_photometry(self): # PhotoCalib on the exposure must be identically 1. self.assertEqual(self.exposure.photoCalib.getCalibrationMean(), 1.0) - # Check that we got reliable magnitudes and fluxes vs. truth. - fitted = SkyCoord(stars['coord_ra'], stars['coord_dec'], unit="radian") + # Check that we got reliable magnitudes and fluxes vs. truth, ignoring + # sky sources. + sky = stars["sky_source"] + fitted = SkyCoord(stars[~sky]['coord_ra'], stars[~sky]['coord_dec'], unit="radian") truth = SkyCoord(self.truth_cat['coord_ra'], self.truth_cat['coord_dec'], unit="radian") idx, _, _ = fitted.match_to_catalog_sky(truth) # Because the input variance image does not include contributions from # the sources, we can't use fluxErr as a bound on the measurement # quality here. - self.assertFloatsAlmostEqual(stars['slot_PsfFlux_flux'], self.truth_cat['truth_flux'][idx], rtol=0.1) - self.assertFloatsAlmostEqual(stars['slot_PsfFlux_mag'], self.truth_cat['truth_mag'][idx], rtol=0.01) + self.assertFloatsAlmostEqual(stars[~sky]['slot_PsfFlux_flux'], + self.truth_cat['truth_flux'][idx], + rtol=0.1) + self.assertFloatsAlmostEqual(stars[~sky]['slot_PsfFlux_mag'], + self.truth_cat['truth_mag'][idx], + rtol=0.01) class CalibrateImageTaskRunQuantumTests(lsst.utils.tests.TestCase): @@ -277,8 +290,14 @@ def setUp(self): self.repo_path = tempfile.TemporaryDirectory() self.repo = butlerTests.makeTestRepo(self.repo_path.name) + # A complete instrument record is necessary for the id generator. + instrumentRecord = self.repo.dimensions["instrument"].RecordClass( + name=instrument, visit_max=1e6, exposure_max=1e6, detector_max=128, + class_name="lsst.obs.base.instrument_tests.DummyCam", + ) + self.repo.registry.syncDimensionData("instrument", instrumentRecord) + # dataIds for fake data - butlerTests.addDataIdValue(self.repo, "instrument", instrument) butlerTests.addDataIdValue(self.repo, "exposure", exposure) butlerTests.addDataIdValue(self.repo, "visit", visit) butlerTests.addDataIdValue(self.repo, "detector", detector) @@ -354,7 +373,7 @@ def test_runQuantum(self): self.assertEqual(task.astrometry.refObjLoader.name, "gaia_dr3_20230707") self.assertEqual(task.photometry.match.refObjLoader.name, "ps1_pv3_3pi_20170110") # Check that the proper kwargs are passed to run(). - self.assertEqual(mock_run.call_args.kwargs.keys(), {"exposure"}) + self.assertEqual(mock_run.call_args.kwargs.keys(), {"exposure", "id_generator"}) def test_runQuantum_no_optional_outputs(self): config = CalibrateImageTask.ConfigClass() @@ -379,7 +398,7 @@ def test_runQuantum_no_optional_outputs(self): self.assertEqual(task.astrometry.refObjLoader.name, "gaia_dr3_20230707") self.assertEqual(task.photometry.match.refObjLoader.name, "ps1_pv3_3pi_20170110") # Check that the proper kwargs are passed to run(). - self.assertEqual(mock_run.call_args.kwargs.keys(), {"exposure"}) + self.assertEqual(mock_run.call_args.kwargs.keys(), {"exposure", "id_generator"}) def test_lintConnections(self): """Check that the connections are self-consistent.