LSSTApplications  10.0+286,10.0+36,10.0+46,10.0-2-g4f67435,10.1+152,10.1+37,11.0,11.0+1,11.0-1-g47edd16,11.0-1-g60db491,11.0-1-g7418c06,11.0-2-g04d2804,11.0-2-g68503cd,11.0-2-g818369d,11.0-2-gb8b8ce7
LSSTDataManagementBasePackage
Public Member Functions | Public Attributes | Private Member Functions | List of all members
lsst.skypix.quadsphere.QuadSpherePixelization Class Reference
Inheritance diagram for lsst.skypix.quadsphere.QuadSpherePixelization:

Public Member Functions

def __init__
 
def getGeometry
 
def getCenter
 
def getNeighbors
 
def intersect
 
def __len__
 
def __iter__
 
def __repr__
 
def pixel
 
def id
 
def coords
 

Public Attributes

 resolution
 
 padding
 
 center
 
 x
 
 y
 
 xrot
 
 yrot
 
 xplane
 
 yplane
 
 cornerNeighbors
 

Private Member Functions

def _fiducialXPlane
 
def _fiducialYPlane
 
def _intersect
 

Detailed Description

A quad-sphere based pixelization of the unit sphere. Each of the 6
root pixels are divided into an R by R grid of sky-pixels, where R
is the resolution (a positive integer). Sky-pixel identifiers are
contiguous integers in range [0, 6 * R**2), and are ordered such that
enumerating the pixels 0, 1, ..., 6 * R**2 results in a roughly spiral
traversal (from south to north pole) of the unit sphere. In particular,
iterating over pixels in this order guarantees that when arriving at
pixel P:

  - all pixels with id P' < P have been visited
  - all neighbors of pixels with id P' < P - 8R have been visited

Each pixel may optionally be expanded by a padding angle A, such that
points on the unit sphere outside of a pixel but within angular separation
A of the fiducial pixel boundaries are also considered to be part of that
pixel. Pixel (and padded pixel) edges are great circles on the unit
sphere.

This class provides methods for computing fiducial or padded geometric
representations of a sky-pixel - spherical convex polygons in both
cases. It also supports finding the list of sky-pixels intersecting
an input spherical convex polygon, as well as determining the neighbors
and center of a sky-pixel.

TODO: Right now, pixels are defined by a simple tangent plane projection
of each cube face onto the sphere, but pixel boundaries could be adjusted
to produce equal area pixels.

Definition at line 44 of file quadsphere.py.

Constructor & Destructor Documentation

def lsst.skypix.quadsphere.QuadSpherePixelization.__init__ (   self,
  resolution,
  paddingRad 
) 
Creates a new quad-sphere sky pixelisation.

resolution: the number of pixels along the x and y axes of each root
    pixel.

paddingRad: the angular separation (rad) by which fiducial sky-pixels
    are to be padded.

Definition at line 73 of file quadsphere.py.

