[PROJ] Pseudo-mercator

[Andrew Patterson]

Even,


> That one is a very odd beast. The
> PARAMETER["latitude_of_center",6.64456744726493] is completely
> inconsistant with the name "WGS_84_Pseudo_Mercator". I see GDAL 2.4
> honoured this latitude of center, as gdalsrsinfo on it returns
> +proj=merc +lon_0=0 +lat_ts=6.64456744726493 +x_0=0 +y_0=0 +ellps=WGS84
> +towgs84=0,0,0,0,0,0,0 +units=m +no_defs.    Newer PROJ versions have a
> heustistics that only look at the name WGS_84_Pseudo_Mercator and ignore
> the projection parameters to decide to use the official EPSG:3857 "WGS
> 84 / Pseudo-Mercator" definition. One could argue it is a bug to ignore
> the latitude_of_center parameter, but this definition is very confusing
> to start with. It is hard to know what the intention of the user was:
> use official EPSG:3857 "WGS 84 / Pseudo-Mercator", or use a modified
> version of it...
>

I should say up front that I have nothing invested in any of these
personally -- if they're wrong, that's an answer I can at least relay back.
It'll be ugly in some cases because we'll have umpteen hundred maps using a
questionable WKT that has to be remedied somehow, but at least I can say
there's no code solution at my end.

That said, to play devil's advocate, what would the right way to specify a
customized pseudo-mercator look like? It sounds like you're saying that the
PROJCS name is the problem? I tried changing the name to 'Unknown' and it
had no effect, so either you're talking about a PROJ version later than
8.1.1 or I'm misunderstanding what you mean.

Of the two cases, this is the one that I feel like there should be a way to
make it work properly. The next map is sounding more & more like a lousy
map to start with. This one is weird, but it should be feasible (if
strange) to tweak a pseudo-mercator like this, shouldn't it?


> As far as I can see, GDAL 2.4 can't transform that to a PROJ.4 string at
> all (it doesn't understand the "Mercator" projection), so I'm not sure
> why you get a result at all with it. With recent PROJ,
>
> $ projinfo
> 'PROJCS["Popular_Visualisation_CRS_Mercator_deprecated",GEOGCS["GCS_Popular
>
> Visualisation
>
> CRS",DATUM["D_Popular_Visualisation_Datum",SPHEROID["Popular_Visualisation_Sphere",6378137,0]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Mercator"],PARAMETER["central_meridian",0],PARAMETER["false_easting",0],PARAMETER["false_northing",0],PARAMETER["standard_parallel_1",0],UNIT["Meter",1]]'
>
> returns
>
> PROJ.4 string:
> +proj=merc +lat_ts=0 +lon_0=0 +x_0=0 +y_0=0 +R=6378137 +units=m +no_defs
> +type=crs
>
> WKT2:2019 string:
> PROJCRS["Popular_Visualisation_CRS_Mercator_deprecated",
>      BASEGEOGCRS["GCS_Popular Visualisation CRS",
>          DATUM["Popular Visualisation Datum",
>              ELLIPSOID["Popular_Visualisation_Sphere",6378137,0,
>                  LENGTHUNIT["metre",1]],
>              ID["EPSG",6055]],
>          PRIMEM["Greenwich",0,
>              ANGLEUNIT["Degree",0.0174532925199433]]],
>      CONVERSION["unnamed",
>          METHOD["Mercator (variant B)",
>              ID["EPSG",9805]],
>          PARAMETER["Latitude of 1st standard parallel",0,
>              ANGLEUNIT["Degree",0.0174532925199433],
>              ID["EPSG",8823]],
>          PARAMETER["Longitude of natural origin",0,
>              ANGLEUNIT["Degree",0.0174532925199433],
>              ID["EPSG",8802]],
>          PARAMETER["False easting",0,
>              LENGTHUNIT["metre",1],
>              ID["EPSG",8806]],
>          PARAMETER["False northing",0,
>              LENGTHUNIT["metre",1],
>              ID["EPSG",8807]]],
>      CS[Cartesian,2],
>          AXIS["(E)",east,
>              ORDER[1],
>              LENGTHUNIT["metre",1,
>                  ID["EPSG",9001]]],
>          AXIS["(N)",north,
>              ORDER[2],
>              LENGTHUNIT["metre",1,
>                  ID["EPSG",9001]]]]
>
> Which is a reasonable enough interpretation of this WKT. It could have
> potentially been better recognized as the deprecated EPSG:3785 (which
> PROJ then interprets as EPSG;3857), but there's no heuristics for that
> particular formulation, so it is recognized as classic Mercator on a
> sphere.  In practice, it should give identical results to official
> EPSG:3857 "WGS 84 / Pseudo-Mercator", except when a datum transformation
> is involved. So if you transform between WGS 84 and this, you should get
> correct results. But results off by 26 km are beyond an ignored datum
> shift (differences should be off by a few hundred of meters max), so I'm
> not sure why you observe this.
>

That actually explains something one of our QA told me. He said there was
an unnamed conversion in there but I wasn't seeing -- I guess he was using
WKT2! Our QA was already converging on the consensus that this particular
map was just bad, so it sounds like your reckoning backs that up. The
strangest part is that it worked in PROJ 6.1 (though I'm told by a
colleague that we were running it in legacy mode, so perhaps it was more
akin to PROJ 4.x).


> Yes that one is the canonical WKT1:GDAL formulation
>

I asked our QA to go through the working maps and so far they've all had
the AUTHORITY["EPSG","3857"] tag. I suspect all of the workings are just
plain-jane vanilla pseudo-mercator. The problem ones all seem to be ones
that wander off the beaten path. The unfortunate thing is that it seems to
be somewhat common in our map store.

I appreciate the attention. The bottom line here seems to be that there was
some latitude (no pun intended) for these questionable versions of
pseudo-mercator that is gone in later versions of PROJ. Ultimately better
for everyone, but we'll have to figure out how to deal with those
problematic maps.

..............................
Andrew Patterson
Lead Software Architect
Avenza Systems Inc.

email: andrew at avenza.com
phone: 416.487.5116
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/proj/attachments/20211201/e8b1c395/attachment.html>