Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Revert extra arguments added to handle zero humidity in thermo computations #26

Merged
merged 2 commits into from
Nov 5, 2024
Merged

Revert extra arguments added to handle zero humidity in thermo computations #26

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
8fa1108
Revert extra arguments added to handle zero humidity in thermo comput…
sandorkertesz Nov 5, 2024
28525b6
Revert extra arguments added to handle zero humidity in thermo comput…
sandorkertesz Nov 5, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions docs/release_notes/index.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,5 +4,6 @@ Release notes
.. toctree::
:maxdepth: 1

version_0.3_updates
version_0.2_updates
version_0.1_updates
14 changes: 14 additions & 0 deletions docs/release_notes/version_0.3_updates.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,14 @@
Version 0.3 Updates
/////////////////////////


Version 0.3.0
===============

- Removed the ``eps`` and ``out`` keywords arguments to control zero humidity in

- :py:meth:`dewpoint_from_specific_humidity <meteo.thermo.array.dewpoint_from_specific_humidity>`
- :py:meth:`dewpoint_from_relative_humidity <meteo.thermo.array.dewpoint_from_relative_humidity>`
- :py:meth:`temperature_from_saturation_vapour_pressure <meteo.thermo.array.temperature_from_saturation_vapour_pressure>`

This was done to revert the changes added in version 0.2.0.
88 changes: 12 additions & 76 deletions src/earthkit/meteo/thermo/array/thermo.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -263,21 +263,9 @@ def compute_slope(self, t, phase):
elif phase == "ice":
return self._es_ice_slope(t)

def t_from_es(self, es, eps=1e-8, out=None):
def _comp(x):
x = x / self.c1
x = np.log(x)
return (x * self.c4w - self.c3w * self.t0) / (x - self.c3w)

if out is not None:
v = np.asarray(es)
z_mask = v < eps
v_mask = ~z_mask
v[v_mask] = _comp(v[v_mask])
v[z_mask] = out
return v
else:
return _comp(es)
def t_from_es(self, es):
v = np.log(es / self.c1)
return (v * self.c4w - self.c3w * self.t0) / (v - self.c3w)

