Coverage for tests/test_assemble_cell_coadd.py: 91%

164 statements  

« prev     ^ index     » next       coverage.py v7.16.0, created at 2026-09-11 11:03 +0000

1# This file is part of drp_tasks. 

2# 

3# LSST Data Management System 

4# This product includes software developed by the 

5# LSST Project (http://www.lsst.org/). 

6# See COPYRIGHT file at the top of the source tree. 

7# 

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

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

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

11# (at your option) any later version. 

12# 

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

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

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

16# GNU General Public License for more details. 

17# 

18# You should have received a copy of the LSST License Statement and 

19# the GNU General Public License along with this program. If not, 

20# see <https://www.lsstcorp.org/LegalNotices/>. 

21# 

22 

23from __future__ import annotations 

24 

25import unittest 

26import warnings 

27from typing import TYPE_CHECKING, Iterable 

28 

29import galsim 

30import hpgeom as hpg 

31import numpy as np 

32from assemble_coadd_test_utils import MockCoaddTestData, makeMockSkyInfo 

33 

34import lsst.afw.geom as afwGeom 

35import lsst.afw.image as afwImage 

36import lsst.pipe.base as pipeBase 

37import lsst.utils.tests 

38from lsst.cell_coadds import CoaddUnits, CommonComponents, PatchIdentifiers 

39from lsst.drp.tasks.assemble_cell_coadd import ( 

40 AssembleCellCoaddConfig, 

41 AssembleCellCoaddTask, 

42 WarpInputs, 

43) 

44from lsst.images.cells import CellCoadd 

45from lsst.skymap import Index2D 

46 

47if TYPE_CHECKING: 

48 from lsst.cell_coadds import ObservationIdentifiers 

49 

50__all__ = ( 

51 "MockAssembleCellCoaddConfig", 

52 "MockAssembleCellCoaddTask", 

53) 

54 

55 

56class MockAssembleCellCoaddConfig(AssembleCellCoaddConfig): 

57 def setDefaults(self): 

58 super().setDefaults() 

59 self.require_artifact_mask = False 

60 

61 

62class MockAssembleCellCoaddTask(AssembleCellCoaddTask): 

63 """Lightly modified version of `AssembleCellCoaddTask` for unit tests. 

64 

65 The modifications bypass the usual middleware for loading data and setting 

66 up the Task, and instead supply in-memory mock data references to the `run` 

67 method so that the coaddition algorithms can be tested without a Butler. 

68 """ 

69 

70 ConfigClass = MockAssembleCellCoaddConfig 

71 

72 def runQuantum( 

73 self, 

74 mockSkyInfo, 

75 warpRefList, 

76 maskedFractionRefList, 

77 noise0RefList, 

78 visitSummaryList=None, 

79 ): 

80 """Modified interface for testing coaddition algorithms without a 

81 Butler. 

82 

83 Parameters 

84 ---------- 

85 mockSkyInfo : `lsst.pipe.base.Struct` 

86 A simple container that supplies a bounding box and WCS in the 

87 same format as the output of 

88 `lsst.pipe.tasks.CoaddBaseTask.getSkyInfo` 

89 warpRefList : `list` of `lsst.pipe.tasks.MockExposureReference` 

90 Data references to the test exposures that will be coadded, 

91 using the Gen 3 API. 

92 

93 Returns 

94 ------- 

95 retStruct : `lsst.pipe.base.Struct` 

96 The coadded exposure and associated metadata. 

97 """ 

98 

99 self.common = CommonComponents( 

100 units=CoaddUnits.nJy, 

101 wcs=mockSkyInfo.wcs, 

102 band="i", 

103 identifiers=PatchIdentifiers(skymap="mock", tract=0, patch=Index2D(0, 0), band="i"), 

104 ) 

105 

106 inputs = {} 

107 for warpInput, maskedFractionInput, noise0Input in zip( 

108 warpRefList, 

109 maskedFractionRefList, 

110 noise0RefList, 

111 ): 

112 inputs[warpInput.dataId] = WarpInputs( 

113 warp=warpInput, 

114 masked_fraction=maskedFractionInput, 

115 noise_warps=[noise0Input], 

116 ) 

