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.geom.geometry.SphericalConvexPolygon Class Reference
Inheritance diagram for lsst.geom.geometry.SphericalConvexPolygon:
lsst.geom.geometry.SphericalRegion

Public Member Functions

def __init__
 
def getVertices
 
def getEdges
 
def getBoundingCircle
 
def getBoundingBox
 
def getZRange
 
def clip
 
def intersect
 
def area
 
def containsPoint
 
def contains
 
def intersects
 
def __repr__
 
def __eq__
 

Public Attributes

 boundingBox
 
 boundingCircle
 
 vertices
 
 edges
 

Private Member Functions

def _computeEdges
 

Detailed Description

A convex polygon on the unit sphere with great circle edges. Points
falling exactly on polygon edges are considered to be inside (contained
by) the polygon.

Definition at line 1075 of file geometry.py.

Constructor & Destructor Documentation

def lsst.geom.geometry.SphericalConvexPolygon.__init__ (   self,
  args 
)
Creates a new polygon. Arguments must be either:

1. a SphericalConvexPolygon
2. a list of vertices
3. a list of vertices and a list of corresponding edges

In the first case, a copy is constructed. In the second case,
the argument must be a sequence of 3 element tuples/lists
(unit cartesian 3-vectors) - a copy of the vertices is stored
and polygon edges are computed. In the last case, copies of the
input vertex and edge lists are stored.

Vertices must be hemispherical and in counter-clockwise order when
viewed from outside the unit sphere in a right handed coordinate
system. They must also form a convex polygon.

When edges are specified, the i-th edge must correspond to the plane
equation of great circle connecting vertices (i - 1, i), that is,
the edge should be a unit cartesian 3-vector parallel to v[i-1] ^ v[i]
(where ^ denotes the vector cross product).

Note that these conditions are NOT verified for performance reasons.
Operations on SphericalConvexPolygon objects constructed with inputs
not satisfying these conditions are undefined. Use the convex()
function to check for convexity and ordering of the vertices. Or,
use the convexHull() function to create a SphericalConvexPolygon
from an arbitrary list of hemispherical input vertices.

Definition at line 1080 of file geometry.py.

1081  def __init__(self, *args):
1082  """Creates a new polygon. Arguments must be either:
1083 
1084  1. a SphericalConvexPolygon
1085  2. a list of vertices
1086  3. a list of vertices and a list of corresponding edges
1087 
1088  In the first case, a copy is constructed. In the second case,
1089  the argument must be a sequence of 3 element tuples/lists
1090  (unit cartesian 3-vectors) - a copy of the vertices is stored
1091  and polygon edges are computed. In the last case, copies of the
1092  input vertex and edge lists are stored.
1093 
1094  Vertices must be hemispherical and in counter-clockwise order when
1095  viewed from outside the unit sphere in a right handed coordinate
1096  system. They must also form a convex polygon.
1097 
1098  When edges are specified, the i-th edge must correspond to the plane
1099  equation of great circle connecting vertices (i - 1, i), that is,
1100  the edge should be a unit cartesian 3-vector parallel to v[i-1] ^ v[i]
1101  (where ^ denotes the vector cross product).
1102 
1103  Note that these conditions are NOT verified for performance reasons.
1104  Operations on SphericalConvexPolygon objects constructed with inputs
1105  not satisfying these conditions are undefined. Use the convex()
1106  function to check for convexity and ordering of the vertices. Or,
1107  use the convexHull() function to create a SphericalConvexPolygon
1108  from an arbitrary list of hemispherical input vertices.
1109  """
1110  self.boundingBox = None
1111  self.boundingCircle = None
1112  if len(args) == 0:
1113  raise RuntimeError('Expecting at least one argument')
1114  elif len(args) == 1:
1115  if isinstance(args[0], SphericalConvexPolygon):
1116  self.vertices = list(args[0].vertices)
1117  self.edges = list(args[0].edges)
1118  else:
1119  self.vertices = list(args[0])
1120  self._computeEdges()
1121  elif len(args) == 2:
1122  self.vertices = list(args[0])
1123  self.edges = list(args[1])
1124  if len(self.edges) != len(self.vertices):
1125  raise RuntimeError(
1126  'number of edges does not match number of vertices')
1127  else:
1128  self.vertices = list(args)
1129  self._computeEdges()
1130  if len(self.vertices) < 3:
1131  raise RuntimeError(
1132  'spherical polygon must contain at least 3 vertices')

