[Proj] EPSG:3994 problem

[pcreso at pcreso.com]

Hi,

I'm unable to get Proj to correctly return lat/long coordinates from Mercator41 data.
http://www.spatialreference.org/ref/epsg/3994/ does not list any proj4 string, Postgis proj4text or Mapserver values for this projection.

I have an ESRI binary grid with a projection definition (prj.adf) of:

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)

gdalinfo fails to interpret this correctly, giving:

gdalinfo w001001.adf 
Driver: AIG/Arc/Info Binary Grid
Files: .
       ./#m41_to_ll#
       ./sta.adf
       ./w001001x.adf
       ./hdr.adf
       ./w001001.adf
       ./m41_to_ll
       ./whbr
       ./prj.adf
       ./dblbnd.adf
Size is 8817, 7901
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 = (6199987.499999998137355,-3709987.500002769753337)
Pixel Size = (25.000000000000000,-25.000000000000000)
Corner Coordinates:
Upper Left  ( 6199987.500,-3709987.500) 
Lower Left  ( 6199987.500,-3907512.500) 
Upper Right ( 6420412.500,-3709987.500) 
Lower Right ( 6420412.500,-3907512.500) 
Center      ( 6310200.000,-3808750.000) 
Band 1 Block=256x4 Type=Float32, ColorInterp=Undefined
  Min=-2570.453 Max=2917.197 
  NoData Value=-3.4028234663852886e+38

With some obvious errors, like latitude of origin = 100... so just extracting the xyz values & trying to reproject with proj, I'm having problems.

I tried: g dal2xyz.py w001001.adf | head 
which returns:
6200000.000 -3710000.000 76.1574
6200025.000 -3710000.000 76.1842
...

which seems to work fine, but when I pass these Mercator41 coordinates into proj I'm unable to get a parameter list which gives the lat/long values for where I know these data sit. There is an offset of many miles. This suggests either my source data has an error in the coordinates (unlikely) or I don't understand enough about proj parameters (more likely)

Based on the contents of the above prj.adf file I have tried the following command:

echo "6200000.000 -3710000.000" | proj -If "%.5f" +proj=merc +lon_0=100 +lat_ts=-41 +k=1 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs +to +proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs
173.69093       -40.48353

I expected values around 174.66  -41.41 (the data is a 25m multibeam bathymetry grid around Cook Strait, between the North & South Islands of New Zealand) 


Any advice appreciated.

  Brent Wood



 

-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://lists.maptools.org/pipermail/proj/attachments/20120415/b444a446/attachment.htm