John Burkhart wrote:
Dear Jeffrey Whitaker,
I've been working with and using matplotlib and your basemap toolbox. Terrific work! It's a lot of fun to use. I'm currently aboard the R/V Knorr and have no internet - only email - hence my correspondence to you directly rather than via a forum.
I've now used your warpimage.py example to transform and create a basemap from a cylindrically projected image, but now I want to use one which is polar stereographic. Unfortunately, I cannot include attachments, but basically, I'm trying to overly images on this file:
http//:polarview.met.no/highres/sarmap2.jpg and another (c_map1.jpg)How do I go the other way? That is, to use m.transform_scalar to convert a polar stereographic image to another projection?
Thanks so much.
Via InMarsat,
johnHere is the geotiff information for the c_map1.jpg file:
Driver: GTiff/GeoTIFF
Files: ice20080416.tif
Size is 4149, 3150
Coordinate System is:
PROJCS["unnamed",
GEOGCS["WGS 84",
DATUM["WGS_1984",
SPHEROID["WGS 84",6378137,298.2572235629972,
AUTHORITY["EPSG","7030"]],
AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich",0],
UNIT["degree",0.0174532925199433],
AUTHORITY["EPSG","4326"]],
PROJECTION["Polar_Stereographic"],
PARAMETER["latitude_of_origin",90],
PARAMETER["central_meridian",0],
PARAMETER["scale_factor",1],
PARAMETER["false_easting",0],
PARAMETER["false_northing",0],
UNIT["metre",1,
AUTHORITY["EPSG","9001"]]]
Origin = (-1700500.241080038715154,-149240.597878495842451)
Pixel Size = (1000.482160077145636,-1000.482160077145522)
Metadata:
AREA_OR_POINT=Area
Image Structure Metadata:
INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left (-1700500.241, -149240.598) ( 84d59'3.89"W, 74d48'10.68"N)
Lower Left (-1700500.241,-3300759.402) ( 27d15'24.57"W, 57d36'46.68"N)
Upper Right ( 2450500.241, -149240.598) ( 86d30'53.54"E, 68d16'24.80"N)
Lower Right ( 2450500.241,-3300759.402) ( 36d35'25.36"E, 54d20'41.99"N)
Center ( 375000.000,-1725000.000) ( 12d15'53.19"E, 74d17'25.72"N)
Band 1 Block=4149x1 Type=Byte, ColorInterp=Red
Band 2 Block=4149x1 Type=Byte, ColorInterp=Green
Band 3 Block=4149x1 Type=Byte, ColorInterp=Blue
Hi John:
If your image had no borders of annotations, you could read it in and 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:
def interp(datain,xin,yin,xout,yout,checkbounds=False,masked=False,order=1):
"""
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 : Jeffrey S. Whitaker: NOAA Physical Sciences Laboratory