LSSTApplications  10.0-2-g4f67435,11.0.rc2+1,11.0.rc2+12,11.0.rc2+3,11.0.rc2+4,11.0.rc2+5,11.0.rc2+6,11.0.rc2+7,11.0.rc2+8
LSSTDataManagementBasePackage
visualization.py
Go to the documentation of this file.
1 #
2 # LSST Data Management System
3 # Copyright 2008, 2009, 2010, 2011, 2012 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 
23 from matplotlib import pyplot
24 
25 import lsst.afw.geom
26 
27 def plotObservations(catalog, wcs):
28  """Plot the bounding boxes of an observation catalog (see MockCoaddTask.buildObservationCatalog)
29  using matplotlib, in the coordinates defined by the given Wcs (usually a skymap Wcs).
30  """
31  for record in catalog:
32  box = lsst.afw.geom.Box2D(record.getBBox())
33  x = []
34  y = []
35  iWcs = record.getWcs()
36  for xi, yi in box.getCorners():
37  try:
38  coord = iWcs.pixelToSky(xi, yi)
39  xo, yo = wcs.skyToPixel(coord)
40  x.append(xo)
41  y.append(yo)
42  except:
43  print "WARNING: point %d, %d failed" % (xi, yi)
44  pyplot.fill(x, y, facecolor='r', alpha=0.1, edgecolor=None)
45 
46 def plotPatches(tractInfo):
47  """Plot the patches in a skymap tract using matplotlib.
48  """
49  nPatchX, nPatchY = tractInfo.getNumPatches()
50  for iPatchX in range(nPatchX):
51  for iPatchY in range(nPatchY):
52  patchInfo = tractInfo.getPatchInfo((iPatchX, iPatchY))
53  xp1, yp1 = zip(*patchInfo.getOuterBBox().getCorners())
54  xp2, yp2 = zip(*patchInfo.getInnerBBox().getCorners())
55  pyplot.fill(xp1, yp1, fill=False, edgecolor='g', linestyle='dashed')
56  pyplot.fill(xp2, yp2, fill=False, edgecolor='g')
57 
58 def plotTruth(catalog, wcs):
59  """Plot the objects in a truth catalog as dots using matplotlib, in the coordinate
60  system defined by the given Wcs.
61  """
62  xp = []
63  yp = []
64  for record in catalog:
65  x, y = wcs.skyToPixel(record.getCoord())
66  xp.append(x)
67  yp.append(y)
68  pyplot.plot(xp, yp, 'k+')
69 
70 def displayImages(root):
71  """Display coadd images with DS9 in different frames, with the bounding boxes of the
72  observations that went into them overlayed.
73  """
76  butler = lsst.daf.persistence.Butler(root=root)
77  skyMap = butler.get("deepCoadd_skyMap")
78  tractInfo = skyMap[0]
79  task = lsst.pipe.tasks.mocks.MockCoaddTask()
80  coadds = [patchRef.get("deepCoadd", immediate=True)
81  for patchRef in task.iterPatchRefs(butler, tractInfo)]
82  for n, coadd in enumerate(coadds):
83  lsst.afw.display.ds9.mtv(coadd, frame=n+1)
84  for n, coadd in enumerate(coadds):
86  return butler
87 
88 def makePlots(root):
89  """Convenience function to make all matplotlib plots.
90  """
93  butler = lsst.daf.persistence.Butler(root=root)
94  skyMap = butler.get("deepCoadd_skyMap")
95  observations = butler.get("observations", tract=0)
96  truth = butler.get("truth", tract=0)
97  tractInfo = skyMap[0]
98  plotPatches(tractInfo)
99  plotObservations(observations, tractInfo.getWcs())
100  plotTruth(truth, tractInfo.getWcs())
101  pyplot.axis("scaled")
102  pyplot.show()
103  return butler
A floating-point coordinate rectangle geometry.
Definition: Box.h:271