73 
74  def __init__(self, resolution, paddingRad):
75  """Creates a new quad-sphere sky pixelisation.
76 
77  resolution: the number of pixels along the x and y axes of each root
78  pixel.
79 
80  paddingRad: the angular separation (rad) by which fiducial sky-pixels
81  are to be padded.
82  """
83  if not isinstance(resolution, (int, long)):
84  raise TypeError(
85  'Quad-sphere resolution parameter must be an int or long')
86  if resolution < 3 or resolution > 16384:
87  raise RuntimeError(
88  'Quad-sphere resolution must be in range [3, 16384]')
89  if not isinstance(paddingRad, float):
90  raise TypeError(
91  'Quad-sphere pixel padding angle must be a float')
92  if paddingRad < 0.0 or paddingRad >= math.pi * 0.25:
93  raise RuntimeError(
94  'Quad-sphere pixel padding radius must be in range [0, 45) deg')
95  R = resolution
96  self.resolution = R
97  self.padding = paddingRad
98  x, y, z = (1.0, 0.0, 0.0), (0.0, 1.0, 0.0), (0.0, 0.0, 1.0)
99  nx, ny, nz = (-1.0, 0.0, 0.0), (0.0, -1.0, 0.0), (0.0, 0.0, -1.0)
100  # center of each root pixel
101  self.center = [z, x, y, nx, ny, nz]
102  # x basis vector of each root pixel
103  self.x = [ny, y, nx, ny, x, ny]
104  # y basis vector of each root pixel
105  self.y = [x, z, z, z, z, nx]
106  # functions for rotating vectors in direction of increasing x (per root pixel)
107  self.xrot = [_rotX, _rotZ, _rotZ, _rotZ, _rotZ, _rotNX]
108  # functions for rotating vectors in direction of increasing y (per root pixel)
109  self.yrot = [_rotY, _rotNY, _rotX, _rotY, _rotNX, _rotY]
110  # Compute fiducial and splitting x/y planes for each root pixel
111  self.xplane = []
112  self.yplane = []
113  sp = math.sin(self.padding)
114  for root in xrange(6):
115  xplanes = []
116  yplanes = []
117  c, x, y = self.center[root], self.x[root], self.y[root]
118  for i in xrange(R + 1):
119  xfp = self._fiducialXPlane(root, i)
120  yfp = self._fiducialYPlane(root, i)
121  f = 2.0 * float(i) / float(R) - 1.0
122  thetaX = math.radians(0.5 * geom.cartesianAngularSep(
123  (c[0] + f * x[0] + y[0],
124  c[1] + f * x[1] + y[1],
125  c[2] + f * x[2] + y[2]),
126  (c[0] + f * x[0] - y[0],
127  c[1] + f * x[1] - y[1],
128  c[2] + f * x[2] - y[2])))
129  thetaY = math.radians(0.5 * geom.cartesianAngularSep(
130  (c[0] + x[0] + f * y[0],
131  c[1] + x[1] + f * y[1],
132  c[2] + x[2] + f * y[2]),
133  (c[0] - x[0] + f * y[0],
134  c[1] - x[1] + f * y[1],
135  c[2] - x[2] + f * y[2])))
136  sinrx = sp / math.cos(thetaX)
137  sinry = sp / math.cos(thetaY)
138  cosrx = math.sqrt(1.0 - sinrx * sinrx)
139  cosry = math.sqrt(1.0 - sinry * sinry)
140  if i == 0:
141  xlp = None
142  ybp = None
143  else:
144  xlp = self.xrot[root](xfp, sinrx, cosrx)
145  ybp = self.yrot[root](yfp, sinry, cosry)
146  xlp = (-xlp[0], -xlp[1], -xlp[2])
147  ybp = (-ybp[0], -ybp[1], -ybp[2])
148  if i == R:
149  xrp = None
150  ytp = None
151  else:
152  xrp = self.xrot[root](xfp, -sinrx, cosrx)
153  ytp = self.yrot[root](yfp, -sinry, cosry)
154  xplanes.append((xfp, xlp, xrp))
155  yplanes.append((yfp, ybp, ytp))
156  self.xplane.append(xplanes)
157  self.yplane.append(yplanes)
158  # Corner pixel neighbors.
159  self.cornerNeighbors = (
160  # root pixel 0
161  (
162  # x, y = 0, 0
163  (self.id(0, 1, 0),
164  self.id(0, 1, 1),
165  self.id(0, 0, 1),
166  self.id(2, R - 1, R - 1),
167  self.id(2, R - 2, R - 1),
168  self.id(3, 0, R - 1),
169  self.id(3, 1, R - 1)),
170  # x, y = 0, R-1
171  (self.id(0, 1, R - 1),
172  self.id(0, 1, R - 2),
173  self.id(0, 0, R - 2),
174  self.id(1, R - 1, R - 1),
175  self.id(1, R - 2, R - 1),
176  self.id(2, 0, R - 1),
177  self.id(2, 1, R - 1)),
178  # x, y = R-1, 0
179  (self.id(0, R - 2, 0),
180  self.id(0, R - 2, 1),
181  self.id(0, R - 1, 1),
182  self.id(3, R - 2, R - 1),
183  self.id(3, R - 1, R - 1),
184  self.id(4, 0, R - 1),
185  self.id(4, 1, R - 1)),
186  # x, y = R-1, R-1
187  (self.id(0, R - 2, R - 1),
188  self.id(0, R - 2, R - 2),
189  self.id(0, R - 1, R - 2),
190  self.id(1, 0, R - 1),
191  self.id(1, 1, R - 1),
192  self.id(4, R - 2, R - 1),
193  self.id(4, R - 1, R - 1))
194  ),
195  # root pixel 1
196  (
197  # x, y = 0, 0
198  (self.id(1, 1, 0),
199  self.id(1, 1, 1),
200  self.id(1, 0, 1),
201  self.id(4, R - 1, 0),
202  self.id(4, R - 1, 1),
203  self.id(5, R - 2, 0),
204  self.id(5, R - 1, 0)),
205  # x, y = 0, R-1
206  (self.id(1, 1, R - 1),
207  self.id(1, 1, R - 2),
208  self.id(1, 0, R - 2),
209  self.id(4, R - 1, R - 2),
210  self.id(4, R - 1, R - 1),
211  self.id(0, R - 2, R - 1),
212  self.id(0, R - 1, R - 1)),
213  # x, y = R-1, 0
214  (self.id(1, R - 2, 0),
215  self.id(1, R - 2, 1),
216  self.id(1, R - 1, 1),
217  self.id(2, 0, 0),
218  self.id(2, 0, 1),
219  self.id(5, 0, 0),
220  self.id(5, 1, 0)),
221  # x, y = R-1, R-1
222  (self.id(1, R - 2, R - 1),
223  self.id(1, R - 2, R - 2),
224  self.id(1, R - 1, R - 2),
225  self.id(2, 0, R - 2),
226  self.id(2, 0, R - 1),
227  self.id(0, 0, R - 1),
228  self.id(0, 1, R - 1))
229  ),
230  # root pixel 2
231  (
232  # x, y = 0, 0
233  (self.id(2, 1, 0),
234  self.id(2, 1, 1),
235  self.id(2, 0, 1),
236  self.id(1, R - 1, 0),
237  self.id(1, R - 1, 1),
238  self.id(5, 0, 0),
239  self.id(5, 0, 1)),
240  # x, y = 0, R-1
241  (self.id(2, 1, R - 1),
242  self.id(2, 1, R - 2),
243  self.id(2, 0, R - 2),
244  self.id(1, R - 1, R - 2),
245  self.id(1, R - 1, R - 1),
246  self.id(0, 0, R - 2),
247  self.id(0, 0, R - 1)),
248  # x, y = R-1, 0
249  (self.id(2, R - 2, 0),
250  self.id(2, R - 2, 1),
251  self.id(2, R - 1, 1),
252  self.id(3, 0, 0),
253  self.id(3, 0, 1),
254  self.id(5, 0, R - 2),
255  self.id(5, 0, R - 1)),
256  # x, y = R-1, R-1
257  (self.id(2, R - 2, R - 1),
258  self.id(2, R - 2, R - 2),
259  self.id(2, R - 1, R - 2),
260  self.id(3, 0, R - 2),
261  self.id(3, 0, R - 1),
262  self.id(0, 0, 0),
263  self.id(0, 0, 1))
264  ),
265  # root pixel 3
266  (
267  # x, y = 0, 0
268  (self.id(3, 1, 0),
269  self.id(3, 1, 1),
270  self.id(3, 0, 1),
271  self.id(2, R - 1, 0),
272  self.id(2, R - 1, 1),
273  self.id(5, 0, R - 1),
274  self.id(5, 1, R - 1)),
275  # x, y = 0, R-1
276  (self.id(3, 1, R - 1),
277  self.id(3, 1, R - 2),
278  self.id(3, 0, R - 2),
279  self.id(2, R - 1, R - 2),
280  self.id(2, R - 1, R - 1),
281  self.id(0, 0, 0),
282  self.id(0, 1, 0)),
283  # x, y = R-1, 0
284  (self.id(3, R - 2, 0),
285  self.id(3, R - 2, 1),
286  self.id(3, R - 1, 1),
287  self.id(4, 0, 0),
288  self.id(4, 0, 1),
289  self.id(5, R - 2, R - 1),
290  self.id(5, R - 1, R - 1)),
291  # x, y = R-1, R-1
292  (self.id(3, R - 2, R - 1),
293  self.id(3, R - 2, R - 2),
294  self.id(3, R - 1, R - 2),
295  self.id(4, 0, R - 2),
296  self.id(4, 0, R - 1),
297  self.id(0, R - 2, 0),
298  self.id(0, R - 1, 0))
299  ),
300  # root pixel 4
301  (
302  # x, y = 0, 0
303  (self.id(4, 1, 0),
304  self.id(4, 1, 1),
305  self.id(4, 0, 1),
306  self.id(3, R - 1, 0),
307  self.id(3, R - 1, 1),
308  self.id(5, R - 1, R - 2),
309  self.id(5, R - 1, R - 1)),
310  # x, y = 0, R-1
311  (self.id(4, 1, R - 1),
312  self.id(4, 1, R - 2),
313  self.id(4, 0, R - 2),
314  self.id(3, R - 1, R - 2),
315  self.id(3, R - 1, R - 1),
316  self.id(0, R - 1, 0),
317  self.id(0, R - 1, 1)),
318  # x, y = R-1, 0
319  (self.id(4, R - 2, 0),
320  self.id(4, R - 2, 1),
321  self.id(4, R - 1, 1),
322  self.id(1, 0, 0),
323  self.id(1, 0, 1),
324  self.id(5, R - 1, 0),
325  self.id(5, R - 1, 1)),
326  # x, y = R-1, R-1
327  (self.id(4, R - 2, R - 1),
328  self.id(4, R - 2, R - 2),
329  self.id(4, R - 1, R - 2),
330  self.id(1, 0, R - 2),
331  self.id(1, 0, R - 1),
332  self.id(0, R - 1, R - 2),
333  self.id(0, R - 1, R - 1))
334  ),
335  # root pixel 5
336  (
337  # x, y = 0, 0
338  (self.id(5, 1, 0),
339  self.id(5, 1, 1),
340  self.id(5, 0, 1),
341  self.id(1, R - 2, 0),
342  self.id(1, R - 1, 0),
343  self.id(2, 0, 0),
344  self.id(2, 1, 0)),
345  # x, y = 0, R-1
346  (self.id(5, 1, R - 1),
347  self.id(5, 1, R - 2),
348  self.id(5, 0, R - 2),
349  self.id(2, R - 2, 0),
350  self.id(2, R - 1, 0),
351  self.id(3, 0, 0),
352  self.id(3, 1, 0)),
353  # x, y = R-1, 0
354  (self.id(5, R - 2, 0),
355  self.id(5, R - 2, 1),
356  self.id(5, R - 1, 1),
357  self.id(1, 0, 0),
358  self.id(1, 1, 0),
359  self.id(4, R - 2, 0),
360  self.id(4, R - 1, 0)),
361  # x, y = R-1, R-1
362  (self.id(5, R - 2, R - 1),
363  self.id(5, R - 2, R - 2),
364  self.id(5, R - 1, R - 2),
365  self.id(3, R - 2, 0),
366  self.id(3, R - 1, 0),
367  self.id(4, 0, 0),
368  self.id(4, 1, 0))
369  ),
370  )

