Hi all
I wanna project air temperature (tmean: [7025 (lon), 3105 (lat)]) over the U.S. from geographic coordinate to Lambert Azimuthal Equal-Area coordinate, and need to know the x,y value of the new variable.
Information of Lambert Azimuthal Equal-Area coordinate:
The central longitude and latitude is -100 and 50, respectively. Both false northing and false easting are 0. The ellipsoid is WGS 84.
I didn't notice much difference among map_proj_image, map_image and map_patch according to the manual, so I used these three routines, and wanna compare them.
I use tmean1, tmean2 and tmean3 to represent results from map_proj_image, map_image and map_patch, respectively.
The spatial pattern of tmean1 and tmean3 are very similar, both of them showed the warped image correctly. Except, the dimension of tmean1 is [7025, 3105], which is the same with original tmean (I can understand this), and the dimension of tmean2 is [1171, 666].
tmean2 is totally wrong, and its dimension is [1, 680].
Could someone help me to clarify the differences among these three routines, please?
It seems all these three routines just warp the image, and won't give any information about the value of x, y. How can I know the x,y of projected image?
Thanks a lot.
Below is my code:
*********************************************
PRO test1
; ------ Read Data ------
data = 'us_tmean.nc'
fid = NCDF_OPEN(data)
NCDF_VARGET, fid, NCDF_VARID(fid, 'tmean'), tmean
NCDF_VARGET, fid, NCDF_VARID(fid, 'lat'), latitude
NCDF_VARGET, fid, NCDF_VARID(fid, 'lon'), longitude
NCDF_CLOSE, fid
tmean = FLOAT(tmean(*,*,0))
tmean(WHERE(tmean LT -9000)) = !values.f_nan
tmean = tmean * 0.01
tmean = REVERSE(tmean, 2)
latitude = REVERSE(latitude, 1)
dims = size(tmean, /dimensions)
llon = dims(0)
llat = dims(1)
lonmin = MIN(longitude)
lonmax = MAX(longitude)
latmin = MIN(latitude)
latmax = MAX(latitude)
limit = [latmin, lonmin, latmax, lonmax]
range = [lonmin, latmin, lonmax, latmax]
; ------ Define the Lambert Azimuthal Equal Area Map Projection ------
sMAP = MAP_PROJ_INIT('LambertAzimuthal', ELLIPSOID='WGS 84', /GCTP, LIMIT=limit, $
CENTER_LONGITUDE=-100, CENTER_LATITUDE=50, $
FALSE_EASTING =0, FALSE_NORTHING = 0)
; ------ warps an image from geographic coordinates (lat/lon) to a specified map projection
tmean1 = MAP_PROJ_IMAGE(tmean, range, MAP_STRUCTURE=sMAP, MASK=mask, MAX_VALUE = 50.0, MIN_VALUE = -50.0, MISSING = -9999, UVRANGE=uvrange, XINDEX=xindex, YINDEX=yindex)
MAP_SET, /LAMBERT, 50, -100, 0, LIMIT = [latmin, lonmin, latmax, lonmax], /GRID, /HORIZON, /USA, /ISOTROPIC
tmean2 = MAP_IMAGE(tmean, startx, starty,latmin=latmin, latmax=latmax, lonmin=lonmin,lonmax=lonmax, $
compress=1, MAP_STRUCTURE=sMAP, mask=mask, max_value = 50.0, min_value=-50.0, missing=-9999)
tmean3 = MAP_PATCH(tmean, longitude, latitude, LAT0=latmin, LAT1=latmax, LON0=lonmin, LON1=lonmax, MAX_VALUE = 50, MISSING = -9999, XSIZE=xsize, XSTART=xstart, YSIZE=ysize, YSTART=ystart)
help, tmean1, tmean2, tmean3
END
*****************************************
Olivia
|
comp.lang.idl-pvwave archive
Members
Search
Help
Login
Home
