Coverage for metadetect / lsst / util.py: 11%

309 statements  

« prev     ^ index     » next       coverage.py v7.12.0, created at 2025-12-17 16:58 -0800

1import logging 

2from contextlib import contextmanager, ExitStack 

3import ngmix 

4import numpy as np 

5 

6logger = logging.getLogger(__name__) 

7 

8 

9class ContextNoiseReplacer(object): 

10 """ 

11 noise replacer that works as a context manager 

12 

13 Parameters 

14 ---------- 

15 mbexp: lsst.afw.image.Exposure 

16 The data 

17 sources: lsst.afw.table.SourceCatalog 

18 Catalog of sources 

19 noise_image: array, optional 

20 Optional noise image to use. If not sent one is generated. 

21 

22 Examples 

23 -------- 

24 with ContextNoiseReplacer(exposure=exp, sources=sources) as replacer: 

25 # do something 

26 """ 

27 

28 def __init__(self, exposure, sources, rng, config=None, noise_image=None): 

29 from lsst.meas.base import NoiseReplacerConfig, NoiseReplacer 

30 import lsst.afw.detection as afw_det 

31 

32 # Notes for metacal. 

33 # 

34 # For metacal we should generate a noise image so that the exact noise 

35 # field is used for all versions of the metacal images. The assumption is 

36 # that, because these noise data should contain no signal, metacal is not 

37 # calibrating it. Thus it doesn't matter whether or not the noise field is 

38 # representative of the full covariance of the true image noise. Rather by 

39 # making the field the same for all metacal images we reduce variance in 

40 # the calculation of the response 

41 

42 if config is None: 

43 config = NoiseReplacerConfig() 

44 

45 # TODO DM needs to fix the crash 

46 # config.noiseSource = 'variance' 

47 config.noiseSeedMultiplier = rng.randint(0, 2**24) 

48 

49 if noise_image is None: 

50 # TODO remove_poisson should be true for real data 

51 noise_image = get_noise_image(exposure, rng=rng, remove_poisson=False) 

52 

53 footprints = {} 

54 for source in sources: 

55 parent_id = source.getParent() 

56 fp = source.getFootprint() 

57 heavy = afw_det.makeHeavyFootprint(fp, exposure.maskedImage) 

58 footprints[source.getId()] = (parent_id, heavy) 

59 # footprints = { 

60 # source.getId(): (source.getParent(), source.getFootprint()) 

61 # for source in sources 

62 # } 

63 

64 # This constructor will replace all detected pixels with noise in the 

65 # image 

66 self.replacer = NoiseReplacer( 

67 config, 

68 exposure=exposure, 

69 footprints=footprints, 

70 noiseImage=noise_image, 

71 ) 

72 

73 @contextmanager 

74 def sourceInserted(self, source_id): 

75 """ 

76 Context manager to insert a source 

77 

78 with replacer.sourceInserted(source_id): 

79 # do something with exposure 

80 """ 

81 self.insertSource(source_id) 

82 

83 try: 

84 yield self 

85 finally: 

86 self.removeSource(source_id) 

87 

88 def insertSource(self, source_id): 

89 """ 

90 Insert a source 

91 """ 

92 self.replacer.insertSource(source_id) 

93 

94 def removeSource(self, source_id): 

95 """ 

96 Remove a source 

97 """ 

98 self.replacer.removeSource(source_id) 

99 

100 def end(self): 

101 self.replacer.end() 

102 

103 def __enter__(self): 

104 return self 

105 

106 def __exit__(self, exception_type, exception_value, traceback): 

107 self.end() 

108 

109 

110class MultibandNoiseReplacer(object): 

111 """ 

112 noise replacer that works on multiple bands 

113 

114 Parameters 

115 ---------- 

116 mbexp: lsst.afw.image.MultibandExposure 

117 The data 

118 sources: lsst.afw.table.SourceCatalog 

119 Catalog of sources 

120 

121 Examples 

122 -------- 

123 with MultibandNoiseReplacer(mbexp=mbexp, sources=sources) as replacer: 

124 with replacer.insertSource(source_id): 

125 # do something with exposures 

126 """ 

127 def __init__(self, mbexp, sources, rng): 

128 self.mbexp = mbexp 

129 self.sources = sources 

130 self.rng = rng 

131 self._set_noise_replacers() 

132 

