Coverage for tests/test_calibrateImage.py: 14%

1# This file is part of pipe_tasks. 

2# 

3# Developed for the LSST Data Management System. 

4# This product includes software developed by the LSST Project 

5# (https://www.lsst.org). 

6# See the COPYRIGHT file at the top-level directory of this distribution 

7# for details of code ownership. 

8# 

9# This program is free software: you can redistribute it and/or modify 

10# it under the terms of the GNU General Public License as published by 

11# the Free Software Foundation, either version 3 of the License, or 

12# (at your option) any later version. 

13# 

14# This program is distributed in the hope that it will be useful, 

15# but WITHOUT ANY WARRANTY; without even the implied warranty of 

16# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 

17# GNU General Public License for more details. 

18# 

19# You should have received a copy of the GNU General Public License 

20# along with this program. If not, see <https://www.gnu.org/licenses/>. 

21 

22import unittest 

23from unittest import mock 

24import tempfile 

25 

26import astropy.units as u 

27from astropy.coordinates import SkyCoord 

28import copy 

29import numpy as np 

30import esutil 

31import os 

32import requests 

33 

34import lsst.afw.image as afwImage 

35import lsst.afw.math as afwMath 

36import lsst.afw.table as afwTable 

37import lsst.daf.base 

38import lsst.daf.butler 

39import lsst.daf.butler.tests as butlerTests 

40import lsst.geom 

41import lsst.meas.algorithms 

42from lsst.meas.algorithms import testUtils 

43import lsst.meas.extensions.psfex 

44import lsst.meas.base 

45import lsst.meas.base.tests 

46import lsst.pipe.base.testUtils 

47from lsst.pipe.tasks.calibrateImage import CalibrateImageTask, \ 

48 NoPsfStarsToStarsMatchError, AllCentroidsFlaggedError 

49import lsst.pex.config as pexConfig 

50import lsst.utils.tests 

51 

52 

53class CalibrateImageTaskTests(lsst.utils.tests.TestCase): 

54 

55 def setUp(self): 

56 # Different x/y dimensions so they're easy to distinguish in a plot, 

57 # and non-zero minimum, to help catch xy0 errors. 

58 bbox = lsst.geom.Box2I(lsst.geom.Point2I(5, 4), lsst.geom.Point2I(205, 184)) 

59 self.sky_center = lsst.geom.SpherePoint(245.0, -45.0, lsst.geom.degrees) 

60 self.photo_calib = 12.3 

61 dataset = lsst.meas.base.tests.TestDataset(bbox, crval=self.sky_center, calibration=self.photo_calib) 

62 # sqrt of area of a normalized 2d gaussian 

63 psf_scale = np.sqrt(4*np.pi*(dataset.psfShape.getDeterminantRadius())**2) 

64 noise = 10.0 # stddev of noise per pixel 

65 # Sources ordered from faintest to brightest. 

66 self.fluxes = np.array((6*noise*psf_scale, 

67 12*noise*psf_scale, 

68 45*noise*psf_scale, 

69 150*noise*psf_scale, 

70 400*noise*psf_scale, 

71 1000*noise*psf_scale)) 

72 self.centroids = np.array(((162, 22), 

73 (40, 70), 

74 (100, 160), 

75 (50, 120), 

76 (92, 35), 

77 (175, 154)), dtype=np.float32) 

78 for flux, centroid in zip(self.fluxes, self.centroids): 

79 dataset.addSource(instFlux=flux, centroid=lsst.geom.Point2D(centroid[0], centroid[1])) 

80 

81 # Bright extended source in the center of the image: should not appear 

82 # in any of the output catalogs. 

83 self.extended_source = lsst.geom.Point2D(100, 100) 

84 shape = lsst.afw.geom.Quadrupole(8, 9, 3) 

85 dataset.addSource(instFlux=500*noise*psf_scale, centroid=self.extended_source, shape=shape) 

86 

87 schema = dataset.makeMinimalSchema() 

88 self.truth_exposure, self.truth_cat = dataset.realize(noise=noise, schema=schema) 

89 # Add in a significant background, so we can test that the output 

90 # background is self-consistent with the calibrated exposure. 

91 self.background_level = 500.0 

92 self.truth_exposure.image += self.background_level 

93 # To make it look like a version=1 (nJy fluxes) refcat 

94 self.truth_cat = self.truth_exposure.photoCalib.calibrateCatalog(self.truth_cat) 

95 self.ref_loader = testUtils.MockReferenceObjectLoaderFromMemory([self.truth_cat]) 

