lenapy.utils.gravity

The gravity module provides functions for manipulating spherical harmonics datasets with respect to different reference frames, tide systems, normal gravity reference and smoothing techniques.

This module includes functions to:

  • Change the reference frame of spherical harmonics datasets.

  • Change the tide system of spherical harmonics datasets.

  • Modify degree 1 Love numbers for different reference frames.

  • Estimate normal gravity values.

  • Generate Gaussian weights for smoothing spherical harmonics datasets.

  • Estimate time variable gravity field from a .gfc temporal file.

Examples

>>> import xarray as xr
>>> from lenapy.utils.gravity import *
# Load a dataset
>>> ds = xr.open_dataset('example_file.nc')
# Change the reference frame of the dataset
>>> ds_new_ref = change_reference(ds, new_radius=6378137.0, new_earth_gravity_constant=3.986004418e14)
# Change the tide system of the dataset
>>> ds_new_tide_system = change_tide_system(ds, new_tide='tide_free')
# Generate Gaussian weights for smoothing spherical harmonics datasets.
>>> ds_gauss_weights = gauss_weights(radius=100000, lmax=60)
# Filter the dataset with Gaussian weights
>>> ds['clm'] *= ds_gauss_weights
>>> ds['slm'] *= ds_gauss_weights
lenapy.utils.gravity.change_reference(ds: Dataset, new_radius: float = 6378137.0, new_earth_gravity_constant: float = 398600441500000.0, old_radius: float | None = None, old_earth_gravity_constant: float | None = None, apply: bool = False) → Dataset[source]

Spherical Harmonics dataset are associated with an earth radius a and µ or GM the earth gravity constant. This function return the dataset updated with the new constants associated to it in input (as a deep copy by default). The current reference frame constants can be given in parameters or can be contained in ds.attrs[‘radius’] and ds.attrs[‘earth_gravity_constant’].

Warning ! This function must be applied before removing the mean values of the dataset through time.

Parameters:

  • ds (xr.Dataset) – xr.Dataset that corresponds to SH data associated with a current reference frame (constants whise) to update.

  • new_radius (float) – New Earth radius constant in meters. Default is LNPY_A_EARTH_GRS80 define in constants.py

  • new_earth_gravity_constant (float) – New gravitational constant of the Earth in m³/s². Default is LNPY_GM_EARTH define in constants.py

  • old_radius (float | None, optional) – Current Earth radius constant of the dataset ds in meters. If not provided, uses ds.attrs[‘radius’].

  • old_earth_gravity_constant (float | None, optional) – Current gravitational constant of the Earth of the dataset ds in m³/s². If not provided, uses ds.attrs[‘earth_gravity_constant’].

  • apply (bool, optional) – If True, apply the update to the current dataset without making a deep copy. Default is False.

Returns:

ds_out – Updated dataset with the new constants.

Return type:

xr.Dataset

Raises:

KeyError – If the current reference frame constants are not provided and not found in the dataset attributes.

Examples

>>> ds_new_ref = change_reference(ds, new_radius=6378137.0, new_earth_gravity_constant=3.986004418e14)
lenapy.utils.gravity.change_tide_system(ds: Dataset, new_tide: Literal['tide_free', 'zero_tide', 'mean_tide'], old_tide: Literal['tide_free', 'zero_tide', 'mean_tide'] | None = None, k20: float | None = None, apply: bool = False) → Dataset[source]

Apply a C20 offset to the dataset to change the tide system. Follows [IERS2010] convention to convert between tide system (‘tide_free’, ‘zero_tide’, ‘mean_tide’). Warning : It should be noted that each center use his one tide offset, it may differ in terms of offset formula that can come from [IERS2010] and [Smith1998]. The value of k20 can also change between centers.

Parameters:

  • ds (xr.Dataset) – xr.Dataset that corresponds to SH data associated with a current tide system to update.

  • new_tide (str) – Output tidal system, either ‘tide_free’, ‘zero_tide’ or ‘mean_tide’.

  • old_tide (str | None, optional) – Input tidal system. If not provided, uses ds.attrs[‘tide_system’].

  • k20 (float | None, optional) – k20 Earth tide external potential Love numbers. If not provided, the default value from [IERS2010] is used.

  • apply (bool, optional) – If True, apply the update to the current dataset without making a deep copy. Default is False.

Returns:

ds_out – Updated dataset with the new tidal system. Warning, the dataset is also updated in place if copy=False (default)

Return type:

xr.Dataset

Raises:

ValueError | KeyError – If the input tidal system is not provided and not found in the dataset attributes.

Examples

>>> ds_new_tide_system = change_tide_system(ds, new_tide='tide_free')

References

[IERS2010] (1,2,3)

G. Petit and B. Luzum, “IERS Conventions (2010)”, IERS Technical Note, 36, 88–89, (2010). doi: 10.1007/s00190-002-0216-2

[Smith1998]

D. A. Smith, “There is no such thing as ‘The’ EGM96 geoid: Subtle points on the use of a global geopotential model”, IGeS Bulletin, International Geoid Service, 8, 17-28, (1998). link to archive

lenapy.utils.gravity.change_love_reference_frame(ds: Dataset, new_frame: Literal['CM', 'CE', 'CF', 'CL', 'CH'], old_frame: Literal['CM', 'CE', 'CF', 'CL', 'CH'], apply: bool = False) → Dataset[source]

Modify degree 1 love numbers of the dataset to change the reference frame. The input dataset need to contain ‘hl’, ‘ll’, ‘kl’ variables that are potential love numbers with a degree dimension named ‘l’ with at least coordinate l=1.

