[postgis-users] ST_Transform Z heights between ODN and WSG

[Paul Ramsey]

The problem *is* in fact at the proj level... you cannot do grid-shifted adjustments (like horiztonal and vertical datum shifts) if you do not have the grid shift files available. You cannot have the grid shift files available (in modern proj) without allowing network access to the grid CDN. To have access to the CDN you have to flip the PROJ_NETWORK environment to ON, so this:

  echo 545068   258591 8.51 | PROJ_NETWORK=ON cs2cs 'EPSG:7405' 'EPSG:4979'

works, while this:

  echo 545068 258591 8.51 | PROJ_NETWORK=OFF cs2cs 'EPSG:7405' 'EPSG:4979'

does not.
How to inject PROJ_NETWORK=ON? Well you can set it in your postgresql startup script environment. It's actually quite finicky much like handling the GDAL_ environment variables, actually. If you're on a DBaaS you may have no control over the value at all. And they may not want to turn it on, since in the limit you could cause GBs of grid shift files to be downloaded into the DB instance container.

P


> On Feb 1, 2023, at 8:51 AM, Paul Ramsey <pramsey at cleverelephant.ca> wrote:
> 
> Before then, if there's a proj wizard around, a command-line call that
> replicates the desired transform would be nice, to take proj out of
> the equation as a source of the problem.
> 
> On Wed, Feb 1, 2023 at 8:48 AM Paul Ramsey <pramsey at cleverelephant.ca> wrote:
>> 
>> I'm seeing something similar in my environment, even going between two
>> systems with explicit vertical datums.
>> 
>> pramsey=# SELECT ST_AsTexT(ST_Transform('SRID=7405;POINT(545068 258591
>> 8.51)'::geometry, 3901));
>>                      st_astext
>> -----------------------------------------------------
>> POINT Z (1681647.7511311213 6134585.208145529 8.51)
>> 
>> Time to put it in the debugger.
>> 
>> On Wed, Feb 1, 2023 at 4:13 AM James Howe <jmh205 at cam.ac.uk> wrote:
>>> 
>>> I can confirm it works in PostGIS 2.5, but not PostGIS 3.0+.
>>> I've been using the official Docker images to test.
>>> 
>>> Time to file a bug then?
>>> 
>>> James
>>> 
>>>> I had postgresql10-postgis2.5 running.
>>>> It didn't have 4979, so I inserted it.
>>>> -- ERROR: GetProj4StringSPI: Cannot find SRID (4979) in spatial_ref_sys
>>>> -- https://epsg.io/4979
>>>> INSERT into spatial_ref_sys (srid, auth_name, auth_srid, proj4text, srtext) values ( 4979, 'EPSG', > 4979, '+proj=longlat +datum=WGS84 +no_defs +type=crs', 'None');
>>>> SELECT ST_Transform('SRID=7405;POINT(545068 258591 8.51)'::geometry, 4979);
>>>> That rendered me:
>>>> st_transform                                                   |
>>>> ---------------------------------------------------------------+
>>>> POINT Z(0.1215563454880247 52.20645574758174 55.02789907436818)|
>>>> 
>>>> 
>>>> Met vriendelijke groet,
>>>> Wouter Scherphof
>>> _______________________________________________
>>> postgis-users mailing list
>>> postgis-users at lists.osgeo.org
>>> https://lists.osgeo.org/mailman/listinfo/postgis-users