This is an additional exploration of the problem described at Reprojecting MODIS Ice Sea Temperature products with gdalwarp which was never answered in a satisfying way (I apologise for the length of the post).
I want to reproject Antarctic MODIS Sea Ice Extent and IST Daily L3 Global 4km EASE-Grid Day data to Antarctic Polar Stereographic using GDAL. I am interested in both AQUA (product MYD29E1D) and TERRA (product MOD29E1D) and I am interested in the full time series available.
Processing of version 6 is currently ongoing and should be completed by the end of 2015 according to https://nsidc.org/data/MYD29E1D and http://nsidc.org/data/MOD29E1D. Therefore I try to use version 5 that is made available via FTP from ftp://n5eil01u.ecs.nsidc.org/SAN/MOSA/MYD29E1D.005/ and ftp://n5eil01u.ecs.nsidc.org/SAN/MOST/MOD29E1D.005/.
These MODIS products use the rather esoteric EASE-Grid projection which exists in several versions described here: http://nsidc.org/data/ease/versions.html. Apparently there is an orignal version of EASE-Grid that was released in 1992 and a version 2.0 released in 2011. The website sited above states that six EPSG codes are available to handle EASE-Grid. For each of the two version there exist three EPSG codes relating to "North", "South" and "Global".
OriginalNorth: 3408, South: 3409, Global 3410
V2.0 North: 6931, South: 6932, Global: 6933
Which adds to the complexity is that there are two deprecated EASE-GRID EPSG codes
Deprecated North: 3973, South: 3974, Global: 3975
The question that I had now were:
gdalsrsinfo 'HDF4_EOS:EOS_GRID:"MOD29E1D.A2009345.005.2009346083234.hdf":MOD_Grid_Seaice_4km_South:Ice_Surface_Temperature_SP'PROJ.4 : '+proj=laea +lat_0=-5156620156.177409 +lon_0=0 +x_0=0 +y_0=0 +ellps=clrk66 +units=m +no_defs 'OGC WKT
ROJCS["unnamed", GEOGCS["Unknown datum based upon the Clarke 1866 ellipsoid", DATUM["Not specified (based on Clarke 1866 spheroid)", SPHEROID["Clarke 1866",6378206.4,294.9786982139006, AUTHORITY["EPSG","7008"]]], PRIMEM["Greenwich",0], UNIT["degree",0.0174532925199433]], PROJECTION["Lambert_Azimuthal_Equal_Area"], PARAMETER["latitude_of_center",-5156620156.177409], PARAMETER["longitude_of_center",0], PARAMETER["false_easting",0], PARAMETER["false_northing",0], UNIT["Meter",1]]As I could not find any conclusive documentation, I tried to convert a dataset from 2009 and a dataset from 2014 using all available "Global" and "South" codes (3409, 3410, 6932, 6933, 3974, 3975). The two datasets I used are ftp://n5eil01u.ecs.nsidc.org/SAN/MOST/MOD29E1D.005/2009.12.11/ and ftp://n5eil01u.ecs.nsidc.org/SAN/MOST/MOD29E1D.005/2014.12.11/.
Regarding the software versions I use:
gdalinfo --versionGDAL 1.11.2, released 2015/02/10projRel. 4.9.0, 27 October 2013I first extract the relevant band (MODIS product) and then convert it to Antarctic Polar Stereographic e.g.
gdal_translate -unscale -srcwin 500 500 3500 3500 -a_srs epsg:3409 HDF4_EOS:EOS_GRID:"MOD29E1D.A2009345.005.2009346083234.hdf":MOD_Grid_Seaice_4km_South:Ice_Surface_Temperature_SP MOD29E1D.A2009345.005.2009346083234_tmp.tif(The -srcwin is necessary to get rid of reprojection errors resulting from values outside the projection area)
followed by e.g.
gdalwarp -overwrite -s_srs EPSG:3409 -t_srs EPSG:3031 MOD29E1D.A2009345.005.2009346083234_tmp.tif MOD29E1D.A2009345.005.2009346083234_3409.tif(Side remark: if I use the two commands without specifying the source CRS, the results are unusable)
The first problem that I encountered is that 6932 and 6933 are not EPSG-codes supported by the GDAL/PROJ version that I use. The errors that I get are:
ERROR 6: EPSG PCS/GCS code 6932 not found in EPSG support files. Is this a valid EPSG coordinate system?ERROR 6: EPSG PCS/GCS code 6933 not found in EPSG support files. Is this a valid EPSG coordinate system?Having a look at http://spatialreference.org/ref/sr-org/6932/ it seems that 6932 indeed is not a registered EPSG code (although stated otherwise on the EASE-GRID website). I therefore switched to using the proj-definition of 6932 and 6933:
+proj=tmerc +lat_0=0 +lon_0=-111 +k=0.9996 +x_0=72449 +y_0=-6259003 +ellps=GRS80 +units=m +no_defsBut if I use this definition I get the following error:
ERROR 1: latitude or longitude exceeded limitsERROR 1: Reprojection failed, err = -14, further errors will be supressed on the transform object.Coordinate reference systems 6932 and 6933 are therefore not usable (for me).
I therefore applied transformatons to the four remaining coordinate reference systems (3409, 3410, 3974, 3975) to the two example files. Of he eight resulting files, only one created an output with a reasonable fit: the "old" EPSG:3409 or the deprecated EPSG:3974 when applied to the 2014 dataset.
So to summarise:
أكثر...
I want to reproject Antarctic MODIS Sea Ice Extent and IST Daily L3 Global 4km EASE-Grid Day data to Antarctic Polar Stereographic using GDAL. I am interested in both AQUA (product MYD29E1D) and TERRA (product MOD29E1D) and I am interested in the full time series available.
Processing of version 6 is currently ongoing and should be completed by the end of 2015 according to https://nsidc.org/data/MYD29E1D and http://nsidc.org/data/MOD29E1D. Therefore I try to use version 5 that is made available via FTP from ftp://n5eil01u.ecs.nsidc.org/SAN/MOSA/MYD29E1D.005/ and ftp://n5eil01u.ecs.nsidc.org/SAN/MOST/MOD29E1D.005/.
These MODIS products use the rather esoteric EASE-Grid projection which exists in several versions described here: http://nsidc.org/data/ease/versions.html. Apparently there is an orignal version of EASE-Grid that was released in 1992 and a version 2.0 released in 2011. The website sited above states that six EPSG codes are available to handle EASE-Grid. For each of the two version there exist three EPSG codes relating to "North", "South" and "Global".
OriginalNorth: 3408, South: 3409, Global 3410
V2.0 North: 6931, South: 6932, Global: 6933
Which adds to the complexity is that there are two deprecated EASE-GRID EPSG codes
Deprecated North: 3973, South: 3974, Global: 3975
The question that I had now were:
- Which EPSG-code should I use in the first place?
- If there are two versions of EASE-GRID, does that mean I have to use a different EPSG-code/CRS for data created pre-2011?
gdalsrsinfo 'HDF4_EOS:EOS_GRID:"MOD29E1D.A2009345.005.2009346083234.hdf":MOD_Grid_Seaice_4km_South:Ice_Surface_Temperature_SP'PROJ.4 : '+proj=laea +lat_0=-5156620156.177409 +lon_0=0 +x_0=0 +y_0=0 +ellps=clrk66 +units=m +no_defs 'OGC WKT
ROJCS["unnamed", GEOGCS["Unknown datum based upon the Clarke 1866 ellipsoid", DATUM["Not specified (based on Clarke 1866 spheroid)", SPHEROID["Clarke 1866",6378206.4,294.9786982139006, AUTHORITY["EPSG","7008"]]], PRIMEM["Greenwich",0], UNIT["degree",0.0174532925199433]], PROJECTION["Lambert_Azimuthal_Equal_Area"], PARAMETER["latitude_of_center",-5156620156.177409], PARAMETER["longitude_of_center",0], PARAMETER["false_easting",0], PARAMETER["false_northing",0], UNIT["Meter",1]]As I could not find any conclusive documentation, I tried to convert a dataset from 2009 and a dataset from 2014 using all available "Global" and "South" codes (3409, 3410, 6932, 6933, 3974, 3975). The two datasets I used are ftp://n5eil01u.ecs.nsidc.org/SAN/MOST/MOD29E1D.005/2009.12.11/ and ftp://n5eil01u.ecs.nsidc.org/SAN/MOST/MOD29E1D.005/2014.12.11/.Regarding the software versions I use:
gdalinfo --versionGDAL 1.11.2, released 2015/02/10projRel. 4.9.0, 27 October 2013I first extract the relevant band (MODIS product) and then convert it to Antarctic Polar Stereographic e.g.
gdal_translate -unscale -srcwin 500 500 3500 3500 -a_srs epsg:3409 HDF4_EOS:EOS_GRID:"MOD29E1D.A2009345.005.2009346083234.hdf":MOD_Grid_Seaice_4km_South:Ice_Surface_Temperature_SP MOD29E1D.A2009345.005.2009346083234_tmp.tif(The -srcwin is necessary to get rid of reprojection errors resulting from values outside the projection area)
followed by e.g.
gdalwarp -overwrite -s_srs EPSG:3409 -t_srs EPSG:3031 MOD29E1D.A2009345.005.2009346083234_tmp.tif MOD29E1D.A2009345.005.2009346083234_3409.tif(Side remark: if I use the two commands without specifying the source CRS, the results are unusable)
The first problem that I encountered is that 6932 and 6933 are not EPSG-codes supported by the GDAL/PROJ version that I use. The errors that I get are:
ERROR 6: EPSG PCS/GCS code 6932 not found in EPSG support files. Is this a valid EPSG coordinate system?ERROR 6: EPSG PCS/GCS code 6933 not found in EPSG support files. Is this a valid EPSG coordinate system?Having a look at http://spatialreference.org/ref/sr-org/6932/ it seems that 6932 indeed is not a registered EPSG code (although stated otherwise on the EASE-GRID website). I therefore switched to using the proj-definition of 6932 and 6933:
+proj=tmerc +lat_0=0 +lon_0=-111 +k=0.9996 +x_0=72449 +y_0=-6259003 +ellps=GRS80 +units=m +no_defsBut if I use this definition I get the following error:
ERROR 1: latitude or longitude exceeded limitsERROR 1: Reprojection failed, err = -14, further errors will be supressed on the transform object.Coordinate reference systems 6932 and 6933 are therefore not usable (for me).
I therefore applied transformatons to the four remaining coordinate reference systems (3409, 3410, 3974, 3975) to the two example files. Of he eight resulting files, only one created an output with a reasonable fit: the "old" EPSG:3409 or the deprecated EPSG:3974 when applied to the 2014 dataset.
So to summarise:
- I could not find any way to reproject pre-2011 MODIS IST data
- I was not able to use the recommended CRS EPSG:6932/EPSG:6933
أكثر...