pull request 6

Author: Jeff Whitaker

Changelog

@@ -1,4 +1,6 @@
 version 1.0.3 (not yet released)
+           * bug fixes for celestial projections, and non-standard sphere
+             radii (https://github.com/matplotlib/basemap/pull/6).
            * fix bug in drawparallels that results in 'KeyError' when drawing
              parallels very close together (0.1 degrees).
            * fix constant in Robinson projection (update PJ_robin.c

lib/mpl_toolkits/basemap/__init__.py

@@ -887,20 +887,26 @@ def __call__(self,x,y,inverse=False):
         Input arguments lon, lat can be either scalar floats,
         sequences, or numpy arrays.
         """
+        if self.celestial:
+            # don't assume center of map is at greenwich
+            # (only relevant for cyl or pseudo-cyl projections)
+            if (self.projection in _pseudocyl) or (self.projection in _cylproj):
+                lon_0=self.projparams['lon_0']
+            else:
+                lon_0 = 0.
         if self.celestial and not inverse:
             try:
-                x = -x
+                x = 2.*lon_0-x
             except TypeError:
-                x = [-xx for xx in x]
+                x = [2*lon_0-xx for xx in x]
         xout,yout = self.projtran(x,y,inverse=inverse)
         if self.celestial and inverse:
             try:
-                xout = -xout
+                xout = -2.*lon_0-xout
             except:
-                xout = [-xx for xx in xout]
+                xout = [-2.*lon_0-xx for xx in xout]
         return xout,yout
 
-
     def makegrid(self,nx,ny,returnxy=False):
         """
         return arrays of shape (ny,nx) containing lon,lat coordinates of
@@ -1892,7 +1898,9 @@ def drawparallels(self,circles,color='k',linewidth=1.,zorder=None, \
                 xd = (x[1:]-x[0:-1])**2
                 yd = (y[1:]-y[0:-1])**2
                 dist = np.sqrt(xd+yd)
-                split = dist > 500000.
+                #split = dist > 500000.
+                # normalize by radius of sphere
+                split = dist > 5000000./6370997.0*self.rmajor
                 if np.sum(split) and self.projection not in _cylproj:
                     ind = (np.compress(split,np.squeeze(split*np.indices(xd.shape)))+1).tolist()
                     xl = []
@@ -2191,7 +2199,9 @@ def addlon(meridians,madd):
                 xd = (x[1:]-x[0:-1])**2
                 yd = (y[1:]-y[0:-1])**2
                 dist = np.sqrt(xd+yd)
-                split = dist > 500000.
+                #split = dist > 500000.
+                # normalize by radius of sphere
+                split = dist > 5000000./6370997.0*self.rmajor
                 if np.sum(split) and self.projection not in _cylproj:
                     ind = (np.compress(split,np.squeeze(split*np.indices(xd.shape)))+1).tolist()
                     xl = []