[GRASS-dev] PROJ6 support in GRASS

[Even Rouault]

On jeudi 7 novembre 2019 19:50:33 CET Markus Metz wrote:
> On Thu, Nov 7, 2019 at 6:47 PM Martin Landa <landa.martin at gmail.com> wrote:
> > Hi,
> > 
> > čt 7. 11. 2019 v 15:17 odesílatel Markus Metz
> > 
> > <markus.metz.giswork at gmail.com> napsal:
> > > PROJ 6 support is not yet finished in GRASS. Currently it sort of works
> 
> in master and relbr78, but I expect troubles because GRASS is heavily
> relying on deprecated proj strings as CRS definitions. This is becoming
> problematic e.g. when creating locations from GDAL/OGR datasets because
> GRASS goes from OGR spatialreference through proj string to GRASS
> definition.
> 
> > OSGeo4W packages are planned to be switched to GDAL3/PROJ6 soon, are
> > you aware of any known bugs?
> 
> yes, EPSG:3857 is causing trouble. This is the WGS 84 / Pseudo-Mercator
> projection on a sphere with radius 6378137 that needs special treatment
> which no longer works with GRASS + PROJ 6 + GDAL 3.
> 
> g.proj epsg=3857 -p gives me
> -PROJ_INFO-------------------------------------------------
> name       : WGS 84 / Pseudo-Mercator
> datum      : wgs84
> ellps      : wgs84
> proj       : merc
> lat_ts     : 0
> lon_0      : 0
> x_0        : 0
> y_0        : 0
> k          : 1
> wktext     : defined
> no_defs    : defined
> -PROJ_EPSG-------------------------------------------------
> epsg       : 3857
> -PROJ_UNITS------------------------------------------------
> unit       : meter
> units      : meters
> meters     : 1
> 
> which is wrong because it should be a Mercator projection on a sphere, not
> an ellipsoid

Well, this is a messy topic. EPSG:3857 uses the WebMercator projection, which 
is not the standard Mercator projection since it uses spherical formulas even 
if the datum and ellipsoid is WGS84

So neither  proj=merc + ellps=wgs84 or proj=merc + ellps=sphere are 
technically correct

Nothing new with PROJ 6 for that case

> projinfo -o PROJ EPSG:3857
> gives me
> PROJ.4 string:
> +proj=merc +a=6378137 +b=6378137 +lat_ts=0 +lon_0=0 +x_0=0 +y_0=0 +k=1
> +units=m +nadgrids=@null +wktext +no_defs +type=crs

Yes, this is a hack of PROJ.4 still supported by PROJ 6 to more or less get 
the correct result of WebMercator using the standard Mercator, by lying on the 
ellipsoid definition (using the Sphere) and applying a nadgrids=null, so that 
no datum transformation between WGS84 and that Sphere is attempted !
> 
> that is correct, but unfortunately not returned by GDAL's function
> OSRExportToProj4() which is now deprecated.

Hum, I do ot agree with this. This *is* still what is returned by GDAL 
OSRExportToProj4() with GDAL + PROJ master :

$ gdalsrsinfo EPSG:3857

PROJ.4 : +proj=merc +a=6378137 +b=6378137 +lat_ts=0 +lon_0=0 +x_0=0 +y_0=0 
+k=1 +units=m +nadgrids=@null +wktext +no_defs

And you also actually use (PROJ 6 required here)

$ gdalsrsinfo "+proj=webmerc +datum=WGS84"
(or $ projinfo "+proj=webmerc +datum=WGS84 +type=crs" -o PROJ )

PROJ.4 : +proj=merc +a=6378137 +b=6378137 +lat_ts=0 +lon_0=0 +x_0=0 +y_0=0 
+k=1 +units=m +nadgrids=@null +wktext +no_defs


> That means GRASS needs to avoid
> OSRExportToProj4() with GDAL 3+:

In the general case, I agree with your conclusion, but not for the above one 
:-)

-- 
Spatialys - Geospatial professional services
http://www.spatialys.com