On Feb 26, 7:43=A0am, David Fanning <n...@[EMAIL PROTECTED]
> 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/
> Se****e 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 =3D read_ascii(file+'.asc',data_start=3D6)
data0 =3D result.field0001
bad =3D where( data0 eq -9999, complement=3Dgood )
data1 =3D data0 & data1[bad] =3D !values.f_nan
mm =3D minmax( data0[good] )
img =3D bytscl( data1, min=3Dmm[0], max=3Dmm[1], top=3D253 ) + 1
img[ bad ] =3D 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 =3D reverse(data,2)
x0 =3D -2713600 ; from data set header
y0 =3D -2304000
xx =3D [x0,x0,-1*x0,-1*x0] ; the four corners
yy =3D [y0,-1*y0,-1*y0,y0]
;; this is the projection the data is distributed on
stereo =3D map_proj_init('Polar Stereographic', /GCTP, DATUM=3D8, $
CENTER_LONGITUDE=3D0, CENTER_LATITUDE=3D-71 )
lonlat =3D MAP_PROJ_INVERSE( xx, yy, MAP_STRUCTURE=3Dstereo )
longitude =3D reform(lonlat[0,*])
latitude =3D reform(lonlat[1,*])
;; output zoom
limit =3D [ -90, -180, max(latitude), 180 ]
;limit =3D [ -80, 150, -70, 180 ]
;; this is the projection I would like it on
cyl =3D map_proj_init('Cylindrical', limit=3Dlimit)
range =3D [ x0, y0, -1*x0, -1*y0 ]
warp =3D MAP_PROJ_IMAGE( data, range, $
image_structure=3D stereo, $ ;; input
map_structure =3D cyl, $ ;; output
missing =3D -2, $
uvrange =3D uvrange, $
min_value =3D 0, $
_EXTRA=3De )
erase
tv, congrid( warp, !d.x_size, !d.y_size )
pos =3D [0,0,1,1]
;; Pick one. Which one?
uv_box =3D cyl.uv_box
uv_box =3D uvrange
Plot, uv_box[[0, 2]], uv_box[[1, 3]], Position=3Dpos, $
/Nodata, XStyle=3D5, YStyle=3D5, /NoErase
MAP_CONTINENTS, Map_Structure=3Dcyl, /HIRES
map_grid, glinest=3D0, color=3D255, /label, map_structure=3Dcyl
end


|