Member Function Documentation

def lsst.geom.geometry.SphericalConvexPolygon.__eq__ (   self,
  other 
)

Definition at line 1417 of file geometry.py.

1418  def __eq__(self, other):
1419  if not isinstance(other, SphericalConvexPolygon):
1420  return False
1421  n = len(self.vertices)
1422  if n != len(other.vertices):
1423  return False
1424  # Test for equality up to cyclic permutation of vertices/edges
1425  try:
1426  offset = other.vertices.index(self.vertices[0])
1427  except ValueError:
1428  return False
1429  for i in xrange(0, n):
1430  j = offset + i
1431  if j >= n:
1432  j -= n
1433  if (self.vertices[i] != self.vertices[j] or
1434  self.edges[i] != self.edges[j]):
1435  return False
1436  return True
1437 
1438 
1439 # -- Finding the median element of an array in linear time ----
1440 #
1441 # This is a necessary primitive for Megiddo's linear time
1442 # 2d linear programming algorithm.
def lsst.geom.geometry.SphericalConvexPolygon.__repr__ (   self)
Returns a string representation of this polygon.

Definition at line 1411 of file geometry.py.

1412  def __repr__(self):
1413  """Returns a string representation of this polygon.
1414  """
1415  return ''.join([self.__class__.__name__ , '(',
1416  repr(map(sphericalCoords, self.vertices)), ')'])
def lsst.geom.geometry.SphericalConvexPolygon._computeEdges (   self)
private
Computes edge plane normals from vertices.

Definition at line 1133 of file geometry.py.

1134  def _computeEdges(self):
1135  """Computes edge plane normals from vertices.
1136  """
1137  v = self.vertices
1138  n = len(v)
1139  edges = []
1140  for i in xrange(n):
1141  edges.append(normalize(cross(v[i - 1], v[i])))
1142  self.edges = edges
def lsst.geom.geometry.SphericalConvexPolygon.area (   self)
Returns the area of this spherical convex polygon.

Definition at line 1275 of file geometry.py.

1276  def area(self):
1277  """Returns the area of this spherical convex polygon.
1278  """
1279  numVerts = len(self.vertices)
1280  angleSum = 0.0
1281  for i in xrange(numVerts):
1282  tmp = cross(self.edges[i - 1], self.edges[i])
1283  sina = math.sqrt(dot(tmp, tmp))
1284  cosa = -dot(self.edges[i - 1], self.edges[i])
1285  a = math.atan2(sina, cosa)
1286  angleSum += a
1287  return angleSum - (numVerts - 2) * math.pi
def lsst.geom.geometry.SphericalConvexPolygon.clip (   self,
  plane 
)
Returns the polygon corresponding to the intersection of this
polygon with the positive half space defined by the given plane.
The plane must be specified with a cartesian unit vector (its
normal) and always passes through the origin.

Clipping is performed using the Sutherland-Hodgman algorithm,
adapted for spherical polygons.

Definition at line 1184 of file geometry.py.

