Source code for meteodatalab.operators.vertical_extrapolation
"""Vertical extrapolation operators."""
# Standard library
from typing import Literal
# Third-party
import numpy as np
import xarray as xr
# Local
from .. import metadata
from .. import physical_constants as pc
LAPSE_RATE = 0.0065 # K m^-1
H1 = 2000.0
H2 = 2500.0
T1 = 298.0
[docs]
def extrapolate_temperature_sfc2p(
t_sfc: xr.DataArray,
h_sfc: xr.DataArray,
p_sfc: xr.DataArray,
p_target: float,
) -> xr.DataArray:
"""Extrapolate temperature to a target pressure level.
Implements the algorithm described in [1]_. The algorithm extrapolates
temperature from the surface to a target pressure level using a
polynomial expression of a dimensionless variable y, which is a function of
the surface temperature, surface pressure, and height. It assumes
a constant lapse rate of 0.0065 K m^-1 and dry air gas constant.
.. caution :
This extrapolation should be used with caution. Its intended use is to
extrapolate temperature to pressure levels below the surface, where
values are undefined. This is useful for applications where no missing values
are allowed, such as when training data-driven models. Results of the
extrapolation are not physically meaningful.
Parameters
----------
t_sfc : xr.DataArray
Surface temperature [K].
h_sfc : xr.DataArray
Surface height [m].
p_sfc : xr.DataArray
Surface pressure [Pa].
p_target : float
Target pressure level [Pa].
Returns
-------
xr.DataArray
Extrapolated temperature at the target pressure level.
References
----------
.. [1] https://www.umr-cnrm.fr/gmapdoc/IMG/pdf/ykfpos46t1r1.pdf
"""
y = _vertical_extrapolation_y_term(t_sfc, p_sfc, h_sfc, p_target)
res = t_sfc * (1 + y + (y**2) / 2 + (y**3) / 6)
res.attrs = metadata.override(
t_sfc.metadata, shortName="T", typeOfLevel="isobaricInPa"
)
res = _assign_vcoord(res, p_target)
return res
[docs]
def extrapolate_geopotential_sfc2p(
h_sfc: xr.DataArray,
t_sfc: xr.DataArray,
p_sfc: xr.DataArray,
p_target: float,
) -> xr.DataArray:
"""Extrapolate geopotential to a target pressure level.
Implements the algorithm described in [1]_. The algorithm extrapolates
geopotential from the surface to a target pressure level using a
polynomial expression of a dimensionless variable y, which is a function of
the surface temperature, surface pressure, and height. It assumes
a constant lapse rate of 0.0065 K m^-1 and dry air gas constant.
.. caution :
This extrapolation should be used with caution. Its intended use is to
extrapolate geopotential to pressure levels below the surface, where
values are undefined. This is useful for applications where no missing values
are allowed, such as when training data-driven models. Results of the
extrapolation are not physically meaningful.
Parameters
----------
h_sfc : xr.DataArray
Surface height [m].
t_sfc : xr.DataArray
Surface temperature [K].
p_sfc : xr.DataArray
Surface pressure [Pa].
p_target : float
Target pressure level [Pa].
Returns
-------
xr.DataArray
Extrapolated geopotential at the target pressure level.
References
----------
.. [1] https://www.umr-cnrm.fr/gmapdoc/IMG/pdf/ykfpos46t1r1.pdf
"""
y = _vertical_extrapolation_y_term(
t_sfc, p_sfc, h_sfc, p_target, lapse_rate=LAPSE_RATE
)
res = h_sfc * pc.g - pc.r_d * t_sfc * np.log(p_target / p_sfc) * (
1 + y / 2 + (y**2) / 6
)
res.attrs = metadata.override(
t_sfc.metadata, shortName="FI", typeOfLevel="isobaricInPa"
)
res = _assign_vcoord(res, p_target)
return res
[docs]
def extrapolate_k2p(
field: xr.DataArray,
p_target: float,
mode: Literal["constant"] = "constant",
) -> xr.DataArray:
"""Extrapolate a field to a target pressure level.
Currently, only the 'constant' extrapolation mode is implemented, where
the extrapolation is done by simply extending the values of the
lowermost model level to the target pressure level.
.. caution :
This extrapolation should be used with caution. Its intended use is to
extrapolate temperature to pressure levels below the surface, where
values are undefined. This is useful for applications where no missing values
are allowed, such as when training data-driven models. Results of the
extrapolation are not physically meaningful.
Parameters
----------
field : xr.DataArray
Field to extrapolate.
p_target : float
Target pressure level [Pa].
mode : str, optional
Extrapolation mode. Currently only 'constant' is implemented.
Returns
-------
xr.DataArray
Extrapolated field at the target pressure level.
References
----------
.. [1] https://www.umr-cnrm.fr/gmapdoc/IMG/pdf/ykfpos46t1r1.pdf
"""
return _assign_vcoord(field[{"z": [-1]}], p_target).assign_attrs(
metadata.override(field.metadata, typeOfLevel="isobaricInPa")
)
def _vertical_extrapolation_t_zero_prime(t_sfc, h_sfc):
t = t_sfc + LAPSE_RATE * h_sfc
t_min = np.minimum(t, T1)
return xr.where(h_sfc > H2, t_min, t_min * 0.5 + t * 0.5)
def _vertical_extrapolation_lapse_rate(h_sfc, t_sfc):
t_zero_prime = _vertical_extrapolation_t_zero_prime(t_sfc, h_sfc)
return xr.where(
h_sfc < H1,
LAPSE_RATE,
1 / h_sfc * np.maximum(t_zero_prime - t_sfc, 0.0),
)
def _vertical_extrapolation_y_term(
t_sfc, p_sfc, h_sfc, p_target, lapse_rate=None
) -> xr.DataArray:
if lapse_rate is None:
lapse_rate = _vertical_extrapolation_lapse_rate(h_sfc, t_sfc)
return lapse_rate * pc.r_d / pc.g * np.log(p_target / p_sfc)
def _assign_vcoord(x: xr.DataArray, p_target: float) -> xr.DataArray:
attrs = {
"units": "Pa",
"positive": "down",
"standard_name": "air_pressure",
"long_name": "pressure",
}
x = x.assign_coords(z=[p_target])
x["z"].attrs = attrs
return x