133 @contextmanager 

134 def sourceInserted(self, source_id): 

135 """ 

136 Context manager to insert a source 

137 

138 with replacer.sourceInserted(source_id): 

139 # do something with exposures 

140 """ 

141 self.insertSource(source_id) 

142 

143 try: 

144 # usually won't use this yielded value 

145 yield self.mbexp 

146 finally: 

147 self.removeSource(source_id) 

148 

149 def insertSource(self, source_id): 

150 """ 

151 Insert a source 

152 """ 

153 for replacer in self.noise_replacers: 

154 replacer.insertSource(source_id) 

155 

156 def removeSource(self, source_id): 

157 """ 

158 Remove a source 

159 """ 

160 for replacer in self.noise_replacers: 

161 replacer.removeSource(source_id) 

162 

163 def end(self): 

164 """ 

165 end the noise replacment. Called automatically upon leaving the 

166 context manager 

167 """ 

168 self.exit_stack.close() 

169 

170 def _set_noise_replacers(self): 

171 """ 

172 set a noise replacer for each exp and put it on the 

173 context lib ExitStack for cleanup later 

174 """ 

175 

176 self.noise_replacers = [] 

177 self.exit_stack = ExitStack() 

178 

179 for exp in self.mbexp.singles: 

180 replacer = ContextNoiseReplacer( 

181 exposure=exp, 

182 sources=self.sources, 

183 rng=self.rng, 

184 ) 

185 self.noise_replacers.append(replacer) 

186 self.exit_stack.enter_context(replacer) 

187 

188 def __enter__(self): 

189 return self 

190 

191 def __exit__(self, exception_type, exception_value, traceback): 

192 self.end() 

193 

194 

195def get_noise_replacer(exposure, sources, noise_image=None): 

196 """ 

197 get a regular noise replacer for the input exposure and source list 

198 

199 Parameters 

200 ---------- 

201 exposure: lsst.afw.image.ExposureF 

202 The exposure data 

203 sources: lsst.afw.table.SourceCatalog 

204 The sources 

205 noise_image: array, optional 

206 Optional array of noise 

207 

208 Returns 

209 ------- 

210 lsst.meas.base.NoiseReplacer 

211 """ 

212 from lsst.meas.base import NoiseReplacerConfig, NoiseReplacer 

213 

214 # Notes for metacal. 

215 # 

216 # For metacal we should generate a noise image so that the exact noise 

217 # field is used for all versions of the metacal images. The assumption is 

218 # that, because these noise data should contain no signal, metacal is not 

219 # calibrating it. Thus it doesn't matter whether or not the noise field is 

220 # representative of the full covariance of the true image noise. Rather by 

221 # making the field the same for all metacal images we reduce variance in 

222 # the calculation of the response 

223 

224 noise_replacer_config = NoiseReplacerConfig() 

225 footprints = { 

226 source.getId(): (source.getParent(), source.getFootprint()) 

227 for source in sources 

228 } 

229 

230 # This constructor will replace all detected pixels with noise in the 

231 # image 

232 return NoiseReplacer( 

233 noise_replacer_config, 

234 exposure=exposure, 

235 footprints=footprints, 

236 noiseImage=noise_image, 

237 ) 

238 

239 

240def get_noise_image(exp, rng, remove_poisson): 

241 """ 

242 get a noise image based on the input exposure 

243 

244 TODO gain correct separately in each amplifier, currently 

245 averaged 

246 

247 Parameters 

248 ---------- 

249 exp: afw.image.ExposureF 

250 The exposure upon which to base the noise 

251 rng: np.random.RandomState 

252 The random number generator for making the noise image 

253 remove_poisson: bool 

254 If True, remove the poisson noise from the variance 

255 estimate. 

256 

257 Returns 

258 ------- 

259 MaskedImage 

260 """ 

261 import lsst.afw.image as afw_image 

262 

263 noise_exp = afw_image.ExposureF(exp, deep=True) 

264 

265 signal = exp.image.array 

266 variance = exp.variance.array 

267 

268 use = np.where(np.isfinite(variance) & np.isfinite(signal)) 

269 

270 if remove_poisson: 

271 gains = [ 

272 amp.getGain() for amp in exp.getDetector().getAmplifiers() 

273 ] 

274 mean_gain = np.mean(gains) 

275 