117 

118 retStruct = self.run( 

119 inputs=inputs, 

120 skyInfo=mockSkyInfo, 

121 visitSummaryList=visitSummaryList, 

122 ) 

123 

124 return retStruct 

125 

126 

127class AssembleCellCoaddTestCase(lsst.utils.tests.TestCase): 

128 """Tests of AssembleCellCoaddTask. 

129 

130 These tests bypass the middleware used for accessing data and managing Task 

131 execution. 

132 """ 

133 

134 @classmethod 

135 def setUpClass(cls) -> None: 

136 np.random.seed(42) 

137 rng = np.random.Generator(np.random.MT19937(42)) 

138 patch = 42 

139 tract = 0 

140 testData = MockCoaddTestData(fluxRange=1e4) 

141 exposures = {} 

142 matchedExposures = {} 

143 masked_fraction_images = {} 

144 noise0_masked_images = {} 

145 for expId in range(100, 110): 

146 exposures[expId], matchedExposures[expId] = testData.makeTestImage(expId) 

147 masked_fraction_images[expId] = afwImage.ImageF(bbox=exposures[expId].getBBox()) 

148 masked_fraction_images[expId].array[:, :] = rng.random(masked_fraction_images[expId].array.shape) 

149 noise0_masked_images[expId] = afwImage.MaskedImageF(bbox=exposures[expId].getBBox()) 

150 noise0_masked_images[expId].image.array[:, :] = rng.normal( 

151 0, 1, noise0_masked_images[expId].image.array.shape 

152 ) * (exposures[expId].variance.array**0.5) 

153 

154 cls.handleList = testData.makeDataRefList(exposures, patch=patch, tract=tract) 

155 cls.maskedFractionRefList = testData.makeDataRefList(masked_fraction_images, patch=patch, tract=tract) 

156 cls.noise0RefList = testData.makeDataRefList(noise0_masked_images, patch=patch, tract=tract) 

157 cls.visitSummaryList = [ 

158 testData.makeVisitSummaryTableHandle(warpHandle) for warpHandle in cls.handleList 

159 ] 

160 cls.skyInfo = makeMockSkyInfo(testData.bbox, testData.wcs, patch=patch) 

161 

162 def runTask( 

163 self, 

164 config=None, 

165 warpRefList=None, 

166 maskedFractionRefList=None, 

167 noise0RefList=None, 

168 visitSummaryList=None, 

169 ) -> None: 

170 if config is None: 

171 config = MockAssembleCellCoaddConfig() 

172 config.do_input_map = True 

173 assembleTask = MockAssembleCellCoaddTask(config=config) 

174 if warpRefList is None: 

175 warpRefList = self.handleList 

176 if maskedFractionRefList is None: 

177 maskedFractionRefList = self.maskedFractionRefList 

178 if noise0RefList is None: 

179 noise0RefList = self.noise0RefList 

180 if visitSummaryList is None: 

181 visitSummaryList = self.visitSummaryList 

182 

183 self.result = assembleTask.runQuantum( 

184 self.skyInfo, 

185 warpRefList=warpRefList, 

186 maskedFractionRefList=maskedFractionRefList, 

187 noise0RefList=noise0RefList, 

188 visitSummaryList=visitSummaryList, 

189 ) 

190 

191 def checkSortOrder(self, inputs: Iterable[ObservationIdentifiers]) -> None: 

192 """Check that the inputs are sorted. 

193 

194 The inputs must be sorted first by visit, and within the same visit, 

195 by detector. 

196 

197 Parameters 

198 ---------- 

199 inputs : `Iterable` [`ObservationIdentifiers`] 

200 The inputs to be checked. 

201 """ 

202 visit, detector = -np.inf, -np.inf # Previous visit, detector IDs. 

203 for _, obsId in enumerate(inputs): 

204 with self.subTest(input_number=repr(obsId)): 

205 self.assertGreaterEqual(obsId.visit, visit) 

206 if visit == obsId.visit: 206 ↛ 207line 206 didn't jump to line 207 because the condition on line 206 was never true

