Coverage for tests/test_createApFakes.py: 96%

227 statements  

« prev     ^ index     » next       coverage.py v7.16.2, created at 2026-10-17 10:11 +0000

1# 

2# This file is part of ap_pipe. 

3# 

4# Developed for the LSST Data Management System. 

5# This product includes software developed by the LSST Project 

6# (http://www.lsst.org). 

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

8# for details of code ownership. 

9# 

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

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

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

13# (at your option) any later version. 

14# 

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

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

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

18# GNU General Public License for more details. 

19# 

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

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

22# 

23 

24import numpy as np 

25import shutil 

26import tempfile 

27import unittest 

28from unittest.mock import MagicMock 

29 

30import lsst.daf.butler.tests as butlerTests 

31import lsst.geom as geom 

32from astropy.table import Table 

33from astropy import units as u 

34from lsst.pipe.base import testUtils 

35import lsst.skymap as skyMap 

36import lsst.utils.tests 

37 

38from lsst.ap.pipe.createApFakes import ( 

39 CreateRandomApFakesTask, 

40 CreateRandomApFakesConfig, 

41 CreateVisitDetectorFakesTask, 

42 CreateVisitDetectorFakesConfig, 

43) 

44 

45 

46class TestCreateApFakes(lsst.utils.tests.TestCase): 

47 

48 def setUp(self): 

49 """ 

50 """ 

51 self.tractId = 0 

52 self.rng = np.random.default_rng(self.tractId) 

53 simpleMapConfig = skyMap.discreteSkyMap.DiscreteSkyMapConfig() 

54 simpleMapConfig.raList = [10] 

55 simpleMapConfig.decList = [-1] 

56 simpleMapConfig.radiusList = [0.1] 

57 self.simpleMap = skyMap.DiscreteSkyMap(simpleMapConfig) 

58 self.tract = self.simpleMap.generateTract(self.tractId) 

59 

60 bBox = self.tract.getOuterSkyPolygon().getBoundingBox() 

61 self.nSources = 50 

62 self.sourceDensity = (self.nSources 

63 / (bBox.getArea() * (180 / np.pi) ** 2)) 

64 self.fraction = 0.5 

65 self.nInVisit = (int(self.nSources * self.fraction) 

66 + int((1 - self.fraction) / 2 * self.nSources)) 

67 self.nInTemplate = (self.nSources - self.nInVisit 

68 + int(self.nSources * self.fraction)) 

69 

70 def testRunQuantum(self): 

71 """Test the run quantum method with a gen3 butler. 

72 """ 

73 root = tempfile.mkdtemp() 

74 dimensions = {"instrument": ["notACam"], 

75 "skymap": ["skyMap"], 

76 "tract": [0, 42], 

77 } 

78 testRepo = butlerTests.makeTestRepo(root, dimensions) 

79 with self.assertWarns(FutureWarning): 

80 fakesTask = CreateRandomApFakesTask() 

81 connections = fakesTask.config.ConnectionsClass( 

82 config=fakesTask.config) 

83 butlerTests.addDatasetType( 

84 testRepo, 

85 connections.skyMap.name, 

86 connections.skyMap.dimensions, 

87 connections.skyMap.storageClass) 

88 butlerTests.addDatasetType( 

89 testRepo, 

90 connections.fakeCat.name, 

91 connections.fakeCat.dimensions, 

92 connections.fakeCat.storageClass) 

93 

94 dataId = {"skymap": "skyMap", "tract": 0} 

95 butler = butlerTests.makeTestCollection(testRepo) 

96 butler.put(self.simpleMap, "skyMap", {"skymap": "skyMap"}) 

97 

98 quantum = testUtils.makeQuantum( 

99 fakesTask, butler, dataId, 

100 {"skyMap": {"skymap": dataId["skymap"]}, "fakeCat": dataId}) 

101 run = testUtils.runTestQuantum(fakesTask, butler, quantum, True) 

102 # Actual input dataset omitted for simplicity 

103 run.assert_called_once_with(tractId=dataId["tract"], skyMap=self.simpleMap) 

104 shutil.rmtree(root, ignore_errors=True) 

105 

106 def testRun(self): 

107 """Test the run method. 

108 """ 

109 with self.assertWarns(FutureWarning): 

110 fakesConfig = CreateRandomApFakesConfig() 

111 fakesConfig.fraction = 0.5 

112 fakesConfig.fakeDensity = self.sourceDensity 

113 

114 with self.assertWarns(FutureWarning): 

115 fakesTask = CreateRandomApFakesTask(config=fakesConfig) 

116 bBox = self.tract.getOuterSkyPolygon().getBoundingBox() 

117 result = fakesTask.run(self.tractId, self.simpleMap) 

118 fakeCat = result.fakeCat 

119 self.assertEqual(len(fakeCat), self.nSources) 

120 self.assertIn("injection_id", fakeCat.columns) 

121 self.assertIn("source_type", fakeCat.columns) 

122 self.assertIn("mag", fakeCat.columns) 

123 self.assertTrue(np.issubdtype(fakeCat["injection_id"].dtype, np.integer)) 

124 self.assertEqual(fakeCat["injection_id"].nunique(), len(fakeCat)) 

125 self.assertTrue(np.all(fakeCat["source_type"] == "Star")) 

126 

127 for idx, row in fakeCat.iterrows(): 

128 self.assertTrue( 

129 bBox.contains( 

130 geom.SpherePoint(row[fakesTask.config.ra_col], 

131 row[fakesTask.config.dec_col], 

132 geom.degrees).getVector())) 

133 self.assertEqual(fakeCat[fakesConfig.visitSourceFlagCol].sum(), 

134 self.nInVisit) 

135 self.assertEqual(fakeCat[fakesConfig.templateSourceFlagCol].sum(), 

136 self.nInTemplate) 

137 for f in fakesConfig.filterSet: 

138 filterMags = fakeCat[fakesConfig.mag_col % f] 

139 self.assertEqual(self.nSources, len(filterMags)) 

140 self.assertTrue( 

141 np.all(fakesConfig.magMin <= filterMags)) 

142 self.assertTrue( 

143 np.all(fakesConfig.magMax > filterMags)) 

144 self.assertTrue(np.allclose(filterMags, fakeCat["mag"])) 

145 

146 def testVisitCoaddSubdivision(self): 

147 """Test that the number of assigned visit to template objects is 

148 correct. 

149 """ 

150 with self.assertWarns(FutureWarning): 

151 fakesConfig = CreateRandomApFakesConfig() 

152 fakesConfig.fraction = 0.5 

153 with self.assertWarns(FutureWarning): 

154 fakesTask = CreateRandomApFakesTask(config=fakesConfig) 

155 subdivision = fakesTask.createVisitCoaddSubdivision(self.nSources) 

156 self.assertEqual( 

157 subdivision[fakesConfig.visitSourceFlagCol].sum(), 

158 self.nInVisit) 

159 self.assertEqual( 

160 subdivision[fakesConfig.templateSourceFlagCol].sum(), 

161 self.nInTemplate) 

162 

163 

164def _make_mock_visit_image(visitId=2024111100094, detId=3, 

165 xmin=0, xmax=4096, ymin=0, ymax=4096, 

166 magLim=25.0, ra_center=10.0, dec_center=-1.0, 

167 return_quantity_angles=False): 

168 """Build a minimal MagicMock that satisfies CreateVisitDetectorFakesTask.run.""" 

169 img = MagicMock() 

170 

171 # visit / detector IDs 

172 img.getInfo().getVisitInfo().id = visitId 

173 img.detector.getId.return_value = detId 

174 

175 # bounding box 

176 bbox = MagicMock() 

177 bbox.getMinX.return_value = xmin 

178 bbox.getMaxX.return_value = xmax 

179 bbox.getMinY.return_value = ymin 

180 bbox.getMaxY.return_value = ymax 

181 img.getBBox.return_value = bbox 

182 

183 # summary stats — magLim drives max_mag 

184 stats = MagicMock() 

185 stats.magLim = magLim 

186 img.getInfo().getSummaryStats.return_value = stats 

187 

188 # convex polygon for density calculation (returns a bbox in steradians) 

189 poly_bbox = MagicMock() 

190 # ~1e-5 sr ≈ small patch, dense enough to yield O(10) fakes at density=1000 

191 poly_bbox.getArea.return_value = 1e-5 

192 img.getConvexPolygon().getBoundingBox.return_value = poly_bbox 

193 

194 # WCS: pixelToSkyArray returns ra/dec arrays of the right length 

195 wcs = MagicMock() 

196 

197 def _pix_to_sky(xs, ys, degrees=True): 

198 ra = ra_center + xs * 1e-4 

199 dec = dec_center + ys * 1e-4 

200 if return_quantity_angles: 200 ↛ 201line 200 didn't jump to line 201 because the condition on line 200 was never true

201 return np.asarray(ra) * u.deg, np.asarray(dec) * u.deg 

202 return np.asarray(ra), np.asarray(dec) 

203 

204 def _sky_to_pix(ra, dec, degrees=False): 

205 ra_arr = np.asarray(ra) 

