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.xmax-m.xmin)/dx)+1;

ny = int((m.ymax-m.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.xmax-m.xmin)/dx)+1; ny =

int((m.ymax-m.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 self-contained 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 9-12, 2009. Register

now!

http://p.sf.net/sfu/devconf

_______________________________________________

Matplotlib-users mailing list

Matplotlib-users@lists.sourceforge.net

https://lists.sourceforge.net/lists/listinfo/matplotlib-users

--

View this message in context: http://www.nabble.com/basemap%2C-transform_scalar-question-tp25649437p25795985.html

Sent from the matplotlib - users mailing list archive at Nabble.com.