#6523 closed defect (fixed)
Incorrect parameters read for Mercator projection in adf-style prj file
| Reported by: | Mike Taves | Owned by: | warmerdam |
|---|---|---|---|
| Priority: | normal | Milestone: | 1.11.5 |
| Component: | default | Version: | unspecified |
| Severity: | normal | Keywords: | |
| Cc: |
Description
Reading "NZ 250m ESRI ASCII grid" from this website (563 MB) with QGIS shows the raster in the wrong location. The .prj file appears to be a copy of the prj.adf file from the "NZ 250m ESRI binary grid" download option, and looks like this:
Projection MERCATOR Datum WGS84 Spheroid WGS84 Units METERS Zunits NO Xshift 0.0 Yshift 0.0 Parameters 100 0 0.0 /* longitude of central meridian -41 0 0.0 /* latitude of true scale 0.0 /* false easting (meters) 0.0 /* false northing (meters)
Further information for the dataset describes the projection as EPSG:3994. However, with GDAL 2.0.2:
$ gdalinfo nzbathymetry_2016_ascii-grid.txt
Driver: AAIGrid/Arc/Info ASCII Grid
Files: nzbathymetry_2016_ascii-grid.txt
nzbathymetry_2016_ascii-grid.prj
Size is 12116, 15391
Coordinate System is:
PROJCS["unnamed",
GEOGCS["WGS 84",
DATUM["WGS_1984",
SPHEROID["WGS 84",6378137,298.257223563,
AUTHORITY["EPSG","7030"]],
TOWGS84[0,0,0,0,0,0,0],
AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich",0,
AUTHORITY["EPSG","8901"]],
UNIT["degree",0.0174532925199433,
AUTHORITY["EPSG","9108"]],
AUTHORITY["EPSG","4326"]],
PROJECTION["Mercator_1SP"],
PARAMETER["latitude_of_origin",100],
PARAMETER["central_meridian",0],
PARAMETER["scale_factor",1],
PARAMETER["false_easting",-41],
PARAMETER["false_northing",0],
UNIT["METERS",1]]
Origin = (4795705.599999999600000,-2067835.799999999800000)
Pixel Size = (250.000000000000000,-250.000000000000000)
Corner Coordinates:
Upper Left ( 4795705.600,-2067835.800)
Lower Left ( 4795705.600,-5915585.800)
Upper Right ( 7824705.600,-2067835.800)
Lower Right ( 7824705.600,-5915585.800)
Center ( 6310205.600,-3991710.800)
Band 1 Block=12116x1 Type=Float32, ColorInterp=Undefined
NoData Value=-9999
Note the discrepancies for latitude_of_origin (latitude 100?), central_meridian, and false_easting. This raster displays correctly in ArcGIS 10.2.
Change History (9)
comment:1 by , 8 years ago
follow-up: 3 comment:2 by , 8 years ago
Could you prepare a small test data with prj.adf and some rows of ascii grid data and attach them to the ticket? I was about to have a look but downloading the 563 MB zip from New Zealand to Finland would have taken 4 hours.
comment:3 by , 8 years ago
Replying to jratike80:
Could you prepare a small test data with prj.adf and some rows of ascii grid data and attach them to the ticket? I was about to have a look but downloading the 563 MB zip from New Zealand to Finland would have taken 4 hours.
The reduced version of this issue was demonstrated in 1, so you don't need to download anything.
comment:8 by , 8 years ago
| Milestone: | → 1.11.5 |
|---|
(As far as I can see, the issue existed since support for MERCATOR was introduced in r19742)
To confirm, the same erroneous happens with the ESRI binary grid file, which uses
prj.adf. From ArcCatalog, here is the Coordinate System info:WGS_1984_Mercator Authority: Custom Projection: Mercator False_Easting: 0.0 False_Northing: 0.0 Central_Meridian: 100.0 Standard_Parallel_1: -41.0 Linear Unit: Meter (1.0) Geographic Coordinate System: GCS_WGS_1984 Angular Unit: Degree (0.0174532925199433) Prime Meridian: Greenwich (0.0) Datum: D_WGS_1984 Spheroid: WGS_1984 Semimajor Axis: 6378137.0 Semiminor Axis: 6356752.314245179 Inverse Flattening: 298.257223563Looking at gdal/ogr/ogr_srs_esri.cpp ~L818, the order of PARAM_0 to PARAM_3 for
SetMercatorlooks correct (as documented by Esri for Mercator > Parameters > Workstation), andSetMercatoris behaving normally with the OSR Python bindings.Here is a reduced version of the issue with OSR through the Python bindings:
shows
2000200 sr1 PROJ.4: +proj=merc +lon_0=100 +lat_ts=-41 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs sr1 WTK: PROJCS["WGS 84 / Mercator 41", GEOGCS["WGS 84", DATUM["WGS_1984", SPHEROID["WGS 84",6378137,298.257223563, AUTHORITY["EPSG","7030"]], AUTHORITY["EPSG","6326"]], PRIMEM["Greenwich",0, AUTHORITY["EPSG","8901"]], UNIT["degree",0.0174532925199433, AUTHORITY["EPSG","9122"]], AUTHORITY["EPSG","4326"]], PROJECTION["Mercator_2SP"], PARAMETER["standard_parallel_1",-41], PARAMETER["central_meridian",100], PARAMETER["false_easting",0], PARAMETER["false_northing",0], UNIT["metre",1, AUTHORITY["EPSG","9001"]], AXIS["X",EAST], AXIS["Y",NORTH], AUTHORITY["EPSG","3994"]] sr2 PROJ.4: +proj=merc +lon_0=100 +lat_ts=-41 +x_0=0 +y_0=0 +ellps=WGS84 +units=m +no_defs sr2 WTK: PROJCS["unnamed", PROJECTION["Mercator_1SP"], PARAMETER["latitude_of_origin",-41], PARAMETER["central_meridian",100], PARAMETER["scale_factor",1], PARAMETER["false_easting",0], PARAMETER["false_northing",0]] sr3 PROJ.4: +proj=merc +lon_0=0 +lat_ts=100 +x_0=-41 +y_0=0 +ellps=WGS84 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs sr3 WTK: PROJCS["unnamed", GEOGCS["WGS 84", DATUM["WGS_1984", SPHEROID["WGS 84",6378137,298.257223563, AUTHORITY["EPSG","7030"]], TOWGS84[0,0,0,0,0,0,0], AUTHORITY["EPSG","6326"]], PRIMEM["Greenwich",0, AUTHORITY["EPSG","8901"]], UNIT["degree",0.0174532925199433, AUTHORITY["EPSG","9108"]], AUTHORITY["EPSG","4326"]], PROJECTION["Mercator_1SP"], PARAMETER["latitude_of_origin",100], PARAMETER["central_meridian",0], PARAMETER["scale_factor",1], PARAMETER["false_easting",-41], PARAMETER["false_northing",0], UNIT["METERS",1]]