Member Function Documentation

def lsst.skypix.quadsphere.QuadSpherePixelization.__iter__ (   self) 
Returns an iterator over the sky-pixel ids of all the pixels
in this pixelization.

Definition at line 583 of file quadsphere.py.

583 
584  def __iter__(self):
585  """Returns an iterator over the sky-pixel ids of all the pixels
586  in this pixelization.
587  """
588  return iter(xrange(len(self)))
iter
int iter
Definition: GaussianCentroid.cc:73
def lsst.skypix.quadsphere.QuadSpherePixelization.__len__ (   self) 
Returns the total number of sky-pixels in this pixelization.

Definition at line 578 of file quadsphere.py.

578 
579  def __len__(self):
580  """Returns the total number of sky-pixels in this pixelization.
581  """
582  return 6 * self.resolution * self.resolution
def lsst.skypix.quadsphere.QuadSpherePixelization.__repr__ (   self)

Definition at line 589 of file quadsphere.py.

589 
590  def __repr__(self):
591  return ''.join([self.__class__.__name__, '(',
592  repr(self.resolution), ', ', repr(self.padding), ')'])
def lsst.skypix.quadsphere.QuadSpherePixelization._fiducialXPlane (   self,
  root,
  ix 
)
private

Definition at line 713 of file quadsphere.py.

