Earth geoid height grid at 10 arc-minute resolution

The grid is grid-node registered and stored in netCDF with CF-compliant metadata. The geoid height is derived from the EIGEN-6C4 spherical harmonic model of the Earth’s gravity field.

Original source: EIGEN-6C4 model

Pre-processing: Source code for preparation of the original dataset for redistribution in Ensaio

import pygmt
import xarray as xr

import ensaio

Download and cache the data and return the path to it on disk.

fname = ensaio.fetch_earth_geoid(version=1)
print(fname)

/home/runner/work/_temp/cache/ensaio/v1/earth-geoid-10arcmin.nc

Load the netCDF grid with xarray.

data = xr.load_dataarray(fname)
data

<xarray.DataArray 'geoid' (latitude: 1081, longitude: 2161)> Size: 19MB
array([[-29.5, -29.5, -29.5, ..., -29.5, -29.5, -29.5],
       [-29.5, -29.5, -29.5, ..., -29.5, -29.5, -29.5],
       [-29.6, -29.6, -29.6, ..., -29.6, -29.6, -29.6],
       ...,
       [ 14.7,  14.7,  14.7, ...,  14.7,  14.7,  14.7],
       [ 15.2,  15.2,  15.2, ...,  15.2,  15.2,  15.2],
       [ 15.4,  15.4,  15.4, ...,  15.4,  15.4,  15.4]],
      shape=(1081, 2161))
Coordinates:
  * longitude  (longitude) float64 17kB -180.0 -179.8 -179.7 ... 179.8 180.0
  * latitude   (latitude) float64 9kB -90.0 -89.83 -89.67 ... 89.67 89.83 90.0
Attributes:
    Conventions:     CF-1.8
    title:           Geoid height (EIGEN-6C4) with respect to WGS84
    crs:             WGS84
    source:          Generated from the EIGEN-6C4 model by the ICGEM Calculat...
    license:         Creative Commons Attribution 4.0 International Licence
    references:      https://doi.org/10.5880/icgem.2015.1
    long_name:       geoid height
    standard_name:   geoid_height_above_reference_ellipsoid
    description:     height of the geoid with respect to the WGS84 ellipsoid
    units:           m
    actual_range:    [-106.5   86. ]
    icgem_metadata:  generating_institute: gfz-potsdam\ngenerating_date: 2021...

xarray.DataArray

'geoid'

• latitude: 1081
• longitude: 2161

• -29.5 -29.5 -29.5 -29.5 -29.5 -29.5 ... 15.4 15.4 15.4 15.4 15.4 15.4

  array([[-29.5, -29.5, -29.5, ..., -29.5, -29.5, -29.5],
         [-29.5, -29.5, -29.5, ..., -29.5, -29.5, -29.5],
         [-29.6, -29.6, -29.6, ..., -29.6, -29.6, -29.6],
         ...,
         [ 14.7,  14.7,  14.7, ...,  14.7,  14.7,  14.7],
         [ 15.2,  15.2,  15.2, ...,  15.2,  15.2,  15.2],
         [ 15.4,  15.4,  15.4, ...,  15.4,  15.4,  15.4]],
        shape=(1081, 2161))

• Coordinates: (2)
  • longitude
    (longitude)
    float64
    -180.0 -179.8 ... 179.8 180.0

    long_name :

    longitude

    standard_name :

    longitude

    units :

    degrees_east

    actual_range :

    [-180. 180.]

    array([-180.      , -179.833333, -179.666667, ...,  179.666667,  179.833333,
            180.      ], shape=(2161,))

  • latitude
    (latitude)
    float64
    -90.0 -89.83 -89.67 ... 89.83 90.0

    long_name :

    latitude

    standard_name :

    latitude

    units :

    degrees_north

    actual_range :

    [-90. 90.]

    array([-90.      , -89.833333, -89.666667, ...,  89.666667,  89.833333,
            90.      ], shape=(1081,))

• Indexes: (2)
  • longitude
    PandasIndex

    PandasIndex(Index([             -180.0, -179.83333333333334, -179.66666666666666,
                        -179.5, -179.33333333333334, -179.16666666666666,
                        -179.0, -178.83333333333334, -178.66666666666666,
                        -178.5,
           ...
                         178.5,  178.66666666666663,  178.83333333333331,
                         179.0,  179.16666666666663,  179.33333333333331,
                         179.5,  179.66666666666663,  179.83333333333331,
                         180.0],
          dtype='float64', name='longitude', length=2161))

  • latitude
    PandasIndex

    PandasIndex(Index([             -90.0, -89.83333333333333, -89.66666666666667,
                        -89.5, -89.33333333333333, -89.16666666666667,
                        -89.0, -88.83333333333333, -88.66666666666667,
                        -88.5,
           ...
                         88.5,  88.66666666666666,  88.83333333333331,
                         89.0,  89.16666666666666,  89.33333333333331,
                         89.5,  89.66666666666666,  89.83333333333331,
                         90.0],
          dtype='float64', name='latitude', length=1081))

• Attributes: (12)

  Conventions :

  CF-1.8

  title :

  Geoid height (EIGEN-6C4) with respect to WGS84

  crs :

  WGS84

  source :

  Generated from the EIGEN-6C4 model by the ICGEM Calculation Service (http://icgem.gfz-potsdam.de)

  license :

  Creative Commons Attribution 4.0 International Licence

  references :

  https://doi.org/10.5880/icgem.2015.1

  long_name :

  geoid height

  standard_name :

  geoid_height_above_reference_ellipsoid

  description :

  height of the geoid with respect to the WGS84 ellipsoid

  units :

  m

  actual_range :

  [-106.5 86. ]

  icgem_metadata :

  generating_institute: gfz-potsdam generating_date: 2021/06/11 product_type: gravity_field body: earth modelname: EIGEN-6C4 max_used_degree: 720 gentlecut: 540 to 720 tide_system: tide_free functional: geoid zero_degree_term: included unit: meter crust_density: 2670.0 kg/m**3 refsysname: WGS84 gmrefpot: 3.98600441800E+14 m**3/s**2 radiusrefpot: 6378137.000 m flatrefpot: 3.352810664747480E-03 (1/298.25722356300) omegarefpot: 7.29211500000E-05 1/s normal_potential: 6.263685171456948E+07 m**2/s**2 long_lat_unit: degree latlimit_north: 90.000000000000 latlimit_south: -90.000000000000 longlimit_west: -180.00000000000 longlimit_east: 180.00000000000 gridstep: 0.16666666666667 latitude_parallels: 1081 longitude_parallels: 2161 number_of_gridpoints: 2336041 gapvalue: 999.0000 weighted_mean: -5.2816814E-02 meter maxvalue: 8.5964166E+01 meter minvalue: -1.0652440E+02 meter signal_wrms: 3.0591615E+01 meter grid_format: long_lat_value attributes: longitude latitude geoid attributes_units: deg. deg. meter

Make a PyGMT pseudo-color map of the grid in a Mollweide projection.

fig = pygmt.Figure()
fig.basemap(
    region="g",
    projection="W15c",
    frame=True,
)
fig.grdimage(data, cmap="polar+h")
fig.colorbar(frame='af+l"geoid height [m]"')
fig.coast(shorelines=True, resolution="c", area_thresh=1e4)
fig.show()