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 

52from utils import makeTestVisitInfo 

53 

54 

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

56 

57 def setUp(self): 

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

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

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

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

62 self.photo_calib = 12.3 

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

64 # sqrt of area of a normalized 2d gaussian 

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

66 noise = 10.0 # stddev of noise per pixel 

67 # Sources ordered from faintest to brightest. 

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

69 12*noise*psf_scale, 

70 45*noise*psf_scale, 

71 150*noise*psf_scale, 

72 400*noise*psf_scale, 

73 1000*noise*psf_scale)) 

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

75 (40, 70), 

76 (100, 160), 

77 (50, 120), 

78 (92, 35), 

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

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

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

82 

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

84 # in any of the output catalogs. 

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

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

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

88 

89 schema = dataset.makeMinimalSchema() 

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

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

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

93 self.background_level = 500.0 

94 self.truth_exposure.image += self.background_level 

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

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

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

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

99 metadata.set("REFCAT_FORMAT_VERSION", 1) 

100 self.truth_cat.setMetadata(metadata) 

101 

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

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

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

105 

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

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

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

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

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

111 self.exposure.info.id = 12345 

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

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

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

115 

116 # Test-specific configuration: 

117 self.config = CalibrateImageTask.ConfigClass() 

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

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

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

121 # Only insert 2 sky sources, for simplicity. 

122 self.config.star_sky_sources.nSources = 2 

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

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

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

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

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

128 self.config.run_sattle = False 

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

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

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

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

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

134 self.config.astrometry.magnitudeOutlierRejectionNSigma = 9.0 

135 

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

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

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

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

140 # tests to fail. 

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

142 instrument="I", 

143 visit=self.truth_exposure.visitInfo.id, 

144 detector=12, 

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

146 ) 

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

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

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

150 self.config.id_generator.n_releases = 8 

151 self.config.id_generator.release_id = 2 

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

153 

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

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

156 

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

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

159 

160 Parameters 

161 ---------- 

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

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

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

165 Result of calling calibrate.run(). 

166 expect_calibrated_pixels : `bool`, optional 

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

168 """ 

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

170 # re-estimation during source detection. 

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

172 

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

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

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

176 

177 # Check that the summary statistics are reasonable. 

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

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

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

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

182 

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

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

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

186 

187 if expect_calibrated_pixels: 

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

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

190 photo_calib = result.applied_photo_calib 

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

192 else: 

193 self.assertIsNone(result.applied_photo_calib) 

194 photo_calib = result.exposure.photoCalib 

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

196 # comparison on just 2-3 stars. 

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

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

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

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

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

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

203 

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

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

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

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

208 

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

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

211 

212 # Check metadata. 

213 key = "LSST CALIB ILLUMCORR APPLIED" 

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

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

216 

217 # Check that the psf_stars cross match worked correctly. 

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

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

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

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

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

223 result.photometry_matches["src_psf_id"]) 

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

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

226 atol=2e-5) 

227 

228 def test_run(self): 

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

230 """ 

231 calibrate = CalibrateImageTask(config=self.config) 

232 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

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

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

235 

236 self._check_run(calibrate, result) 

237 

238 def test_run_downsample(self): 

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

240 """ 

241 config = copy.copy(self.config) 

242 config.do_downsample_footprints = True 

243 config.downsample_max_footprints = 5 

244 

245 calibrate = CalibrateImageTask(config=config) 

246 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

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

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

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

250 self.assertIn( 

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

252 cm.output, 

253 ) 

254 

255 self._check_run(calibrate, result) 

256 

257 def test_run_2_snaps(self): 

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

259 passed two exposures to combine as snaps. 

260 """ 

261 calibrate = CalibrateImageTask(config=self.config) 

262 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

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

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

265 self.exposure.image /= 2 

266 self.exposure.variance /= 2 

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

268 

269 self._check_run(calibrate, result) 

270 

271 def test_run_no_optionals(self): 

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

273 struct, as appropriate. 

274 """ 

275 self.config.optional_outputs = [] 

276 calibrate = CalibrateImageTask(config=self.config) 

277 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

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

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

280 

281 self._check_run(calibrate, result) 

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

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

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

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

286 

287 def test_run_no_calibrate_pixels(self): 

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

289 do_calibrate_pixels=False. 

290 """ 

291 self.config.do_calibrate_pixels = False 

292 calibrate = CalibrateImageTask(config=self.config) 

293 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

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

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

296 

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

