On Feb 26, 7:43 am, David Fanning <n...@dfanning.com> wrote:
> mankoff writes:
>> It does help. Image is better aligned. But still not accurate :(.
>
> I just got back to my office and I'm doing the usual
> up-at-3AM-thing for a week or so. Do you mean "not accurate"
> in the way using the UV-BOX from the map structure, rather
> from the MAP_PROJ_IMAGE UV-BOX, is not accurate? This wouldn't
> surprise me. Did you try using MAP_PROJ_IMAGE for creating
> the UV-BOX, as I outlined in my article?
>
> If you make the data available, I'll schedule an appointment
> for tomorrow at 4AM. :-)
>
> Cheers,
>
> David
> --
> David Fanning, Ph.D.
> Fanning Software Consulting, Inc.
> Coyote's Guide to IDL Programming:http://www.dfanning.com/
> Sepore ma de ni thui. ("Perhaps thou speakest truth.")
Here is the data set website:
http://www.antarctica.ac.uk/bas_research/data/access/bedmap/ download/
And the actual one I've started with is:
http://www.antarctica.ac.uk/bas_research/data/access/bedmap/ download/surface.asc.gz
Code to read in this file (once un-gzipped) is:
pro load_asc, file, data0, data1, img
if not keyword_set(file) then begin
print, 'bathy, bedelev, groundbed, icethic, surface, water'
return
end
result = read_ascii(file+'.asc',data_start=6)
data0 = result.field0001
bad = where( data0 eq -9999, complement=good )
data1 = data0 & data1[bad] = !values.f_nan
mm = minmax( data0[good] )
img = bytscl( data1, min=mm[0], max=mm[1], top=253 ) + 1
img[ bad ] = 0
end
My code to attempt to 'unroll' this data is above in this thread, and
re-pasted here (slightly different than above perhaps... 2 days
later). Note that I have uv_box from both map_proj_init and
map_proj_image. I think the map_proj_image code provides slightly
better match. It appears to mach East/West perfectly (?) but there is
still a north/south error.
pro unroll_foo
end
;; load the data
load_asc,'surface', data & save, data
;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=-71 )
lonlat = MAP_PROJ_INVERSE( xx, yy, MAP_STRUCTURE=stereo )
longitude = reform(lonlat[0,*])
latitude = reform(lonlat[1,*])
;; output zoom
limit = [ -90, -180, max(latitude), 180 ]
;limit = [ -80, 150, -70, 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, $
uvrange = uvrange, $
min_value = 0, $
_EXTRA=e )
erase
tv, congrid( warp, !d.x_size, !d.y_size )
pos = [0,0,1,1]
;; Pick one. Which one?
uv_box = cyl.uv_box
uv_box = uvrange
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
end
|
comp.lang.idl-pvwave archive
Members
Search
Help
Login
Home
