Coverage for tests/test_edgeBleed.py: 98%

1# This file is part of ip_isr. 

2# 

3# Developed for the LSST Data Management System. 

4# This product includes software developed by the LSST Project 

5# (https://www.lsst.org). 

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

7# for details of code ownership. 

8# 

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

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

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

12# (at your option) any later version. 

13# 

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

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

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

17# GNU General Public License for more details. 

18# 

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

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

21import unittest 

22 

23import numpy as np 

24 

25import lsst.geom as geom 

26import lsst.afw.cameraGeom as cameraGeom 

27import lsst.afw.image as afwImage 

28import lsst.utils.tests 

29import lsst.ip.isr.isrFunctions as ipIsrFunctions 

30from lsst.ip.isr.isrTask import IsrTaskConfig 

31from lsst.ip.isr.masking import DECamEdgeBleedMaskTask 

32 

33AMP_WIDTH = 200 

34AMP_HEIGHT = 800 

35SKY = 1000.0 

36NOISE = 10.0 

37DIP = 200.0 

38 

39 

40def makeDetector(): 

41 """Two side-by-side amplifiers: A reads out at the top, B at the bottom.""" 

42 camBuilder = cameraGeom.Camera.Builder("testCam") 

43 detBuilder = camBuilder.add("testDet", 0) 

44 detBuilder.setBBox(geom.Box2I(geom.Point2I(0, 0), geom.Extent2I(2*AMP_WIDTH, AMP_HEIGHT))) 

45 detBuilder.setPixelSize(geom.Extent2D(0.015, 0.015)) 

46 detBuilder.setOrientation(cameraGeom.Orientation()) 

47 for name, x0, corner in (("A", 0, cameraGeom.ReadoutCorner.UL), 

48 ("B", AMP_WIDTH, cameraGeom.ReadoutCorner.LL)): 

49 bbox = geom.Box2I(geom.Point2I(x0, 0), geom.Extent2I(AMP_WIDTH, AMP_HEIGHT)) 

50 amp = cameraGeom.Amplifier.Builder() 

51 amp.setName(name) 

52 amp.setBBox(bbox) 

53 amp.setRawBBox(bbox) 

54 amp.setRawDataBBox(bbox) 

55 amp.setReadoutCorner(corner) 

56 amp.setGain(1.0) 

57 amp.setReadNoise(5.0) 

58 amp.setSaturation(60000.0) 

59 detBuilder.append(amp) 

60 return camBuilder.finish()["testDet"] 

61 

62 

63class MaskDECamEdgeBleedTestCase(lsst.utils.tests.TestCase): 

64 def setUp(self): 

65 self.detector = makeDetector() 

66 self.exposure = afwImage.ExposureF(self.detector.getBBox()) 

67 self.exposure.setDetector(self.detector) 

68 rng = np.random.default_rng(12345) 

69 self.exposure.image.array[:] = SKY + rng.normal(0.0, NOISE, self.exposure.image.array.shape) 

70 self.exposure.variance.array[:] = NOISE**2 

71 self.satBit = self.exposure.mask.getPlaneBitMask("SAT") 

72 

73 def addSatBlock(self, x0, y0, width, height): 

74 """Mark a width x height block as SAT (image value irrelevant).""" 

75 self.exposure.mask.array[y0:y0 + height, x0:x0 + width] |= self.satBit 

76 self.exposure.image.array[y0:y0 + height, x0:x0 + width] = 60000.0 

77 

78 def addDip(self, x0, y0, width, height): 

79 """Depress a width x height block below sky.""" 

80 self.exposure.image.array[y0:y0 + height, x0:x0 + width] -= DIP 

81 

82 def runMasking(self): 

83 before = self.exposure.mask.array.copy() 

84 ipIsrFunctions.maskDECamEdgeBleed(self.exposure) 

85 return before, self.exposure.mask.array 

86 

87 def test_topReadEdgeBleedIsMasked(self): 

88 # 100 x 200 = 20000 SAT pixels touching the top edge of amp A. 

89 self.addSatBlock(50, AMP_HEIGHT - 200, 100, 200) 

90 # Dip in the 100 rows nearest the top edge, full width of amp A. 

91 self.addDip(0, AMP_HEIGHT - 100, AMP_WIDTH, 100) 

92 

93 before, after = self.runMasking() 

94 

95 # height 100 -> margin int(100*0.125) + 1 = 13 -> 113 rows masked. 

96 expected = before.copy() 

97 expected[AMP_HEIGHT - 113:, 0:AMP_WIDTH] |= self.satBit 

98 np.testing.assert_array_equal(after, expected) 

99 

100 def test_noDipIsNotMasked(self): 

101 self.addSatBlock(50, AMP_HEIGHT - 200, 100, 200) 

102 

103 before, after = self.runMasking() 

104 

105 np.testing.assert_array_equal(after, before) 

106 

107 def test_smallFootprintIsNotMasked(self): 

108 # 50 x 100 = 5000 SAT pixels: below satMinArea. 