713 
714  def _fiducialXPlane(self, root, ix):
715  assert isinstance(ix, (int,long)) and ix >= 0 and ix <= self.resolution
716  x = 2.0 * float(ix) / float(self.resolution) - 1.0
717  c, b = self.center[root], self.x[root]
718  v = (c[0] + x * b[0], c[1] + x * b[1], c[2] + x * b[2])
719  return geom.normalize(self.xrot[root](v, 1.0, 0.0))
def lsst.skypix.quadsphere.QuadSpherePixelization._fiducialYPlane (   self,
  root,
  iy 
)
private

Definition at line 720 of file quadsphere.py.

720 
721  def _fiducialYPlane(self, root, iy):
722  assert isinstance(iy, (int,long)) and iy >= 0 and iy <= self.resolution
723  y = 2.0 * float(iy) / float(self.resolution) - 1.0
724  c, b = self.center[root], self.y[root]
725  v = (c[0] + y * b[0], c[1] + y * b[1], c[2] + y * b[2])
726  return geom.normalize(self.yrot[root](v, 1.0, 0.0))
def lsst.skypix.quadsphere.QuadSpherePixelization._intersect (   self,
  pixels,
  poly,
  root,
  box 
)
private

Definition at line 727 of file quadsphere.py.

727 
728  def _intersect(self, pixels, poly, root, box):
729  dx = box[1] - box[0]
730  dy = box[3] - box[2]
731  if dx == 1 and dy == 1:
732  i = self.id(root, box[0], box[2])
733  if poly.intersects(self.getGeometry(i)):
734  pixels.append(i)
735  return
736  if dx >= dy:
737  # Split x range and recurse
738  xsplit = box[0] + dx / 2
739  p = poly.clip(self.xplane[root][xsplit][1])
740  if p != None:
741  self._intersect(pixels, p, root, (box[0], xsplit, box[2], box[3]))
742  p = poly.clip(self.xplane[root][xsplit][2])
743  if p != None:
744  self._intersect(pixels, p, root, (xsplit, box[1], box[2], box[3]))
745  else:
746  # Split y range and recurse
747  ysplit = box[2] + dy / 2
748  p = poly.clip(self.yplane[root][ysplit][1])
749  if p != None:
750  self._intersect(pixels, p, root, (box[0], box[1], box[2], ysplit))
751  p = poly.clip(self.yplane[root][ysplit][2])
752  if p != None:
753  self._intersect(pixels, p, root, (box[0], box[1], ysplit, box[3]))
754 
def lsst.skypix.quadsphere.QuadSpherePixelization.coords (   self,
  pixelId 
) 
Maps from a sky-pixel id to a root pixel number and x, y 
pixel coordinates.