96 metadata = lsst.daf.base.PropertyList() 

97 metadata.set("REFCAT_FORMAT_VERSION", 1) 

98 self.truth_cat.setMetadata(metadata) 

99 

100 # TODO: a cosmic ray (need to figure out how to insert a fake-CR) 

101 # self.truth_exposure.image.array[10, 10] = 100000 

102 # self.truth_exposure.variance.array[10, 10] = 100000/noise 

103 

104 # Copy the truth exposure, because CalibrateImage modifies the input. 

105 # Post-ISR ccds only contain: initial WCS, VisitInfo, filter 

106 self.exposure = afwImage.ExposureF(self.truth_exposure, deep=True) 

107 self.exposure.setWcs(self.truth_exposure.wcs) 

108 self.exposure.info.setVisitInfo(self.truth_exposure.visitInfo) 

109 self.exposure.info.id = 12345 

110 # "truth" filter, to match the "truth" refcat. 

111 self.exposure.setFilter(lsst.afw.image.FilterLabel(physical='truth', band="truth")) 

112 self.exposure.metadata["LSST ISR FLAT APPLIED"] = True 

113 

114 # Test-specific configuration: 

115 self.config = CalibrateImageTask.ConfigClass() 

116 # We don't have many sources, so have to fit simpler models. 

117 self.config.psf_detection.background.approxOrderX = 1 

118 self.config.star_detection.background.approxOrderX = 1 

119 # Only insert 2 sky sources, for simplicity. 

120 self.config.star_sky_sources.nSources = 2 

121 # Use PCA psf fitter, as psfex fails if there are only 4 stars. 

122 self.config.psf_measure_psf.psfDeterminer = 'pca' 

123 # We don't have many test points, so can't match on complicated shapes. 

124 self.config.astrometry.sourceSelector["science"].flags.good = [] 

125 self.config.astrometry.matcher.numPointsForShape = 3 

126 self.config.run_sattle = False 

127 # ApFlux has more noise than PsfFlux (the latter unrealistically small 

128 # in this test data), so we need to do magnitude rejection at higher 

129 # sigma, otherwise we can lose otherwise good sources. 

130 # TODO DM-39203: Once we are using Compensated Gaussian Fluxes, we 

131 # will use those fluxes here, and hopefully can remove this. 

132 self.config.astrometry.magnitudeOutlierRejectionNSigma = 9.0 

133 

134 # Make a realistic id generator so that output catalog ids are useful. 

135 # NOTE: The id generator is used to seed the noise replacer during 

136 # measurement, so changes to values here can have subtle effects on 

137 # the centroids and fluxes mesaured on the image, which might cause 

138 # tests to fail. 

139 data_id = lsst.daf.butler.DataCoordinate.standardize( 

140 instrument="I", 

141 visit=self.truth_exposure.visitInfo.id, 

142 detector=12, 

143 universe=lsst.daf.butler.DimensionUniverse(), 

144 ) 

145 self.config.id_generator.packer.name = "observation" 

146 self.config.id_generator.packer["observation"].n_observations = 10000 

147 self.config.id_generator.packer["observation"].n_detectors = 99 

148 self.config.id_generator.n_releases = 8 

149 self.config.id_generator.release_id = 2 

150 self.id_generator = self.config.id_generator.apply(data_id) 

151 

152 # Something about this test dataset prefers a larger threshold here. 

153 self.config.star_selector["science"].unresolved.maximum = 0.2 

154 

155 def _check_run(self, calibrate, result, expect_calibrated_pixels: bool = True): 

156 """Test the result of CalibrateImage.run(). 

157 

158 Parameters 

159 ---------- 

160 calibrate : `lsst.pipe.tasks.calibrateImage.CalibrateImageTask` 

161 Configured task that had `run` called on it. 

162 result : `lsst.pipe.base.Struct` 

163 Result of calling calibrate.run(). 

164 expect_calibrated_pixels : `bool`, optional 

165 Whether to expect image and background pixels to be calibrated. 

166 """ 

167 # Background should have 4 elements: 3 from compute_psf and one from 

168 # re-estimation during source detection. 

169 self.assertEqual(len(result.background), 4) 

170 

171 # Both afw and astropy psf_stars catalogs should be populated. 

172 self.assertEqual(result.psf_stars["calib_psf_used"].sum(), 3) 

173 self.assertEqual(result.psf_stars_footprints["calib_psf_used"].sum(), 3) 

174 

175 # Check that the summary statistics are reasonable. 

