Coverage for metadetect / tests / test_metadetect.py: 9%

473 statements  

« prev     ^ index     » next       coverage.py v7.13.1, created at 2026-01-22 02:30 -0800

1""" 

2test with super simple sim. The purpose here is not 

3to make sure it gets the right answer or anything, just 

4to test all the moving parts 

5""" 

6import time 

7import copy 

8import itertools 

9 

10import pytest 

11 

12from packaging import version 

13import ngmix 

14import numpy as np 

15 

16from metadetect import metadetect 

17from metadetect import fitting 

18from metadetect import procflags 

19from .sim import Sim 

20 

21 

22TEST_METADETECT_CONFIG = { 

23 "model": "wmom", 

24 

25 "weight": { 

26 "fwhm": 1.2, # arcsec 

27 }, 

28 

29 "metacal": { 

30 "psf": "fitgauss", 

31 "types": ["noshear", "1p", "1m", "2p", "2m"], 

32 }, 

33 

34 "sx": { 

35 # in sky sigma 

36 # DETECT_THRESH 

37 "detect_thresh": 0.8, 

38 

39 # Minimum contrast parameter for deblending 

40 # DEBLEND_MINCONT 

41 "deblend_cont": 0.00001, 

42 

43 # minimum number of pixels above threshold 

44 # DETECT_MINAREA: 6 

45 "minarea": 4, 

46 

47 "filter_type": "conv", 

48 

49 # 7x7 convolution mask of a gaussian PSF with FWHM = 3.0 pixels. 

50 "filter_kernel": [ 

51 [0.004963, 0.021388, 0.051328, 0.068707, 0.051328, 0.021388, 0.004963], # noqa 

52 [0.021388, 0.092163, 0.221178, 0.296069, 0.221178, 0.092163, 0.021388], # noqa 

53 [0.051328, 0.221178, 0.530797, 0.710525, 0.530797, 0.221178, 0.051328], # noqa 

54 [0.068707, 0.296069, 0.710525, 0.951108, 0.710525, 0.296069, 0.068707], # noqa 

55 [0.051328, 0.221178, 0.530797, 0.710525, 0.530797, 0.221178, 0.051328], # noqa 

56 [0.021388, 0.092163, 0.221178, 0.296069, 0.221178, 0.092163, 0.021388], # noqa 

57 [0.004963, 0.021388, 0.051328, 0.068707, 0.051328, 0.021388, 0.004963], # noqa 

58 ] 

59 }, 

60 

61 "meds": { 

62 "min_box_size": 32, 

63 "max_box_size": 256, 

64 

65 "box_type": "iso_radius", 

66 

67 "rad_min": 4, 

68 "rad_fac": 2, 

69 "box_padding": 2, 

70 }, 

71 

72 # check for an edge hit 

73 "bmask_flags": 2**30, 

74 

75 "nodet_flags": 2**0, 

76} 

77 

78 

79def _show_mbobs(mer): 

80 import images 

81 

82 mbobs = mer.mbobs 

83 

84 rgb = images.get_color_image( 

85 mbobs[2][0].image.transpose(), 

86 mbobs[1][0].image.transpose(), 

87 mbobs[0][0].image.transpose(), 

88 nonlinear=0.1, 

89 ) 

90 rgb *= 1.0/rgb.max() 

91 

92 images.view_mosaic( 

93 [rgb, 

94 mer.seg, 

95 mer.detim], 

96 titles=["image", "seg", "detim"], 

97 ) 

98 

99 

100def test_detect(ntrial=1, show=False): 

101 """ 

102 just test the detection 

103 """ 

104 pytest.importorskip("meds") 

105 pytest.importorskip("sxdes") 

106 from .. import detect 

107 

108 rng = np.random.RandomState(seed=45) 

109 

110 tm0 = time.time() 

111 nobj_meas = 0 

112 

113 sim = Sim(rng) 

114 

115 config = {} 

116 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

117 

118 for trial in range(ntrial): 

119 print("trial: %d/%d" % (trial+1, ntrial)) 

120 

121 mbobs = sim.get_mbobs() 

122 mer = detect.MEDSifier( 

123 mbobs=mbobs, 

124 sx_config=config["sx"], 

125 meds_config=config["meds"], 

126 ) 

127 

128 mbm = mer.get_multiband_meds() 

129 

130 nobj = mbm.size 

