shadedrelief basemap question

Hi Jeff,

I'm continuing to enjoy using basemap, but have a question about the shaded relief background. I frequently use greyscale shaded relief on the continents, but blue or white for the oceans. the shadedrelief() function is really convenient, but it includes shading for the oceans. Is there a way to afterwards just shade oceans? I've included my modified shaded_relief function in this email that uses some of the other natural earth products in case you need to see it.

Example:
bmap.drawmapboundary(fill_color='aqua')
bmap.shadedrelief()

Or something like:
bmap.shadedrelief()
**bmap.filloceans('aqua')

Thanks,
Scott

Scott: You could overlay the shaded relief image on a land-sea mask, where the land part of the mask is transparent. Like this:

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np
bmap = Basemap(projection='ortho',lat_0=45,lon_0=-100,resolution='l')
bmap.drawmapboundary()
bmap.drawmeridians(np.arange(0,360,30))
bmap.drawparallels(np.arange(-90,90,30))
bmap.shadedrelief()
bmap.drawlsmask(ocean_color='aqua',land_color=(255,255,255,1))
plt.show()

PS. is there a basemap mailing list I should send questions like this to?

No - just send questions to matplotlib-users (I've cc'ed the list)

-Jeff

···

On 8/1/12 4:09 PM, Scott Henderson wrote:

from PIL import Image
from matplotlib.image import pil_to_array

     def shadedrelief_highres(self,style='hypsometric', region=(-75.0, -35.0, -65.0, -15.0),ax=None,scale=None,**kwargs):
         ''' can avoid Memory Error with warpimage() by pre-cutting high-res arrays with GDAL?
             region = (minLon, minLat, maxLon, maxLat)'''
         print 'NOTE: only set up for cyl coordinates with extent defined by region'
         basedir = '/Users/scotthenderson/data/natural_earth/'
         if style == 'hypsometric':
             path = os.path.join(basedir,'HYP_HR_SR_OB_DR.tif')
         elif style == 'natural':
             path = os.path.join(basedir,'NE1_HR_LC_SR_W_DR.tif')
         elif style == 'hillshade':
             path = os.path.join(basedir,'SR_HR.tif')

         outfile = os.path.join(basedir,'tmp.tif')
         if os.path.isfile(outfile): os.remove(outfile)
         os.system('gdalwarp -te {0} {1} {2} {3} {infile} {out}'.format(*region, infile=path, out=outfile))

         #copy only relevant commands from warpimage imshow() overrides colorbar for rgb or PIL arrays
         pilImage = Image.open(outfile)
         if scale is not None:
             w, h = pilImage.size
             width = int(np.round(w*scale))
             height = int(np.round(h*scale))
             pilImage = pilImage.resize((width,height),Image.ANTIALIAS)
         self._bm_rgba = pil_to_array(pilImage)
         # if pil_to_array returns a 2D array, it's a grayscale image.
         # create an RGB image, with R==G==B.
         if self._bm_rgba.ndim == 2:
             tmp = np.empty(self._bm_rgba.shape+(3,),np.uint8)
             for k in range(3):
                 tmp[:,:,k] = self._bm_rgba
             self._bm_rgba = tmp

         im = self.imshow(self._bm_rgba, ax=ax, **kwargs)
         return im

--
Jeffrey S. Whitaker Phone : (303)497-6313
Meteorologist FAX : (303)497-6449
NOAA/OAR/PSD R/PSD1 Email : Jeffrey.S.Whitaker@...259...
325 Broadway Office : Skaggs Research Cntr 1D-113
Boulder, CO, USA 80303-3328 Web : http://tinyurl.com/5telg