LSSTApplications  11.0-24-g0a022a1,14.0+77,15.0,15.0+1
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__ (self, resolution, paddingRad)
 
def getGeometry (self, pixelId, fiducial=False)
 
def getCenter (self, pixelId)
 
def getNeighbors (self, pixelId)
 
def intersect (self, polygon)
 
def __len__ (self)
 
def __iter__ (self)
 
def __repr__ (self)
 
def pixel (self, theta, phi)
 
def id (self, root, x, y)
 
def coords (self, pixelId)
 

Public Attributes

 resolution
 
 padding
 
 center
 
 x
 
 y
 
 xrot
 
 yrot
 
 xplane
 
 yplane
 
 cornerNeighbors
 

Private Member Functions

def _fiducialXPlane (self, root, ix)
 
def _fiducialYPlane (self, root, iy)
 
def _intersect (self, pixels, poly, root, box)
 

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 59 of file quadsphere.py.

Constructor & Destructor Documentation

◆ __init__()

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 89 of file quadsphere.py.

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

Member Function Documentation

◆ __iter__()

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 601 of file quadsphere.py.

601  def __iter__(self):
602  """Returns an iterator over the sky-pixel ids of all the pixels
603  in this pixelization.
604  """
605  return iter(range(len(self)))
606 

◆ __len__()

def lsst.skypix.quadsphere.QuadSpherePixelization.__len__ (   self)
Returns the total number of sky-pixels in this pixelization.

Definition at line 596 of file quadsphere.py.

596  def __len__(self):
597  """Returns the total number of sky-pixels in this pixelization.
598  """
599  return 6 * self.resolution * self.resolution
600 

◆ __repr__()

def lsst.skypix.quadsphere.QuadSpherePixelization.__repr__ (   self)

Definition at line 607 of file quadsphere.py.

607  def __repr__(self):
608  return ''.join([self.__class__.__name__, '(',
609  repr(self.resolution), ', ', repr(self.padding), ')'])
610 

◆ _fiducialXPlane()

def lsst.skypix.quadsphere.QuadSpherePixelization._fiducialXPlane (   self,
  root,
  ix 
)
private

Definition at line 731 of file quadsphere.py.

731  def _fiducialXPlane(self, root, ix):
732  assert isinstance(ix, numbers.Integral) and ix >= 0 and ix <= self.resolution
733  x = 2.0 * float(ix) / float(self.resolution) - 1.0
734  c, b = self.center[root], self.x[root]
735  v = (c[0] + x * b[0], c[1] + x * b[1], c[2] + x * b[2])
736  return geom.normalize(self.xrot[root](v, 1.0, 0.0))
737 

◆ _fiducialYPlane()

def lsst.skypix.quadsphere.QuadSpherePixelization._fiducialYPlane (   self,
  root,
  iy 
)
private

Definition at line 738 of file quadsphere.py.

738  def _fiducialYPlane(self, root, iy):
739  assert isinstance(iy, numbers.Integral) and iy >= 0 and iy <= self.resolution
740  y = 2.0 * float(iy) / float(self.resolution) - 1.0
741  c, b = self.center[root], self.y[root]
742  v = (c[0] + y * b[0], c[1] + y * b[1], c[2] + y * b[2])
743  return geom.normalize(self.yrot[root](v, 1.0, 0.0))
744 

◆ _intersect()

def lsst.skypix.quadsphere.QuadSpherePixelization._intersect (   self,
  pixels,
  poly,
  root,
  box 
)
private

Definition at line 745 of file quadsphere.py.

745  def _intersect(self, pixels, poly, root, box):
746  dx = box[1] - box[0]
747  dy = box[3] - box[2]
748  if dx == 1 and dy == 1:
749  i = self.id(root, box[0], box[2])
750  if poly.intersects(self.getGeometry(i)):
751  pixels.append(i)
752  return
753  if dx >= dy:
754  # Split x range and recurse
755  xsplit = box[0] + dx//2
756  p = poly.clip(self.xplane[root][xsplit][1])
757  if p is not None:
758  self._intersect(pixels, p, root, (box[0], xsplit, box[2], box[3]))
759  p = poly.clip(self.xplane[root][xsplit][2])
760  if p is not None:
761  self._intersect(pixels, p, root, (xsplit, box[1], box[2], box[3]))
762  else:
763  # Split y range and recurse
764  ysplit = box[2] + dy//2
765  p = poly.clip(self.yplane[root][ysplit][1])
766  if p is not None:
767  self._intersect(pixels, p, root, (box[0], box[1], box[2], ysplit))
768  p = poly.clip(self.yplane[root][ysplit][2])
769  if p is not None:
770  self._intersect(pixels, p, root, (box[0], box[1], ysplit, box[3]))
771 
772 