131 nobj_meas += nobj 

132 

133 if show: 

134 _show_mbobs(mer) 

135 if ntrial > 1 and trial != (ntrial-1): 

136 if "q" == input("hit a key: "): 

137 return 

138 

139 total_time = time.time()-tm0 

140 print("found", nobj_meas, "objects") 

141 print("time per group:", total_time/ntrial) 

142 print("time per object:", total_time/nobj_meas) 

143 

144 

145def test_detect_masking(ntrial=1, show=False): 

146 """ 

147 just test the detection 

148 """ 

149 pytest.importorskip("meds") 

150 pytest.importorskip("sxdes") 

151 from .. import detect 

152 

153 rng = np.random.RandomState(seed=45) 

154 

155 sim = Sim(rng) 

156 

157 config = {} 

158 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

159 

160 for trial in range(ntrial): 

161 print("trial: %d/%d" % (trial+1, ntrial)) 

162 

163 mbobs = sim.get_mbobs() 

164 for obslist in mbobs: 

165 for obs in obslist: 

166 obs.bmask = obs.bmask | config["nodet_flags"] 

167 

168 mer = detect.MEDSifier( 

169 mbobs=mbobs, 

170 sx_config=config["sx"], 

171 meds_config=config["meds"], 

172 nodet_flags=config["nodet_flags"], 

173 ) 

174 assert mer.cat.size == 0 

175 

176 

177def _check_result_array(res, shear, msk, model): 

178 for col in res[shear].dtype.names: 

179 if col == "shear_bands": 

180 assert np.all(res[shear][msk][col] == "012") 

181 elif col == "det_bands": 

182 assert np.all(res[shear][msk][col] == "012") 

183 else: 

184 # admom doesn't make band fluxes 

185 if model in ["admom", "am"] and "band_flux" in col: 

186 if col.endswith("band_flux_flags"): 

187 assert np.array_equal( 

188 res[shear][msk][col], 

189 np.zeros_like(res[shear][msk][col]) 

190 + procflags.NO_ATTEMPT, 

191 ), ( 

192 "result column '%s' is not NO_ATTEMPT: %s" % ( 

193 col, res[shear][msk][col] 

194 ) 

195 ) 

196 elif any( 

197 col.endswith(s) for s in ["band_flux", "band_flux_err"] 

198 ): 

199 assert np.all(np.isnan( 

200 res[shear][msk][col], 

201 )), ( 

202 "result column '%s' is not NaN: %s" % ( 

203 col, res[shear][msk][col] 

204 ) 

205 ) 

206 else: 

207 assert np.all(np.isfinite(res[shear][msk][col])), ( 

208 "result column '%s' has NaNs: %s" % ( 

209 col, res[shear][msk][col] 

210 ) 

211 ) 

212 

213 

214@pytest.mark.parametrize("model", ["gauss", "wmom"]) 

215def test_metadetect_coadd_faster(model): 

216 """ 

217 test coadding is faster 

218 """ 

219 pytest.importorskip("sxdes") 

220 

221 config = {} 

222 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

223 del config["model"] 

224 del config["weight"] 

225 config["fitters"] = [ 

226 {"model": model, "coadd": False, "weight": {"fwhm": 1.2}} 

227 ] 

228 

229 # warm up once 

230 mbobs = Sim( 

231 np.random.RandomState(seed=116), config={"nband": 6} 

232 ).get_mbobs() 

233 metadetect.do_metadetect( 

234 config, mbobs, np.random.RandomState(seed=116) 

235 ) 

236 

237 mbobs = Sim( 

238 np.random.RandomState(seed=116), config={"nband": 6} 

239 ).get_mbobs() 

240 tm0 = time.time() 

241 metadetect.do_metadetect( 

242 config, mbobs, np.random.RandomState(seed=116) 

243 ) 

244 no_coadd_time = time.time() - tm0 

245 

246 config["fitters"][0]["coadd"] = True 

247 mbobs = Sim( 

248 np.random.RandomState(seed=116), config={"nband": 6} 

249 ).get_mbobs() 

250 tm0 = time.time() 

251 metadetect.do_metadetect( 

252 config, mbobs, np.random.RandomState(seed=116) 

253 ) 

254 coadd_time = time.time() - tm0 

255 

256 print("coadd|nocoadd: %f|%f" % (coadd_time, no_coadd_time)) 

257 

258 if model == "gauss": 

259 assert coadd_time < no_coadd_time*0.7, (coadd_time, no_coadd_time) 

260 else: 

261 assert np.allclose(coadd_time, no_coadd_time, atol=0, rtol=0.3) 

262 

263 

264@pytest.mark.parametrize("model", ["wmom", "pgauss", "ksigma", "am", "gauss"]) 

265def test_metadetect_smoke(model): 

266 """ 

267 test full metadetection 

268 """ 

269 pytest.importorskip("sxdes") 

270 

271 ntrial = 1 

272 rng = np.random.RandomState(seed=116) 

273 

274 tm0 = time.time() 

275 

276 sim = Sim(rng) 

277 config = {} 

278 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

279 config["model"] = model 

280 

281 for trial in range(ntrial): 

282 print("trial: %d/%d" % (trial+1, ntrial)) 

283 

284 mbobs = sim.get_mbobs() 

285 res = metadetect.do_metadetect(config, mbobs, rng) 

286 for shear in ["noshear", "1p", "1m", "2p", "2m"]: 

287 assert np.all(res[shear]["mfrac"] == 0) 

288 assert any(c.endswith("band_flux") for c in res[shear].dtype.names) 

289 assert np.any(res[shear]["psfrec_g"] != 0) 

290 assert np.any(res[shear]["psfrec_T"] != 0) 

291 msk = res[shear][model + '_flags'] == 0 

292 _check_result_array(res, shear, msk, model) 

293 

294 total_time = time.time()-tm0 

295 print("time per:", total_time/ntrial) 

296 

297 

298@pytest.mark.parametrize("model", ["wmom", "pgauss", "ksigma", "am", "gauss"]) 

299def test_metadetect_uberseg(model): 

300 """ 

301 test full metadetection 

302 """ 

303 pytest.importorskip("sxdes") 

304 

305 ntrial = 1 

306 rng = np.random.RandomState(seed=116) 

307 

308 tm0 = time.time() 

309 

310 sim = Sim(rng) 

311 config = {} 

312 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

313 config["model"] = model 

314 config["meds"]["weight_type"] = "uberseg" 

315 

316 mbobs = sim.get_mbobs() 

317 res = metadetect.do_metadetect(config, mbobs, rng) 

318 

319 for trial in range(ntrial): 

320 print("trial: %d/%d" % (trial+1, ntrial)) 

321 

322 res = metadetect.do_metadetect(config, mbobs, rng) 

323 for shear in ["noshear", "1p", "1m", "2p", "2m"]: 

324 assert np.all(res[shear]["mfrac"] == 0) 

325 assert any(c.endswith("band_flux") for c in res[shear].dtype.names) 

326 assert np.any(res[shear]["psfrec_g"] != 0) 

327 assert np.any(res[shear]["psfrec_T"] != 0) 

328 msk = res[shear][model + '_flags'] == 0 

329 _check_result_array(res, shear, msk, model) 

330 

331 total_time = time.time()-tm0 

332 print("time per:", total_time/ntrial) 

333 

334 

335@pytest.mark.parametrize("model", ["wmom", "pgauss", "ksigma", "am", "gauss"]) 

336def test_metadetect_mfrac(model): 

337 """ 

338 test full metadetection w/ mfrac 

339 """ 

340 pytest.importorskip("sxdes") 

341 

342 ntrial = 1 

343 rng = np.random.RandomState(seed=53341) 

344 

345 tm0 = time.time() 

346 

347 sim = Sim(rng) 

348 config = {} 

349 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

350 config["model"] = model 

351 

352 for trial in range(ntrial): 

353 print("trial: %d/%d" % (trial+1, ntrial)) 

354 

355 mbobs = sim.get_mbobs() 

356 for band in range(len(mbobs)): 

357 mbobs[band][0].mfrac = rng.uniform( 

358 size=mbobs[band][0].image.shape, low=0.2, high=0.8 

359 ) 

360 res = metadetect.do_metadetect(config, mbobs, rng) 

361 for shear in ["noshear", "1p", "1m", "2p", "2m"]: 

362 assert np.all( 

363 (res[shear]["mfrac"] > 0.45) 

364 & (res[shear]["mfrac"] < 0.55) 

365 ) 

366 assert np.all( 

367 (res[shear]["mfrac_img"] > 0.45) 

368 & (res[shear]["mfrac_img"] < 0.55) 

369 ) 