276 corrected_var = variance[use] - signal[use] / mean_gain 

277 

278 var = np.median(corrected_var) 

279 else: 

280 var = np.median(variance[use]) 

281 

282 noise = rng.normal(scale=np.sqrt(var), size=signal.shape) 

283 

284 noise_exp.image.array[:, :] = noise 

285 noise_exp.variance.array[:, :] = var 

286 

287 return noise_exp.getMaskedImage() 

288 

289 

290def get_mbexp(exposures): 

291 """ 

292 convert a list of exposures into an MultibandExposure 

293 

294 Parameters 

295 ---------- 

296 exposures: [lsst.afw.image.ExposureF] 

297 List of exposures of type F, D etc. 

298 

299 Returns 

300 ------- 

301 lsst.afw.image.MultibandExposure 

302 """ 

303 from lsst.afw.image import MultibandExposure 

304 

305 filters = [exp.getFilter().bandLabel for exp in exposures] 

306 mbexp = MultibandExposure.fromExposures(filters, exposures) 

307 

308 for exp, sexp in zip(exposures, mbexp.singles): 

309 sexp.setFilter(exp.getFilter()) 

310 sexp.setWcs(exp.getWcs()) 

311 

312 return mbexp 

313 

314 

315def copy_mbexp(mbexp, clear=False): 

316 """ 

317 copy a MultibandExposure with psfs 

318 

319 clone() does not copy the psfs 

320 

321 Parameters 

322 ---------- 

323 mbexp: lsst.afw.image.MultibandExposure 

324 The exposures to copy 

325 clear: bool, optional 

326 If set to True, clear the data 

327 

328 Returns 

329 ------- 

330 lsst.afw.image.MultibandExposure 

331 """ 

332 

333 new_mbexp = mbexp.clone() 

334 

335 # clone does not copy the psfs 

336 for band in mbexp.bands: 

337 psf = try_clone_psf(mbexp[band].getPsf()) 

338 new_mbexp[band].setPsf(psf) 

339 

340 for new_exp, exp in zip(new_mbexp.singles, mbexp.singles): 

341 new_exp.setFilter(exp.getFilter()) 

342 new_exp.setWcs(exp.getWcs()) 

343 

344 if clear: 

345 new_mbexp.image.array[:, :, :] = 0 

346 

347 return new_mbexp 

348 

349 

350def try_clone_psf(psf): 

351 """ 

352 try to clone it, if not return a reference 

353 """ 

354 try: 

355 new_psf = psf.clone() 

356 except RuntimeError: 

357 # this happens on some psfs, not sure why 

358 # proceed without a copy 

359 new_psf = psf 

360 

361 return new_psf 

362 

363 

364def coadd_exposures(exposures): 

365 """ 

366 coadd a set of exposures, assuming they share the same wcs 

367 

368 Parameters 

369 ---------- 

370 exposures: [lsst.afw.image.Exposure] 

371 List of exposures to coadd 

372 

373 Returns 

374 -------- 

375 lsst.afw.image.ExposureF 

376 """ 

377 import lsst.geom as geom 

378 import lsst.afw.image as afw_image 

379 

380 wsum = 0.0 

381 

382 for i, exp in enumerate(exposures): 

383 

384 shape = exp.image.array.shape 

385 

386 ycen, xcen = (np.array(shape) - 1)/2 

387 cen = geom.Point2D(xcen, ycen) 

388 

389 psfobj = exp.getPsf() 

390 this_psfim = psfobj.computeKernelImage(cen).array 

391 

392 if i == 0: 

393 coadd_exp = afw_image.ExposureF(exp, deep=True) 

394 coadd_exp.image.array[:, :] = 0.0 

395 

396 weight = np.zeros(shape, dtype='f4') 

397 

398 psf_im = this_psfim.copy() 

399 psf_im[:, :] = 0 

400 

401 coadd_exp.mask.array |= exp.mask.array 

402 

403 w = np.where(exp.variance.array > 0) 

404 medvar = np.median(exp.variance.array[w]) 

405 this_weight = 1.0/medvar 

406 

407 coadd_exp.image.array[w] += exp.image.array[w] * this_weight 

408 psf_im += this_psfim * this_weight 

409 

410 weight[w] += 1.0/exp.variance.array[w] 

411 

412 wsum += this_weight 

413 

414 if wsum <= 0: 

