Vertical grid shift with PROJ 9

Question

I want to perform a vertical grid shift using PROJ 9.I have never used PROJ syntax and I am having some issues to have it right.

I have checked documentation and this previous post from which I got a sample syntax:

echo lat lng z | proj +proj=vgridshift +grids=egm08_25.gtx

For my case, I wan to transform (long, lat, ellipsoidal heigh to X&Y in EPSG25831,Z in epsg 3855). So following the recommendation I got this syntax:

echo 41.75758799 2.32152788  1235.577 | proj +proj=vgridshift +grids=egm08_25.gtx

But it produces an error:

C:\Program Files\QGIS 3.26.2>echo 2.32152788 41.75758799 1235.577 | proj +proj=vgridshift +grids=egm08_25.gtx
Rel. 9.0.1, June 15th, 2022
<proj>:
can't initialize operations that produce angular output coordinates
program abnormally terminated

I perfomed a test with a 2 D coordinates and it works just fine:

echo 2.32152788 41.75758799 | cs2cs +init=epsg:4326 +to +init=epsg:25831 

[

I also tried to follow the example code (from proj vgrid shift link) and I get an error:

+proj=vgridshift +grids=egm96_15.gtx
C:\Program Files\QGIS 3.26.2>proj=vgridshift +grids=egm96_15.gtx
proj_create: unrecognized format / unknown name
Rel. 9.0.1, June 15th, 2022
<proj>:
projection initialization failure
cause: Invalid PROJ string syntax
program abnormally terminated

, tested a modification but got another error:

C:\Program Files\QGIS 3.26.2>proj +proj=vgridshift +grids=egm96_15.gtx
Rel. 9.0.1, June 15th, 2022
<proj>:
can't initialize operations that produce angular output coordinates
program abnormally terminated

What am I missing?

Answer 1

I got help to solve this and the solution is here to anyone looking for this type of operation.

Using from WGS84 3D to EPSG:25830 + EGM2008

$echo 43.37551925 -2.69713764 102.635 | PROJ_NETWORK=ON cs2cs EPSG:4979
EPSG:25830+3855

Note:PROJ_NETWORK=ON is not needed if proj-data is properly installed.

$ echo 41.75758799 2.32152788  1235.577| cs2cs -d 3 EPSG:4979 EPSG:25831+3855

using Legacy syntax with:

$ echo -2.69713764 43.37551925 102.635| cs2cs -d 3 +init=epsg:4979 +to
+init=epsg:25830 +geoidgrids=egm08_25.gtx

And to transform coordinates in a file I create a small script in python:

import os
import pandas as pd
#command
os.system("echo 43.37551925 -2.69713764 102.635| cs2cs -d 3 EPSG:4326 EPSG:25830+3855 >> transformed.txt")
#import csv
data= pd.read_csv("tgernika0.csv")
data.columns
dft = data[['Longitude','Latitude', 'Ellipsoidal height']]
dft
#transform coordinates and save
for index, row in dft.iterrows():
       os.system("echo " + str(row['Latitude']) +" " + str(row['Longitude']) + " " + str(row['Ellipsoidal height']) + "| cs2cs -d 3 EPSG:4326 EPSG:25830+3855 >> transformed.txt")