From: Jeff Whitaker <jswhit@fa...>  20071212 16:54:09

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 Dec 11, 2007, at 4:57 PM, Jeff Whitaker wrote: > >> mbauer wrote: >>> Matplotlib users, I looking to tap your wealth of ideas and >>> experience to help solve a problem I'm working on. >>> >>> The problem: I have a series of 2d scalar arrays representing a >>> fixed width/height lon/lat box centered on an arbitrary lon/lat. I >>> need to average these composites on a common basis that >>> accommodates the scale changes due to latitude, preferably by >>> shifting everything to a common central lon/lat (a polar/radial >>> distance basis would work too). I want a plot of the end result too >>> and I'm like to do everything with matplotlib and python so that it >>> folds into the rest of my program. >>> >>> Something similar can be seen at >>> http://www.atmos.washington.edu/~robwood/topic_cyclones.htm >>> >>> I've been looking at transform_scalar from basemap but I'm not >>> quite sure this is what I should use. >>> >> Mike: >> >> transform_scalar does simple bilinear interpolation from a lat/lon >> grid to a regular grid in map projection coordinates. If your map >> projection is just a lat/lon projection, then this amounts to >> interpolating from one lat/lon grid to another. >>> If anyone can offer a solution, a point in the right direction, or >>> just wave me off this path I'd be most appreciative. >>> >> I'm sure numpy/matplotlib can do what you need to do. Matplotlib >> can certainly make a plot similar to the one given in your link. I >> think you question relates more to the processing of your arrays >> though, and not specifically the plotting. Are all your 2d arrays >> the same shape (the same number of lats and lons)? Are they just >> centered on different regions? If so, I think you can just multiply >> each grid point by the cosine of latitude to get the proper area >> weighting before summing them together. But perhaps I'm missing the >> essence of your question .... >> >> 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 1D124 >> 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 1D124 Boulder, CO, USA 803033328 Web : http://tinyurl.com/5telg 