4

I want to convert wgs84 coordinates with ellipsoidal height to wgs84 coordinates with height above goid. Therefore I have set up a new conda environment and installed pyproj, proj_data, see conda list below. I found the EPSG codes 'EPSG:4979' for ellipsoidal height and 'EPSG:9705' for mean sea level.

Then I run:

import pyproj

# define coordinates
lat = 54.44
lon = 9.88
z = 100

# define CRSs
wgs84_ellipse = pyproj.crs.CRS.from_epsg(4979)
wgs84_msl = pyproj.crs.CRS.from_epsg(9705)

# create transformer
tform = pyproj.transformer.Transformer.from_crs(crs_from=wgs84_ellipse, crs_to=wgs84_msl, always_xy=True)

# transform coordinates
print(tform.transform(lon,lat,z))

Which gives me the same z-component back as I put in:

(9.88, 54.44, 100)

However, here I see that the goid height should be approx. 40m at this location:

https://geographiclib.sourceforge.io/cgi-bin/GeoidEval?input=54.44+9.88&option=Submit

Therefore, I would expect my code to return a value of about 100m-40m = 60m for the z-component.

Why is pyproj not giving me the height I expect?


Edit:
When I further run:

from pyproj.transformer import TransformerGroup
t = TransformerGroup(wgs84_ellipse, wgs84_msl, always_xy=True)

I get the following warning:

UserWarning: Best transformation is not available due to missing Grid(short_name=Und_min1x1_egm2008_isw=82_WGS84_TideFree.gz, full_name=, package_name=, url=, direct_download=False, open_license=False, available=False)

Then I went to

https://earth-info.nga.mil/index.php?dir=wgs84&action=wgs84

and downloaded under Data/Apps -> EGM2008, this file:

https://earth-info.nga.mil/php/download.php?file=egm-08spherical

The readme mentions the file I am looking for called "Und_min1x1_egm2008_isw=82_WGS84_TideFree.gz". However, after extracting the archive, I can not find the file.

Now the questions is: How do I get the missing grid file?


To create the conda environment I ran:

conda create -n "test_env" python=3.11.0
conda install -c conda-forge pyproj
conda install -c conda-forge proj-data

conda list returns:

# Name                    Version                   Build  Channel
bzip2                     1.0.8                h8ffe710_4    conda-forge
ca-certificates           2022.12.7            h5b45459_0    conda-forge
certifi                   2022.12.7          pyhd8ed1ab_0    conda-forge
jpeg                      9e                   hcfcfb64_3    conda-forge
krb5                      1.20.1               heb0366b_0    conda-forge
lerc                      4.0.0                h63175ca_0    conda-forge
libcurl                   7.88.1               h68f0423_0    conda-forge
libdeflate                1.17                 hcfcfb64_0    conda-forge
libffi                    3.4.2                h8ffe710_5    conda-forge
libsqlite                 3.40.0               hcfcfb64_0    conda-forge
libssh2                   1.10.0               h9a1e1f7_3    conda-forge
libtiff                   4.5.0                hf8721a0_2    conda-forge
libzlib                   1.2.13               hcfcfb64_4    conda-forge
openssl                   3.0.8                hcfcfb64_0    conda-forge
pip                       23.0.1             pyhd8ed1ab_0    conda-forge
proj                      9.1.1                heca977f_2    conda-forge
proj-data                 1.13                 hd8ed1ab_0    conda-forge
pyproj                    3.4.1           py311h36482e4_1    conda-forge
python                    3.11.0          hcf16a7b_0_cpython    conda-forge
python_abi                3.11                    3_cp311    conda-forge
setuptools                67.4.0             pyhd8ed1ab_0    conda-forge
sqlite                    3.40.0               hcfcfb64_0    conda-forge
tk                        8.6.12               h8ffe710_0    conda-forge
tzdata                    2022g                h191b570_0    conda-forge
ucrt                      10.0.22621.0         h57928b3_0    conda-forge
vc                        14.3                hb6edc58_10    conda-forge
vs2015_runtime            14.34.31931         h4c5c07a_10    conda-forge
wheel                     0.38.4             pyhd8ed1ab_0    conda-forge
xz                        5.2.6                h8d14728_0    conda-forge
zstd                      1.5.2                h12be248_6    conda-forge
4
  • Check this thread, please. It may help. Commented Mar 3, 2023 at 11:30
  • Thank you that did help to a degree! However, now I am looking for the missing grid file. Do you have an idea where I can download/find the file?
    – sehan2
    Commented Mar 3, 2023 at 12:44
  • Do you mean egmgrid? Commented Mar 3, 2023 at 13:08
  • pyproj is refering to this pyproj4.github.io/pyproj/stable/transformation_grids.html Or is there an alternative to the file pyproj suggests?
    – sehan2
    Commented Mar 3, 2023 at 13:11