176 summary = result.exposure.info.getSummaryStats() 

177 self.assertFloatsAlmostEqual(summary.psfSigma, 2.0, rtol=1e-2) 

178 self.assertFloatsAlmostEqual(summary.ra, self.sky_center.getRa().asDegrees(), rtol=1e-7) 

179 self.assertFloatsAlmostEqual(summary.dec, self.sky_center.getDec().asDegrees(), rtol=1e-7) 

180 

181 # Should have finite sky coordinates in the afw and astropy catalogs. 

182 self.assertTrue(np.isfinite(result.stars_footprints["coord_ra"]).all()) 

183 self.assertTrue(np.isfinite(result.stars["coord_ra"]).all()) 

184 

185 if expect_calibrated_pixels: 

186 # Fit photoCalib should be the applied value if we calibrated 

187 # pixels, not the ==1 one on the exposure. 

188 photo_calib = result.applied_photo_calib 

189 self.assertEqual(result.exposure.photoCalib.getCalibrationMean(), 1.0) 

190 else: 

191 self.assertIsNone(result.applied_photo_calib) 

192 photo_calib = result.exposure.photoCalib 

193 # PhotoCalib comparison is very approximate because we are basing this 

194 # comparison on just 2-3 stars. 

195 self.assertFloatsAlmostEqual(photo_calib.getCalibrationMean(), self.photo_calib, rtol=1e-2) 

196 # Should have calibrated flux/magnitudes in the afw and astropy catalogs 

197 self.assertIn("slot_PsfFlux_flux", result.stars_footprints.schema) 

198 self.assertIn("slot_PsfFlux_mag", result.stars_footprints.schema) 

199 self.assertEqual(result.stars["slot_PsfFlux_flux"].unit, u.nJy) 

200 self.assertEqual(result.stars["slot_PsfFlux_mag"].unit, u.ABmag) 

201 

202 # Should have detected all S/N >= 10 sources plus 2 sky sources, whether 1 or 2 snaps. 

203 self.assertEqual(len(result.stars), 6) 

204 # Did the psf flags get propagated from the psf_stars catalog? 

205 self.assertEqual(result.stars["calib_psf_used"].sum(), 3) 

206 

207 # Check that all necessary fields are in the output. 

208 lsst.pipe.base.testUtils.assertValidOutput(calibrate, result) 

209 

210 # Check metadata. 

211 key = "LSST CALIB ILLUMCORR APPLIED" 

212 self.assertIn(key, result.exposure.metadata) 

213 self.assertEqual(result.exposure.metadata[key], False) 

214 

215 # Check that the psf_stars cross match worked correctly. 

216 matches = esutil.numpy_util.match(result.psf_stars["id"], result.stars["psf_id"]) 

217 self.assertFloatsAlmostEqual(result.psf_stars["slot_Centroid_x"][matches[0]], 

218 result.stars["slot_Centroid_x"][matches[1]], atol=2e-5) 

219 if "astrometry_matches" in self.config.optional_outputs: 

220 matches = esutil.numpy_util.match(result.astrometry_matches["src_id"], 

221 result.photometry_matches["src_psf_id"]) 

222 self.assertFloatsAlmostEqual(result.astrometry_matches["src_slot_Centroid_x"][matches[0]], 

223 result.photometry_matches["src_slot_Centroid_x"][matches[1]], 

224 atol=2e-5) 

225 

226 def test_run(self): 

227 """Test that run() returns reasonable values to be butler put. 

228 """ 

229 calibrate = CalibrateImageTask(config=self.config) 

230 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

231 calibrate.photometry.match.setRefObjLoader(self.ref_loader) 

232 result = calibrate.run(exposures=self.exposure) 

233 

234 self._check_run(calibrate, result) 

235 

236 def test_run_downsample(self): 

237 """Test that run() runs with downsample. 

238 """ 

239 config = copy.copy(self.config) 

240 config.do_downsample_footprints = True 

241 config.downsample_max_footprints = 5 

242 

243 calibrate = CalibrateImageTask(config=config) 

244 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

245 calibrate.photometry.match.setRefObjLoader(self.ref_loader) 

246 with self.assertLogs("lsst.calibrateImage", level="INFO") as cm: 

247 result = calibrate.run(exposures=self.exposure) 

248 self.assertIn( 

249 "INFO:lsst.calibrateImage:Downsampling from 8 to 7 non-sky-source footprints.", 

250 cm.output, 

251 ) 

252 

253 self._check_run(calibrate, result) 

254 

255 def test_run_2_snaps(self): 

