SF.net SVN: matplotlib: [4968] trunk/toolkits/basemap

[<js...@us...>]

Revision: 4968
          http://matplotlib.svn.sourceforge.net/matplotlib/?rev=4968&view=rev
Author:   jswhit
Date:     2008-02-15 05:16:54 -0800 (Fri, 15 Feb 2008)

Log Message:
-----------
add warpimage method.

Modified Paths:
--------------
    trunk/toolkits/basemap/MANIFEST.in
    trunk/toolkits/basemap/examples/README
    trunk/toolkits/basemap/examples/warpimage.py
    trunk/toolkits/basemap/lib/mpl_toolkits/basemap/basemap.py
    trunk/toolkits/basemap/setup.py

Added Paths:
-----------
    trunk/toolkits/basemap/examples/earth_lights_lrg.jpg
    trunk/toolkits/basemap/lib/mpl_toolkits/basemap/data/bmng.jpg

Removed Paths:
-------------
    trunk/toolkits/basemap/examples/bmng.jpg

Modified: trunk/toolkits/basemap/MANIFEST.in
===================================================================
--- trunk/toolkits/basemap/MANIFEST.in	2008-02-14 23:18:09 UTC (rev 4967)
+++ trunk/toolkits/basemap/MANIFEST.in	2008-02-15 13:16:54 UTC (rev 4968)
@@ -49,7 +49,7 @@
 include examples/run_all.py
 include examples/polarmaps.py
 include examples/warpimage.py
-include examples/bmng.jpg
+include examples/earth_lights_lrg.jpg
 include examples/pnganim.py
 include examples/garp.py
 include examples/setwh.py

Modified: trunk/toolkits/basemap/examples/README
===================================================================
--- trunk/toolkits/basemap/examples/README	2008-02-14 23:18:09 UTC (rev 4967)
+++ trunk/toolkits/basemap/examples/README	2008-02-15 13:16:54 UTC (rev 4968)
@@ -81,8 +81,8 @@
 projections (prefixed by 'np' and 'sp').
 
 warpimage.py shows how to interpolate (warp) an image from one
-map projection to another. Requires PIL, and an image from
-http://www.space-graphics.com/earth_topo-bathy.htm, 
+map projection to another using the 'warpimage' Basemap method.
+Requires PIL.
 
 garp.py makes a 'World According to Garp' map - an azimuthal equidistant
 projection centered on a specified location.  Straight lines from that

Deleted: trunk/toolkits/basemap/examples/bmng.jpg
===================================================================
(Binary files differ)

Added: trunk/toolkits/basemap/examples/earth_lights_lrg.jpg
===================================================================
(Binary files differ)


Property changes on: trunk/toolkits/basemap/examples/earth_lights_lrg.jpg
___________________________________________________________________
Name: svn:mime-type
   + application/octet-stream

Modified: trunk/toolkits/basemap/examples/warpimage.py
===================================================================
--- trunk/toolkits/basemap/examples/warpimage.py	2008-02-14 23:18:09 UTC (rev 4967)
+++ trunk/toolkits/basemap/examples/warpimage.py	2008-02-15 13:16:54 UTC (rev 4968)
@@ -1,53 +1,35 @@
 import pylab as P
 import numpy
-from mpl_toolkits.basemap import Basemap as Basemap1 
-from numpy import ma
+from mpl_toolkits.basemap import Basemap
 
-class Basemap(Basemap1):
-    # subclass Basemap and add bluemarble method.
-    def bluemarble(self,masked=False):
-        """display 'blue marble next generation' image from http://visibleearth.nasa.gov/"""
-        try:
-            from PIL import Image
-        except ImportError:
-            raise ImportError('bluemarble method requires PIL (http://www.pythonware.com/products/pil/)')
-        from matplotlib.image import pil_to_array
-        
-        # read in jpeg image to rgba array of normalized floats.
-        pilImage = Image.open('bmng.jpg')
-        rgba = pil_to_array(pilImage)
-        rgba = rgba.astype(numpy.float32)/255. # convert to normalized floats.
-        
-        # define lat/lon grid that image spans (projection='cyl').
-        nlons = rgba.shape[1]; nlats = rgba.shape[0]
-        delta = 360./float(nlons)
-        lons = numpy.arange(-180.+0.5*delta,180.,delta)
-        lats = numpy.arange(-90.+0.5*delta,90.,delta)
+# illustrate use of warpimage method to display an image background
+# on the map projection region.  Default background is the 'blue
+# marble' image from NASA (http://visibleearth.nasa.gov).
 