def _es_water(self, t):
return self.c1 * np.exp(self.c3w * (t - self.t0) / (t - self.c4w))
Expand Down Expand Up @@ -572,35 +560,26 @@ def saturation_specific_humidity_slope(t, p, es=None, es_slope=None, phase="mixe
return constants.epsilon * es_slope * p / v


def temperature_from_saturation_vapour_pressure(es, eps=1e-8, out=None):
def temperature_from_saturation_vapour_pressure(es):
r"""Compute the temperature from saturation vapour pressure.

Parameters
----------
es: ndarray
:func:`saturation_vapour_pressure` (Pa)
eps: number
If ``out`` is not None, return ``out`` when ``es`` < ``eps``. If out
is None, ``eps`` is ignored and return np.nan for ``es`` values very
close to zero. *New in version 0.2.0*
out: number or None
If not None, return ``out`` when ``es`` < ``eps``. If None, ``eps`` is
ignored and return np.nan for ``es`` values very close to zero.
*New in version 0.2.0*


Returns
-------
ndarray
Temperature (K)
Temperature (K). For zero ``es`` values returns np.nan.


The computation is always based on the "water" phase of
the :func:`saturation_vapour_pressure` formulation irrespective of the
phase ``es`` was computed to.

"""
return _EsComp().t_from_es(es, eps=eps, out=out)
return _EsComp().t_from_es(es)


def relative_humidity_from_dewpoint(t, td):
Expand Down Expand Up @@ -775,7 +754,7 @@ def specific_humidity_from_relative_humidity(t, r, p):
return specific_humidity_from_vapour_pressure(e, p)


def dewpoint_from_relative_humidity(t, r, eps=1e-8, out=None):
def dewpoint_from_relative_humidity(t, r):
r"""Compute the dewpoint temperature from relative humidity.

Parameters
Expand All @@ -784,19 +763,11 @@ def dewpoint_from_relative_humidity(t, r, eps=1e-8, out=None):
Temperature (K)
r: ndarray
Relative humidity (%)
eps: number
If ``out`` is not None, return ``out`` when ``r`` < ``eps``.
If out is None, ``eps`` is ignored and return np.nan for ``r``
values very close to zero. *New in version 0.2.0*
out: number or None
If not None, return ``out`` when ``r`` < ``eps``. If None, ``eps`` is
ignored and return np.nan for ``r`` values very close to zero.
*New in version 0.2.0*

Returns
-------
ndarray
Dewpoint temperature (K)
Dewpoint temperature (K). For zero ``r`` values returns np.nan.


The computation starts with determining the the saturation vapour pressure over
Expand All @@ -815,24 +786,11 @@ def dewpoint_from_relative_humidity(t, r, eps=1e-8, out=None):
equations used in :func:`saturation_vapour_pressure`.

"""
# by masking upfront we avoid RuntimeWarning when calling log() in
# the computation of td when r is very small
if out is not None:
r = np.asarray(r)
mask = r < eps
r[mask] = np.nan

es = saturation_vapour_pressure(t, phase="water") * r / 100.0
td = temperature_from_saturation_vapour_pressure(es)
return temperature_from_saturation_vapour_pressure(es)

if out is not None and not np.isnan(out):
td = np.asarray(td)
td[mask] = out

return td


def dewpoint_from_specific_humidity(q, p, eps=1e-8, out=None):
def dewpoint_from_specific_humidity(q, p):
r"""Compute the dewpoint temperature from specific humidity.

Parameters
Expand All @@ -841,19 +799,11 @@ def dewpoint_from_specific_humidity(q, p, eps=1e-8, out=None):
Specific humidity (kg/kg)
p: ndarray
Pressure (Pa)
eps: number
If ``out`` is not None, return ``out`` when ``q`` < ``eps``.
If out is None, ``eps`` is ignored and return np.nan for ``q``
values very close to zero. *New in version 0.2.0*
out: number or None
If not None, return ``out`` when ``q`` < ``eps``. If None, ``eps`` is
ignored and return np.nan for ``q`` values very close to zero.
*New in version 0.2.0*

Returns
-------
ndarray
Dewpoint temperature (K)
Dewpoint temperature (K). For zero ``q`` values returns np.nan.


The computation starts with determining the the saturation vapour pressure over
Expand All @@ -873,21 +823,7 @@ def dewpoint_from_specific_humidity(q, p, eps=1e-8, out=None):
used in :func:`saturation_vapour_pressure`.

"""
# by masking upfront we avoid RuntimeWarning when calling log() in
# the computation of td when q is very small
if out is not None:
q = np.asarray(q)
mask = q < eps
q[mask] = np.nan

es = vapour_pressure_from_specific_humidity(q, p)
td = temperature_from_saturation_vapour_pressure(es)

if out is not None and not np.isnan(out):
td = np.asarray(td)
td[mask] = out

return td
return temperature_from_saturation_vapour_pressure(vapour_pressure_from_specific_humidity(q, p))


def virtual_temperature(t, q):
Expand Down
80 changes: 13 additions & 67 deletions tests/thermo/test_thermo_array.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -334,24 +334,20 @@ def test_temperature_from_saturation_vapour_pressure_1():


@pytest.mark.parametrize(
"es,kwargs,expected_values",
"es,expected_values",
[
(4.2, {}, 219.7796336743947),
([4.2, 0, 0.001, np.nan], {"eps": 1e-2, "out": 100}, [219.7796336743947, 100, 100, np.nan]),
([4.2, 0, 0.001, np.nan], {"eps": 1e-2, "out": np.nan}, [219.7796336743947, np.nan, np.nan, np.nan]),
(0, {}, np.nan),
(0.001, {"eps": 1e-2, "out": 100}, 100.0),
(0.001, {"eps": 1e-2, "out": np.nan}, np.nan),
(4.2, 219.7796336743947),
(0, np.nan),
],
)
def test_temperature_from_saturation_vapour_pressure_2(es, kwargs, expected_values):
def test_temperature_from_saturation_vapour_pressure_2(es, expected_values):

multi = isinstance(es, list)
if multi:
es = np.array(es)
expected_values = np.array(expected_values)

t = thermo.array.temperature_from_saturation_vapour_pressure(es, **kwargs)
t = thermo.array.temperature_from_saturation_vapour_pressure(es)
if multi:
np.testing.assert_allclose(t, expected_values, equal_nan=True)
else:
Expand Down Expand Up @@ -448,7 +444,7 @@ def test_specific_humidity_from_relative_humidity():


@pytest.mark.parametrize(
"t,r,kwargs,expected_values",
"t,r,expected_values",
[
(
[20.0, 20, 0, 35, 5, -15, 25, 25],
Expand All @@ -462,37 +458,12 @@ def test_specific_humidity_from_relative_humidity():
15.4779832381,
0,
],
{},
[20.0, 10, -10, 32, -15, -24, -3, np.nan],
),
(
[20.0, 20.0, 20.0],
[
52.5224541378,
0.0,
0.000001,
],
{"eps": 1e-3, "out": thermo.array.celsius_to_kelvin(100)},
[10, 100, 100],
),
(
[20.0, 20.0, 20.0],
[
52.5224541378,
0.0,
0.000001,
],
{"eps": 1e-3, "out": np.nan},
[10, np.nan, np.nan],
),
(20.0, 52.5224541378, {}, 10.0),
(20.0, 0.0, {"eps": 1e-3, "out": thermo.array.celsius_to_kelvin(100)}, 100),
(20.0, 0.000001, {"eps": 1e-3, "out": thermo.array.celsius_to_kelvin(100)}, 100),
(20.0, 0.0, {"eps": 1e-3, "out": np.nan}, np.nan),
(20, 0.000001, {"eps": 1e-3, "out": np.nan}, np.nan),
(20.0, 52.5224541378, 10.0),
],
)
def test_dewpoint_from_relative_humidity(t, r, kwargs, expected_values):
def test_dewpoint_from_relative_humidity(t, r, expected_values):
# reference was tested with an online relhum calculator at:
# https://bmcnoldy.rsmas.miami.edu/Humidity.html