256 """Test that run() returns reasonable values to be butler put, when 

257 passed two exposures to combine as snaps. 

258 """ 

259 calibrate = CalibrateImageTask(config=self.config) 

260 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

261 calibrate.photometry.match.setRefObjLoader(self.ref_loader) 

262 # Halve the flux in each exposure to get the expected visit sum. 

263 self.exposure.image /= 2 

264 self.exposure.variance /= 2 

265 result = calibrate.run(exposures=[self.exposure, self.exposure]) 

266 

267 self._check_run(calibrate, result) 

268 

269 def test_run_no_optionals(self): 

270 """Test that disabling optional outputs removes them from the output 

271 struct, as appropriate. 

272 """ 

273 self.config.optional_outputs = [] 

274 calibrate = CalibrateImageTask(config=self.config) 

275 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

276 calibrate.photometry.match.setRefObjLoader(self.ref_loader) 

277 result = calibrate.run(exposures=self.exposure) 

278 

279 self._check_run(calibrate, result) 

280 # These are the only optional outputs that require extra computation, 

281 # the others are included in the output struct regardless. 

282 self.assertNotIn("astrometry_matches", result.getDict()) 

283 self.assertNotIn("photometry_matches", result.getDict()) 

284 

285 def test_run_no_calibrate_pixels(self): 

286 """Test that run() returns reasonable values when 

287 do_calibrate_pixels=False. 

288 """ 

289 self.config.do_calibrate_pixels = False 

290 calibrate = CalibrateImageTask(config=self.config) 

291 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

292 calibrate.photometry.match.setRefObjLoader(self.ref_loader) 

293 result = calibrate.run(exposures=self.exposure) 

294 

295 self._check_run(calibrate, result, expect_calibrated_pixels=False) 

296 

297 def test_run_no_astrom_errors(self): 

298 """Test that run() returns reasonable values when 

299 do_calibrate_pixels=False. 

300 """ 

301 self.config.do_include_astrometric_errors = False 

302 calibrate = CalibrateImageTask(config=self.config) 

303 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

304 calibrate.photometry.match.setRefObjLoader(self.ref_loader) 

305 result = calibrate.run(exposures=self.exposure) 

306 

307 self._check_run(calibrate, result) 

308 

309 def test_compute_psf(self): 

310 """Test that our brightest sources are found by _compute_psf(), 

311 that a PSF is assigned to the expopsure. 

312 """ 

313 calibrate = CalibrateImageTask(config=self.config) 

314 psf_stars, background, candidates = calibrate._compute_psf(self.exposure, self.id_generator) 

315 

316 # Catalog ids should be very large from this id generator. 

317 self.assertTrue(all(psf_stars['id'] > 1000000000)) 

318 

319 # Background should have 3 elements: initial subtraction, and two from 

320 # re-estimation during the two detection passes. 

321 self.assertEqual(len(background), 3) 

322 

323 # Only the point-sources with S/N > 50 should be in this output. 

324 self.assertEqual(psf_stars["calib_psf_used"].sum(), 3) 

325 # Sort in order of brightness, to easily compare with expected positions. 

326 psf_stars.sort(psf_stars.getPsfFluxSlot().getMeasKey()) 

327 for record, flux, center in zip(psf_stars[::-1], self.fluxes, self.centroids[self.fluxes > 50]): 

328 self.assertFloatsAlmostEqual(record.getX(), center[0], rtol=0.01) 

329 self.assertFloatsAlmostEqual(record.getY(), center[1], rtol=0.01) 

330 # PsfFlux should match the values inserted. 

331 self.assertFloatsAlmostEqual(record["slot_PsfFlux_instFlux"], flux, rtol=0.01) 

332 

333 # TODO: While debugging DM-32701, we're using PCA instead of psfex. 

334 # Check that we got a useable PSF. 

335 # self.assertIsInstance(self.exposure.psf, lsst.meas.extensions.psfex.PsfexPsf) 

336 self.assertIsInstance(self.exposure.psf, lsst.meas.algorithms.PcaPsf) 

337 # TestDataset sources have PSF radius=2 pixels. 

338 radius = self.exposure.psf.computeShape(self.exposure.psf.getAveragePosition()).getDeterminantRadius() 

339 self.assertFloatsAlmostEqual(radius, 2.0, rtol=1e-2) 

340 

341 # To look at images for debugging (`setup display_ds9` and run ds9): 

342 # import lsst.afw.display 

343 # display = lsst.afw.display.getDisplay() 

