# plotting / overlay. help via satellite

John Burkhart wrote:

Jeff,

Apologies for dual emails... my earlier one wasn't meant to be sent...

I was hoping you could expand slightly (or at least provide an example) on #1 below. As I said earlier, I now have the GEOTIFF information which should help, and a customized version of the image without the white space. Next, is just to get the raw data! Probably easier - but it seems this is a good exercise nevertheless.

Regards,
john

John: It would definitely be easier to have the raw data. Regarding step number 1, you can do something like this

# m1 is the Basemap instance for map projection 1 (the original
# projection that the GEOTIFF image is in)
# m2 is the Basemap instance for map projection 2 (what you want
# to interpolate the image to)
# x2,y2 are the map projection coordinates of the m2 grid.
lons2, lats2 = m2(x2, y2, inverse=True) # x2, y2 must be 2d arrays
# x1, y1 is the m2 grid in m1 coordinates.
x1, y1 = m1(lons2, lat2)
# if x,y are 1d arrays defining the m1 grid, and data1 is the data on
# the m1 grid, here's how to interpolate to the m2 grid.
# data2 is now data1 interpolated to the m2 grid.
data2 = interp(data1,x,y,x1,y1)

-Jeff

···

get the RGB values of each pixel using PIL (following the example in
warpimage.py). You could then

1) compute the polar stereographic coordinates of the rectilinear grid
you want to interpolate to

2) use the interp function to interpolate the RGB values from the
original polar stereographic grid to the new grid.

Here's part of the docstring for the interp function:

"""
dataout = interp(datain,xin,yin,xout,yout,order=1)

interpolate data (datain) on a rectilinear grid (with x=xin
y=yin) to a grid with x=xout, y=yout.

datain is a rank-2 array with 1st dimension corresponding to y,
2nd dimension x.

xin, yin are rank-1 arrays containing x and y of
datain grid in increasing order.

xout, yout are rank-2 arrays containing x and y of desired output grid."""

Here xin and yin would be the (1d) polar stereographic coords of the
original image grid. xout, yout would be the (2d) coordinates of the
new grid (in the same polar stereographic coordinates as the original
image grid, even though the new grid is a different map projection).

You can use the Basemap instances defined for each projection to
compute the coordinates of each grid, and to transform the new grid
into the projection coordinates of the original grid.

It's tricky, but should be possible if the image doesn't have any
whitespace or annotations around the edges. Unfortunately, the image
you point to doesn't appear to be that simple.

I'm copying the matplotlib-users list just in case anyone has a better
suggestion ...

-Jeff

--
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-124
Boulder, CO, USA 80303-3328 Web : http://tinyurl.com/5telg

----------------------------------------------------------------
This message was sent using IMP, the Internet Messaging Program.

--
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-124
Boulder, CO, USA 80303-3328 Web : http://tinyurl.com/5telg