415 logger.info('found wsum <= 0') 

416 return None 

417 

418 fac = 1.0/wsum 

419 

420 coadd_exp.image.array[:, :] *= fac 

421 

422 # the psf is always normalized 

423 psf_im *= 1.0/psf_im.sum() 

424 

425 coadd_exp.variance.array[:, :] = np.inf 

426 w = np.where(weight > 0) 

427 coadd_exp.variance.array[w] = 1/weight[w] 

428 

429 coadd_psf = get_stack_kernel_psf(psf_im) 

430 coadd_exp.setPsf(coadd_psf) 

431 

432 coadd_exp.setWcs(exposures[0].getWcs()) 

433 

434 return coadd_exp 

435 

436 

437def get_stack_kernel_psf(psf_image): 

438 """ 

439 create a KernelPsf from the input image 

440 

441 Parameters 

442 ---------- 

443 psf_image: array 

444 The array data 

445 

446 Returns 

447 ------- 

448 lsst.meas.algorithms.KernelPsf 

449 """ 

450 import lsst.afw.image as afw_image 

451 from lsst.meas.algorithms import KernelPsf 

452 from lsst.afw.math import FixedKernel 

453 

454 return KernelPsf( 

455 FixedKernel( 

456 afw_image.ImageD(psf_image.astype(float)) 

457 ) 

458 ) 

459 

460 

461def coadd_mbobs(mbobs): 

462 """ 

463 coadd a set of exposures, assuming they share the same wcs 

464 

465 Parameters 

466 ---------- 

467 exposures: [lsst.afw.image.Exposure] 

468 List of exposures to coadd 

469 

470 Returns 

471 -------- 

472 lsst.afw.image.ExposureF 

473 """ 

474 

475 wsum = 0.0 

476 

477 meta = {} 

478 for i, obslist in enumerate(mbobs): 

479 obs = obslist[0] 

480 meta.update(obs.meta) 

481 

482 shape = obs.image.shape 

483 

484 if i == 0: 

485 coadd_image = np.zeros(shape, dtype='f4') 

486 coadd_weight = np.zeros(shape, dtype='f4') 

487 

488 coadd_psf = np.zeros(obs.psf.image.shape) 

489 

490 w = np.where(obs.weight > 0) 

491 this_weight = np.median(obs.weight[w]) 

492 

493 coadd_image[:, :] += obs.image * this_weight 

494 coadd_weight[:, :] += obs.weight 

495 coadd_psf[:, :] += obs.psf.image * this_weight 

496 

497 wsum += this_weight 

498 

499 fac = 1.0/wsum 

500 

501 coadd_image[:, :] *= fac 

502 

503 # the psf is always normalized 

504 coadd_psf *= 1.0/coadd_psf.sum() 

505 

506 # use the jacobians from the last obs 

507 jac = obs.jacobian 

508 psf_jac = obs.psf.jacobian 

509 

510 fake_noise = coadd_psf.max() / 1.e6 

511 psf_weight = np.zeros(coadd_psf.shape) 

512 psf_weight += (1.0/fake_noise**2) 

513 

514 psf_obs = ngmix.Observation( 

515 image=coadd_psf, 

516 weight=psf_weight, 

517 jacobian=psf_jac, 

518 ) 

519 output_obs = ngmix.Observation( 

520 image=coadd_image, 

521 weight=coadd_weight, 

522 jacobian=jac, 

523 psf=psf_obs, 

524 meta=meta, 

525 ) 

526 

527 return output_obs 

528 

529 

530def trim_odd_image(im): 

531 """ 

532 trim an odd dimension image to by square and with equal distance from 

533 canonical center to all edges 

534 """ 

535 

536 dims = im.shape 

537 if dims[0] != dims[1]: 

538 logger.debug('original dims: %s' % str(dims)) 

539 assert dims[0] % 2 != 0, 'image must have odd dims' 

540 assert dims[1] % 2 != 0, 'image must have odd dims' 

541 

542 dims = np.array(dims) 

543 cen = (dims-1)//2 

544 cen = cen.astype('i4') 

545 

546 distances = ( 

547 cen[0]-0, 

548 dims[0]-cen[0]-1, 

549 cen[1]-0, 

550 dims[1]-cen[1]-1, 

551 ) 

552 logger.debug('distances: %s' % str(distances)) 

553 min_dist = min(distances) 