344 # display.mtv(self.exposure) 

345 

346 def test_compute_psf_bad_centroids(self): 

347 """Test that we raise an appropriate error if all measured centroids 

348 are bad in compute_psf. 

349 The root cause of this is likely an interaction between the PSF fit 

350 and cosmic ray repair. An ugly true PSF can result in the centers of 

351 sources being masked and repaired as CRs, thus removing the center of 

352 the sources and causing them all to be flagged. 

353 """ 

354 # make the sources have a (not too deep) hole in the middle 

355 for point in self.centroids: 

356 for i in (-1, 0, 1): 

357 for j in (-1, 0, 1): 

358 self.exposure.image[point[0]+i, point[1]+j] -= \ 

359 (self.exposure.image[point[0]+i, point[1]+j] - self.background_level)/1.1 

360 # and the extended central source 

361 point = self.extended_source 

362 self.exposure.image[point[0], point[1]] -= \ 

363 (self.exposure.image[point[0], point[1]] - self.background_level) 

364 

365 # add a diagonal gradient to each source 

366 size = 5 

367 for point, flux in zip(self.centroids, self.fluxes): 

368 box = lsst.geom.Box2I.makeCenteredBox(lsst.geom.Point2D(point[0], point[1]), 

369 lsst.geom.Extent2I(size*2, size*2)) 

370 X, Y = np.ogrid[point[0] - size:point[0] + size, point[1] - size:point[1] + size] 

371 distance = np.sqrt((X - point[0])**2 + (Y - point[1])**2) 

372 gradient = (X - point[0] + size)*10 + (Y - point[1] + size)*10 

373 gradient[distance > size] = 0 

374 cutout = self.exposure.image.subset(box) 

375 cutout.array += gradient / flux 

376 

377 calibrate = CalibrateImageTask(config=self.config) 

378 with self.assertRaisesRegex(AllCentroidsFlaggedError, r"source centroids \(out of 4\) flagged"): 

379 calibrate.run(exposures=[self.exposure], id_generator=self.id_generator) 

380 

381 def test_measure_aperture_correction(self): 

382 """Test that _measure_aperture_correction() assigns an ApCorrMap to the 

383 exposure. 

384 """ 

385 calibrate = CalibrateImageTask(config=self.config) 

386 psf_stars, background, candidates = calibrate._compute_psf(self.exposure, self.id_generator) 

387 

388 # First check that the exposure doesn't have an ApCorrMap. 

389 self.assertIsNone(self.exposure.apCorrMap) 

390 calibrate._measure_aperture_correction(self.exposure, psf_stars) 

391 self.assertIsInstance(self.exposure.apCorrMap, afwImage.ApCorrMap) 

392 # We know that there are 2 fields from the normalization, plus more from 

393 # other configured plugins. 

394 self.assertGreater(len(self.exposure.apCorrMap), 2) 

395 

396 def test_find_stars(self): 

397 """Test that _find_stars() correctly identifies the S/N>10 stars 

398 in the image and returns them in the output catalog. 

399 """ 

400 calibrate = CalibrateImageTask(config=self.config) 

401 psf_stars, background, candidates = calibrate._compute_psf(self.exposure, self.id_generator) 

402 calibrate._measure_aperture_correction(self.exposure, psf_stars) 

403 

404 stars = calibrate._find_stars(self.exposure, background, self.id_generator) 

405 

406 # Catalog ids should be very large from this id generator. 

407 self.assertTrue(all(stars['id'] > 1000000000)) 

408 

409 # Background should have 4 elements: 3 from compute_psf and one from 

410 # re-estimation during source detection. 

411 self.assertEqual(len(background), 4) 

412 

413 # Only 5 psf-like sources with S/N>10 should be in the output catalog, 

414 # plus two sky sources. 

415 self.assertEqual(len(stars), 6) 

416 self.assertTrue(stars.isContiguous()) 

417 # Sort in order of brightness, to easily compare with expected positions. 

418 stars.sort(stars.getPsfFluxSlot().getMeasKey()) 

419 for record, flux, center in zip(stars[::-1], self.fluxes, self.centroids[self.fluxes > 50]): 

420 self.assertFloatsAlmostEqual(record.getX(), center[0], rtol=0.01) 

421 self.assertFloatsAlmostEqual(record.getY(), center[1], rtol=0.01) 

422 self.assertFloatsAlmostEqual(record["slot_PsfFlux_instFlux"], flux, rtol=0.01) 

423 

424 def test_astrometry(self): 

425 """Test that the fitted WCS gives good catalog coordinates. 

426 """ 

