From: <js...@us...> - 2007-11-07 16:42:54
|
Revision: 4145 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4145&view=rev Author: jswhit Date: 2007-11-07 08:42:46 -0800 (Wed, 07 Nov 2007) Log Message: ----------- remove debug code Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-07 16:25:58 UTC (rev 4144) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-07 16:42:46 UTC (rev 4145) @@ -18,7 +18,6 @@ from shapely.geometry import LineString as LineShape from shapely.geometry import Point as PointShape from shapely import wkb -import time # basemap data files now installed in lib/matplotlib/toolkits/basemap/data basemap_datadir = os.sep.join([os.path.dirname(__file__), 'data']) @@ -744,7 +743,6 @@ bdatmetafile = open(os.path.join(basemap_datadir,name+'meta_'+self.resolution+'.dat'),'r') except: raise IOError, msg - timetot = 0. polygons = [] polygon_types = [] # see if map projection region polygon contains a pole. This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-08 22:06:41
|
Revision: 4164 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4164&view=rev Author: jswhit Date: 2007-11-08 14:06:40 -0800 (Thu, 08 Nov 2007) Log Message: ----------- fix typo (byteswapped() --> byteswap()) Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-08 21:53:39 UTC (rev 4163) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-08 22:06:40 UTC (rev 4164) @@ -792,7 +792,7 @@ # numpy array (first column is lons, second is lats). polystring = bdatfile.read(bytecount) if not npy.little_endian: - b = npy.reshape(npy.fromstring(polystring,dtype=npy.float64).byteswapped(),(npts,2)) + b = npy.reshape(npy.fromstring(polystring,dtype=npy.float64).byteswap(),(npts,2)) else: b = npy.reshape(npy.fromstring(polystring,dtype=npy.float64),(npts,2)) # if map boundary polygon is a valid one in lat/lon @@ -853,8 +853,8 @@ # pylab.show() if poly.is_valid: poly = poly.intersection(boundarypolyll) - else: - print 'warning, invalid ',name,' geometry',poly.area + #else: + # print 'warning, invalid ',name,' geometry',poly.area # create iterable object with geometries # that intersect map region. if hasattr(poly,'geoms'): @@ -919,8 +919,8 @@ #try: if poly.is_valid: poly = boundarypolyxy.intersection(poly) - else: - print 'warning, invalid ',name,' geometry',poly.area + #else: + # print 'warning, invalid ',name,' geometry',poly.area #except: # continue # create iterable object with geometries This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-08 22:12:02
|
Revision: 4165 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4165&view=rev Author: jswhit Date: 2007-11-08 14:11:50 -0800 (Thu, 08 Nov 2007) Log Message: ----------- specify endian-ness in dtype, instead of using byteswap() Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-08 22:06:40 UTC (rev 4164) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-08 22:11:50 UTC (rev 4165) @@ -791,10 +791,8 @@ # read in binary string convert into an npts by 2 # numpy array (first column is lons, second is lats). polystring = bdatfile.read(bytecount) - if not npy.little_endian: - b = npy.reshape(npy.fromstring(polystring,dtype=npy.float64).byteswap(),(npts,2)) - else: - b = npy.reshape(npy.fromstring(polystring,dtype=npy.float64),(npts,2)) + # binary data is little endian. + b = npy.reshape(npy.fromstring(polystring,dtype='<f8'),(npts,2)) # if map boundary polygon is a valid one in lat/lon # coordinates (i.e. it does not contain either pole), # the intersections of the boundary geometries This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-09 20:00:43
|
Revision: 4193 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4193&view=rev Author: jswhit Date: 2007-11-09 12:00:39 -0800 (Fri, 09 Nov 2007) Log Message: ----------- if poly intersection fails, attempt to use the polygon anyway Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-09 19:48:42 UTC (rev 4192) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-09 20:00:39 UTC (rev 4193) @@ -841,18 +841,12 @@ # if polygon instersects map projection # region, process it. if poly.intersects(boundarypolyll): - #if not poly.is_valid: - # print poly.geom_type, poly.is_ring, boundarypolyll.is_valid - # import pylab - # a = npy.asarray(boundarypolyll.boundary) - # b = npy.asarray(poly.boundary) - # pylab.plot(a[:,0],a[:,1],'b') - # pylab.plot(b[:,0],b[:,1],'b') - # pylab.show() - if poly.is_valid: + # if geometry intersection calculation fails, + # just move on. + try: poly = poly.intersection(boundarypolyll) - #else: - # print 'warning, invalid ',name,' geometry',poly.area + except: + pass # create iterable object with geometries # that intersect map region. if hasattr(poly,'geoms'): @@ -912,15 +906,22 @@ # if geometry instersects map projection # region, and doesn't have any invalid points, process it. if not badmask.any() and boundarypolyxy.intersects(poly): + #if not poly.is_valid: + # print poly.geom_type, poly.is_ring, boundarypolyll.is_valid + # import pylab + # a = npy.asarray(boundarypolyxy.boundary) + # b = npy.asarray(poly.boundary) + # poly2 = boundarypolyxy.intersection(poly) + # c = npy.asarray(poly2.boundary) + # pylab.plot(a[:,0],a[:,1],'b') + # pylab.fill(c[:,0],c[:,1],'r') + # pylab.show() # if geometry intersection calculation fails, - # just skip this polygon. - #try: - if poly.is_valid: + # just move on. + try: poly = boundarypolyxy.intersection(poly) - #else: - # print 'warning, invalid ',name,' geometry',poly.area - #except: - # continue + except: + pass # create iterable object with geometries # that intersect map region. if hasattr(poly,'geoms'): This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-09 22:11:18
|
Revision: 4201 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4201&view=rev Author: jswhit Date: 2007-11-09 14:11:14 -0800 (Fri, 09 Nov 2007) Log Message: ----------- check to see if Greenwich meridian is in map region. If not, some speedups can be had. Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-09 21:28:58 UTC (rev 4200) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-09 22:11:14 UTC (rev 4201) @@ -775,6 +775,15 @@ b = npy.asarray(self._boundarypolyll.boundary) blons = b[:,0]; blats = b[:,1] boundarypolyxy = PolygonShape(zip(*maptran(blons,blats))) + # for map regions that don't contain Pole, see if Greenwich + # meridian is contained. + if not containsPole: + lon_0,lat_0 = self(self.xmin+0.5*(self.xmax-self.xmin),\ + self.ymin+0.5*(self.ymax-self.ymin),inverse=True) + GM0 = PointShape((0,lat_0)) + GM360 = PointShape((360,lat_0)) + hasGM0 = GM0.within(boundarypolyll) + hasGM360 = GM360.within(boundarypolyll) for line in bdatmetafile: linesplit = line.split() area = float(linesplit[1]) @@ -830,11 +839,16 @@ # (so as to properly treat polygons that cross # Greenwich meridian). if not antart: - blons = b[:,0]-360 - poly1 = Shape(zip(blons,blats)) - blons = b[:,0]+360 - poly2 = Shape(zip(blons,blats)) - polys = [poly1,poly,poly2] + if hasGM0: + blons = b[:,0]-360 + poly1 = Shape(zip(blons,blats)) + polys = [poly1,poly] + elif hasGM360: + blons = b[:,0]+360 + poly2 = Shape(zip(blons,blats)) + polys = [poly,poly2] + else: + polys = [poly] else: # Antartica already extends from -360 to +720. polys = [poly] for poly in polys: @@ -906,16 +920,6 @@ # if geometry instersects map projection # region, and doesn't have any invalid points, process it. if not badmask.any() and boundarypolyxy.intersects(poly): - #if not poly.is_valid: - # print poly.geom_type, poly.is_ring, boundarypolyll.is_valid - # import pylab - # a = npy.asarray(boundarypolyxy.boundary) - # b = npy.asarray(poly.boundary) - # poly2 = boundarypolyxy.intersection(poly) - # c = npy.asarray(poly2.boundary) - # pylab.plot(a[:,0],a[:,1],'b') - # pylab.fill(c[:,0],c[:,1],'r') - # pylab.show() # if geometry intersection calculation fails, # just move on. try: This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-09 22:29:50
|
Revision: 4202 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4202&view=rev Author: jswhit Date: 2007-11-09 14:29:49 -0800 (Fri, 09 Nov 2007) Log Message: ----------- back off last commit, doesn't do the right thing in all cases. Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-09 22:11:14 UTC (rev 4201) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-09 22:29:49 UTC (rev 4202) @@ -775,15 +775,6 @@ b = npy.asarray(self._boundarypolyll.boundary) blons = b[:,0]; blats = b[:,1] boundarypolyxy = PolygonShape(zip(*maptran(blons,blats))) - # for map regions that don't contain Pole, see if Greenwich - # meridian is contained. - if not containsPole: - lon_0,lat_0 = self(self.xmin+0.5*(self.xmax-self.xmin),\ - self.ymin+0.5*(self.ymax-self.ymin),inverse=True) - GM0 = PointShape((0,lat_0)) - GM360 = PointShape((360,lat_0)) - hasGM0 = GM0.within(boundarypolyll) - hasGM360 = GM360.within(boundarypolyll) for line in bdatmetafile: linesplit = line.split() area = float(linesplit[1]) @@ -839,16 +830,11 @@ # (so as to properly treat polygons that cross # Greenwich meridian). if not antart: - if hasGM0: - blons = b[:,0]-360 - poly1 = Shape(zip(blons,blats)) - polys = [poly1,poly] - elif hasGM360: - blons = b[:,0]+360 - poly2 = Shape(zip(blons,blats)) - polys = [poly,poly2] - else: - polys = [poly] + blons = b[:,0]-360 + poly1 = Shape(zip(blons,blats)) + blons = b[:,0]+360 + poly2 = Shape(zip(blons,blats)) + polys = [poly1,poly,poly2] else: # Antartica already extends from -360 to +720. polys = [poly] for poly in polys: This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-11 21:12:00
|
Revision: 4219 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4219&view=rev Author: jswhit Date: 2007-11-11 13:11:35 -0800 (Sun, 11 Nov 2007) Log Message: ----------- use numpy arrays to create Shapely instances (10% speedup), allow for the possibility to use full res coastlines. Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-10 23:26:36 UTC (rev 4218) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-11 21:11:35 UTC (rev 4219) @@ -670,6 +670,8 @@ area_thresh = 100. elif resolution == 'h': area_thresh = 10. + elif resolution == 'f': + area_thresh = 1. else: raise ValueError, "boundary resolution must be one of 'c','l','i' or 'h'" self.area_thresh = area_thresh @@ -792,7 +794,9 @@ # numpy array (first column is lons, second is lats). polystring = bdatfile.read(bytecount) # binary data is little endian. - b = npy.reshape(npy.fromstring(polystring,dtype='<f4'),(npts,2)) + b = npy.array(npy.fromstring(polystring,dtype='<f4'),'f8') + b.shape = (npts,2) + b2 = b.copy() # if map boundary polygon is a valid one in lat/lon # coordinates (i.e. it does not contain either pole), # the intersections of the boundary geometries @@ -819,21 +823,20 @@ lats.append(-90.) poly = PolygonShape(zip(lons,lats)) antart = True - b = npy.empty((len(lons),2),npy.float32) + b = npy.empty((len(lons),2),npy.float64) b[:,0] = lons; b[:,1] = lats else: antart = False # create Shapely geometry from lons/lons array. - blons = b[:,0]; blats = b[:,1] - poly = Shape(zip(blons,blats)) + poly = Shape(b) # create duplicate polygons shifted by -360 and +360 # (so as to properly treat polygons that cross # Greenwich meridian). if not antart: - blons = b[:,0]-360 - poly1 = Shape(zip(blons,blats)) - blons = b[:,0]+360 - poly2 = Shape(zip(blons,blats)) + b2[:,0] = b[:,0]-360 + poly1 = Shape(b2) + b2[:,0] = b[:,0]+360 + poly2 = Shape(b2) polys = [poly1,poly,poly2] else: # Antartica already extends from -360 to +720. polys = [poly] @@ -874,17 +877,15 @@ # projection region and boundary geometries in map # projection coordinates. else: - blons = b[:,0]; blats = b[:,1] # transform coordinates from lat/lon # to map projection coordinates. # special case for ortho, compute coastline polygon # vertices in stereographic coords. if name == 'gshhs' and self.projection == 'ortho': - bx, by = maptran(blons, blats) + b[:,0], b[:,1] = maptran(b[:,0], b[:,1]) else: - bx, by = self(blons, blats) - goodmask = npy.logical_and(bx<1.e20,by<1.e20) - badmask = npy.logical_or(bx>1.e20,by>1.e20) + b[:,0], b[:,1] = self(b[:,0], b[:,1]) + goodmask = npy.logical_and(b[:,0]<1.e20,b[:,1]<1.e20) # if less than two points are valid in # map proj coords, skip this geometry. if npy.sum(goodmask) <= 1: continue @@ -892,8 +893,8 @@ # if not a polygon, # just remove parts of geometry that are undefined # in this map projection. - bx = npy.compress(goodmask, bx) - by = npy.compress(goodmask, by) + bx = npy.compress(goodmask, b[:,0]) + by = npy.compress(goodmask, b[:,1]) # for orthographic projection, all points # outside map projection region are eliminated # with the above step, so we're done. @@ -902,10 +903,10 @@ polygon_types.append(type) continue # create a Shapely geometry object. - poly = Shape(zip(bx,by)) + poly = Shape(b) # if geometry instersects map projection # region, and doesn't have any invalid points, process it. - if not badmask.any() and boundarypolyxy.intersects(poly): + if goodmask.all() and boundarypolyxy.intersects(poly): # if geometry intersection calculation fails, # just move on. try: @@ -935,14 +936,11 @@ if name == 'gshhs' and\ self._fulldisk and\ areafrac > 0.99: continue - bx = b[:,0]; by = b[:,1] # inverse transform from stereographic # to lat/lon. - blons, blats = maptran(bx, by, inverse=True) - # forward transform from lat/lon to + b[:,0], b[:,1] = maptran(b[:,0], b[:,1], inverse=True) # orthographic. - bx, by = self(blons, blats) - b[:,0] = bx; b[:,1] = by + b[:,0], b[:,1]= self(b[:,0], b[:,1]) polygons.append(zip(b[:,0],b[:,1])) polygon_types.append(type) return polygons, polygon_types This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-13 20:55:34
|
Revision: 4253 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4253&view=rev Author: jswhit Date: 2007-11-13 12:55:29 -0800 (Tue, 13 Nov 2007) Log Message: ----------- now use custom libgeos pyrex wrapper, Shapely no longer needed. Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-13 20:52:54 UTC (rev 4252) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-13 20:55:29 UTC (rev 4253) @@ -14,9 +14,9 @@ from numpy import linspace, squeeze, ma from matplotlib.cbook import popd, is_scalar from shapelib import ShapeFile -from shapely.geometry import Polygon as PolygonShape -from shapely.geometry import LineString as LineShape -from shapely.geometry import Point as PointShape +from geos import Polygon as PolygonShape +from geos import LineString as LineShape +from geos import Point as PointShape # basemap data files now installed in lib/matplotlib/toolkits/basemap/data basemap_datadir = os.sep.join([os.path.dirname(__file__), 'data']) @@ -774,9 +774,10 @@ maptran = pyproj.Proj(proj='stere',lon_0=lon_0,lat_0=lat_0) # boundary polygon for orthographic projection # in stereographic coorindates. - b = npy.asarray(self._boundarypolyll.boundary) + b = self._boundarypolyll.get_coords() blons = b[:,0]; blats = b[:,1] - boundarypolyxy = PolygonShape(zip(*maptran(blons,blats))) + b[:,0], b[:,1] = maptran(blons, blats) + boundarypolyxy = PolygonShape(b) for line in bdatmetafile: linesplit = line.split() area = float(linesplit[1]) @@ -821,14 +822,13 @@ lats.insert(0,-90.) lons.append(lonend) lats.append(-90.) - poly = PolygonShape(zip(lons,lats)) - antart = True b = npy.empty((len(lons),2),npy.float64) b[:,0] = lons; b[:,1] = lats + poly = PolygonShape(b) + antart = True else: + poly = Shape(b) antart = False - # create Shapely geometry from lons/lons array. - poly = Shape(b) # create duplicate polygons shifted by -360 and +360 # (so as to properly treat polygons that cross # Greenwich meridian). @@ -847,25 +847,16 @@ # if geometry intersection calculation fails, # just move on. try: - poly = poly.intersection(boundarypolyll) + geoms = poly.intersection(boundarypolyll) except: pass - # create iterable object with geometries - # that intersect map region. - if hasattr(poly,'geoms'): - geoms = poly.geoms - else: - geoms = [poly] # iterate over geometries in intersection. for psub in geoms: # only coastlines are polygons, # which have a 'boundary' attribute. # otherwise, use 'coords' attribute # to extract coordinates. - if name == 'gshhs': - b = npy.asarray(psub.boundary) - else: - b = npy.asarray(psub.coords) + b = psub.get_coords() blons = b[:,0]; blats = b[:,1] # transformation from lat/lon to # map projection coordinates. @@ -906,25 +897,16 @@ poly = Shape(b) # if geometry instersects map projection # region, and doesn't have any invalid points, process it. - if goodmask.all() and boundarypolyxy.intersects(poly): + if goodmask.all() and poly.intersects(boundarypolyxy): # if geometry intersection calculation fails, # just move on. try: - poly = boundarypolyxy.intersection(poly) + geoms = poly.intersection(boundarypolyxy) except: pass - # create iterable object with geometries - # that intersect map region. - if hasattr(poly,'geoms'): - geoms = poly.geoms - else: - geoms = [poly] # iterate over geometries in intersection. for psub in geoms: - if name == 'gshhs': - b = npy.asarray(psub.boundary) - else: - b = npy.asarray(psub.coords) + b = psub.get_coords() # if projection == 'ortho', # transform polygon from stereographic # to orthographic coordinates. @@ -932,7 +914,7 @@ # if coastline polygon covers more than 99% # of map region for fulldisk projection, # it's probably bogus, so skip it. - areafrac = psub.area/boundarypolyxy.area + areafrac = psub.area()/boundarypolyxy.area() if name == 'gshhs' and\ self._fulldisk and\ areafrac > 0.99: continue @@ -963,7 +945,9 @@ rmajor = self._width x = rmajor*npy.cos(thetas) + rmajor y = rminor*npy.sin(thetas) + rminor - boundaryxy = PolygonShape(zip(x,y)) + b = npy.empty((len(x),2),npy.float64) + b[:,0]=x; b[:,1]=y + boundaryxy = PolygonShape(b) # compute proj instance for full disk, if necessary. if not self._fulldisk: projparms = self.projparams @@ -1000,7 +984,9 @@ lons = npy.array(lons1+lons2+lons3+lons4,npy.float64) lats = npy.array(lats1+lats2+lats3+lats4,npy.float64) x, y = maptran(lons,lats) - boundaryxy = PolygonShape(zip(x,y)) + b = npy.empty((len(x),2),npy.float64) + b[:,0]=x; b[:,1]=y + boundaryxy = PolygonShape(b) else: # all other projections are rectangular. # left side (x = xmin, ymin <= y <= ymax) yy = linspace(self.ymin, self.ymax, ny)[:-1] @@ -1018,10 +1004,10 @@ y = y + len(xx)*[self.ymin] x = npy.array(x,npy.float64) y = npy.array(y,npy.float64) - boundaryxy = PolygonShape([(self.xmin,self.ymin),\ - (self.xmin,self.ymax),\ - (self.xmax,self.ymax),\ - (self.xmax,self.ymin)]) + b = npy.empty((4,2),npy.float64) + b[:,0]=[self.xmin,self.xmin,self.xmax,self.xmax] + b[:,1]=[self.ymin,self.ymax,self.ymax,self.ymin] + boundaryxy = PolygonShape(b) if self.projection in ['mill','merc','cyl']: # make sure map boundary doesn't quite include pole. if self.urcrnrlat > 89.9999: @@ -1035,7 +1021,9 @@ lons = [self.llcrnrlon, self.llcrnrlon, self.urcrnrlon, self.urcrnrlon] lats = [llcrnrlat, urcrnrlat, urcrnrlat, llcrnrlat] x, y = self(lons, lats) - boundaryxy = PolygonShape(zip(x,y)) + b = npy.empty((len(x),2),npy.float64) + b[:,0]=x; b[:,1]=y + boundaryxy = PolygonShape(b) else: if self.projection not in ['moll','robin','sinu']: lons, lats = maptran(x,y,inverse=True) @@ -1051,7 +1039,9 @@ lons[n] = lon lonprev = lon n = n + 1 - boundaryll = PolygonShape(zip(lons,lats)) + b = npy.empty((len(lons),2),npy.float64) + b[:,0]=lons; b[:,1]=lats + boundaryll = PolygonShape(b) return boundaryll, boundaryxy This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-13 21:05:04
|
Revision: 4254 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4254&view=rev Author: jswhit Date: 2007-11-13 13:04:59 -0800 (Tue, 13 Nov 2007) Log Message: ----------- remove docstring reference to Shapely Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-13 20:55:29 UTC (rev 4253) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-13 21:04:59 UTC (rev 4254) @@ -893,7 +893,7 @@ polygons.append(zip(bx,by)) polygon_types.append(type) continue - # create a Shapely geometry object. + # create a GEOS geometry object. poly = Shape(b) # if geometry instersects map projection # region, and doesn't have any invalid points, process it. This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-13 23:29:43
|
Revision: 4262 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4262&view=rev Author: jswhit Date: 2007-11-13 15:29:42 -0800 (Tue, 13 Nov 2007) Log Message: ----------- fix deleterious effects of global search and replace Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-13 23:11:28 UTC (rev 4261) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-13 23:29:42 UTC (rev 4262) @@ -29,7 +29,7 @@ transverse mercator, miller cylindrical, lambert conformal conic, azimuthal equidistant, equidistant conic, lambert azimuthal equal area, albers equal area conic, gnomonic, orthographic, sinusoidal, mollweide, - _geostationary, robinson, cassini-soldner or stereographic). + geostationary, robinson, cassini-soldner or stereographic). Doesn't actually draw anything, but sets up the map projection class and creates the coastline, lake river and political boundary data structures in native map projection coordinates. @@ -40,7 +40,7 @@ projection - map projection ('cyl','merc','mill','lcc','eqdc','aea', 'laea', 'nplaea', 'splaea', 'tmerc', 'omerc', 'cass', 'gnom', 'poly', 'sinu', 'moll', 'ortho', 'robin', 'aeqd', 'npaeqd', 'spaeqd', 'stere', - '_geos', 'npstere' or 'spstere') + 'geos', 'npstere' or 'spstere') (projections prefixed with 'np' or 'sp' are special case polar-centric versions of the parent projection) aspect - map aspect ratio (size of y dimension / size of x dimension). @@ -110,7 +110,7 @@ 'cass' - cassini-soldner (transverse cylindrical equidistant), 'poly' - polyconic, 'omerc' - oblique mercator, 'ortho' - orthographic, 'sinu' - sinusoidal, 'moll' - mollweide, 'robin' - robinson, - '_geos' - _geostationary, and 'gnom' - gnomonic are currently available. + 'geos' - geostationary, and 'gnom' - gnomonic are currently available. Default is 'cyl'. The map projection region can either be specified by setting these keywords: @@ -133,9 +133,9 @@ either they are computed internally, or entire globe is always plotted). For the cylindrical projections ('cyl','merc' and 'mill'), the default is to use llcrnrlon=-180,llcrnrlat=-90, urcrnrlon=180 and urcrnrlat=90). For all other - projections except 'ortho' and '_geos', either the lat/lon values of the + projections except 'ortho' and 'geos', either the lat/lon values of the corners or width and height must be specified by the user. - For 'ortho' and '_geos', the lat/lon values of the corners may be specified, + For 'ortho' and 'geos', the lat/lon values of the corners may be specified, but if they are not, the entire globe is plotted. resolution - resolution of boundary database to use. Can be 'c' (crude), @@ -208,7 +208,7 @@ on the north or south pole. The longitude lon_0 is at 6-o'clock, and the latitude circle boundinglat is tangent to the edge of the map at lon_0. satellite_height - height of satellite (in m) above equator - - only relevant for _geostationary projections ('_geos'). + only relevant for geostationary projections ('geos'). """ @@ -502,7 +502,7 @@ # FIXME: won't work for points exactly on equator?? if npy.abs(lat_0) < 1.e-2: lat_0 = 1.e-2 projparams['lat_0'] = lat_0 - elif projection == '_geos': + elif projection == 'geos': if lon_0 is None and satellite_height is None: raise ValueError, 'must specify lon_0 and satellite_height for Geostationary basemap' if width is not None or height is not None: @@ -591,7 +591,7 @@ 'splaea' - lambert azimuthal, special case centered on south pole, 'cass' - cassini-soldner (transverse cylindrical equidistant), 'poly' - polyconic, 'omerc' - oblique mercator, 'ortho' - orthographic, - '_geos' - _geostationary, 'sinu' - sinusoidal, 'moll' - mollweide, + 'geos' - geostationary, 'sinu' - sinusoidal, 'moll' - mollweide, 'robin' - robinson, or 'gnom' - gnomonic. You tried '%s'""" % projection # initialize proj4 @@ -606,7 +606,7 @@ atts = ['rmajor','rminor','esq','flattening','ellipsoid','projparams'] for att in atts: self.__dict__[att] = proj.__dict__[att] - # these only exist for _geostationary projection. + # these only exist for geostationary projection. if hasattr(proj,'_width'): self.__dict__['_width'] = proj.__dict__['_width'] if hasattr(proj,'_height'): @@ -645,7 +645,7 @@ if projection in ['mill','cyl','merc']: self.latmin = self.llcrnrlat self.latmax = self.urcrnrlat - elif projection in ['ortho','_geos','moll','robin','sinu']: + elif projection in ['ortho','geos','moll','robin','sinu']: self.latmin = -90. self.latmax = 90. else: @@ -759,7 +759,7 @@ containsPole = hasNP or hasSP # these projections cannot cross pole. if containsPole and\ - self.projection in ['tmerc','cass','omerc','merc','mill','cyl','robin','moll','sinu','_geos']: + self.projection in ['tmerc','cass','omerc','merc','mill','cyl','robin','moll','sinu','geos']: raise ValueError('%s projection cannot cross pole'%(self.projection)) # make sure orthographic projection has containsPole=True @@ -932,7 +932,7 @@ nx = 100 ny = 100 maptran = self - if self.projection in ['ortho','_geos']: + if self.projection in ['ortho','geos']: # circular region. thetas = linspace(0.,2.*npy.pi,2*nx*ny)[:-1] if self.projection == 'ortho': @@ -1065,7 +1065,7 @@ circle.set_edgecolor(color) circle.set_linewidth(linewidth) circle.set_clip_on(False) - elif self.projection == '_geos' and self._fulldisk: # elliptical region + elif self.projection == 'geos' and self._fulldisk: # elliptical region # define an Ellipse patch, add it to axes instance. ellps = Ellipse((self._width,self._height),2.*self._width,2.*self._height) ax.add_patch(ellps) @@ -1507,8 +1507,8 @@ l.set_zorder(zorder) ax.add_line(l) # draw labels for parallels - # parallels not labelled for fulldisk orthographic or _geostationary - if self.projection in ['ortho','_geos'] and max(labels): + # parallels not labelled for fulldisk orthographic or geostationary + if self.projection in ['ortho','geos'] and max(labels): if self._fulldisk: print 'Warning: Cannot label parallels on full-disk Orthographic or Geostationary basemap' labels = [0,0,0,0] @@ -1714,11 +1714,11 @@ ax.add_line(l) # draw labels for meridians. # meridians not labelled for sinusoidal, mollweide, or - # or full-disk orthographic/_geostationary. + # or full-disk orthographic/geostationary. if self.projection in ['sinu','moll'] and max(labels): print 'Warning: Cannot label meridians on Sinusoidal or Mollweide basemap' labels = [0,0,0,0] - if self.projection in ['ortho','_geos'] and max(labels): + if self.projection in ['ortho','geos'] and max(labels): if self._fulldisk: print 'Warning: Cannot label meridians on full-disk Geostationary or Orthographic basemap' labels = [0,0,0,0] @@ -2072,7 +2072,7 @@ # turn off axes frame for non-rectangular projections. if self.projection in ['moll','robin','sinu']: ax.set_frame_on(False) - if self.projection in ['ortho','_geos'] and self._fulldisk: + if self.projection in ['ortho','geos'] and self._fulldisk: ax.set_frame_on(False) # make sure aspect ratio of map preserved. # plot is re-centered in bounding rectangle. This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-14 16:28:00
|
Revision: 4273 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4273&view=rev Author: jswhit Date: 2007-11-14 08:27:55 -0800 (Wed, 14 Nov 2007) Log Message: ----------- support for full resolution datasets. Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-14 16:17:55 UTC (rev 4272) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-14 16:27:55 UTC (rev 4273) @@ -139,7 +139,7 @@ but if they are not, the entire globe is plotted. resolution - resolution of boundary database to use. Can be 'c' (crude), - 'l' (low), 'i' (intermediate), 'h' (high), or None. Default is 'c'. + 'l' (low), 'i' (intermediate), 'h' (high), 'f' (full) or None. If None, no boundary data will be read in (and class methods such as drawcoastlines will raise an exception if invoked). Resolution drops off by roughly 80% @@ -150,8 +150,8 @@ Tools (http://gmt.soest.hawaii.edu). area_thresh - coastline or lake with an area smaller than area_thresh - in km^2 will not be plotted. Default 10000,1000,100,10 for resolution - 'c','l','i','h'. + in km^2 will not be plotted. Default 10000,1000,100,10,1 for resolution + 'c','l','i','h','f'. rsphere - radius of the sphere used to define map projection (default 6370997 meters, close to the arithmetic mean radius of the earth). If @@ -671,7 +671,7 @@ elif resolution == 'f': area_thresh = 1. else: - raise ValueError, "boundary resolution must be one of 'c','l','i' or 'h'" + raise ValueError, "boundary resolution must be one of 'c','l','i','h' or 'f'" self.area_thresh = area_thresh # define map boundary polygon (in lat/lon coordinates) self._boundarypolyll, self._boundarypolyxy = self._getmapboundary() @@ -733,10 +733,10 @@ def _readboundarydata(self,name): msg = """ -Unable to open boundary dataset file. Only the 'crude', 'low' -and 'intermediate' resolution datasets are installed by default. If you -are requesting a 'high' resolution dataset, you need to download -and install those files manually (see the basemap README for details).""" +Unable to open boundary dataset file. Only the 'crude', 'low', +'intermediate' and 'high' resolution datasets are installed by default. If you +are requesting a 'full' resolution dataset, you need to download +and install those files separately(see the basemap README for details).""" try: bdatfile = open(os.path.join(basemap_datadir,name+'_'+self.resolution+'.dat'),'rb') bdatmetafile = open(os.path.join(basemap_datadir,name+'meta_'+self.resolution+'.dat'),'r') This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-16 00:12:03
|
Revision: 4327 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4327&view=rev Author: jswhit Date: 2007-11-15 16:11:50 -0800 (Thu, 15 Nov 2007) Log Message: ----------- refine definition of Antarctica Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-15 21:33:13 UTC (rev 4326) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-16 00:11:50 UTC (rev 4327) @@ -805,7 +805,7 @@ # regions and high-resolution boundary geometries). if not containsPole: # close Antarctica. - if name == 'gshhs' and south < -68: + if name == 'gshhs' and south < -68 and area > 10000: lons = b[:,0] lats = b[:,1] lons2 = lons[:-2][::-1] This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-17 14:46:53
|
Revision: 4350 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4350&view=rev Author: jswhit Date: 2007-11-17 05:22:47 -0800 (Sat, 17 Nov 2007) Log Message: ----------- raise exception if width, height keywords used with omerc. Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-17 13:08:05 UTC (rev 4349) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-17 13:22:47 UTC (rev 4350) @@ -541,7 +541,7 @@ projparams['lon_2'] = lon_2 if None in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]: if width is None or height is None: - raise ValueError, 'must either specify lat/lon values of corners (llcrnrlon,llcrnrlat,ucrnrlon,urcrnrlat) in degrees or width and height in meters' + raise ValueError, 'cannot specify map region with width and height keywords for this projection, please specify lat/lon values of corners' else: if lon_0 is None or lat_0 is None: raise ValueError, 'must specify lon_0 and lat_0 when using width, height to specify projection region' This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-17 14:55:37
|
Revision: 4351 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4351&view=rev Author: jswhit Date: 2007-11-17 05:29:15 -0800 (Sat, 17 Nov 2007) Log Message: ----------- fix omerc projection parameter checking Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-17 13:22:47 UTC (rev 4350) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-17 13:29:15 UTC (rev 4351) @@ -540,14 +540,7 @@ projparams['lat_2'] = lat_2 projparams['lon_2'] = lon_2 if None in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]: - if width is None or height is None: - raise ValueError, 'cannot specify map region with width and height keywords for this projection, please specify lat/lon values of corners' - else: - if lon_0 is None or lat_0 is None: - raise ValueError, 'must specify lon_0 and lat_0 when using width, height to specify projection region' - llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat = _choosecorners(width,height,**projparams) - self.llcrnrlon = llcrnrlon; self.llcrnrlat = llcrnrlat - self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat + raise ValueError, 'cannot specify map region with width and height keywords for this projection, please specify lat/lon values of corners' elif projection == 'aeqd': if lat_0 is None or lon_0 is None: raise ValueError, 'must specify lat_0 and lon_0 for Azimuthal Equidistant basemap' This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-17 15:17:31
|
Revision: 4356 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4356&view=rev Author: jswhit Date: 2007-11-17 07:17:30 -0800 (Sat, 17 Nov 2007) Log Message: ----------- fix orthographic bug Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-17 14:40:01 UTC (rev 4355) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-17 15:17:30 UTC (rev 4356) @@ -766,7 +766,8 @@ containsPole = True lon_0=self.projparams['lon_0'] lat_0=self.projparams['lat_0'] - maptran = pyproj.Proj(proj='stere',lon_0=lon_0,lat_0=lat_0) + re = self.projparams['R'] + maptran = pyproj.Proj(proj='stere',lon_0=lon_0,lat_0=lat_0,R=re) # boundary polygon for orthographic projection # in stereographic coorindates. b = self._boundarypolyll.boundary This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-18 13:32:30
|
Revision: 4370 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4370&view=rev Author: jswhit Date: 2007-11-18 05:32:27 -0800 (Sun, 18 Nov 2007) Log Message: ----------- port over some of Eric's refactoring from trunk Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-17 22:11:19 UTC (rev 4369) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-18 13:32:27 UTC (rev 4370) @@ -275,8 +275,8 @@ projparams['lat_ts'] = lat_ts if satellite_height is not None: projparams['h'] = satellite_height - - if None not in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]: + using_corners = (None not in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]) + if using_corners: # make sure lat/lon limits are converted to floats. self.llcrnrlon = float(llcrnrlon) self.llcrnrlat = float(llcrnrlat) @@ -303,7 +303,7 @@ raise ValueError, 'must specify lat_1 or lat_0 and lon_0 for Lambert Conformal basemap (lat_2 is optional)' if lat_2 is None: projparams['lat_2'] = lat_1 - if None in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]: + if not using_corners: if width is None or height is None: raise ValueError, 'must either specify lat/lon values of corners (llcrnrlon,llcrnrlat,ucrnrlon,urcrnrlat) in degrees or width and height in meters' else: @@ -323,7 +323,7 @@ raise ValueError, 'must specify lat_1 or lat_0 and lon_0 for Equidistant Conic basemap (lat_2 is optional)' if lat_2 is None: projparams['lat_2'] = lat_1 - if None in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]: + if not using_corners: if width is None or height is None: raise ValueError, 'must either specify lat/lon values of corners (llcrnrlon,llcrnrlat,ucrnrlon,urcrnrlat) in degrees or width and height in meters' else: @@ -342,7 +342,7 @@ raise ValueError, 'must specify lat_1 or lat_0 and lon_0 for Albers Equal Area basemap (lat_2 is optional)' if lat_2 is None: projparams['lat_2'] = lat_1 - if None in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]: + if not using_corners: if width is None or height is None: raise ValueError, 'must either specify lat/lon values of corners (llcrnrlon,llcrnrlat,ucrnrlon,urcrnrlat) in degrees or width and height in meters' else: @@ -354,7 +354,7 @@ elif projection == 'stere': if lat_0 is None or lon_0 is None: raise ValueError, 'must specify lat_0 and lon_0 for Stereographic basemap (lat_ts is optional)' - if None in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]: + if not using_corners: if width is None or height is None: raise ValueError, 'must either specify lat/lon values of corners (llcrnrlon,llcrnrlat,ucrnrlon,urcrnrlat) in degrees or width and height in meters' else: @@ -446,7 +446,7 @@ elif projection == 'laea': if lat_0 is None or lon_0 is None: raise ValueError, 'must specify lat_0 and lon_0 for Lambert Azimuthal basemap' - if None in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]: + if not using_corners: if width is None or height is None: raise ValueError, 'must either specify lat/lon values of corners (llcrnrlon,llcrnrlat,ucrnrlon,urcrnrlat) in degrees or width and height in meters' else: @@ -460,7 +460,7 @@ raise ValueError, 'must specify lat_ts for Mercator basemap' # clip plot region to be within -89.99S to 89.99N # (mercator is singular at poles) - if None in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]: + if not using_corners: llcrnrlon = -180. llcrnrlat = -90. urcrnrlon = 180 @@ -476,7 +476,7 @@ elif projection in ['tmerc','gnom','cass','poly'] : if lat_0 is None or lon_0 is None: raise ValueError, 'must specify lat_0 and lon_0 for Transverse Mercator, Gnomonic, Cassini-Soldnerr Polyconic basemap' - if None in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]: + if not using_corners: if width is None or height is None: raise ValueError, 'must either specify lat/lon values of corners (llcrnrlon,llcrnrlat,ucrnrlon,urcrnrlat) in degrees or width and height in meters' else: @@ -493,7 +493,7 @@ raise ValueError, 'must specify lat_0 and lon_0 for Orthographic basemap' if width is not None or height is not None: print 'warning: width and height keywords ignored for %s projection' % self.projection - if None in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]: + if not using_corners: llcrnrlon = -180. llcrnrlat = -90. urcrnrlon = 180 @@ -511,7 +511,7 @@ raise ValueError, 'must specify lon_0 and satellite_height for Geostationary basemap' if width is not None or height is not None: print 'warning: width and height keywords ignored for %s projection' % self.projection - if None in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]: + if not using_corners: llcrnrlon = -180. llcrnrlat = -90. urcrnrlon = 180 @@ -539,12 +539,12 @@ projparams['lon_1'] = lon_1 projparams['lat_2'] = lat_2 projparams['lon_2'] = lon_2 - if None in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]: + if not using_corners: raise ValueError, 'cannot specify map region with width and height keywords for this projection, please specify lat/lon values of corners' elif projection == 'aeqd': if lat_0 is None or lon_0 is None: raise ValueError, 'must specify lat_0 and lon_0 for Azimuthal Equidistant basemap' - if None in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]: + if not using_corners: if width is None or height is None: raise ValueError, 'must either specify lat/lon values of corners (llcrnrlon,llcrnrlat,ucrnrlon,urcrnrlat) in degrees or width and height in meters' else: @@ -554,7 +554,7 @@ self.llcrnrlon = llcrnrlon; self.llcrnrlat = llcrnrlat self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat elif projection == 'mill': - if None in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]: + if not using_corners: llcrnrlon = -180. llcrnrlat = -90. urcrnrlon = 180 @@ -564,7 +564,7 @@ if width is not None or height is not None: print 'warning: width and height keywords ignored for %s projection' % self.projection elif projection == 'cyl': - if None in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]: + if not using_corners: llcrnrlon = -180. llcrnrlat = -90. urcrnrlon = 180 This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-18 14:31:10
|
Revision: 4371 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4371&view=rev Author: jswhit Date: 2007-11-18 06:31:06 -0800 (Sun, 18 Nov 2007) Log Message: ----------- port over the rest of Eric's changes from trunk Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-18 13:32:27 UTC (rev 4370) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-18 14:31:06 UTC (rev 4371) @@ -12,7 +12,7 @@ from proj import Proj import numpy as npy from numpy import linspace, squeeze, ma -from matplotlib.cbook import popd, is_scalar +from matplotlib.cbook import is_scalar, dedent from shapelib import ShapeFile import _geos @@ -21,77 +21,10 @@ __version__ = '0.9.7' -class Basemap(object): - - """ - Set up a basemap with one of 19 supported map projections - (cylindrical equidistant, mercator, polyconic, oblique mercator, - transverse mercator, miller cylindrical, lambert conformal conic, - azimuthal equidistant, equidistant conic, lambert azimuthal equal area, - albers equal area conic, gnomonic, orthographic, sinusoidal, mollweide, - geostationary, robinson, cassini-soldner or stereographic). - Doesn't actually draw anything, but sets up the map projection class and - creates the coastline, lake river and political boundary data - structures in native map projection coordinates. - Uses a pyrex interface to C-code from proj.4 (http://proj.maptools.org). - - Useful instance variables: - - projection - map projection ('cyl','merc','mill','lcc','eqdc','aea', - 'laea', 'nplaea', 'splaea', 'tmerc', 'omerc', 'cass', 'gnom', 'poly', - 'sinu', 'moll', 'ortho', 'robin', 'aeqd', 'npaeqd', 'spaeqd', 'stere', - 'geos', 'npstere' or 'spstere') - (projections prefixed with 'np' or 'sp' are special case polar-centric - versions of the parent projection) - aspect - map aspect ratio (size of y dimension / size of x dimension). - llcrnrlon - longitude of lower left hand corner of the desired map domain. - llcrnrlon - latitude of lower left hand corner of the desired map domain. - urcrnrlon - longitude of upper right hand corner of the desired map domain. - urcrnrlon - latitude of upper right hand corner of the desired map domain. - llcrnrx,llcrnry,urcrnrx,urcrnry - corners of map domain in projection coordinates. - rmajor,rminor - equatorial and polar radii of ellipsoid used (in meters). - resolution - resolution of boundary dataset being used ('c' for crude, - 'l' for low, etc.). If None, no boundary dataset is associated with the - Basemap instance. - srs - a string representing the 'spatial reference system' for the map - projection as defined by PROJ.4. - - Example Usage: - ->>> from matplotlib.toolkits.basemap import Basemap ->>> from pylab import load, meshgrid, title, arange, show ->>> # read in topo data (on a regular lat/lon grid) ->>> etopo = load('etopo20data.gz') ->>> lons = load('etopo20lons.gz') ->>> lats = load('etopo20lats.gz') ->>> # create Basemap instance for Robinson projection. ->>> m = Basemap(projection='robin',lon_0=0.5*(lons[0]+lons[-1])) ->>> # compute native map projection coordinates for lat/lon grid. ->>> x, y = m(*meshgrid(lons,lats)) ->>> # make filled contour plot. ->>> cs = m.contourf(x,y,etopo,30,cmap=cm.jet) ->>> m.drawcoastlines() # draw coastlines ->>> m.drawmapboundary() # draw a line around the map region ->>> m.drawparallels(arange(-90.,120.,30.),labels=[1,0,0,0]) # draw parallels ->>> m.drawmeridians(arange(0.,420.,60.),labels=[0,0,0,1]) # draw meridians ->>> title('Robinson Projection') # add a title ->>> show() - - [this example (simpletest.py) plus many others can be found in the - examples directory of source distribution. The "OO" version of this - example (which does not use pylab) is called "simpletest_oo.py".] - - Contact: Jeff Whitaker <jef...@no...> - """ - - def __init__(self,llcrnrlon=None,llcrnrlat=None, - urcrnrlon=None,urcrnrlat=None,\ - width=None,height=None,\ - projection='cyl',resolution='c',area_thresh=None,rsphere=6370997.0,\ - lat_ts=None,lat_1=None,lat_2=None,lat_0=None,lon_0=None,\ - lon_1=None,lon_2=None,suppress_ticks=True,\ - satellite_height=None,boundinglat=None,anchor='C',ax=None): - """ +# The __init__ docstring is pulled out here because it is so long; +# Having it in the usual place makes it hard to get from the +# __init__ argument list to the code that uses the arguments. +_Basemap_init_doc = """ create a Basemap instance. arguments: @@ -127,19 +60,19 @@ lon_0 - center of desired map domain (in degrees). lat_0 - center of desired map domain (in degrees). - For 'sinu', 'moll', 'npstere', 'spstere', 'nplaea', 'splaea', 'nplaea', - 'splaea', 'npaeqd', 'spaeqd' or 'robin', the values of - llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat,width and height are ignored (because - either they are computed internally, or entire globe is always plotted). For the - cylindrical projections ('cyl','merc' and 'mill'), the default is to use - llcrnrlon=-180,llcrnrlat=-90, urcrnrlon=180 and urcrnrlat=90). For all other + For 'sinu', 'moll', 'npstere', 'spstere', 'nplaea', 'splaea', 'nplaea', + 'splaea', 'npaeqd', 'spaeqd' or 'robin', the values of + llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat,width and height are ignored (because + either they are computed internally, or entire globe is always plotted). For the + cylindrical projections ('cyl','merc' and 'mill'), the default is to use + llcrnrlon=-180,llcrnrlat=-90, urcrnrlon=180 and urcrnrlat=90). For all other projections except 'ortho' and 'geos', either the lat/lon values of the corners or width and height must be specified by the user. For 'ortho' and 'geos', the lat/lon values of the corners may be specified, but if they are not, the entire globe is plotted. resolution - resolution of boundary database to use. Can be 'c' (crude), - 'l' (low), 'i' (intermediate), 'h' (high), 'f' (full) or None. + 'l' (low), 'i' (intermediate), 'h' (high), or None. Default is 'c'. If None, no boundary data will be read in (and class methods such as drawcoastlines will raise an exception if invoked). Resolution drops off by roughly 80% @@ -150,8 +83,8 @@ Tools (http://gmt.soest.hawaii.edu). area_thresh - coastline or lake with an area smaller than area_thresh - in km^2 will not be plotted. Default 10000,1000,100,10,1 for resolution - 'c','l','i','h','f'. + in km^2 will not be plotted. Default 10000,1000,100,10 for resolution + 'c','l','i','h'. rsphere - radius of the sphere used to define map projection (default 6370997 meters, close to the arithmetic mean radius of the earth). If @@ -185,7 +118,7 @@ The following parameters are map projection parameters which all default to None. Not all parameters are used by all projections, some are ignored. - lat_ts - latitude of natural origin (used for mercator, and + lat_ts - latitude of natural origin (used for mercator, and optionally for stereographic projection). lat_1 - first standard parallel for lambert conformal, albers equal area projection and equidistant conic projections. Latitude of one @@ -195,13 +128,13 @@ lat_2 - second standard parallel for lambert conformal, albers equal area projection and equidistant conic projections. Latitude of one of the two points on the projection centerline for oblique mercator. - If lat_2 is not given, it is set to lat_1 for + If lat_2 is not given, it is set to lat_1 for lambert conformal, albers equal area and equidistant conic. lon_1 - longitude of one of the two points on the projection centerline for oblique mercator. lon_2 - longitude of one of the two points on the projection centerline for oblique mercator. - lat_0 - central latitude (y-axis origin) - used by all projections, + lat_0 - central latitude (y-axis origin) - used by all projections, lon_0 - central meridian (x-axis origin) - used by all projections, boundinglat - bounding latitude for pole-centered projections (npstere,spstere, nplaea,splaea,npaeqd,spaeqd). These projections are square regions centered @@ -209,10 +142,146 @@ latitude circle boundinglat is tangent to the edge of the map at lon_0. satellite_height - height of satellite (in m) above equator - only relevant for geostationary projections ('geos'). - + """ +_unsupported_projection = """ + unsupported projection, use 'cyl' - cylindrical equidistant, 'merc' - + mercator, 'lcc' - lambert conformal conic, 'stere' - stereographic, + 'npstere' - stereographic, special case centered on north pole. + 'spstere' - stereographic, special case centered on south pole, + 'aea' - albers equal area conic, 'tmerc' - transverse mercator, + 'aeqd' - azimuthal equidistant, 'mill' - miller cylindrical, + 'npaeqd' - azimuthal equidistant, special case centered on north pole, + 'spaeqd' - azimuthal equidistant, special case centered on south pole, + 'eqdc' - equidistant conic, 'laea' - lambert azimuthal equal area, + 'nplaea' - lambert azimuthal, special case centered on north pole, + 'splaea' - lambert azimuthal, special case centered on south pole, + 'cass' - cassini-soldner (transverse cylindrical equidistant), + 'poly' - polyconic, 'omerc' - oblique mercator, 'ortho' - orthographic, + 'geos' - geostationary, 'sinu' - sinusoidal, 'moll' - mollweide, + 'robin' - robinson, or 'gnom' - gnomonic. You tried '%s' + """ + +# This allows substitution of longer names into error messages. +projnames = {'cyl' : 'Cylindrical Equidistant', + 'merc' : 'Mercator', + 'tmerc' : 'Transverse Mercator', + 'omerc' : 'Oblique Mercator', + 'mill' : 'Miller Cylindrical', + 'llc' : 'Lambert Conformal', + 'laea' : 'Lambert Azimuthal Equal Area', + 'nplaea' : 'North-Polar Lambert Azimuthal', + 'splaea' : 'South-Polar Lambert Azimuthal', + 'eqdc' : 'Equidistant Conic', + 'eaqd' : 'Azimuthal Equidistant', + 'npaeqd' : 'North-Polar Azimuthal Equidistant', + 'spaeqd' : 'South-Polar Azimuthal Equidistant', + 'aea' : 'Albers Equal Area', + 'stere' : 'Stereographic', + 'npstere' : 'Nouth-Polar Stereographic', + 'spstere' : 'South-Polar Stereographic', + 'cass' : 'Cassini-Soldner', + 'poly' : 'Polyconic', + 'ortho' : 'Orthographic', + 'geos' : 'Geostationary', + 'sinu' : 'Sinusoidal', + 'moll' : 'Mollweide', + 'robin' : 'Robinson', + 'gnom' : 'Gnomonic', + } + +def _validated_ll(param, name, minval, maxval): + param = float(param) + if param > maxval or param < minval: + raise ValueError('%s must be between %f and %f degrees' % + (name, minval, maxval)) + return param + +def _insert_validated(d, param, name, minval, maxval): + if param is not None: + d[name] = _validated_ll(param, name, minval, maxval) + +class Basemap(object): + """ + Set up a basemap with one of 19 supported map projections + (cylindrical equidistant, mercator, polyconic, oblique mercator, + transverse mercator, miller cylindrical, lambert conformal conic, + azimuthal equidistant, equidistant conic, lambert azimuthal equal area, + albers equal area conic, gnomonic, orthographic, sinusoidal, mollweide, + geostationary, robinson, cassini-soldner or stereographic). + Doesn't actually draw anything, but sets up the map projection class and + creates the coastline, lake river and political boundary data + structures in native map projection coordinates. + Uses a pyrex interface to C-code from proj.4 (http://proj.maptools.org). + + Useful instance variables: + + projection - map projection ('cyl','merc','mill','lcc','eqdc','aea', + 'laea', 'nplaea', 'splaea', 'tmerc', 'omerc', 'cass', 'gnom', 'poly', + 'sinu', 'moll', 'ortho', 'robin', 'aeqd', 'npaeqd', 'spaeqd', 'stere', + 'geos', 'npstere' or 'spstere') + (projections prefixed with 'np' or 'sp' are special case polar-centric + versions of the parent projection) + aspect - map aspect ratio (size of y dimension / size of x dimension). + llcrnrlon - longitude of lower left hand corner of the desired map domain. + llcrnrlon - latitude of lower left hand corner of the desired map domain. + urcrnrlon - longitude of upper right hand corner of the desired map domain. + urcrnrlon - latitude of upper right hand corner of the desired map domain. + llcrnrx,llcrnry,urcrnrx,urcrnry - corners of map domain in projection coordinates. + rmajor,rminor - equatorial and polar radii of ellipsoid used (in meters). + resolution - resolution of boundary dataset being used ('c' for crude, + 'l' for low, etc.). If None, no boundary dataset is associated with the + Basemap instance. + srs - a string representing the 'spatial reference system' for the map + projection as defined by PROJ.4. + + Example Usage: + + >>> from matplotlib.toolkits.basemap import Basemap + >>> from pylab import load, meshgrid, title, arange, show + >>> # read in topo data (on a regular lat/lon grid) + >>> etopo = load('etopo20data.gz') + >>> lons = load('etopo20lons.gz') + >>> lats = load('etopo20lats.gz') + >>> # create Basemap instance for Robinson projection. + >>> m = Basemap(projection='robin',lon_0=0.5*(lons[0]+lons[-1])) + >>> # compute native map projection coordinates for lat/lon grid. + >>> x, y = m(*meshgrid(lons,lats)) + >>> # make filled contour plot. + >>> cs = m.contourf(x,y,etopo,30,cmap=cm.jet) + >>> m.drawcoastlines() # draw coastlines + >>> m.drawmapboundary() # draw a line around the map region + >>> m.drawparallels(arange(-90.,120.,30.),labels=[1,0,0,0]) # draw parallels + >>> m.drawmeridians(arange(0.,420.,60.),labels=[0,0,0,1]) # draw meridians + >>> title('Robinson Projection') # add a title + >>> show() + + [this example (simpletest.py) plus many others can be found in the + examples directory of source distribution. The "OO" version of this + example (which does not use pylab) is called "simpletest_oo.py".] + + Contact: Jeff Whitaker <jef...@no...> + """ + + + def __init__(self, llcrnrlon=None, llcrnrlat=None, + urcrnrlon=None, urcrnrlat=None, + width=None, height=None, + projection='cyl', resolution='c', + area_thresh=None, rsphere=6370997.0, + lat_ts=None, + lat_1=None, lat_2=None, + lat_0=None, lon_0=None, + lon_1=None, lon_2=None, + suppress_ticks=True, + satellite_height=None, + boundinglat=None, + anchor='C', + ax=None): + # docstring is added after definition + #print "starting: ", time.clock() # where to put plot in figure (default is 'C' or center) self.anchor = anchor # map projection. @@ -238,223 +307,91 @@ # set units to meters. projparams['units']='m' # check for sane values of lon_0, lat_0, lat_ts, lat_1, lat_2 - if lat_0 is not None: - if lat_0 > 90. or lat_0 < -90.: - raise ValueError, 'lat_0 must be between -90 and +90 degrees' - else: - projparams['lat_0'] = lat_0 - if lon_0 is not None: - if lon_0 < -720. or lon_0 > 720.: - raise ValueError, 'lon_0 must be between -720 and +720 degrees' - else: - projparams['lon_0'] = lon_0 - if lon_1 is not None: - if lon_1 < -720. or lon_1 > 720.: - raise ValueError, 'lon_1 must be between -720 and +720 degrees' - else: - projparams['lon_1'] = lon_1 - if lon_2 is not None: - if lon_2 < -720. or lon_2 > 720.: - raise ValueError, 'lon_2 must be between -720 and +720 degrees' - else: - projparams['lon_2'] = lon_2 - if lat_1 is not None: - if lat_1 > 90. or lat_1 < -90.: - raise ValueError, 'lat_1 must be between -90 and +90 degrees' - else: - projparams['lat_1'] = lat_1 - if lat_2 is not None: - if lat_2 > 90. or lat_2 < -90.: - raise ValueError, 'lat_2 must be between -90 and +90 degrees' - else: - projparams['lat_2'] = lat_2 - if lat_ts is not None: - if lat_ts > 90. or lat_ts < -90.: - raise ValueError, 'lat_ts must be between -90 and +90 degrees' - else: - projparams['lat_ts'] = lat_ts + _insert_validated(projparams, lat_0, 'lat_0', -90, 90) + _insert_validated(projparams, lat_1, 'lat_1', -90, 90) + _insert_validated(projparams, lat_2, 'lat_2', -90, 90) + _insert_validated(projparams, lat_ts, 'lat_ts', -90, 90) + _insert_validated(projparams, lon_0, 'lon_0', -360, 720) + _insert_validated(projparams, lon_1, 'lon_1', -360, 720) + _insert_validated(projparams, lon_2, 'lon_2', -360, 720) if satellite_height is not None: projparams['h'] = satellite_height + using_corners = (None not in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]) if using_corners: - # make sure lat/lon limits are converted to floats. - self.llcrnrlon = float(llcrnrlon) - self.llcrnrlat = float(llcrnrlat) - self.urcrnrlon = float(urcrnrlon) - self.urcrnrlat = float(urcrnrlat) - # check values of urcrnrlon,urcrnrlat and llcrnrlon,llcrnrlat - if self.urcrnrlat > 90.0 or self.llcrnrlat > 90.0: - raise ValueError, 'urcrnrlat and llcrnrlat must be less than 90' - if self.urcrnrlat < -90.0 or self.llcrnrlat < -90.0: - raise ValueError, 'urcrnrlat and llcrnrlat must be greater than -90' - if self.urcrnrlon > 720. or self.llcrnrlon > 720.: - raise ValueError, 'urcrnrlon and llcrnrlon must be less than 720' - if self.urcrnrlon < -360. or self.llcrnrlon < -360.: - raise ValueError, 'urcrnrlon and llcrnrlon must be greater than -360' - # for each of the supported projections, compute lat/lon of domain corners + self.llcrnrlon = _validated_ll(llcrnrlon, 'llcrnrlon', -360, 720) + self.urcrnrlon = _validated_ll(urcrnrlon, 'urcrnrlon', -360, 720) + self.llcrnrlat = _validated_ll(llcrnrlat, 'llcrnrlat', -90, 90) + self.urcrnrlat = _validated_ll(urcrnrlat, 'urcrnrlat', -90, 90) + + # for each of the supported projections, + # compute lat/lon of domain corners # and set values in projparams dict as needed. - if projection == 'lcc': + if projection in ['lcc', 'eqdc', 'aea']: # if lat_0 is given, but not lat_1, # set lat_1=lat_0 if lat_1 is None and lat_0 is not None: lat_1 = lat_0 projparams['lat_1'] = lat_1 if lat_1 is None or lon_0 is None: - raise ValueError, 'must specify lat_1 or lat_0 and lon_0 for Lambert Conformal basemap (lat_2 is optional)' + raise ValueError('must specify lat_1 or lat_0 and lon_0 for %(projection)s basemap (lat_2 is optional)' % projnames) if lat_2 is None: projparams['lat_2'] = lat_1 if not using_corners: if width is None or height is None: raise ValueError, 'must either specify lat/lon values of corners (llcrnrlon,llcrnrlat,ucrnrlon,urcrnrlat) in degrees or width and height in meters' - else: - if lon_0 is None or lat_0 is None: - raise ValueError, 'must specify lon_0 and lat_0 when using width, height to specify projection region' - llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat = _choosecorners(width,height,**projparams) - self.llcrnrlon = llcrnrlon; self.llcrnrlat = llcrnrlat - self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat - - elif projection == 'eqdc': - # if lat_0 is given, but not lat_1, - # set lat_1=lat_0 - if lat_1 is None and lat_0 is not None: - lat_1 = lat_0 - projparams['lat_1'] = lat_1 - if lat_1 is None or lon_0 is None: - raise ValueError, 'must specify lat_1 or lat_0 and lon_0 for Equidistant Conic basemap (lat_2 is optional)' - if lat_2 is None: - projparams['lat_2'] = lat_1 - if not using_corners: - if width is None or height is None: - raise ValueError, 'must either specify lat/lon values of corners (llcrnrlon,llcrnrlat,ucrnrlon,urcrnrlat) in degrees or width and height in meters' - else: - if lon_0 is None or lat_0 is None: - raise ValueError, 'must specify lon_0 and lat_0 when using width, height to specify projection region' - llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat = _choosecorners(width,height,**projparams) - self.llcrnrlon = llcrnrlon; self.llcrnrlat = llcrnrlat - self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat - elif projection == 'aea': - # if lat_0 is given, but not lat_1, - # set lat_1=lat_0 - if lat_1 is None and lat_0 is not None: - lat_1 = lat_0 - projparams['lat_1'] = lat_1 - if lat_1 is None or lon_0 is None: - raise ValueError, 'must specify lat_1 or lat_0 and lon_0 for Albers Equal Area basemap (lat_2 is optional)' - if lat_2 is None: - projparams['lat_2'] = lat_1 - if not using_corners: - if width is None or height is None: - raise ValueError, 'must either specify lat/lon values of corners (llcrnrlon,llcrnrlat,ucrnrlon,urcrnrlat) in degrees or width and height in meters' - else: - if lon_0 is None or lat_0 is None: - raise ValueError, 'must specify lon_0 and lat_0 when using width, height to specify projection region' - llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat = _choosecorners(width,height,**projparams) - self.llcrnrlon = llcrnrlon; self.llcrnrlat = llcrnrlat - self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat + if lon_0 is None or lat_0 is None: + raise ValueError, 'must specify lon_0 and lat_0 when using width, height to specify projection region' + llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat = _choosecorners(width,height,**projparams) + self.llcrnrlon = llcrnrlon; self.llcrnrlat = llcrnrlat + self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat + + # skipping over the following for now; it can be beautified and + # consolidated later elif projection == 'stere': if lat_0 is None or lon_0 is None: raise ValueError, 'must specify lat_0 and lon_0 for Stereographic basemap (lat_ts is optional)' if not using_corners: if width is None or height is None: raise ValueError, 'must either specify lat/lon values of corners (llcrnrlon,llcrnrlat,ucrnrlon,urcrnrlat) in degrees or width and height in meters' - else: - if lon_0 is None or lat_0 is None: - raise ValueError, 'must specify lon_0 and lat_0 when using width, height to specify projection region' - llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat = _choosecorners(width,height,**projparams) - self.llcrnrlon = llcrnrlon; self.llcrnrlat = llcrnrlat - self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat - elif projection == 'spstere': + if lon_0 is None or lat_0 is None: + raise ValueError, 'must specify lon_0 and lat_0 when using width, height to specify projection region' + llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat = _choosecorners(width,height,**projparams) + self.llcrnrlon = llcrnrlon; self.llcrnrlat = llcrnrlat + self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat + elif projection in ['spstere', 'npstere', + 'splaea', 'nplaea', + 'spaeqd', 'npaeqd']: if boundinglat is None or lon_0 is None: - raise ValueError, 'must specify boundinglat and lon_0 for South-Polar Stereographic basemap' - projparams['lat_ts'] = -90. - projparams['lat_0'] = -90. - projparams['proj'] = 'stere' - self.llcrnrlon = lon_0+45. - self.urcrnrlon = lon_0-135. + raise ValueError('must specify boundinglat and lon_0 for %(projection) basemap' % projnames) + if projection[0] == 's': + sgn = -1 + else: + sgn = 1 + rootproj = projection[2:] + projparams['proj'] = rootproj + if rootproj == 'stere': + projparams['lat_ts'] = sgn * 90. + projparams['lat_0'] = sgn * 90. + self.llcrnrlon = lon_0 - sgn*45. + self.urcrnrlon = lon_0 + sgn*135. proj = pyproj.Proj(projparams) x,y = proj(lon_0,boundinglat) lon,self.llcrnrlat = proj(math.sqrt(2.)*y,0.,inverse=True) self.urcrnrlat = self.llcrnrlat if width is not None or height is not None: print 'warning: width and height keywords ignored for %s projection' % self.projection - elif projection == 'npstere': - if boundinglat is None or lon_0 is None: - raise ValueError, 'must specify boundinglat and lon_0 for North-Polar Stereographic basemap' - projparams['lat_ts'] = 90. - projparams['lat_0'] = 90. - projparams['proj'] = 'stere' - self.llcrnrlon = lon_0-45. - self.urcrnrlon = lon_0+135. - proj = pyproj.Proj(projparams) - x,y = proj(lon_0,boundinglat) - lon,self.llcrnrlat = proj(math.sqrt(2.)*y,0.,inverse=True) - self.urcrnrlat = self.llcrnrlat - if width is not None or height is not None: - print 'warning: width and height keywords ignored for %s projection' % self.projection - elif projection == 'splaea': - if boundinglat is None or lon_0 is None: - raise ValueError, 'must specify boundinglat and lon_0 for South-Polar Lambert Azimuthal basemap' - projparams['lat_0'] = -90. - projparams['proj'] = 'laea' - self.llcrnrlon = lon_0+45. - self.urcrnrlon = lon_0-135. - proj = pyproj.Proj(projparams) - x,y = proj(lon_0,boundinglat) - lon,self.llcrnrlat = proj(math.sqrt(2.)*y,0.,inverse=True) - self.urcrnrlat = self.llcrnrlat - if width is not None or height is not None: - print 'warning: width and height keywords ignored for %s projection' % self.projection - elif projection == 'nplaea': - if boundinglat is None or lon_0 is None: - raise ValueError, 'must specify boundinglat and lon_0 for North-Polar Lambert Azimuthal basemap' - projparams['lat_0'] = 90. - projparams['proj'] = 'laea' - self.llcrnrlon = lon_0-45. - self.urcrnrlon = lon_0+135. - proj = pyproj.Proj(projparams) - x,y = proj(lon_0,boundinglat) - lon,self.llcrnrlat = proj(math.sqrt(2.)*y,0.,inverse=True) - self.urcrnrlat = self.llcrnrlat - if width is not None or height is not None: - print 'warning: width and height keywords ignored for %s projection' % self.projection - elif projection == 'spaeqd': - if boundinglat is None or lon_0 is None: - raise ValueError, 'must specify boundinglat and lon_0 for South-Polar Azimuthal Equidistant basemap' - projparams['lat_0'] = -90. - projparams['proj'] = 'aeqd' - self.llcrnrlon = lon_0+45. - self.urcrnrlon = lon_0-135. - proj = pyproj.Proj(projparams) - x,y = proj(lon_0,boundinglat) - lon,self.llcrnrlat = proj(math.sqrt(2.)*y,0.,inverse=True) - self.urcrnrlat = self.llcrnrlat - if width is not None or height is not None: - print 'warning: width and height keywords ignored for %s projection' % self.projection - elif projection == 'npaeqd': - if boundinglat is None or lon_0 is None: - raise ValueError, 'must specify boundinglat and lon_0 for North-Polar Azimuthal Equidistant basemap' - projparams['lat_0'] = 90. - projparams['proj'] = 'aeqd' - self.llcrnrlon = lon_0-45. - self.urcrnrlon = lon_0+135. - proj = pyproj.Proj(projparams) - x,y = proj(lon_0,boundinglat) - lon,self.llcrnrlat = proj(math.sqrt(2.)*y,0.,inverse=True) - self.urcrnrlat = self.llcrnrlat - if width is not None or height is not None: - print 'warning: width and height keywords ignored for %s projection' % self.projection elif projection == 'laea': if lat_0 is None or lon_0 is None: raise ValueError, 'must specify lat_0 and lon_0 for Lambert Azimuthal basemap' if not using_corners: if width is None or height is None: raise ValueError, 'must either specify lat/lon values of corners (llcrnrlon,llcrnrlat,ucrnrlon,urcrnrlat) in degrees or width and height in meters' - else: - if lon_0 is None or lat_0 is None: - raise ValueError, 'must specify lon_0 and lat_0 when using width, height to specify projection region' - llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat = _choosecorners(width,height,**projparams) - self.llcrnrlon = llcrnrlon; self.llcrnrlat = llcrnrlat - self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat + if lon_0 is None or lat_0 is None: + raise ValueError, 'must specify lon_0 and lat_0 when using width, height to specify projection region' + llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat = _choosecorners(width,height,**projparams) + self.llcrnrlon = llcrnrlon; self.llcrnrlat = llcrnrlat + self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat elif projection == 'merc': if lat_ts is None: raise ValueError, 'must specify lat_ts for Mercator basemap' @@ -479,13 +416,12 @@ if not using_corners: if width is None or height is None: raise ValueError, 'must either specify lat/lon values of corners (llcrnrlon,llcrnrlat,ucrnrlon,urcrnrlat) in degrees or width and height in meters' - else: - if lon_0 is None or lat_0 is None: - raise ValueError, 'must specify lon_0 and lat_0 when using width, height to specify projection region' - llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat = _choosecorners(width,height,**projparams) - self.llcrnrlon = llcrnrlon; self.llcrnrlat = llcrnrlat - self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat - + if lon_0 is None or lat_0 is None: + raise ValueError, 'must specify lon_0 and lat_0 when using width, height to specify projection region' + llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat = _choosecorners(width,height,**projparams) + self.llcrnrlon = llcrnrlon; self.llcrnrlat = llcrnrlat + self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat + elif projection == 'ortho': if not projparams.has_key('R'): raise ValueError, 'orthographic projection only works for perfect spheres - not ellipsoids' @@ -547,12 +483,11 @@ if not using_corners: if width is None or height is None: raise ValueError, 'must either specify lat/lon values of corners (llcrnrlon,llcrnrlat,ucrnrlon,urcrnrlat) in degrees or width and height in meters' - else: - if lon_0 is None or lat_0 is None: - raise ValueError, 'must specify lon_0 and lat_0 when using width, height to specify projection region' - llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat = _choosecorners(width,height,**projparams) - self.llcrnrlon = llcrnrlon; self.llcrnrlat = llcrnrlat - self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat + if lon_0 is None or lat_0 is None: + raise ValueError, 'must specify lon_0 and lat_0 when using width, height to specify projection region' + llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat = _choosecorners(width,height,**projparams) + self.llcrnrlon = llcrnrlon; self.llcrnrlat = llcrnrlat + self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat elif projection == 'mill': if not using_corners: llcrnrlon = -180. @@ -574,22 +509,7 @@ if width is not None or height is not None: print 'warning: width and height keywords ignored for %s projection' % self.projection else: - raise ValueError, """ - unsupported projection, use 'cyl' - cylindrical equidistant, 'merc' - - mercator, 'lcc' - lambert conformal conic, 'stere' - stereographic, - 'npstere' - stereographic, special case centered on north pole. - 'spstere' - stereographic, special case centered on south pole, - 'aea' - albers equal area conic, 'tmerc' - transverse mercator, - 'aeqd' - azimuthal equidistant, 'mill' - miller cylindrical, - 'npaeqd' - azimuthal equidistant, special case centered on north pole, - 'spaeqd' - azimuthal equidistant, special case centered on south pole, - 'eqdc' - equidistant conic, 'laea' - lambert azimuthal equal area, - 'nplaea' - lambert azimuthal, special case centered on north pole, - 'splaea' - lambert azimuthal, special case centered on south pole, - 'cass' - cassini-soldner (transverse cylindrical equidistant), - 'poly' - polyconic, 'omerc' - oblique mercator, 'ortho' - orthographic, - 'geos' - geostationary, 'sinu' - sinusoidal, 'moll' - mollweide, - 'robin' - robinson, or 'gnom' - gnomonic. You tried '%s'""" % projection + raise ValueError(_unsupported_projection % projection) # initialize proj4 proj = Proj(projparams,self.llcrnrlon,self.llcrnrlat,self.urcrnrlon,self.urcrnrlat) @@ -701,34 +621,37 @@ def __call__(self,x,y,inverse=False): """ - Calling a Basemap class instance with the arguments lon, lat will - convert lon/lat (in degrees) to x/y native map projection - coordinates (in meters). If optional keyword 'inverse' is - True (default is False), the inverse transformation from x/y - to lon/lat is performed. + Calling a Basemap class instance with the arguments lon, lat will + convert lon/lat (in degrees) to x/y native map projection + coordinates (in meters). If optional keyword 'inverse' is + True (default is False), the inverse transformation from x/y + to lon/lat is performed. - For cylindrical equidistant projection ('cyl'), this - does nothing (i.e. x,y == lon,lat). + For cylindrical equidistant projection ('cyl'), this + does nothing (i.e. x,y == lon,lat). - For non-cylindrical projections, the inverse transformation - always returns longitudes between -180 and 180 degrees. For - cylindrical projections (self.projection == 'cyl','mill' or 'merc') - the inverse transformation will return longitudes between - self.llcrnrlon and self.llcrnrlat. + For non-cylindrical projections, the inverse transformation + always returns longitudes between -180 and 180 degrees. For + cylindrical projections (self.projection == 'cyl','mill' or 'merc') + the inverse transformation will return longitudes between + self.llcrnrlon and self.llcrnrlat. - input arguments lon, lat can be either scalar floats or N arrays. + input arguments lon, lat can be either scalar floats or N arrays. """ return self.projtran(x,y,inverse=inverse) def makegrid(self,nx,ny,returnxy=False): """ - return arrays of shape (ny,nx) containing lon,lat coordinates of - an equally spaced native projection grid. - if returnxy = True, the x,y values of the grid are returned also. + return arrays of shape (ny,nx) containing lon,lat coordinates of + an equally spaced native projection grid. + if returnxy = True, the x,y values of the grid are returned also. """ return self.projtran.makegrid(nx,ny,returnxy=returnxy) def _readboundarydata(self,name): + """ + read boundary data, clip to map projection region + """ msg = """ Unable to open boundary dataset file. Only the 'crude', 'low', 'intermediate' and 'high' resolution datasets are installed by default. If you @@ -925,7 +848,7 @@ def _getmapboundary(self): """ - create map boundary polygon (in lat/lon and x/y coordinates) + create map boundary polygon (in lat/lon and x/y coordinates) """ nx = 100 ny = 100 @@ -1043,8 +966,9 @@ def drawmapboundary(self,color='k',linewidth=1.0,ax=None): """ - draw boundary around map projection region. If ax=None (default), - default axis instance is used, otherwise specified axis instance is used. + draw boundary around map projection region. If ax=None (default), + default axis instance is used, otherwise specified axis + instance is used. """ # get current axes instance (if none specified). if ax is None and self.ax is None: @@ -1105,15 +1029,15 @@ def fillcontinents(self,color='0.8',ax=None,zorder=None): """ - Fill continents. + Fill continents. - color - color to fill continents (default gray). - ax - axes instance (overrides default axes instance). - zorder - sets the zorder for the continent polygons (if not specified, - uses default zorder for a Polygon patch). Set to zero if you want to paint - over the filled continents). + color - color to fill continents (default gray). + ax - axes instance (overrides default axes instance). + zorder - sets the zorder for the continent polygons (if not specified, + uses default zorder for a Polygon patch). Set to zero if you want to paint + over the filled continents). - After filling continents, lakes are re-filled with axis background color. + After filling continents, lakes are re-filled with axis background color. """ if self.resolution is None: raise AttributeError, 'there are no boundary datasets associated with this Basemap instance' @@ -1161,14 +1085,14 @@ def drawcoastlines(self,linewidth=1.,color='k',antialiased=1,ax=None,zorder=None): """ - Draw coastlines. + Draw coastlines. - linewidth - coastline width (default 1.) - color - coastline color (default black) - antialiased - antialiasing switch for coastlines (default True). - ax - axes instance (overrides default axes instance) - zorder - sets the zorder for the coastlines (if not specified, - uses default zorder for LineCollections). + linewidth - coastline width (default 1.) + color - coastline color (default black) + antialiased - antialiasing switch for coastlines (default True). + ax - axes instance (overrides default axes instance) + zorder - sets the zorder for the coastlines (if not specified, + uses default zorder for LineCollections). """ if self.resolution is None: raise AttributeError, 'there are no boundary datasets associated with this Basemap instance' @@ -1192,14 +1116,14 @@ def drawcountries(self,linewidth=0.5,color='k',antialiased=1,ax=None,zorder=None): """ - Draw country boundaries. + Draw country boundaries. - linewidth - country boundary line width (default 0.5) - color - country boundary line color (default black) - antialiased - antialiasing switch for country boundaries (default True). - ax - axes instance (overrides default axes instance) - zorder - sets the zorder for the country boundaries (if not specified, - uses default zorder for LineCollections). + linewidth - country boundary line width (default 0.5) + color - country boundary line color (default black) + antialiased - antialiasing switch for country boundaries (default True). + ax - axes instance (overrides default axes instance) + zorder - sets the zorder for the country boundaries (if not specified, + uses default zorder for LineCollections). """ if self.resolution is None: raise AttributeError, 'there are no boundary datasets associated with this Basemap instance' @@ -1227,14 +1151,14 @@ def drawstates(self,linewidth=0.5,color='k',antialiased=1,ax=None,zorder=None): """ - Draw state boundaries in Americas. + Draw state boundaries in Americas. - linewidth - state boundary line width (default 0.5) - color - state boundary line color (default black) - antialiased - antialiasing switch for state boundaries (default True). - ax - axes instance (overrides default axes instance) - zorder - sets the zorder for the state boundaries (if not specified, - uses default zorder for LineCollections). + linewidth - state boundary line width (default 0.5) + color - state boundary line color (default black) + antialiased - antialiasing switch for state boundaries (default True). + ax - axes instance (overrides default axes instance) + zorder - sets the zorder for the state boundaries (if not specified, + uses default zorder for LineCollections). """ if self.resolution is None: raise AttributeError, 'there are no boundary datasets associated with this Basemap instance' @@ -1262,14 +1186,14 @@ def drawrivers(self,linewidth=0.5,color='k',antialiased=1,ax=None,zorder=None): """ - Draw major rivers. + Draw major rivers. - linewidth - river boundary line width (default 0.5) - color - river boundary line color (default black) - antialiased - antialiasing switch for river boundaries (default True). - ax - axes instance (overrides default axes instance) - zorder - sets the zorder for the rivers (if not specified, - uses default zorder for LineCollections). + linewidth - river boundary line width (default 0.5) + color - river boundary line color (default black) + antialiased - antialiasing switch for river boundaries (default True). + ax - axes instance (overrides default axes instance) + zorder - sets the zorder for the rivers (if not specified, + uses default zorder for LineCollections). """ if self.resolution is None: raise AttributeError, 'there are no boundary datasets associated with this Basemap instance' @@ -1298,40 +1222,40 @@ def readshapefile(self,shapefile,name,drawbounds=True,zorder=None, linewidth=0.5,color='k',antialiased=1,ax=None): """ - read in shape file, draw boundaries on map. + read in shape file, draw boundaries on map. - Restrictions: - - Assumes shapes are 2D - - vertices must be in geographic (lat/lon) coordinates. + Restrictions: + - Assumes shapes are 2D + - vertices must be in geographic (lat/lon) coordinates. - shapefile - path to shapefile components. Example: - shapefile='/home/jeff/esri/world_borders' assumes that - world_borders.shp, world_borders.shx and world_borders.dbf - live in /home/jeff/esri. - name - name for Basemap attribute to hold the shapefile - vertices in native map projection coordinates. - Class attribute name+'_info' is a list of dictionaries, one - for each shape, containing attributes of each shape from dbf file. - For example, if name='counties', self.counties - will be a list of vertices for each shape in map projection - coordinates and self.counties_info will be a list of dictionaries - with shape attributes. Rings in individual shapes are split out - into separate polygons. Additional keys - 'RINGNUM' and 'SHAPENUM' are added to shape attribute dictionary. - drawbounds - draw boundaries of shapes (default True) - zorder = shape boundary zorder (if not specified, default for LineCollection - is used). - linewidth - shape boundary line width (default 0.5) - color - shape boundary line color (default black) - antialiased - antialiasing switch for shape boundaries (default True). - ax - axes instance (overrides default axes instance) + shapefile - path to shapefile components. Example: + shapefile='/home/jeff/esri/world_borders' assumes that + world_borders.shp, world_borders.shx and world_borders.dbf + live in /home/jeff/esri. + name - name for Basemap attribute to hold the shapefile + vertices in native map projection coordinates. + Class attribute name+'_info' is a list of dictionaries, one + for each shape, containing attributes of each shape from dbf file. + For example, if name='counties', self.counties + will be a list of vertices for each shape in map projection + coordinates and self.counties_info will be a list of dictionaries + with shape attributes. Rings in individual shapes are split out + into separate polygons. Additional keys + 'RINGNUM' and 'SHAPENUM' are added to shape attribute dictionary. + drawbounds - draw boundaries of shapes (default True) + zorder = shape boundary zorder (if not specified, default for LineCollection + is used). + linewidth - shape boundary line width (default 0.5) + color - shape boundary line color (default black) + antialiased - antialiasing switch for shape boundaries (default True). + ax - axes instance (overrides default axes instance) - returns a tuple (num_shapes, type, min, max) containing shape file info. - num_shapes is the number of shapes, type is the type code (one of - the SHPT* constants defined in the shapelib module, see - http://shapelib.maptools.org/shp_api.html) and min and - max are 4-element lists with the minimum and maximum values of the - vertices. + returns a tuple (num_shapes, type, min, max) containing shape file info. + num_shapes is the number of shapes, type is the type code (one of + the SHPT* constants defined in the shapelib module, see + http://shapelib.maptools.org/shp_api.html) and min and + max are 4-element lists with the minimum and maximum values of the + vertices. """ # open shapefile, read vertices for each object, convert # to map projection coordinates (only works for 2D shape types). @@ -1355,11 +1279,11 @@ for ring in range(rings): lons, lats = zip(*verts[ring]) if max(lons) > 721. or min(lons) < -721. or max(lats) > 91. or min(lats) < -91: - msg=""" -shapefile must have lat/lon vertices - it looks like this one has vertices -in map projection coordinates. You can convert the shapefile to geographic -coordinates using the shpproj utility from the shapelib tools -(http://shapelib.maptools.org/shapelib-tools.html)""" + msg=dedent(""" + shapefile must have lat/lon vertices - it looks like this one has vertices + in map projection coordinates. You can convert the shapefile to geographic + coordinates using the shpproj utility from the shapelib tools + (http://shapelib.maptools.org/shapelib-tools.html)""") raise ValueError,msg x, y = self(lons, lats) shpsegs.append(zip(x,y)) @@ -1400,32 +1324,32 @@ linestyle='--',dashes=[1,1],labels=[0,0,0,0],labelstyle=None, \ fmt='%g',xoffset=None,yoffset=None,ax=None,**kwargs): """ - draw parallels (latitude lines). + draw parallels (latitude lines). - circles - list containing latitude values to draw (in degrees). - color - color to draw parallels (default black). - linewidth - line width for parallels (default 1.) - zorder - sets the zorder for parallels (if not specified, - uses default zorder for Line2D class). - linestyle - line style for parallels (default '--', i.e. dashed). - dashes - dash pattern for parallels (default [1,1], i.e. 1 pixel on, - 1 pixel off). - labels - list of 4 values (default [0,0,0,0]) that control whether - parallels are labelled where they intersect the left, right, top or - bottom of the plot. For example labels=[1,0,0,1] will cause parallels - to be labelled where they intersect the left and bottom of the plot, - but not the right and top. - labelstyle - if set to "+/-", north and south latitudes are labelled - with "+" and "-", otherwise they are labelled with "N" and "S". - fmt is a format string to format the parallel labels (default '%g'). - xoffset - label offset from edge of map in x-direction - (default is 0.01 times width of map in map projection coordinates). - yoffset - label offset from edge of map in y-direction - (default is 0.01 times height of map in map projection coordinates). - ax - axes instance (overrides default axes instance) + circles - list containing latitude values to draw (in degrees). + color - color to draw parallels (default black). + linewidth - line width for parallels (default 1.) + zorder - sets the zorder for parallels (if not specified, + uses default zorder for Line2D class). + linestyle - line style for parallels (default '--', i.e. dashed). + dashes - dash pattern for parallels (default [1,1], i.e. 1 pixel on, + 1 pixel off). + labels - list of 4 values (default [0,0,0,0]) that control whether + parallels are labelled where they intersect the left, right, top or + bottom of the plot. For example labels=[1,0,0,1] will cause parallels + to be labelled where they intersect the left and bottom of the plot, + but not the right and top. + labelstyle - if set to "+/-", north and south latitudes are labelled + with "+" and "-", otherwise they are labelled with "N" and "S". + fmt is a format string to format the parallel labels (default '%g'). + xoffset - label offset from edge of map in x-direction + (default is 0.01 times width of map in map projection coordinates). + yoffset - label offset from edge of map in y-direction + (default is 0.01 times height of map in map projection coordinates). + ax - axes instance (overrides default axes instance) - additional keyword arguments control text properties for labels (see - pylab.text documentation) + additional keyword arguments control text properties for labels (see + pylab.text documentation) """ # get current axes instance (if none specified). if ax is None and self.ax is None: @@ -1620,32 +1544,32 @@ linestyle='--',dashes=[1,1],labels=[0,0,0,0],labelstyle=None,\ fmt='%g',xoffset=None,yoffset=None,ax=None,**kwargs): """ - draw meridians (longitude lines). + draw meridians (longitude lines). - meridians - list containing longitude values to draw (in degrees). - color - color to draw meridians (default black). - linewidth - line width for meridians (default 1.) - zorder - sets the zorder for meridians (if not specified, - uses default zorder for Line2D class). - linestyle - line style for meridians (default '--', i.e. dashed). - dashes - dash pattern for meridians (default [1,1], i.e. 1 pixel on, - 1 pixel off). - labels - list of 4 values (default [0,0,0,0]) that control whether - meridians are labelled where they intersect the left, right, top or - bottom of the plot. For example labels=[1,0,0,1] will cause meridians - to be labelled where they intersect the left and bottom of the plot, - but not the right and top. - labelstyle - if set to "+/-", east and west longitudes are labelled - with "+" and "-", otherwise they are labelled with "E" and "W". - fmt is a format string to format the meridian labels (default '%g'). - xoffset - label offset from edge of map in x-direction - (default is 0.01 times width of map in map projection coordinates). - yoffset - label offset from edge of map in y-direction - (default is 0.01 times height of map in map projection coordinates). - ax - axes instance (overrides default axes instance) + meridians - list containing longitude values to draw (in degrees). + color - color to draw meridians (default black). + linewidth - line width for meridians (default 1.) + zorder - sets the zorder for meridians (if not specified, + uses default zorder for Line2D class). + linestyle - line style for meridians (default '--', i.e. dashed). + dashes - dash pattern for meridians (default [1,1], i.e. 1 pixel on, + 1 pixel off). + labels - list of 4 values (default [0,0,0,0]) that control whether + meridians are labelled where they intersect the left, right, top or + bottom of the plot. For example labels=[1,0,0,1] will cause meridians + to be labelled where they intersect the left and bottom of the plot, + but not the right and top. + labelstyle - if set to "+/-", east and west longitudes are labelled + with "+" and "-", otherwise they are labelled with "E" and "W". + fmt is a format string to format the meridian labels (default '%g'). + xoffset - label offset from edge of map in x-direction + (default is 0.01 times width of map in map projection coordinates). + yoffset - label offset from edge of map in y-direction + (default is 0.01 times height of map in map projection coordinates). + ax - axes instance (overrides default axes instance) - additional keyword arguments control text properties for labels (see - pylab.text documentation) + additional keyword arguments control text properties for labels (see + pylab.text documentation) """ # get current axes instance (if none specified). if ax is None and self.ax is None: @@ -1824,9 +1748,9 @@ def gcpoints(self,lon1,lat1,lon2,lat2,npoints): """ - compute npoints points along a great circle with endpoints - (lon1,lat1) and (lon2,lat2). Returns arrays x,y - with map projection coordinates. + compute npoints points along a great circle with endpoints + (lon1,lat1) and (lon2,lat2). Returns arrays x,y + with map projection coordinates. """ gc = pyproj.Geod(a=self.rmajor,b=self.rminor) lonlats = gc.npts(lon1,lat1,lon2,lat2,npoints-2) @@ -1839,17 +1763,17 @@ def drawgreatcircle(self,lon1,lat1,lon2,lat2,del_s=100.,**kwargs): """ - draw a great circle on the map. + draw a great circle on the map. - lon1,lat1 - longitude,latitude of one endpoint of the great circle. - lon2,lat2 - longitude,latitude of the other endpoint of the great circle. - del_s - points on great circle computed every delta kilometers (default 100). + lon1,lat1 - longitude,latitude of one endpoint of the great circle. + lon2,lat2 - longitude,latitude of the other endpoint of the great circle. + del_s - points on great circle computed every delta kilometers (default 100). - Other keyword arguments (**kwargs) control plotting of great circle line, - see pylab.plot documentation for details. + Other keyword arguments (**kwargs) control plotting of great circle line, + see pylab.plot documentation for details. - Note: cannot handle situations in which the great circle intersects - the edge of the map projection domain, and then re-enters the domain. + Note: cannot handle situations in which the great circle intersects + the edge of the map projection domain, and then re-enters the domain. """ # use great circle formula for a perfect sphere. gc = pyproj.Geod(a=self.rmajor,b=self.rminor) @@ -1866,29 +1790,29 @@ def transform_scalar(self,datin,lons,lats,nx,ny,returnxy=False,checkbounds=False,order=1,masked=False): """ - interpolate a scalar field (datin) from a lat/lon grid with longitudes = - lons and latitudes = lats to a (ny,nx) native map projection grid. - Typically used to transform data to map projection coordinates - so it can be plotted on the map with imshow. + interpolate a scalar field (datin) from a lat/lon grid with longitudes = + lons and latitudes = lats to a (ny,nx) native map projection grid. + Typically used to transform data to map projection coordinates + so it can be plotted on the map with imshow. - lons, lats must be rank-1 arrays containing longitudes and latitudes - (in degrees) of datin grid in increasing order - (i.e. from dateline eastward, and South Pole northward). - For non-cylindrical projections (those other than - cylindrical equidistant, mercator and miller) - lons must fit within range -180 to 180. + lons, lats must be rank-1 arrays containing longitudes and latitudes + (in degrees) of datin grid in increasing order + (i.e. from dateline eastward, and South Pole northward). + For non-cylindrical projections (those other than + cylindrical equidistant, mercator and miller) + lons must fit within range -180 to 180. - if returnxy=True, the x and y values of the native map projection grid - are also returned. + if returnxy=True, the x and y values of the native map projection grid + are also returned. - If checkbounds=True, values of lons and lats are checked to see that - they lie within the map projection region. Default is False. - If checkbounds=False, points outside map projection region will - be clipped to values on the boundary if masked=False. If masked=True, - the return value will be a masked array with those points masked. + If checkbounds=True, values of lons and lats are checked to see that + they lie within the map projection region. Default is False. + If checkbounds=False, points outside map projection region will + be clipped to values on the boundary if masked=False. If masked=True, + the return value will be a masked array with those points masked. - The order keyword can be 0 for nearest-neighbor interpolation, - or 1 for bilinear interpolation (default 1). + The order keyword can be 0 for nearest-neighbor interpolation, + or 1 for bilinear interpolation (default 1). """ # check that lons, lats increasing delon = lons[1:]-lons[0:-1] @@ -1916,33 +1840,33 @@ def transform_vector(self,uin,vin,lons,lats,nx,ny,returnxy=False,checkbounds=False,order=1,masked=False):... [truncated message content] |
From: <js...@us...> - 2007-11-18 16:16:17
|
Revision: 4372 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4372&view=rev Author: jswhit Date: 2007-11-18 08:16:15 -0800 (Sun, 18 Nov 2007) Log Message: ----------- add 'f' to coastline resolutions. Reformat unsupported projection message. Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-18 14:31:06 UTC (rev 4371) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-18 16:16:15 UTC (rev 4372) @@ -72,7 +72,7 @@ but if they are not, the entire globe is plotted. resolution - resolution of boundary database to use. Can be 'c' (crude), - 'l' (low), 'i' (intermediate), 'h' (high), or None. Default is 'c'. + 'l' (low), 'i' (intermediate), 'h' (high), 'f' (full) or None. If None, no boundary data will be read in (and class methods such as drawcoastlines will raise an exception if invoked). Resolution drops off by roughly 80% @@ -83,8 +83,8 @@ Tools (http://gmt.soest.hawaii.edu). area_thresh - coastline or lake with an area smaller than area_thresh - in km^2 will not be plotted. Default 10000,1000,100,10 for resolution - 'c','l','i','h'. + in km^2 will not be plotted. Default 10000,1000,100,10,1 for resolution + 'c','l','i','h','f'. rsphere - radius of the sphere used to define map projection (default 6370997 meters, close to the arithmetic mean radius of the earth). If @@ -142,28 +142,8 @@ latitude circle boundinglat is tangent to the edge of the map at lon_0. satellite_height - height of satellite (in m) above equator - only relevant for geostationary projections ('geos'). - - """ -_unsupported_projection = """ - unsupported projection, use 'cyl' - cylindrical equidistant, 'merc' - - mercator, 'lcc' - lambert conformal conic, 'stere' - stereographic, - 'npstere' - stereographic, special case centered on north pole. - 'spstere' - stereographic, special case centered on south pole, - 'aea' - albers equal area conic, 'tmerc' - transverse mercator, - 'aeqd' - azimuthal equidistant, 'mill' - miller cylindrical, - 'npaeqd' - azimuthal equidistant, special case centered on north pole, - 'spaeqd' - azimuthal equidistant, special case centered on south pole, - 'eqdc' - equidistant conic, 'laea' - lambert azimuthal equal area, - 'nplaea' - lambert azimuthal, special case centered on north pole, - 'splaea' - lambert azimuthal, special case centered on south pole, - 'cass' - cassini-soldner (transverse cylindrical equidistant), - 'poly' - polyconic, 'omerc' - oblique mercator, 'ortho' - orthographic, - 'geos' - geostationary, 'sinu' - sinusoidal, 'moll' - mollweide, - 'robin' - robinson, or 'gnom' - gnomonic. You tried '%s' - """ - # This allows substitution of longer names into error messages. projnames = {'cyl' : 'Cylindrical Equidistant', 'merc' : 'Mercator', @@ -180,7 +160,7 @@ 'spaeqd' : 'South-Polar Azimuthal Equidistant', 'aea' : 'Albers Equal Area', 'stere' : 'Stereographic', - 'npstere' : 'Nouth-Polar Stereographic', + 'npstere' : 'North-Polar Stereographic', 'spstere' : 'South-Polar Stereographic', 'cass' : 'Cassini-Soldner', 'poly' : 'Polyconic', @@ -192,6 +172,12 @@ 'gnom' : 'Gnomonic', } +_unsupported_projection = ["'%s' is an unsupported projection.\n"] +_unsupported_projection.append("The supported projections are:\n") +for k,v in projnames.iteritems(): + _unsupported_projection.append("'%s' = %s\n" % (k,v)) +_unsupported_projection = ''.join(_unsupported_projection) + def _validated_ll(param, name, minval, maxval): param = float(param) if param > maxval or param < minval: @@ -650,13 +636,14 @@ def _readboundarydata(self,name): """ - read boundary data, clip to map projection region + read boundary data, clip to map projection region. """ - msg = """ -Unable to open boundary dataset file. Only the 'crude', 'low', -'intermediate' and 'high' resolution datasets are installed by default. If you -are requesting a 'full' resolution dataset, you need to download -and install those files separately(see the basemap README for details).""" + msg = dedent(""" + Unable to open boundary dataset file. Only the 'crude', 'low', + 'intermediate' and 'high' resolution datasets are installed by default. + If you are requesting a 'full' resolution dataset, you may need to + download and install those files separately + (see the basemap README for details).""") try: bdatfile = open(os.path.join(basemap_datadir,name+'_'+self.resolution+'.dat'),'rb') bdatmetafile = open(os.path.join(basemap_datadir,name+'meta_'+self.resolution+'.dat'),'r') This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-18 16:45:29
|
Revision: 4373 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4373&view=rev Author: jswhit Date: 2007-11-18 08:45:28 -0800 (Sun, 18 Nov 2007) Log Message: ----------- further streamlining of docstrings Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-18 16:16:15 UTC (rev 4372) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-18 16:45:28 UTC (rev 4373) @@ -21,6 +21,38 @@ __version__ = '0.9.7' +# supported map projections. +projnames = {'cyl' : 'Cylindrical Equidistant', + 'merc' : 'Mercator', + 'tmerc' : 'Transverse Mercator', + 'omerc' : 'Oblique Mercator', + 'mill' : 'Miller Cylindrical', + 'lcc' : 'Lambert Conformal', + 'laea' : 'Lambert Azimuthal Equal Area', + 'nplaea' : 'North-Polar Lambert Azimuthal', + 'splaea' : 'South-Polar Lambert Azimuthal', + 'eqdc' : 'Equidistant Conic', + 'aeqd' : 'Azimuthal Equidistant', + 'npaeqd' : 'North-Polar Azimuthal Equidistant', + 'spaeqd' : 'South-Polar Azimuthal Equidistant', + 'aea' : 'Albers Equal Area', + 'stere' : 'Stereographic', + 'npstere' : 'North-Polar Stereographic', + 'spstere' : 'South-Polar Stereographic', + 'cass' : 'Cassini-Soldner', + 'poly' : 'Polyconic', + 'ortho' : 'Orthographic', + 'geos' : 'Geostationary', + 'sinu' : 'Sinusoidal', + 'moll' : 'Mollweide', + 'robin' : 'Robinson', + 'gnom' : 'Gnomonic', + } +supported_projections = [] +for k,v in projnames.iteritems(): + supported_projections.append("'%s' = %s\n" % (k,v)) +supported_projections = ''.join(supported_projections) + # The __init__ docstring is pulled out here because it is so long; # Having it in the usual place makes it hard to get from the # __init__ argument list to the code that uses the arguments. @@ -29,21 +61,8 @@ arguments: - projection - map projection. 'cyl' - cylindrical equidistant, 'merc' - - mercator, 'lcc' - lambert conformal conic, 'stere' - stereographic, - 'npstere' - stereographic, special case centered on north pole. - 'spstere' - stereographic, special case centered on south pole, - 'aea' - albers equal area conic, 'tmerc' - transverse mercator, - 'aeqd' - azimuthal equidistant, 'mill' - miller cylindrical, - 'npaeqd' - azimuthal equidistant, special case centered on north pole, - 'spaeqd' - azimuthal equidistant, special case centered on south pole, - 'eqdc' - equidistant conic, 'laea' - lambert azimuthal equal area, - 'nplaea' - lambert azimuthal, special case centered on north pole, - 'splaea' - lambert azimuthal, special case centered on south pole, - 'cass' - cassini-soldner (transverse cylindrical equidistant), - 'poly' - polyconic, 'omerc' - oblique mercator, 'ortho' - orthographic, - 'sinu' - sinusoidal, 'moll' - mollweide, 'robin' - robinson, - 'geos' - geostationary, and 'gnom' - gnomonic are currently available. + projection - map projection. Supported projections are:\n"""+\ +supported_projections+""" Default is 'cyl'. The map projection region can either be specified by setting these keywords: @@ -142,40 +161,12 @@ latitude circle boundinglat is tangent to the edge of the map at lon_0. satellite_height - height of satellite (in m) above equator - only relevant for geostationary projections ('geos'). - """ +""" -# This allows substitution of longer names into error messages. -projnames = {'cyl' : 'Cylindrical Equidistant', - 'merc' : 'Mercator', - 'tmerc' : 'Transverse Mercator', - 'omerc' : 'Oblique Mercator', - 'mill' : 'Miller Cylindrical', - 'llc' : 'Lambert Conformal', - 'laea' : 'Lambert Azimuthal Equal Area', - 'nplaea' : 'North-Polar Lambert Azimuthal', - 'splaea' : 'South-Polar Lambert Azimuthal', - 'eqdc' : 'Equidistant Conic', - 'eaqd' : 'Azimuthal Equidistant', - 'npaeqd' : 'North-Polar Azimuthal Equidistant', - 'spaeqd' : 'South-Polar Azimuthal Equidistant', - 'aea' : 'Albers Equal Area', - 'stere' : 'Stereographic', - 'npstere' : 'North-Polar Stereographic', - 'spstere' : 'South-Polar Stereographic', - 'cass' : 'Cassini-Soldner', - 'poly' : 'Polyconic', - 'ortho' : 'Orthographic', - 'geos' : 'Geostationary', - 'sinu' : 'Sinusoidal', - 'moll' : 'Mollweide', - 'robin' : 'Robinson', - 'gnom' : 'Gnomonic', - } - +# unsupported projection error message. _unsupported_projection = ["'%s' is an unsupported projection.\n"] _unsupported_projection.append("The supported projections are:\n") -for k,v in projnames.iteritems(): - _unsupported_projection.append("'%s' = %s\n" % (k,v)) +_unsupported_projection.append(supported_projections) _unsupported_projection = ''.join(_unsupported_projection) def _validated_ll(param, name, minval, maxval): @@ -191,12 +182,7 @@ class Basemap(object): """ - Set up a basemap with one of 19 supported map projections - (cylindrical equidistant, mercator, polyconic, oblique mercator, - transverse mercator, miller cylindrical, lambert conformal conic, - azimuthal equidistant, equidistant conic, lambert azimuthal equal area, - albers equal area conic, gnomonic, orthographic, sinusoidal, mollweide, - geostationary, robinson, cassini-soldner or stereographic). + Set up a basemap with specified map projection. Doesn't actually draw anything, but sets up the map projection class and creates the coastline, lake river and political boundary data structures in native map projection coordinates. @@ -204,10 +190,8 @@ Useful instance variables: - projection - map projection ('cyl','merc','mill','lcc','eqdc','aea', - 'laea', 'nplaea', 'splaea', 'tmerc', 'omerc', 'cass', 'gnom', 'poly', - 'sinu', 'moll', 'ortho', 'robin', 'aeqd', 'npaeqd', 'spaeqd', 'stere', - 'geos', 'npstere' or 'spstere') + projection - map projection. Print the module variable + "supported_projections" to see a list. (projections prefixed with 'np' or 'sp' are special case polar-centric versions of the parent projection) aspect - map aspect ratio (size of y dimension / size of x dimension). @@ -247,11 +231,8 @@ [this example (simpletest.py) plus many others can be found in the examples directory of source distribution. The "OO" version of this example (which does not use pylab) is called "simpletest_oo.py".] - - Contact: Jeff Whitaker <jef...@no...> """ - def __init__(self, llcrnrlon=None, llcrnrlat=None, urcrnrlon=None, urcrnrlat=None, width=None, height=None, @@ -267,7 +248,7 @@ anchor='C', ax=None): # docstring is added after definition - #print "starting: ", time.clock() + # where to put plot in figure (default is 'C' or center) self.anchor = anchor # map projection. @@ -604,6 +585,7 @@ else: coastsegs.append(seg) self.coastsegs = coastsegs + __init__.__doc__ = _Basemap_init_doc def __call__(self,x,y,inverse=False): """ This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-18 19:28:36
|
Revision: 4378 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4378&view=rev Author: jswhit Date: 2007-11-18 11:28:13 -0800 (Sun, 18 Nov 2007) Log Message: ----------- more docstring tweaks Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-18 19:06:49 UTC (rev 4377) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-18 19:28:13 UTC (rev 4378) @@ -22,7 +22,7 @@ __version__ = '0.9.7' # supported map projections. -projnames = {'cyl' : 'Cylindrical Equidistant', +_projnames = {'cyl' : 'Cylindrical Equidistant', 'merc' : 'Mercator', 'tmerc' : 'Transverse Mercator', 'omerc' : 'Oblique Mercator', @@ -49,7 +49,7 @@ 'gnom' : 'Gnomonic', } supported_projections = [] -for k,v in projnames.iteritems(): +for k,v in _projnames.iteritems(): supported_projections.append("'%s' = %s\n" % (k,v)) supported_projections = ''.join(supported_projections) @@ -191,9 +191,7 @@ Useful instance variables: projection - map projection. Print the module variable - "supported_projections" to see a list. - (projections prefixed with 'np' or 'sp' are special case polar-centric - versions of the parent projection) + "supported_projections" to see a list of supported projections. aspect - map aspect ratio (size of y dimension / size of x dimension). llcrnrlon - longitude of lower left hand corner of the desired map domain. llcrnrlon - latitude of lower left hand corner of the desired map domain. @@ -301,7 +299,7 @@ lat_1 = lat_0 projparams['lat_1'] = lat_1 if lat_1 is None or lon_0 is None: - raise ValueError('must specify lat_1 or lat_0 and lon_0 for %(projection)s basemap (lat_2 is optional)' % projnames) + raise ValueError('must specify lat_1 or lat_0 and lon_0 for %(projection)s basemap (lat_2 is optional)' % _projnames) if lat_2 is None: projparams['lat_2'] = lat_1 if not using_corners: @@ -330,7 +328,7 @@ 'splaea', 'nplaea', 'spaeqd', 'npaeqd']: if boundinglat is None or lon_0 is None: - raise ValueError('must specify boundinglat and lon_0 for %(projection) basemap' % projnames) + raise ValueError('must specify boundinglat and lon_0 for %(projection) basemap' % _projnames) if projection[0] == 's': sgn = -1 else: This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-19 13:10:52
|
Revision: 4380 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4380&view=rev Author: jswhit Date: 2007-11-19 05:10:49 -0800 (Mon, 19 Nov 2007) Log Message: ----------- more minor reformatting Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-18 20:54:22 UTC (rev 4379) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-19 13:10:49 UTC (rev 4380) @@ -1,10 +1,10 @@ from matplotlib import rcParams -from matplotlib import __version__ as matplotlib_version +from matplotlib import __version__ as _matplotlib_version # check to make sure matplotlib is not too old. -mpl_required_version = '0.90' -if matplotlib_version < mpl_required_version: +_mpl_required_version = '0.90' +if _matplotlib_version < _mpl_required_version: raise ImportError('your matplotlib is too old - basemap ' - 'requires at least %s, you have %s'%(mpl_required_version,matplotlib_version)) + 'requires version %s or higher'% _matplotlib_version) from matplotlib.collections import LineCollection from matplotlib.patches import Ellipse, Circle, Polygon from matplotlib.lines import Line2D @@ -49,8 +49,8 @@ 'gnom' : 'Gnomonic', } supported_projections = [] -for k,v in _projnames.iteritems(): - supported_projections.append("'%s' = %s\n" % (k,v)) +for _items in _projnames.iteritems(): + supported_projections.append("'%s' = %s\n" % (_items)) supported_projections = ''.join(supported_projections) # The __init__ docstring is pulled out here because it is so long; @@ -61,8 +61,8 @@ arguments: - projection - map projection. Supported projections are:\n"""+\ -supported_projections+""" + projection - map projection. Supported projections are: +%(supported_projections)s Default is 'cyl'. The map projection region can either be specified by setting these keywords: @@ -94,7 +94,7 @@ 'l' (low), 'i' (intermediate), 'h' (high), 'f' (full) or None. If None, no boundary data will be read in (and class methods such as drawcoastlines will raise an exception if invoked). - Resolution drops off by roughly 80% + Resolution drops off by roughly 80%% between datasets. Higher res datasets are much slower to draw. Default 'c'. Coastline data is from the GSHHS (http://www.soest.hawaii.edu/wessel/gshhs/gshhs.html). @@ -161,7 +161,7 @@ latitude circle boundinglat is tangent to the edge of the map at lon_0. satellite_height - height of satellite (in m) above equator - only relevant for geostationary projections ('geos'). -""" +""" % locals() # unsupported projection error message. _unsupported_projection = ["'%s' is an unsupported projection.\n"] @@ -191,7 +191,7 @@ Useful instance variables: projection - map projection. Print the module variable - "supported_projections" to see a list of supported projections. + "supported_projections" to see a list. aspect - map aspect ratio (size of y dimension / size of x dimension). llcrnrlon - longitude of lower left hand corner of the desired map domain. llcrnrlon - latitude of lower left hand corner of the desired map domain. @@ -245,7 +245,7 @@ boundinglat=None, anchor='C', ax=None): - # docstring is added after definition + # docstring is added after __init__ method definition # where to put plot in figure (default is 'C' or center) self.anchor = anchor @@ -281,7 +281,7 @@ _insert_validated(projparams, lon_2, 'lon_2', -360, 720) if satellite_height is not None: projparams['h'] = satellite_height - + # check for sane values of projection corners. using_corners = (None not in [llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat]) if using_corners: self.llcrnrlon = _validated_ll(llcrnrlon, 'llcrnrlon', -360, 720) @@ -292,6 +292,7 @@ # for each of the supported projections, # compute lat/lon of domain corners # and set values in projparams dict as needed. + if projection in ['lcc', 'eqdc', 'aea']: # if lat_0 is given, but not lat_1, # set lat_1=lat_0 @@ -310,9 +311,6 @@ llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat = _choosecorners(width,height,**projparams) self.llcrnrlon = llcrnrlon; self.llcrnrlat = llcrnrlat self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat - - # skipping over the following for now; it can be beautified and - # consolidated later elif projection == 'stere': if lat_0 is None or lon_0 is None: raise ValueError, 'must specify lat_0 and lon_0 for Stereographic basemap (lat_ts is optional)' @@ -386,7 +384,6 @@ llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat = _choosecorners(width,height,**projparams) self.llcrnrlon = llcrnrlon; self.llcrnrlat = llcrnrlat self.urcrnrlon = urcrnrlon; self.urcrnrlat = urcrnrlat - elif projection == 'ortho': if not projparams.has_key('R'): raise ValueError, 'orthographic projection only works for perfect spheres - not ellipsoids' @@ -583,6 +580,7 @@ else: coastsegs.append(seg) self.coastsegs = coastsegs + # set __init__'s docstring __init__.__doc__ = _Basemap_init_doc def __call__(self,x,y,inverse=False): @@ -648,7 +646,6 @@ if containsPole and\ self.projection in ['tmerc','cass','omerc','merc','mill','cyl','robin','moll','sinu','geos']: raise ValueError('%s projection cannot cross pole'%(self.projection)) - # make sure orthographic projection has containsPole=True # we will compute the intersections in stereographic # coordinates, then transform to orthographic. @@ -2141,9 +2138,6 @@ """ Make a pseudo-color plot over the map. see pylab.pcolormesh documentation for definition of **kwargs - Unlike pcolor, pcolormesh cannot handle masked arrays, and so - cannot be used to plot data when the grid lies partially outside - the projection limb (use pcolor or contourf instead). extra keyword 'ax' can be used to override the default axis instance. """ if not kwargs.has_key('ax') and self.ax is None: This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-19 13:14:55
|
Revision: 4381 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4381&view=rev Author: jswhit Date: 2007-11-19 05:14:53 -0800 (Mon, 19 Nov 2007) Log Message: ----------- run reindent.py Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-19 13:10:49 UTC (rev 4380) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-19 13:14:53 UTC (rev 4381) @@ -49,8 +49,8 @@ 'gnom' : 'Gnomonic', } supported_projections = [] -for _items in _projnames.iteritems(): - supported_projections.append("'%s' = %s\n" % (_items)) +for _items in _projnames.iteritems(): + supported_projections.append("'%s' = %s\n" % (_items)) supported_projections = ''.join(supported_projections) # The __init__ docstring is pulled out here because it is so long; @@ -91,7 +91,7 @@ but if they are not, the entire globe is plotted. resolution - resolution of boundary database to use. Can be 'c' (crude), - 'l' (low), 'i' (intermediate), 'h' (high), 'f' (full) or None. + 'l' (low), 'i' (intermediate), 'h' (high), 'f' (full) or None. If None, no boundary data will be read in (and class methods such as drawcoastlines will raise an exception if invoked). Resolution drops off by roughly 80%% @@ -289,8 +289,8 @@ self.llcrnrlat = _validated_ll(llcrnrlat, 'llcrnrlat', -90, 90) self.urcrnrlat = _validated_ll(urcrnrlat, 'urcrnrlat', -90, 90) - # for each of the supported projections, - # compute lat/lon of domain corners + # for each of the supported projections, + # compute lat/lon of domain corners # and set values in projparams dict as needed. if projection in ['lcc', 'eqdc', 'aea']: @@ -302,7 +302,7 @@ if lat_1 is None or lon_0 is None: raise ValueError('must specify lat_1 or lat_0 and lon_0 for %(projection)s basemap (lat_2 is optional)' % _projnames) if lat_2 is None: - projparams['lat_2'] = lat_1 + projparams['lat_2'] = lat_1 if not using_corners: if width is None or height is None: raise ValueError, 'must either specify lat/lon values of corners (llcrnrlon,llcrnrlat,ucrnrlon,urcrnrlat) in degrees or width and height in meters' @@ -475,7 +475,7 @@ # initialize proj4 proj = Proj(projparams,self.llcrnrlon,self.llcrnrlat,self.urcrnrlon,self.urcrnrlat) - + # make sure axis ticks are suppressed. self.noticks = suppress_ticks @@ -571,12 +571,12 @@ dist = npy.sqrt(xd+yd) split = dist > 5000000. if npy.sum(split) and self.projection not in ['merc','cyl','mill']: - ind = (npy.compress(split,squeeze(split*npy.indices(xd.shape)))+1).tolist() - iprev = 0 - ind.append(len(xd)) - for i in ind: - coastsegs.append(zip(x[iprev:i],y[iprev:i])) - iprev = i + ind = (npy.compress(split,squeeze(split*npy.indices(xd.shape)))+1).tolist() + iprev = 0 + ind.append(len(xd)) + for i in ind: + coastsegs.append(zip(x[iprev:i],y[iprev:i])) + iprev = i else: coastsegs.append(seg) self.coastsegs = coastsegs @@ -619,7 +619,7 @@ msg = dedent(""" Unable to open boundary dataset file. Only the 'crude', 'low', 'intermediate' and 'high' resolution datasets are installed by default. - If you are requesting a 'full' resolution dataset, you may need to + If you are requesting a 'full' resolution dataset, you may need to download and install those files separately (see the basemap README for details).""") try: @@ -649,7 +649,7 @@ # make sure orthographic projection has containsPole=True # we will compute the intersections in stereographic # coordinates, then transform to orthographic. - if self.projection == 'ortho' and name == 'gshhs': + if self.projection == 'ortho' and name == 'gshhs': containsPole = True lon_0=self.projparams['lon_0'] lat_0=self.projparams['lat_0'] @@ -674,7 +674,7 @@ offsetbytes = int(linesplit[5]) bytecount = int(linesplit[6]) bdatfile.seek(offsetbytes,0) - # read in binary string convert into an npts by 2 + # read in binary string convert into an npts by 2 # numpy array (first column is lons, second is lats). polystring = bdatfile.read(bytecount) # binary data is little endian. @@ -685,7 +685,7 @@ # coordinates (i.e. it does not contain either pole), # the intersections of the boundary geometries # and the map projection region can be computed before - # transforming the boundary geometry to map projection + # transforming the boundary geometry to map projection # coordinates (this saves time, especially for small map # regions and high-resolution boundary geometries). if not containsPole: @@ -713,7 +713,7 @@ poly = Shape(b) antart = False # create duplicate polygons shifted by -360 and +360 - # (so as to properly treat polygons that cross + # (so as to properly treat polygons that cross # Greenwich meridian). if not antart: b2[:,0] = b[:,0]-360 @@ -760,7 +760,7 @@ else: b[:,0], b[:,1] = self(b[:,0], b[:,1]) goodmask = npy.logical_and(b[:,0]<1.e20,b[:,1]<1.e20) - # if less than two points are valid in + # if less than two points are valid in # map proj coords, skip this geometry. if npy.sum(goodmask) <= 1: continue if name != 'gshhs': @@ -849,7 +849,7 @@ projparms['x_0']=-llcrnrx projparms['y_0']=-llcrnry maptran = pyproj.Proj(projparms) - elif self.projection in ['moll','robin','sinu']: + elif self.projection in ['moll','robin','sinu']: # quasi-elliptical region. lon_0 = self.projparams['lon_0'] # left side @@ -893,7 +893,7 @@ boundaryxy = _geos.Polygon(b) if self.projection in ['mill','merc','cyl']: # make sure map boundary doesn't quite include pole. - if self.urcrnrlat > 89.9999: + if self.urcrnrlat > 89.9999: urcrnrlat = 89.9999 else: urcrnrlat = self.urcrnrlat @@ -908,14 +908,14 @@ b[:,0]=x; b[:,1]=y boundaryxy = _geos.Polygon(b) else: - if self.projection not in ['moll','robin','sinu']: + if self.projection not in ['moll','robin','sinu']: lons, lats = maptran(x,y,inverse=True) # fix lons so there are no jumps. n = 1 lonprev = lons[0] for lon,lat in zip(lons[1:],lats[1:]): if npy.abs(lon-lonprev) > 90.: - if lonprev < 0: + if lonprev < 0: lon = lon - 360. else: lon = lon + 360 @@ -1004,7 +1004,7 @@ After filling continents, lakes are re-filled with axis background color. """ if self.resolution is None: - raise AttributeError, 'there are no boundary datasets associated with this Basemap instance' + raise AttributeError, 'there are no boundary datasets associated with this Basemap instance' # get current axes instance (if none specified). if ax is None and self.ax is None: try: @@ -1059,7 +1059,7 @@ uses default zorder for LineCollections). """ if self.resolution is None: - raise AttributeError, 'there are no boundary datasets associated with this Basemap instance' + raise AttributeError, 'there are no boundary datasets associated with this Basemap instance' # get current axes instance (if none specified). if ax is None and self.ax is None: try: @@ -1090,7 +1090,7 @@ uses default zorder for LineCollections). """ if self.resolution is None: - raise AttributeError, 'there are no boundary datasets associated with this Basemap instance' + raise AttributeError, 'there are no boundary datasets associated with this Basemap instance' # read in country line segments, only keeping those that # intersect map boundary polygon. if not hasattr(self,'cntrysegs'): @@ -1125,7 +1125,7 @@ uses default zorder for LineCollections). """ if self.resolution is None: - raise AttributeError, 'there are no boundary datasets associated with this Basemap instance' + raise AttributeError, 'there are no boundary datasets associated with this Basemap instance' # read in state line segments, only keeping those that # intersect map boundary polygon. if not hasattr(self,'statesegs'): @@ -1160,7 +1160,7 @@ uses default zorder for LineCollections). """ if self.resolution is None: - raise AttributeError, 'there are no boundary datasets associated with this Basemap instance' + raise AttributeError, 'there are no boundary datasets associated with this Basemap instance' # read in river line segments, only keeping those that # intersect map boundary polygon. if not hasattr(self,'riversegs'): @@ -1374,15 +1374,15 @@ dist = npy.sqrt(xd+yd) split = dist > 500000. if npy.sum(split) and self.projection not in ['merc','cyl','mill','moll','robin','sinu']: - ind = (npy.compress(split,squeeze(split*npy.indices(xd.shape)))+1).tolist() - xl = [] - yl = [] - iprev = 0 - ind.append(len(xd)) - for i in ind: - xl.append(x[iprev:i]) - yl.append(y[iprev:i]) - iprev = i + ind = (npy.compress(split,squeeze(split*npy.indices(xd.shape)))+1).tolist() + xl = [] + yl = [] + iprev = 0 + ind.append(len(xd)) + for i in ind: + xl.append(x[iprev:i]) + yl.append(y[iprev:i]) + iprev = i else: xl = [x] yl = [y] @@ -1580,15 +1580,15 @@ dist = npy.sqrt(xd+yd) split = dist > 500000. if npy.sum(split) and self.projection not in ['merc','cyl','mill','moll','robin','sinu']: - ind = (npy.compress(split,squeeze(split*npy.indices(xd.shape)))+1).tolist() - xl = [] - yl = [] - iprev = 0 - ind.append(len(xd)) - for i in ind: - xl.append(x[iprev:i]) - yl.append(y[iprev:i]) - iprev = i + ind = (npy.compress(split,squeeze(split*npy.indices(xd.shape)))+1).tolist() + xl = [] + yl = [] + iprev = 0 + ind.append(len(xd)) + for i in ind: + xl.append(x[iprev:i]) + yl.append(y[iprev:i]) + iprev = i else: xl = [x] yl = [y] @@ -2189,7 +2189,7 @@ ax = self.ax else: ax = kwargs.pop('ax') - # make sure x is monotonically increasing - if not, + # make sure x is monotonically increasing - if not, # print warning suggesting that the data be shifted in longitude # with the shiftgrid function. # only do this check for global projections. @@ -2233,7 +2233,7 @@ # reset current active image (only if pylab is imported). try: try: # new contour. - if CS._A is not None: pylab.gci._current = CS + if CS._A is not None: pylab.gci._current = CS except: # old contour. if CS[1].mappable is not None: pylab.gci._current = CS[1].mappable except: @@ -2257,7 +2257,7 @@ ax = self.ax else: ax = kwargs.pop('ax') - # make sure x is monotonically increasing - if not, + # make sure x is monotonically increasing - if not, # print warning suggesting that the data be shifted in longitude # with the shiftgrid function. # only do this check for global projections. @@ -2559,7 +2559,7 @@ if masked: xmask = npy.logical_or(npy.less(xcoords,0),npy.greater(xcoords,len(xin)-1)) ymask = npy.logical_or(npy.less(ycoords,0),npy.greater(ycoords,len(yin)-1)) - xymask = npy.logical_or(xmask,ymask) + xymask = npy.logical_or(xmask,ymask) xcoords = npy.clip(xcoords,0,len(xin)-1) ycoords = npy.clip(ycoords,0,len(yin)-1) # interpolate to output grid using bilinear interpolation. @@ -2657,6 +2657,6 @@ corners = llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat # test for invalid projection points on output if llcrnrlon > 1.e20 or urcrnrlon > 1.e20: - raise ValueError, 'width and/or height too large for this projection, try smaller values' + raise ValueError, 'width and/or height too large for this projection, try smaller values' else: - return corners + return corners This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-20 02:37:13
|
Revision: 4384 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4384&view=rev Author: jswhit Date: 2007-11-19 18:37:12 -0800 (Mon, 19 Nov 2007) Log Message: ----------- lake_color keyword arg to fill continents added Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-19 21:43:24 UTC (rev 4383) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-20 02:37:12 UTC (rev 4384) @@ -991,11 +991,12 @@ # set axes limits to fit map region. self.set_axes_limits(ax=ax) - def fillcontinents(self,color='0.8',ax=None,zorder=None): + def fillcontinents(self,color='0.8',lake_color=None,ax=None,zorder=None): """ Fill continents. color - color to fill continents (default gray). + lake_color - color to fill inland lakes (default axes background). ax - axes instance (overrides default axes instance). zorder - sets the zorder for the continent polygons (if not specified, uses default zorder for a Polygon patch). Set to zero if you want to paint @@ -1038,8 +1039,11 @@ xy = zip(xa.tolist(),ya.tolist()) if self.coastpolygontypes[np] != 2: poly = Polygon(xy,facecolor=color,edgecolor=color,linewidth=0) - else: # lakes filled with background color. - poly = Polygon(xy,facecolor=axisbgc,edgecolor=axisbgc,linewidth=0) + else: # lakes filled with background color by default + if lake_color is None: + poly = Polygon(xy,facecolor=axisbgc,edgecolor=axisbgc,linewidth=0) + else: + poly = Polygon(xy,facecolor=lake_color,edgecolor=lake_color,linewidth=0) if zorder is not None: poly.set_zorder(zorder) ax.add_patch(poly) This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-20 13:18:23
|
Revision: 4389 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4389&view=rev Author: jswhit Date: 2007-11-20 05:18:22 -0800 (Tue, 20 Nov 2007) Log Message: ----------- added module variable 'projection_params' - a dict containing the relevant projection parameters for each projection. Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-20 13:14:34 UTC (rev 4388) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-20 13:18:22 UTC (rev 4389) @@ -53,6 +53,34 @@ supported_projections.append("'%s' = %s\n" % (_items)) supported_projections = ''.join(supported_projections) +# projection specific parameters. +projection_params = {'cyl' : 'corners only (no width/height)', + 'merc' : 'corners plus lat_ts (no width/height)', + 'tmerc' : 'lon_0,lat_0', + 'omerc' : 'lon_0,lat_0,lat_1,lat_2,lon_1,lon_2,no width/height', + 'mill' : 'corners only (no width/height)', + 'lcc' : 'lon_0,lat_0,lat_1,lat_2', + 'laea' : 'lon_0,lat_0', + 'nplaea' : 'bounding_lat,lon_0,lat_0,no corners or width/height', + 'splaea' : 'bounding_lat,lon_0,lat_0,no corners or width/height', + 'eqdc' : 'lon_0,lat_0,lat_1,lat_2', + 'aeqd' : 'lon_0,lat_0', + 'npaeqd' : 'bounding_lat,lon_0,lat_0,no corners or width/height', + 'spaeqd' : 'bounding_lat,lon_0,lat_0,no corners or width/height', + 'aea' : 'lon_0,lat_0,lat_1', + 'stere' : 'lon_0,lat_0,lat_ts', + 'npstere' : 'bounding_lat,lon_0,lat_0,no corners or width/height', + 'spstere' : 'bounding_lat,lon_0,lat_0,no corners or width/height', + 'cass' : 'lon_0,lat_0', + 'poly' : 'lon_0,lat_0', + 'ortho' : 'lon_0,lat_0', + 'geos' : 'lon_0,lat_0,satellite_height', + 'sinu' : 'lon_0,lat_0,no corners or width/height', + 'moll' : 'lon_0,lat_0,no corners or width/height', + 'robin' : 'lon_0,lat_0,no corners or width/height', + 'gnom' : 'lon_0,lat_0', + } + # The __init__ docstring is pulled out here because it is so long; # Having it in the usual place makes it hard to get from the # __init__ argument list to the code that uses the arguments. @@ -65,7 +93,8 @@ %(supported_projections)s Default is 'cyl'. - The map projection region can either be specified by setting these keywords: + For most map projections, the map projection region can either be + specified by setting these keywords: llcrnrlon - longitude of lower left hand corner of the desired map domain (degrees). llcrnrlat - latitude of lower left hand corner of the desired map domain (degrees). @@ -81,9 +110,10 @@ For 'sinu', 'moll', 'npstere', 'spstere', 'nplaea', 'splaea', 'nplaea', 'splaea', 'npaeqd', 'spaeqd' or 'robin', the values of - llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat,width and height are ignored (because - either they are computed internally, or entire globe is always plotted). For the - cylindrical projections ('cyl','merc' and 'mill'), the default is to use + llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat,width and height are ignored + (because either they are computed internally, or entire globe is + always plotted). For the cylindrical projections + ('cyl','merc' and 'mill'), the default is to use llcrnrlon=-180,llcrnrlat=-90, urcrnrlon=180 and urcrnrlat=90). For all other projections except 'ortho' and 'geos', either the lat/lon values of the corners or width and height must be specified by the user. @@ -117,9 +147,10 @@ in map projection coordinates. Default False, so parallels and meridians can be labelled instead. If parallel or meridian labelling is requested (using drawparallels and drawmeridians methods), automatic tick labelling - will be supressed even is suppress_ticks=False. Typically, you will - only want to override the default if you want to label the axes in meters - using native map projection coordinates. + will be supressed even is suppress_ticks=False. suppress_ticks=False + is useful if you want to use your own custom tick formatter, or + if you want to let matplotlib label the axes in meters + using native map projection coordinates anchor - determines how map is placed in axes rectangle (passed to axes.set_aspect). Default is 'C', which means map is centered. @@ -136,9 +167,11 @@ The following parameters are map projection parameters which all default to None. Not all parameters are used by all projections, some are ignored. + The module variable 'projection_params' lists which parameters apply + to which projections. - lat_ts - latitude of natural origin (used for mercator, and - optionally for stereographic projection). + lat_ts - latitude of true scale for mercator projection, optional + for stereographic projection. lat_1 - first standard parallel for lambert conformal, albers equal area projection and equidistant conic projections. Latitude of one of the two points on the projection centerline for oblique mercator. @@ -154,6 +187,7 @@ lon_2 - longitude of one of the two points on the projection centerline for oblique mercator. lat_0 - central latitude (y-axis origin) - used by all projections, + Must be equator for mercator projection. lon_0 - central meridian (x-axis origin) - used by all projections, boundinglat - bounding latitude for pole-centered projections (npstere,spstere, nplaea,splaea,npaeqd,spaeqd). These projections are square regions centered This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <js...@us...> - 2007-11-20 13:20:57
|
Revision: 4390 http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4390&view=rev Author: jswhit Date: 2007-11-20 05:20:54 -0800 (Tue, 20 Nov 2007) Log Message: ----------- docstring tweak. Modified Paths: -------------- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py Modified: trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py =================================================================== --- trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-20 13:18:22 UTC (rev 4389) +++ trunk/toolkits/basemap-testing/lib/matplotlib/toolkits/basemap/basemap.py 2007-11-20 13:20:54 UTC (rev 4390) @@ -167,8 +167,8 @@ The following parameters are map projection parameters which all default to None. Not all parameters are used by all projections, some are ignored. - The module variable 'projection_params' lists which parameters apply - to which projections. + The module variable 'projection_params' is a dictionary which + lists which parameters apply to which projections. lat_ts - latitude of true scale for mercator projection, optional for stereographic projection. This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |