On Mon, Jul 20, 2020 at 11:25 AM Shira Bezalel <[email protected]> wrote:
> We recently recompiled PostGIS in one of our Postgres 12.3 environments > with PROJ 6 (6.3.2) from PROJ 4 (4.9.3). > > postgis_full_version() output: > POSTGIS="3.0.1 ec2a9aa" [EXTENSION] PGSQL="120" GEOS="3.7.1-CAPI-1.11.1 > 27a5e771" PROJ="6.3.2" LIBXML="2.9.4" LIBJSON="0.12.1" > > We're now noticing that using ST_Transform(geometry, srid) to transform > some points in California from SRID 4326 to 900913 is causing a shift along > the Y axis of about 20-30 meters. The difference is evidenced by a > comparison of the old and new ST_Y() values. > > If I force the PROJ 4 definition to be used, the transformation works fine > with no shift. (To force the PROJ 4 definition, I used the method described > in http://blog.cleverelephant.ca/2019/02/proj4-postgis.html and set > authname, authsrid and srtext to null so that the proj4text column is used.) > > Is this shift expected with PROJ 6? I realize 900913 is deprecated and we > plan to migrate to 3857 as a longer term solution, but looking for > a shorter term fix if possible since much data and code is impacted by > this. Note we are seeing no shift issues with SRID 3857. > > Let me know if I can provide more information. > > Thanks much, > Shira > > Hello again. Wanted to follow up on this. We have found a workaround that seems somewhat better than just setting the auth_name, auth_srid and srtext columns to null especially since external applications rely on the srtext column. If we replace the record for srid 900913 in the spatial_ref_sys table with the one at epsg.io (https://epsg.io/900913) then the points are correctly transforming without a shift with st_transform(geometry,900913). As a related consequence of this, I am noticing a shift when exporting these 900913 points to shapefile format with pgsql2shp (it's fine with the old srtext definition), but we don't use pgsql2shp much so I think that's okay. Our long-term plan is to move away from 900913, but this might be suitable for the short-term. If I'm missing any pitfalls with this approach (apart from the pgsql2shp one), please feel free to let me know. New srtext value from epsg.io: PROJCS["Google Maps Global Mercator",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.01745329251994328,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Mercator_2SP"],PARAMETER["standard_parallel_1",0],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",0],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["Meter",1],EXTENSION["PROJ4","+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +wktext +no_defs"],AUTHORITY["EPSG","900913"]] Original srtext value in spatial_ref_sys: PROJCS["Popular Visualisation CRS / Mercator (deprecated)",GEOGCS["Popular Visualisation CRS",DATUM["Popular_Visualisation_Datum",SPHEROID["Popular Visualisation Sphere",6378137,0,AUTHORITY["EPSG","7059"]],TOWGS84[0,0,0,0,0,0,0],AUTHORITY["EPSG","6055"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.01745329251994328,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4055"]],UNIT["metre",1,AUTHORITY["EPSG","9001"]],PROJECTION["Mercator_1SP"],PARAMETER["central_meridian",0],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],AUTHORITY["EPSG","3785"],AXIS["X",EAST],AXIS["Y",NORTH]] Thanks, Shira
_______________________________________________ postgis-users mailing list [email protected] https://lists.osgeo.org/mailman/listinfo/postgis-users
