Coverage for python / lsst / obs / base / yamlCamera.py: 9%

141 statements  

« prev     ^ index     » next       coverage.py v7.13.1, created at 2026-01-20 01:41 -0800

1# This file is part of obs_base. 

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/>. 

21 

22from functools import lru_cache 

23 

24import numpy as np 

25import yaml 

26 

27import lsst.afw.cameraGeom as cameraGeom 

28import lsst.afw.geom as afwGeom 

29import lsst.geom as geom 

30from lsst.afw.cameraGeom import Amplifier, Camera, ReadoutCorner 

31 

32__all__ = ["makeCamera"] 

33 

34 

35@lru_cache 

36def makeCamera(cameraFile): 

37 """Construct an imaging camera (e.g. the LSST 3Gpix camera). 

38 

39 Parameters 

40 ---------- 

41 cameraFile : `str` 

42 Camera description YAML file. 

43 

44 Returns 

45 ------- 

46 camera : `lsst.afw.cameraGeom.Camera` 

47 The desired Camera 

48 """ 

49 with open(cameraFile) as fd: 

50 cameraParams = yaml.load(fd, Loader=yaml.CLoader) 

51 

52 cameraName = cameraParams["name"] 

53 

54 # 

55 # Handle distortion models. 

56 # 

57 plateScale = geom.Angle(cameraParams["plateScale"], geom.arcseconds) 

58 nativeSys = cameraGeom.CameraSys(cameraParams["transforms"].pop("nativeSys")) 

59 transforms = makeTransformDict(nativeSys, cameraParams["transforms"], plateScale) 

60 

61 ccdParams = cameraParams["CCDs"] 

62 detectorConfigList = makeDetectorConfigList(ccdParams) 

63 focalPlaneParity = cameraParams.get("focalPlaneParity", False) 

64 

65 amplifierDict = {} 

66 for ccdName, ccdValues in ccdParams.items(): 

67 amplifierDict[ccdName] = makeAmplifierList(ccdValues) 

68 

69 return makeCameraFromCatalogs( 

70 cameraName, 

71 detectorConfigList, 

72 nativeSys, 

73 transforms, 

74 amplifierDict, 

75 focalPlaneParity=focalPlaneParity, 

76 ) 

77 

78 

79def makeDetectorConfigList(ccdParams): 

80 """Make a list of detector configs. 

81 

82 Returns 

83 ------- 

84 detectorConfig : `list` of `lsst.afw.cameraGeom.DetectorConfig` 

85 A list of detector configs. 

86 """ 

87 detectorConfigs = [] 

88 for name, ccd in ccdParams.items(): 

89 detectorConfig = cameraGeom.DetectorConfig() 

90 detectorConfigs.append(detectorConfig) 

91 

92 detectorConfig.name = name 

93 detectorConfig.id = ccd["id"] 

94 detectorConfig.serial = ccd["serial"] 

95 detectorConfig.detectorType = ccd["detectorType"] 

96 if "physicalType" in ccd: 

97 detectorConfig.physicalType = ccd["physicalType"] 

98 # This is the orientation we need to put the serial direction along 

99 # the x-axis 

100 detectorConfig.bbox_x0, detectorConfig.bbox_y0 = ccd["bbox"][0] 

101 detectorConfig.bbox_x1, detectorConfig.bbox_y1 = ccd["bbox"][1] 

102 detectorConfig.pixelSize_x, detectorConfig.pixelSize_y = ccd["pixelSize"] 

103 detectorConfig.transformDict.nativeSys = ccd["transformDict"]["nativeSys"] 

104 transforms = ccd["transformDict"]["transforms"] 

105 detectorConfig.transformDict.transforms = None if transforms == "None" else transforms 

106 detectorConfig.refpos_x, detectorConfig.refpos_y = ccd["refpos"] 

107 if len(ccd["offset"]) == 2: 

108 detectorConfig.offset_x, detectorConfig.offset_y = ccd["offset"] 

109 detectorConfig.offset_z = 0.0 

110 else: 

111 detectorConfig.offset_x, detectorConfig.offset_y, detectorConfig.offset_z = ccd["offset"] 

112 detectorConfig.transposeDetector = ccd["transposeDetector"] 