370 assert np.all( 

371 (res[shear]["mfrac_noshear"] > 0.45) 

372 & (res[shear]["mfrac_noshear"] < 0.55) 

373 ) 

374 msk = res[shear][model + '_flags'] == 0 

375 _check_result_array(res, shear, msk, model) 

376 

377 total_time = time.time()-tm0 

378 print("time per:", total_time/ntrial) 

379 

380 

381@pytest.mark.parametrize("model", ["wmom", "pgauss", "ksigma", "am", "gauss"]) 

382def test_metadetect_mfrac_all(model): 

383 """ 

384 test full metadetection w/ mfrac all 1 

385 """ 

386 ntrial = 1 

387 rng = np.random.RandomState(seed=53341) 

388 

389 sim = Sim(rng) 

390 config = {} 

391 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

392 config["model"] = model 

393 

394 for trial in range(ntrial): 

395 print("trial: %d/%d" % (trial+1, ntrial)) 

396 

397 mbobs = sim.get_mbobs() 

398 for band in range(len(mbobs)): 

399 mbobs[band][0].mfrac = np.ones_like(mbobs[band][0].image) 

400 

401 res = metadetect.do_metadetect(config, mbobs, rng) 

402 assert res is None 

403 

404 

405@pytest.mark.parametrize("model", ["wmom", "pgauss", "ksigma", "am", "gauss"]) 

406def test_metadetect_zero_weight_all(model): 

407 """ 

408 test full metadetection w/ all zero weight 

409 """ 

410 

411 ntrial = 1 

412 rng = np.random.RandomState(seed=53341) 

413 

414 sim = Sim(rng) 

415 config = {} 

416 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

417 config["model"] = model 

418 

419 for trial in range(ntrial): 

420 print("trial: %d/%d" % (trial+1, ntrial)) 

421 

422 mbobs = sim.get_mbobs() 

423 for band in range(len(mbobs)): 

424 mbobs[band][0].weight = np.zeros_like(mbobs[band][0].image) 

425 

426 res = metadetect.do_metadetect(config, mbobs, rng) 

427 assert res is None 

428 

429 

430@pytest.mark.parametrize("model", ["wmom", "pgauss", "ksigma", "am", "gauss"]) 

431def test_metadetect_zero_weight_some(model): 

432 """ 

433 test full metadetection w/ some zero weight 

434 """ 

435 ntrial = 1 

436 rng = np.random.RandomState(seed=53341) 

437 

438 sim = Sim(rng) 

439 config = {} 

440 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

441 config["model"] = model 

442 

443 for trial in range(ntrial): 

444 print("trial: %d/%d" % (trial+1, ntrial)) 

445 

446 mbobs = sim.get_mbobs() 

447 for band in range(len(mbobs)): 

448 if band == 1: 

449 mbobs[band][0].weight = np.zeros_like(mbobs[band][0].image) 

450 

451 res = metadetect.do_metadetect(config, mbobs, rng) 

452 assert res is None 

453 

454 

455@pytest.mark.parametrize("model", ["wmom", "pgauss", "ksigma", "am", "gauss"]) 

456def test_metadetect_nodet_flags_all(model): 

457 """ 

458 test full metadetection w/ all bmask all nodet_flags 

459 """ 

460 ntrial = 1 

461 rng = np.random.RandomState(seed=53341) 

462 

463 sim = Sim(rng) 

464 config = {} 

465 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

466 config["model"] = model 

467 

468 for trial in range(ntrial): 

469 print("trial: %d/%d" % (trial+1, ntrial)) 

470 

471 mbobs = sim.get_mbobs() 

472 for band in range(len(mbobs)): 

473 mbobs[band][0].bmask = np.ones_like( 

474 mbobs[band][0].image, dtype=np.int32 

475 ) 

476 

477 res = metadetect.do_metadetect(config, mbobs, rng) 

478 assert res is None 

479 

480 

481@pytest.mark.parametrize("model", ["wmom", "pgauss", "ksigma", "am", "gauss"]) 

482def test_metadetect_nodet_flags_some(model): 

483 """ 

484 test full metadetection w/ some bmask nodet_flags 

485 """ 

486 ntrial = 1 

487 rng = np.random.RandomState(seed=53341) 

488 

489 sim = Sim(rng) 

490 config = {} 

491 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

492 config["model"] = model 

