title: Python Work out area of a polygon on a spherical surface
body: <p>I have a series of points, of right ascension and declination values.
  These points correspond to the vertices of a polygon on the surface of a sphere.</p>
  
  <p>What would be the best way to calculate the area enclosed by these points? I would assume that converting the points with an equal-area projection, and then carrying out typical polygonal area calculating on a flat surface would be an appropriate solution.</p>
  
  <p>note: I cannot use custom python libraries. eg pyproj or shapely</p>
  
  <p>Example code (works for latitude longitude, what modifications would be required to enure this works with sky coordinates?)</p>
  
  <pre><code>def reproject(latitude, longitude):
      """Returns the x &amp; y coordinates in metres using a sinusoidal projection"""
      from math import pi, cos, radians
      earth_radius = 6371009
      lat_dist = pi * earth_radius / 180.0
  
      y = [lat * lat_dist for lat in latitude]
      x = [long * lat_dist * cos(radians(lat)) 
                  for lat, long in zip(latitude, longitude)]
      return x, y
  
  def area_of_polygon(x, y):
      """Calculates the area of an arbitrary polygon given its vertices"""
      area = 0.0
      for i in xrange(-1, len(x)-1):
          area += x[i] * (y[i+1] - y[i-1])
      return abs(area) / 2.0
  
  dec = [-15.,89.,89.,-15.,-15.]    
  ra = [105.,105.,285.,285.,105.]
  
  x,y = reproject(dec, ra)
  print area_of_polygon(x,y)
  </code></pre>