[GRASSLIST:6319] Goode in Grass

Hi

Could somebody please share the PROJ_INFO for Goode's projection? Or the proj4 definition?

I'd like to work on MODIS imagery of Poland, Eastern Europe. I couldn't figure out the appropriate projection definition by myself.

Thanks.

Maciek

On Sun, Apr 03, 2005 at 11:26:52AM +0200, Maciek Sieczka wrote:

Hi

Could somebody please share the PROJ_INFO for Goode's projection? Or the
proj4 definition?

I'd like to work on MODIS imagery of Poland, Eastern Europe. I couldn't
figure out the appropriate projection definition by myself.

Most MODIS maps are in Sinusoidal (on MODIS sphere).
Where did you get the map?

Markus

From: "Markus Neteler" <neteler@itc.it>

On Sun, Apr 03, 2005 at 11:26:52AM +0200, Maciek Sieczka wrote:

Hi

Could somebody please share the PROJ_INFO for Goode's projection? Or the
proj4 definition?

I'd like to work on MODIS imagery of Poland, Eastern Europe. I couldn't
figure out the appropriate projection definition by myself.

Most MODIS maps are in Sinusoidal (on MODIS sphere).
Where did you get the map?

say ftp://ftp.glcf.umiacs.umd.edu/modis/500m/Eurasia/Goodes.EUAS.2000321/
via the GLCF UMIACS www

Maciek

On Tue, Apr 05, 2005 at 10:13:58AM +0200, Maciek Sieczka wrote:

From: "Markus Neteler" <neteler@itc.it>
>On Sun, Apr 03, 2005 at 11:26:52AM +0200, Maciek Sieczka wrote:
>>Hi
>>
>>Could somebody please share the PROJ_INFO for Goode's projection? Or the
>>proj4 definition?
>>
>>I'd like to work on MODIS imagery of Poland, Eastern Europe. I couldn't
>>figure out the appropriate projection definition by myself.

>Most MODIS maps are in Sinusoidal (on MODIS sphere).
>Where did you get the map?

say ftp://ftp.glcf.umiacs.umd.edu/modis/500m/Eurasia/Goodes.EUAS.2000321/
via the GLCF UMIACS www

Looking at
http://xserve.flids.com/pipermail/gdal-dev/2004-March/005095.html
"... cutting the image into pieces which I seperately reprojected using either a
  sinusoidal or a mollweide projection ..."

which contains the link to
http://gis.esri.com/library/userconf/proc98/PROCEED/TO850/PAP844/P844.HTM