113 detectorConfig.pitchDeg = ccd["pitch"] 

114 detectorConfig.yawDeg = ccd["yaw"] 

115 detectorConfig.rollDeg = ccd["roll"] 

116 if "crosstalk" in ccd: 

117 detectorConfig.crosstalk = ccd["crosstalk"] 

118 

119 return detectorConfigs 

120 

121 

122def makeAmplifierList(ccd): 

123 """Construct a list of AmplifierBuilder objects.""" 

124 # Much of this will need to be filled in when we know it. 

125 assert len(ccd) > 0 

126 amp = list(ccd["amplifiers"].values())[0] 

127 

128 rawBBox = makeBBoxFromList(amp["rawBBox"]) # total in file 

129 xRawExtent, yRawExtent = rawBBox.getDimensions() 

130 

131 readCorners = { 

132 "LL": ReadoutCorner.LL, 

133 "LR": ReadoutCorner.LR, 

134 "UL": ReadoutCorner.UL, 

135 "UR": ReadoutCorner.UR, 

136 } 

137 

138 amplifierList = [] 

139 for name, amp in sorted(ccd["amplifiers"].items(), key=lambda x: x[1]["hdu"]): 

140 amplifier = Amplifier.Builder() 

141 amplifier.setName(name) 

142 

143 ix, iy = amp["ixy"] 

144 perAmpData = amp["perAmpData"] 

145 if perAmpData: 

146 x0, y0 = 0, 0 # origin of data within each amp image 

147 else: 

148 x0, y0 = ix * xRawExtent, iy * yRawExtent 

149 

150 rawDataBBox = makeBBoxFromList(amp["rawDataBBox"]) # Photosensitive area 

151 xDataExtent, yDataExtent = rawDataBBox.getDimensions() 

152 amplifier.setBBox( 

153 geom.BoxI(geom.PointI(ix * xDataExtent, iy * yDataExtent), rawDataBBox.getDimensions()) 

154 ) 

155 

156 rawBBox = makeBBoxFromList(amp["rawBBox"]) 

157 rawBBox.shift(geom.ExtentI(x0, y0)) 

158 amplifier.setRawBBox(rawBBox) 

159 

160 rawDataBBox = makeBBoxFromList(amp["rawDataBBox"]) 

161 rawDataBBox.shift(geom.ExtentI(x0, y0)) 

162 amplifier.setRawDataBBox(rawDataBBox) 

163 

164 rawSerialOverscanBBox = makeBBoxFromList(amp["rawSerialOverscanBBox"]) 

165 rawSerialOverscanBBox.shift(geom.ExtentI(x0, y0)) 

166 amplifier.setRawHorizontalOverscanBBox(rawSerialOverscanBBox) 

167 

168 rawParallelOverscanBBox = makeBBoxFromList(amp["rawParallelOverscanBBox"]) 

169 rawParallelOverscanBBox.shift(geom.ExtentI(x0, y0)) 

170 amplifier.setRawVerticalOverscanBBox(rawParallelOverscanBBox) 

171 

172 rawSerialPrescanBBox = makeBBoxFromList(amp["rawSerialPrescanBBox"]) 

173 rawSerialPrescanBBox.shift(geom.ExtentI(x0, y0)) 

174 amplifier.setRawPrescanBBox(rawSerialPrescanBBox) 

175 

176 if perAmpData: 

177 amplifier.setRawXYOffset(geom.Extent2I(ix * xRawExtent, iy * yRawExtent)) 

178 else: 

179 amplifier.setRawXYOffset(geom.Extent2I(0, 0)) 

180 

181 amplifier.setReadoutCorner(readCorners[amp["readCorner"]]) 

182 amplifier.setGain(amp["gain"]) 

183 amplifier.setReadNoise(amp["readNoise"]) 

184 amplifier.setSaturation(amp["saturation"]) 

185 amplifier.setSuspectLevel(amp.get("suspect", np.nan)) 

186 

187 # flip data when assembling if needs be (e.g. data from the serial at 

188 # the top of a CCD) 

189 flipX, flipY = amp.get("flipXY") 

190 

191 amplifier.setRawFlipX(flipX) 

192 amplifier.setRawFlipY(flipY) 

193 # linearity placeholder stuff 

194 amplifier.setLinearityCoeffs([float(val) for val in amp["linearityCoeffs"]]) 

195 amplifier.setLinearityType(amp["linearityType"]) 

196 amplifier.setLinearityThreshold(float(amp["linearityThreshold"])) 

197 amplifier.setLinearityMaximum(float(amp["linearityMax"])) 

198 amplifier.setLinearityUnits("DN") 

199 amplifierList.append(amplifier) 

200 return amplifierList 

201 

202 

203def makeAmpInfoCatalog(ccd): 

204 """Backward compatible name.""" 

205 return makeAmplifierList(ccd) 

206 

207 

208def makeBBoxFromList(ylist): 

209 """Given a list [(x0, y0), (xsize, ysize)], probably from a yaml file, 

210 return a BoxI. 

211 """ 

212 (x0, y0), (xsize, ysize) = ylist 

213 return geom.BoxI(geom.PointI(x0, y0), geom.ExtentI(xsize, ysize)) 

214 

215 

216def makeTransformDict(nativeSys, transformDict, plateScale): 

217 """Make a dictionary of TransformPoint2ToPoint2s from yaml, mapping from 

218 nativeSys. 

219 

220 Parameters 

221 ---------- 

222 nativeSys : `lsst.afw.cameraGeom.CameraSys` 

223 transformDict : `dict` 

224 A dict specifying parameters of transforms; keys are camera system 

225 names. 

226 plateScale : `lsst.geom.Angle` 

227 The size of a pixel in angular units/mm (e.g. 20 arcsec/mm for LSST) 

228 

229 Returns 

230 ------- 

231 transforms : `dict` 

232 A dict of `lsst.afw.cameraGeom.CameraSys` : 

233 `lsst.afw.geom.TransformPoint2ToPoint2` 

234 

235 The resulting dict's keys are `~lsst.afw.cameraGeom.CameraSys`, 

236 and the values are Transforms *from* NativeSys *to* CameraSys 

237 """ 

238 # As other comments note this is required, and this is one function where 

239 # it's assumed 

240 assert nativeSys == cameraGeom.FOCAL_PLANE, "Cameras with nativeSys != FOCAL_PLANE are not supported." 

241 

242 resMap = {} 

243 

244 for key, transform in transformDict.items(): 

245 transformType = transform["transformType"] 

246 knownTransformTypes = ["affine", "radial"] 

247 if transformType not in knownTransformTypes: 

248 raise RuntimeError( 

249 "Saw unknown transform type for {}: {} (known types are: [{}])".format( 

250 key, transform["transformType"], ", ".join(knownTransformTypes) 

251 ) 

252 ) 

253 

254 if transformType == "affine": 

255 affine = geom.AffineTransform(np.array(transform["linear"]), np.array(transform["translation"])) 

256 

257 transform = afwGeom.makeTransform(affine) 

258 elif transformType == "radial": 

259 # radial coefficients of the form 

260 # [0, 1 (no units), C2 (mm^-1), usually 0, C4 (mm^-3), ...] 

261 # Radial distortion is modeled as a radial polynomial that converts 

262 # from focal plane radius (in mm) to field angle (in radians). 

263 # The provided coefficients are divided by the plate 

264 # scale (in radians/mm) meaning that C1 is always 1. 

265 radialCoeffs = np.array(transform["coeffs"]) 

266 

267 radialCoeffs *= plateScale.asRadians() 

268 transform = afwGeom.makeRadialTransform(radialCoeffs) 

269 else: 

270 raise RuntimeError( 

271 'Impossible condition "{}" is not in: [{}])'.format( 

272 transform["transformType"], ", ".join(knownTransformTypes) 

273 ) 

274 ) 

275 

276 resMap[cameraGeom.CameraSys(key)] = transform 

277 

278 return resMap 

279 

280 

281def makeCameraFromCatalogs( 

282 cameraName, 

283 detectorConfigList, 

284 nativeSys, 

285 transformDict, 

286 amplifierDict, 

287 pupilFactoryClass=cameraGeom.pupil.PupilFactory, 

288 focalPlaneParity=False, 

289): 