206 dec_arr = np.asarray(dec) 

207 if not degrees: 207 ↛ 208line 207 didn't jump to line 208 because the condition on line 207 was never true

208 ra_arr = np.rad2deg(ra_arr) 

209 dec_arr = np.rad2deg(dec_arr) 

210 xs = (ra_arr - ra_center) / 1e-4 

211 ys = (dec_arr - dec_center) / 1e-4 

212 return np.asarray(xs), np.asarray(ys) 

213 

214 wcs.pixelToSkyArray.side_effect = _pix_to_sky 

215 wcs.skyToPixelArray.side_effect = _sky_to_pix 

216 img.getWcs.return_value = wcs 

217 

218 # photoCalib — not used in non-hosted paths, but must exist 

219 img.getPhotoCalib.return_value = MagicMock() 

220 

221 return img 

222 

223 

224def _make_calibrated_source_table(n_sources, host_kind="galaxy", ra_deg=10.0, dec_deg=-1.0): 

225 """Build a calibrated source table that can pass host selection cuts.""" 

226 if host_kind == "galaxy": 

227 size_ext = np.ones(n_sources) 

228 ext = np.ones(n_sources, dtype=int) 

229 elif host_kind == "star": 229 ↛ 233line 229 didn't jump to line 233 because the condition on line 229 was always true

230 size_ext = np.full(n_sources, 0.2) 

231 ext = np.zeros(n_sources, dtype=int) 

232 else: 

233 raise ValueError(f"Unknown host_kind={host_kind}") 

234 

235 return Table({ 

236 "slot_Centroid_x": np.linspace(1000.0, 3000.0, n_sources), 

237 "slot_Centroid_y": np.linspace(1100.0, 3100.0, n_sources), 

238 "slot_ModelFlux_mag": np.full(n_sources, 20.0), 

239 "slot_ModelFlux_flux": np.full(n_sources, 1e4), 

240 "slot_ModelFlux_fluxErr": np.full(n_sources, 100.0), 

241 "slot_PsfFlux_mag": np.full(n_sources, 20.0), 

242 "slot_PsfFlux_flux": np.full(n_sources, 1e4), 

243 "slot_PsfFlux_fluxErr": np.full(n_sources, 100.0), 

244 "slot_Shape_xx": np.full(n_sources, 4.0), 

245 "slot_Shape_xy": np.zeros(n_sources), 

246 "slot_Shape_yy": np.full(n_sources, 4.0), 

247 "id": np.arange(n_sources, dtype=np.int64), 

248 "coord_ra": np.deg2rad(np.full(n_sources, ra_deg)), 

249 "coord_dec": np.deg2rad(np.full(n_sources, dec_deg)), 

250 "sky_source": np.zeros(n_sources, dtype=bool), 

251 "base_ClassificationSizeExtendedness_flag": np.zeros(n_sources, dtype=bool), 

252 "base_ClassificationExtendedness_flag": np.zeros(n_sources, dtype=bool), 

253 "slot_Shape_flag": np.zeros(n_sources, dtype=bool), 

254 "slot_Centroid_flag": np.zeros(n_sources, dtype=bool), 

255 "base_PixelFlags_flag": np.zeros(n_sources, dtype=bool), 

256 "base_ClassificationSizeExtendedness_value": size_ext, 

257 "base_ClassificationExtendedness_value": ext, 

258 }) 

259 

260 

261class TestCreateVisitDetectorFakesTask(lsst.utils.tests.TestCase): 

262 

263 def setUp(self): 

264 self.visit_image = _make_mock_visit_image() 

265 self.source_cat = MagicMock() # not used unless doAddHostedFakes 

266 

267 def _make_task(self, **config_overrides): 

268 cfg = CreateVisitDetectorFakesConfig() 

269 cfg.doAddRandomVisitFakes = True 

270 cfg.doAddRandomTemplateFakes = False 

271 cfg.doAddHostedFakes = False 

272 cfg.doAddVariableFakes = False 

273 cfg.nRandomFakes = 50 

274 for k, v in config_overrides.items(): 

275 setattr(cfg, k, v) 

276 return CreateVisitDetectorFakesTask(config=cfg) 

277 

278 # ------------------------------------------------------------------ 

279 # Basic smoke test: random-only path 

280 # ------------------------------------------------------------------ 

281 def testRunRandomOnly(self): 

282 task = self._make_task() 

283 result = task.run(self.source_cat, self.visit_image) 

284 cat = result.outputCat 

285 

286 self.assertIsInstance(cat, Table) 

287 self.assertEqual(len(cat), 50) 

288 # Required columns 

289 for col in ("ra", "dec", "mag", "x", "y", "source_type", 

290 "injection_id", "visit", "detector"): 

291 self.assertIn(col, cat.colnames) 

292 # All sources are stars 

293 self.assertTrue(np.all(cat["source_type"] == "Star")) 

294 # IDs are unique 

295 self.assertEqual(len(np.unique(cat["injection_id"])), len(cat)) 

296 # visit/detector are set correctly 

297 self.assertTrue(np.all(cat["visit"] == 2024111100094)) 

298 self.assertTrue(np.all(cat["detector"] == 3)) 

299 

300 # ------------------------------------------------------------------ 

301 # Template-fraction split 

302 # ------------------------------------------------------------------ 

303 def testTemplateFakeFraction(self): 

304 task = self._make_task( 

305 doAddRandomTemplateFakes=True, 

306 templateFakeFraction=0.25, 

307 nRandomFakes=200, 

308 ) 

309 result = task.run(self.source_cat, self.visit_image) 

310 cat = result.outputCat 

311 

312 n_template = int(np.sum(cat["isTemplateSource"])) 

313 # No source should be both visit AND template (mutually exclusive after split) 

314 self.assertFalse(np.any(cat["isVisitSource"] & cat["isTemplateSource"])) 

315 # template fraction should be roughly 25 % 

316 frac = n_template / len(cat) 

317 self.assertAlmostEqual(frac, 0.25, delta=0.10) 

318 

319 # ------------------------------------------------------------------ 

320 # Variable fakes: twin_id bookkeeping 

321 # ------------------------------------------------------------------ 

322 def testVariableFakesTwinTracking(self): 

323 task = self._make_task( 

324 doAddVariableFakes=True, 

325 variableFakeFraction=0.2, 

326 nRandomFakes=100, 

327 ) 

328 result = task.run(self.source_cat, self.visit_image) 

329 cat = result.outputCat 

330 

331 # Catalog is larger than the base 100 fakes 

332 self.assertGreater(len(cat), 100) 

333 # twin_id column exists 

334 self.assertIn("twin_id", cat.colnames) 

335 # Every injection_id is unique 

336 self.assertEqual(len(np.unique(cat["injection_id"])), len(cat)) 

337 # mag_offset column is present on variable rows 

338 self.assertIn("mag_offset", cat.colnames) 

339 

340 # ------------------------------------------------------------------ 

341 # Empty-mode guard raises RuntimeError 

342 # ------------------------------------------------------------------ 

343 def testEmptyCatalogRaises(self): 

344 task = self._make_task( 

345 doAddRandomVisitFakes=False, 

346 doAddHostedFakes=False, 

347 ) 

348 with self.assertRaises(RuntimeError): 

349 task.run(self.source_cat, self.visit_image) 

350 

351 # ------------------------------------------------------------------ 

352 # Hosted fakes: zero valid hosts — warning issued, no crash when 

353 # random fakes are also enabled as fallback 

354 # ------------------------------------------------------------------ 

355 def testHostedFakesNoHostsWarning(self): 

356 task = self._make_task( 

357 doAddRandomVisitFakes=True, 

358 nRandomFakes=20, 

359 doAddHostedFakes=True, 

360 fracHostedFakes=0.1, 

361 ) 

362 empty_table = Table({ 

363 "slot_Centroid_x": np.array([]), 

364 "slot_Centroid_y": np.array([]), 

365 "slot_ModelFlux_mag": np.array([]), 

366 "slot_ModelFlux_flux": np.array([]), 

367 "slot_ModelFlux_fluxErr": np.array([]), 

368 "slot_Shape_xx": np.array([]), 

369 "slot_Shape_xy": np.array([]), 

370 "slot_Shape_yy": np.array([]), 

371 "id": np.array([], dtype=np.int64), 

372 "coord_ra": np.array([]), 

373 "coord_dec": np.array([]), 

374 "sky_source": np.array([], dtype=bool), 

375 "base_ClassificationSizeExtendedness_flag": np.array([], dtype=bool), 

376 "base_ClassificationExtendedness_flag": np.array([], dtype=bool), 

377 "slot_Shape_flag": np.array([], dtype=bool), 

378 "slot_Centroid_flag": np.array([], dtype=bool), 

379 "base_PixelFlags_flag": np.array([], dtype=bool), 

380 "base_ClassificationSizeExtendedness_value": np.array([]), 

381 "base_ClassificationExtendedness_value": np.array([], dtype=int), 

382 }) 

383 self.visit_image.getPhotoCalib().calibrateCatalog.return_value.asAstropy.return_value = ( 

384 empty_table 

385 ) 