109 self.addSatBlock(50, AMP_HEIGHT - 100, 50, 100) 

110 self.addDip(0, AMP_HEIGHT - 100, AMP_WIDTH, 100) 

111 

112 before, after = self.runMasking() 

113 

114 np.testing.assert_array_equal(after, before) 

115 

116 def test_largeFootprintIsNotMasked(self): 

117 # 300 x 400 = 120000 SAT pixels: above satMaxArea. Only x 0-49 of the 

118 # dip rows are unsaturated, still enough low pixels per row to be 

119 # detectable, so only the area cut prevents masking. 

120 self.addSatBlock(50, 400, 300, 400) 

121 self.addDip(0, AMP_HEIGHT - 100, AMP_WIDTH, 100) 

122 

123 before, after = self.runMasking() 

124 

125 np.testing.assert_array_equal(after, before) 

126 

127 def test_footprintFarFromEdgeIsNotMasked(self): 

128 # Footprint top at row 599, 200 rows short of the read edge. 

129 self.addSatBlock(50, 400, 100, 200) 

130 self.addDip(0, AMP_HEIGHT - 100, AMP_WIDTH, 100) 

131 

132 before, after = self.runMasking() 

133 

134 np.testing.assert_array_equal(after, before) 

135 

136 def test_bottomReadEdgeBleedIsMasked(self): 

137 # Amp B reads out at the bottom. 

138 self.addSatBlock(AMP_WIDTH + 50, 0, 100, 200) 

139 self.addDip(AMP_WIDTH, 0, AMP_WIDTH, 100) 

140 

141 before, after = self.runMasking() 

142 

143 expected = before.copy() 

144 # height 100 -> margin int(100*0.125) + 1 = 13 -> 113 rows masked. 

145 expected[:113, AMP_WIDTH:2*AMP_WIDTH] |= self.satBit 

146 np.testing.assert_array_equal(after, expected) 

147 

148 def test_dipInOtherAmpIsNotMasked(self): 

149 # Footprint in amp A (top read edge); dip only in amp B's top rows, 

150 # which is not B's read edge and has no footprint. 

151 self.addSatBlock(50, AMP_HEIGHT - 200, 100, 200) 

152 self.addDip(AMP_WIDTH, AMP_HEIGHT - 100, AMP_WIDTH, 100) 

153 

154 before, after = self.runMasking() 

155 

156 np.testing.assert_array_equal(after, before) 

157 

158 def test_unusableEdgeRowsAreSkippedButMasked(self): 

159 # The 20 rows nearest the read edge are NO_DATA (as after trimming 

160 # on DECam), so the dip can only be seen from row 20 inward. The 

161 # measured height still counts from the physical edge. 

162 noDataBit = self.exposure.mask.getPlaneBitMask("NO_DATA") 

163 self.exposure.mask.array[AMP_HEIGHT - 20:, 0:AMP_WIDTH] |= noDataBit 

164 self.addSatBlock(50, AMP_HEIGHT - 200, 100, 200) 

165 self.addDip(0, AMP_HEIGHT - 100, AMP_WIDTH, 80) 

166 

167 before, after = self.runMasking() 

168 

169 expected = before.copy() 

170 expected[AMP_HEIGHT - 113:, 0:AMP_WIDTH] |= self.satBit 

171 np.testing.assert_array_equal(after, expected) 

172 

173 def test_suspectPixelsCountAsLow(self): 

174 # DECam flags the rows nearest the read edge SUSPECT; those pixels 

175 # must still be counted when confirming the dip. 

176 suspectBit = self.exposure.mask.getPlaneBitMask("SUSPECT") 

177 self.exposure.mask.array[AMP_HEIGHT - 35:, 0:AMP_WIDTH] |= suspectBit 

178 self.addSatBlock(50, AMP_HEIGHT - 200, 100, 200) 

179 self.addDip(0, AMP_HEIGHT - 100, AMP_WIDTH, 100) 

180 

181 before, after = self.runMasking() 

182 

183 expected = before.copy() 

184 expected[AMP_HEIGHT - 113:, 0:AMP_WIDTH] |= self.satBit 

185 np.testing.assert_array_equal(after, expected) 

186 

187 def test_nonDefaultParameters(self): 

188 # A 100-row footprint ending 150 rows short of the edge is only a 

189 # candidate with a larger approachRows; a 50 ADU dip is only seen 

190 # with a smaller nSigma; a zero marginFraction masks height + 1 rows. 

191 self.addSatBlock(50, AMP_HEIGHT - 250, 200, 100) 

192 self.exposure.image.array[AMP_HEIGHT - 100:, 0:AMP_WIDTH] -= 50.0 

193 

194 before = self.exposure.mask.array.copy() 

195 ipIsrFunctions.maskDECamEdgeBleed(self.exposure, approachRows=150, nSigma=3.0, 

196 marginFraction=0.0) 

197 

198 expected = before.copy() 

199 expected[AMP_HEIGHT - 101:, 0:AMP_WIDTH] |= self.satBit 