290 """Construct a Camera instance from a dictionary of 

291 detector name and `lsst.afw.cameraGeom.Amplifier`. 

292 

293 Parameters 

294 ---------- 

295 cameraName : `str` 

296 The name of the camera 

297 detectorConfigList : `list` 

298 A list of `lsst.afw.cameraGeom.cameraConfig.DetectorConfig` 

299 nativeSys : `lsst.afw.cameraGeom.CameraSys` 

300 The native transformation type; must be 

301 `lsst.afw.cameraGeom.FOCAL_PLANE` 

302 transformDict : `dict` 

303 A dict of lsst.afw.cameraGeom.CameraSys : 

304 `lsst.afw.geom.TransformPoint2ToPoint2` 

305 amplifierDict : `dict` [`str`, `lsst.afw.cameraGeom.Amplifier.Builder` ] 

306 A dictionary of detector name and amplifier builders. 

307 pupilFactoryClass : `type` [ `lsst.default afw.cameraGeom.PupilFactory`], \ 

308 optional 

309 Class to attach to camera. 

310 focalPlaneParity : `bool`, optional 

311 If `True`, the X axis is flipped between the FOCAL_PLANE and 

312 FIELD_ANGLE coordinate systems. 

313 

314 Returns 

315 ------- 

316 camera : `lsst.afw.cameraGeom.Camera` 

317 New Camera instance. 

318 

319 Notes 

320 ----- 

321 Copied from `lsst.afw.cameraGeom.cameraFactory` with permission and 

322 encouragement from Jim Bosch. 

323 """ 

324 # nativeSys=FOCAL_PLANE seems to be assumed in various places in this file 

325 # (e.g. the definition of TAN_PIXELS), despite CameraConfig providing the 

326 # illusion that it's configurable. 

327 # Note that we can't actually get rid of the nativeSys config option 

328 # without breaking lots of on-disk camera configs. 

329 assert nativeSys == cameraGeom.FOCAL_PLANE, "Cameras with nativeSys != FOCAL_PLANE are not supported." 

330 

331 focalPlaneToField = transformDict[cameraGeom.FIELD_ANGLE] 

332 

333 cameraBuilder = Camera.Builder(cameraName) 

334 cameraBuilder.setPupilFactoryClass(pupilFactoryClass) 

335 cameraBuilder.setFocalPlaneParity(focalPlaneParity) 

336 

337 # Ensure all transforms in the camera transform dict are included. 

338 for toSys, transform in transformDict.items(): 

339 cameraBuilder.setTransformFromFocalPlaneTo(toSys, transform) 

340 

341 for detectorConfig in detectorConfigList: 

342 # This should build all detector pixel -> focalPlane transforms. 

343 cameraGeom.addDetectorBuilderFromConfig( 

344 cameraBuilder, detectorConfig, amplifierDict[detectorConfig.name], focalPlaneToField 

345 ) 

346 

347 # For reasons I don't understand, some obs_ packages (e.g. HSC) set 

348 # nativeSys to None for their detectors (which doesn't seem to be 

349 # permitted by the config class!), but they really mean PIXELS. For 

350 # backwards compatibility we use that as the default... 

351 detectorNativeSys = detectorConfig.transformDict.nativeSys 

352 detectorNativeSys = ( 

353 cameraGeom.PIXELS if detectorNativeSys is None else cameraGeom.CameraSysPrefix(detectorNativeSys) 

354 ) 

355 

356 # ...well, actually, it seems that we've always assumed down in C++ 

357 # that the answer is always PIXELS without ever checking that it is. 

358 # So let's assert that it is, since there are hints all over this file 

359 # (e.g. the definition of TAN_PIXELS) that other parts of the codebase 

360 # have regularly made that assumption as well. Note that we can't 

361 # actually get rid of the nativeSys config option without breaking 

362 # lots of on-disk camera configs. 

363 assert detectorNativeSys == cameraGeom.PIXELS, "Detectors with nativeSys != PIXELS are not supported." 

364 detectorNativeSys = cameraGeom.CameraSys(detectorNativeSys, detectorConfig.name) 

365 

366 return cameraBuilder.finish()