◆ coords()

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 679 of file quadsphere.py.

679  def coords(self, pixelId):
680  """Maps from a sky-pixel id to a root pixel number and x, y
681  pixel coordinates.
682  """
683  if not isinstance(pixelId, numbers.Integral):
684  raise TypeError(
685  'Sky-pixel id must be an integer')
686  R = self.resolution
687  R2 = R ** 2
688  if pixelId < 0 or pixelId >= 6 * R2:
689  raise RuntimeError(
690  'Invalid sky-pixel id %d - value must be in range [0, %d)' %
691  (pixelId, 6 * R2))
692  r = pixelId//R2
693  if r > 0 and r < 5:
694  # equatorial pixel case
695  pixelId -= R2
696  y = pixelId//(4 * R)
697  pixelId -= 4 * R * y
698  root = 1 + pixelId//R
699  x = pixelId - (root - 1) * R
700  return (root, x, y)
701  # polar root pixels
702  if r == 5:
703  root = 0
704  i = 6 * R2 - pixelId - 1
705  else:
706  root = 5
707  i = pixelId
708  s = int(math.floor(math.sqrt(i)))
709  if R & 1 == 0:
710  ring = 0.5 * (s + ((s & 1) ^ 1))
711  else:
712  ring = 0.5 * (s + (s & 1))
713  r = 2.0 * ring - 1.0
714  if r < 0.0:
715  return (root, R//2, R//2)
716  i -= r * r
717  if i <= 2.0 * ring:
718  dx = ring - i
719  dy = -ring
720  elif i <= 4.0 * ring:
721  dx = -ring
722  dy = i - 3.0 * ring
723  elif i <= 6.0 * ring:
724  dx = i - 5.0 * ring
725  dy = ring
726  else:
727  dx = ring
728  dy = 7.0 * ring - i
729  return (root, int(dx + 0.5 * (R - 1)), int(dy + 0.5 * (R - 1)))
730 

◆ getCenter()

def lsst.skypix.quadsphere.QuadSpherePixelization.getCenter (   self,
  pixelId 
)
Returns the center of a sky-pixel as a unit cartesian 3-vector.

Definition at line 433 of file quadsphere.py.

433  def getCenter(self, pixelId):
434  """Returns the center of a sky-pixel as a unit cartesian 3-vector.
435  """
436  root, ix, iy = self.coords(pixelId)
437  xc = 2.0 * (float(ix) + 0.5) / float(self.resolution) - 1.0
438  yc = 2.0 * (float(iy) + 0.5) / float(self.resolution) - 1.0
439  c, x, y = self.center[root], self.x[root], self.y[root]
440  return geom.normalize((c[0] + xc * x[0] + yc * y[0],
441  c[1] + xc * x[1] + yc * y[1],
442  c[2] + xc * x[2] + yc * y[2]))
443 

◆ getGeometry()

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 387 of file quadsphere.py.

387  def getGeometry(self, pixelId, fiducial=False):
388  """Returns a spherical convex polygon corresponding to the fiducial
389  or padded boundaries of the sky-pixel with the specified id.
390  """
391  root, ix, iy = self.coords(pixelId)
392  xl = 2.0 * float(ix) / float(self.resolution) - 1.0
393  xr = 2.0 * float(ix + 1) / float(self.resolution) - 1.0
394  yb = 2.0 * float(iy) / float(self.resolution) - 1.0
395  yt = 2.0 * float(iy + 1) / float(self.resolution) - 1.0
396  c, x, y = self.center[root], self.x[root], self.y[root]
397  v = list(map(geom.normalize, [(c[0] + xl * x[0] + yb * y[0],
398  c[1] + xl * x[1] + yb * y[1],
399  c[2] + xl * x[2] + yb * y[2]),
400  (c[0] + xr * x[0] + yb * y[0],
401  c[1] + xr * x[1] + yb * y[1],
402  c[2] + xr * x[2] + yb * y[2]),
403  (c[0] + xr * x[0] + yt * y[0],
404  c[1] + xr * x[1] + yt * y[1],
405  c[2] + xr * x[2] + yt * y[2]),
406  (c[0] + xl * x[0] + yt * y[0],
407  c[1] + xl * x[1] + yt * y[1],
408  c[2] + xl * x[2] + yt * y[2])]))
409  if not fiducial and self.padding > 0.0:
410  # Determine angles by which edge planes must be rotated outwards
411  sp = math.sin(self.padding)
412  theta = [0.5 * geom.cartesianAngularSep(x[0], x[1]) for
413  x in ((v[0], v[3]), (v[1], v[0]), (v[2], v[1]), (v[3], v[2]))]
414  sina = [sp / math.cos(math.radians(x)) for x in theta]
415  cosa = [math.sqrt(1.0 - x * x) for x in sina]
416  # find plane equations of fiducial pixel boundaries
417  xlp = self.xplane[root][ix][0]
418  ybp = self.yplane[root][iy][0]
419  xrp = self.xplane[root][ix + 1][0]
420  ytp = self.yplane[root][iy + 1][0]
421  # rotate edge planes outwards
422  xlp = self.xrot[root](xlp, -sina[0], cosa[0])
423  ybp = self.yrot[root](ybp, -sina[1], cosa[1])
424  xrp = self.xrot[root](xrp, sina[2], cosa[2])
425  ytp = self.yrot[root](ytp, sina[3], cosa[3])
426  # intersect rotated planes to find vertices of padded sky-pixel
427  v = list(map(geom.normalize, [geom.cross(xlp, ybp),
428  geom.cross(xrp, ybp),
429  geom.cross(xrp, ytp),
430  geom.cross(xlp, ytp)]))
431  return geom.SphericalConvexPolygon(v)
432 

◆ getNeighbors()

def lsst.skypix.quadsphere.QuadSpherePixelization.getNeighbors (   self,
  pixelId 
)
Returns a list of ids for sky-pixels adjacent to specified pixel.

Definition at line 444 of file quadsphere.py.

444  def getNeighbors(self, pixelId):
445  """Returns a list of ids for sky-pixels adjacent to specified pixel.
446  """
447  root, ix, iy = self.coords(pixelId)
448  R = self.resolution
449  if ix == 0:
450  # left edge
451  if iy == 0:
452  return self.cornerNeighbors[root][0]
453  elif iy == R - 1:
454  return self.cornerNeighbors[root][1]
455  elif root == 0:
456  neighbors = (self.id(2, R - iy - 2, R - 1),
457  self.id(2, R - iy - 1, R - 1),
458  self.id(2, R - iy, R - 1))
459  elif root == 5:
460  neighbors = (self.id(2, iy - 1, 0),
461  self.id(2, iy, 0),
462  self.id(2, iy + 1, 0))
463  else:
464  r = root - 1
465  if r == 0:
466  r = 4
467  neighbors = (self.id(r, R - 1, iy - 1),
468  self.id(r, R - 1, iy),
469  self.id(r, R - 1, iy + 1))
470  return neighbors + (self.id(root, 1, iy - 1),
471  self.id(root, 1, iy),
472  self.id(root, 1, iy + 1),
473  self.id(root, 0, iy - 1),
474  self.id(root, 0, iy + 1))
475  elif ix == R - 1:
476  # right edge
477  if iy == 0:
478  return self.cornerNeighbors[root][2]
479  elif iy == R - 1:
480  return self.cornerNeighbors[root][3]
481  elif root == 0:
482  neighbors = (self.id(4, iy - 1, R - 1),
483  self.id(4, iy, R - 1),
484  self.id(4, iy + 1, R - 1))
485  elif root == 5:
486  neighbors = (self.id(4, R - iy - 2, 0),
487  self.id(4, R - iy - 1, 0),
488  self.id(4, R - iy, 0))
489  else:
490  r = root + 1
491  if r == 5:
492  r = 1
493  neighbors = (self.id(r, 0, iy - 1),
494  self.id(r, 0, iy),
495  self.id(r, 0, iy + 1))
496  return neighbors + (self.id(root, R - 2, iy - 1),
497  self.id(root, R - 2, iy),
498  self.id(root, R - 2, iy + 1),
499  self.id(root, R - 1, iy - 1),
500  self.id(root, R - 1, iy + 1))
501  elif iy == 0:
502  # bottom edge
503  if root == 0:
504  neighbors = (self.id(3, ix - 1, R - 1),
505  self.id(3, ix, R - 1),
506  self.id(3, ix + 1, R - 1))
507  elif root == 1:
508  neighbors = (self.id(5, R - ix - 2, 0),
509  self.id(5, R - ix - 1, 0),
510  self.id(5, R - ix, 0))
511  elif root == 2:
512  neighbors = (self.id(5, 0, ix - 1),
513  self.id(5, 0, ix),
514  self.id(5, 0, ix + 1))
515  elif root == 3:
516  neighbors = (self.id(5, ix - 1, R - 1),
517  self.id(5, ix, R - 1),
518  self.id(5, ix + 1, R - 1))
519  elif root == 4:
520  neighbors = (self.id(5, R - 1, R - ix - 2),
521  self.id(5, R - 1, R - ix - 1),
522  self.id(5, R - 1, R - ix))
523  else:
524  neighbors = (self.id(1, R - ix - 2, 0),
525  self.id(1, R - ix - 1, 0),
526  self.id(1, R - ix, 0))
527  return neighbors + (self.id(root, ix - 1, 1),
528  self.id(root, ix, 1),
529  self.id(root, ix + 1, 1),
530  self.id(root, ix - 1, 0),
531  self.id(root, ix + 1, 0))
532  elif iy == R - 1:
533  # top edge
534  if root == 0:
535  neighbors = (self.id(1, R - ix - 2, R - 1),
536  self.id(1, R - ix - 1, R - 1),
537  self.id(1, R - ix, R - 1))
538  elif root == 1:
539  neighbors = (self.id(0, R - ix - 2, R - 1),
540  self.id(0, R - ix - 1, R - 1),
541  self.id(0, R - ix, R - 1))
542  elif root == 2:
543  neighbors = (self.id(0, 0, R - ix - 2),
544  self.id(0, 0, R - ix - 1),
545  self.id(0, 0, R - ix))
546  elif root == 3:
547  neighbors = (self.id(0, ix - 1, 0),
548  self.id(0, ix, 0),
549  self.id(0, ix + 1, 0))
550  elif root == 4:
551  neighbors = (self.id(0, R - 1, ix - 1),
552  self.id(0, R - 1, ix),
553  self.id(0, R - 1, ix + 1))
554  else:
555  neighbors = (self.id(3, ix - 1, 0),
556  self.id(3, ix, 0),
557  self.id(3, ix + 1, 0))
558  return neighbors + (self.id(root, ix - 1, R - 2),
559  self.id(root, ix, R - 2),
560  self.id(root, ix + 1, R - 2),
561  self.id(root, ix - 1, R - 1),
562  self.id(root, ix + 1, R - 1))
563  # interior pixel
564  return (self.id(root, ix - 1, iy - 1),
565  self.id(root, ix, iy - 1),
566  self.id(root, ix + 1, iy - 1),
567  self.id(root, ix - 1, iy),
568  self.id(root, ix + 1, iy),
569  self.id(root, ix - 1, iy + 1),
570  self.id(root, ix, iy + 1),
571  self.id(root, ix + 1, iy + 1))
572 

◆ id()

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 639 of file quadsphere.py.

639  def id(self, root, x, y):
640  """Maps from a root pixel number and x, y pixel coordinates
641  to a sky-pixel id.
642  """
643  if not all(isinstance(p, numbers.Integral) for p in (root, x, y)):
644  raise TypeError('Pixel coordinates must all be an integer')
645  if root < 0 or root > 5:
646  raise RuntimeError('root pixel number must be in range [0,6)')
647  R = self.resolution
648  if x < 0 or x >= R or y < 0 or y >= R:
649  raise RuntimeError('x and y pixel coordinates must ' +
650  'be in range [0,%d)' % R)
651  if root > 0 and root < 5:
652  # easy equatorial pixel case
653  return R ** 2 + 4 * R * y + (root - 1) * R + x
654  # assign ids to pixels in expanding concentric squares
655  dx = x - 0.5 * (R - 1)
656  dy = y - 0.5 * (R - 1)
657  ring = max(abs(dx), abs(dy))
658  r = 2.0 * ring - 1.0
659  if r < 0:
660  if root == 5:
661  return 0
662  else:
663  return 6 * R ** 2 - 1
664  if dy == -ring:
665  i = int(r * r + ring - dx)
666  elif dx == -ring:
667  i = int(r * r + 3.0 * ring + dy)
668  elif dy == ring:
669  i = int(r * r + 5.0 * ring + dx)
670  else: # dx == ring
671  i = int(r * r + 7.0 * ring - dy)
672  if root == 5:
673  # south pole - done
674  return i
675  # north pole: use shrinking concentric squares and wind
676  # in the opposite direction
677  return 6 * R ** 2 - i - 1
678 
bool all(CoordinateExpr< N > const &expr)
Return true if all elements are true.
int max
Definition: BoundedField.cc:99
int id
Definition: CR.cc:155

◆ intersect()

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 573 of file quadsphere.py.

573  def intersect(self, polygon):
574  """Returns a list of ids for sky-pixels intersecting the given
575  polygon. Intersection is computed against padded sky-pixels.
576  """
577  if not isinstance(polygon, geom.SphericalConvexPolygon):
578  raise TypeError('Expecting a SphericalConvexPolygon as input ' +
579  'to the sky-pixel intersection computation')
580  pixels = []
581  R = self.resolution
582  for root in range(6):
583  # clip against padded root pixel edge planes
584  poly = polygon
585  for p in (self.xplane[root][0][2], self.xplane[root][R][1],
586  self.yplane[root][0][2], self.yplane[root][R][1]):
587  poly = poly.clip(p)
588  if poly is None:
589  break
590  if poly is None:
591  # polygon doesn't intersect this root triangle
592  continue
593  self._intersect(pixels, poly, root, (0, R, 0, R))
594  return pixels
595 

◆ pixel()

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 611 of file quadsphere.py.

611  def pixel(self, theta, phi):
612  """Maps from spherical coordinates (longitude and latitude
613  angles theta and phi, both in radians) to a sky-pixel id.
614  """
615  R = self.resolution
616  root = int(math.fmod(0.5 + 2.0 * theta / math.pi, 4.0))
617  theta1 = theta - 0.5 * math.pi * root
618  root += 1
619  tanPhi = math.tan(phi)
620  y = tanPhi / math.cos(theta1)
621  if y > 1:
622  root = 0
623  x = -math.sin(theta) / tanPhi
624  y = math.cos(theta) / tanPhi
625  elif y < -1:
626  root = 5
627  x = math.sin(theta) / tanPhi
628  y = math.cos(theta) / tanPhi
629  else:
630  x = math.tan(theta1)
631  x = int(R * 0.5 * (x + 1.0))
632  y = int(R * 0.5 * (y + 1.0))
633  if x >= R:
634  x = R - 1
635  if y >= R:
636  y = R - 1
637  return self.id(root, x, y)
638 
table::PointKey< int > pixel

Member Data Documentation

◆ center

lsst.skypix.quadsphere.QuadSpherePixelization.center

Definition at line 116 of file quadsphere.py.

◆ cornerNeighbors

lsst.skypix.quadsphere.QuadSpherePixelization.cornerNeighbors

Definition at line 174 of file quadsphere.py.

◆ padding

lsst.skypix.quadsphere.QuadSpherePixelization.padding

Definition at line 112 of file quadsphere.py.

◆ resolution

lsst.skypix.quadsphere.QuadSpherePixelization.resolution

Definition at line 111 of file quadsphere.py.

◆ x

lsst.skypix.quadsphere.QuadSpherePixelization.x

Definition at line 118 of file quadsphere.py.

◆ xplane

lsst.skypix.quadsphere.QuadSpherePixelization.xplane

Definition at line 126 of file quadsphere.py.

◆ xrot

lsst.skypix.quadsphere.QuadSpherePixelization.xrot

Definition at line 122 of file quadsphere.py.

◆ y

lsst.skypix.quadsphere.QuadSpherePixelization.y

Definition at line 120 of file quadsphere.py.

◆ yplane

lsst.skypix.quadsphere.QuadSpherePixelization.yplane

Definition at line 127 of file quadsphere.py.

◆ yrot

lsst.skypix.quadsphere.QuadSpherePixelization.yrot

Definition at line 124 of file quadsphere.py.


The documentation for this class was generated from the following file: