From: John [H2O] <washakie@gm...>  20091007 23:10:21

Jeff, Here's a quick snippet. I've looked at the test.py file provided with the basemap examples. What I am unclear on are the different ways in which nx and ny are defined. I would like to have this 'automatically' defined, based solely on variables from my input object.. say for example a netcdf file that has len and lon dimensions defined. Below is my crude stab at it, but I am clearly having some problems. I guess the point is, maybe it's not possible to have a Basemap instance with extents beyond the imshow object. Then perhaps I need to make sure that when I set up the Basemap instance, I pass the H.outlon0 to llcrnrlon for example. But is that necessary? Thanks! #!/usr/bin/env python import matplotlib.pyplot as plt from mpl_toolkits.basemap import Basemap import numpy as np def plot_imshow_custom(H,transform=True ): """ function to automagically plot an mxn array of arbitrary lats/lons """ data = H.data print data.shape m = Basemap(projection='npstere',boundinglat=10,lon_0=270,resolution='l') fig = plt.figure() ax = fig.gca() print "Preparing to plot %s with dimensions:" % H.name print "lon0, numx, dx:" print H.outlon0, H.numxgrid, H.dxout print "lat0, numy, dy:" print H.outlat0, H.numygrid, H.dyout ## set up transformations for the data array ## THIS IS WHERE I NEED SOME HELP: if m.projection not in ['cyl','merc','mill']: lats = np.arange( H.outlat0, ( H.outlat0 + ( H.numygrid*H.dyout ) ), H.dyout )[:1] lons = np.arange( H.outlon0, ( H.outlon0 + ( H.numxgrid*H.dxout ) ), H.dxout )[:1] data = data[:1,:1] else: lats = np.arange( H.outlat0, ( H.outlat0 + ( H.numygrid*H.dyout ) ), H.dyout ) lons = np.arange( H.outlon0, ( H.outlon0 + ( H.numxgrid*H.dxout ) ), H.dxout ) print data.shape ## transform to nx x ny regularly spaced native projection grid if transform: if m.projection not in ['cyl','merc','mill']: dx = 2.*np.pi*m.rmajor/len(lons) dy = 2.*np.pi*m.rminor/len(lats) else: dx = len(lons) dy = len(lats) nx = int((m.xmaxm.xmin)/dx)+1; ny = int((m.ymaxm.ymin)/dy)+1 print nx if nx is 1: topodat = data else: topodat = m.transform_scalar(data,lons,lats,nx,ny) else: topodat = data ## Get the current axes, and properties for use later pos = ax.get_position() l, b, w, h = pos.bounds ## Set up the IMAGE colmap = plt.get_cmap('gist_ncar') im = m.imshow(topodat,cmap=colmap) m.drawcoastlines() return fig class SuperDict(dict): """just so I can use . notation""" def __getattr__(self, attr): return self[attr] def __setattr__(self, attr, value): self[attr] = value if __name__ == "__main__": H = SuperDict() H.name = 'working example' H.outlat0 = 90 H.numygrid = 180 H.dyout = 1. H.outlon0 = 179 H.numxgrid = 360 H.dxout = 1.0 H.data = np.random.rand(H.numygrid,H.numxgrid) print H.data.shape fig = plot_imshow_custom(H,transform=True) plt.show() print 'it worked' try: H.name = 'Not working example' H.outlat0 = 40 H.numygrid = 100 H.dyout = 0.5 H.outlon0 = 179 H.numxgrid = 110 H.dxout = 0.5 H.data = np.random.rand(H.numygrid,H.numxgrid) fig = plot_imshow_custom(H) print 'huh?' plt.show() except: print "As I said, it's not working..." Jeff Whitaker wrote: > > John [H2O] wrote: >> I'm trying to 'automate' a few components within basemap. I have a pretty >> complicated, and assuredly poorly written, set of functions that allow me >> to >> 'dynamically' plot a grid of data (lon,lat). >> >> Here is one section where I try to deal with transforming the data based >> on >> the projection. 'data' is a grid, often of size 720x360 or 720x180, >> representing full globe or hemisphere at 0.5 degree resolution. >> 'outlon0', >> outlat0', and 'd*out' are the llcrnr coordinates and step. 'transform' is >> an >> option, that is set to True by default: >> >> 1680 ## set up transformations for the data array >> 1681 if m.projection not in ['cyl','merc','mill']: >> 1682 lats = np.arange( outlat0, ( outlat0 + ( numygrid*dyout ) ), >> dyout )[:1] >> 1683 lons = np.arange( outlon0, ( outlon0 + ( numxgrid*dxout ) ), >> dxout )[:1] >> 1684 data = data[:1,:1] >> 1685 else: >> 1686 lats = np.arange( outlat0, ( outlat0 + ( numygrid*dyout ) ), >> dyout ) >> 1687 lons = np.arange( outlon0, ( outlon0 + ( numxgrid*dxout ) ), >> dxout ) >> 1688 >> 1689 ## transform to nx x ny regularly spaced native projection grid >> 1690 if transform: >> 1691 dx = 2.*np.pi*m.rmajor/len(lons) >> 1692 nx = int((m.xmaxm.xmin)/dx)+1; ny = >> int((m.ymaxm.ymin)/dx)+1 >> 1693 if nx is 1: >> 1694 topodat = data >> 1695 else: >> 1696 topodat = m.transform_scalar(data,lons,lats,nx,ny) >> 1697 else: >> 1698 topodat = data >> >> The problem is, when I use the approach with a 'cyl' grid, then >> subsequently >> try to draw the lsmask, I get a failure. Is this approach incorrect? I >> had >> to use the if nx is 1 statement because it was crashing with zero >> division >> error in some cases. >> >> Thanks. >> > John: Please supply us with a selfcontained example triggering the > error that we can run. > > Jeff > >  > Come build with us! The BlackBerry® Developer Conference in SF, CA > is the only developer event you need to attend this year. Jumpstart your > developing skills, take BlackBerry mobile applications to market and stay > ahead of the curve. Join us from November 912, 2009. Register > now! > http://p.sf.net/sfu/devconf > _______________________________________________ > Matplotlibusers mailing list > Matplotlibusers@... > https://lists.sourceforge.net/lists/listinfo/matplotlibusers > >  View this message in context: http://www.nabble.com/basemap%2Ctransform_scalarquestiontp25649437p25795985.html Sent from the matplotlib  users mailing list archive at Nabble.com. 