LSSTApplications  18.1.0
LSSTDataManagementBasePackage
yamlCamera.py
Go to the documentation of this file.
1 #
2 # LSST Data Management System
3 # Copyright 2017 LSST Corporation.
4 #
5 # This product includes software developed by the
6 # LSST Project (http://www.lsst.org/).
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 <http://www.lsstcorp.org/LegalNotices/>.
21 #
22 import yaml
23 
24 import numpy as np
25 import lsst.afw.cameraGeom as cameraGeom
26 import lsst.geom as geom
27 import lsst.afw.geom as afwGeom
28 from lsst.afw.table import AmpInfoCatalog, AmpInfoTable
29 from lsst.afw.cameraGeom.cameraFactory import makeDetectorData
30 
31 __all__ = ["makeCamera"]
32 
33 
34 def makeCamera(cameraFile):
35  """An imaging camera (e.g. the LSST 3Gpix camera)
36 
37  Parameters
38  ----------
39  cameraFile : `str`
40  Camera description YAML file.
41 
42  Returns
43  -------
44  camera : `lsst.afw.cameraGeom.Camera`
45  The desired Camera
46  """
47 
48  with open(cameraFile) as fd:
49  cameraParams = yaml.load(fd, Loader=yaml.CLoader)
50 
51  cameraName = cameraParams["name"]
52 
53  #
54  # Handle distortion models.
55  #
56  plateScale = geom.Angle(cameraParams["plateScale"], geom.arcseconds)
57  nativeSys = cameraGeom.CameraSys(cameraParams["transforms"].pop("nativeSys"))
58  transforms = makeTransformDict(nativeSys, cameraParams["transforms"], plateScale)
59 
60  ccdParams = cameraParams["CCDs"]
61  detectorConfigList = makeDetectorConfigList(ccdParams)
62 
63  ampInfoCatDict = {}
64  for ccdName, ccdValues in ccdParams.items():
65  ampInfoCatDict[ccdName] = makeAmpInfoCatalog(ccdValues)
66 
67  return makeCameraFromCatalogs(cameraName, detectorConfigList, nativeSys, transforms, ampInfoCatDict)
68 
69 
70 def makeDetectorConfigList(ccdParams):
71  """Make a list of detector configs
72 
73  Returns
74  -------
75  detectorConfig : `list` of `lsst.afw.cameraGeom.DetectorConfig`
76  A list of detector configs.
77  """
78  detectorConfigs = []
79  for name, ccd in ccdParams.items():
80  detectorConfig = cameraGeom.DetectorConfig()
81  detectorConfigs.append(detectorConfig)
82 
83  detectorConfig.name = name
84  detectorConfig.id = ccd['id']
85  detectorConfig.serial = ccd['serial']
86  detectorConfig.detectorType = ccd['detectorType']
87  if 'physicalType' in ccd:
88  detectorConfig.physicalType = ccd['physicalType']
89  # This is the orientation we need to put the serial direction along the x-axis
90  detectorConfig.bbox_x0, detectorConfig.bbox_y0 = ccd['bbox'][0]
91  detectorConfig.bbox_x1, detectorConfig.bbox_y1 = ccd['bbox'][1]
92  detectorConfig.pixelSize_x, detectorConfig.pixelSize_y = ccd['pixelSize']
93  detectorConfig.transformDict.nativeSys = ccd['transformDict']['nativeSys']
94  transforms = ccd['transformDict']['transforms']
95  detectorConfig.transformDict.transforms = None if transforms == 'None' else transforms
96  detectorConfig.refpos_x, detectorConfig.refpos_y = ccd['refpos']
97  detectorConfig.offset_x, detectorConfig.offset_y = ccd['offset']
98  detectorConfig.transposeDetector = ccd['transposeDetector']
99  detectorConfig.pitchDeg = ccd['pitch']
100  detectorConfig.yawDeg = ccd['yaw']
101  detectorConfig.rollDeg = ccd['roll']
102  if 'crosstalk' in ccd:
103  detectorConfig.crosstalk = ccd['crosstalk']
104 
105  return detectorConfigs
106 
107 
109  """Construct an amplifier info catalog
110  """
111  # Much of this will need to be filled in when we know it.
112  assert len(ccd['amplifiers']) > 0
113  amp = list(ccd['amplifiers'].values())[0]
114 
115  rawBBox = makeBBoxFromList(amp['rawBBox']) # total in file
116  xRawExtent, yRawExtent = rawBBox.getDimensions()
117 
118  from lsst.afw.table import LL, LR, UL, UR
119  readCorners = dict(LL=LL, LR=LR, UL=UL, UR=UR)
120 
121  schema = AmpInfoTable.makeMinimalSchema()
122 
123  linThreshKey = schema.addField('linearityThreshold', type=float)
124  linMaxKey = schema.addField('linearityMaximum', type=float)
125  linUnitsKey = schema.addField('linearityUnits', type=str, size=9)
126  hduKey = schema.addField('hdu', type=np.int32)
127  # end placeholder
128 
129  ampCatalog = AmpInfoCatalog(schema)
130  for name, amp in sorted(ccd['amplifiers'].items(), key=lambda x: x[1]['hdu']):
131  record = ampCatalog.addNew()
132  record.setName(name)
133  record.set(hduKey, amp['hdu'])
134 
135  ix, iy = amp['ixy']
136  perAmpData = amp['perAmpData']
137  if perAmpData:
138  x0, y0 = 0, 0 # origin of data within each amp image
139  else:
140  x0, y0 = ix*xRawExtent, iy*yRawExtent
141 
142  rawDataBBox = makeBBoxFromList(amp['rawDataBBox']) # Photosensitive area
143  xDataExtent, yDataExtent = rawDataBBox.getDimensions()
144  record.setBBox(afwGeom.BoxI(
145  afwGeom.PointI(ix*xDataExtent, iy*yDataExtent), rawDataBBox.getDimensions()))
146 
147  rawBBox = makeBBoxFromList(amp['rawBBox'])
148  rawBBox.shift(afwGeom.ExtentI(x0, y0))
149  record.setRawBBox(rawBBox)
150 
151  rawDataBBox = makeBBoxFromList(amp['rawDataBBox'])
152  rawDataBBox.shift(afwGeom.ExtentI(x0, y0))
153  record.setRawDataBBox(rawDataBBox)
154 
155  rawSerialOverscanBBox = makeBBoxFromList(amp['rawSerialOverscanBBox'])
156  rawSerialOverscanBBox.shift(afwGeom.ExtentI(x0, y0))
157  record.setRawHorizontalOverscanBBox(rawSerialOverscanBBox)
158 
159  rawParallelOverscanBBox = makeBBoxFromList(amp['rawParallelOverscanBBox'])
160  rawParallelOverscanBBox.shift(afwGeom.ExtentI(x0, y0))
161  record.setRawVerticalOverscanBBox(rawParallelOverscanBBox)
162 
163  rawSerialPrescanBBox = makeBBoxFromList(amp['rawSerialPrescanBBox'])
164  rawSerialPrescanBBox.shift(afwGeom.ExtentI(x0, y0))
165  record.setRawPrescanBBox(rawSerialPrescanBBox)
166 
167  if perAmpData:
168  record.setRawXYOffset(afwGeom.Extent2I(ix*xRawExtent, iy*yRawExtent))
169  else:
170  record.setRawXYOffset(afwGeom.Extent2I(0, 0))
171 
172  record.setReadoutCorner(readCorners[amp['readCorner']])
173  record.setGain(amp['gain'])
174  record.setReadNoise(amp['readNoise'])
175  record.setSaturation(amp['saturation'])
176  record.setHasRawInfo(True)
177  # flip data when assembling if needs be (e.g. data from the serial at the top of a CCD)
178  flipX, flipY = amp.get("flipXY")
179 
180  record.setRawFlipX(flipX)
181  record.setRawFlipY(flipY)
182  # linearity placeholder stuff
183  record.setLinearityCoeffs([float(val) for val in amp['linearityCoeffs']])
184  record.setLinearityType(amp['linearityType'])
185  record.set(linThreshKey, float(amp['linearityThreshold']))
186  record.set(linMaxKey, float(amp['linearityMax']))
187  record.set(linUnitsKey, "DN")
188  return ampCatalog
189 
190 
191 def makeBBoxFromList(ylist):
192  """Given a list [(x0, y0), (xsize, ysize)], probably from a yaml file,
193  return a BoxI
194  """
195  (x0, y0), (xsize, ysize) = ylist
196  return afwGeom.BoxI(afwGeom.PointI(x0, y0), afwGeom.ExtentI(xsize, ysize))
197 
198 
199 def makeTransformDict(nativeSys, transformDict, plateScale):
200  """Make a dictionary of TransformPoint2ToPoint2s from yaml, mapping from nativeSys
201 
202  Parameters
203  ----------
204  nativeSys : `lsst.afw.cameraGeom.CameraSys`
205  transformDict : `dict`
206  A dict specifying parameters of transforms; keys are camera system names.
207  plateScale : `lsst.geom.Angle`
208  The size of a pixel in angular units/mm (e.g. 20 arcsec/mm for LSST)
209 
210  Returns
211  -------
212  transforms : `dict`
213  A dict of `lsst.afw.cameraGeom.CameraSys` : `lsst.afw.geom.TransformPoint2ToPoint2`
214 
215  The resulting dict's keys are `~lsst.afw.cameraGeom.CameraSys`,
216  and the values are Transforms *from* NativeSys *to* CameraSys
217  """
218  # As other comments note this is required, and this is one function where it's assumed
219  assert nativeSys == cameraGeom.FOCAL_PLANE, "Cameras with nativeSys != FOCAL_PLANE are not supported."
220 
221  resMap = dict()
222 
223  for key, transform in transformDict.items():
224  transformType = transform["transformType"]
225  knownTransformTypes = ["affine", "radial"]
226  if transformType not in knownTransformTypes:
227  raise RuntimeError("Saw unknown transform type for %s: %s (known types are: [%s])" % (
228  key, transform["transformType"], ", ".join(knownTransformTypes)))
229 
230  if transformType == "affine":
231  affine = geom.AffineTransform(np.array(transform["linear"]),
232  np.array(transform["translation"]))
233 
234  transform = afwGeom.makeTransform(affine)
235  elif transformType == "radial":
236  # radial coefficients of the form [0, 1 (no units), C2 (rad), usually 0, C3 (rad^2), ...]
237  # Radial distortion is modeled as a radial polynomial that converts from focal plane radius
238  # (in mm) to field angle (in radians). The provided coefficients are divided by the plate
239  # scale (in radians/mm) meaning that C1 is always 1.
240  radialCoeffs = np.array(transform["coeffs"])
241 
242  radialCoeffs *= plateScale.asRadians()
243  transform = afwGeom.makeRadialTransform(radialCoeffs)
244  else:
245  raise RuntimeError("Impossible condition \"%s\" is not in: [%s])" % (
246  transform["transformType"], ", ".join(knownTransformTypes)))
247 
248  resMap[cameraGeom.CameraSys(key)] = transform
249 
250  return resMap
251 
252 
253 def makeCameraFromCatalogs(cameraName, detectorConfigList, nativeSys, transformDict, ampInfoCatDict,
254  pupilFactoryClass=cameraGeom.pupil.PupilFactory):
255  """Construct a Camera instance from a dictionary of
256  detector name : `lsst.afw.table.ampInfo.AmpInfoCatalog`
257 
258  Parameters
259  ----------
260  cameraName : `str`
261  The name of the camera
262  detectorConfig : `list`
263  A list of `lsst.afw.cameraGeom.cameraConfig.DetectorConfig`
264  nativeSys : `lsst.afw.cameraGeom.CameraSys`
265  The native transformation type; must be `lsst.afw.cameraGeom.FOCAL_PLANE`
266  transformDict : `dict`
267  A dict of lsst.afw.cameraGeom.CameraSys : `lsst.afw.geom.TransformPoint2ToPoint2`
268  ampInfoCatDict : `dict`
269  A dictionary of detector name :
270  `lsst.afw.table.ampInfo.AmpInfoCatalog`
271  pupilFactoryClass : `type`, optional
272  Class to attach to camera;
273  `lsst.default afw.cameraGeom.PupilFactory`
274 
275  Returns
276  -------
277  camera : `lsst.afw.cameraGeom.Camera`
278  New Camera instance.
279 
280  Notes
281  ------
282  Copied from `lsst.afw.cameraGeom.cameraFactory` with permission and encouragement
283  from Jim Bosch
284  """
285 
286  # nativeSys=FOCAL_PLANE seems to be assumed in various places in this file
287  # (e.g. the definition of TAN_PIXELS), despite CameraConfig providing the
288  # illusion that it's configurable.
289  # Note that we can't actually get rid of the nativeSys config option
290  # without breaking lots of on-disk camera configs.
291  assert nativeSys == cameraGeom.FOCAL_PLANE, "Cameras with nativeSys != FOCAL_PLANE are not supported."
292 
293  focalPlaneToField = transformDict[cameraGeom.FIELD_ANGLE]
294  transformMapBuilder = cameraGeom.TransformMap.Builder(nativeSys)
295  transformMapBuilder.connect(transformDict)
296 
297  # First pass: build a list of all Detector ctor kwargs, minus the
298  # transformMap (which needs information from all Detectors).
299  detectorData = []
300  for detectorConfig in detectorConfigList:
301 
302  # Get kwargs that could be used to construct each Detector
303  # if we didn't care about giving each of them access to
304  # all of the transforms.
305  thisDetectorData = makeDetectorData(
306  detectorConfig=detectorConfig,
307  ampInfoCatalog=ampInfoCatDict[detectorConfig.name],
308  focalPlaneToField=focalPlaneToField,
309  )
310 
311  # Pull the transforms dictionary out of the data dict; we'll replace
312  # it with a TransformMap argument later.
313  thisDetectorTransforms = thisDetectorData.pop("transforms")
314 
315  # Save the rest of the Detector data dictionary for later
316  detectorData.append(thisDetectorData)
317 
318  # For reasons I don't understand, some obs_ packages (e.g. HSC) set
319  # nativeSys to None for their detectors (which doesn't seem to be
320  # permitted by the config class!), but they really mean PIXELS. For
321  # backwards compatibility we use that as the default...
322  detectorNativeSys = detectorConfig.transformDict.nativeSys
323  detectorNativeSys = (cameraGeom.PIXELS if detectorNativeSys is None else
324  cameraGeom.CameraSysPrefix(detectorNativeSys))
325 
326  # ...well, actually, it seems that we've always assumed down in C++
327  # that the answer is always PIXELS without ever checking that it is.
328  # So let's assert that it is, since there are hints all over this file
329  # (e.g. the definition of TAN_PIXELS) that other parts of the codebase
330  # have regularly made that assumption as well. Note that we can't
331  # actually get rid of the nativeSys config option without breaking
332  # lots of on-disk camera configs.
333  assert detectorNativeSys == cameraGeom.PIXELS, \
334  "Detectors with nativeSys != PIXELS are not supported."
335  detectorNativeSys = cameraGeom.CameraSys(detectorNativeSys, detectorConfig.name)
336 
337  # Add this detector's transform dict to the shared TransformMapBuilder
338  transformMapBuilder.connect(detectorNativeSys, thisDetectorTransforms)
339 
340  # Now that we've collected all of the Transforms, we can finally build the
341  # (immutable) TransformMap.
342  transformMap = transformMapBuilder.build()
343 
344  # Second pass through the detectorConfigs: actually make Detector instances
345  detectorList = [cameraGeom.Detector(transformMap=transformMap, **kw) for kw in detectorData]
346 
347  return cameraGeom.Camera(cameraName, detectorList, transformMap, pupilFactoryClass)
def makeBBoxFromList(ylist)
Definition: yamlCamera.py:191
def makeDetectorData(detectorConfig, ampInfoCatalog, focalPlaneToField)
An affine coordinate transformation consisting of a linear transformation and an offset.
def makeCamera(cameraFile)
Definition: yamlCamera.py:34
A class representing an angle.
Definition: Angle.h:127
def makeDetectorConfigList(ccdParams)
Definition: yamlCamera.py:70
def makeTransformDict(nativeSys, transformDict, plateScale)
Definition: yamlCamera.py:199
std::shared_ptr< TransformPoint2ToPoint2 > makeRadialTransform(std::vector< double > const &forwardCoeffs, std::vector< double > const &inverseCoeffs)
A purely radial polynomial distortion.
CatalogT< AmpInfoRecord > AmpInfoCatalog
Definition: fwd.h:97
def makeCameraFromCatalogs(cameraName, detectorConfigList, nativeSys, transformDict, ampInfoCatDict, pupilFactoryClass=cameraGeom.pupil.PupilFactory)
Definition: yamlCamera.py:254
std::vector< SchemaItem< Flag > > * items
An integer coordinate rectangle.
Definition: Box.h:54
daf::base::PropertyList * list
Definition: fits.cc:885
std::shared_ptr< TransformPoint2ToPoint2 > makeTransform(lsst::geom::AffineTransform const &affine)
Wrap an lsst::geom::AffineTransform as a Transform.