493 

494 for trial in range(ntrial): 

495 print("trial: %d/%d" % (trial+1, ntrial)) 

496 

497 mbobs = sim.get_mbobs() 

498 for band in range(len(mbobs)): 

499 if band == 1: 

500 mbobs[band][0].bmask = np.ones_like( 

501 mbobs[band][0].image, dtype=np.int32 

502 ) 

503 

504 res = metadetect.do_metadetect(config, mbobs, rng) 

505 assert res is None 

506 

507 

508@pytest.mark.skipif( 

509 version.parse(ngmix.__version__) < version.parse("2.1.0"), 

510 reason="ngmix version 2.1.0 or greater is needed for smoothing prepsf moments", 

511) 

512@pytest.mark.parametrize("model", ["pgauss", "ksigma"]) 

513def test_metadetect_fitter_fwhm_smooth(model): 

514 pytest.importorskip("sxdes") 

515 

516 nband = 3 

517 rng = np.random.RandomState(seed=116) 

518 

519 sim = Sim(rng, config={"nband": nband}) 

520 mbobs = sim.get_mbobs() 

521 

522 config = {} 

523 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

524 config["model"] = model 

525 config["weight"]["fwhm"] = 1.2 

526 

527 md = metadetect.Metadetect( 

528 config, mbobs, rng, 

529 ) 

530 md.go() 

531 res = md.result 

532 assert md._fitters[0].fwhm_smooth == 0 

533 

534 config["weight"]["fwhm_smooth"] = 0.8 

535 md = metadetect.Metadetect( 

536 config, mbobs, rng, 

537 ) 

538 md.go() 

539 res_smooth = md.result 

540 assert md._fitters[0].fwhm_smooth == 0.8 

541 

542 for shear in ["noshear", "1p", "1m", "2p", "2m"]: 

543 msk = res[shear][model + "_flags"] == 0 

544 msk_smooth = res_smooth[shear][model + "_flags"] == 0 

545 assert ( 

546 np.mean(res_smooth[shear][model + "_T"][msk_smooth]) 

547 > np.mean(res[shear][model + "_T"][msk]) 

548 ) 

549 assert ( 

550 np.mean(res_smooth[shear][model + "_g_cov"][msk_smooth, 0, 0]) 

551 < np.mean(res[shear][model + "_g_cov"][msk, 0, 0]) 

552 ) 

553 

554 

555@pytest.mark.parametrize("model", ["pgauss", "ksigma"]) 

556def test_metadetect_fitter_fwhm_reg(model): 

557 pytest.importorskip("sxdes") 

558 

559 nband = 3 

560 rng = np.random.RandomState(seed=116) 

561 

562 sim = Sim(rng, config={"nband": nband}) 

563 mbobs = sim.get_mbobs() 

564 

565 config = {} 

566 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

567 config["model"] = model 

568 config["weight"]["fwhm"] = 1.2 

569 

570 rng = np.random.RandomState(seed=116) 

571 md = metadetect.Metadetect( 

572 config, mbobs, rng, 

573 ) 

574 md.go() 

575 res = md.result 

576 

577 config["weight"]["fwhm_reg"] = 0.8 

578 rng = np.random.RandomState(seed=116) 

579 md = metadetect.Metadetect( 

580 config, mbobs, rng, 

581 ) 

582 md.go() 

583 res_reg = md.result 

584 

585 for shear in ["noshear", "1p", "1m", "2p", "2m"]: 

586 msk = res[shear][model + "_flags"] == 0 

587 msk_reg = res_reg[shear][model + "_reg0.80_" + "flags"] == 0 

588 assert np.allclose( 

589 res_reg[shear][model + "_reg0.80" + "_T"][msk_reg], 

590 res[shear][model + "_T"][msk] 

591 ) 

592 assert not np.allclose( 

593 res_reg[shear][model + "_reg0.80" + "_g"][msk_reg], 

594 res[shear][model + "_g"][msk] 

595 ) 

596 

597 

598def test_metadetect_fitter_multi_meas(): 

599 pytest.importorskip("sxdes") 

600 

601 nband = 3 

602 rng = np.random.RandomState(seed=116) 

603 

604 sim = Sim(rng, config={"nband": nband}) 

605 mbobs = sim.get_mbobs() 

606 

607 config = {} 

608 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

609 del config["model"] 

610 del config["weight"] 