From this I tried:

     /*
     /* Northern Eurasia
     /* REGION = "02 "
4, 30.0000000, 65.0000000
     -40.0000000, 40.7366111
     -40.0000000, 60.0000000
     -10.0000000, 60.0000000
     -10.0000000, 90.0000000
     180.0000000, 90.0000000
     180.0000000, 40.7366111
end

#####################################
dd.r02:

input
projection GEOGRAPHIC
units DD
parameters
output
projection MOLLWEIDE
units METERS
spheroid SPHERE
xshift 3335846.22854
yshift -336410.83237
parameters
30 00 00
end
#####################################

#note that we cannot cut out on the fly:
gdalwarp -s_srs '+proj=moll +ellps=sphere +lon_0=30E \
+y_0=-336410.83237 +x_0=3335846.22854' \
Goodes.EUAS.2000321.band4.tif d02_europe_EUAS.2000321.band4.moll.tif

#Since it's only defined for a subregion, extract the valid part
#(not completely true for part of greenland, I think):
gdalwarp -s_srs '+proj=moll +ellps=sphere +lon_0=30E \
+y_0=-336410.83237 +x_0=3335846.22854' \
-t_srs '+init=epsg:4326' -te -40 40 180 90 \
d02_europe_EUAS.2000321.band4.moll.tif d02_europe_EUAS.2000321.band4.LL.tif

... I didn't fully try, was running out of disk space.

gdalinfo d02_europe_EUAS.2000321.band4.LL.tif
qgis d02_europe_EUAS.2000321.band4.LL.tif
r.in.gdal ....

Maybe this gives you some ideas.
Please post a solution if you find anything.

Markus

Thanks Markus

I'll look into this in some spare time soon. Will let know for sure.

Maciek

----- Original Message ----- From: "Markus Neteler" <neteler@itc.it>
To: "Maciek Sieczka" <werchowyna@epf.pl>
Cc: <grasslist@baylor.edu>
Sent: Wednesday, April 06, 2005 6:17 PM
Subject: Re: [GRASSLIST:6332] Re: Goode in Grass

On Tue, Apr 05, 2005 at 10:13:58AM +0200, Maciek Sieczka wrote:

From: "Markus Neteler" <neteler@itc.it>
>On Sun, Apr 03, 2005 at 11:26:52AM +0200, Maciek Sieczka wrote:
>>Hi
>>
>>Could somebody please share the PROJ_INFO for Goode's projection? Or >>the
>>proj4 definition?
>>
>>I'd like to work on MODIS imagery of Poland, Eastern Europe. I couldn't
>>figure out the appropriate projection definition by myself.

>Most MODIS maps are in Sinusoidal (on MODIS sphere).
>Where did you get the map?

say ftp://ftp.glcf.umiacs.umd.edu/modis/500m/Eurasia/Goodes.EUAS.2000321/
via the GLCF UMIACS www

Looking at
http://xserve.flids.com/pipermail/gdal-dev/2004-March/005095.html
"... cutting the image into pieces which I seperately reprojected using either a
sinusoidal or a mollweide projection ..."

which contains the link to
http://gis.esri.com/library/userconf/proc98/PROCEED/TO850/PAP844/P844.HTM

From this I tried:

    /*
    /* Northern Eurasia
    /* REGION = "02 "
4, 30.0000000, 65.0000000
    -40.0000000, 40.7366111
    -40.0000000, 60.0000000
    -10.0000000, 60.0000000
    -10.0000000, 90.0000000
    180.0000000, 90.0000000
    180.0000000, 40.7366111
end

#####################################
dd.r02:

input
projection GEOGRAPHIC
units DD
parameters
output
projection MOLLWEIDE
units METERS
spheroid SPHERE
xshift 3335846.22854
yshift -336410.83237
parameters
30 00 00
end
#####################################

#note that we cannot cut out on the fly:
gdalwarp -s_srs '+proj=moll +ellps=sphere +lon_0=30E \
+y_0=-336410.83237 +x_0=3335846.22854' \
Goodes.EUAS.2000321.band4.tif d02_europe_EUAS.2000321.band4.moll.tif

#Since it's only defined for a subregion, extract the valid part
#(not completely true for part of greenland, I think):
gdalwarp -s_srs '+proj=moll +ellps=sphere +lon_0=30E \
+y_0=-336410.83237 +x_0=3335846.22854' \
-t_srs '+init=epsg:4326' -te -40 40 180 90 \
d02_europe_EUAS.2000321.band4.moll.tif d02_europe_EUAS.2000321.band4.LL.tif

... I didn't fully try, was running out of disk space.

gdalinfo d02_europe_EUAS.2000321.band4.LL.tif
qgis d02_europe_EUAS.2000321.band4.LL.tif
r.in.gdal ....

Maybe this gives you some ideas.
Please post a solution if you find anything.

Markus

Markus

Thank you for the links. I did some trial and error. Sorry it took so long. Please comment.

>>I'd like to work on MODIS imagery of Poland, Eastern Europe. I couldn't
>>figure out the appropriate projection definition by myself.

>Most MODIS maps are in Sinusoidal (on MODIS sphere).
>Where did you get the map?

say ftp://ftp.glcf.umiacs.umd.edu/modis/500m/Eurasia/Goodes.EUAS.2000321/
via the GLCF UMIACS www

Looking at
http://xserve.flids.com/pipermail/gdal-dev/2004-March/005095.html
"... cutting the image into pieces which I seperately reprojected using
either a
sinusoidal or a mollweide projection ..."

which contains the link to
http://gis.esri.com/library/userconf/proc98/PROCEED/TO850/PAP844/P844.HTM

From this I tried:

    /*
    /* Northern Eurasia
    /* REGION = "02 "
4, 30.0000000, 65.0000000
    -40.0000000, 40.7366111
    -40.0000000, 60.0000000
    -10.0000000, 60.0000000
    -10.0000000, 90.0000000
    180.0000000, 90.0000000
    180.0000000, 40.7366111
end

#####################################
dd.r02:

input
projection GEOGRAPHIC
units DD
parameters
output
projection MOLLWEIDE
units METERS
spheroid SPHERE
xshift 3335846.22854
yshift -336410.83237
parameters
30 00 00
end
#####################################

#note that we cannot cut out on the fly:
gdalwarp -s_srs '+proj=moll +ellps=sphere +lon_0=30E \
+y_0=-336410.83237 +x_0=3335846.22854' \
Goodes.EUAS.2000321.band4.tif d02_europe_EUAS.2000321.band4.moll.tif

This step seems redundant to me. Is there a reason you don't reproject to
wgs84 at once? Whether I apply the first step or not, the final wgs84
layer is the same.

#Since it's only defined for a subregion, extract the valid part
#(not completely true for part of greenland, I think):
gdalwarp -s_srs '+proj=moll +ellps=sphere +lon_0=30E \
+y_0=-336410.83237 +x_0=3335846.22854' \
-t_srs '+init=epsg:4326' -te -40 40 180 90 \
d02_europe_EUAS.2000321.band4.moll.tif
d02_europe_EUAS.2000321.band4.LL.tif

1. These settings you found give more-less good result. But, according to
that ESRI document, do you understand how the x,y shifts were derived?

2. As we move east or west from the central meridian in Mollweide's projection, there is an growing E distortion. How to calculate an optimal x_0 for Mollweide's, so the E distortion is minimalized in the area of my interest (a part of Poland)?

3. I noticed that +y_o=-316410.83237 (20000 more than in this ESRI document) gives a better match in the y axis (when the VMAP0 "COASTLINE" is overlayed, go figure). Is it me, the ESRI document or this MODIS raster in error ? Or maybe the VMAP0 (20 km?!). I don't have any other global data to compare to. But also when some wgs84 Landsat TM of Poland is overlayed on the MODIS reprojected to wgs84, they fit better when "my" +y_o=-316410.83237 is used.

4. There is Goode's projection in proj 4.4.9 (+proj=goode) which takes care of all the "zones" in some automated way. In case of the Eurasian MODIS we discuss here it works, as expected, similarly to +proj=moll, but:

a. It handles the N and S utmost parts of my MODIS image properly while Mollweide's introduces a vast distorion there.
b. y_0 is intended not to be provided in the proj definition for +proj=goode (I guess), but, following my observation of the 20km vertical shift, if y_0=20000 is given, the resulting wgs84 MODIS fits the VMAP0 data *much* better than when y_0 is omitted. Example:

gdalwarp -s_srs '+proj=goode +ellps=sphere +lon_0=30E +towgs84=0,0,0,0,0,0,0 +x_0=3335846.22854 +y_0=20000' -t_srs '+proj=latlong +ellps=wgs84 +towgs84=0,0,0,0,0,0,0' -te 7 33 34 59 Goodes.EUAS.2000321.band3.tif modis_goode_esri_y0+20000.tif

Yours badly puzzled,
Maciek

(attachments)

modis&vmap0_y_no_offs.jpg
modis&vmap0_y_20km_offs.jpg