Datum recognition problem Alber's Equal Area Conic Conformal

[Markus Neteler]

Frank,

(this is based on a GRASS user request, Steve in CC:

​http://grass.itc.it/pipermail/grassuser/2007-April/039354.html

"...the data come in UTM about 19.3km north of where they should be..." )

downloading the US MODIS NDVI data, it appears that the datum isn't read correctly:

Should be: ​ftp://ftp.glcf.umiacs.umd.edu/modis/NDVI/United_States/2000/US.Albers.2000081/USVI.UserGuide.txt " Projection: Alber's Equal Area Conic Conformal

• reprojected from native MODIS projection: Sinusoidal
• Central Meridian used was 96d W
• std Parallels of 29d 30' 0" and 45d 30' 0"
• Latitude of Origin is 0d 0' 0"
• Easting & Northing are 0.0000000000 m
• Datum: sphere (a sphere of radius 6370997 m)
• Resolution is 250 meters

"

Reported (recent GDAL SVN HEAD):

gdalinfo US.Albers.2000081.ndvi.tif
Driver: GTiff/GeoTIFF
Size is 21000, 13000
Coordinate System is:
PROJCS["unnamed",
    GEOGCS["Unknown datum based upon the Authalic Sphere",
        DATUM["Not_specified_based_on_Authalic_Sphere",
            SPHEROID["Sphere",6370997,0,
                AUTHORITY["EPSG","7035"]],
            AUTHORITY["EPSG","6035"]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433],
        AUTHORITY["EPSG","4035"]],
    PROJECTION["Albers_Conic_Equal_Area"],
    PARAMETER["standard_parallel_1",29.5],
    PARAMETER["standard_parallel_2",45.5],
    PARAMETER["latitude_of_center",0],
    PARAMETER["longitude_of_center",-96],
    PARAMETER["false_easting",0],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]]]
Origin = (-2356750.000000000000000,5769500.000000000000000)
Pixel Size = (250.000000000000000,-250.000000000000000)
Metadata:
  AREA_OR_POINT=Point
  TIFFTAG_SOFTWARE=hdf2geotiff v0.1.5
Corner Coordinates:
Upper Left  (-2356750.000, 5769500.000) (128d59'15.31"W, 49d50'21.28"N)
Lower Left  (-2356750.000, 2519500.000) (118d30'25.87"W, 21d47'47.16"N)
Upper Right ( 2893250.000, 5769500.000) ( 56d16'50.11"W, 48d 0'54.53"N)
Lower Right ( 2893250.000, 2519500.000) ( 68d37'26.12"W, 20d30'41.05"N)
Center      (  268250.000, 4144500.000) ( 92d52'8.97"W, 38d54'34.96"N)
Band 1 Block=21000x1 Type=Byte, ColorInterp=Gray

Dataset: ​ftp://ftp.glcf.umiacs.umd.edu/modis/NDVI/United_States/2000/US.Albers.2000081/ -> US.Albers.2000081.ndvi.tif.gz

I wonder if this could be fixed if there is anything to fix (probably PROJCS?). Myself, I didn't try to reproject to UTM. Steve, the reporter, is in CC.

Best regards, Markus

[Even Rouault]

I've reprojected it into EPSG:4326 and displayed on top of it VMAP1 coastlines. There's obviously something wrong with that image : I can also observe a major shift towards north.

Now, here's the associated proj4 parameters that are used by gdalinfo : +proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=0 +lon_0=-96 +x_0=0 +y_0=0 +a=6370997 +b=6370997 +units=m +no_defs. They look coherent with the gdalinfo and the information of USVI.UserGuide.txt.

This file contains also the following indication (with a correction on the extra '-' for the easting of the lower right corner. the technical guide ​http://glcf.umiacs.umd.edu/data/guide/technical/techguide_ndvi.pdf has not this typo) :

* Geographic Upper Left Corner: 128d 58' 17.49" W  49d 59' 50.83" N
* Geographic Lower Right Corner: 68d 38'  9.69" W  20d 36' 24.18" N

* Geocoded Upper Left Corner (m):  -2356875.000 E   5769625.000 N
* Geocoded Lower Right Corner (m): 2893125.000 E   2519625.000 N

The geocoded corners are coherent with the gdalinfo output, but the geographic coordinates are not !

Now when I do :

cs2cs +proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=0 +lon_0=-96 +x_0=0 +y_0=0 +a=6370997 +b=6370997 +units=m +no_defs

with

-2356875.000 5769625.000

I get : 128d59'23.212"W 49d50'23.75"N 0.000

So it's either a bug in PROJ4 or something wrong with the georeferencing of the GeoTiff and the associated text file.

[Even Rouault]

The metadata text file is very confusing. The geographic coordinates that are mentionned for the corners (ie upper left corner : 128d 58' 17.49" W 49d 59' 50.83" N) are expressed on ellipsoid WGS84, not on the 6370997 sphere !!!

cs2cs +proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=0 +lon_0=-96 +x_0=0 +y_0=0 +ellps=WGS84 +units=m +no_defs

with

-2356875.000 5769625.000

gives

128d58'17.486"W 49d59'50.827"N 0.000

So, I've tried

gdalwarp -s_srs '+proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=0 +lon_0=-96 +x_0=0 +y_0=0 +ellps=WGS84 +units=m +no_defs' -t_srs EPSG:4326 US.Albers.2000081.ndvi.tif output.tif

The shift towards north has a bit decreased, but it's still there...

[Even Rouault]

Apparently, working with this data source is a nightmare for many people. Just look at : ​http://grass.itc.it/pipermail/grassuser/2005-May/028996.html.

I'm thinking more and more that there was something wrong when they reprojected from the original SRS to Albers's Equal Area Conic Conformal

[warmerdam]

I'm going to close this based on the assumption that the data is wrong, not GDAL's reprojection support. Please reopen if you find there really is a GDAL problem, and can provide some support.