200 np.testing.assert_array_equal(self.exposure.mask.array, expected) 

201 

202 def test_alternateSaturatedMaskName(self): 

203 self.exposure.mask.addMaskPlane("MYSAT") 

204 mySatBit = self.exposure.mask.getPlaneBitMask("MYSAT") 

205 self.exposure.mask.array[AMP_HEIGHT - 200:, 50:150] |= mySatBit 

206 self.addDip(0, AMP_HEIGHT - 100, AMP_WIDTH, 100) 

207 

208 before = self.exposure.mask.array.copy() 

209 ipIsrFunctions.maskDECamEdgeBleed(self.exposure, saturatedMaskName="MYSAT") 

210 

211 expected = before.copy() 

212 expected[AMP_HEIGHT - 113:, 0:AMP_WIDTH] |= mySatBit 

213 np.testing.assert_array_equal(self.exposure.mask.array, expected) 

214 

215 def test_partialAmpIsSkipped(self): 

216 # An exposure covering only the bottom 600 rows of amp B has a 

217 # qualifying footprint and dip, but the amp is not fully contained 

218 # so nothing is masked (rather than failing on the sub-image). 

219 self.addSatBlock(AMP_WIDTH + 50, 0, 100, 200) 

220 self.addDip(AMP_WIDTH, 0, AMP_WIDTH, 100) 

221 partial = geom.Box2I(geom.Point2I(AMP_WIDTH, 0), geom.Extent2I(AMP_WIDTH, 600)) 

222 subExposure = self.exposure[partial] 

223 

224 before = subExposure.mask.array.copy() 

225 ipIsrFunctions.maskDECamEdgeBleed(subExposure) 

226 

227 np.testing.assert_array_equal(subExposure.mask.array, before) 

228 

229 def test_noDetectorRaises(self): 

230 exposure = afwImage.ExposureF(geom.Box2I(geom.Point2I(0, 0), geom.Extent2I(50, 50))) 

231 with self.assertRaises(RuntimeError): 

232 ipIsrFunctions.maskDECamEdgeBleed(exposure) 

233 

234 

235class DECamEdgeBleedMaskTaskTestCase(lsst.utils.tests.TestCase): 

236 def setUp(self): 

237 self.detector = makeDetector() 

238 self.exposure = afwImage.ExposureF(self.detector.getBBox()) 

239 self.exposure.setDetector(self.detector) 

240 rng = np.random.default_rng(12345) 

241 self.exposure.image.array[:] = SKY + rng.normal(0.0, NOISE, self.exposure.image.array.shape) 

242 self.exposure.variance.array[:] = NOISE**2 

243 self.satBit = self.exposure.mask.getPlaneBitMask("SAT") 

244 

245 def test_defaults(self): 

246 config = DECamEdgeBleedMaskTask.ConfigClass() 

247 self.assertEqual(config.satMinArea, 10000) 

248 self.assertEqual(config.satMaxArea, 100000) 

249 self.assertEqual(config.approachRows, 20) 

250 self.assertEqual(config.nSigma, 5.0) 

251 self.assertEqual(config.nRowsCheck, 20) 

252 self.assertEqual(config.minLowPixelsPerRow, 30) 

253 self.assertEqual(config.minLowPixelsExtent, 10) 

254 self.assertEqual(config.marginFraction, 0.125) 

255 self.assertEqual(config.saturatedMaskName, "SAT") 

256 config.validate() 

257 

258 def test_validateExtentBelowPerRow(self): 

259 config = DECamEdgeBleedMaskTask.ConfigClass() 

260 config.minLowPixelsExtent = config.minLowPixelsPerRow 

261 with self.assertRaises(ValueError): 

262 config.validate() 

263 

264 def test_runMasksEdgeBleed(self): 

265 self.exposure.mask.array[AMP_HEIGHT - 200:, 50:150] |= self.satBit 

266 self.exposure.image.array[AMP_HEIGHT - 100:, 0:AMP_WIDTH] -= DIP 

267 before = self.exposure.mask.array.copy() 

268 

269 task = DECamEdgeBleedMaskTask() 

270 task.run(self.exposure) 

271 

272 expected = before.copy() 

273 expected[AMP_HEIGHT - 113:, 0:AMP_WIDTH] |= self.satBit 

274 np.testing.assert_array_equal(self.exposure.mask.array, expected) 

275 

276 def test_retargetIntoIsrTaskConfig(self): 

277 config = IsrTaskConfig() 

278 config.doCameraSpecificMasking = True 

279 config.masking.retarget(DECamEdgeBleedMaskTask) 

280 config.masking.nSigma = 3.0 

281 config.validate() 

282 self.assertEqual(config.masking.nSigma, 3.0) 

283 self.assertIs(config.masking.target, DECamEdgeBleedMaskTask) 

284 

285 

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

287 pass 

288 

289 

290def setup_module(module): 

291 lsst.utils.tests.init() 

292 

293 

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

295 lsst.utils.tests.init() 

296 unittest.main()