427 calibrate = CalibrateImageTask(config=self.config) 

428 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

429 psf_stars, background, candidates = calibrate._compute_psf(self.exposure, self.id_generator) 

430 calibrate._measure_aperture_correction(self.exposure, psf_stars) 

431 stars = calibrate._find_stars(self.exposure, background, self.id_generator) 

432 

433 calibrate._fit_astrometry(self.exposure, stars) 

434 

435 # Check that we got reliable matches with the truth coordinates. 

436 sky = stars["sky_source"] 

437 fitted = SkyCoord(stars[~sky]['coord_ra'], stars[~sky]['coord_dec'], unit="radian") 

438 truth = SkyCoord(self.truth_cat['coord_ra'], self.truth_cat['coord_dec'], unit="radian") 

439 idx, d2d, _ = fitted.match_to_catalog_sky(truth) 

440 np.testing.assert_array_less(d2d.to_value(u.milliarcsecond), 35.0) 

441 

442 def test_photometry(self): 

443 """Test that the fitted photoCalib matches the one we generated, 

444 and that the exposure is calibrated. 

445 """ 

446 calibrate = CalibrateImageTask(config=self.config) 

447 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

448 calibrate.photometry.match.setRefObjLoader(self.ref_loader) 

449 psf_stars, background, candidates = calibrate._compute_psf(self.exposure, self.id_generator) 

450 calibrate._measure_aperture_correction(self.exposure, psf_stars) 

451 stars = calibrate._find_stars(self.exposure, background, self.id_generator) 

452 calibrate._fit_astrometry(self.exposure, stars) 

453 

454 stars, matches, meta, photoCalib = calibrate._fit_photometry(self.exposure, stars) 

455 calibrate._apply_photometry(self.exposure, background) 

456 

457 # NOTE: With this test data, PhotoCalTask returns calibrationErr==0, 

458 # so we can't check that the photoCal error has been set. 

459 self.assertFloatsAlmostEqual(photoCalib.getCalibrationMean(), self.photo_calib, rtol=1e-2) 

460 # The exposure should be calibrated by the applied photoCalib, 

461 # and the background should be calibrated to match. 

462 uncalibrated = self.exposure.image.clone() 

463 uncalibrated += background.getImage() 

464 uncalibrated /= self.photo_calib 

465 self.assertFloatsAlmostEqual(uncalibrated.array, self.truth_exposure.image.array, rtol=1e-2) 

466 # PhotoCalib on the exposure must be identically 1. 

467 self.assertEqual(self.exposure.photoCalib.getCalibrationMean(), 1.0) 

468 

469 # Check that we got reliable magnitudes and fluxes vs. truth, ignoring 

470 # sky sources. 

471 sky = stars["sky_source"] 

472 fitted = SkyCoord(stars[~sky]['coord_ra'], stars[~sky]['coord_dec'], unit="radian") 

473 truth = SkyCoord(self.truth_cat['coord_ra'], self.truth_cat['coord_dec'], unit="radian") 

474 idx, _, _ = fitted.match_to_catalog_sky(truth) 

475 # Because the input variance image does not include contributions from 

476 # the sources, we can't use fluxErr as a bound on the measurement 

477 # quality here. 

478 self.assertFloatsAlmostEqual(stars[~sky]['slot_PsfFlux_flux'], 

479 self.truth_cat['truth_flux'][idx], 

480 rtol=0.1) 

481 self.assertFloatsAlmostEqual(stars[~sky]['slot_PsfFlux_mag'], 

482 self.truth_cat['truth_mag'][idx], 

483 rtol=0.01) 

484 

485 def test_match_psf_stars(self): 

486 """Test that _match_psf_stars() flags the correct stars as psf stars 

487 and candidates. 

488 """ 

489 calibrate = CalibrateImageTask(config=self.config) 

490 psf_stars, background, candidates = calibrate._compute_psf(self.exposure, self.id_generator) 

491 calibrate._measure_aperture_correction(self.exposure, psf_stars) 

492 stars = calibrate._find_stars(self.exposure, background, self.id_generator) 

493 

494 # There should be no psf-related flags set at first. 

495 self.assertEqual(stars["calib_psf_candidate"].sum(), 0) 

496 self.assertEqual(stars["calib_psf_used"].sum(), 0) 

497 self.assertEqual(stars["calib_psf_reserved"].sum(), 0) 

498 

499 # Reorder stars to be out of order with psf_stars (putting the sky 

500 # sources in front); this tests that I get the indexing right. 