207 with self.subTest(detector_number=repr(obsId.detector)): 

208 self.assertGreaterEqual(obsId.detector, detector) 

209 

210 visit, detector = obsId.visit, obsId.detector 

211 

212 def checkRun(self, assembleTask): 

213 """Check that the task runs successfully.""" 

214 result = assembleTask.runQuantum(self.skyInfo, self.handleList) 

215 

216 # Check that we produced an exposure. 

217 self.assertTrue(result.multipleCellCoadd is not None) 

218 # Check that the visit_count method returns a number less than or equal 

219 # to the total number of input exposures available. 

220 max_visit_count = len(self.handleList) 

221 for cellId, singleCellCoadd in result.multipleCellCoadd.cells.items(): 

222 with self.subTest(x=repr(cellId.x), y=repr(cellId.y)): 

223 self.assertLessEqual(singleCellCoadd.visit_count, max_visit_count) 

224 # Check that the aperture correction maps are not None. 

225 with self.subTest(x=repr(cellId.x), y=repr(cellId.y)): 

226 self.assertTrue(singleCellCoadd.aperture_correction_map is not None) 

227 # Check that the inputs are sorted. 

228 with self.subTest(x=repr(cellId.x), y=repr(cellId.y)): 

229 self.checkSortOrder(singleCellCoadd.inputs) 

230 

231 def test_assemble_basic(self): 

232 """Test that AssembleCellCoaddTask runs successfully without errors. 

233 

234 This test does not check the correctness of the coaddition algorithms. 

235 This is intended to prevent the code from bit rotting. 

236 """ 

237 self.runTask() 

238 # Check that we produced an exposure. 

239 self.assertTrue(self.result.multipleCellCoadd is not None) 

240 self.assertTrue(self.result.inputMap is not None) 

241 

242 # Check the input map. 

243 inputMap = self.result.inputMap 

244 

245 warp_input_list = [handle.get() for handle in self.handleList] 

246 visit_detectors = [] 

247 for warp_input in warp_input_list: 

248 for row in warp_input.getInfo().getCoaddInputs().ccds: 

249 visit_detectors.append((int(row["visit"]), int(row["ccd"]))) 

250 

251 for bit, (visit, detector) in enumerate(visit_detectors): 

252 self.assertEqual(inputMap.metadata[f"B{bit:04d}VIS"], visit) 

253 self.assertEqual(inputMap.metadata[f"B{bit:04d}CCD"], detector) 

254 self.assertGreater(inputMap.metadata[f"B{bit:04d}WT"], 0.0) 

255 

256 coadd_poly = afwGeom.Polygon(lsst.geom.Box2D(self.result.multipleCellCoadd.outer_bbox)) 

257 sph_pts = self.result.multipleCellCoadd.wcs.pixelToSky(coadd_poly) 

258 radec = np.asarray([(sph.getRa().asDegrees(), sph.getDec().asDegrees()) for sph in sph_pts]) 

259 pixels = hpg.query_polygon(inputMap.nside_sparse, radec[:-1, 0], radec[:-1, 1]) 

260 np.testing.assert_array_equal(pixels, inputMap.valid_pixels) 

261 

262 def test_assemble_future_output(self): 

263 """Test that AssembleCellCoaddTask runs successfully and produces 

264 an lsst.images.cells.CellCoadd when configured to. 

265 """ 

266 config = MockAssembleCellCoaddConfig() 

267 config.output_image_type = "future" 

268 self.runTask(config=config) 

269 # Check that we produced an exposure. 

270 self.assertIsInstance(self.result.multipleCellCoadd, CellCoadd) 

271 

272 def test_assemble_empty(self): 

273 """Test that AssembleCellCoaddTask runs successfully without errors 

274 when no input exposures are provided.""" 

275 self.result = None # so tearDown has something. 

276 with self.assertRaises(pipeBase.NoWorkFound, msg="No cells could be populated for the cell coadd."): 

277 self.runTask(warpRefList=[], maskedFractionRefList=[], noise0RefList=[], visitSummaryList=[]) 