-        if self.projection != 'cyl':
-            # transform to nx x ny regularly spaced native projection grid
-            # nx and ny chosen to have roughly the same horizontal res as original image.
-            dx = 2.*numpy.pi*m.rmajor/float(nlons)
-            nx = int((self.xmax-self.xmin)/dx)+1; ny = int((self.ymax-self.ymin)/dx)+1
-            rgba_warped = ma.zeros((ny,nx,4),numpy.float64)
-            # interpolate rgba values from proj='cyl' (geographic coords) to 'lcc'
-            # if masked=True, values outside of projection limb will be masked.
-            for k in range(4):
-                rgba_warped[:,:,k] = self.transform_scalar(rgba[:,:,k],lons,lats,nx,ny,masked=masked)
-            # make points outside projection limb transparent.
-            rgba_warped = rgba_warped.filled(0.)
-            # plot warped rgba image.
-            im = self.imshow(rgba_warped)
-        else:
-            im = self.imshow(rgba)
-        return im
+# create new figure
+fig=P.figure()
+# define orthographic projection centered on North America.
+m = Basemap(projection='ortho',lat_0=40,lon_0=-100,resolution='l')
+# display a non-default image.
+m.warpimage(file='earth_lights_lrg.jpg')
+# draw coastlines.
+m.drawcoastlines(linewidth=0.5,color='0.5')
+# draw lat/lon grid lines every 30 degrees.
+m.drawmeridians(numpy.arange(0,360,30),color='0.5')
+m.drawparallels(numpy.arange(-90,90,30),color='0.5')
+P.title("Lights at Night image warped from 'cyl' to 'ortho' projection",fontsize=12)
+print 'warp to orthographic map ...'
 
+# redisplay should be fast.
+fig = P.figure()
+m.warpimage()
+
 # create new figure
 fig=P.figure()
 # define cylindrical equidistant projection.
 m = Basemap(projection='cyl',llcrnrlon=-180,llcrnrlat=-90,urcrnrlon=180,urcrnrlat=90,resolution='l')
 # plot (unwarped) rgba image.
-im = m.bluemarble()
+im = m.warpimage()
 # draw coastlines.
 m.drawcoastlines(linewidth=0.5,color='0.5')
 # draw lat/lon grid lines.
@@ -58,10 +40,10 @@
 
 # create new figure
 fig=P.figure()
-# define orthographic projection centered on North America.
+# define orthographic projection centered on Europe.
 m = Basemap(projection='ortho',lat_0=40,lon_0=40,resolution='l')
 # plot warped rgba image.
-im = m.bluemarble(masked=True)
+im = m.warpimage(masked=True)
 # draw coastlines.
 m.drawcoastlines(linewidth=0.5,color='0.5')
 # draw lat/lon grid lines every 30 degrees.
@@ -76,7 +58,7 @@
 m = Basemap(llcrnrlon=-145.5,llcrnrlat=1.,urcrnrlon=-2.566,urcrnrlat=46.352,\
             rsphere=(6378137.00,6356752.3142),lat_1=50.,lon_0=-107.,\
             resolution='i',area_thresh=1000.,projection='lcc')
-im = m.bluemarble()
+im = m.warpimage()
 # draw coastlines.
 m.drawcoastlines(linewidth=0.5,color='0.5')
 # draw parallels and meridians.
@@ -95,7 +77,7 @@
             resolution=None,projection='omerc',\
             lon_0=-100,lat_0=15,lon_2=-120,lat_2=65,lon_1=-50,lat_1=-55)
 # plot warped rgba image.
-im = m.bluemarble()
+im = m.warpimage()
 # draw lat/lon grid lines every 20 degrees.
 m.drawmeridians(numpy.arange(0,360,20),color='0.5')
 m.drawparallels(numpy.arange(-80,81,20),color='0.5')

Modified: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/basemap.py
===================================================================
--- trunk/toolkits/basemap/lib/mpl_toolkits/basemap/basemap.py	2008-02-14 23:18:09 UTC (rev 4967)
+++ trunk/toolkits/basemap/lib/mpl_toolkits/basemap/basemap.py	2008-02-15 13:16:54 UTC (rev 4968)
@@ -2674,6 +2674,84 @@
         im = self.imshow(rgba,interpolation='nearest',ax=ax,**kwargs)
         return im
 
+    def warpimage(self,file=None,**kwargs):
+        """
+        display an image (given by file keyword) as a background.
+        Default (if file not specified) is to display 
+        'blue marble next generation' image from http://visibleearth.nasa.gov/.
+        Specified image must have pixels covering the whole globe in a regular
+        lat/lon grid, starting and -180W and the South Pole.
+        extra keyword 'ax' can be used to override the default axis instance.
+        """
+        try:
+            from PIL import Image
+        except ImportError:
+            raise ImportError('bluemarble method requires PIL (http://www.pythonware.com/products/pil/)')
+        from matplotlib.image import pil_to_array
+        if not kwargs.has_key('ax') and self.ax is None:
+            try:
+                ax = pylab.gca()
+            except:
+                import pylab
+                ax = pylab.gca()
+        elif not kwargs.has_key('ax') and self.ax is not None:
+            ax = self.ax
+        else:
+            ax = kwargs.pop('ax')
+        # default image file is blue marble next generation
+        # from NASA (http://visibleearth.nasa.gov).
+        if file is None:
+            file = os.path.join(basemap_datadir,'bmng.jpg')
+            newfile = False
+        else:
+            newfile = True
+        # read in jpeg image to rgba array of normalized floats.
+        if not hasattr(self,'bm_rgba') or newfile:
+            pilImage = Image.open(file)
+            self.bm_rgba = pil_to_array(pilImage)
+            # convert to normalized floats.
+            self.bm_rgba = self.bm_rgba.astype(npy.float32)/255.
+            # define lat/lon grid that image spans.
+            nlons = self.bm_rgba.shape[1]; nlats = self.bm_rgba.shape[0]
+            delta = 360./float(nlons)
+            self.bm_lons = npy.arange(-180.+0.5*delta,180.,delta)
+            self.bm_lats = npy.arange(-90.+0.5*delta,90.,delta)
+
+        if self.projection != 'cyl':
+            if newfile or not hasattr(self,'bm_rgba_warped'):
+                # transform to nx x ny regularly spaced native
+                # projection grid.
+                # nx and ny chosen to have roughly the 
+                # same horizontal res as original image.
+                dx = 2.*npy.pi*self.rmajor/float(nlons)
+                nx = int((self.xmax-self.xmin)/dx)+1
+                ny = int((self.ymax-self.ymin)/dx)+1
+                self.bm_rgba_warped = npy.zeros((ny,nx,4),npy.float64)
+                # interpolate rgba values from geographic coords (proj='cyl')
+                # to map projection coords.
+                # if masked=True, values outside of
+                # projection limb will be masked.
+                for k in range(4):
+                    self.bm_rgba_warped[:,:,k],x,y = \
+                    self.transform_scalar(self.bm_rgba[:,:,k],\
+                    self.bm_lons,self.bm_lats,nx,ny,returnxy=True)
+                # for ortho,geos mask pixels outside projection limb.
+                if self.projection in ['geos','ortho']:
+                    lonsr,latsr = self(x,y,inverse=True)
+                    mask = ma.zeros((nx,ny,4),npy.int8)
+                    mask[:,:,0] = npy.logical_or(lonsr>1.e20,latsr>1.e30)
+                    for k in range(1,4):
+                        mask[:,:,k] = mask[:,:,0]
+                    self.bm_rgba_warped = \
+                    ma.masked_array(self.bm_rgba_warped,mask=mask)
+                    # make points outside projection limb transparent.
+                    self.bm_rgba_warped = self.bm_rgba_warped.filled(0.)
+            # plot warped rgba image.
+            im = self.imshow(self.bm_rgba_warped,ax=ax)
+        else:
+            im = self.imshow(self.bm_rgba,ax=ax)
+        return im
+
 ### End of Basemap class
 
 def _searchlist(a,x):

Added: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/data/bmng.jpg
===================================================================
(Binary files differ)


Property changes on: trunk/toolkits/basemap/lib/mpl_toolkits/basemap/data/bmng.jpg
___________________________________________________________________
Name: svn:mime-type
   + application/octet-stream

Modified: trunk/toolkits/basemap/setup.py
===================================================================
--- trunk/toolkits/basemap/setup.py	2008-02-14 23:18:09 UTC (rev 4967)
+++ trunk/toolkits/basemap/setup.py	2008-02-15 13:16:54 UTC (rev 4968)
@@ -137,7 +137,7 @@
     package_dirs['httlib2'] = os.path.join('lib','httplib2')
 
 # Specify all the required mpl data
-pyproj_datafiles = ['data/epsg', 'data/esri', 'data/esri.extra', 'data/GL27', 'data/nad.lst', 'data/nad27', 'data/nad83', 'data/ntv2_out.dist', 'data/other.extra', 'data/pj_out27.dist', 'data/pj_out83.dist', 'data/proj_def.dat', 'data/README', 'data/td_out.dist', 'data/test27', 'data/test83', 'data/testntv2', 'data/testvarious', 'data/world']
+pyproj_datafiles = ['data/epsg', 'data/esri', 'data/esri.extra', 'data/GL27', 'data/nad.lst', 'data/nad27', 'data/nad83', 'data/ntv2_out.dist', 'data/other.extra', 'data/pj_out27.dist', 'data/pj_out83.dist', 'data/proj_def.dat', 'data/README', 'data/td_out.dist', 'data/test27', 'data/test83', 'data/testntv2', 'data/testvarious', 'data/world','data/bmng.jpg']
 boundaryfiles = []
 for resolution in ['c','l','i','h','f']:
     boundaryfiles = boundaryfiles + glob.glob("lib/mpl_toolkits/basemap/data/*_"+resolution+".dat")


This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.