Definition at line 661 of file quadsphere.py.

661 
662  def coords(self, pixelId):
663  """Maps from a sky-pixel id to a root pixel number and x, y
664  pixel coordinates.
665  """
666  if not isinstance(pixelId, (int, long)):
667  raise TypeError(
668  'Sky-pixel id must be an int or long')
669  R = self.resolution
670  R2 = R ** 2
671  if pixelId < 0 or pixelId >= 6 * R2:
672  raise RuntimeError(
673  'Invalid sky-pixel id %d - value must be in range [0, %d)' %
674  (pixelId, 6 * R2))
675  r = pixelId / R2
676  if r > 0 and r < 5:
677  # equatorial pixel case
678  pixelId -= R2
679  y = pixelId / (4 * R)
680  pixelId -= 4 * R * y
681  root = 1 + pixelId / R
682  x = pixelId - (root - 1) * R
683  return (root, x, y)
684  # polar root pixels
685  if r == 5:
686  root = 0
687  i = 6 * R2 - pixelId - 1
688  else:
689  root = 5
690  i = pixelId
691  s = int(math.floor(math.sqrt(i)))
692  if R & 1 == 0:
693  ring = 0.5 * (s + ((s & 1) ^ 1))
694  else:
695  ring = 0.5 * (s + (s & 1))
696  r = 2.0 * ring - 1.0
697  if r < 0.0:
698  return (root, R / 2, R / 2)
699  i -= r * r
700  if i <= 2.0 * ring:
701  dx = ring - i
702  dy = -ring
703  elif i <= 4.0 * ring:
704  dx = -ring
705  dy = i - 3.0 * ring
706  elif i <= 6.0 * ring:
707  dx = i - 5.0 * ring
708  dy = ring
709  else:
710  dx = ring
711  dy = 7.0 * ring - i
712  return (root, int(dx + 0.5 * (R - 1)), int(dy + 0.5 * (R - 1)))
def lsst.skypix.quadsphere.QuadSpherePixelization.getCenter (   self,
  pixelId 
) 
Returns the center of a sky-pixel as a unit cartesian 3-vector.

Definition at line 417 of file quadsphere.py.

417 
418  def getCenter(self, pixelId):
419  """Returns the center of a sky-pixel as a unit cartesian 3-vector.
420  """
421  root, ix, iy = self.coords(pixelId)
422  xc = 2.0 * (float(ix) + 0.5) / float(self.resolution) - 1.0
423  yc = 2.0 * (float(iy) + 0.5) / float(self.resolution) - 1.0
424  c, x, y = self.center[root], self.x[root], self.y[root]
425  return geom.normalize((c[0] + xc * x[0] + yc * y[0],
426  c[1] + xc * x[1] + yc * y[1],
427  c[2] + xc * x[2] + yc * y[2]))
def lsst.skypix.quadsphere.QuadSpherePixelization.getGeometry (   self,
  pixelId,
  fiducial = False 
) 
Returns a spherical convex polygon corresponding to the fiducial
or padded boundaries of the sky-pixel with the specified id.

Definition at line 371 of file quadsphere.py.

371 
372  def getGeometry(self, pixelId, fiducial=False):
373  """Returns a spherical convex polygon corresponding to the fiducial
374  or padded boundaries of the sky-pixel with the specified id.
375  """
376  root, ix, iy = self.coords(pixelId)
377  xl = 2.0 * float(ix) / float(self.resolution) - 1.0
378  xr = 2.0 * float(ix + 1) / float(self.resolution) - 1.0
379  yb = 2.0 * float(iy) / float(self.resolution) - 1.0
380  yt = 2.0 * float(iy + 1) / float(self.resolution) - 1.0
381  c, x, y = self.center[root], self.x[root], self.y[root]
382  v = map(geom.normalize, [(c[0] + xl * x[0] + yb * y[0],
383  c[1] + xl * x[1] + yb * y[1],
384  c[2] + xl * x[2] + yb * y[2]),
385  (c[0] + xr * x[0] + yb * y[0],
386  c[1] + xr * x[1] + yb * y[1],
387  c[2] + xr * x[2] + yb * y[2]),
388  (c[0] + xr * x[0] + yt * y[0],
389  c[1] + xr * x[1] + yt * y[1],
390  c[2] + xr * x[2] + yt * y[2]),
391  (c[0] + xl * x[0] + yt * y[0],
392  c[1] + xl * x[1] + yt * y[1],
393  c[2] + xl * x[2] + yt * y[2])])
394  if not fiducial and self.padding > 0.0:
395  # Determine angles by which edge planes must be rotated outwards
396  sp = math.sin(self.padding)
397  theta = map(lambda x: 0.5 * geom.cartesianAngularSep(x[0], x[1]),
398  ((v[0],v[3]), (v[1],v[0]), (v[2],v[1]), (v[3],v[2])))
399  sina = map(lambda x: sp / math.cos(math.radians(x)), theta)
400  cosa = map(lambda x: math.sqrt(1.0 - x * x), sina)
401  # find plane equations of fiducial pixel boundaries
402  xlp = self.xplane[root][ix][0]
403  ybp = self.yplane[root][iy][0]
404  xrp = self.xplane[root][ix + 1][0]
405  ytp = self.yplane[root][iy + 1][0]
406  # rotate edge planes outwards
407  xlp = self.xrot[root](xlp, -sina[0], cosa[0])
408  ybp = self.yrot[root](ybp, -sina[1], cosa[1])
409  xrp = self.xrot[root](xrp, sina[2], cosa[2])
410  ytp = self.yrot[root](ytp, sina[3], cosa[3])
411  # intersect rotated planes to find vertices of padded sky-pixel
412  v = map(geom.normalize, [geom.cross(xlp, ybp),
413  geom.cross(xrp, ybp),
414  geom.cross(xrp, ytp),
415  geom.cross(xlp, ytp)])
416  return geom.SphericalConvexPolygon(v)
def lsst.skypix.quadsphere.QuadSpherePixelization.getNeighbors (   self,
  pixelId 
) 
Returns a list of ids for sky-pixels adjacent to specified pixel.

Definition at line 428 of file quadsphere.py.