611 config["fitters"] = [ 

612 {"model": "wmom", "weight": {"fwhm": 1.2}}, 

613 {"model": "pgauss", "weight": {"fwhm": 2.0}}, 

614 {"model": "pgauss", "weight": {"fwhm": 2.0, "fwhm_reg": 0.8}}, 

615 {"model": "am"}, 

616 {"model": "gauss"}, 

617 ] 

618 

619 rng = np.random.RandomState(seed=116) 

620 md = metadetect.Metadetect( 

621 config, mbobs, rng, 

622 ) 

623 md.go() 

624 res = md.result 

625 

626 model = "pgauss" 

627 for shear in ["noshear", "1p", "1m", "2p", "2m"]: 

628 msk = res[shear][model + "_flags"] == 0 

629 msk_reg = res[shear][model + "_reg0.80_" + "flags"] == 0 

630 msk_wmom = res[shear]["wmom_flags"] == 0 

631 msk_admom = res[shear]["am_flags"] == 0 

632 msk_gauss = res[shear]["gauss_flags"] == 0 

633 assert np.allclose( 

634 res[shear][model + "_reg0.80" + "_T"][msk_reg], 

635 res[shear][model + "_T"][msk] 

636 ) 

637 assert not np.allclose( 

638 res[shear][model + "_reg0.80" + "_g"][msk_reg], 

639 res[shear][model + "_g"][msk] 

640 ) 

641 assert not np.allclose( 

642 res[shear]["wmom_T"][msk_wmom], 

643 res[shear][model + "_T"][msk] 

644 ) 

645 assert not np.allclose( 

646 res[shear]["wmom_g"][msk_wmom], 

647 res[shear][model + "_g"][msk] 

648 ) 

649 

650 # admom can fail so look at intersection 

651 assert not np.allclose( 

652 res[shear]["am_T"][msk_admom & msk], 

653 res[shear][model + "_T"][msk_admom & msk] 

654 ) 

655 assert not np.allclose( 

656 res[shear]["am_g"][msk_admom & msk], 

657 res[shear][model + "_g"][msk_admom & msk] 

658 ) 

659 

660 # gauss can fail so look at intersection 

661 assert not np.allclose( 

662 res[shear]["gauss_T"][msk_gauss & msk], 

663 res[shear][model + "_T"][msk_gauss & msk] 

664 ) 

665 assert not np.allclose( 

666 res[shear]["gauss_g"][msk_gauss & msk], 

667 res[shear][model + "_g"][msk_gauss & msk] 

668 ) 

669 

670 

671@pytest.mark.parametrize("model", ["wmom", "pgauss", "ksigma"]) 

672@pytest.mark.parametrize("nband,nshear", [(3, 2), (1, 1), (4, 2), (3, 1)]) 

673def test_metadetect_flux(model, nband, nshear): 

674 """ 

675 test full metadetection w/ fluxes 

676 """ 

677 pytest.importorskip("sxdes") 

678 

679 ntrial = 1 

680 rng = np.random.RandomState(seed=116) 

681 

682 tm0 = time.time() 

683 

684 sim = Sim(rng, config={"nband": nband}) 

685 config = {} 

686 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

687 config["model"] = model 

688 

689 for trial in range(ntrial): 

690 print("trial: %d/%d" % (trial+1, ntrial)) 

691 

692 mbobs = sim.get_mbobs() 

693 for shear_bands in itertools.combinations(list(range(nband)), nshear): 

694 res = metadetect.do_metadetect( 

695 config, mbobs, rng, shear_band_combs=[shear_bands], 

696 det_band_combs="shear_bands", 

697 ) 

698 for shear in ["noshear", "1p", "1m", "2p", "2m"]: 

699 assert np.all(res[shear]["mfrac"] == 0) 

700 assert np.all( 

701 res[shear]["shear_bands"] == "".join("%s" % b for b in shear_bands) 

702 ) 

703 assert np.all( 

704 res[shear]["det_bands"] == "".join("%s" % b for b in shear_bands) 

705 ) 

706 for c in res[shear].dtype.names: 

707 if c.endswith("band_flux"): 

708 if nband > 1: 

709 assert res[shear][c][0].shape == (nband,) 

710 else: 

711 assert res[shear][c][0].shape == tuple() 

712 

713 total_time = time.time()-tm0 

714 print("time per:", total_time/ntrial) 

715 

