From: Mike Bauer <mbauer@gi...>  20100225 16:01:38

Howdy All, I'm hoping someone can give me a quick solution to a couple of problems. I think I'm just missing an idea or two. Problem 1: I'm creating a map using the 'llc' lambert conformal projection and pcolormesh. Here is a sampling of the source. self.m = Basemap(lat_0=self.lat_0,lon_0=self.lon_0,projection='lcc', width=store_comp.base_width,height=store_comp.base_height, resolution='i',area_thresh=100000) self.fig = plt.figure(figsize=(width,hieght),frameon=True) self.ax = self.fig.add_subplot(111) self.xx, self.yy = self.m(*numpy.meshgrid(self.x,self.y)) self.the_image = self.m.pcolormesh(self.xx,self.yy,z,edgecolors='None',cmap=self.color_scheme) The problem I have is two fold: 1) the map segment isn't fully shown unless I drive up the width and size in basemap so that the map floats in a lot of whitespace. 2) when I use drawparallels etc. the lines extend beyond the map in a way that I wish they wouldn't. See map1 in http://gallery.me.com/ohtinsel#100149 Problem 2: This problem comes up when I use contourf on the same data, which occupies only a limited domain (i.e., there is no data outside the lat/lon bounds shown in map1). Here the contours spill out onto the plot in a way that I wish they wound't (some of the source is below). See map2 in http://gallery.me.com/ohtinsel#100149 self.x = numpy.where(self.x < 180.0,self.x,self.x360.0) scale = 1 dx = self.width/((len(self.x)1)*scale) dy = self.height/((len(self.y)1)*scale) nx = int((self.m.xmaxself.m.xmin)/dx)+1 ny = int((self.m.ymaxself.m.ymin)/dy)+1 self.z,self.xx,self.yy = self.m.transform_scalar( self.z,self.x,self.y,nx,ny,returnxy=True) self.the_image = self.m.contour(self.xx,self.yy,self.z,colors='k') self.the_image = self.m.contourf(self.xx,self.yy,self.z, cmap=self.color_scheme,extend=self.extend) No doubt I'm doing something wrong and probably obvious, but I can't figure it out. Suggestions are much appreciated. Mike 
From: Jeff Whitaker <jswhit@fa...>  20100225 17:16:38

Mike Bauer wrote: > Howdy All, > > I'm hoping someone can give me a quick solution to a couple of > problems. I think I'm just missing an idea or two. > > Problem 1: I'm creating a map using the 'llc' lambert conformal > projection and pcolormesh. Here is a sampling of the source. > > self.m = Basemap(lat_0=self.lat_0,lon_0=self.lon_0,projection='lcc', > > width=store_comp.base_width,height=store_comp.base_height, > resolution='i',area_thresh=100000) > self.fig = plt.figure(figsize=(width,hieght),frameon=True) > self.ax = self.fig.add_subplot(111) > self.xx, self.yy = self.m(*numpy.meshgrid(self.x,self.y)) > self.the_image = > self.m.pcolormesh(self.xx,self.yy,z,edgecolors='None',cmap=self.color_scheme) > > The problem I have is two fold: 1) the map segment isn't > fully shown unless I drive up the width and size in basemap so that > the map floats in a lot of whitespace. 2) when I use drawparallels > etc. the lines extend beyond the map in a way that I wish they > wouldn't. See map1 in http://gallery.me.com/ohtinsel#100149 > > Problem 2: This problem comes up when I use contourf on the same > data, which occupies only a limited domain (i.e., there is no data > outside the lat/lon bounds shown in map1). Here the contours spill out > onto the plot in a way that I wish they wound't (some of the source > is below). See map2 in http://gallery.me.com/ohtinsel#100149 > > self.x = numpy.where(self.x < 180.0,self.x,self.x360.0) > scale = 1 > dx = self.width/((len(self.x)1)*scale) > dy = self.height/((len(self.y)1)*scale) > nx = int((self.m.xmaxself.m.xmin)/dx)+1 > ny = int((self.m.ymaxself.m.ymin)/dy)+1 > self.z,self.xx,self.yy = self.m.transform_scalar( > self.z,self.x,self.y,nx,ny,returnxy=True) > self.the_image = > self.m.contour(self.xx,self.yy,self.z,colors='k') > self.the_image = self.m.contourf(self.xx,self.yy,self.z, > cmap=self.color_scheme,extend=self.extend) > > No doubt I'm doing something wrong and probably obvious, but I can't > figure it out. Suggestions are much appreciated. > > Mike > > > Mike: It's difficult to tell what's going on (or even what you're really trying to do) from the code snippets you posted. If you can send a complete, selfcontained example of the problem you're having that anybody can run, then we can probably help. Jeff  Jeffrey S. Whitaker Phone : (303)4976313 Meteorologist FAX : (303)4976449 NOAA/OAR/PSD R/PSD1 Email : Jeffrey.S.Whitaker@... 325 Broadway Office : Skaggs Research Cntr 1D113 Boulder, CO, USA 803033328 Web : http://tinyurl.com/5telg 
From: Jeff Whitaker <jswhit@fa...>  20100225 20:19:05