501 stars.sort(stars.getCentroidSlot().getMeasKey().getX()) 

502 stars = stars.copy(deep=True) 

503 # Re-number the ids: the matcher requires sorted ids: this is always 

504 # true in the code itself, but we've permuted them by sorting on 

505 # flux. We don't care what the actual ids themselves are here. 

506 stars["id"] = np.arange(len(stars)) 

507 

508 calibrate._match_psf_stars(psf_stars, stars) 

509 

510 # Check that the three brightest stars have the psf flags transfered 

511 # from the psf_stars catalog by sorting in order of brightness. 

512 stars.sort(stars.getPsfFluxSlot().getMeasKey()) 

513 # sort() above leaves the catalog non-contiguous. 

514 stars = stars.copy(deep=True) 

515 np.testing.assert_array_equal(stars["calib_psf_candidate"], 

516 [False, False, False, True, True, True]) 

517 np.testing.assert_array_equal(stars["calib_psf_used"], 

518 [False, False, False, True, True, True]) 

519 # Too few sources to reserve any in these tests. 

520 self.assertEqual(stars["calib_psf_reserved"].sum(), 0) 

521 

522 def test_match_psf_stars_no_matches(self): 

523 """Check that _match_psf_stars handles the case of no cross-matches. 

524 """ 

525 calibrate = CalibrateImageTask(config=self.config) 

526 # Make two catalogs that cannot have matches. 

527 stars = self.truth_cat[2:].copy(deep=True) 

528 psf_stars = self.truth_cat[:2].copy(deep=True) 

529 

530 with self.assertRaisesRegex(NoPsfStarsToStarsMatchError, 

531 "No psf stars out of 2 matched 5 calib stars") as cm: 

532 calibrate._match_psf_stars(psf_stars, stars) 

533 self.assertEqual(cm.exception.metadata["n_psf_stars"], 2) 

534 self.assertEqual(cm.exception.metadata["n_stars"], 5) 

535 

536 def test_calibrate_image_illumcorr(self): 

537 """Test running through with an illumination correction.""" 

538 config = copy.copy(self.config) 

539 config.do_illumination_correction = True 

540 config.psf_subtract_background.doApplyFlatBackgroundRatio = True 

541 config.psf_detection.doApplyFlatBackgroundRatio = True 

542 config.star_detection.doApplyFlatBackgroundRatio = True 

543 

544 calibrate = CalibrateImageTask(config=config) 

545 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

546 calibrate.photometry.match.setRefObjLoader(self.ref_loader) 

547 

548 # Assume that the exposure has been flattened by a flat-flat. 

549 background_flat = self.exposure.clone() 

550 background_flat.image.array[:, :] = 1.0 

551 background_flat.mask.array[:, :] = 0 

552 background_flat.variance.array[:, :] = 0.0 

553 

554 # And create an illumination correction of 1.1. 

555 illum_corr_value = 1.1 

556 illumination_correction = self.exposure.clone() 

557 illumination_correction.image.array[:, :] = illum_corr_value 

558 illumination_correction.mask.array[:, :] = 0 

559 illumination_correction.variance.array[:, :] = 0.0 

560 

561 result = calibrate.run( 

562 exposures=self.exposure, 

563 id_generator=self.id_generator, 

564 background_flat=background_flat, 

565 illumination_correction=illumination_correction, 

566 ) 

567 

568 # We divide the image by the illumination correction, but the reference 

569 # sources stay the same, so the applied photocalib will increase by 

570 # the illum_corr_value. 

571 # Tolerance is the same as the direct test with no illumination 

572 # correction. 

573 self.assertFloatsAlmostEqual( 

574 result.applied_photo_calib.getCalibrationMean(), 

575 self.photo_calib * illum_corr_value, 

576 rtol=1e-2, 

577 ) 

578 

579 self.assertEqual(len(result.background), 4) 

580 self.assertFloatsAlmostEqual( 

581 np.median(result.background.getImage().array), 

582 result.applied_photo_calib.getCalibrationMean() * self.background_level, 

583 rtol=1e-3, 

584 ) 

585 

586 # Check metadata. 

587 key = "LSST CALIB ILLUMCORR APPLIED" 

588 self.assertIn(key, result.exposure.metadata) 

589 self.assertEqual(result.exposure.metadata[key], True) 

590 

591 @mock.patch.dict(os.environ, {"SATTLE_URI_BASE": ""}) 

592 def test_fail_on_sattle_miconfiguration(self): 