716 

717@pytest.mark.parametrize("coadd", [True, False]) 

718@pytest.mark.parametrize("det_bands", [None, "shear_bands", "single"]) 

719@pytest.mark.parametrize("model", ["wmom", "pgauss", "am", "gauss"]) 

720def test_metadetect_multiband(model, det_bands, coadd): 

721 """ 

722 test full metadetection w/ multiple bands 

723 """ 

724 pytest.importorskip("sxdes") 

725 

726 nband = 3 

727 ntrial = 1 

728 rng = np.random.RandomState(seed=116) 

729 

730 tm0 = time.time() 

731 

732 sim = Sim(rng, config={"nband": nband}) 

733 config = {} 

734 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

735 config["model"] = model 

736 config["coadd"] = coadd 

737 

738 for trial in range(ntrial): 

739 print("trial: %d/%d" % (trial+1, ntrial)) 

740 

741 mbobs = sim.get_mbobs() 

742 shear_band_combs = [list(range(nband))] 

743 shear_band_combs += [ 

744 list(shear_bands) 

745 for shear_bands in itertools.combinations(list(range(nband)), 2) 

746 ] 

747 shear_band_combs += [ 

748 list(shear_bands) 

749 for shear_bands in itertools.combinations(list(range(nband)), 1) 

750 ] 

751 det_band_combs = ( 

752 det_bands 

753 if det_bands != "single" 

754 else [[0]] * len(shear_band_combs) 

755 ) 

756 res = metadetect.do_metadetect( 

757 config, mbobs, rng, shear_band_combs=shear_band_combs, 

758 det_band_combs=det_band_combs, 

759 ) 

760 if det_band_combs is None: 

761 det_band_combs = [list(range(nband))] * len(shear_band_combs) 

762 elif det_band_combs == "shear_bands": 

763 det_band_combs = shear_band_combs 

764 

765 for shear in ["noshear", "1p", "1m", "2p", "2m"]: 

766 assert np.all(res[shear]["mfrac"] == 0) 

767 for det_bands, shear_bands in zip(det_band_combs, shear_band_combs): 

768 assert np.any( 

769 ( 

770 res[shear]["shear_bands"] 

771 == "".join("%s" % b for b in shear_bands) 

772 ) 

773 & ( 

774 res[shear]["det_bands"] 

775 == "".join("%s" % b for b in det_bands) 

776 ) 

777 ) 

778 for c in res[shear].dtype.names: 

779 if c.endswith("band_flux"): 

780 if nband > 1: 

781 assert res[shear][c][0].shape == (nband,) 

782 else: 

783 assert res[shear][c][0].shape == tuple() 

784 

785 total_time = time.time()-tm0 

786 print("time per:", total_time/ntrial) 

787 

788 

789def test_metadetect_with_color_is_same(): 

790 pytest.importorskip("sxdes") 

791 

792 model = "wmom" 

793 nband = 3 

794 ntrial = 1 

795 

796 tm0 = time.time() 

797 

798 config = {} 

799 config.update(copy.deepcopy(TEST_METADETECT_CONFIG)) 

800 config["model"] = model 

801 

802 for trial in range(ntrial): 

803 print("trial: %d/%d" % (trial+1, ntrial)) 

804 

805 shear_band_combs = [list(range(nband))] 

806 shear_band_combs += [ 

807 list(shear_bands) 

808 for shear_bands in itertools.combinations(list(range(nband)), 2) 

809 ] 

810 shear_band_combs += [ 

811 list(shear_bands) 

812 for shear_bands in itertools.combinations(list(range(nband)), 1) 

813 ] 

814 

815 rng = np.random.RandomState(seed=116) 

816 sim = Sim(rng, config={"nband": nband}) 

817 mbobs = sim.get_mbobs() 

818 rng = np.random.RandomState(seed=11) 

819 res = metadetect.do_metadetect( 

820 config, mbobs, rng, shear_band_combs=shear_band_combs, 

821 ) 

822 

823 rng = np.random.RandomState(seed=116) 

824 sim = Sim(rng, config={"nband": nband}) 

825 mbobs = sim.get_mbobs() 

826 rng = np.random.RandomState(seed=11) 

827 res_color = metadetect.do_metadetect( 

828 config, mbobs, rng, shear_band_combs=shear_band_combs, 

829 color_key_func=lambda x: "blah", color_dep_mbobs={"blah": mbobs}, 

830 ) 