1185  def clip(self, plane):
1186  """Returns the polygon corresponding to the intersection of this
1187  polygon with the positive half space defined by the given plane.
1188  The plane must be specified with a cartesian unit vector (its
1189  normal) and always passes through the origin.
1190 
1191  Clipping is performed using the Sutherland-Hodgman algorithm,
1192  adapted for spherical polygons.
1193  """
1194  vertices, edges, classification = [], [], []
1195  vin, vout = False, False
1196  for i in xrange(len(self.vertices)):
1197  d = dot(plane, self.vertices[i])
1198  if d > SIN_MIN: vin = True
1199  elif d < -SIN_MIN: vout = True
1200  else: d = 0.0
1201  classification.append(d)
1202  if not vin and not vout:
1203  # polygon and clipping plane are coplanar
1204  return None
1205  if not vout:
1206  return self
1207  elif not vin:
1208  return None
1209  v0 = self.vertices[-1]
1210  d0 = classification[-1]
1211  for i in xrange(len(self.vertices)):
1212  v1 = self.vertices[i]
1213  d1 = classification[i]
1214  if d0 > 0.0:
1215  if d1 >= 0.0:
1216  # positive to positive, positive to zero: no intersection,
1217  # add current input vertex to output polygon
1218  vertices.append(v1)
1219  edges.append(self.edges[i])
1220  else:
1221  # positive to negative: add intersection point to polygon
1222  f = d0 / (d0 - d1)
1223  v = normalize((v0[0] + (v1[0] - v0[0]) * f,
1224  v0[1] + (v1[1] - v0[1]) * f,
1225  v0[2] + (v1[2] - v0[2]) * f))
1226  vertices.append(v)
1227  edges.append(self.edges[i])
1228  elif d0 == 0.0:
1229  if d1 >= 0.0:
1230  # zero to positive: no intersection, add current input
1231  # vertex to output polygon.
1232  vertices.append(v1)
1233  edges.append(self.edges[i])
1234  # zero to zero: coplanar edge - since the polygon has vertices
1235  # on both sides of plane, this implies concavity or a
1236  # duplicate vertex. Under the convexity assumption, this
1237  # must be caused by a near duplicate vertex, so skip the
1238  # vertex.
1239  #
1240  # zero to negative: no intersection, skip the vertex
1241  else:
1242  if d1 > 0.0:
1243  # negative to positive: add intersection point to output
1244  # polygon followed by the current input vertex
1245  f = d0 / (d0 - d1)
1246  v = normalize((v0[0] + (v1[0] - v0[0]) * f,
1247  v0[1] + (v1[1] - v0[1]) * f,
1248  v0[2] + (v1[2] - v0[2]) * f))
1249  vertices.append(v)
1250  edges.append(tuple(plane))
1251  vertices.append(v1)
1252  edges.append(self.edges[i])
1253  elif d1 == 0.0:
1254  # negative to zero: add current input vertex to output
1255  # polygon
1256  vertices.append(v1)
1257  edges.append(tuple(plane))
1258  # negative to negative: no intersection, skip vertex
1259  v0 = v1
1260  d0 = d1
1261  return SphericalConvexPolygon(vertices, edges)
def lsst.geom.geometry.SphericalConvexPolygon.contains (   self,
  pointOrRegion 
)
Returns True if the specified point or spherical region is
completely contained in this polygon. Note that the implementation
is conservative where ellipses are concerned: False may be returned
for an ellipse that is actually completely contained by this polygon.

Definition at line 1294 of file geometry.py.

