| trouble with map projections [message #84986] |
Thu, 20 June 2013 11:37  |
chris.orphanides
Messages: 3
Registered: June 2013
|
Junior Member |
|
|
Hello,
I am trying to take a satellite image of sea surface temperature (SST), subset it, and then project it, and I can't seem to get it right. I am able to read the image and subset the area I am interested in without a problem. Getting it projected and having things line up is another story.
What makes sense to me is first using the map_proj_init function to create the map projection I want to put the image into (Lambert Conformal Conic), then use the map_proj_image function to warp the image to the proper projection. However, when I do this, the resulting array has lost its SST values and is all 0.0s. Can anyone tell me what I am doing wrong? I have experimented with many of the mapping capabilities in IDL, but I just can't get it right. The code I described is below Thanks in advance for your help.
range = [-78.1900, 34.0300, -61.8100, 45.4900]
nec_prj = MAP_PROJ_INIT('Lambert Conformal Conic', /GCTP, $
ELLIPSOID='WGS 84', $
LIMIT=range, $
CENTER_LATITUDE=40.00, $
CENTER_LONGITUDE=-70.00, $
STANDARD_PAR1 = 36.1667, $
STANDARD_PAR2 = 43.8333 )
necprj_sst = MAP_PROJ_IMAGE(nec_region, range, MAP_STRUCTURE = nec_prj)
; nec_region is the SST data for the region I am interested in,
; subset to fit the range in the lines above
A little additional information: The main input image before I subset it is described as being in a Cylindrical Lat-Lon projection with a regular 0.01 degree grid and a WGS 84 Ellipsoid. Since Cylindrical is the default IDL projection I didn't set a map projection for this image. I would prefer to set the ellipsoid to WGS 84 for the Cylindrical projection, though it doesn't appear possible in IDL (I would love it if I was wrong about this). Also, I am working with the type of images described here: ( http://podaac.jpl.nasa.gov/dataset/JPL_OUROCEAN-L4UHfnd-GLOB -G1SST)
Thanks.
Chris
|
|
|
|
| Re: trouble with map projections [message #84988 is a reply to message #84986] |
Thu, 20 June 2013 12:25  |
David Fanning
Messages: 11724
Registered: August 2001
|
Senior Member |
|
|
chris.orphanides@noaa.gov writes:
> I am trying to take a satellite image of sea surface temperature
(SST), subset it, and then project it, and I can't seem to get it right.
I am able to read the image and subset the area I am interested in
without a problem. Getting it projected and having things line up is
another story.
>
> What makes sense to me is first using the map_proj_init function to create the map projection I want to put the image into (Lambert Conformal Conic), then use the map_proj_image function to warp the image to the proper projection. However, when I do this, the resulting array has lost its SST values and is all 0.0s. Can anyone tell me what I am doing wrong? I have experimented with many of the mapping capabilities in IDL, but I just can't get it right. The code I
described is below Thanks in advance for your help.
>
> range = [-78.1900, 34.0300, -61.8100, 45.4900]
>
> nec_prj = MAP_PROJ_INIT('Lambert Conformal Conic', /GCTP, $
> ELLIPSOID='WGS 84', $
> LIMIT=range, $
> CENTER_LATITUDE=40.00, $
> CENTER_LONGITUDE=-70.00, $
> STANDARD_PAR1 = 36.1667, $
> STANDARD_PAR2 = 43.8333 )
>
> necprj_sst = MAP_PROJ_IMAGE(nec_region, range, MAP_STRUCTURE = nec_prj)
> ; nec_region is the SST data for the region I am interested in,
> ; subset to fit the range in the lines above
>
> A little additional information: The main input image before I subset it is described as being in a Cylindrical Lat-Lon projection with a regular 0.01 degree grid and a WGS 84 Ellipsoid. Since Cylindrical is the default IDL projection I didn't set a map projection for this image. I would prefer to set the ellipsoid to WGS 84 for the Cylindrical projection, though it doesn't appear possible in IDL (I would love it if I was wrong about this). Also, I am working with the
type of images described here: ( http://podaac.jpl.nasa.gov/dataset/JPL_OUROCEAN-L4UHfnd-GLOB -G1SST)
As the Stooges would say, "No, no, no. You're doing it all wrong!"
You need to create a map projection that describes your image as you
downloaded it. I'm not sure why you think IDL can't do a Cylindrical map
projection with a WGS-84 ellipsoid, but this is a VERY common projection
for satellite images and IDL handles it perfectly. Then, you create a
map projection for what you want the image to end up as. Finally, you
use Map_Image to warp the image from one map projection to the other.
Here is an article that describes the process:
http://www.idlcoyote.com/map_tips/warpimage.html
Cheers,
David
--
David Fanning, Ph.D.
Fanning Software Consulting, Inc.
Coyote's Guide to IDL Programming: http://www.idlcoyote.com/
Sepore ma de ni thue. ("Perhaps thou speakest truth.")
|
|
|
|
| Re: trouble with map projections [message #84992 is a reply to message #84988] |
Thu, 20 June 2013 13:50 |
chris.orphanides
Messages: 3
Registered: June 2013
|
Junior Member |
|
|
On Thursday, June 20, 2013 3:25:20 PM UTC-4, David Fanning wrote:
>
>
>
>> I am trying to take a satellite image of sea surface temperature
>
> (SST), subset it, and then project it, and I can't seem to get it right.
>
> I am able to read the image and subset the area I am interested in
>
> without a problem. Getting it projected and having things line up is
>
> another story.
>
>>
>
>> What makes sense to me is first using the map_proj_init function to create the map projection I want to put the image into (Lambert Conformal Conic), then use the map_proj_image function to warp the image to the proper projection. However, when I do this, the resulting array has lost its SST values and is all 0.0s. Can anyone tell me what I am doing wrong? I have experimented with many of the mapping capabilities in IDL, but I just can't get it right. The code I
>
> described is below Thanks in advance for your help.
>
>>
>
>> range = [-78.1900, 34.0300, -61.8100, 45.4900]
>
>>
>
>> nec_prj = MAP_PROJ_INIT('Lambert Conformal Conic', /GCTP, $
>
>> ELLIPSOID='WGS 84', $
>
>> LIMIT=range, $
>
>> CENTER_LATITUDE=40.00, $
>
>> CENTER_LONGITUDE=-70.00, $
>
>> STANDARD_PAR1 = 36.1667, $
>
>> STANDARD_PAR2 = 43.8333 )
>
>>
>
>> necprj_sst = MAP_PROJ_IMAGE(nec_region, range, MAP_STRUCTURE = nec_prj)
>
>> ; nec_region is the SST data for the region I am interested in,
>
>> ; subset to fit the range in the lines above
>
>>
>
>> A little additional information: The main input image before I subset it is described as being in a Cylindrical Lat-Lon projection with a regular 0.01 degree grid and a WGS 84 Ellipsoid. Since Cylindrical is the default IDL projection I didn't set a map projection for this image. I would prefer to set the ellipsoid to WGS 84 for the Cylindrical projection, though it doesn't appear possible in IDL (I would love it if I was wrong about this). Also, I am working with the
>
> type of images described here: ( http://podaac.jpl.nasa.gov/dataset/JPL_OUROCEAN-L4UHfnd-GLOB -G1SST)
>
>
>
> As the Stooges would say, "No, no, no. You're doing it all wrong!"
>
>
>
> You need to create a map projection that describes your image as you
>
> downloaded it. I'm not sure why you think IDL can't do a Cylindrical map
>
> projection with a WGS-84 ellipsoid, but this is a VERY common projection
>
> for satellite images and IDL handles it perfectly. Then, you create a
>
> map projection for what you want the image to end up as. Finally, you
>
> use Map_Image to warp the image from one map projection to the other.
>
> Here is an article that describes the process:
>
>
>
> http://www.idlcoyote.com/map_tips/warpimage.html
>
>
>
> Cheers,
>
>
>
> David
>
>
>
>
>
> --
>
> David Fanning, Ph.D.
>
> Fanning Software Consulting, Inc.
>
> Coyote's Guide to IDL Programming: http://www.idlcoyote.com/
>
> Sepore ma de ni thue. ("Perhaps thou speakest truth.")
David, thank you for your quick response. I didn't think that I could do a Cylindrical map projection with a WGS 84 Ellipsoid because in the map_proj_init() help page it lists Sphere as the only available ellipsoid when using IDL's own map projections. In the GCTP map projections it says that Equirectangular only takes a sphere as well and doesn't say you can specify the semimajor or semiminor axes. What am I missing here? Does the below work even though it doesn't seem like it should?
g1_prj = MAP_PROJ_INIT('Equirectangular', ELLIPSOID='WGS 84', /GCTP, LIMIT=[-80, -180, 80, 180])
It runs successfully, and when peeking at the result some of it looks right, but I am hesitant.
|
|
|
|
| Re: trouble with map projections [message #84993 is a reply to message #84992] |
Thu, 20 June 2013 14:20 |
David Fanning
Messages: 11724
Registered: August 2001
|
Senior Member |
|
|
chris.orphanides@noaa.gov writes:
> David, thank you for your quick response. I didn't think that I could do a Cylindrical map projection with a WGS 84 Ellipsoid because in the map_proj_init() help page it lists Sphere as the only available ellipsoid when using IDL's own map projections. In the GCTP map projections it says that Equirectangular only takes a sphere as well and doesn't say you can specify the semimajor or semiminor axes. What am I missing here? Does the below work even though it doesn't seem
like it should?
>
> g1_prj = MAP_PROJ_INIT('Equirectangular', ELLIPSOID='WGS 84', /GCTP, LIMIT=[-80, -180, 80, 180])
>
> It runs successfully, and when peeking at the result some of it looks right, but I am hesitant.
Ah, yes, I guess I was thinking of a Cylindrical Equal Area projection,
which was introduced in IDL 8.0.
Yeah, you're probably screwed. :-)
You probably have to use ENVI to get your map projections right. I have
NO idea my MAP_PROJ_INIT allows that ellipsoid, although in the back on
my mind I seem to remember a change that allowed any ellipsoid with map
projections. But, I can't find any mention of it anywhere. Sorry!
I guess your only solace is that on a map with those limits, the
difference between a sphere and a WGS84 ellipsoid are going to be very
small. I've seen a hell of a lot worse in scientific papers. :-)
Cheers,
David
--
David Fanning, Ph.D.
Fanning Software Consulting, Inc.
Coyote's Guide to IDL Programming: http://www.idlcoyote.com/
Sepore ma de ni thue. ("Perhaps thou speakest truth.")
|
|
|
|
| Re: trouble with map projections [message #84994 is a reply to message #84993] |
Thu, 20 June 2013 17:05 |
chris.orphanides
Messages: 3
Registered: June 2013
|
Junior Member |
|
|
On Thursday, June 20, 2013 5:20:27 PM UTC-4, David Fanning wrote:
> chris writes:
>
>
>
>> David, thank you for your quick response. I didn't think that I could do a Cylindrical map projection with a WGS 84 Ellipsoid because in the map_proj_init() help page it lists Sphere as the only available ellipsoid when using IDL's own map projections. In the GCTP map projections it says that Equirectangular only takes a sphere as well and doesn't say you can specify the semimajor or semiminor axes. What am I missing here? Does the below work even though it doesn't seem
>
> like it should?
>
>>
>
>> g1_prj = MAP_PROJ_INIT('Equirectangular', ELLIPSOID='WGS 84', /GCTP, LIMIT=[-80, -180, 80, 180])
>
>>
>
>> It runs successfully, and when peeking at the result some of it looks right, but I am hesitant.
>
>
>
> Ah, yes, I guess I was thinking of a Cylindrical Equal Area projection,
>
> which was introduced in IDL 8.0.
>
>
>
> Yeah, you're probably screwed. :-)
>
>
>
> You probably have to use ENVI to get your map projections right. I have
>
> NO idea my MAP_PROJ_INIT allows that ellipsoid, although in the back on
>
> my mind I seem to remember a change that allowed any ellipsoid with map
>
> projections. But, I can't find any mention of it anywhere. Sorry!
>
>
>
> I guess your only solace is that on a map with those limits, the
>
> difference between a sphere and a WGS84 ellipsoid are going to be very
>
> small. I've seen a hell of a lot worse in scientific papers. :-)
>
>
>
> Cheers,
>
>
>
> David
>
>
>
>
>
>
>
> --
>
> David Fanning, Ph.D.
>
> Fanning Software Consulting, Inc.
>
> Coyote's Guide to IDL Programming: http://www.idlcoyote.com/
>
> Sepore ma de ni thue. ("Perhaps thou speakest truth.")
OK. Thanks. If I print the result of the map_proj_init() and look at some of the !MAP fields, the array in P (whatever that is ) appears to have the semi-major and semi-minor ellipse axes for WGS 84 in the first two values. (Officially the help menu lists P as: "A 16-element, double-precision floating point array indicating additional projection parameters"). Not too helpful. Whether IDL actually uses these additional parameters for anything, I don't know. When I ran it I was surprised that it didn't give me an error and tell me I couldn't use that ellipse. Maybe it stores those ellipse numbers in that field but doesn't do anything with them, I don't know.
I'll try to actually get the mapping to work with a sphere tomorrow, hopefully I can figure that out.
|
|
|
|
comp.lang.idl-pvwave archive
Members
Search
Help
Login
Home