278 

279 def test_assemble_without_visitSummary(self): 

280 """Test that AssembleCellCoaddTask calculates detector weights and 

281 runs successfully without errors when no visit summaries are provided. 

282 """ 

283 self.runTask(visitSummaryList=[]) 

284 # Check that we produced an exposure. 

285 self.assertTrue(self.result.multipleCellCoadd is not None) 

286 self.assertTrue(self.result.inputMap is not None) 

287 

288 # TODO: Remove this test in DM-49401 

289 @lsst.utils.tests.methodParameters(do_scale_zero_point=[False, True]) 

290 def test_do_scale_zero_point(self, do_scale_zero_point): 

291 config = MockAssembleCellCoaddConfig() 

292 with warnings.catch_warnings(): 

293 warnings.filterwarnings("ignore", category=FutureWarning) 

294 config.do_scale_zero_point = do_scale_zero_point 

295 self.runTask(config) 

296 # Check that we produced an exposure. 

297 self.assertTrue(self.result.multipleCellCoadd is not None) 

298 self.assertTrue(self.result.inputMap is not None) 

299 

300 @lsst.utils.tests.methodParameters(do_calculate_weight_from_warp=[False, True]) 

301 def test_do_calculate_weight_from_warp(self, do_calculate_weight_from_warp): 

302 config = MockAssembleCellCoaddConfig() 

303 config.do_calculate_weight_from_warp = do_calculate_weight_from_warp 

304 self.runTask(config) 

305 # Check that we produced an exposure. 

306 self.assertTrue(self.result.multipleCellCoadd is not None) 

307 self.assertTrue(self.result.inputMap is not None) 

308 

309 def test_visit_count(self): 

310 """Check that the visit_count method returns a number less than or 

311 equal to the total number of input exposures available. 

312 """ 

313 self.runTask() 

314 max_visit_count = len(self.handleList) 

315 for cellId, singleCellCoadd in self.result.multipleCellCoadd.cells.items(): 

316 with self.subTest(x=repr(cellId.x), y=repr(cellId.y)): 

317 self.assertLessEqual(singleCellCoadd.visit_count, max_visit_count) 

318 

319 def test_inputs_sorted(self): 

320 """Check that the inputs are sorted. 

321 

322 The ordering is that inputs are sorted first by visit, and within the 

323 same visit, they are ordered by detector. 

324 """ 

325 self.runTask() 

326 for _, singleCellCoadd in self.result.multipleCellCoadd.cells.items(): 

327 self.checkSortOrder(singleCellCoadd.inputs) 

328 

329 def test_psf_normalization(self): 

330 """Check that the sum of PSF images is close to 1.""" 

331 self.runTask() 

332 for cellId, singleCellCoadd in self.result.multipleCellCoadd.cells.items(): 

333 with self.subTest(x=repr(cellId.x), y=repr(cellId.y)): 

334 self.assertFloatsAlmostEqual(singleCellCoadd.psf_image.array.sum(), 1.0, rtol=None, atol=1e-7) 

335 

336 def test_psf_centering(self): 

337 """Check that the PSF images are centered.""" 

338 self.runTask() 

339 for cellId, singleCellCoadd in self.result.multipleCellCoadd.cells.items(): 

340 with self.subTest(x=repr(cellId.x), y=repr(cellId.y)): 

341 shape = galsim.hsm.FindAdaptiveMom(galsim.Image(singleCellCoadd.psf_image.array)) 

342 self.assertFloatsAlmostEqual( 

343 shape.moments_centroid.x, round(shape.moments_centroid.x), rtol=None, atol=0.01 

344 ) 

345 self.assertFloatsAlmostEqual( 

346 shape.moments_centroid.y, round(shape.moments_centroid.y), rtol=None, atol=0.01 

347 ) 

348 

349 

350class MyMemoryTestCase(lsst.utils.tests.MemoryTestCase): 

351 pass 

352 

353 

354def setup_module(module): 

355 lsst.utils.tests.init() 

356 

357 

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

359 lsst.utils.tests.init() 

360 unittest.main()