831 for shear in ["noshear", "1p", "1m", "2p", "2m"]: 

832 for col in res[shear].dtype.names: 

833 assert col in res_color[shear].dtype.names 

834 if col == "shear_bands" or col == "det_bands": 

835 assert np.array_equal( 

836 res[shear][col], 

837 res_color[shear][col], 

838 ) 

839 else: 

840 np.testing.assert_allclose( 

841 res[shear][col], 

842 res_color[shear][col], 

843 atol=0, 

844 rtol=0, 

845 equal_nan=True, 

846 ) 

847 

848 for shear_bands in shear_band_combs: 

849 assert np.any( 

850 res[shear]["shear_bands"] == "".join("%s" % b for b in shear_bands) 

851 ) 

852 assert np.any( 

853 res_color[shear]["shear_bands"] 

854 == "".join("%s" % b for b in shear_bands) 

855 ) 

856 

857 total_time = time.time()-tm0 

858 print("time per:", total_time/ntrial) 

859 

860 

861@pytest.mark.parametrize("mask_region", [1, 7]) 

862def test_fill_in_mask_col(mask_region): 

863 rng = np.random.RandomState(seed=10) 

864 

865 rows = np.array([31.2]) 

866 cols = np.array([51.7]) 

867 mask = rng.randint(low=0, high=64, size=(100, 100)) 

868 

869 vals = metadetect._fill_in_mask_col( 

870 mask_region=mask_region, 

871 rows=rows, 

872 cols=cols, 

873 mask=mask) 

874 

875 row = 31 

876 col = 52 

877 if mask_region == 1: 

878 assert vals[0] == mask[row, col] 

879 else: 

880 assert vals[0] == np.bitwise_or.reduce( 

881 mask[ 

882 row-mask_region:row+mask_region+1, 

883 col-mask_region:col+mask_region+1 

884 ] 

885 )[0] 

886 

887 

888def test_get_psf_stats(): 

889 rng = np.random.RandomState(seed=10) 

890 sim = Sim(rng) 

891 mbobs = sim.get_mbobs() 

892 fitting.fit_all_psfs(mbobs, rng) 

893 

894 psf_stats = metadetect._get_psf_stats(mbobs, 0) 

895 assert psf_stats["flags"] == 0 

896 assert np.isfinite(psf_stats["g1"]) 

897 assert np.isfinite(psf_stats["g2"]) 

898 assert np.isfinite(psf_stats["T"]) 

899 

900 psf_stats = metadetect._get_psf_stats(mbobs, 2) 

901 assert psf_stats["flags"] == (procflags.PSF_FAILURE | 2) 

902 assert not np.isfinite(psf_stats["g1"]) 

903 assert not np.isfinite(psf_stats["g2"]) 

904 assert not np.isfinite(psf_stats["T"]) 

905 

906 for obslist in mbobs: 

907 for obs in obslist: 

908 obs.weight = -1.0*obs.weight 

909 psf_stats = metadetect._get_psf_stats(mbobs, 0) 

910 assert psf_stats["flags"] == procflags.PSF_FAILURE 

911 assert not np.isfinite(psf_stats["g1"]) 

912 assert not np.isfinite(psf_stats["g2"]) 

913 assert not np.isfinite(psf_stats["T"]) 

914 

915 e1s = np.arange(len(mbobs)) + 0.1 

916 e2s = 2*np.arange(len(mbobs)) + 0.1 

917 Ts = 3*np.arange(len(mbobs)) + 0.1 

918 wgts = np.arange(len(mbobs)) + 1 

919 for i, obslist in enumerate(mbobs): 

920 for obs in obslist: 

921 obs.weight = 0*obs.weight + wgts[i] 

922 obs.psf.meta["result"]["e"] = (e1s[i], e2s[i]) 

923 obs.psf.meta["result"]["T"] = Ts[i] 

924 

925 psf_stats = metadetect._get_psf_stats(mbobs, 0) 

926 assert psf_stats["flags"] == 0 

927 assert np.allclose(psf_stats["g1"], np.sum(wgts * e1s)/np.sum(wgts)) 

928 assert np.allclose(psf_stats["g2"], np.sum(wgts * e2s)/np.sum(wgts)) 

929 assert np.allclose(psf_stats["T"], np.sum(wgts * Ts)/np.sum(wgts))