593 """Test for failure if sattle is requested without appropriate configurations. 

594 """ 

595 self.config.run_sattle = True 

596 with self.assertRaises(pexConfig.FieldValidationError): 

597 CalibrateImageTask(config=self.config) 

598 

599 @mock.patch.dict(os.environ, {"SATTLE_URI_BASE": "fake_host:1234"}) 

600 def test_continue_on_sattle_failure(self): 

601 """Processing should continue when sattle returns status codes other than 200. 

602 """ 

603 response = MockResponse({}, 500, "internal sattle error") 

604 

605 self.config.run_sattle = True 

606 calibrate = CalibrateImageTask(config=self.config) 

607 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

608 calibrate.photometry.match.setRefObjLoader(self.ref_loader) 

609 with mock.patch('requests.put', return_value=response) as mock_put: 

610 calibrate.run(exposures=self.exposure) 

611 mock_put.assert_called_once() 

612 

613 @mock.patch.dict(os.environ, {"SATTLE_URI_BASE": "fake_host:1234"}) 

614 def test_sattle(self): 

615 """Test for successful completion when sattle call returns successfully. 

616 """ 

617 response = MockResponse({}, 200, "success") 

618 

619 self.config.run_sattle = True 

620 calibrate = CalibrateImageTask(config=self.config) 

621 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

622 calibrate.photometry.match.setRefObjLoader(self.ref_loader) 

623 with mock.patch('requests.put', return_value=response) as mock_put: 

624 calibrate.run(exposures=self.exposure) 

625 mock_put.assert_called_once() 

626 

627 

628class CalibrateImageTaskRunQuantumTests(lsst.utils.tests.TestCase): 

629 """Tests of ``CalibrateImageTask.runQuantum``, which need a test butler, 

630 but do not need real images. 

631 """ 

632 def setUp(self): 

633 instrument = "testCam" 

634 exposure0 = 101 

635 exposure1 = 102 

636 visit = 100101 

637 detector = 42 

638 physical_filter = "r" 

639 

640 # Create a and populate a test butler for runQuantum tests. 

641 self.repo_path = tempfile.TemporaryDirectory(ignore_cleanup_errors=True) 

642 self.repo = butlerTests.makeTestRepo(self.repo_path.name) 

643 

644 # A complete instrument record is necessary for the id generator. 

645 instrumentRecord = self.repo.dimensions["instrument"].RecordClass( 

646 name=instrument, visit_max=1e6, exposure_max=1e6, detector_max=128, 

647 class_name="lsst.obs.base.instrument_tests.DummyCam", 

648 ) 

649 self.repo.registry.syncDimensionData("instrument", instrumentRecord) 

650 

651 # dataIds for fake data 

652 butlerTests.addDataIdValue(self.repo, "detector", detector) 

653 butlerTests.addDataIdValue(self.repo, "exposure", exposure0) 

654 butlerTests.addDataIdValue(self.repo, "exposure", exposure1) 

655 butlerTests.addDataIdValue(self.repo, "visit", visit) 

656 butlerTests.addDataIdValue(self.repo, "physical_filter", physical_filter) 

657 

658 # inputs 

659 butlerTests.addDatasetType(self.repo, "postISRCCD", {"instrument", "exposure", "detector"}, 

660 "ExposureF") 

661 butlerTests.addDatasetType(self.repo, "gaia_dr3_20230707", {"htm7"}, "SimpleCatalog") 

662 butlerTests.addDatasetType(self.repo, "ps1_pv3_3pi_20170110", {"htm7"}, "SimpleCatalog") 

663 butlerTests.addDatasetType(self.repo, "flat", {"instrument", "detector", "physical_filter"}, 

664 "Exposure") 

665 butlerTests.addDatasetType(self.repo, 

666 "illuminationCorrection", 

667 {"instrument", "detector", "physical_filter"}, 

668 "Exposure") 

669 

670 # outputs 

671 butlerTests.addDatasetType(self.repo, "initial_pvi", {"instrument", "visit", "detector"}, 

672 "ExposureF") 

673 butlerTests.addDatasetType(self.repo, "initial_stars_footprints_detector", 

674 {"instrument", "visit", "detector"}, 

675 "SourceCatalog") 

676 butlerTests.addDatasetType(self.repo, "initial_stars_detector", 

677 {"instrument", "visit", "detector"}, 

678 "ArrowAstropy") 

679 butlerTests.addDatasetType(self.repo, "initial_photoCalib_detector", 

680 {"instrument", "visit", "detector"}, 

681 "PhotoCalib")