554 

555 start_row = cen[0] - min_dist 

556 end_row = cen[0] + min_dist 

557 start_col = cen[1] - min_dist 

558 end_col = cen[1] + min_dist 

559 

560 # adding +1 for slices 

561 new_im = im[ 

562 start_row:end_row+1, 

563 start_col:end_col+1, 

564 ].copy() 

565 

566 logger.debug('new dims: %s' % str(new_im.shape)) 

567 

568 else: 

569 new_im = im 

570 

571 return new_im 

572 

573 

574def obs2exp(obs, exp=None, copy_mask_from='bmask', copy_psf=False): 

575 """ 

576 copy ngmix.Observation data into an Exposure. The pixel data are copied, 

577 and the weight is converted to variance. 

578 

579 If a complex WCS or PSF is needed, it is better to send an existing one. 

580 The wcs will not be altered by this code. 

581 

582 The psf image can be optionally copied in as a KernelPsf 

583 

584 Parameters 

585 ---------- 

586 obs: ngmix.Observation 

587 The observation data 

588 exp: lsst.afw.image.ExposureF, optional 

589 Optional exposure into which the data are to copied. Note the psf is 

590 not copied unless copy_psf is set to True. The wcs is not modified. 

591 copy_mask_from: bool, optional 

592 Copy from this mask into the exp mask. 

593 copy_psf: bool, optional 

594 If set to True, the psf is copied in as a KernelPsf(FixedKernel) 

595 

596 Returns 

597 ------- 

598 exp: lsst.afw.image.ExposureF 

599 Either a reference to the input exposure, or a new one if it was 

600 created internally. 

601 """ 

602 from lsst.meas.algorithms import KernelPsf 

603 from lsst.afw.math import FixedKernel 

604 import lsst.afw.image as afw_image 

605 

606 if exp is None: 

607 ny, nx = obs.image.shape 

608 exp = afw_image.ExposureF(nx, ny) 

609 

610 exp.image.array[:, :] = obs.image 

611 

612 exp.variance.array[:, :] = np.inf 

613 

614 if copy_mask_from == 'bmask': 

615 mask = obs.bmask 

616 elif copy_mask_from == 'ormask': 

617 mask = obs.ormask 

618 else: 

619 raise RuntimeError( 

620 f"copy_mask_from '{copy_mask_from}' should be 'bmask' or 'ormask'" 

621 ) 

622 exp.mask.array[:, :] = mask 

623 

624 w = np.where(obs.weight > 0) 

625 exp.variance.array[w] = 1.0/obs.weight[w] 

626 

627 if copy_psf: 

628 psf_image = obs.psf.image 

629 stack_psf = KernelPsf( 

630 FixedKernel( 

631 afw_image.ImageD(psf_image.astype(float)) 

632 ) 

633 ) 

634 exp.setPsf(stack_psf) 

635 

636 return exp 

637 

638 

639def exp2obs(exp, copy_mask_to='ormask', store_exp=False): 

640 """ 

641 convert an exposure to an observation. The wcs is linearized 

642 and the psf is reconstructed as a kernel image. 

643 

644 The origin for the jacobian is set to the image center, relative to the 

645 bounding box corner. The psf image is reconstructed at the center of the 

646 image. 

647 

648 TODO do we ever need to do computeImage in this code? 

649 

650 Optionally, the original exposure can be stored in .meta so we can reverse 

651 the process and go back to an exposure, copying in any changes to the image 

652 data 

653 

654 Parameters 

655 ---------- 

656 exp: lsst.afw.image.ExposureF 

657 The exposure 

658 copy_mask_to: str 

659 Copy from the exp.mask into this attribute of the Observation. For 

660 example, if set to 'ormask', the exposure mask is copied to the 

661 obs.ormask and obs.bmask is set to zero. 

662 store_exp: bool, optional 

663 If set to True, store the original exposure in .meta['exposure'] 

664 

665 Returns 

666 ------- 

667 ngmix.Observation 

668 """ 

669 

670 wcs = exp.getWcs() 

671 bbox = exp.getBBox() 

672 

673 # this is to be consistent with the coadd code 

674 cen, _ = get_integer_center( 

675 wcs=wcs, 

676 bbox=bbox, 

677 as_double=True, 

678 ) 

679 

680 jac = get_jacobian(exp=exp, cen=cen) 

