|
Re: map: integerized sinusoidal [message #65224 is a reply to message #65221] |
Fri, 20 February 2009 17:12  |
ben.bighair
Messages: 221 Registered: April 2007
|
Senior Member |
|
|
On Feb 20, 5:18 pm, David Fanning <n...@dfanning.com> wrote:
> Roman writes:
>> After trying out for several days to map (hundreds of) MODIS MOD14A2
>> products in IDL, I hope someone can help me here ! I went through
>> older comments on that but I really don't know where is my mistake !
>> The dataset and the boundaries of the continents/coast lines dont fitt
>> to each other. In the upper right corner yes but than it gets worse to
>> the lower right corner !
>
>> Thank you in advance for any help - cheers Roman
>> --------------------
>> Datasets are HDF files, 1200 x 1200 rows/columns, resoultion 1 km in
>> integerized sinusoidal projection.
>
> This lousy problem has taken over my life! :-(
>
> I don't know what it is, but I hate it when I think
> I know something and I don't. (I don't mind being dumb,
> but I don't like to be stupid.) So, because I thought
> this was a simple problem, I got most annoyed when
> my simple solutions wouldn't work with it. And when
> I get annoyed, I do tend to get a little over-the-top
> tenacious.
>
> But, as with so much in life, there is always more to learn.
> I *finally* figured out how to do this in a simple way.
> I know a little bit more about this now then I knew a couple
> of days ago. And the price I paid to gain this knowledge was
> pretty dear, so I'm not likely to forget it. That's a plus,
> too, I guess.
>
> Anyway, I wrote an article about this in case anyone else
> wants to avoid banging their head against a wall. :-)
>
> http://www.dfanning.con/map_tips/modis_overlay.html
>
> And a special thanks to Garrison Keillor and Ana Ivanovic
> for inspiration. And, of course, to Roman and to all the
> others on this newsgroup whose simple questions are sometimes
> the hardest to answer.
>
Wow! Gives a fellah a hankerin' for a corned beef sandwich with a cup
of black coffee - just like Carson Wiler would have!
Nice!
|
|
|
Re: map: integerized sinusoidal [message #65225 is a reply to message #65224] |
Fri, 20 February 2009 14:18  |
David Fanning
Messages: 11724 Registered: August 2001
|
Senior Member |
|
|
Roman writes:
> After trying out for several days to map (hundreds of) MODIS MOD14A2
> products in IDL, I hope someone can help me here ! I went through
> older comments on that but I really don't know where is my mistake !
> The dataset and the boundaries of the continents/coast lines dont fitt
> to each other. In the upper right corner yes but than it gets worse to
> the lower right corner !
>
> Thank you in advance for any help - cheers Roman
> --------------------
> Datasets are HDF files, 1200 x 1200 rows/columns, resoultion 1 km in
> integerized sinusoidal projection.
This lousy problem has taken over my life! :-(
I don't know what it is, but I hate it when I think
I know something and I don't. (I don't mind being dumb,
but I don't like to be stupid.) So, because I thought
this was a simple problem, I got most annoyed when
my simple solutions wouldn't work with it. And when
I get annoyed, I do tend to get a little over-the-top
tenacious.
But, as with so much in life, there is always more to learn.
I *finally* figured out how to do this in a simple way.
I know a little bit more about this now then I knew a couple
of days ago. And the price I paid to gain this knowledge was
pretty dear, so I'm not likely to forget it. That's a plus,
too, I guess.
Anyway, I wrote an article about this in case anyone else
wants to avoid banging their head against a wall. :-)
http://www.dfanning.con/map_tips/modis_overlay.html
And a special thanks to Garrison Keillor and Ana Ivanovic
for inspiration. And, of course, to Roman and to all the
others on this newsgroup whose simple questions are sometimes
the hardest to answer.
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.")
|
|
|
Re: map: integerized sinusoidal [message #65250 is a reply to message #65225] |
Wed, 18 February 2009 08:47  |
Roman
Messages: 3 Registered: February 2009
|
Junior Member |
|
|
On 18 Feb., 17:31, David Fanning <n...@dfanning.com> wrote:
> Roman writes:
>> Thanks for all your comments ! Below you find the way I did it
>> succesfully in IDL:
>
> Well, I have to admit, it appears to be a miracle to me.
> But if it works, it works. :-)
>
> One question:
>
>> ; transform into UTM coordinates
>> mapStruct = map_proj_init('integerized sinusoidal', $
>> center_longitude=0., false_northing = 0., $
>> false_easting = 0., is_zones=86400.0, $
>> sphere_radius = 6371007.181)
>
> Where did you come up with that number, 86400, for the
> IS_ZONES keyword? I thought there were a total of 648
> "tiles" or "zones" in an integerized sinusoidal projection.
>
> Cheers,
>
> David
> --
> David Fanning, Ph.D.
> Coyote's Guide to IDL Programming (www.dfanning.com)
> Sepore ma de ni thui. ("Perhaps thou speakest truth.")
the details I got from here:
http://modland.nascom.nasa.gov/developers/is_tiles/is_gctp.h tml
I will try to install the MODIS RTool to cross check the reprojection.
But as you say: it works in IDL (IDL 6.4) :) !
Cheers, Roman
|
|
|
Re: map: integerized sinusoidal [message #65251 is a reply to message #65250] |
Wed, 18 February 2009 08:31  |
David Fanning
Messages: 11724 Registered: August 2001
|
Senior Member |
|
|
Roman writes:
> Thanks for all your comments ! Below you find the way I did it
> succesfully in IDL:
Well, I have to admit, it appears to be a miracle to me.
But if it works, it works. :-)
One question:
> ; transform into UTM coordinates
> mapStruct = map_proj_init('integerized sinusoidal', $
> center_longitude=0., false_northing = 0., $
> false_easting = 0., is_zones=86400.0, $
> sphere_radius = 6371007.181)
Where did you come up with that number, 86400, for the
IS_ZONES keyword? I thought there were a total of 648
"tiles" or "zones" in an integerized sinusoidal projection.
Cheers,
David
--
David Fanning, Ph.D.
Coyote's Guide to IDL Programming (www.dfanning.com)
Sepore ma de ni thui. ("Perhaps thou speakest truth.")
|
|
|
Re: map: integerized sinusoidal [message #65252 is a reply to message #65251] |
Wed, 18 February 2009 07:40  |
Roman
Messages: 3 Registered: February 2009
|
Junior Member |
|
|
On 17 Feb., 21:59, David Fanning <n...@dfanning.com> wrote:
> Roman writes:
>> After trying out for several days to map (hundreds of) MODIS MOD14A2
>> products in IDL, I hope someone can help me here ! I went through
>> older comments on that but I really don't know where is my mistake !
>> The dataset and the boundaries of the continents/coast lines dont fitt
>> to each other. In the upper right corner yes but than it gets worse to
>> the lower right corner !
>
>> Thank you in advance for any help - cheers Roman
>> --------------------
>> Datasets are HDF files, 1200 x 1200 rows/columns, resoultion 1 km in
>> integerized sinusoidal projection.
>
> Humm. Well, you seem surprised by this result, but nowhere
> in your code do you set up an integerized sinusoidal projection.
> You appear to be using a cylindrical projection. I'm pretty
> sure that if the data is in one projection, and the overlays
> are in another, that the chances of good alignment are pretty
> poor, indeed. :-)
>
> IDL does have an integerized sinusoidal projection, but
> you will have to use MAP_PROJ_INIT to set it up, and you
> will have to pass the map structure created by MAP_PROJ_INIT
> to Map_Continents and Map_Grid, so they will be able to
> draw on the map. I've never used this projection, but it has
> a number of 131. I am not sure it will be available in your
> version of IDL. It may have been added *after* IDL 5.6.
>
> Your code will look something like this:
>
> window, xsize=800, ysize=800
> position = [0.1, 0.1, 0.9, 0.9]
> mapStruct = Map_Proj_Init(131, ...)
> TVScale, image, Position=position, /KEEP
> Plot, mapStruct.uv_box[[0,2]], mapStruct.uv_box[[1,3]], $
> Pos=position, /NODATA, /NOERASE
> Map_Continents, MAP_STRUCTURE=mapStruct
> Map_Grid, MAP_STRUCTURE=mapStruct
>
> Cheers,
>
> David
>
> --
> David Fanning, Ph.D.
> Coyote's Guide to IDL Programming (www.dfanning.com)
> Sepore ma de ni thui. ("Perhaps thou speakest truth.")
Thanks for all your comments ! Below you find the way I did it
succesfully in IDL: Cheers, Roman
....
....
; create geolocation for specific tiles (lat/lon for each grid point)
lats=replicate(1.d,nx)#((90.d)-double(ty)/(18.d)*(180.d)-(di ndgen(ny)+.
5d)/double(wy)*(180.d))
lons=((dindgen(nx)+tx*nx)/double(wx)*(360.d)-(180.d))#((1.d) /(sin
((dindgen(ny)+ty*ny)*!dpi/double(wy))>1E-8))
; transform into UTM coordinates
mapStruct = map_proj_init('integerized sinusoidal', center_longitude
= 0., false_northing = 0., false_easting = 0., is_zones=86400.0,
sphere_radius = 6371007.181)
limit = map_proj_inverse(lons,lats, map_structure=mapStruct)
xy=map_proj_forward(lons,lats,map_structure=mapStruct)
xutm=dblarr(nx,ny)
yutm=dblarr(nx,ny)
xutm[*]=xy[0,*]
yutm[*]=xy[1,*]
xy=map_proj_forward(lonmin,latmin,map_structure=mapStruct)
xmin=xy[0]
ymin=xy[1]
xy=map_proj_forward(lonmax,latmax,map_structure=mapStruct)
xmax=xy[0]
ymax=xy[1]
print,'UTM boundary coordinates for the chosen area:
',xmin,xmax,ymin,ymax
;; the actual triangulation and reprojection goes here
triangulate,xutm,yutm,tr
dx=out_resolution*1000.
dy=out_resolution*1000.
....
....
Map_Set,/grid,/continent,/HIRES, Position=[0.05, 0.05, 0.95, 0.80],$
charsize=1,LIMIT=[x2,y1,x1,y2]
imagemap,congrid(dummy,mag*nx,mag*ny), congrid(lats,mag*nx,mag*ny,/
interp), congrid(lons,mag*nx,mag*ny,/interp),$
PIXRANGE=[nodata,maxdata],colrange=[0,254],/
nomap,thick=2,xthick=2,ythick=2
Map_Continents,/COUNTRIES, /COASTS,/HIRES,mlinethick=2
Map_Grid,/box_axes,/clip_text,charsize=1.5,charthick=2,/HORI ZON
|
|
|
Re: map: integerized sinusoidal [message #65264 is a reply to message #65252] |
Tue, 17 February 2009 12:59  |
David Fanning
Messages: 11724 Registered: August 2001
|
Senior Member |
|
|
Roman writes:
> After trying out for several days to map (hundreds of) MODIS MOD14A2
> products in IDL, I hope someone can help me here ! I went through
> older comments on that but I really don't know where is my mistake !
> The dataset and the boundaries of the continents/coast lines dont fitt
> to each other. In the upper right corner yes but than it gets worse to
> the lower right corner !
>
> Thank you in advance for any help - cheers Roman
> --------------------
> Datasets are HDF files, 1200 x 1200 rows/columns, resoultion 1 km in
> integerized sinusoidal projection.
Humm. Well, you seem surprised by this result, but nowhere
in your code do you set up an integerized sinusoidal projection.
You appear to be using a cylindrical projection. I'm pretty
sure that if the data is in one projection, and the overlays
are in another, that the chances of good alignment are pretty
poor, indeed. :-)
IDL does have an integerized sinusoidal projection, but
you will have to use MAP_PROJ_INIT to set it up, and you
will have to pass the map structure created by MAP_PROJ_INIT
to Map_Continents and Map_Grid, so they will be able to
draw on the map. I've never used this projection, but it has
a number of 131. I am not sure it will be available in your
version of IDL. It may have been added *after* IDL 5.6.
Your code will look something like this:
window, xsize=800, ysize=800
position = [0.1, 0.1, 0.9, 0.9]
mapStruct = Map_Proj_Init(131, ...)
TVScale, image, Position=position, /KEEP
Plot, mapStruct.uv_box[[0,2]], mapStruct.uv_box[[1,3]], $
Pos=position, /NODATA, /NOERASE
Map_Continents, MAP_STRUCTURE=mapStruct
Map_Grid, MAP_STRUCTURE=mapStruct
Cheers,
David
--
David Fanning, Ph.D.
Coyote's Guide to IDL Programming (www.dfanning.com)
Sepore ma de ni thui. ("Perhaps thou speakest truth.")
|
|
|
Re: map: integerized sinusoidal [message #65266 is a reply to message #65264] |
Tue, 17 February 2009 11:09  |
liamgumley
Messages: 74 Registered: June 2005
|
Member |
|
|
On Feb 17, 10:59 am, Roman <sinus_oi...@gmx.ch> wrote:
> Hi there
>
> After trying out for several days to map (hundreds of) MODIS MOD14A2
> products in IDL, I hope someone can help me here ! I went through
> older comments on that but I really don't know where is my mistake !
> The dataset and the boundaries of the continents/coast lines dont fitt
> to each other. In the upper right corner yes but than it gets worse to
> the lower right corner !
>
> Thank you in advance for any help - cheers Roman
> --------------------
> Datasets are HDF files, 1200 x 1200 rows/columns, resoultion 1 km in
> integerized sinusoidal projection.
You will save yourself a lot of time and aggravation if you use the
MODIS Reprojection Tool to convert the MOD14A2 ISIN tiles to a more
workable projection.
http://igskmncnwb001.cr.usgs.gov/landdaac/tools/modis/
Cheers,
Liam.
Practical IDL Programming
http://www.gumley.com/
|
|
|