428 
429  def getNeighbors(self, pixelId):
430  """Returns a list of ids for sky-pixels adjacent to specified pixel.
431  """
432  root, ix, iy = self.coords(pixelId)
433  R = self.resolution
434  if ix == 0:
435  # left edge
436  if iy == 0:
437  return self.cornerNeighbors[root][0]
438  elif iy == R - 1:
439  return self.cornerNeighbors[root][1]
440  elif root == 0:
441  neighbors = (self.id(2, R - iy - 2, R - 1),
442  self.id(2, R - iy - 1, R - 1),
443  self.id(2, R - iy, R - 1))
444  elif root == 5:
445  neighbors = (self.id(2, iy - 1, 0),
446  self.id(2, iy, 0),
447  self.id(2, iy + 1, 0))
448  else:
449  r = root - 1
450  if r == 0: r = 4
451  neighbors = (self.id(r, R - 1, iy - 1),
452  self.id(r, R - 1, iy ),
453  self.id(r, R - 1, iy + 1))
454  return neighbors + (self.id(root, 1, iy - 1),
455  self.id(root, 1, iy ),
456  self.id(root, 1, iy + 1),
457  self.id(root, 0, iy - 1),
458  self.id(root, 0, iy + 1))
459  elif ix == R - 1:
460  # right edge
461  if iy == 0:
462  return self.cornerNeighbors[root][2]
463  elif iy == R - 1:
464  return self.cornerNeighbors[root][3]
465  elif root == 0:
466  neighbors = (self.id(4, iy - 1, R - 1),
467  self.id(4, iy, R - 1),
468  self.id(4, iy + 1, R - 1))
469  elif root == 5:
470  neighbors = (self.id(4, R - iy - 2, 0),
471  self.id(4, R - iy - 1, 0),
472  self.id(4, R - iy, 0))
473  else:
474  r = root + 1
475  if r == 5: r = 1
476  neighbors = (self.id(r, 0, iy - 1),
477  self.id(r, 0, iy ),
478  self.id(r, 0, iy + 1))
479  return neighbors + (self.id(root, R - 2, iy - 1),
480  self.id(root, R - 2, iy),
481  self.id(root, R - 2, iy + 1),
482  self.id(root, R - 1, iy - 1),
483  self.id(root, R - 1, iy + 1))
484  elif iy == 0:
485  # bottom edge
486  if root == 0:
487  neighbors = (self.id(3, ix - 1, R - 1),
488  self.id(3, ix , R - 1),
489  self.id(3, ix + 1, R - 1))
490  elif root == 1:
491  neighbors = (self.id(5, R - ix - 2, 0),
492  self.id(5, R - ix - 1, 0),
493  self.id(5, R - ix, 0))
494  elif root == 2:
495  neighbors = (self.id(5, 0, ix - 1),
496  self.id(5, 0, ix ),
497  self.id(5, 0, ix + 1))
498  elif root == 3:
499  neighbors = (self.id(5, ix - 1, R - 1),
500  self.id(5, ix , R - 1),
501  self.id(5, ix + 1, R - 1))
502  elif root == 4:
503  neighbors = (self.id(5, R - 1, R - ix - 2),
504  self.id(5, R - 1, R - ix - 1),
505  self.id(5, R - 1, R - ix ))
506  else:
507  neighbors = (self.id(1, R - ix - 2, 0),
508  self.id(1, R - ix - 1, 0),
509  self.id(1, R - ix, 0))
510  return neighbors + (self.id(root, ix - 1, 1),
511  self.id(root, ix, 1),
512  self.id(root, ix + 1, 1),
513  self.id(root, ix - 1, 0),
514  self.id(root, ix + 1, 0))
515  elif iy == R - 1:
516  # top edge
517  if root == 0:
518  neighbors = (self.id(1, R - ix - 2, R - 1),
519  self.id(1, R - ix - 1, R - 1),
520  self.id(1, R - ix, R - 1))
521  elif root == 1:
522  neighbors = (self.id(0, R - ix - 2, R - 1),
523  self.id(0, R - ix - 1, R - 1),
524  self.id(0, R - ix, R - 1))
525  elif root == 2:
526  neighbors = (self.id(0, 0, R - ix - 2),
527  self.id(0, 0, R - ix - 1),
528  self.id(0, 0, R - ix ))
529  elif root == 3:
530  neighbors = (self.id(0, ix - 1, 0),
531  self.id(0, ix, 0),
532  self.id(0, ix + 1, 0))
533  elif root == 4:
534  neighbors = (self.id(0, R - 1, ix - 1),
535  self.id(0, R - 1, ix ),
536  self.id(0, R - 1, ix + 1))
537  else:
538  neighbors = (self.id(3, ix - 1, 0),
539  self.id(3, ix, 0),
540  self.id(3, ix + 1, 0))
541  return neighbors + (self.id(root, ix - 1, R - 2),
542  self.id(root, ix, R - 2),
543  self.id(root, ix + 1, R - 2),
544  self.id(root, ix - 1, R - 1),
545  self.id(root, ix + 1, R - 1))
546  # interior pixel
547  return (self.id(root, ix - 1, iy - 1),
548  self.id(root, ix , iy - 1),
549  self.id(root, ix + 1, iy - 1),
550  self.id(root, ix - 1, iy),
551  self.id(root, ix + 1, iy),
552  self.id(root, ix - 1, iy + 1),
553  self.id(root, ix , iy + 1),
554  self.id(root, ix + 1, iy + 1))
def lsst.skypix.quadsphere.QuadSpherePixelization.id (   self,
  root,
  x,
  y 
) 
Maps from a root pixel number and x, y pixel coordinates
to a sky-pixel id.

Definition at line 621 of file quadsphere.py.