1295  def contains(self, pointOrRegion):
1296  """Returns True if the specified point or spherical region is
1297  completely contained in this polygon. Note that the implementation
1298  is conservative where ellipses are concerned: False may be returned
1299  for an ellipse that is actually completely contained by this polygon.
1300  """
1301  pr = pointOrRegion
1302  if isinstance(pr, SphericalConvexPolygon):
1303  for v in pr.getVertices():
1304  if not self.containsPoint(v):
1305  return False
1306  return True
1307  elif isinstance(pr, SphericalEllipse):
1308  return self.contains(pr.getBoundingCircle())
1309  elif isinstance(pr, SphericalCircle):
1310  cv = cartesianUnitVector(pr.getCenter())
1311  if not self.containsPoint(cv):
1312  return False
1313  else:
1314  minSep = INF
1315  for i in xrange(len(self.vertices)):
1316  s = minEdgeSep(cv, self.edges[i],
1317  self.vertices[i - 1], self.vertices[i])
1318  minSep = min(minSep, s)
1319  return minSep >= pr.getRadius()
1320  elif isinstance(pr, SphericalBox):
1321  # check that all box vertices are inside the polygon
1322  bMin, bMax = pr.getMin(), pr.getMax()
1323  bzMin = math.sin(math.radians(bMin[1]))
1324  bzMax = math.sin(math.radians(bMax[1]))
1325  verts = map(cartesianUnitVector,
1326  (bMin, bMax, (bMin[0], bMax[1]), (bMax[0], bMin[1])))
1327  for v in verts:
1328  if not self.containsPoint(v):
1329  return False
1330  # check that intersections of box small circles with polygon
1331  # edges either don't exist or fall outside the box.
1332  for i in xrange(len(self.vertices)):
1333  v0 = self.vertices[i - 1]
1334  v1 = self.vertices[i]
1335  e = self.edges[i]
1336  d = e[0] * e[0] + e[1] * e[1]
1337  if abs(e[2]) >= COS_MAX or d < MIN_FLOAT:
1338  # polygon edge is approximately described by z = +/-1.
1339  # box vertices are inside the polygon, so they
1340  # cannot intersect the edge.
1341  continue
1342  D = d - bzMin * bzMin
1343  if D >= 0.0:
1344  # polygon edge intersects z = bzMin
1345  D = math.sqrt(D)
1346  xr = -e[0] * e[2] * bzMin
1347  yr = -e[1] * e[2] * bzMin
1348  i1 = ((xr + e[1] * D) / d, (yr - e[0] * D) / d, bzMin)
1349  i2 = ((xr - e[1] * D) / d, (yr + e[0] * D) / d, bzMin)
1350  if (pr.containsPoint(sphericalCoords(i1)) or
1351  pr.containsPoint(sphericalCoords(i2))):
1352  return False
1353  D = d - bzMax * bzMax
1354  if D >= 0.0:
1355  # polygon edge intersects z = bzMax
1356  D = math.sqrt(D)
1357  xr = -e[0] * e[2] * bzMax
1358  yr = -e[1] * e[2] * bzMax
1359  i1 = ((xr + e[1] * D) / d, (yr - e[0] * D) / d, bzMax)
1360  i2 = ((xr - e[1] * D) / d, (yr + e[0] * D) / d, bzMax)
1361  if (pr.containsPoint(sphericalCoords(i1)) or
1362  pr.containsPoint(sphericalCoords(i2))):
1363  return False
1364  return True
1365  else:
1366  return self.containsPoint(cartesianUnitVector(pr))
def lsst.geom.geometry.SphericalConvexPolygon.containsPoint (   self,
  v 
)

Definition at line 1288 of file geometry.py.

1289  def containsPoint(self, v):
1290  for e in self.edges:
1291  if dot(v, e) < 0.0:
1292  return False
1293  return True
def lsst.geom.geometry.SphericalConvexPolygon.getBoundingBox (   self)
Returns a bounding box for this spherical convex polygon.

Definition at line 1166 of file geometry.py.

1167  def getBoundingBox(self):
1168  """Returns a bounding box for this spherical convex polygon.
1169  """
1170  if self.boundingBox == None:
1171  self.boundingBox = SphericalBox()
1172  for i in xrange(0, len(self.vertices)):
1173  self.boundingBox.extend(SphericalBox.edge(
1174  self.vertices[i - 1], self.vertices[i], self.edges[i]))
1175  return self.boundingBox
def lsst.geom.geometry.SphericalConvexPolygon.getBoundingCircle (   self)
Returns a bounding circle (not necessarily minimal) for this
spherical convex polygon.

Definition at line 1154 of file geometry.py.

1155  def getBoundingCircle(self):
1156  """Returns a bounding circle (not necessarily minimal) for this
1157  spherical convex polygon.
1158  """
1159  if self.boundingCircle == None:
1160  center = centroid(self.vertices)
1161  radius = 0.0
1162  for v in self.vertices:
1163  radius = max(radius, cartesianAngularSep(center, v))
1164  self.boundingCircle = SphericalCircle(center, radius)
1165  return self.boundingCircle
def lsst.geom.geometry.SphericalConvexPolygon.getEdges (   self)
Returns the list of polygon edges. The i-th edge is the plane
equation for the great circle edge formed by vertices i-1 and i.

Definition at line 1148 of file geometry.py.

1149  def getEdges(self):
1150  """Returns the list of polygon edges. The i-th edge is the plane
1151  equation for the great circle edge formed by vertices i-1 and i.
1152  """
1153  return self.edges
def lsst.geom.geometry.SphericalConvexPolygon.getVertices (   self)
Returns the list of polygon vertices.

Definition at line 1143 of file geometry.py.