298 

299 def test_run_no_astrom_errors(self): 

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

301 do_calibrate_pixels=False. 

302 """ 

303 self.config.do_include_astrometric_errors = False 

304 calibrate = CalibrateImageTask(config=self.config) 

305 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

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

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

308 

309 self._check_run(calibrate, result) 

310 

311 def test_compute_psf(self): 

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

313 that a PSF is assigned to the expopsure. 

314 """ 

315 calibrate = CalibrateImageTask(config=self.config) 

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

317 

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

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

320 

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

322 # re-estimation during the two detection passes. 

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

324 

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

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

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

328 psf_stars.sort(psf_stars.getPsfFluxSlot().getMeasKey()) 

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

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

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

332 # PsfFlux should match the values inserted. 

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

334 

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

336 # Check that we got a useable PSF. 

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

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

339 # TestDataset sources have PSF radius=2 pixels. 

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

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

342 

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

344 # import lsst.afw.display 

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

346 # display.mtv(self.exposure) 

347 

348 def test_compute_psf_bad_centroids(self): 

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

350 are bad in compute_psf. 

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

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

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

354 the sources and causing them all to be flagged. 

355 """ 

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

357 for point in self.centroids: 

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

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

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

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

362 # and the extended central source 

363 point = self.extended_source 

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

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

366 

367 # add a diagonal gradient to each source 

368 size = 5 

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

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

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

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

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

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

375 gradient[distance > size] = 0 

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

377 cutout.array += gradient / flux 

378 

379 calibrate = CalibrateImageTask(config=self.config) 

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

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

382 

383 def test_measure_aperture_correction(self): 

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

385 exposure. 

386 """ 

387 calibrate = CalibrateImageTask(config=self.config) 

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

389 

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

391 self.assertIsNone(self.exposure.apCorrMap) 

392 calibrate._measure_aperture_correction(self.exposure, psf_stars) 

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

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

395 # other configured plugins. 

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

397 

398 def test_find_stars(self): 

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

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

401 """ 

402 calibrate = CalibrateImageTask(config=self.config) 

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

404 calibrate._measure_aperture_correction(self.exposure, psf_stars) 

405 

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

407 

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

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

410 

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

412 # re-estimation during source detection. 

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

414 

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

416 # plus two sky sources. 

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

418 self.assertTrue(stars.isContiguous()) 

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

420 stars.sort(stars.getPsfFluxSlot().getMeasKey()) 

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

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

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

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

425 

426 def test_astrometry(self): 

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

428 """ 

429 calibrate = CalibrateImageTask(config=self.config) 

430 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

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

432 calibrate._measure_aperture_correction(self.exposure, psf_stars) 

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

434 

435 calibrate._fit_astrometry(self.exposure, stars) 

436 

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

438 sky = stars["sky_source"] 

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

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

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

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

443 

444 def test_photometry(self): 

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

446 and that the exposure is calibrated. 

447 """ 

448 calibrate = CalibrateImageTask(config=self.config) 

449 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

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

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

452 calibrate._measure_aperture_correction(self.exposure, psf_stars) 

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

454 calibrate._fit_astrometry(self.exposure, stars) 

455 

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

457 calibrate._apply_photometry(self.exposure, background) 

458 

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

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

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

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

463 # and the background should be calibrated to match. 

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

465 uncalibrated += background.getImage() 

466 uncalibrated /= self.photo_calib 

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

468 # PhotoCalib on the exposure must be identically 1. 

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

470 

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

472 # sky sources. 

473 sky = stars["sky_source"] 

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

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

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

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

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

479 # quality here. 

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

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

482 rtol=0.1) 

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

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

485 rtol=0.01) 

486 

487 def test_match_psf_stars(self): 

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

489 and candidates. 

490 """ 

491 calibrate = CalibrateImageTask(config=self.config) 

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

493 calibrate._measure_aperture_correction(self.exposure, psf_stars) 

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

495 

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

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

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

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

500 

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

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

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

504 stars = stars.copy(deep=True) 

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

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

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

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

509 

510 calibrate._match_psf_stars(psf_stars, stars) 

511 

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

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

514 stars.sort(stars.getPsfFluxSlot().getMeasKey()) 

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

516 stars = stars.copy(deep=True) 

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

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

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

520 [False, False, False, True, True, True]) 

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

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

523 

524 def test_match_psf_stars_no_matches(self): 

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

526 """ 

527 calibrate = CalibrateImageTask(config=self.config) 

528 # Make two catalogs that cannot have matches. 

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

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

531 

532 with self.assertRaisesRegex(NoPsfStarsToStarsMatchError, 

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

534 calibrate._match_psf_stars(psf_stars, stars) 

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

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

537 

538 def test_calibrate_image_illumcorr(self): 

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

540 config = copy.copy(self.config) 

541 config.do_illumination_correction = True 

542 config.psf_subtract_background.doApplyFlatBackgroundRatio = True 

543 config.psf_detection.doApplyFlatBackgroundRatio = True 

544 config.star_detection.doApplyFlatBackgroundRatio = True 

545 

546 calibrate = CalibrateImageTask(config=config) 

547 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

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

549 

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

551 background_flat = self.exposure.clone() 

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

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

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

555 

556 # And create an illumination correction of 1.1. 

557 illum_corr_value = 1.1 

558 illumination_correction = self.exposure.clone() 

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

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

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

562 

563 result = calibrate.run( 

564 exposures=self.exposure, 

565 id_generator=self.id_generator, 

566 background_flat=background_flat, 

567 illumination_correction=illumination_correction, 

568 ) 

569 

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

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

572 # the illum_corr_value. 

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

574 # correction. 

575 self.assertFloatsAlmostEqual( 

576 result.applied_photo_calib.getCalibrationMean(), 

577 self.photo_calib * illum_corr_value, 

578 rtol=1e-2, 

579 ) 

580 

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

582 self.assertFloatsAlmostEqual( 

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

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

585 rtol=1e-3, 

586 ) 

587 

588 # Check metadata. 

589 key = "LSST CALIB ILLUMCORR APPLIED" 

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

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

592 

593 def test_run_with_diffraction_spike_mask(self): 

594 """Test that the diffraction spike mask subtask runs. 

595 """ 

596 config = self.config 

597 config.doMaskDiffractionSpikes = True 

598 config.diffractionSpikeMask.magnitudeThreshold = 19 

599 # Define a fake SATURATED mask plane for use by the diffraction spike 

600 # task, so that it does not affect the rest of calibrate 

601 config.diffractionSpikeMask.saturatedMaskPlane = "FAKESATURATED" 

602 calibrate = CalibrateImageTask(config=config) 

603 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

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

605 calibrate.diffractionSpikeMask.setRefObjLoader(self.ref_loader) 

606 

607 exposure = self.exposure.clone() 

608 

609 exposure.info.setVisitInfo(makeTestVisitInfo()) 

610 exposure.mask.addMaskPlane(config.diffractionSpikeMask.saturatedMaskPlane) 

611 

612 # Set the saturated mask plane in half of the image 

613 saturatedMaskBit = exposure.mask.getPlaneBitMask(config.diffractionSpikeMask.saturatedMaskPlane) 

614 bbox = exposure.getBBox() 

615 bbox.grow(-lsst.geom.Extent2I(0, bbox.height//4)) 

616 exposure[bbox].mask.array |= saturatedMaskBit 

617 

618 result = calibrate.run(exposures=exposure) 

619 

620 self._check_run(calibrate, result) 

621 

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

623 def test_fail_on_sattle_miconfiguration(self): 

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

625 """ 

626 self.config.run_sattle = True 

627 with self.assertRaises(pexConfig.FieldValidationError): 

628 CalibrateImageTask(config=self.config) 

629 

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

631 def test_continue_on_sattle_failure(self): 

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

633 """ 

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

635 

636 self.config.run_sattle = True 

637 calibrate = CalibrateImageTask(config=self.config) 

638 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

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

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

641 calibrate.run(exposures=self.exposure) 

642 mock_put.assert_called_once() 

643 

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

645 def test_sattle(self): 

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

647 """ 

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

649 

650 self.config.run_sattle = True 

651 calibrate = CalibrateImageTask(config=self.config) 

652 calibrate.astrometry.setRefObjLoader(self.ref_loader) 

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

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

655 calibrate.run(exposures=self.exposure) 

656 mock_put.assert_called_once() 

657 

658 

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

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

661 but do not need real images. 

662 """ 

663 def setUp(self): 

664 instrument = "testCam" 

665 exposure0 = 101 

666 exposure1 = 102 

667 visit = 100101 

668 detector = 42 

669 physical_filter = "r" 

670 

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

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

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

674 

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

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

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

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

679 ) 

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

681 

682 # dataIds for fake data 

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

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

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