| Re: Inverse Map Projection Help [message #58810] |
Mon, 18 February 2008 22:30  |
David Fanning
Messages: 11724
Registered: August 2001
|
Senior Member |
|
|
Ken Mankoff writes:
> I've read the previous posts on inverse map projections and the
> lengthy tutorial by David Fanning, but still cannot get things to line
> up quite right. So I'm posting here for help...
>
> I have a data set (BEDMAP) with this information in the header:
>
> ncols 1371
> nrows 1371
> xllcorner -3426225.75
> yllcorner -3426225.75
> cellsize 5000
> NODATA_value -9999
>
> And this information on the website:
>
> Polar Stereographic projection with 71=B0S as the latitude of true scale
> and 0=B0E as the central meridian.
>
> I've managed to load the data, and inverse project it approximately
> such that things roughly line up. But I cannot get it accurate where
> my reference for 'accurate' is the /MAP_CONTINENTS, /HIRES keywords.
Oh, oh. This leads me to think you are using MAP_SET
to set up your map projections, instead of the more
accurate MAP_PROJ_INIT. After my mapping "epiphany" of
a couple of weeks ago, I have given up on the MAP_SET
projections entirely, except--possibly--for cartoon
maps.
And, of course, there is apparently a bug in the
"more accurate" MAP_PROJ_INIT routines, in that if
you use the UVBOX information in the map structure coming
directly from MAP_PROJ_INIT to set up your "data
coordinate space" for map overlays, you will still
be "slightly off". You need to use the UVBOX information
coming from MAP_PROJ_IMAGE for completely accurate
results. I've had a call into ITTVIS for three weeks
about this, but so far without results.
http://www.dfanning.com/map_tips/tiffoverlay.html
There could also be some confusion about whether the
reported corner pixel coordinates are in the center
of the pixel (likely) or on the edge of the pixel.
If it is the center, then you are going to have to
move the coordinates to the edge of the pixel, which
is what IDL needs.
http://www.dfanning.com/map_tips/precipmap.html
Let's see what you are doing. And can your provide
a link to an image?
Cheers,
David
--
David Fanning, Ph.D.
Fanning Software Consulting, Inc.
Coyote's Guide to IDL Programming (www.dfanning.com)
Sepore ma de ni thui. ("Perhaps thou speakest truth.")
|
|
|
|
| Re: Inverse Map Projection Help [message #58848 is a reply to message #58810] |
Fri, 22 February 2008 09:48  |
mankoff
Messages: 131
Registered: March 2004
|
Senior Member |
|
|
On Feb 19, 1:30 am, David Fanning <n...@dfanning.com> wrote:
> Ken Mankoff writes:
>> I've read the previous posts on inverse map projections and the
>> lengthy tutorial by David Fanning, but still cannot get things to line
>> up quite right. So I'm posting here for help...
>
>> I have a data set (BEDMAP) with this information in the header:
>
>> ncols 1371
>> nrows 1371
>> xllcorner -3426225.75
>> yllcorner -3426225.75
>> cellsize 5000
>> NODATA_value -9999
>
>> And this information on the website:
>
>> Polar Stereographic projection with 71=B0S as the latitude of true scale
>> and 0=B0E as the central meridian.
>
>> I've managed to load the data, and inverse project it approximately
>> such that things roughly line up. But I cannot get it accurate where
>> my reference for 'accurate' is the /MAP_CONTINENTS, /HIRES keywords.
>
> Oh, oh. This leads me to think you are using MAP_SET
> to set up your map projections, instead of the more
> accurate MAP_PROJ_INIT. After my mapping "epiphany" of
> a couple of weeks ago, I have given up on the MAP_SET
> projections entirely, except--possibly--for cartoon
> maps.
>
> And, of course, there is apparently a bug in the
> "more accurate" MAP_PROJ_INIT routines, in that if
> you use the UVBOX information in the map structure coming
> directly from MAP_PROJ_INIT to set up your "data
> coordinate space" for map overlays, you will still
> be "slightly off". You need to use the UVBOX information
> coming from MAP_PROJ_IMAGE for completely accurate
> results. I've had a call into ITTVIS for three weeks
> about this, but so far without results.
>
> http://www.dfanning.com/map_tips/tiffoverlay.html
>
> There could also be some confusion about whether the
> reported corner pixel coordinates are in the center
> of the pixel (likely) or on the edge of the pixel.
> If it is the center, then you are going to have to
> move the coordinates to the edge of the pixel, which
> is what IDL needs.
>
> http://www.dfanning.com/map_tips/precipmap.html
>
> Let's see what you are doing. And can your provide
> a link to an image?
>
> Cheers,
>
> David
Hi David,
Thanks for the offer to help. I went over your pages again and began
coding from scratch but am in the same place. I realize there in on
piece of info I have (71°S as the latitude of true scale) that I am
not using.
Here is the result:
http://edgcm.columbia.edu/~mankoff/tmp/unroll.png
And the code looks like this:
restore
data = reverse(data,2)
x0 = -2713600 ; from data set header
y0 = -2304000
xx = [x0,x0,-1*x0,-1*x0] ; the four corners
yy = [y0,-1*y0,-1*y0,y0]
;; this is the projection the data is distributed on
stereo = map_proj_init('Polar Stereographic', /GCTP, DATUM=8, $
CENTER_LONGITUDE=0, CENTER_LATITUDE=-90)
lonlat = MAP_PROJ_INVERSE( xx, yy, MAP_STRUCTURE=stereo )
limit = [ -90, -180, max(latitude), 180 ]
;; this is the projection I would like it on
cyl = map_proj_init('Cylindrical', limit=limit)
range = [ x0, y0, -1*x0, -1*y0 ]
warp = MAP_PROJ_IMAGE( data, range, $
image_structure= stereo, $ ;; input
map_structure = cyl, $ ;; output
missing = -2, $
min_value = 0, $
_EXTRA=e )
erase
tv, congrid( warp, !d.x_size, !d.y_size )
pos = [0,0,1,1]
uv_box = cyl.uv_box
Plot, uv_box[[0, 2]], uv_box[[1, 3]], Position=pos, $
/Nodata, XStyle=5, YStyle=5, /NoErase
MAP_CONTINENTS, Map_Structure=cyl, /HIRES
map_grid, glinest=0, color=255, /label, map_structure=cyl
|
|
|
|
comp.lang.idl-pvwave archive
Members
Search
Help
Login
Home