621 
622  def id(self, root, x, y):
623  """Maps from a root pixel number and x, y pixel coordinates
624  to a sky-pixel id.
625  """
626  if not all(isinstance(p, (int, long)) for p in (root, x, y)):
627  raise TypeError('Pixel coordinates must all be of type int or long')
628  if root < 0 or root > 5:
629  raise RuntimeError('root pixel number must be in range [0,6)')
630  R = self.resolution
631  if x < 0 or x >= R or y < 0 or y >= R:
632  raise RuntimeError('x and y pixel coordinates must ' +
633  'be in range [0,%d)' % R)
634  if root > 0 and root < 5:
635  # easy equatorial pixel case
636  return R ** 2 + 4 * R * y + (root - 1) * R + x
637  # assign ids to pixels in expanding concentric squares
638  dx = x - 0.5 * (R - 1)
639  dy = y - 0.5 * (R - 1)
640  ring = max(abs(dx), abs(dy))
641  r = 2.0 * ring - 1.0
642  if r < 0:
643  if root == 5:
644  return 0
645  else:
646  return 6 * R ** 2 - 1
647  if dy == -ring:
648  i = int(r * r + ring - dx)
649  elif dx == -ring:
650  i = int(r * r + 3.0 * ring + dy)
651  elif dy == ring:
652  i = int(r * r + 5.0 * ring + dx)
653  else: # dx == ring
654  i = int(r * r + 7.0 * ring - dy)
655  if root == 5:
656  # south pole - done
657  return i
658  # north pole: use shrinking concentric squares and wind
659  # in the opposite direction
660  return 6 * R ** 2 - i - 1
lsst::afw::geom::all
bool all(CoordinateExpr< N > const &expr)
Return true if all elements are true.
Definition: CoordinateExpr.h:80
def lsst.skypix.quadsphere.QuadSpherePixelization.intersect (   self,
  polygon 
) 
Returns a list of ids for sky-pixels intersecting the given
polygon. Intersection is computed against padded sky-pixels.

Definition at line 555 of file quadsphere.py.

555 
556  def intersect(self, polygon):
557  """Returns a list of ids for sky-pixels intersecting the given
558  polygon. Intersection is computed against padded sky-pixels.
559  """
560  if not isinstance(polygon, geom.SphericalConvexPolygon):
561  raise TypeError('Expecting a SphericalConvexPolygon as input ' +
562  'to the sky-pixel intersection computation')
563  pixels = []
564  R = self.resolution
565  for root in xrange(6):
566  # clip against padded root pixel edge planes
567  poly = polygon
568  for p in (self.xplane[root][0][2], self.xplane[root][R][1],
569  self.yplane[root][0][2], self.yplane[root][R][1]):
570  poly = poly.clip(p)
571  if poly == None:
572  break
573  if poly == None:
574  # polygon doesn't intersect this root triangle
575  continue
576  self._intersect(pixels, poly, root, (0, R, 0, R))
577  return pixels
def lsst.skypix.quadsphere.QuadSpherePixelization.pixel (   self,
  theta,
  phi 
) 
Maps from spherical coordinates (longitude and latitude
angles theta and phi, both in radians) to a sky-pixel id.

Definition at line 593 of file quadsphere.py.

593 
594  def pixel(self, theta, phi):
595  """Maps from spherical coordinates (longitude and latitude
596  angles theta and phi, both in radians) to a sky-pixel id.
597  """
598  R = self.resolution
599  root = int(math.fmod(0.5 + 2.0 * theta / math.pi, 4.0))
600  theta1 = theta - 0.5 * math.pi * root
601  root += 1
602  tanPhi = math.tan(phi)
603  y = tanPhi / math.cos(theta1)
604  if y > 1:
605  root = 0
606  x = -math.sin(theta) / tanPhi
607  y = math.cos(theta) / tanPhi
608  elif y < -1:
609  root = 5
610  x = math.sin(theta) / tanPhi
611  y = math.cos(theta) / tanPhi
612  else:
613  x = math.tan(theta1)
614  x = int(R * 0.5 * (x + 1.0))
615  y = int(R * 0.5 * (y + 1.0))
616  if x >= R:
617  x = R - 1
618  if y >= R:
619  y = R - 1
620  return self.id(root, x, y)

Member Data Documentation

lsst.skypix.quadsphere.QuadSpherePixelization.center

Definition at line 100 of file quadsphere.py.

lsst.skypix.quadsphere.QuadSpherePixelization.cornerNeighbors

Definition at line 158 of file quadsphere.py.

lsst.skypix.quadsphere.QuadSpherePixelization.padding

Definition at line 96 of file quadsphere.py.

lsst.skypix.quadsphere.QuadSpherePixelization.resolution

Definition at line 95 of file quadsphere.py.

lsst.skypix.quadsphere.QuadSpherePixelization.x

Definition at line 102 of file quadsphere.py.

lsst.skypix.quadsphere.QuadSpherePixelization.xplane

Definition at line 110 of file quadsphere.py.

lsst.skypix.quadsphere.QuadSpherePixelization.xrot

Definition at line 106 of file quadsphere.py.

lsst.skypix.quadsphere.QuadSpherePixelization.y

Definition at line 104 of file quadsphere.py.

lsst.skypix.quadsphere.QuadSpherePixelization.yplane

Definition at line 111 of file quadsphere.py.

lsst.skypix.quadsphere.QuadSpherePixelization.yrot

Definition at line 108 of file quadsphere.py.

The documentation for this class was generated from the following file:
Generated on Thu Sep 24 2015 02:29:43 for LSSTApplications by   doxygen 1.8.5