1144  def getVertices(self):
1145  """Returns the list of polygon vertices.
1146  """
1147  return self.vertices
def lsst.geom.geometry.SphericalConvexPolygon.getZRange (   self)
Returns the z coordinate range spanned by this spherical
convex polygon.

Definition at line 1176 of file geometry.py.

1177  def getZRange(self):
1178  """Returns the z coordinate range spanned by this spherical
1179  convex polygon.
1180  """
1181  bbox = self.getBoundingBox()
1182  return (math.sin(math.radians(bbox.getMin()[1])),
1183  math.sin(math.radians(bbox.getMax()[1])))
def lsst.geom.geometry.SphericalConvexPolygon.intersect (   self,
  poly 
)
Returns the intersection of poly with this polygon, or
None if the intersection does not exist.

Definition at line 1262 of file geometry.py.

1263  def intersect(self, poly):
1264  """Returns the intersection of poly with this polygon, or
1265  None if the intersection does not exist.
1266  """
1267  if not isinstance(poly, SphericalConvexPolygon):
1268  raise TypeError('Expecting a SphericalConvexPolygon object')
1269  p = self
1270  for edge in poly.getEdges():
1271  p = p.clip(edge)
1272  if p == None:
1273  break
1274  return p
def lsst.geom.geometry.SphericalConvexPolygon.intersects (   self,
  pointOrRegion 
)
Returns True if the given point or spherical region intersects
this polygon. Note that the implementation is conservative where
ellipses are concerned: True may be returned for an ellipse that
is actually disjoint from this polygon.

Definition at line 1367 of file geometry.py.

1368  def intersects(self, pointOrRegion):
1369  """Returns True if the given point or spherical region intersects
1370  this polygon. Note that the implementation is conservative where
1371  ellipses are concerned: True may be returned for an ellipse that
1372  is actually disjoint from this polygon.
1373  """
1374  pr = pointOrRegion
1375  if isinstance(pr, SphericalConvexPolygon):
1376  return self.intersect(pr) != None
1377  elif isinstance(pr, SphericalEllipse):
1378  return self.intersects(pr.getBoundingCircle())
1379  elif isinstance(pr, SphericalCircle):
1380  cv = cartesianUnitVector(pr.getCenter())
1381  if self.containsPoint(cv):
1382  return True
1383  else:
1384  minSep = INF
1385  for i in xrange(len(self.vertices)):
1386  s = minEdgeSep(cv, self.edges[i],
1387  self.vertices[i - 1], self.vertices[i])
1388  minSep = min(minSep, s)
1389  return minSep <= pr.getRadius()
1390  elif isinstance(pr, SphericalBox):
1391  minTheta = math.radians(pr.getMin()[0])
1392  maxTheta = math.radians(pr.getMax()[0])
1393  bzMin = math.sin(math.radians(pr.getMin()[1]))
1394  bzMax = math.sin(math.radians(pr.getMax()[1]))
1395  p = self.clip((-math.sin(minTheta), math.cos(minTheta), 0.0))
1396  if pr.getThetaExtent() > 180.0:
1397  if p != None:
1398  zMin, zMax = p.getZRange()
1399  if zMin <= bzMax and zMax >= bzMin:
1400  return True
1401  p = self.clip((math.sin(maxTheta), -math.cos(maxTheta), 0.0))
1402  else:
1403  if p != None:
1404  p = p.clip((math.sin(maxTheta), -math.cos(maxTheta), 0.0))
1405  if p == None:
1406  return False
1407  zMin, zMax = p.getZRange()
1408  return zMin <= bzMax and zMax >= bzMin
1409  else:
1410  return self.containsPoint(cartesianUnitVector(pr))

Member Data Documentation

lsst.geom.geometry.SphericalConvexPolygon.boundingBox

Definition at line 1109 of file geometry.py.

lsst.geom.geometry.SphericalConvexPolygon.boundingCircle

Definition at line 1110 of file geometry.py.

lsst.geom.geometry.SphericalConvexPolygon.edges

Definition at line 1116 of file geometry.py.

lsst.geom.geometry.SphericalConvexPolygon.vertices

Definition at line 1115 of file geometry.py.


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