Coverage for metadetect / lsst / tests / test_lsst_noise_replacers.py: 15%

99 statements  

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

1import sys 

2import numpy as np 

3import logging 

4import pytest 

5 

6from metadetect.lsst import util 

7from lsst.utils import getPackageDir 

8 

9try: 

10 getPackageDir('descwl_shear_sims') 

11 skip_tests_on_simulations = False 

12except LookupError: 

13 skip_tests_on_simulations = True 

14 

15logging.basicConfig(stream=sys.stdout, level=logging.WARN) 

16 

17 

18def make_lsst_sim(rng): 

19 import descwl_shear_sims 

20 coadd_dim = 200 

21 buff = 20 

22 # coadd_dim = 351 

23 # buff = 50 

24 

25 # the EDGE region is 5 pixels wide but, give a bit more space because the 

26 # sky sub seems to fail with a lsst.pex.exceptions.wrappers.LengthError, 

27 # presumably due to an object near the edge 

28 

29 bands = ['r', 'i', 'z'] 

30 

31 galaxy_catalog = descwl_shear_sims.galaxies.FixedGalaxyCatalog( 

32 rng=rng, 

33 coadd_dim=coadd_dim, 

34 buff=buff, 

35 layout='random', 

36 mag=22, 

37 hlr=0.5, 

38 ) 

39 

40 psf = descwl_shear_sims.psfs.make_fixed_psf(psf_type='gauss') 

41 

42 sim_data = descwl_shear_sims.make_sim( 

43 rng=rng, 

44 galaxy_catalog=galaxy_catalog, 

45 coadd_dim=coadd_dim, 

46 bands=bands, 

47 g1=0, 

48 g2=0, 

49 psf=psf, 

50 se_dim=coadd_dim, 

51 ) 

52 return sim_data 

53 

54 

55def detect_and_deblend(mbexp): 

56 import lsst.afw.table as afw_table 

57 from lsst.meas.base import ( 

58 SingleFrameMeasurementConfig, 

59 SingleFrameMeasurementTask, 

60 ) 

61 from lsst.meas.algorithms import SourceDetectionTask, SourceDetectionConfig 

62 from lsst.meas.deblender import SourceDeblendTask, SourceDeblendConfig 

63 

64 detexp = util.coadd_exposures(mbexp.singles) 

65 

66 schema = afw_table.SourceTable.makeMinimalSchema() 

67 

68 # note we won't run any of these measurements, but it is needed so that 

69 # getCentroid will return the peak position rather than NaN. 

70 # I think it modifies the schema and sets defaults 

71 

72 meas_config = SingleFrameMeasurementConfig() 

73 meas_config.plugins.names = [ 

74 "base_SdssCentroid", 

75 "base_PsfFlux", 

76 "base_SkyCoord", 

77 ] 

78 

79 meas_config.slots.apFlux = None 

80 meas_config.slots.gaussianFlux = None 

81 meas_config.slots.calibFlux = None 

82 meas_config.slots.modelFlux = None 

83 

84 # goes with SdssShape above 

85 meas_config.slots.shape = None 

86 

87 # fix odd issue where it thinks things are near the edge 

88 meas_config.plugins['base_SdssCentroid'].binmax = 1 

89 

90 _ = SingleFrameMeasurementTask( 

91 config=meas_config, 

92 schema=schema, 

93 ) 

94 

95 detection_config = SourceDetectionConfig() 

96 detection_config.reEstimateBackground = False 

97 detection_config.thresholdValue = 5 

98 detection_task = SourceDetectionTask(config=detection_config) 

99 

100 deblend_task = SourceDeblendTask( 

101 config=SourceDeblendConfig(), 

102 schema=schema, 

103 ) 

104 

105 table = afw_table.SourceTable.make(schema) 

106 result = detection_task.run(table, detexp) 

107 

108 if result is not None: 

109 sources = result.sources 

110 deblend_task.run(detexp, sources) 

111 else: 

112 sources = [] 

113 

114 return sources 

115 

116 

117@pytest.mark.skipif( 

118 skip_tests_on_simulations, 

119 reason='descwl_shear_sims not available' 

120) 

121def test_noise_replacer(): 

122 import lsst.afw.image as afw_image 

123 seed = 981 

124 rng = np.random.RandomState(seed) 

125 

126 sim = make_lsst_sim(rng) 

127 

128 exps = [texps[0] for _, texps in sim['band_data'].items()] 

129 mbexp = util.get_mbexp(exps) 

130 sources = detect_and_deblend(mbexp) 

131 

132 exposure = mbexp.singles[0] 

133 exp_copy = afw_image.ExposureF(exposure, deep=True) 

134 with util.ContextNoiseReplacer(exposure, sources, rng) as replacer: 

135 

136 assert np.any(exp_copy.image.array != exposure.image.array) 

137 

138 replaced_copy = afw_image.ExposureF(exposure, deep=True) 

139 for source in sources[:5]: 

140 with replacer.sourceInserted(source.getId()): 

141 assert np.any(exp_copy.image.array != exposure.image.array) 

142 

143 assert np.all(replaced_copy.image.array == exposure.image.array) 

144 

145 assert np.all(exp_copy.image.array == exposure.image.array) 

146 

147 

148@pytest.mark.skipif( 

149 skip_tests_on_simulations, 

150 reason='descwl_shear_sims not available' 

151) 

152def test_multiband_noise_replacer(show=False): 

153 import lsst.afw.image as afw_image 

154 from metadetect.lsst import vis 

155 

156 seed = 981 

157 rng = np.random.RandomState(seed) 

158 

159 sim = make_lsst_sim(rng) 

160 

161 exps = [texps[0] for _, texps in sim['band_data'].items()] 

162 mbexp = util.get_mbexp(exps) 

163 

164 sources = detect_and_deblend(mbexp) 

165 

166 exp_copys = [afw_image.ExposureF(exp, deep=True) for exp in exps] 

167 mbexp_copy = util.get_mbexp(exp_copys) 

168 

169 bands = list(sim['band_data'].keys()) 

170 

171 with util.MultibandNoiseReplacer(mbexp, sources, rng) as replacer: 

172 if show: 

173 vis.compare_mbexp(mbexp_copy, replacer.mbexp) 

174 

175 for exp_copy, exp in zip(exp_copys, replacer.mbexp.singles): 

176 assert np.any(exp_copy.image.array != exp.image.array) 

177 

178 replaced_copys = [ 

179 afw_image.ExposureF(exp, deep=True) for exp in replacer.mbexp.singles 

180 ] 

181 rmbexp_copy = util.get_mbexp(replaced_copys) 

182 for source in sources[:5]: 

183 source_id = source.getId() 

184 with replacer.sourceInserted(source_id): 

185 

186 if show: 

187 vis.compare_mbexp(rmbexp_copy, replacer.mbexp) 

188 

189 for band, rexp, exp in zip( 

190 bands, replaced_copys, replacer.mbexp.singles 

191 ): 

192 

193 assert np.any(rexp.image.array != exp.image.array), ( 

194 f'after inserting source {source_id}, band {band} does not ' 

195 f'differ from replaced exp' 

196 ) 

197 

198 for cexp, exp in zip(replaced_copys, replacer.mbexp.singles): 

199 assert np.all(cexp.image.array == exp.image.array), ( 

200 f'band {band} did not get restored to noise' 

201 ) 

202 

203 for cexp, exp in zip(exp_copys, mbexp.singles): 

204 assert np.all(cexp.image.array == exp.image.array) 

205 

206 

207if __name__ == '__main__': 

208 test_multiband_noise_replacer(show=True)