This function follows [Blewitt2003] to convert between reference frames. Reference frames can be : * ‘CM’ for Center of Mass of the Earth System * ‘CE’ for Center of Mass of the Solid Earth * ‘CF’ for Center of Surface Figure * ‘CL’ for Center of Surface Lateral Figure * ‘CH’ for Center of Surface Height Figure

Parameters:

  • ds (xr.Dataset) – Dataset with ‘hl’, ‘ll’, ‘kl’ variables and ‘l’ dimension.

  • new_frame (str) – Reference frame of the output dataset. Either ‘CM’, ‘CE’, ‘CF’, ‘CL’ or ‘CH’.

  • old_frame (str) – Reference frame of the input dataset. Either ‘CM’, ‘CE’, ‘CF’, ‘CL’ or ‘CH’.

  • apply (bool) – If True, apply the update to the current dataset without making a deep copy. Default is False.

Returns:

ds_out – Dataset with the updated degree 1 love numbers. Warning, the dataset is also updated in place if copy=False (default)

Return type:

xr.Dataset

References

[Blewitt2003]

G. Blewitt, “Self-consistency in reference frames, geocenter definition, and surface loading of the solid Earth”, Journal of Geophysical Research: Solid Earth, 108, 2103, (2003). doi: 10.1029/2002JB002082 <https://doi.org/10.1029/2002JB002082>`_

lenapy.utils.gravity.estimate_normal_gravity(geographic_latitude: DataArray | ndarray | None = None, a_earth: float | None = 6378137.0, earth_gravity_constant: float | None = 398600441500000.0, f_earth: float = 0.003352810681182319, omega_earth: float = 7.292115e-05) → DataArray[source]

Estimate normal acceleration with the Somigliana equation at given geographic latitude for a group of Earth’s parameters that represent the ellipsoid.

Parameters:

  • geographic_latitude (xr.DataArray | np.ndarray, optional) – Geographic / geodetic latitude for ellipsoidal Earth radius computation in radians.

  • a_earth (float, optional) – Earth semi-major axis [m]. Default is LNPY_A_EARTH_GRS80.

  • earth_gravity_constant (float, optional) – Standard gravitational parameter for Earth [m³.s⁻²]. Default is LNPY_GM_EARTH.

  • f_earth (float, optional) – Earth flattening. Default is LNPY_F_EARTH_GRS80.

  • omega_earth (float, optional) – Earth’s rotation rate [rad.s⁻¹]. Default is LNPY_OMEGA_EARTH.

Returns:

gamma_0 – Normal gravity at given geographic colatitude.

Return type:

xr.DataArray | np.ndarray

lenapy.utils.gravity.apply_normal_zonal_correction(ds: Dataset, radius: float | None = None, earth_gravity_constant: float | None = None, f_earth: float = 0.003352810681182319, omega_earth: float = 7.292115e-05, reverse: bool = False, apply: bool = False) → Dataset[source]

Apply a correction of the normal gravity field on zonal coefficients on a SH dataset for a specified ellipsoid.

Parameters:

  • ds (xr.Dataset) – xr.Dataset that corresponds to SH data to be correct for the normal gravity field.

  • radius (float | None, optional) – Earth radius constant of the dataset ds in meters. If not provided, uses ds.attrs[‘radius’].

  • earth_gravity_constant (float | None, optional) – Current gravitational constant of the Earth of the dataset ds in m³/s². If not provided, uses ds.attrs[‘earth_gravity_constant’].

  • f_earth (float, optional) – Earth flattening. Default is LNPY_F_EARTH_GRS80.

  • omega_earth (float, optional) – Earth rotation rate. Default is LNPY_OMEGA_EARTH_GRS80.

  • reverse (bool, optional) – False to apply the correction, True to remove the correction.

  • apply (bool, optional) – If True, apply the update to the current dataset without making a deep copy. Default is False.

Returns:

ds_out – Updated dataset with the correction.

Return type:

xr.Dataset

Raises:

KeyError – If the current reference frame constants are not provided and not found in the dataset attributes.

lenapy.utils.gravity.gauss_weights(radius: float, lmax: int, a_earth: float = 6378137.0, cutoff: float = 1e-10) → DataArray[source]

Generate a xr.DataArray with Gaussian weights as a function of degree.

This function uses the method described in [Jekeli1981] to generate Gaussian weights for spherical harmonics coefficients.

Parameters:

  • radius (float) – Gaussian smoothing radius in meters.

  • lmax (int) – Maximum degree of spherical harmonic coefficients.

  • a_earth (float, optional) – Radius of the Earth in meters. Default is LNPY_A_EARTH_GRS80 define in constants.py

  • cutoff (float, optional) – Minimum value for the tail of the Gaussian averaging function (see [Jekeli1981] p.18), default is 1e-10.

Returns:

gaussian_weights – Degree-dependent Gaussian weights xr.DataArray.

Return type:

xr.DataArray

Examples

>>> ds_gauss_weights = gauss_weights(radius=100000, lmax=60)

References

[Jekeli1981] (1,2)

C. Jekeli, “Alternative Methods to Smooth the Earth’s Gravity Field”, NASA Grant No. NGR 36-008-161, OSURF Proj. No. 783210, 48 pp., (1981).

lenapy.utils.gravity.gfct_field_estimation(ds: Dataset, time: datetime64 | ndarray | list | tuple | DataArray) → Dataset[source]

Compute time-variable gravity field from variation coefficients contain in ‘.gfc’ file with icgem1.0 or icgem2.0 format.

Parameters:

  • ds (xr.Dataset) – Output of the reader from the ‘.gfc’ file with gfct information.

  • time (np.datetime64 | list, tuple, np.ndarray of np.datetime64 | xr.DataArray of np.datetime64) – Time information to compute the coefficients variations on.

Returns:

ds_out – Estimated coefficients time variations.

Return type:

xr.Dataset