2 Answers 2

1

One would think that this is most common and most basic height transformation there is, however, as of now, I was not able to compute it using pyproj.

Alternative using rasterio:

So here is an alternative using rasterio and the egm08_25.gtx-file (EGM2008) that gives me the result I am looking for:

import rasterio
from pathlib import Path

# define coordinates
lat = 54.44
lon = 9.88
z = 100

# open .gtx file using rasterio
path = Path(r"C:\<path_to_gtx_file>\egm08_25.gtx")
src = rasterio.open(path)

# compute undulation
undulation = next(src.sample([(lon, lat)]))[0]
# compute height over goid, i.e. mean sea level
orthometric_height = z - undulation
print('undulation:', undulation, ', orthometric_height:', orthometric_height)

which returns:

undulation: 39.903717 , orthometric_height: 60.096282958984375

This is close enough to the value I get from https://geographiclib.sourceforge.io/cgi-bin/GeoidEval?input=54.44%209.88&option=Submit.


Used conda environment:

Since in my experience GIS related python packages do not tend to work reliably when using different python versions (and in conjunction with various other python packages) here is how I set up my conda environment:

conda create -n "test_env2" python=3.11.0
conda install -c conda-forge rasterio

conda list returns:

# Name                    Version                   Build  Channel
affine                    2.4.0              pyhd8ed1ab_0    conda-forge
attrs                     22.2.0             pyh71513ae_0    conda-forge
blosc                     1.21.3               hdccc3a2_0    conda-forge
boost-cpp                 1.78.0               h9f4b32c_1    conda-forge
bzip2                     1.0.8                h8ffe710_4    conda-forge
ca-certificates           2022.12.7            h5b45459_0    conda-forge
cairo                     1.16.0            hd694305_1014    conda-forge
certifi                   2022.12.7          pyhd8ed1ab_0    conda-forge
cfitsio                   4.2.0                h9ebe7e4_0    conda-forge
click                     8.1.3           win_pyhd8ed1ab_2    conda-forge
click-plugins             1.1.1                      py_0    conda-forge
cligj                     0.7.2              pyhd8ed1ab_1    conda-forge
colorama                  0.4.6              pyhd8ed1ab_0    conda-forge
curl                      7.88.1               h68f0423_0    conda-forge
expat                     2.5.0                h1537add_0    conda-forge
font-ttf-dejavu-sans-mono 2.37                 hab24e00_0    conda-forge
font-ttf-inconsolata      3.000                h77eed37_0    conda-forge
font-ttf-source-code-pro  2.038                h77eed37_0    conda-forge
font-ttf-ubuntu           0.83                 hab24e00_0    conda-forge
fontconfig                2.14.2               hbde0cde_0    conda-forge
fonts-conda-ecosystem     1                             0    conda-forge
fonts-conda-forge         1                             0    conda-forge
freetype                  2.12.1               h546665d_1    conda-forge
freexl                    1.0.6                h67ca5e6_1    conda-forge
geos                      3.11.1               h1537add_0    conda-forge
geotiff                   1.7.1                hc256dc0_6    conda-forge
gettext                   0.21.1               h5728263_0    conda-forge
hdf4                      4.2.15               h1b1b6ef_5    conda-forge
hdf5                      1.12.2          nompi_h57737ce_101    conda-forge
icu                       70.1                 h0e60522_0    conda-forge
intel-openmp              2023.0.0         h57928b3_25922    conda-forge
jpeg                      9e                   hcfcfb64_3    conda-forge
kealib                    1.5.0                h61be68b_0    conda-forge
krb5                      1.20.1               heb0366b_0    conda-forge
lcms2                     2.14                 ha5c8aab_1    conda-forge
lerc                      4.0.0                h63175ca_0    conda-forge
libaec                    1.0.6                h63175ca_1    conda-forge
libblas                   3.9.0              16_win64_mkl    conda-forge
libcblas                  3.9.0              16_win64_mkl    conda-forge
libcurl                   7.88.1               h68f0423_0    conda-forge
libdeflate                1.17                 hcfcfb64_0    conda-forge
libffi                    3.4.2                h8ffe710_5    conda-forge
libgdal                   3.6.2                h329cd9e_7    conda-forge
libglib                   2.74.1               he8f3873_1    conda-forge
libhwloc                  2.9.0                h51c2c0f_0    conda-forge
libiconv                  1.17                 h8ffe710_0    conda-forge
libkml                    1.3.0             hf2ab4e4_1015    conda-forge
liblapack                 3.9.0              16_win64_mkl    conda-forge
libnetcdf                 4.9.1           nompi_hc41bf00_101    conda-forge
libpng                    1.6.39               h19919ed_0    conda-forge
libpq                     15.2                 ha9684e8_0    conda-forge
librttopo                 1.1.0               he22b5cd_12    conda-forge
libspatialite             5.0.1               hfdcade0_23    conda-forge
libsqlite                 3.40.0               hcfcfb64_0    conda-forge
libssh2                   1.10.0               h9a1e1f7_3    conda-forge
libtiff                   4.5.0                hf8721a0_2    conda-forge
libwebp-base              1.2.4                h8ffe710_0    conda-forge
libxml2                   2.10.3               hc3477c8_0    conda-forge
libzip                    1.9.2                h519de47_1    conda-forge
libzlib                   1.2.13               hcfcfb64_4    conda-forge
lz4-c                     1.9.4                hcfcfb64_0    conda-forge
mkl                       2022.1.0           h6a75c08_874    conda-forge
numpy                     1.24.2          py311h0b4df5a_0    conda-forge
openjpeg                  2.5.0                ha2aaf27_2    conda-forge
openssl                   3.0.8                hcfcfb64_0    conda-forge
pcre2                     10.40                h17e33f8_0    conda-forge
pip                       23.0.1             pyhd8ed1ab_0    conda-forge
pixman                    0.40.0               h8ffe710_0    conda-forge
poppler                   23.01.0              h183ae7b_0    conda-forge
poppler-data              0.4.12               hd8ed1ab_0    conda-forge
postgresql                15.2                 hd87cd2b_0    conda-forge
proj                      9.1.1                heca977f_2    conda-forge
pthreads-win32            2.9.1                hfa6e2cd_3    conda-forge
pyparsing                 3.0.9              pyhd8ed1ab_0    conda-forge
python                    3.11.0          hcf16a7b_0_cpython    conda-forge
python_abi                3.11                    3_cp311    conda-forge
rasterio                  1.3.6           py311hc06ee40_0    conda-forge
setuptools                67.4.0             pyhd8ed1ab_0    conda-forge
snappy                    1.1.9                hfb803bf_2    conda-forge
snuggs                    1.4.7                      py_0    conda-forge
sqlite                    3.40.0               hcfcfb64_0    conda-forge
tbb                       2021.8.0             h91493d7_0    conda-forge
tiledb                    2.13.2               h3132609_0    conda-forge
tk                        8.6.12               h8ffe710_0    conda-forge
tzdata                    2022g                h191b570_0    conda-forge
ucrt                      10.0.22621.0         h57928b3_0    conda-forge
vc                        14.3                hb6edc58_10    conda-forge
vs2015_runtime            14.34.31931         h4c5c07a_10    conda-forge
wheel                     0.38.4             pyhd8ed1ab_0    conda-forge
xerces-c                  3.2.4                h63175ca_1    conda-forge
xz                        5.2.6                h8d14728_0    conda-forge
zlib                      1.2.13               hcfcfb64_4    conda-forge
zstd                      1.5.2                h12be248_6    conda-forge
1