Mike Bauer wrote: > Jeff, > > > Attached you'll find an example and sample data. Sample output is > here: http://gallery.me.com/ohtinsel#100149 > > Below you'll also find an old email about the project that you've > answered before (might help you understand what I'm doing). > > Thank you so much. This is a huge help. > > Mike > Mike: OK, I figure out why I have a black plot. In the latest version of netCDF4, scale_factor and add_offset are applied automatically, so in your script they were effectively being applied twice. Concerning the contours spilling out over the domain, just set the keyword masked to True in the call to transform_scalar. Concerning the "floating whitespace"  your data is on a lat/lon grid, so just plot the data with projection='cyl' and specify the lat/lon values of the lower left and upper right corners of the grid. Then the data will exactly fit the map. Your are using a lambert conformal projection, which covers a much larger area than your data. Jeff > mbauer wrote: > >> Thanks Jeff, >> >> To clarify, I'm sampling a numpy array (regular lon/lat grid) and >> extracting a series of same size frames (say 60 longitude grids and 30 >> latitude grids) around a feature of interest, which can be centered >> somewhere on the map. What I want to do is accumulate statistics with >> these frames such that the relative size/distances are persevered, >> which of course means that I can't just add a frame centered on 30N >> with one centered on 80N. Ideally, I'd like to interpolate each frame >> to a common point (lon/lat) and display the results either in the >> common grid space or as radial distances from the common point. >> >> Since you're a meteorologist I can simply say I'm creating an ensemble >> average of extra tropical cyclones from a dozen or so computer models >> (each with very different resolutions). I want to see how cloud and >> precipitation features in each model's cyclones compare to a similar >> product I'm producing from satellite data using weather model output >> to locate the cyclones. Much the same thing as the link I provided. >> >> Thanks for your suggests as transform_scalar sounds like a good place >> to begin. >> >> Mike >> > > Mike: Thanks for the explanation, I get it now  and I think I have > just the thing for you. First, define a Basemap instance for a Lambert > Conformal projection centered on each of you frames that is 5000 km wide > and 5000 km tall. > > # for frame centered on lon_0, lat_0. > # resolution=None skips processing of boundary datasets to save time. > m = Basemap(lon_0=lon_0, lat_0=lat_0, projection='lcc', \ > width=5000000, height=5000000, resolution=None) > > > Now interpolate your data to a nx by ny grid that is regular in the map > projection region. > > # data is the lat/lon gridded data, lons and lats are 1D arrays in degrees. > data2 = m.transform_scalar(data, lons, lats, nx, ny) > > The data2 grids will be approximately equally spaced on the surface of > the earth, regardless on lat_0 and lon_0. Therefore, you should be able > to just ensemble average all the data2 grids and preserve the relative > shapes and sizes of the features (this is helped by the fact that the > projection is conformal, or shapepreserving). > > HTH, > > Jeff > > ================================================================================ > > On Thu, Feb 25, 2010 at 12:16 PM, Jeff Whitaker <jswhit@...> wrote: > >> Mike Bauer wrote: >> >>> Howdy All, >>> >>> I'm hoping someone can give me a quick solution to a couple of >>> problems. I think I'm just missing an idea or two. >>> >>> Problem 1: I'm creating a map using the 'llc' lambert conformal >>> projection and pcolormesh. Here is a sampling of the source. >>> >>> self.m = >>> Basemap(lat_0=self.lat_0,lon_0=self.lon_0,projection='lcc', >>> >>> width=store_comp.base_width,height=store_comp.base_height, >>> resolution='i',area_thresh=100000) >>> self.fig = plt.figure(figsize=(width,hieght),frameon=True) >>> self.ax = self.fig.add_subplot(111) >>> self.xx, self.yy = self.m(*numpy.meshgrid(self.x,self.y)) >>> self.the_image = >>> >>> self.m.pcolormesh(self.xx,self.yy,z,edgecolors='None',cmap=self.color_scheme) >>> >>> The problem I have is two fold: 1) the map segment isn't >>> fully shown unless I drive up the width and size in basemap so that >>> the map floats in a lot of whitespace. 2) when I use drawparallels >>> etc. the lines extend beyond the map in a way that I wish they >>> wouldn't. See map1 in http://gallery.me.com/ohtinsel#100149 >>> >>> Problem 2: This problem comes up when I use contourf on the same >>> data, which occupies only a limited domain (i.e., there is no data >>> outside the lat/lon bounds shown in map1). Here the contours spill out >>> onto the plot in a way that I wish they wound't (some of the source >>> is below). See map2 in http://gallery.me.com/ohtinsel#100149 >>> >>> self.x = numpy.where(self.x < 180.0,self.x,self.x360.0) >>> scale = 1 >>> dx = self.width/((len(self.x)1)*scale) >>> dy = self.height/((len(self.y)1)*scale) >>> nx = int((self.m.xmaxself.m.xmin)/dx)+1 >>> ny = int((self.m.ymaxself.m.ymin)/dy)+1 >>> self.z,self.xx,self.yy = self.m.transform_scalar( >>> self.z,self.x,self.y,nx,ny,returnxy=True) >>> self.the_image = >>> self.m.contour(self.xx,self.yy,self.z,colors='k') >>> self.the_image = self.m.contourf(self.xx,self.yy,self.z, >>> cmap=self.color_scheme,extend=self.extend) >>> >>> No doubt I'm doing something wrong and probably obvious, but I can't >>> figure it out. Suggestions are much appreciated. >>> >>> Mike >>> >>> >>> >>> >> Mike: It's difficult to tell what's going on (or even what you're really >> trying to do) from the code snippets you posted. If you can send a >> complete, selfcontained example of the problem you're having that anybody >> can run, then we can probably help. >> >> Jeff >> >> >>  >> Jeffrey S. Whitaker Phone : (303)4976313 >> Meteorologist FAX : (303)4976449 >> NOAA/OAR/PSD R/PSD1 Email : Jeffrey.S.Whitaker@... >> 325 Broadway Office : Skaggs Research Cntr 1D113 >> Boulder, CO, USA 803033328 Web : http://tinyurl.com/5telg >> >> >> >> >> > >  Jeffrey S. Whitaker Phone : (303)4976313 Meteorologist FAX : (303)4976449 NOAA/OAR/PSD R/PSD1 Email : Jeffrey.S.Whitaker@... 325 Broadway Office : Skaggs Research Cntr 1D113 Boulder, CO, USA 803033328 Web : http://tinyurl.com/5telg 