681 

682 dims = exp.image.array.shape 

683 weight = np.zeros(dims) 

684 w = np.where(exp.variance.array > 0) 

685 weight[w] = 1.0/exp.variance.array[w] 

686 

687 psf_obj = exp.getPsf() 

688 psf_image = psf_obj.computeKernelImage(cen).array 

689 

690 assert psf_image.shape[0] == psf_image.shape[1], 'psf is not square' 

691 assert psf_image.shape[0] % 2 != 0, 'psf dims are not odd' 

692 psf_cen = (np.array(psf_image.shape)-1.0)/2.0 

693 

694 psf_jac = jac.copy() 

695 psf_jac.set_cen(row=psf_cen[0], col=psf_cen[1]) 

696 

697 psf_err = psf_image.max()*0.0001 

698 psf_weight = psf_image*0 + 1.0/psf_err**2 

699 

700 psf_obs = ngmix.Observation( 

701 image=psf_image, 

702 weight=psf_weight, 

703 jacobian=psf_jac, 

704 ) 

705 

706 if copy_mask_to == 'ormask': 

707 ormask = exp.mask.array 

708 bmask = ormask * 0 

709 elif copy_mask_to == 'bmask': 

710 bmask = exp.mask.array 

711 ormask = bmask * 0 

712 else: 

713 raise RuntimeError( 

714 f"copy_mask_to '{copy_mask_to}' should be 'bmask' or 'ormask'" 

715 ) 

716 

717 meta = {} 

718 if store_exp: 

719 meta['exposure'] = exp 

720 

721 return ngmix.Observation( 

722 image=exp.image.array, 

723 weight=weight, 

724 bmask=bmask, 

725 ormask=ormask, 

726 jacobian=jac, 

727 psf=psf_obs, 

728 meta=meta, 

729 ) 

730 

731 

732def get_jacobian(exp, cen): 

733 """ 

734 get an ngmix jacobian at the specified location 

735 

736 Parameters 

737 ---------- 

738 exp: lsst.afw.image.Exposure 

739 The exposure data 

740 cen: lsst.geom.Point2D 

741 The position at which to get the jacobian. It is also used as the 

742 center for the ngmix jacobian as relative to the corner of the bounding 

743 box. 

744 

745 Returns 

746 ------- 

747 ngmix.Jacobian 

748 """ 

749 import lsst.geom as geom 

750 

751 wcs = exp.getWcs() 

752 xy0 = exp.getXY0() 

753 

754 dm_jac = wcs.linearizePixelToSky(cen, geom.arcseconds) 

755 matrix = dm_jac.getLinear().getMatrix() 

756 

757 relative_cen = cen - geom.Extent2D(xy0) 

758 

759 # ESS reverse engineered this convention mismatch. No documentation 

760 # was found. Don't change this unless you know what you are doing. 

761 

762 return ngmix.Jacobian( 

763 x=relative_cen.getX(), 

764 y=relative_cen.getY(), 

765 dudx=matrix[1, 1], 

766 dudy=-matrix[1, 0], 

767 dvdx=matrix[0, 1], 

768 dvdy=-matrix[0, 0], 

769 ) 

770 

771 

772def get_integer_center(wcs, bbox, as_double=False): 

773 """ 

774 get the pixel and sky center of the coadd within the bbox 

775 The pixel center is forced to be integer 

776 

777 Parameters 

778 ----------- 

779 wcs: DM wcs 

780 The wcs for the exposure 

781 bbox: geom.Box2I 

782 The bounding box for the exposure, possibly within a larger wcs system 

783 

784 Returns 

785 ------- 

786 pixcen, skycen: 

787 pixcen as Point2I, skycen as SpherePoint 

788 """ 

789 from lsst.geom import Point2D, Point2I 

790 

791 # force integer location 

792 pixcen = Point2I(bbox.getCenter()) 

793 skycen = wcs.pixelToSky(Point2D(pixcen)) 

794 

795 if as_double: 

796 pixcen = Point2D(pixcen) 

797 

798 return pixcen, skycen 

799 

800 

801def get_stats_mask(*args, **kwargs): 

802 """ 

803 Get a stats mask for use in getting image statistics. 

804 

805 ['BAD', 'EDGE', 'DETECTED', 'DETECTED_NEGATIVE', 'NO_DATA', 'BRIGHT'] 

806 

807 Returns 

808 ------- 

809 A list of mask planes 

810 """ 

811 

812 # these will be ignored when finding the image standard deviation 

813 # stats_mask = ['BAD', 'SAT', 'EDGE', 'NO_DATA'] 

814 

815 import lsst.afw.image as afw_image 

816 afw_image.Mask.addMaskPlane('BRIGHT') 

817 stats_mask = ['BAD', 'EDGE', 'DETECTED', 'DETECTED_NEGATIVE', 'NO_DATA', 'BRIGHT'] 

818 

819 return stats_mask 

820 

821 

822def get_detection_mask(*args, **kwargs): 

823 """ 

824 Get a mask for detection. Regions with these flags set will not be searched 

825 for objects. 

826 

827 Default is ['EDGE', 'NO_DATA', 'BRIGHT] 

828 

829 BRIGHT is defined in descwl_coadd, so it is not guaranteed to be present in 

830 the global mask bit space. We add BRIGHT to detection mask. If it already 

831 exists, it does not add a new plane and is harmless. 

832 

833 Returns 

834 ------- 

835 A list of mask planes 

836 """ 

837 

838 import lsst.afw.image as afw_image 

839 afw_image.Mask.addMaskPlane('BRIGHT') 

840 mask = ['EDGE', 'NO_DATA', 'BRIGHT'] 

841 

842 return mask 

843 

844 

845def extract_multiband_coadd_data(coadd_data_list): 

846 """ 

847 Convert a list of coadd data for a set of bands into MultibandExposure data 

848 

849 Side Effects 

850 ------------ 

851 The masks are zerod, but the original ormasks are returned in the 'ormask' 

852 item. 

853 

854 Parameters 

855 ---------- 

856 coadd_data_list: list of dict 

857 Each should be the output of descwl_coadd.make_coadd for a given band. 

858 

859 Returns 

860 ------- 

861 data: dict 

862 With items mbexp, noise_mbexp, mfrac_mbexp, ormasks 

863 

864 ormasks is a copy of .mask.array for each exp; the mask.array 

865 are zerod for the exposures 

866 """ 

867 

868 exps = [] 

869 noise_exps = [] 

870 mfrac_exps = [] 

871 ormasks = [] 

872 for band, coadd_data in enumerate(coadd_data_list): 

873 exp = coadd_data['coadd_exp'] 

874 noise_exp = coadd_data['coadd_noise_exp'] 

875 mfrac_exp = coadd_data['coadd_mfrac_exp'] 

876 

877 ormasks.append(exp.mask.array.copy()) 

878 exp.mask.array[:, :] = 0 

879 noise_exp.mask.array[:, :] = 0 

880 mfrac_exp.mask.array[:, :] = 0 

881 

882 exps.append(exp) 

883 noise_exps.append(noise_exp) 

884 mfrac_exps.append(mfrac_exp) 

885 

886 mbexp = get_mbexp(exps) 

887 noise_mbexp = get_mbexp(noise_exps) 

888 mfrac_mbexp = get_mbexp(mfrac_exps) 

889 

890 return { 

891 'mbexp': mbexp, 

892 'noise_mbexp': noise_mbexp, 

893 'mfrac_mbexp': mfrac_mbexp, 

894 'ormasks': ormasks, 

895 } 

896 

897 

898def override_config(config, config_override: dict): 

899 """Override the values in a Config instance with those in a dict. 

900 

901 The `config` object can have arbitrary level of nesting, and to override 

902 the values, `config_override` should have the same nested structure. 

903 To allow for arbitrary nesting, this function calls itself recursively. 

904 

905 Parameters 

906 ---------- 

907 config : `lsst.pex.config.Config` 

908 The configuration object to override. 

909 config_override : `dict` 

910 The dictionary of values to override the configuration with. 

911 

912 Raises 

913 ------ 

914 AttributeError 

915 Raised if a key in `config_override` is not an attribute of `config`. 

916 """ 

917 for key, value in config_override.items(): 

918 if isinstance(value, dict): 

919 subconfig = getattr(config, key) 

920 # Use recursion to allow for arbitrarily nested config objects. 

921 override_config(subconfig, value) 

922 else: 

923 if key == 'stamp_size': 

924 config.toDict()[key] = value 

925 else: 

926 setattr(config, key, value)