The solution (derived from Scott Staniewicz's answer to this question) is to use a composite transformation:

from pyproj import Transformer
t=Transformer.from_crs("epsg:4979", "epsg:4326+3855", always_xy=True)
t.transform(9.88, 54.44, 100)
(9.88, 54.44, 60.110880718994146)

This builds a transformer that converts from the 3D ellipsoidal coordinate system epsg:4979 with elevation as height above ellipsoid in meters to a compound CRS constructed from epsg:4326 (2D ellipsoidal coordinate system) and epsg:3855 (the EGM2008 geoid model). It returns the orthometric elevation in meters above the geoid ("elevation MSL").

If the transformer returns the unchanged input elevation, you are missing the required EGM2008 geoid model. I found the filename us_nga_egm08_25.tif via TransformerGroup error message, but you can also download all available transformation grids:

PROJ_DOWNLOAD_DIR="$(python -c 'import pyproj; print(pyproj.datadir.get_user_data_dir())')"
cd "$PROJ_DOWNLOAD_DIR"
curl -O 'https://cdn.proj.org/us_nga_egm08_25.tif' 

I spot-checked the results against a version of GeoTrans and they match within a meter, at least for coordinates near the geoid. For example, Geotrans returns a geoid height of 39.603 m @ (lat=54.44, lon=9.88), so 100 m above the ellipsoid would translate to an orthometric height of (100 m - 39.603 m) or 60.397 m. It would appear pyproj is using a tiff raster for the geoid definition, rather than a polynominal model like Geotrans.

Thanks for asking this question! My data originates from GPS/GNSS systems, which compute height above ellipsoid, and I often need to plot it using tools like Google Earth, which (best as I can tell) uses MSL. I could not figure out how to do this from the pyproj documentation.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.