Convert raster from ellispoidal height to orthometric altitude

Differences between ellipsoidal height and orthometric altitude (or geoid height) are explained here https://spatialthoughts.com/2019/10/26/convert-between-orthometric-and-ellipsoidal-elevations-using-gdal/

Note: Geoid files format used to be GTX but switched to TIFF.

We use GDAL function gdalwarp to proceed :

gdalwarp -s_srs EPSG:XXXX -t_srs EPSG:XXXX Input_DEM.tif Output_DEM.tif

The spatial reference (SRS) can be set as EPSG code, or PROJ4 string. Prefer EPSG for simplicity. See https://spatialreference.org/ for list of SRS in EPSG code or other formats

if you want to use local geoid file, it must be downloaded to :

/usr/share/proj

to download a known geoid file, you can use projsync (ex. RAF20) :

projsync --source-id fr_ign --file fr_ign_RAF20.tif

or alternatively use wget from CDN

wget https://cdn.proj.org/fr_ign_RAF20.tif -O /usr/share/proj/fr_ign_RAF20.tif

Examples

  • from Lambert93 planimetric and NGF-IGN69 (with RAF20) to WGS84 ellipsoidal heights
gdalwarp -s_srs EPSG:2154+5720 -t_srs EPSG:4979 \
  MNT_LidarHD_5m_2154.tif MNT_LidarHD_5m_4326_He.tif

EPSG:2154+5720 = RGF93 Lambert-93 (horizontal) + NGF-IGN69 height (orthometric, gravity-related)

EPSG:4979 = WGS84 3D geographic (ellipsoidal height).

  • From UTM33N-WGS84 ellipsoidal height to UTM33N-EGM96 geoid height

gdalwarp -s_srs "+proj=utm +zone=33 +datum=WGS84 +units=m +no_defs" -t_srs "+proj=utm +zone=33 +datum=WGS84 +units=m +no_defs +geoidgrids=egm96_15.gtx" input.tif output.tif

Note the add of the geoid file with +geoidgrids=egm96_15.gtx in the command.

  • from altitude to ellipsoidal height (change srs and geoid grid if needed):

gdalwarp -s_srs "+proj=utm +zone=33 +datum=WGS84 +units=m +no_defs +geoidgrids=egm96_15.gtx" -t_srs "+proj=utm +zone=33 +datum=WGS84 +units=m +no_defs" input.tif output.tif

Examples of proj4 strings definitions

WGS84 ellipsoid with EGM96 elevation

+proj=longlat +datum=WGS84 +no_defs +geoidgrids=egm96_15.gtx

WGS84 ellipsoid with WGS84 elevation

+proj=longlat +datum=WGS84 +no_def

RGF93/Lambert 93 with NGF-IGN69 elevation

"+proj=lcc +lat_0=46.5 +lon_0=3 +lat_1=49 +lat_2=44 +x_0=700000 +y_0=6600000 +ellps=GRS80 +units=m +no_defs +type=crs +geoidgrids=fr_ign_RAF20.tif"