386 

387 with self.assertLogs("lsst.ap.pipe.createApFakes", level="WARNING") as cm: 

388 result = task.run(self.source_cat, self.visit_image) 

389 

390 self.assertTrue(any("no valid hosts" in msg.lower() for msg in cm.output)) 

391 # Random fakes still produced 

392 self.assertEqual(len(result.outputCat), 18) 

393 

394 def testHostedBlendedFakesPerHostCap(self): 

395 task = self._make_task( 

396 doAddRandomVisitFakes=True, 

397 nRandomFakes=100, 

398 doAddBlendedFakes=True, 

399 fracBlendedFakes=0.25, 

400 fracHostedBlendedFakes=0.5, 

401 maxHostedBlendedFakesPerHost=2, 

402 maxHostedBlendedFakesTotal=-1, 

403 ) 

404 star_hosts = _make_calibrated_source_table(5, host_kind="star") 

405 self.visit_image.getPhotoCalib().calibrateCatalog.return_value.asAstropy.return_value = star_hosts 

406 

407 cat = task.run(self.source_cat, self.visit_image).outputCat 

408 hosted_blended = cat[cat["isBlended"] & cat["hosted_fake"]] 

409 hosted_mask = ~np.ma.getmaskarray(hosted_blended["host_flux"]) 

410 hosted_star_blended = hosted_blended[hosted_mask] 

411 

412 # Requested 20, capped to n_hosts * per_host = 5 * 2 = 6 

413 self.assertLessEqual(len(hosted_star_blended), 8) 

414 host_ids, counts = np.unique(hosted_star_blended["host_id"], return_counts=True) 

415 self.assertLessEqual(len(host_ids), 5) 

416 self.assertTrue(np.all(counts[host_ids > 0] <= 2)) 

417 

418 def testHostedBlendedFakesTotalCap(self): 

419 task = self._make_task( 

420 doAddRandomVisitFakes=True, 

421 nRandomFakes=300, 

422 doAddBlendedFakes=True, 

423 fracBlendedFakes=0.25, 

424 fracHostedBlendedFakes=0.5, 

425 maxHostedBlendedFakesPerHost=10, 

426 maxHostedBlendedFakesTotal=4, 

427 ) 

428 star_hosts = _make_calibrated_source_table(10, host_kind="star") 

429 self.visit_image.getPhotoCalib().calibrateCatalog.return_value.asAstropy.return_value = star_hosts 

430 

431 cat = task.run(self.source_cat, self.visit_image).outputCat 

432 hosted_blended = cat[cat["isBlended"] & cat["hosted_fake"]] 

433 hosted_mask = ~np.ma.getmaskarray(hosted_blended["host_flux"]) 

434 hosted_star_blended = hosted_blended[hosted_mask] 

435 

436 self.assertEqual(len(hosted_star_blended), 4) 

437 

438 def testAngleColumnsAreUnitless(self): 

439 cfg = CreateVisitDetectorFakesConfig() 

440 cfg.doAddRandomVisitFakes = True 

441 cfg.doAddRandomTemplateFakes = False 

442 cfg.doAddHostedFakes = True 

443 cfg.doAddVariableFakes = False 

444 cfg.fracHostedFakes = 0.5 

445 task = CreateVisitDetectorFakesTask(config=cfg) 

446 

447 img = _make_mock_visit_image() 

448 host_table = _make_calibrated_source_table(5, host_kind="galaxy") 

449 img.getPhotoCalib().calibrateCatalog.return_value.asAstropy.return_value = host_table 

450 

451 cat = task.run(self.source_cat, img).outputCat 

452 for col in ("ra", "dec", "delta_ra", "delta_dec", "host_ra", "host_dec"): 

453 self.assertIn(col, cat.colnames) 

454 self.assertIsNone(cat[col].unit) 

455 

456 # Sanity check that RA values are in degrees-like scale, not radians. 

457 self.assertTrue(np.all(np.asarray(cat["ra"], dtype=float) > 2.0 * np.pi)) 

458 hosted_cat = cat[cat['hosted_fake']] 

459 self.assertTrue(np.all(np.asarray(hosted_cat["host_ra"], dtype=float) > 2.0 * np.pi)) 

460 

461 

462class MemoryTester(lsst.utils.tests.MemoryTestCase): 

463 pass 

464 

465 

466def setup_module(module): 

467 lsst.utils.tests.init() 

468 

469 

470if __name__ == "__main__": 470 ↛ 471line 470 didn't jump to line 471 because the condition on line 470 was never true

471 lsst.utils.tests.init() 

472 unittest.main()