Expand All @@ -505,50 +476,25 @@ def test_dewpoint_from_relative_humidity(t, r, kwargs, expected_values):
t = thermo.array.celsius_to_kelvin(t)
v_ref = thermo.array.celsius_to_kelvin(expected_values)

td = thermo.array.dewpoint_from_relative_humidity(t, r, **kwargs)
td = thermo.array.dewpoint_from_relative_humidity(t, r)
if multi:
assert np.allclose(td, v_ref, equal_nan=True)
else:
assert np.isclose(td, v_ref, equal_nan=True)


@pytest.mark.parametrize(
"q,p,kwargs,expected_values",
"q,p,expected_values",
[
(
[0.0169461501, 0.0155840075, 0.0134912382, 0.0083409720, 0.0057268584, 0.0025150791, 0],
[967.5085, 936.3775, 872.248, 756.1647, 649.157, 422.4207, 422.4207],
{},
[21.78907, 19.90885, 16.50236, 7.104064, -0.3548709, -16.37916, np.nan],
),
(
[
0.0169461501,
0.0,
0.000001,
],
[967.5085, 967.5085, 967.5085],
{"eps": 1e-3, "out": thermo.array.celsius_to_kelvin(100)},
[21.78907, 100, 100],
),
(
[
0.0169461501,
0.0,
0.000001,
],
[967.5085, 967.5085, 967.5085],
{"eps": 1e-3, "out": np.nan},
[21.78907, np.nan, np.nan],
),
(0.0169461501, 967.508, {}, 21.78907),
(0.0, 967.5085, {"eps": 1e-3, "out": thermo.array.celsius_to_kelvin(100)}, 100),
(0.000001, 967.5085, {"eps": 1e-3, "out": thermo.array.celsius_to_kelvin(100)}, 100),
(0.0, 967.5085, {"eps": 1e-3, "out": np.nan}, np.nan),
(0.000001, 967.5085, {"eps": 1e-3, "out": np.nan}, np.nan),
(0.0169461501, 967.508, 21.78907),
],
)
def test_dewpoint_from_specific_humidity(q, p, kwargs, expected_values):
def test_dewpoint_from_specific_humidity(q, p, expected_values):
multi = isinstance(q, list)
if multi:
q = np.array(q)
Expand All @@ -558,7 +504,7 @@ def test_dewpoint_from_specific_humidity(q, p, kwargs, expected_values):
p = p * 100.0
v_ref = thermo.array.celsius_to_kelvin(expected_values)

td = thermo.array.dewpoint_from_specific_humidity(q, p, **kwargs)
td = thermo.array.dewpoint_from_specific_humidity(q, p)
if multi:
assert np.allclose(td, v_ref, equal_nan=True)
else:
Expand Down
Loading