From 9e232eee2faf865f8ed4ed540c961cadbe1dea0e Mon Sep 17 00:00:00 2001 From: daflack Date: Fri, 2 Jan 2026 11:37:50 +0000 Subject: [PATCH 01/91] Basic conversion operators for humidity, temperature and pressure Fixes #1852 --- src/CSET/operators/humidity.py | 35 ++++++++++++++++++++++++++++++++++ 1 file changed, 35 insertions(+) create mode 100644 src/CSET/operators/humidity.py diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py new file mode 100644 index 000000000..c54590061 --- /dev/null +++ b/src/CSET/operators/humidity.py @@ -0,0 +1,35 @@ +# © Crown copyright, Met Office (2022-2026) and CSET contributors. +# +# Licensed under the Apache License, Version 2.0 (the "License"); +# you may not use this file except in compliance with the License. +# You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software +# distributed under the License is distributed on an "AS IS" BASIS, +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +# See the License for the specific language governing permissions and +# limitations under the License. + +"""Operators for humidity conversions.""" + +import iris.cube + +from CSET._common import iter_maybe + + +def specific_humidity_to_mixing_ratio( + cubes: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Convert specific humidity to mixing ratio.""" + w = iris.cube.CubeList([]) + for cube in iter_maybe(cubes): + mr = cube.copy() + mr = cube / (1 - cube) + mr.rename("mixing_ratio") + w.append(mr) + if len(w) == 1: + return w[0] + else: + return w From c072a26c2942543fff238f7f4ddf68d06ac6e3c7 Mon Sep 17 00:00:00 2001 From: daflack Date: Fri, 2 Jan 2026 11:54:22 +0000 Subject: [PATCH 02/91] Add mixing ratio to specific humidity operator --- src/CSET/operators/humidity.py | 16 ++++++++++++++++ 1 file changed, 16 insertions(+) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index c54590061..8f37371fb 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -33,3 +33,19 @@ def specific_humidity_to_mixing_ratio( return w[0] else: return w + + +def mixing_ratio_to_specific_humidity( + cubes: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Convert mixing ratio to specific humidity.""" + q = iris.cube.CubeList([]) + for cube in iter_maybe(cubes): + sh = cube.copy() + sh = cube / (1 + cube) + sh.rename("specific_humidity") + q.append(sh) + if len(q) == 1: + return q[0] + else: + return q From 191ac7cfd090745fb4470de874192663a6a1eb02 Mon Sep 17 00:00:00 2001 From: daflack Date: Fri, 2 Jan 2026 12:43:54 +0000 Subject: [PATCH 03/91] Add standard atmospheric constants --- src/CSET/operators/atmospheric_constants.py | 38 +++++++++++++++++++++ 1 file changed, 38 insertions(+) create mode 100644 src/CSET/operators/atmospheric_constants.py diff --git a/src/CSET/operators/atmospheric_constants.py b/src/CSET/operators/atmospheric_constants.py new file mode 100644 index 000000000..2815c642a --- /dev/null +++ b/src/CSET/operators/atmospheric_constants.py @@ -0,0 +1,38 @@ +# © Crown copyright, Met Office (2022-2026) and CSET contributors. +# +# Licensed under the Apache License, Version 2.0 (the "License"); +# you may not use this file except in compliance with the License. +# You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software +# distributed under the License is distributed on an "AS IS" BASIS, +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +# See the License for the specific language governing permissions and +# limitations under the License. +"""Constants for the atmosphere.""" + +# Reference pressure. +P0 = 1000.0 # hPa + +# Specific gas constant for dry air. +RD = 287.0 + +# Specific gas constant for water vapour. +RV = 461.0 + +# Specific heat capacity for dry air. +CPD = 1005.7 + +# Latent heat of vaporization. +LV = 2.501e6 + +# Reference vapour pressure. +E0 = 6.1078 # hPa + +# Reference temperature. +T0 = 273.15 # K + +# Ratio between mixing ratio of dry and moist air. +EPSILON = 0.622 From 2dc4795434caa97d2b926792647199ddeec42915 Mon Sep 17 00:00:00 2001 From: daflack Date: Fri, 2 Jan 2026 12:45:21 +0000 Subject: [PATCH 04/91] Add vapour pressure conversion --- src/CSET/operators/pressure.py | 37 ++++++++++++++++++++++++++++++++++ 1 file changed, 37 insertions(+) create mode 100644 src/CSET/operators/pressure.py diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py new file mode 100644 index 000000000..9d1bdc58b --- /dev/null +++ b/src/CSET/operators/pressure.py @@ -0,0 +1,37 @@ +# © Crown copyright, Met Office (2022-2026) and CSET contributors. +# +# Licensed under the Apache License, Version 2.0 (the "License"); +# you may not use this file except in compliance with the License. +# You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software +# distributed under the License is distributed on an "AS IS" BASIS, +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +# See the License for the specific language governing permissions and +# limitations under the License. + +"""Operators for pressure conversions.""" + +import iris.cube +import numpy as np + +from CSET._common import iter_maybe +from CSET.operators.atmospheric_constants import E0 + + +def vapour_pressure( + temperature: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Calculate the vapour pressure of the atmosphere.""" + vapour_pressure = iris.cube.CubeList([]) + for T in iter_maybe(temperature): + es = T.copy() + exponent = 17.27 * (T - 273.16) / (T - 35.86) + es.data[:] = E0 * np.exp(exponent.core_data()) + vapour_pressure.append(es) + if len(vapour_pressure) == 1: + return vapour_pressure[0] + else: + return vapour_pressure From f6ca91294170c2dd02678ab2a2196623d0f886d4 Mon Sep 17 00:00:00 2001 From: daflack Date: Fri, 2 Jan 2026 13:17:08 +0000 Subject: [PATCH 05/91] Add vapour pressure if dewpoint unknown --- src/CSET/operators/pressure.py | 22 ++++++++++++++++++++++ 1 file changed, 22 insertions(+) diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index 9d1bdc58b..92f996836 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -30,8 +30,30 @@ def vapour_pressure( es = T.copy() exponent = 17.27 * (T - 273.16) / (T - 35.86) es.data[:] = E0 * np.exp(exponent.core_data()) + es.units = "hPa" + es.rename("vapour_pressure") vapour_pressure.append(es) if len(vapour_pressure) == 1: return vapour_pressure[0] else: return vapour_pressure + + +def vapour_pressure_if_Td_unknown( + temperature: iris.cube.Cube | iris.cube.CubeList, + relative_humidity: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Calculate the vapour pressure using RH.""" + vapour_pressure = iris.cube.CubeList([]) + for T, RH in zip( + iter_maybe(temperature), iter_maybe(relative_humidity), strict=True + ): + if RH.units == "%": + RH /= 100.0 + RH.units = "1" + vp = vapour_pressure(T) * RH + vapour_pressure.append(vp) + if len(vapour_pressure) == 1: + return vapour_pressure[0] + else: + return vapour_pressure From b49f776e3351584ba0963ab626cd66c6b781bc79 Mon Sep 17 00:00:00 2001 From: daflack Date: Fri, 2 Jan 2026 13:21:49 +0000 Subject: [PATCH 06/91] Update naming convections for vapour pressures --- src/CSET/operators/pressure.py | 22 +++++++++++----------- 1 file changed, 11 insertions(+), 11 deletions(-) diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index 92f996836..7e59f0dfa 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -25,26 +25,26 @@ def vapour_pressure( temperature: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: """Calculate the vapour pressure of the atmosphere.""" - vapour_pressure = iris.cube.CubeList([]) + v_pressure = iris.cube.CubeList([]) for T in iter_maybe(temperature): es = T.copy() exponent = 17.27 * (T - 273.16) / (T - 35.86) es.data[:] = E0 * np.exp(exponent.core_data()) es.units = "hPa" es.rename("vapour_pressure") - vapour_pressure.append(es) - if len(vapour_pressure) == 1: - return vapour_pressure[0] + v_pressure.append(es) + if len(v_pressure) == 1: + return v_pressure[0] else: - return vapour_pressure + return v_pressure -def vapour_pressure_if_Td_unknown( +def vapour_pressure_if_dewpoint_unknown( temperature: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: """Calculate the vapour pressure using RH.""" - vapour_pressure = iris.cube.CubeList([]) + v_pressure = iris.cube.CubeList([]) for T, RH in zip( iter_maybe(temperature), iter_maybe(relative_humidity), strict=True ): @@ -52,8 +52,8 @@ def vapour_pressure_if_Td_unknown( RH /= 100.0 RH.units = "1" vp = vapour_pressure(T) * RH - vapour_pressure.append(vp) - if len(vapour_pressure) == 1: - return vapour_pressure[0] + v_pressure.append(vp) + if len(v_pressure) == 1: + return v_pressure[0] else: - return vapour_pressure + return v_pressure From cef08a9c45ef88071b29d580cd667f9e6d2bd971 Mon Sep 17 00:00:00 2001 From: daflack Date: Fri, 2 Jan 2026 13:32:53 +0000 Subject: [PATCH 07/91] Calculate saturation mixing ratio and specific humidity --- src/CSET/operators/humidity.py | 36 ++++++++++++++++++++++++++++++++++ 1 file changed, 36 insertions(+) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index 8f37371fb..e39b76345 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -17,6 +17,8 @@ import iris.cube from CSET._common import iter_maybe +from CSET.operators.constants import EPSILON +from CSET.operators.pressure import vapour_pressure def specific_humidity_to_mixing_ratio( @@ -49,3 +51,37 @@ def mixing_ratio_to_specific_humidity( return q[0] else: return q + + +def saturation_mixing_ratio( + temperature: iris.cube.Cube | iris.cube.CubeList, + pressure: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Calculate saturation mixing ratio.""" + w = iris.cube.CubeList([]) + for T, P in zip(iter_maybe(temperature), iter_maybe(pressure), strict=True): + mr = (EPSILON * vapour_pressure(T)) / ((P / 100.0) - vapour_pressure(T)) + mr.units("kg/kg") + mr.rename("mixing_ratio") + w.append(mr) + if len(w) == 1: + return w[0] + else: + return w + + +def saturation_specific_humidity( + temperature: iris.cube.Cube | iris.cube.CubeList, + pressure: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Calculate saturation specific humidity.""" + q = iris.cube.CubeList([]) + for T, P in zip(iter_maybe(temperature), iter_maybe(pressure), strict=True): + sh = (EPSILON * vapour_pressure(T)) / (P / 100.0) + sh.units("kg/kg") + sh.rename("specific_humidity") + q.append(sh) + if len(q) == 1: + return q[0] + else: + return q From 9cf07b257f84ae45aa41288430cde19c44c330b6 Mon Sep 17 00:00:00 2001 From: daflack Date: Fri, 2 Jan 2026 13:34:48 +0000 Subject: [PATCH 08/91] Rename vapour_pressure_if_dewpoint_unknown to vapour_pressure_from_RH --- src/CSET/operators/pressure.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index 7e59f0dfa..f25c33c1e 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -39,7 +39,7 @@ def vapour_pressure( return v_pressure -def vapour_pressure_if_dewpoint_unknown( +def vapour_pressure_from_RH( temperature: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: From 80b08f669b6201c7e7edebe3bca6e2899519c4f8 Mon Sep 17 00:00:00 2001 From: daflack Date: Fri, 2 Jan 2026 14:02:32 +0000 Subject: [PATCH 09/91] Add mixing ratio from relative humidity conversion --- src/CSET/operators/humidity.py | 24 ++++++++++++++++++++++++ 1 file changed, 24 insertions(+) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index e39b76345..d5272937b 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -85,3 +85,27 @@ def saturation_specific_humidity( return q[0] else: return q + + +def mixing_ratio_from_RH( + temperature: iris.cube.Cube | iris.cube.CubeList, + pressure: iris.cube.Cube | iris.cube.CubeList, + relative_humidity: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Calculate the mixing ratio from RH.""" + w = iris.cube.CubeList([]) + for T, P, RH in zip( + iter_maybe(temperature), + iter_maybe(pressure), + iter_maybe(relative_humidity), + strict=True, + ): + if RH.units == "%": + RH /= 100.0 + RH.units = "1" + mr = saturation_mixing_ratio(T, P) * RH + w.append(mr) + if len(w) == 1: + return w[0] + else: + return w From 3d87c7a2257a47b85c09fdcf9b94a51d30370bab Mon Sep 17 00:00:00 2001 From: daflack Date: Fri, 2 Jan 2026 14:06:45 +0000 Subject: [PATCH 10/91] Add specific humidity from RH conversion --- src/CSET/operators/humidity.py | 24 ++++++++++++++++++++++++ 1 file changed, 24 insertions(+) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index d5272937b..b48f452d6 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -109,3 +109,27 @@ def mixing_ratio_from_RH( return w[0] else: return w + + +def specific_humidity_from_RH( + temperature: iris.cube.Cube | iris.cube.CubeList, + pressure: iris.cube.Cube | iris.cube.CubeList, + relative_humidity: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Calculate the mixing ratio from RH.""" + q = iris.cube.CubeList([]) + for T, P, RH in zip( + iter_maybe(temperature), + iter_maybe(pressure), + iter_maybe(relative_humidity), + strict=True, + ): + if RH.units == "%": + RH /= 100.0 + RH.units = "1" + sh = saturation_specific_humidity(T, P) * RH + q.append(sh) + if len(q) == 1: + return q[0] + else: + return q From bff14ae19d93625405f28ab83e31d0602d106d02 Mon Sep 17 00:00:00 2001 From: daflack Date: Fri, 2 Jan 2026 14:19:29 +0000 Subject: [PATCH 11/91] Add relative humidity conversions from mixing ratio and specific humidity --- src/CSET/operators/humidity.py | 44 ++++++++++++++++++++++++++++++++++ 1 file changed, 44 insertions(+) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index b48f452d6..d83497649 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -133,3 +133,47 @@ def specific_humidity_from_RH( return q[0] else: return q + + +def relative_humidity_from_mixing_ratio( + mixing_ratio: iris.cube.Cube | iris.cube.CubeList, + temperature: iris.cube.Cube | iris.cube.CubeList, + pressure: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Convert mixing ratio to relative humidity.""" + RH = iris.cube.CubeList([]) + for W, T, P in zip( + iter_maybe(mixing_ratio), + iter_maybe(temperature), + iter_maybe(pressure), + strict=True, + ): + rel_h = W / saturation_mixing_ratio(T, P) + rel_h.rename("relative_humidity") + RH.append(rel_h) + if len(RH) == 1: + return RH[0] + else: + return RH + + +def relative_humidity_from_specific_humidity( + specific_humidity: iris.cube.Cube | iris.cube.CubeList, + temperature: iris.cube.Cube | iris.cube.CubeList, + pressure: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Convert specific humidity to relative humidity.""" + RH = iris.cube.CubeList([]) + for Q, T, P in zip( + iter_maybe(specific_humidity), + iter_maybe(temperature), + iter_maybe(pressure), + strict=True, + ): + rel_h = Q / saturation_specific_humidity(T, P) + rel_h.rename("relative_humidity") + RH.append(rel_h) + if len(RH) == 1: + return RH[0] + else: + return RH From 7be43e7484d3c5be49e01d6f9fd1f973a54b1e4d Mon Sep 17 00:00:00 2001 From: daflack Date: Fri, 2 Jan 2026 14:32:46 +0000 Subject: [PATCH 12/91] Add dewpoint temperature calculation --- src/CSET/operators/temperature.py | 47 +++++++++++++++++++++++++++++++ 1 file changed, 47 insertions(+) create mode 100644 src/CSET/operators/temperature.py diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py new file mode 100644 index 000000000..ed69f6d8c --- /dev/null +++ b/src/CSET/operators/temperature.py @@ -0,0 +1,47 @@ +# © Crown copyright, Met Office (2022-2026) and CSET contributors. +# +# Licensed under the Apache License, Version 2.0 (the "License"); +# you may not use this file except in compliance with the License. +# You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software +# distributed under the License is distributed on an "AS IS" BASIS, +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +# See the License for the specific language governing permissions and +# limitations under the License. + +"""Operators for temperature conversions.""" + +import iris.cube +import numpy as np + +from CSET._common import iter_maybe +from CSET.operators.constants import T0 +from CSET.operators.pressure import vapour_pressure + + +def dewpoint( + temperature: iris.cube.Cube | iris.cube.CubeList, + relative_humidity: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Calculate the dewpoint temperature.""" + Td = iris.cube.CubeList([]) + for T, RH in zip( + iter_maybe(temperature), iter_maybe(relative_humidity), strict=True + ): + vp = vapour_pressure(T, RH) + td = vp.copy() + td.data = (243.5 * np.log(vp.core_data()) - 440.8) / ( + 19.48 - np.log(vp.core_data()) + ) + td.data[td.data - T0 < -35.0] = np.nan + td.data[td.data - T0 > 35.0] = np.nan + td.units("K") + td.rename("dewpoint temperature") + Td.append(td) + if len(Td) == 1: + return Td[0] + else: + return Td From 678d89c009395909a97b2fdbaf14b3923abcff1d Mon Sep 17 00:00:00 2001 From: daflack Date: Fri, 2 Jan 2026 14:55:00 +0000 Subject: [PATCH 13/91] Add virtual temperature calculation --- src/CSET/operators/temperature.py | 20 ++++++++++++++++++-- 1 file changed, 18 insertions(+), 2 deletions(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index ed69f6d8c..7d88e9f47 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -18,11 +18,11 @@ import numpy as np from CSET._common import iter_maybe -from CSET.operators.constants import T0 +from CSET.operators.constants import EPSILON, T0 from CSET.operators.pressure import vapour_pressure -def dewpoint( +def dewpoint_temperature( temperature: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: @@ -45,3 +45,19 @@ def dewpoint( return Td[0] else: return Td + + +def virtual_temperature( + temperature: iris.cube.Cube | iris.cube.CubeList, + mixing_ratio: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Calculate the virtual temperature.""" + Tv = iris.cube.CubeList([]) + for T, W in zip(iter_maybe(temperature), iter_maybe(mixing_ratio), strict=True): + virT = T * ((W + EPSILON) / (EPSILON * (1 + W))) + virT.rename("virtual_temperature") + Tv.append(virT) + if len(Tv) == 1: + return Tv[0] + else: + return Tv From d258924ba6f1740708454e16f1f7c3cb29f0a1a2 Mon Sep 17 00:00:00 2001 From: daflack Date: Fri, 2 Jan 2026 14:59:17 +0000 Subject: [PATCH 14/91] Update atmospheric constants with kappa --- src/CSET/operators/atmospheric_constants.py | 3 +++ 1 file changed, 3 insertions(+) diff --git a/src/CSET/operators/atmospheric_constants.py b/src/CSET/operators/atmospheric_constants.py index 2815c642a..85390b730 100644 --- a/src/CSET/operators/atmospheric_constants.py +++ b/src/CSET/operators/atmospheric_constants.py @@ -36,3 +36,6 @@ # Ratio between mixing ratio of dry and moist air. EPSILON = 0.622 + +# Ratio between specific gas constant and specific heat capacity. +KAPPA = RD / CPD From a92510320b1066851e17488ef61e220cf609791c Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 5 Jan 2026 09:34:57 +0000 Subject: [PATCH 15/91] Add potential temperature and exner pressure convertors --- src/CSET/operators/pressure.py | 19 ++++++++++++++++++- src/CSET/operators/temperature.py | 19 ++++++++++++++++++- 2 files changed, 36 insertions(+), 2 deletions(-) diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index f25c33c1e..f18d048cc 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -18,7 +18,7 @@ import numpy as np from CSET._common import iter_maybe -from CSET.operators.atmospheric_constants import E0 +from CSET.operators.atmospheric_constants import E0, KAPPA, P0 def vapour_pressure( @@ -57,3 +57,20 @@ def vapour_pressure_from_RH( return v_pressure[0] else: return v_pressure + + +def exner_pressure( + pressure: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Calculate the exner pressure.""" + pi = iris.cube.CubeList([]) + for P in iter_maybe(pressure): + PI = P.copy() + PI.data[:] = (P.core_data() / P0) ** KAPPA + PI.rename("exner_pressure") + PI.units("1") + pi.append(PI) + if len(pi) == 1: + return pi[0] + else: + return pi diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 7d88e9f47..50fb8bb2d 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -19,7 +19,7 @@ from CSET._common import iter_maybe from CSET.operators.constants import EPSILON, T0 -from CSET.operators.pressure import vapour_pressure +from CSET.operators.pressure import exner_pressure, vapour_pressure def dewpoint_temperature( @@ -61,3 +61,20 @@ def virtual_temperature( return Tv[0] else: return Tv + + +def potential_temperature( + temperature: iris.cube.Cube | iris.cube.CubeList, + pressure: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Calculate the potenital temperature.""" + theta = iris.cube.CubeList([]) + for T, P in zip(iter_maybe(temperature), iter_maybe(pressure), strict=True): + TH = T / exner_pressure(P) + TH.rename("potential_temperature") + TH.units("K") + theta.append(TH) + if len(theta) == 1: + return theta[0] + else: + return theta From 68ea4385e7a20c940e205fc4dd44faf6fd8abd22 Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 5 Jan 2026 09:58:48 +0000 Subject: [PATCH 16/91] Adds virtual potential temperature convertor --- src/CSET/operators/temperature.py | 23 ++++++++++++++++++++++- 1 file changed, 22 insertions(+), 1 deletion(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 50fb8bb2d..49ce02e8a 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -72,9 +72,30 @@ def potential_temperature( for T, P in zip(iter_maybe(temperature), iter_maybe(pressure), strict=True): TH = T / exner_pressure(P) TH.rename("potential_temperature") - TH.units("K") theta.append(TH) if len(theta) == 1: return theta[0] else: return theta + + +def virtual_potential_temperature( + temperature: iris.cube.Cube | iris.cube.CubeList, + mixing_ratio: iris.cube.Cube | iris.cube.CubeList, + pressure: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Calculate the virtual potential temperature.""" + theta_v = iris.cube.CubeList([]) + for T, W, P in zip( + iter_maybe(temperature), + iter_maybe(mixing_ratio), + iter_maybe(pressure), + strict=True, + ): + TH_V = virtual_temperature(T, W) / exner_pressure(P) + TH_V.rename("virtual_potential_temperature") + theta_v.append(TH_V) + if len(theta_v) == 1: + return theta_v[0] + else: + return theta_v From 4084822607fe12fc1c59a91c1abb1d738a6c20dc Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 5 Jan 2026 11:03:39 +0000 Subject: [PATCH 17/91] Adds equivalent potential temperature conversion and fixes unit assignment --- src/CSET/operators/humidity.py | 4 ++-- src/CSET/operators/pressure.py | 2 +- src/CSET/operators/temperature.py | 37 +++++++++++++++++++++++++++++-- 3 files changed, 38 insertions(+), 5 deletions(-) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index d83497649..6af3917f2 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -61,7 +61,7 @@ def saturation_mixing_ratio( w = iris.cube.CubeList([]) for T, P in zip(iter_maybe(temperature), iter_maybe(pressure), strict=True): mr = (EPSILON * vapour_pressure(T)) / ((P / 100.0) - vapour_pressure(T)) - mr.units("kg/kg") + mr.units = "kg/kg" mr.rename("mixing_ratio") w.append(mr) if len(w) == 1: @@ -78,7 +78,7 @@ def saturation_specific_humidity( q = iris.cube.CubeList([]) for T, P in zip(iter_maybe(temperature), iter_maybe(pressure), strict=True): sh = (EPSILON * vapour_pressure(T)) / (P / 100.0) - sh.units("kg/kg") + sh.units = "kg/kg" sh.rename("specific_humidity") q.append(sh) if len(q) == 1: diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index f18d048cc..a0a30832c 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -68,7 +68,7 @@ def exner_pressure( PI = P.copy() PI.data[:] = (P.core_data() / P0) ** KAPPA PI.rename("exner_pressure") - PI.units("1") + PI.units = "1" pi.append(PI) if len(pi) == 1: return pi[0] diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 49ce02e8a..3cb87e502 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -18,7 +18,8 @@ import numpy as np from CSET._common import iter_maybe -from CSET.operators.constants import EPSILON, T0 +from CSET.operators.constants import CPD, EPSILON, LV, RV, T0 +from CSET.operators.humidity import mixing_ratio_from_RH from CSET.operators.pressure import exner_pressure, vapour_pressure @@ -38,7 +39,7 @@ def dewpoint_temperature( ) td.data[td.data - T0 < -35.0] = np.nan td.data[td.data - T0 > 35.0] = np.nan - td.units("K") + td.units = "K" td.rename("dewpoint temperature") Td.append(td) if len(Td) == 1: @@ -99,3 +100,35 @@ def virtual_potential_temperature( return theta_v[0] else: return theta_v + + +def equivalent_potential_temperature( + temperature: iris.cube.Cube | iris.cube.CubeList, + relative_humidity: iris.cube.Cube | iris.cube.CubeList, + pressure: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Calculate the equivalent potenital temperature.""" + theta_e = iris.cube.CubeList([]) + for T, RH, P in zip( + iter_maybe(temperature), + iter_maybe(relative_humidity), + iter_maybe(pressure), + strict=True, + ): + if RH.units == "%": + RH /= 100.0 + RH.units = "1" + theta = potential_temperature(T, P) + w = mixing_ratio_from_RH(T, P, RH) + second_term_power = -(w * RV) / CPD + second_term = RH.core_data() ** second_term_power.core_data() + third_term_power = LV * w / (CPD * T) + third_term = np.exp(third_term_power.core_data()) + TH_E = theta * second_term * third_term + TH_E.rename("equivalent_potential_temperature") + TH_E.units = "K" + theta_e.append(TH_E) + if len(theta_e) == 1: + return theta_e[0] + else: + return theta_e From d9f1d0c38d28cdaf05ed10009c9b20527a7620d0 Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 5 Jan 2026 11:38:37 +0000 Subject: [PATCH 18/91] Remove if loop for RH and switch to convert units --- src/CSET/operators/humidity.py | 9 +++------ src/CSET/operators/pressure.py | 5 ++--- src/CSET/operators/temperature.py | 5 ++--- 3 files changed, 7 insertions(+), 12 deletions(-) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index 6af3917f2..d15ed319e 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -18,6 +18,7 @@ from CSET._common import iter_maybe from CSET.operators.constants import EPSILON +from CSET.operators.misc import convert_units from CSET.operators.pressure import vapour_pressure @@ -100,9 +101,7 @@ def mixing_ratio_from_RH( iter_maybe(relative_humidity), strict=True, ): - if RH.units == "%": - RH /= 100.0 - RH.units = "1" + RH = convert_units(RH, "1") mr = saturation_mixing_ratio(T, P) * RH w.append(mr) if len(w) == 1: @@ -124,9 +123,7 @@ def specific_humidity_from_RH( iter_maybe(relative_humidity), strict=True, ): - if RH.units == "%": - RH /= 100.0 - RH.units = "1" + RH = convert_units(RH, "1") sh = saturation_specific_humidity(T, P) * RH q.append(sh) if len(q) == 1: diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index a0a30832c..5f5f92c4a 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -19,6 +19,7 @@ from CSET._common import iter_maybe from CSET.operators.atmospheric_constants import E0, KAPPA, P0 +from CSET.operators.misc import convert_units def vapour_pressure( @@ -48,9 +49,7 @@ def vapour_pressure_from_RH( for T, RH in zip( iter_maybe(temperature), iter_maybe(relative_humidity), strict=True ): - if RH.units == "%": - RH /= 100.0 - RH.units = "1" + RH = convert_units(RH, "1") vp = vapour_pressure(T) * RH v_pressure.append(vp) if len(v_pressure) == 1: diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 3cb87e502..738d264c5 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -20,6 +20,7 @@ from CSET._common import iter_maybe from CSET.operators.constants import CPD, EPSILON, LV, RV, T0 from CSET.operators.humidity import mixing_ratio_from_RH +from CSET.operators.misc import convert_units from CSET.operators.pressure import exner_pressure, vapour_pressure @@ -115,9 +116,7 @@ def equivalent_potential_temperature( iter_maybe(pressure), strict=True, ): - if RH.units == "%": - RH /= 100.0 - RH.units = "1" + RH = convert_units(RH, "1") theta = potential_temperature(T, P) w = mixing_ratio_from_RH(T, P, RH) second_term_power = -(w * RV) / CPD From ae72d4e03a401ad2cf6032e40b52a3b01b96f6ae Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 5 Jan 2026 11:57:55 +0000 Subject: [PATCH 19/91] Adds wet-bulb temperature convertor --- src/CSET/operators/temperature.py | 29 +++++++++++++++++++++++++++++ 1 file changed, 29 insertions(+) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 738d264c5..9117e59ca 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -65,6 +65,35 @@ def virtual_temperature( return Tv +def wet_bulb_temperature( + temperature: iris.cube.Cube | iris.cube.CubeList, + relative_humidity: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Calculate the wet-bulb temperature.""" + Tw = iris.cube.CubeList([]) + for T, RH in zip( + iter_maybe(temperature), iter_maybe(relative_humidity), strict=True + ): + RH = convert_units(RH, "1") + T = convert_units(T, "Celsius") + wetT = ( + T * np.arctan(0.151977 * (RH.core_data() + 8.313659) ** 0.5) + + np.arctan(T.core_data() + RH.core_data()) + - np.arctan(RH.core_data() - 1.676331) + + 0.00391838 + * (RH.core_data()) ** (3.0 / 2.0) + * np.arctan(0.023101 * RH.core_data()) + - 4.686035 + ) + wetT.rename("wet_bulb_temperature") + wetT = convert_units(wetT, "K") + Tw.append(wetT) + if len(Tw) == 1: + return Tw[0] + else: + return Tw + + def potential_temperature( temperature: iris.cube.Cube | iris.cube.CubeList, pressure: iris.cube.Cube | iris.cube.CubeList, From 4fe1e39621d7b21053683ea910d88e9673f1f8f5 Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 5 Jan 2026 12:20:31 +0000 Subject: [PATCH 20/91] Adds relevant references for temperature calculations where required --- src/CSET/operators/temperature.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 9117e59ca..ea53b091d 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -28,7 +28,7 @@ def dewpoint_temperature( temperature: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate the dewpoint temperature.""" + """Calculate the dewpoint temperature following Bolton (1980).""" Td = iris.cube.CubeList([]) for T, RH in zip( iter_maybe(temperature), iter_maybe(relative_humidity), strict=True @@ -69,7 +69,7 @@ def wet_bulb_temperature( temperature: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate the wet-bulb temperature.""" + """Calculate the wet-bulb temperature following Stull (2011).""" Tw = iris.cube.CubeList([]) for T, RH in zip( iter_maybe(temperature), iter_maybe(relative_humidity), strict=True From e7d7297177ee978887cc0b90087316530187bbdb Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 5 Jan 2026 12:38:28 +0000 Subject: [PATCH 21/91] Adds wet-bulb potential temperature convertor --- src/CSET/operators/temperature.py | 40 +++++++++++++++++++++++++++++++ 1 file changed, 40 insertions(+) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index ea53b091d..3d8b3dc4d 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -160,3 +160,43 @@ def equivalent_potential_temperature( return theta_e[0] else: return theta_e + + +def wet_bulb_potential_temperature( + temperature: iris.cube.Cube | iris.cube.CubeList, + relative_humidity: iris.cube.Cube | iris.cube.CubeList, + pressure: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Calculate wet-bulb potential temperature after Davies-Jones (2008).""" + theta_w = iris.cube.CubeList([]) + for T, RH, P in zip( + iter_maybe(temperature), + iter_maybe(relative_humidity), + iter_maybe(pressure), + strict=True, + ): + TH_E = equivalent_potential_temperature(T, P, RH) + X = TH_E / T0 + X.units = "1" + A0 = 7.101574 + A1 = -20.68208 + A2 = 16.11182 + A3 = 2.574631 + A4 = -5.205688 + B1 = -3.552497 + B2 = 3.781782 + B3 = -0.6899655 + B4 = -0.5929340 + exponent = (A0 + A1 * X + A2 * X**2 + A3 * X**3 + A4 * X**4) / ( + 1.0 + B1 * X + B2 * X**2 + B3 * X**3 + B4 * X**4 + ) + TH_W = TH_E.copy() + TH_W.data[:] = TH_E.core_data() - np.exp(exponent.core_data()) + TH_W.rename("wet_bulb_potential_temperature") + TH_W.data[TH_W.data - T0 < -30.0] = np.nan + TH_W.data[TH_W.data - T0 > 50.0] = np.nan + theta_w.append(TH_W) + if len(theta_w) == 1: + return theta_w[0] + else: + return theta_w From dfe7657656d527dd71697e69586636ec7419e741 Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 5 Jan 2026 13:36:30 +0000 Subject: [PATCH 22/91] Adds saturation equivalent potential temperature convertor --- src/CSET/operators/temperature.py | 27 ++++++++++++++++++++++++++- 1 file changed, 26 insertions(+), 1 deletion(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 3d8b3dc4d..087ff24ca 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -19,7 +19,7 @@ from CSET._common import iter_maybe from CSET.operators.constants import CPD, EPSILON, LV, RV, T0 -from CSET.operators.humidity import mixing_ratio_from_RH +from CSET.operators.humidity import mixing_ratio_from_RH, saturation_mixing_ratio from CSET.operators.misc import convert_units from CSET.operators.pressure import exner_pressure, vapour_pressure @@ -200,3 +200,28 @@ def wet_bulb_potential_temperature( return theta_w[0] else: return theta_w + + +def saturation_equivalent_potential_temperature( + temperature: iris.cube.Cube | iris.cube.CubeList, + pressure: iris.cube.Cube | iris.cube.CubeList, +) -> iris.cube.Cube | iris.cube.CubeList: + """Calculate the saturation equivalent potenital temperature.""" + theta_s = iris.cube.CubeList([]) + for T, P in zip( + iter_maybe(temperature), + iter_maybe(pressure), + strict=True, + ): + theta = potential_temperature(T, P) + ws = saturation_mixing_ratio(T, P) + second_term_power = LV * ws / (CPD * T) + second_term = np.exp(second_term_power.core_data()) + TH_S = theta * second_term + TH_S.rename("saturation_equivalent_potential_temperature") + TH_S.units = "K" + theta_s.append(TH_S) + if len(theta_s) == 1: + return theta_s[0] + else: + return theta_s From 1c6fef4894f80ba173f63a2aedde08d7ec0c25eb Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 5 Jan 2026 13:54:59 +0000 Subject: [PATCH 23/91] Adds to init file and updates name for atmospheric constants --- src/CSET/operators/__init__.py | 6 ++++++ .../{atmospheric_constants.py => _atmospheric_constants.py} | 0 src/CSET/operators/humidity.py | 2 +- src/CSET/operators/pressure.py | 2 +- src/CSET/operators/temperature.py | 2 +- 5 files changed, 9 insertions(+), 3 deletions(-) rename src/CSET/operators/{atmospheric_constants.py => _atmospheric_constants.py} (100%) diff --git a/src/CSET/operators/__init__.py b/src/CSET/operators/__init__.py index ef40fdf22..bb9247c75 100644 --- a/src/CSET/operators/__init__.py +++ b/src/CSET/operators/__init__.py @@ -34,12 +34,15 @@ convection, ensembles, filters, + humidity, imageprocessing, mesoscale, misc, plot, + pressure, read, regrid, + temperature, transect, wind, write, @@ -56,13 +59,16 @@ "ensembles", "execute_recipe", "filters", + "humidity", "get_operator", "imageprocessing", "mesoscale", "misc", "plot", + "pressure", "read", "regrid", + "temperature", "transect", "wind", "write", diff --git a/src/CSET/operators/atmospheric_constants.py b/src/CSET/operators/_atmospheric_constants.py similarity index 100% rename from src/CSET/operators/atmospheric_constants.py rename to src/CSET/operators/_atmospheric_constants.py diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index d15ed319e..25b4782c7 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -17,7 +17,7 @@ import iris.cube from CSET._common import iter_maybe -from CSET.operators.constants import EPSILON +from CSET.operators._atmospheric_constants import EPSILON from CSET.operators.misc import convert_units from CSET.operators.pressure import vapour_pressure diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index 5f5f92c4a..9210fc3b4 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -18,7 +18,7 @@ import numpy as np from CSET._common import iter_maybe -from CSET.operators.atmospheric_constants import E0, KAPPA, P0 +from CSET.operators._atmospheric_constants import E0, KAPPA, P0 from CSET.operators.misc import convert_units diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 087ff24ca..5b66f1b78 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -18,7 +18,7 @@ import numpy as np from CSET._common import iter_maybe -from CSET.operators.constants import CPD, EPSILON, LV, RV, T0 +from CSET.operators._atmospheric_constants import CPD, EPSILON, LV, RV, T0 from CSET.operators.humidity import mixing_ratio_from_RH, saturation_mixing_ratio from CSET.operators.misc import convert_units from CSET.operators.pressure import exner_pressure, vapour_pressure From 2e0bd497e5276d05f6cf54f2f906ea27af311d1f Mon Sep 17 00:00:00 2001 From: daflack Date: Tue, 6 Jan 2026 09:51:11 +0000 Subject: [PATCH 24/91] Update argument names in specific humidity to mixing ratio and vice versa --- src/CSET/operators/humidity.py | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index 25b4782c7..cf0ba2ae0 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -23,13 +23,13 @@ def specific_humidity_to_mixing_ratio( - cubes: iris.cube.Cube | iris.cube.CubeList, + specific_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: """Convert specific humidity to mixing ratio.""" w = iris.cube.CubeList([]) - for cube in iter_maybe(cubes): - mr = cube.copy() - mr = cube / (1 - cube) + for q in iter_maybe(specific_humidity): + mr = q.copy() + mr = q / (1 - q) mr.rename("mixing_ratio") w.append(mr) if len(w) == 1: @@ -39,13 +39,13 @@ def specific_humidity_to_mixing_ratio( def mixing_ratio_to_specific_humidity( - cubes: iris.cube.Cube | iris.cube.CubeList, + mixing_ratio: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: """Convert mixing ratio to specific humidity.""" q = iris.cube.CubeList([]) - for cube in iter_maybe(cubes): - sh = cube.copy() - sh = cube / (1 + cube) + for w in iter_maybe(mixing_ratio): + sh = w.copy() + sh = w / (1 + w) sh.rename("specific_humidity") q.append(sh) if len(q) == 1: From a1a9fd0e78e9a487f6c66de8750cc83f526d62d3 Mon Sep 17 00:00:00 2001 From: daflack Date: Tue, 6 Jan 2026 11:47:55 +0000 Subject: [PATCH 25/91] Uses convert_units for pressure consistency --- src/CSET/operators/humidity.py | 6 ++++-- src/CSET/operators/pressure.py | 1 + 2 files changed, 5 insertions(+), 2 deletions(-) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index cf0ba2ae0..072f2d170 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -61,7 +61,8 @@ def saturation_mixing_ratio( """Calculate saturation mixing ratio.""" w = iris.cube.CubeList([]) for T, P in zip(iter_maybe(temperature), iter_maybe(pressure), strict=True): - mr = (EPSILON * vapour_pressure(T)) / ((P / 100.0) - vapour_pressure(T)) + P = convert_units(P, "hPa") + mr = (EPSILON * vapour_pressure(T)) / (P - vapour_pressure(T)) mr.units = "kg/kg" mr.rename("mixing_ratio") w.append(mr) @@ -78,7 +79,8 @@ def saturation_specific_humidity( """Calculate saturation specific humidity.""" q = iris.cube.CubeList([]) for T, P in zip(iter_maybe(temperature), iter_maybe(pressure), strict=True): - sh = (EPSILON * vapour_pressure(T)) / (P / 100.0) + P = convert_units(P, "hPa") + sh = (EPSILON * vapour_pressure(T)) / P sh.units = "kg/kg" sh.rename("specific_humidity") q.append(sh) diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index 9210fc3b4..fa394b321 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -65,6 +65,7 @@ def exner_pressure( pi = iris.cube.CubeList([]) for P in iter_maybe(pressure): PI = P.copy() + P = convert_units(P, "hPa") PI.data[:] = (P.core_data() / P0) ** KAPPA PI.rename("exner_pressure") PI.units = "1" From be0eb2d4bcbb4fa81f5106c36b457f7ffe87020f Mon Sep 17 00:00:00 2001 From: daflack Date: Wed, 7 Jan 2026 10:12:13 +0000 Subject: [PATCH 26/91] Adds names to humidity cubes where missing --- src/CSET/operators/humidity.py | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index 072f2d170..b68e2c3a2 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -64,7 +64,7 @@ def saturation_mixing_ratio( P = convert_units(P, "hPa") mr = (EPSILON * vapour_pressure(T)) / (P - vapour_pressure(T)) mr.units = "kg/kg" - mr.rename("mixing_ratio") + mr.rename("saturation_mixing_ratio") w.append(mr) if len(w) == 1: return w[0] @@ -82,7 +82,7 @@ def saturation_specific_humidity( P = convert_units(P, "hPa") sh = (EPSILON * vapour_pressure(T)) / P sh.units = "kg/kg" - sh.rename("specific_humidity") + sh.rename("saturation_specific_humidity") q.append(sh) if len(q) == 1: return q[0] @@ -105,6 +105,8 @@ def mixing_ratio_from_RH( ): RH = convert_units(RH, "1") mr = saturation_mixing_ratio(T, P) * RH + mr.rename("mixing_ratio") + mr.units = "kg/kg" w.append(mr) if len(w) == 1: return w[0] @@ -127,6 +129,8 @@ def specific_humidity_from_RH( ): RH = convert_units(RH, "1") sh = saturation_specific_humidity(T, P) * RH + sh.rename("specific_humidity") + sh.units = "kg/kg" q.append(sh) if len(q) == 1: return q[0] From cdecfc584015c99d4b7da25f26fa88168245a8c7 Mon Sep 17 00:00:00 2001 From: daflack Date: Wed, 7 Jan 2026 10:57:21 +0000 Subject: [PATCH 27/91] Correct units --- src/CSET/operators/pressure.py | 2 ++ src/CSET/operators/temperature.py | 2 +- 2 files changed, 3 insertions(+), 1 deletion(-) diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index fa394b321..38d4d5ee9 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -51,6 +51,8 @@ def vapour_pressure_from_RH( ): RH = convert_units(RH, "1") vp = vapour_pressure(T) * RH + vp.units = "hPa" + vp.rename("vapour_pressure") v_pressure.append(vp) if len(v_pressure) == 1: return v_pressure[0] diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 5b66f1b78..a4eaa555f 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -74,7 +74,7 @@ def wet_bulb_temperature( for T, RH in zip( iter_maybe(temperature), iter_maybe(relative_humidity), strict=True ): - RH = convert_units(RH, "1") + RH = convert_units(RH, "%") T = convert_units(T, "Celsius") wetT = ( T * np.arctan(0.151977 * (RH.core_data() + 8.313659) ** 0.5) From 127fecbcec5f0493aa6e30a730edf215e6f2c67b Mon Sep 17 00:00:00 2001 From: daflack Date: Wed, 7 Jan 2026 11:34:02 +0000 Subject: [PATCH 28/91] Adds tests data and fixtures --- tests/conftest.py | 38 ++++++++++++++++++ tests/test_data/pressure/air_temperature.nc | Bin 0 -> 15681 bytes tests/test_data/pressure/pressure.nc | Bin 0 -> 15516 bytes tests/test_data/pressure/relative_humidity.nc | Bin 0 -> 15759 bytes 4 files changed, 38 insertions(+) create mode 100644 tests/test_data/pressure/air_temperature.nc create mode 100644 tests/test_data/pressure/pressure.nc create mode 100644 tests/test_data/pressure/relative_humidity.nc diff --git a/tests/conftest.py b/tests/conftest.py index 7bc1c40bb..b372970fa 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -323,3 +323,41 @@ def orography_4D_cube_read_only(): def orography_4D_cube(orography_4D_cube_read_only): """Get 4D orography cube to run tests on. It is safe to modify.""" return orography_4D_cube_read_only.copy() + + +@pytest.fixture() +def temperature_for_conversions_cube_read_only(): + """Get temperature cube for conversions to run tests on. It is NOT safe to modify.""" + return read.read_cube("tests/test_data/pressure/air_temperature.nc") + + +@pytest.fixture() +def temperature_for_conversions_cube(temperature_for_conversions_cube_read_only): + """Get temperature cube for conversions to run tests on. It is safe to modify.""" + return temperature_for_conversions_cube_read_only.copy() + + +@pytest.fixture() +def pressure_for_conversions_cube_read_only(): + """Get pressure cube for conversions to run tests on. It is NOT safe to modify.""" + return read.read_cube("tests/test_data/pressure/pressure.nc") + + +@pytest.fixture() +def pressure_for_conversions_cube(pressure_for_conversions_cube_read_only): + """Get pressure cube for conversions to run tests on. It is safe to modify.""" + return pressure_for_conversions_cube_read_only.copy() + + +@pytest.fixture() +def relative_humidity_for_conversions_cube_read_only(): + """Get relative humidity cube for conversions to run tests on. It is NOT safe to modify.""" + return read.read_cube("tests/test_data/pressure/relative_humidity.nc") + + +@pytest.fixture() +def relative_humidity_for_conversions_cube( + relative_humidity_for_conversions_cube_read_only, +): + """Get relative humidity cube for conversions to run tests on. It is safe to modify.""" + return relative_humidity_for_conversions_cube_read_only.copy() diff --git a/tests/test_data/pressure/air_temperature.nc b/tests/test_data/pressure/air_temperature.nc new file mode 100644 index 0000000000000000000000000000000000000000..c49d20d9d5387ed0cb904ab6e6c1aa59af2ee587 GIT binary patch literal 15681 zcmeHOZ)_aJ6@PbYU(U`taYC9vh?>PoX$Z;X9J@g_HTZnKi>>%i@HsS4l+AIwXRqqp zJ@(d?KhOeA)RqqjwP_zRf=k@YNZ(lGT2#~x$3c7@-=8$Pn0Up>te;H7u&;FtG{QfQzI?+j*s`P@@~ zXX-a~ghQLDDjW`lgj}eW);Ch0l-qe{E{dRpu<=(zm#ffB^6a8VYQOF0eijZ zO)A1KQCyUC^d zn~4U(Y_@U(fk&gONum;xbBpuMT17OLBG1vbDM{6KkP=}ulRx>mu|3nq%TvNl@a-PjgUH?%K@_Rn!tai)Z% z*~pZ!;E$T7#uaGdqRp)9jq7gKuqI2j$m)J?aDPG_AtO@YCN7Kh)^(V4*YKr9Tk z$OU*b(XD1^Y!n2Fb$^~*xz5xk_$`wOuAdcLp(|^njB{;*5mn{WKON;B!%UZU5;EQ1{=w8xrhj-y&mFu zkXD8{=9*?nZ3S4yddv7qUe-3Ut}(kF_(77!T^fk}fSQ;4vs^S+<{j6}WlNrpJGA8( zo7lb#@aMIOlR~qdhwUqLYZx^|QbY%S(lDn@G+}M8#E6r-<=#DV8+9qpQO0}cAt|3hTVrnvOM-x!61?| zEB)8ORe(T%K!8AiK!8AiK!8AiK!8AiK!8AiK!8Ai!2b^dI6Bl z6<7&zarmn_GIDAjYfvbLbR z_9E8%>G+P*iS$U{UM4Bb8W0T9{VW|{(y1XDkPVf^Nl#6=ph^by<0@i3ccN{WNit=+)cdZT+id2S>}rp7#M-uv z8~LwwwZ>Z8W4j`^ARBt3WVc{$Z{hVxNilM$$L%z9IElk%tP!$;^_Fs+1pk<@>cnmzGQIyZzcno?IO(!w`v5LuG zazF-C2w&Gzgif`9jQ#(+6Edm8>8GbtFQ@x@d&)6;toY@Fn;R57r7t-khVUVhDdnk= z@&R{#1c!umPUz$Uu^}Xg1n*Js1kd9}3OWibjq4?6(#`P*F*>0WS0Ub7inm!i+jc3q zZ}9MuS*{kr)r)=xE$=!7ZP4_z;e6i8nc6srgmD} zVsQ!yV=XY-5J78<>x<*n$P>~`o{-tNAm$07_$BgmNCybT8&P{iG9^Tu23EviuDm9m z8p4DW9AwsoTt&PSIUw5bOS0yK&S>Rp*UX!)Y3EGcLm!PcZdzmGo<@1RE1>HK0lw%% zs3P1`1l`RKldNAs2B&P5aXSFVcFSO|w1`KNb#kc>R}G?4H9~2DDNkN0GQM$~xg{Ra zyM(2rY??gE3$~K31hydS}*``;l-L_VY>+ zW6XW15r%T8PKF(?E02P|g&i?|^I!JA`e0=g+t3`VQEw$h6+9L2WhpyhQtTgFRmrlb z6S8~Bqu)JqwIe32_ue^rNb6l>%<#Cg;F!q`tDxc|#iFA^ps~?qiziKoimn%(lI2-b ze0w7<`c4*E6d~s3&~xWc{*iE%dRL!y+(vuwd`I*UPpBiDMsIm!-?akuK<%Jj?L@#p z!bKsM`!{V=9Ell$M7uPYx5 zjA12Gg+^kZO*;0JX?y6?=Z=190?G9D*4Q1VPQRd2gSE_ppTBtg+0l_yivCMOzSodx z3)i9Q??&-2Fp@zUflTE~1z8`Y0RjO60RjO60RjO60RjO6*Ea&w2mUtw>Z;?b9@^CX d?r)YIyYu2<{qm15e)!?{wXquk63;gu{5K9!(9Hk< literal 0 HcmV?d00001 diff --git a/tests/test_data/pressure/pressure.nc b/tests/test_data/pressure/pressure.nc new file mode 100644 index 0000000000000000000000000000000000000000..9d3f1864c5265df0016842e827e69868e13bd365 GIT binary patch literal 15516 zcmeHNeQaA-6+iZKJjZ=Wvr<-C3iU#vP1j~lmUUW0(r#&O3J}2NNyfO<{j1M4Dbz^x(T?OH_6I@J9;+W68mB zzj@li!Uh*U4V<-X9}Lj9IO zzi%sg{eGWM_tq<$!dQWEja2>o?Q5>WZ@(W4^ATE6%35dxq#AxB*^P>v8iK+0=FOro z!t}vS=p%YTp%-qJzm4@Qhpvol5sX9;+7O}*A3^6lcPwv}!y5ZGc~c?mwUhZaUjNF& z&1!3CyP$jnn};jI{8T}pPsSlP23OT zf+F?Li+Ke8(1={3 zqy78L6!ykopv3Z*$<=P9F~N(r5ZqK1T%)U)qk?l|f*DogEscCW-Z8lc%RE+H z_fczsMCWgpj`DeB{>LGZNjs>}MeT1rnJ^q{)GlYTxt>^|P|P}c_Rb}j;ny=1^BQS? zL21UZ;r>)gm@fP`Om}c-BtDuN8XL7zvHsz>>!5!7WBWIb#ggxpYnWFqlr36PSXX3u z^KJimPnCOdqip@8VluZ&Yi*K>jP{Sjg&wWk4?X0zX*+}3u*X(jJ)b_k`ax&L$-39A z#lKY=ni8Kr;*^AV92$s!KT-Vj*^`M>a**hY>L}?GrAsA?XFNbhh!s}e`^)O0!E*@G zU#{N|cI}iNP;yM!O=dgLH~1r5&RpIp%SBwl z3|K_go-EpS$?PJ!LVzRqO67Dun=WRnZpn?QK4aVex%eS96P^9zBcHDB)nC;`@&VNy zjG*ERU`|krYjaL$r@sb$e0QK5V%xqqF)+aZhoJ41+?m}z{|+z_h51j?-<~-AqTp%< zR|_{X*E?ljmFvu8_Sdyo&A(Qx_=ajlte^na>0071p6=^Nq*8IZ6VNmH*kle0IDV4< z{^z;j&=dDOx<58HEKEn-3u$8^UcMxaix|gEx?S6YVxz~`>Y7)22zUs12zUs12zUs1 z2zUs12zUs12zUs12zUtmzadZ!phQ)G=nkc#IgQY`L2D_B;N5sU5g(6_#^ToSP-=p9 zb3gzggKWB-F4<)wyvjguAj|RKe90-~DTIC*PKsPG2&`Y-Fq@)8!n!jto;)BV`NagP zM~Aj`*h>dB)C77WEvOQysWjB7uzot?r2=K!xxA}GvyXr2nPe5CkE@qN9BGy$LZVVP zljDvE__ea~EYcc{9VixN?PA&CFed`-?fgL)<*;1@K6|=%=0!YoPfzEKo!9N0%I2=$ z)79D4)45y5Du1lvHvDcxl88Q4x&HCz-dDKe@z2fv`?2J+FF!{6PQ~U2PQ|vwPsQYb z^LtU6{hRwywy6?zk30lC1Uv*h1Uv*h1TJ3$1}26Re7|uER-h6!zH{W-ae5vNDnH`q zzdLn|&f+BxzmUW&nr-&oOP>LU&dhKgI!UL|pz<)@7yjwP6vq&i^>F67g|AZxcP}GE zxKTV~efUYb8XQ|?2F1`znKN_@4IjYUDx%E2pEyTXafQT<5jQ-wc<$EUe7A%qk1h=Y zV*i_6FuRPwz0I{->@pG8PP#VUwWTgAa+yQxJ+I7PniY9c8bLs2d7Qh1SsnrPhF7GK zu_|v!mIvHM!&x3CKcpddXlXXaS{4RRiXXAr#`*^~n3FmyU7M8~20^btEnH zjoqI}4h~e30!8z`cK_g*!B6R{PO||%AYG!4Nh*)8%M&t?3j_ALnV1(y`T##Z2hs*V z%*%|B8H;6?rmRw7u9%T2Xn3#%+f1~p)V15`?(R1D5+V~rs$A`gYmh?@kzBDb6B)6~ zk+EFP$=H#Dc_-)C*~myCYfl@rrZSm&+7?gXc+$r_UDtI+k#il6XFmQ zX+uz?ns_9FB7g`%4;>3ZLKEf#=@KX}s!$8+wZFmBK)#T{`z@A`sY$3pqyfpZkd-SG z?M%8zkWhr1PEQ=pf_F(G z4_j|psK(}A;(J%n&9YEQfHF20Y)7Ln05w5tm&|;ZNFFEpeeciKy@GcpZGY!k4 zav3%1W7`YX3t+p2n=9BxN)8{}NY{-wM0Z8I_&*Yj_KJU#$FGPwyb|P8*6#KKOGBn_ zANqTsitcgF*8_2+S5i&nGbr}W7f|`3w4u%n&#CC?GYh}%g=B4qUs+;@=DVP%!Ot1e zyQQaXBnJxT3|qw_^!3RZx~=2&?wHo&ktKTQ^YB01vuzjIGkwK^hn9x5o(~tOM;5i| z>${%z4**xKcl~o5ge~xSP!VPdK>nwc$(l+jb#aLm@`vqwnIqORjE!B!4#X0>x;k%o z^u#y(n6r^R`|8(@pPU$v^Yf)|*=g7e#a5uB{6JbuOOu*~ORQz^>K+0f0v-Y$0v-Y$ o0v-Y$0+%%cM-FAC$QS%*u(`Q;sAH(ZIq|LNuUZXy)@kYcFWO4AbpQYW literal 0 HcmV?d00001 diff --git a/tests/test_data/pressure/relative_humidity.nc b/tests/test_data/pressure/relative_humidity.nc new file mode 100644 index 0000000000000000000000000000000000000000..908b6e32c96aec869838440042045006d3ea5d8a GIT binary patch literal 15759 zcmeHOeQX@X6@PnszFeI3M<697iOAN0CIoUm$0jsszRu^n*uo!)&!7<6ET`K&dnMoP zvA2#~HEojy1^P!7rD^jIDvAVcRTb2x0i>!QL~8yhjSwJ#1hlCHQVBtYrYe6_NHl$K zKJL6{pW`EG)kr(`-R#c1dGqEsZ{F^_dH0?EzI03Tisn!>3Nl?0>DV@25*HOOwJv^d zq`z-8q&*kvU}frJQ6-zy=FK5ch$M15eJf%5rnJz%nML3btu`U;8_+fXD{A zP{OPkUafE)FY8&RLG<2~X&IAFahzRkQhA&9cFJR70Qd7kfF<%mRnv)%LkX-LieE%4NL(lrJP?yms3}2NZs7c`y{6HA+pb{VwsBk zR`A**LgJr_TwR8IV5z>E5d!G2JPOOdMZ`z3VshW_CJZ-g+Ee95cVTXcwF$(h8}FU&wAd;tz6c4`@# z8x?_K!=J@hT9w8aziKhYO;yG-WMzF+aBhq-qN;xKmy=x64k!xF7@vb#hN`aBn5lsP zmk)9~>D)@30KUhwid~9cRsWWOMv7w$TVBr0_aq(1HSHoA=dw@n>nVzM7uUZCSL5Vh zB9mdNOXizYw{Ku5HJlk388$M>#9#`vbzG~Oc!S}@P>RvF3VIH&)CXcudee$h;wq6+ z_SyJ+@yjO$uC3$>A=qc9?J{Fe%1z=t2yo%U7oJFG`uhN$6)G<^U*>`r_P86S0RO5?jlgUGP@+61g-UTLd4CSA*6;=1j#8 zGYwNk+v9^+gakAQ(;9aztE_baEMnj=y0VuonpxMZnh!swaxIhxV~>lTVB@33Pe~N6 zMJL2F>dP z+uS3+WmJoas)bYCQA>61AGrr-Vio(WSmFE5yjYO{G)~_Tr*{`mrZbrobP)8I6jYfE zgX5$3WAqTuy#4R@cO^##nd-1t7O5M%>a-$ruT*7DD#=DEmZ@H}}wmH^)_k@|h zXIodSt0%UJd$c^^x4|w_w_t#i)_UY26FZRFl$&0Y_O!8RondIu!nIvy; zc?acMVpW3Q|BKzXSMGUf@Esr!AP^uBAP^uBAaG+L&^tDm#`lefNa`Mbdh>_r&QW-q z3Pc>kpFTKq7+$AC9z2kwS2Ro8egcjWg$$supF9meqXK-_pyS*<&2K#qTd9eNPB{PW zx!=GQDi{Nxlf?7JH(r9zQ^D;#(qZ7E+yO{NJ|haaG=cz+#<+AH zqcNTmyABrRg`&tPvHxl`2APePBJ7f}vr`|tq|4F*PueTDj$BH}GzIj~5UHXE5JWyJ zP;!o74`M~Mcc{3B=Yb#v16h_QjIuN3=6Hk{?PrOq5br9-yX}oV_o(>f;PD_;st!Tb zhn@yK?>Ys2$nx}&eBRDk`u#;aZ(F85EwGqK5L8*kd+Q$JnZu{ecb3()sAvknB^PXrL2oq9qn%NL? zGvb}d!O)E#;+1z?E0-;M26`WC$J8tLu9df3tC+J454|+{gk_IUcsgBUt03Jd2ymkh z%8X!75pp*>NV2JwR83tk<8}a^+$CeT(jp#7)=5(zsvmf#6Gn3^_$DkS)tkU|Th=y` ze$Qh{q=l|1CRbEgPuv%t_8q=&OVu~`U3_X_N0-f;5|&}i%PZ*mdpIn-rz_W;e~vgy zJshslL+y(ssEYdxZsxL+76k_~D^0pegOJas_RSD1i-~sOU3QVsVLhRg-@C&d>i$xJtR}uVb#G z-^X=I=@FjMKzRDwo;S{2FHrY44(hcogtb?JT)q#`gO4N8M2}8PxP|6P8VZ0PnN#}+ z)Ch}GR>I=-7$W@X*&XK%BDA{D!bn>hCgzWeMWYFc#B)386c1QM56=Mj!GKhOWO_qa z?A~QZw)RtpjjV!S{OZ){vC&kDe)l0yN@RW;7qIU4BjQIjl0Yhf%qzzuk|4+f1OfyC z1OfyC1OfyC1OfzZXauJB=N>+!JQjH@a`b5A=u^8|mDi5{e$fvf26+18C8z%bqG;Lw literal 0 HcmV?d00001 From db884234f17dfdcbe5d3d39b49f5e33252bf639b Mon Sep 17 00:00:00 2001 From: daflack Date: Wed, 7 Jan 2026 11:36:10 +0000 Subject: [PATCH 29/91] Update copyright year --- tests/conftest.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/conftest.py b/tests/conftest.py index b372970fa..dbc08c7c2 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -1,4 +1,4 @@ -# © Crown copyright, Met Office (2022-2024) and CSET contributors. +# © Crown copyright, Met Office (2022-2026) and CSET contributors. # # Licensed under the Apache License, Version 2.0 (the "License"); # you may not use this file except in compliance with the License. From 19d63be6befc127aac4d93686746c734f80c3140 Mon Sep 17 00:00:00 2001 From: daflack Date: Wed, 7 Jan 2026 11:55:40 +0000 Subject: [PATCH 30/91] Add tests for vapour pressure --- tests/operators/test_pressure.py | 74 ++++++++++++++++++++++++++++++++ 1 file changed, 74 insertions(+) create mode 100644 tests/operators/test_pressure.py diff --git a/tests/operators/test_pressure.py b/tests/operators/test_pressure.py new file mode 100644 index 000000000..ec83bb182 --- /dev/null +++ b/tests/operators/test_pressure.py @@ -0,0 +1,74 @@ +# © Crown copyright, Met Office (2022-2026) and CSET contributors. +# +# Licensed under the Apache License, Version 2.0 (the "License"); +# you may not use this file except in compliance with the License. +# You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software +# distributed under the License is distributed on an "AS IS" BASIS, +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +# See the License for the specific language governing permissions and +# limitations under the License. + +"""Test pressure operators.""" + +import cf_units +import iris.cube +import numpy as np + +from CSET.operators import pressure + + +def test_vapour_pressure(temperature_for_conversions_cube): + """Test calculation of vapour pressure for a cube.""" + expected_data = temperature_for_conversions_cube.copy() + exponent = ( + 17.27 + * (temperature_for_conversions_cube - 273.16) + / (temperature_for_conversions_cube - 35.86) + ) + expected_data.data = 6.1078 * np.exp(exponent.core_data()) + assert np.allclose( + expected_data.data, + pressure.vapour_pressure(temperature_for_conversions_cube).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_vapour_pressure_name(temperature_for_conversions_cube): + """Test renaming of vapour pressure cube.""" + expected_name = "vapour_pressure" + assert ( + expected_name + == pressure.vapour_pressure(temperature_for_conversions_cube).name() + ) + + +def test_vapour_pressure_units(temperature_for_conversions_cube): + """Test units for vapour pressure cube.""" + expected_units = cf_units.Unit("hPa") + assert ( + expected_units + == pressure.vapour_pressure(temperature_for_conversions_cube).units + ) + + +def test_vapour_pressure_cube_list(temperature_for_conversions_cube): + """Test calculation of vapour pressure for a CubeList.""" + expected_data = temperature_for_conversions_cube.copy() + exponent = ( + 17.27 + * (temperature_for_conversions_cube - 273.16) + / (temperature_for_conversions_cube - 35.86) + ) + expected_data.data = 6.1078 * np.exp(exponent.core_data()) + expected_list = iris.cube.CubeList([expected_data, expected_data]) + input_cubelist = iris.cube.CubeList( + [temperature_for_conversions_cube, temperature_for_conversions_cube] + ) + actual_cubelist = pressure.vapour_pressure(input_cubelist) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From 51da41f9253eb2f4ecaffcb37f0b05b4dad2e01d Mon Sep 17 00:00:00 2001 From: daflack Date: Wed, 7 Jan 2026 12:02:11 +0000 Subject: [PATCH 31/91] Update to assignment in pressure convertors --- src/CSET/operators/pressure.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index 38d4d5ee9..b381ae1c9 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -30,7 +30,7 @@ def vapour_pressure( for T in iter_maybe(temperature): es = T.copy() exponent = 17.27 * (T - 273.16) / (T - 35.86) - es.data[:] = E0 * np.exp(exponent.core_data()) + es.data = E0 * np.exp(exponent.core_data()) es.units = "hPa" es.rename("vapour_pressure") v_pressure.append(es) @@ -68,7 +68,7 @@ def exner_pressure( for P in iter_maybe(pressure): PI = P.copy() P = convert_units(P, "hPa") - PI.data[:] = (P.core_data() / P0) ** KAPPA + PI.data = (P.core_data() / P0) ** KAPPA PI.rename("exner_pressure") PI.units = "1" pi.append(PI) From 25ebb9b3cc31ffda75c8db393deb8d587a1e6dc0 Mon Sep 17 00:00:00 2001 From: daflack Date: Wed, 7 Jan 2026 12:12:57 +0000 Subject: [PATCH 32/91] Rename vapour_pressure_from_RH to relative_humidity_to_vapour_pressure --- src/CSET/operators/pressure.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index b381ae1c9..14fa734c6 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -40,7 +40,7 @@ def vapour_pressure( return v_pressure -def vapour_pressure_from_RH( +def relative_humidity_to_vapour_pressure( temperature: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: From 5251814e792d0de7293da971c204e54ade0bb6aa Mon Sep 17 00:00:00 2001 From: daflack Date: Wed, 7 Jan 2026 12:25:59 +0000 Subject: [PATCH 33/91] Adds tests for relative_humidity_to_vapour_pressure --- tests/operators/test_pressure.py | 64 ++++++++++++++++++++++++++++++++ 1 file changed, 64 insertions(+) diff --git a/tests/operators/test_pressure.py b/tests/operators/test_pressure.py index ec83bb182..cd838fc32 100644 --- a/tests/operators/test_pressure.py +++ b/tests/operators/test_pressure.py @@ -72,3 +72,67 @@ def test_vapour_pressure_cube_list(temperature_for_conversions_cube): actual_cubelist = pressure.vapour_pressure(input_cubelist) for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) + + +def test_relative_humidity_to_vapour_pressure( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube +): + """Test calculation of vapour pressure from relative humidity.""" + expected_data = pressure.vapour_pressure(temperature_for_conversions_cube) * ( + relative_humidity_for_conversions_cube / 100.0 + ) + assert np.allclose( + expected_data.data, + pressure.relative_humidity_to_vapour_pressure( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_relative_humidity_to_vapour_pressure_name( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube +): + """Test naming of vapour pressure cube.""" + expected_name = "vapour_pressure" + assert ( + expected_name + == pressure.relative_humidity_to_vapour_pressure( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube + ).name() + ) + + +def test_relative_humidity_to_vapour_pressure_units( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube +): + """Test units of vapour pressure cube.""" + expected_units = cf_units.Unit("hPa") + assert ( + expected_units + == pressure.relative_humidity_to_vapour_pressure( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube + ).units + ) + + +def test_relative_humidity_to_vapour_pressure_cubelist( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube +): + """Test calculation of vapour pressure from relative humidity as a CubeList.""" + expected_data = pressure.vapour_pressure(temperature_for_conversions_cube) * ( + relative_humidity_for_conversions_cube / 100.0 + ) + expected_list = iris.cube.CubeList([expected_data, expected_data]) + temperature_input_list = iris.cube.CubeList( + [temperature_for_conversions_cube, temperature_for_conversions_cube] + ) + relative_humidity_input_list = iris.cube.CubeList( + [relative_humidity_for_conversions_cube, relative_humidity_for_conversions_cube] + ) + actual_cubelist = pressure.relative_humidity_to_vapour_pressure( + temperature_input_list, relative_humidity_input_list + ) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From dea81b03b27fb29b064b73f44270b8228e191f0c Mon Sep 17 00:00:00 2001 From: daflack Date: Wed, 7 Jan 2026 12:57:20 +0000 Subject: [PATCH 34/91] Adds tests for exner pressure --- tests/operators/test_pressure.py | 45 ++++++++++++++++++++++++++++++++ 1 file changed, 45 insertions(+) diff --git a/tests/operators/test_pressure.py b/tests/operators/test_pressure.py index cd838fc32..4ac249b83 100644 --- a/tests/operators/test_pressure.py +++ b/tests/operators/test_pressure.py @@ -136,3 +136,48 @@ def test_relative_humidity_to_vapour_pressure_cubelist( ) for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) + + +def test_exner_pressure(pressure_for_conversions_cube): + """Test calculation of exner pressure.""" + expected_data = pressure_for_conversions_cube.copy() + expected_data.data = ( + (pressure_for_conversions_cube / 100).core_data() / 1000.0 + ) ** (1005.7 / 287.0) + assert np.allclose( + expected_data.data, + pressure.exner_pressure(pressure_for_conversions_cube).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_exner_pressure_name(pressure_for_conversions_cube): + """Test naming of exner pressure.""" + expected_name = "exner_pressure" + assert ( + expected_name == pressure.exner_pressure(pressure_for_conversions_cube).name() + ) + + +def test_exner_pressure_units(pressure_for_conversions_cube): + """Test units of exner pressure.""" + expected_units = cf_units.Unit("1") + assert ( + expected_units == pressure.exner_pressure(pressure_for_conversions_cube).units + ) + + +def test_exner_pressure_cubelist(pressure_for_conversions_cube): + """Test calculation of exner pressure for a CubeList.""" + expected_data = pressure_for_conversions_cube.copy() + expected_data.data = ( + (pressure_for_conversions_cube / 100).core_data() / 1000.0 + ) ** (1005.7 / 287.0) + expected_list = iris.cube.CubeList([expected_data, expected_data]) + input_list = iris.cube.CubeList( + [pressure_for_conversions_cube, pressure_for_conversions_cube] + ) + actual_cubelist = pressure.exner_pressure(input_list) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From 9041869d82f286e0fa0e154860387899e3abc990 Mon Sep 17 00:00:00 2001 From: daflack Date: Wed, 7 Jan 2026 13:11:33 +0000 Subject: [PATCH 35/91] Uses atmospheric constants in tests --- tests/operators/test_pressure.py | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/tests/operators/test_pressure.py b/tests/operators/test_pressure.py index 4ac249b83..ee194b161 100644 --- a/tests/operators/test_pressure.py +++ b/tests/operators/test_pressure.py @@ -18,7 +18,7 @@ import iris.cube import numpy as np -from CSET.operators import pressure +from CSET.operators import _atmospheric_constants, pressure def test_vapour_pressure(temperature_for_conversions_cube): @@ -29,7 +29,7 @@ def test_vapour_pressure(temperature_for_conversions_cube): * (temperature_for_conversions_cube - 273.16) / (temperature_for_conversions_cube - 35.86) ) - expected_data.data = 6.1078 * np.exp(exponent.core_data()) + expected_data.data = _atmospheric_constants.E0 * np.exp(exponent.core_data()) assert np.allclose( expected_data.data, pressure.vapour_pressure(temperature_for_conversions_cube).data, @@ -64,7 +64,7 @@ def test_vapour_pressure_cube_list(temperature_for_conversions_cube): * (temperature_for_conversions_cube - 273.16) / (temperature_for_conversions_cube - 35.86) ) - expected_data.data = 6.1078 * np.exp(exponent.core_data()) + expected_data.data = _atmospheric_constants.E0 * np.exp(exponent.core_data()) expected_list = iris.cube.CubeList([expected_data, expected_data]) input_cubelist = iris.cube.CubeList( [temperature_for_conversions_cube, temperature_for_conversions_cube] @@ -142,8 +142,8 @@ def test_exner_pressure(pressure_for_conversions_cube): """Test calculation of exner pressure.""" expected_data = pressure_for_conversions_cube.copy() expected_data.data = ( - (pressure_for_conversions_cube / 100).core_data() / 1000.0 - ) ** (1005.7 / 287.0) + (pressure_for_conversions_cube / 100).core_data() / _atmospheric_constants.P0 + ) ** _atmospheric_constants.KAPPA assert np.allclose( expected_data.data, pressure.exner_pressure(pressure_for_conversions_cube).data, @@ -172,8 +172,8 @@ def test_exner_pressure_cubelist(pressure_for_conversions_cube): """Test calculation of exner pressure for a CubeList.""" expected_data = pressure_for_conversions_cube.copy() expected_data.data = ( - (pressure_for_conversions_cube / 100).core_data() / 1000.0 - ) ** (1005.7 / 287.0) + (pressure_for_conversions_cube / 100).core_data() / _atmospheric_constants.P0 + ) ** _atmospheric_constants.KAPPA expected_list = iris.cube.CubeList([expected_data, expected_data]) input_list = iris.cube.CubeList( [pressure_for_conversions_cube, pressure_for_conversions_cube] From 26dccaf3c5a9e2259b8dd2e115a28b4270f667f6 Mon Sep 17 00:00:00 2001 From: daflack Date: Wed, 7 Jan 2026 14:13:25 +0000 Subject: [PATCH 36/91] Change to vapour pressure from relative humidity as more intuitive --- src/CSET/operators/pressure.py | 2 +- tests/operators/test_pressure.py | 16 ++++++++-------- 2 files changed, 9 insertions(+), 9 deletions(-) diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index 14fa734c6..fe03c3b8c 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -40,7 +40,7 @@ def vapour_pressure( return v_pressure -def relative_humidity_to_vapour_pressure( +def vapour_pressure_from_relative_humidity( temperature: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: diff --git a/tests/operators/test_pressure.py b/tests/operators/test_pressure.py index ee194b161..42e7d78ff 100644 --- a/tests/operators/test_pressure.py +++ b/tests/operators/test_pressure.py @@ -74,7 +74,7 @@ def test_vapour_pressure_cube_list(temperature_for_conversions_cube): assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) -def test_relative_humidity_to_vapour_pressure( +def test_vapour_pressure_from_relative_humidity( temperature_for_conversions_cube, relative_humidity_for_conversions_cube ): """Test calculation of vapour pressure from relative humidity.""" @@ -83,7 +83,7 @@ def test_relative_humidity_to_vapour_pressure( ) assert np.allclose( expected_data.data, - pressure.relative_humidity_to_vapour_pressure( + pressure.vapour_pressure_from_relative_humidity( temperature_for_conversions_cube, relative_humidity_for_conversions_cube ).data, rtol=1e-6, @@ -91,33 +91,33 @@ def test_relative_humidity_to_vapour_pressure( ) -def test_relative_humidity_to_vapour_pressure_name( +def test_vapour_pressure_from_relative_humidity_name( temperature_for_conversions_cube, relative_humidity_for_conversions_cube ): """Test naming of vapour pressure cube.""" expected_name = "vapour_pressure" assert ( expected_name - == pressure.relative_humidity_to_vapour_pressure( + == pressure.vapour_pressure_from_relative_humidity( temperature_for_conversions_cube, relative_humidity_for_conversions_cube ).name() ) -def test_relative_humidity_to_vapour_pressure_units( +def test_vapour_pressure_from_relative_humidity_units( temperature_for_conversions_cube, relative_humidity_for_conversions_cube ): """Test units of vapour pressure cube.""" expected_units = cf_units.Unit("hPa") assert ( expected_units - == pressure.relative_humidity_to_vapour_pressure( + == pressure.vapour_pressure_from_relative_humidity( temperature_for_conversions_cube, relative_humidity_for_conversions_cube ).units ) -def test_relative_humidity_to_vapour_pressure_cubelist( +def test_vapour_pressure_from_relative_humidity_cubelist( temperature_for_conversions_cube, relative_humidity_for_conversions_cube ): """Test calculation of vapour pressure from relative humidity as a CubeList.""" @@ -131,7 +131,7 @@ def test_relative_humidity_to_vapour_pressure_cubelist( relative_humidity_input_list = iris.cube.CubeList( [relative_humidity_for_conversions_cube, relative_humidity_for_conversions_cube] ) - actual_cubelist = pressure.relative_humidity_to_vapour_pressure( + actual_cubelist = pressure.vapour_pressure_from_relative_humidity( temperature_input_list, relative_humidity_input_list ) for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): From 3a7a197037d59aca32304cab7d1418d58760b54d Mon Sep 17 00:00:00 2001 From: daflack Date: Wed, 7 Jan 2026 14:18:47 +0000 Subject: [PATCH 37/91] Update aming convention to be x_from_y rather than y_from_x --- src/CSET/operators/humidity.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index b68e2c3a2..a0305beaf 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -22,7 +22,7 @@ from CSET.operators.pressure import vapour_pressure -def specific_humidity_to_mixing_ratio( +def mixing_ratio_from_specific_humidity( specific_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: """Convert specific humidity to mixing ratio.""" @@ -38,7 +38,7 @@ def specific_humidity_to_mixing_ratio( return w -def mixing_ratio_to_specific_humidity( +def specific_humidity_from_mixing_ratio( mixing_ratio: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: """Convert mixing ratio to specific humidity.""" @@ -90,7 +90,7 @@ def saturation_specific_humidity( return q -def mixing_ratio_from_RH( +def mixing_ratio_from_relative_humidity( temperature: iris.cube.Cube | iris.cube.CubeList, pressure: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, @@ -114,7 +114,7 @@ def mixing_ratio_from_RH( return w -def specific_humidity_from_RH( +def specific_humidity_from_relative_humidity( temperature: iris.cube.Cube | iris.cube.CubeList, pressure: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, From 3a46b54808b7550fe8ff7aa12dc5bbe833538be7 Mon Sep 17 00:00:00 2001 From: daflack Date: Wed, 7 Jan 2026 15:32:09 +0000 Subject: [PATCH 38/91] Update operators in temperature for correct conversion name --- src/CSET/operators/temperature.py | 7 +++++-- 1 file changed, 5 insertions(+), 2 deletions(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index a4eaa555f..c1221bc5b 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -19,7 +19,10 @@ from CSET._common import iter_maybe from CSET.operators._atmospheric_constants import CPD, EPSILON, LV, RV, T0 -from CSET.operators.humidity import mixing_ratio_from_RH, saturation_mixing_ratio +from CSET.operators.humidity import ( + mixing_ratio_from_relative_humidity, + saturation_mixing_ratio, +) from CSET.operators.misc import convert_units from CSET.operators.pressure import exner_pressure, vapour_pressure @@ -147,7 +150,7 @@ def equivalent_potential_temperature( ): RH = convert_units(RH, "1") theta = potential_temperature(T, P) - w = mixing_ratio_from_RH(T, P, RH) + w = mixing_ratio_from_relative_humidity(T, P, RH) second_term_power = -(w * RV) / CPD second_term = RH.core_data() ** second_term_power.core_data() third_term_power = LV * w / (CPD * T) From 7d782603f7ffd0594f595fd6858911b15f8d63e0 Mon Sep 17 00:00:00 2001 From: daflack Date: Wed, 7 Jan 2026 15:50:30 +0000 Subject: [PATCH 39/91] Adds test data and fixtures for humidity convertors --- tests/conftest.py | 28 ++++++++++++++++++ tests/test_data/humidity/mixing_ratio.nc | Bin 0 -> 20369 bytes tests/test_data/humidity/specific_humidity.nc | Bin 0 -> 20369 bytes 3 files changed, 28 insertions(+) create mode 100644 tests/test_data/humidity/mixing_ratio.nc create mode 100644 tests/test_data/humidity/specific_humidity.nc diff --git a/tests/conftest.py b/tests/conftest.py index dbc08c7c2..a04207353 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -361,3 +361,31 @@ def relative_humidity_for_conversions_cube( ): """Get relative humidity cube for conversions to run tests on. It is safe to modify.""" return relative_humidity_for_conversions_cube_read_only.copy() + + +@pytest.fixture() +def specific_humidity_for_conversions_cube_read_only(): + """Get specific humidity cube for conversions to run tests on. It is NOT safe to modify.""" + return read.read_cube("tests/test_data/humidity/specific_humidity.nc") + + +@pytest.fixture() +def specific_humidity_for_conversions_cube( + specific_humidity_for_conversions_cube_read_only, +): + """Get specific humidity cube for conversions to run tests on. It is safe to modify.""" + return specific_humidity_for_conversions_cube_read_only.copy() + + +@pytest.fixture() +def mixing_ratio_for_conversions_cube_read_only(): + """Get mixing ratio cube for conversions to run tests on. It is NOT safe to modify.""" + return read.read_cube("tests/test_data/humidity/mixing_ratio.nc") + + +@pytest.fixture() +def mixing_ratio_for_conversions_cube( + mixing_ratio_for_conversions_cube_read_only, +): + """Get mixing ratio cube for conversions to run tests on. It is safe to modify.""" + return mixing_ratio_for_conversions_cube_read_only.copy() diff --git a/tests/test_data/humidity/mixing_ratio.nc b/tests/test_data/humidity/mixing_ratio.nc new file mode 100644 index 0000000000000000000000000000000000000000..a30b3751b8a4857c3040c02a6c752fde1e3d7187 GIT binary patch literal 20369 zcmeGkZEPIH_3iE1P3%kJ1b;vXVG{*7oW!y73B(Ehyg0=lscj&Xs>k_m&snhVoV&9L zDB&X@rKD9;f*KW6C=DO2C_;i(1q#H_($bKCP+LGP6(Ud-f=c*M#a9dIdo%O)ZZG!b zF2+%R>?G^ioq6+i=FOWo@6CIAx3skcrj^Vt@l2TloQ|q=YnLpkRTax#xhc`r*3#{n z)~R{si$Z<6s=KaL-u%5<=BMj2FY$mzcNKHMi<#~L#-5$~1O^?Z!xRDyb~bNF#D~m8 z$}*F13bv&;dhi-C6G%d0_)=0l!_dQvsZ=0672c?(h4%n%gd#GK(s7 zEp@5{V7~U@e0Y#4rb?y=hX%}8B$N=WMKIKjWGX~uk=}}hr(P`vQB37Cuqxv3Ym*3l z$X7C@-l%;K@6-U`YQ^(Bf}#})PGyHbKYAP9-t)xJR=JFRr)dkQ&g0QE)fPJ_(=mBD zQD&xbrdUS+22@Rfs_PNBsDx_tYFb|aA8h5NRs!3I(=qK=J(6q#z$5}iBGCx~Tta}n zFzr&RDKz=Mb?i(QiX%7tB_DloVQz?%Q+sgpE$af^Z5xoSh!x7xZ)HLYO z^?!+0)DMO#>won{d40u=74yLyFz~y)Y%5@mN?bgFkNd>L+2O#{yONXDoe6Mgo$rAVPL@6NVf2xDb zj!AG#vx7~X^PR~r4qkE0v4bTX?6CK>ZHFd07Po!j6p(=k-}Y$2sY+{E_mRJ|3>U5R zvl@;71FX7%!w1;EE#cBJKE-wj0Yh-^Yos3PFbib@w(3tFrxXvm872VN`Y#S~ni!yG8be_%VDgYjr^C}9nn!9f#ITcH>ROkAjZjSaB86qX-F35U@Z z8jeFG&ebFMQ07sT$TW%C3?)*Vf{9SX8cr61tKsdAjg7u;8mkVcNnB#Iu{a?diVvWl zTsh1O!Hd*yZ1DBa+~_cS9K5#wy`aIrWD3MZX~Y?Q6@9nqTv(5hwEE|-?~$DLfTrP$ z$vNN<0y{r}5SM^*Sl1HAV7f%*%%k-x#D)FKa5l}j&5xoD2gIm{d@>(yeA%p%vz^Z` znK~zsT}kMfDt2M68y|)o@G=||&Q~M%@-ZvNDZFyWm928t&dZoJMPRW+si7`sZQdAx zNh0=7FDb*=0v*~_e9B(+HV#uNvS)*5g>#|;wJes>qPe}Xr$tU)W{keJjtNvHkyJ-1iVaA1sJFA0;L zYx_CxmB4(yK=HZC?08FkeCi6VN=F&u!^<8_b3C4iSTPLgNjiswGLyTqM|s<2`qFi4bAN9&f=)RHc6jJA@eqGE?2hci|c#WK^GC%qsVPnh9QG8M!?XDc2t z?B|3TH4|nmYz9*(ql`_a)xRlaB+Sr&buCdhW>7g|aMS+}1(Y%UGWAGp$UB~DS_PL< z0qz~dw9mf1ryic9yIN@o|8dKbZuo#8)E2$;pZ8n?+o`}WnI!n#8+XD3RDd{0EfD7R z9fSi}U91u^x93v9R%uLE{`D8#Cu>+P-aIKF6}<;Dlx0@p?L2SAj@9E52&oS-7U1fqbmRh= zv+8$6wu$2>L0~J<5?ICjtA~6tQn?6G-pO~a6jymRApA_0yqX-}#=&n#!=9Z|3TCQ6m*(DQZv3<1-%q83eJ!U8ftL>kB(fhIrHlqWn>VxHx% zZ}@q^$en#cLFXd_0e$vW?P1?pOeL+WaB#&75OsF88 zT%inYv80|&eIt~Y+1hbmd^i@7mGpVkG4BiZVfAFr4u2bkFBgeti=S;0rmY1Md(_d0 z*`F{??49~DkU!E&n;-No+c868%x>;dfGLpuU1D^NtVvt)@pta8mLCD}r#ZaPIcHqi z9|7H`08sH4iSveAJ@+fa;1#OKX&CGyB>1r;ndQPT_>h{YeJyPmsGU|A2EQN(wU-LR z;KzhxeuiOSP(d!kU`yiJyAi}qlSP623?r(?2) z<>Jkg0#d;QEyK&E41r_4Q$^EBh0&{7bHW43b*$kn^F?Lrzi`#_f zW=Y>t(wlPZu=MpE?XdRwlsvBj^|e2}BcF8V4|-EV2VVH?L%!srYzt+BAf?xbm9)MHAlBg>}@BA$pA=2T>POpMyiT?z;a zOjH(-jr$KvB~Y3GPh$r`qu`{X0Ab4p;RGBtKpE$nZ; z&y}O~Y`2sWz`zmpeVyBE|K8e*rS8jbW;ww*LrUG(z(-l*LH(1`Y?;~lySvZz3*D#a zh4X%fjx?;E2mVj~a|%3;w`alrFLzZ#7;n#nz~3&b0WYGP4;TGpXeGRkrDwy1r{1&{ z{X_tZ;jT|-G{Fu8sDghCzY&JLSXu=S&1&8Q7M9M1nce!2gnPdl-g+!}uQ+Tq!29Zf zKfwkBI1N@*?CZovZaS$H5c6MbGoGAn``}N`tvuJkZf(cQNFUq3fS$Z`@%@{;}I5ES=OEgZq3R&oTLT!|z(h|Up2l5Ez!;*pJ9;^rqx0io~l&qEQC zzs_sgS_6-GJE+lJ1RDJ!Td6`qIYzXK6k!X=Kpm|HxR1{0E zaG%-(H7~PsX@;LTrmR2B%iQfH&+)Yo~Sgjvapv5o(MbA4P%K2?xn# zk#kQc4j|r=l=d0)WdDTJP&{d+tPvS>Ku~lI4~ZY(!n>^W=XU=`W6&b+Xv>?`_?eek z1NcIpR3qmH$;po2KoevJd_Xmu)}Q!Cn5d(OT0yxq4+mR*wI zpLdhp``*3ho_Fs#=bn4cz31(>8*1w(=N9E=O_%_jK2zzg7FklJDwe!*OQfZtwl!<= zMOxMzQK(N*bvM?@i#%QCOI4nm1savim<1lp^ekZP$+=fx&|wNpAW&a(%@vVwpBaf+ zW)uprEj=*{KLchYYK2297cU5zu|TkMWku=5rAroc2RkoZSzcORQM%N7KHyVz@JduQ z;n08?ieVdxC`Tcxsb5%Lx(vGU7LqW?PD9RlRF_uDnVlBcs*jJpRNvNC2UpQhFlp1X z;Gcp-zQ6qLqq^p{H7(7))$I+n*mne(h9JS>rbPl4CgPt=11(oS)ZT4H4gB{TNQmSn zU`W7R$W=|wd3txdfSO2%Cu3=_SP&Jj;!uGX`(}#{PW^&xr`@W-X3^QTgPXS{BWNm@ zQ1PB^xbf*=eH->4aAYEh(cu;6vJ{ZZSHMOnwrgSv;EM`?b8F@`3tMU$E}xK36}lEX zO9HUOdgqZ~=tvtOIRSr<84CI%LTCtvx)F`}sVvw&WByqi#Q3JKkXW>7S#VrD`QY~| z?i1t7rY4iHI#}FMFQL3BOfofAA3r9batLawglb5JiXAC?Pdp}1A<(`Cu!kP=deNbQ zEpi$i)U>(OFe^*b)O^8V$TTFMPV}2*oGa!LfB}UQpl}-kXXjFlUPAK=;KUZ5YPqnL zBpuW9^ju7Gm;g*jP=pknAb^VqkSC^HN;R1#@3lkuohgpg@K=3sVpVF0GhDlU_pNu< zw>DgXQb;nF7k@lclqUqweZ@+e368{*+)+5p4K}k|N){KRC>&eLpbzdn{7fOu92N=9 zW-fP$^D)%Ed1B;$FOj{d;LZp&6g02H$x3PO#9zLjnp&uLL}Eo`kJTQttY`fu??(O} zo?1=wxy5yBo_kj+E`GA=08KAPN>*Gml%ivW#d%cAeJk}CQTK?%%4no9j!q&XRowlV z&;9WjUyJ#45RUF6cuaEF&RadJTHc?uck%M< zGFjt(yA*IK;8MV)fJ*_F0^>u0)$NV-xEcRCVaKByI(7WR#zn2Llnw>(=#tI@j=b|o zD{P^HO91E~;BW8#`X+c3F8~iu=mf*`Ua|vb5`>|Sq^x@%hS#W|QttI1c>US^P(-ln zK^$Tq@NIn__7TKf06J=3b1ZNGzKR#1jt&Dm?s^AiQh|sJ00_Kr=$$u@M(&{^H;oMi z1pg;G$n0Q(gUt>xaq>G+UL3OG;ITs_9OAJ5!B1}(n-RHKGNFJBMEH7^CJt3;e(GL2 z@65x4IC;)$I06i?YCDGy@Y%71i^urnjC{&mO?}fLA5{*v>WWxFU$4Kf&kA+FWuHy2O?DN*UiX3#Q{A*yX-(r>`$!^_sq36H#0K z5C%-FQfG}-P%#NAK0^(M(H0sh9!i|cI{+{0Jc=cmEK!^ONUYlz@dvH`XePL--fmi3 z?QNy8>TtHiB}N+x6T-f54-Q9-D-*m({pu=j8}UYmq7m>inKDCs$P|dPPGFFCKWK(R#V+r2i|L>UdWaWnyZ`Y1bmf(cLdDU)Ycn8)D$ghbb5ghSc^oi}=HeTC3LzIy5GAPSWZ5AN~G*NyDCDg1sb6 zoY(e!+QD8d=W~^utIUqK)JLYSu&Q*_Azu90gQ*FJBS9;Knm6hk66#Fu9_KhpF-#fk zrDHH$9-)e}rEyJLJE}iiPViHR9}4dv^~)GI7jsF)e$b3gyeEq@Bg!K@G|l?e3W|Kn znMZe;N8yG7rzZ1tgo4rR3ke3`u3cBXBqI9~A}>2v$hvrZo3Jc+fvn)~qT4rL zSl>{)+77xC6?1pLSTbr76EVh$xI~9U;!>PWngUG&vcq#Zo-)2{Q=}@fGumectWGQ7 z>+Z*CVZ}Dt6SzNQ#bg97v6x`5?^<+ymku+e2nM2N%ohmv_TkivhC{v%f7BE*D4AC^ z?7}SugQ#?LB7{*94fjU^(vp~?9#iX&`HjvQk&;E%?981Za{ zgNA({F+0tO848%b7%D5H+qAm6V@AaE_gI_JfI!tn<%q#epF0JVGyNy(h+2?udt&u6 zxS9%-#qjLgJ1b#7No|>QgoV3(f*{lqz4)y?o8V?DD3-F}`_4_b!(J*-W`S_GcM%Sh za4|W^-9Cp3wn%4s=D*(it#swi!^--NJHP3t^VGaWvTP}OVrr-Kx7pIsMfo(*m z-+d<}nOuY@?c_y^h?cT66X73Y9|zx(AKmb4XB)*G;eK>dF(IvkG^dU(gCN%s-;(Xikf&ss6~P%Eu)41PusYAqFx!99dyPJ&}#P(dokU}NOj`w+xU<3)k=90RtW z5;+QqG6RXy0Euh`2h$U!==kFEo5pKcF5WaLAPr2=5>DCF4)4&@sq86aTkn^VOCe(g6qcZvZ!}4S#&8@SbF=8^qajcCAX?;^%u8T_`vRr zSz**;N^B#`y96ShjuqxKWO+=C+Ra@G2nvi=7s}Vmf1<1{-t=NknBu%kY@66>A}-=%O5TFp=>^~BKTd}kd_8E2Cfe!8F!q4+- zc0vdN&VgyI`uAZzUY5WckNF`Q zZh&O7+}FW=NuJr8c#_?804$Oe8nqqD;>5O)5+3rLM3Q1<|2Tvc4BEO9NlJu0bEE{( zK7vS6;@&vY@l!5w^V6k(u=n`qp$N&F`Rgy!Jk0XRj(sow`3TL!0yAyMc`|P$X&#oM zqLD+zD7uuA=3$}o(2+m?yoBaq;mUXPLs(+GN2PgKKj2{&U9u{#Aof>a$me|Jid$j94 z`#xJn@O@IO__JtW!>stY8i Date: Wed, 7 Jan 2026 15:57:56 +0000 Subject: [PATCH 40/91] Adds tests for specific humidity and mixing ratio conversions --- tests/operators/test_humidity.py | 160 +++++++++++++++++++++++++++++++ 1 file changed, 160 insertions(+) create mode 100644 tests/operators/test_humidity.py diff --git a/tests/operators/test_humidity.py b/tests/operators/test_humidity.py new file mode 100644 index 000000000..d208d3960 --- /dev/null +++ b/tests/operators/test_humidity.py @@ -0,0 +1,160 @@ +# © Crown copyright, Met Office (2022-2026) and CSET contributors. +# +# Licensed under the Apache License, Version 2.0 (the "License"); +# you may not use this file except in compliance with the License. +# You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software +# distributed under the License is distributed on an "AS IS" BASIS, +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +# See the License for the specific language governing permissions and +# limitations under the License. + +"""Test humidity operators.""" + +import cf_units +import iris.cube +import numpy as np + +from CSET.operators import humidity + + +def test_mixing_ratio_from_specific_humidity(specific_humidity_for_conversions_cube): + """Test calculation of mixing ratio from specific humidity.""" + expected_data = specific_humidity_for_conversions_cube / ( + 1 - specific_humidity_for_conversions_cube + ) + assert np.allclose( + expected_data.data, + humidity.mixing_ratio_from_specific_humidity( + specific_humidity_for_conversions_cube + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_mixing_ratio_from_specific_humidity_name( + specific_humidity_for_conversions_cube, +): + """Test naming of mixing ratio cube.""" + expected_name = "mixing_ratio" + assert ( + expected_name + == humidity.mixing_ratio_from_specific_humidity( + specific_humidity_for_conversions_cube + ).name() + ) + + +def test_mixing_ratio_from_specific_humidity_unit( + specific_humidity_for_conversions_cube, +): + """Test units for mixing ratio.""" + expected_units = cf_units.Unit("kg/kg") + assert ( + expected_units + == humidity.mixing_ratio_from_specific_humidity( + specific_humidity_for_conversions_cube + ).units + ) + + +def test_mixing_ratio_from_specific_humidity_cubelist( + specific_humidity_for_conversions_cube, +): + """Test calculation of mixing ratio from specific humidity for a CubeList.""" + expected_data = specific_humidity_for_conversions_cube / ( + 1 - specific_humidity_for_conversions_cube + ) + expected_list = iris.cube.CubeList([expected_data, expected_data]) + input_list = iris.cube.CubeList( + [specific_humidity_for_conversions_cube, specific_humidity_for_conversions_cube] + ) + actual_cubelist = humidity.mixing_ratio_from_specific_humidity(input_list) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) + + +def test_mixing_ratio_from_specific_humidity_using_calculated( + mixing_ratio_for_conversions_cube, +): + """Test mixing ratio calculation using calculated q from mixing ratio.""" + calculated_w = humidity.mixing_ratio_from_specific_humidity( + humidity.specific_humidity_from_mixing_ratio(mixing_ratio_for_conversions_cube) + ) + assert np.allclose( + mixing_ratio_for_conversions_cube.data, calculated_w.data, rtol=1e-6, atol=1e-2 + ) + + +def test_specific_humidity_from_mixing_ratio(mixing_ratio_for_conversions_cube): + """Test specific humidity calculation from mixing ratio.""" + expected_data = mixing_ratio_for_conversions_cube / ( + 1 + mixing_ratio_for_conversions_cube + ) + assert np.allclose( + expected_data.data, + humidity.specific_humidity_from_mixing_ratio( + mixing_ratio_for_conversions_cube + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_specific_humidity_from_mixing_ratio_name(mixing_ratio_for_conversions_cube): + """Test specific humidity cube is named.""" + expected_name = "specific_humidity" + assert ( + expected_name + == humidity.specific_humidity_from_mixing_ratio( + mixing_ratio_for_conversions_cube + ).name() + ) + + +def test_specific_humidity_from_mixing_ratio_units(mixing_ratio_for_conversions_cube): + """Test units of specific humidity.""" + expected_units = cf_units.Unit("kg/kg") + assert ( + expected_units + == humidity.specific_humidity_from_mixing_ratio( + mixing_ratio_for_conversions_cube + ).units + ) + + +def test_specific_humidity_from_mixing_ratio_cubelist( + mixing_ratio_for_conversions_cube, +): + """Test specific humidity calculation from mixing ratio with a CuebList.""" + expected_data = mixing_ratio_for_conversions_cube / ( + 1 + mixing_ratio_for_conversions_cube + ) + expected_list = iris.cube.CubeList([expected_data, expected_data]) + input_list = iris.cube.CubeList( + [mixing_ratio_for_conversions_cube, mixing_ratio_for_conversions_cube] + ) + actual_cubelist = humidity.specific_humidity_from_mixing_ratio(input_list) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) + + +def test_specific_humidity_from_mixing_ratio_using_calculated( + specific_humidity_for_conversions_cube, +): + """Test specific humidity calculation using calculated w from specific humidity.""" + calculated_q = humidity.specific_humidity_from_mixing_ratio( + humidity.mixing_ratio_from_specific_humidity( + specific_humidity_for_conversions_cube + ) + ) + assert np.allclose( + specific_humidity_for_conversions_cube.data, + calculated_q.data, + rtol=1e-6, + atol=1e-2, + ) From 95386b0c9a557f03caf26e648ae91444b427bcff Mon Sep 17 00:00:00 2001 From: daflack Date: Thu, 8 Jan 2026 08:54:33 +0000 Subject: [PATCH 41/91] Adds tests for saturation mixing ratio --- tests/operators/test_humidity.py | 72 +++++++++++++++++++++++++++++++- 1 file changed, 71 insertions(+), 1 deletion(-) diff --git a/tests/operators/test_humidity.py b/tests/operators/test_humidity.py index d208d3960..587cd354b 100644 --- a/tests/operators/test_humidity.py +++ b/tests/operators/test_humidity.py @@ -18,7 +18,7 @@ import iris.cube import numpy as np -from CSET.operators import humidity +from CSET.operators import _atmospheric_constants, humidity, misc, pressure def test_mixing_ratio_from_specific_humidity(specific_humidity_for_conversions_cube): @@ -158,3 +158,73 @@ def test_specific_humidity_from_mixing_ratio_using_calculated( rtol=1e-6, atol=1e-2, ) + + +def test_saturation_mixing_ratio( + temperature_for_conversions_cube, pressure_for_conversions_cube +): + """Test calculation of saturation mixing ratio.""" + expected_data = ( + _atmospheric_constants.EPSILON + * pressure.vapour_pressure(temperature_for_conversions_cube) + ) / ( + (misc.convert_units(pressure_for_conversions_cube, "hPa")) + - pressure.vapour_pressure(temperature_for_conversions_cube) + ) + assert np.allclose( + expected_data.data, + humidity.saturation_mixing_ratio( + temperature_for_conversions_cube, pressure_for_conversions_cube + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_saturation_mixing_ratio_name( + temperature_for_conversions_cube, pressure_for_conversions_cube +): + """Test naming of saturation mixing ratio cube.""" + expected_name = "saturation_mixing_ratio" + assert ( + expected_name + == humidity.saturation_mixing_ratio( + temperature_for_conversions_cube, pressure_for_conversions_cube + ).name() + ) + + +def test_saturation_mixing_ratio_units( + temperature_for_conversions_cube, pressure_for_conversions_cube +): + """Test units of saturation mixing ratio cube.""" + expected_units = cf_units.Unit("kg/kg") + assert ( + expected_units + == humidity.saturation_mixing_ratio( + temperature_for_conversions_cube, pressure_for_conversions_cube + ).units + ) + + +def test_saturation_mixing_ratio_cubelist( + temperature_for_conversions_cube, pressure_for_conversions_cube +): + """Test calculation of saturation mixing ratio as CubeList.""" + expected_data = ( + _atmospheric_constants.EPSILON + * pressure.vapour_pressure(temperature_for_conversions_cube) + ) / ( + (misc.convert_units(pressure_for_conversions_cube, "hPa")) + - pressure.vapour_pressure(temperature_for_conversions_cube) + ) + expected_list = iris.cube.CubeList([expected_data, expected_data]) + temperature_list = iris.cube.CubeList( + [temperature_for_conversions_cube, temperature_for_conversions_cube] + ) + pressure_list = iris.cube.CubeList( + [pressure_for_conversions_cube, pressure_for_conversions_cube] + ) + actual_cubelist = humidity.saturation_mixing_ratio(temperature_list, pressure_list) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From 00e4039cc1aae18c584532c2146c63e0355fa187 Mon Sep 17 00:00:00 2001 From: daflack Date: Thu, 8 Jan 2026 09:08:32 +0000 Subject: [PATCH 42/91] Adds tests for saturation specific humidity --- tests/operators/test_humidity.py | 66 ++++++++++++++++++++++++++++++++ 1 file changed, 66 insertions(+) diff --git a/tests/operators/test_humidity.py b/tests/operators/test_humidity.py index 587cd354b..d2a0a6a29 100644 --- a/tests/operators/test_humidity.py +++ b/tests/operators/test_humidity.py @@ -228,3 +228,69 @@ def test_saturation_mixing_ratio_cubelist( actual_cubelist = humidity.saturation_mixing_ratio(temperature_list, pressure_list) for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) + + +def test_saturation_specific_humidity( + temperature_for_conversions_cube, pressure_for_conversions_cube +): + """Test calculation of saturation specific humidity.""" + expected_data = ( + _atmospheric_constants.EPSILON + * pressure.vapour_pressure(temperature_for_conversions_cube) + ) / misc.convert_units(pressure_for_conversions_cube, "hPa") + assert np.allclose( + expected_data.data, + humidity.saturation_specific_humidity( + temperature_for_conversions_cube, pressure_for_conversions_cube + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_saturation_specific_humidity_name( + temperature_for_conversions_cube, pressure_for_conversions_cube +): + """Test naming of saturation specific humidity cube.""" + expected_name = "saturation_specific_humidity" + assert ( + expected_name + == humidity.saturation_specific_humidity( + temperature_for_conversions_cube, pressure_for_conversions_cube + ).name() + ) + + +def test_saturation_specific_humidity_units( + temperature_for_conversions_cube, pressure_for_conversions_cube +): + """Test units of saturation specific humidity.""" + expected_units = cf_units.Unit("kg/kg") + assert ( + expected_units + == humidity.saturation_specific_humidity( + temperature_for_conversions_cube, pressure_for_conversions_cube + ).units + ) + + +def test_saturation_specific_humidity_cubelist( + temperature_for_conversions_cube, pressure_for_conversions_cube +): + """Test calculation of saturation specific humidity for a CubeList.""" + expected_data = ( + _atmospheric_constants.EPSILON + * pressure.vapour_pressure(temperature_for_conversions_cube) + ) / misc.convert_units(pressure_for_conversions_cube, "hPa") + expected_list = iris.cube.CubeList([expected_data, expected_data]) + temperature_list = iris.cube.CubeList( + [temperature_for_conversions_cube, temperature_for_conversions_cube] + ) + pressure_list = iris.cube.CubeList( + [pressure_for_conversions_cube, pressure_for_conversions_cube] + ) + actual_cubelist = humidity.saturation_specific_humidity( + temperature_list, pressure_list + ) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From 2c7c2671a355f51138f16d2b2babf7212e43175b Mon Sep 17 00:00:00 2001 From: daflack Date: Thu, 8 Jan 2026 11:53:05 +0000 Subject: [PATCH 43/91] Adds tests for mixing ratio from relative humidity --- tests/operators/test_humidity.py | 81 ++++++++++++++++++++++++++++++++ 1 file changed, 81 insertions(+) diff --git a/tests/operators/test_humidity.py b/tests/operators/test_humidity.py index d2a0a6a29..f235f4320 100644 --- a/tests/operators/test_humidity.py +++ b/tests/operators/test_humidity.py @@ -294,3 +294,84 @@ def test_saturation_specific_humidity_cubelist( ) for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) + + +def test_mixing_ratio_from_relative_humidity( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, +): + """Test calculation of mixing ratio from relative humidity.""" + expected_data = humidity.saturation_mixing_ratio( + temperature_for_conversions_cube, pressure_for_conversions_cube + ) * misc.convert_units(relative_humidity_for_conversions_cube, "1") + assert np.allclose( + expected_data.data, + humidity.mixing_ratio_from_relative_humidity( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_mixing_ratio_from_relative_humidity_name( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, +): + """Test naming of mixing ratio cube.""" + expected_name = "mixing_ratio" + assert ( + expected_name + == humidity.mixing_ratio_from_relative_humidity( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, + ).name() + ) + + +def test_mixing_ratio_from_relative_humidity_units( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, +): + """Test units of mixing ratio cube.""" + expected_units = cf_units.Unit("kg/kg") + assert ( + expected_units + == humidity.mixing_ratio_from_relative_humidity( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, + ).units + ) + + +def test_mixing_ratio_from_relative_humidity_cubelist( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, +): + """Test calculation of mixing ration from relative humidity for a CubeList.""" + expected_data = humidity.saturation_mixing_ratio( + temperature_for_conversions_cube, pressure_for_conversions_cube + ) * misc.convert_units(relative_humidity_for_conversions_cube, "1") + expected_list = iris.cube.CubeList([expected_data, expected_data]) + temperature_list = iris.cube.CubeList( + [temperature_for_conversions_cube, temperature_for_conversions_cube] + ) + pressure_list = iris.cube.CubeList( + [pressure_for_conversions_cube, pressure_for_conversions_cube] + ) + rh_list = iris.cube.CubeList( + [relative_humidity_for_conversions_cube, relative_humidity_for_conversions_cube] + ) + actual_cubelist = humidity.mixing_ratio_from_relative_humidity( + temperature_list, pressure_list, rh_list + ) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From bdd90eb8f5d041db3530c8a5a51b77cc0d8bcfa2 Mon Sep 17 00:00:00 2001 From: daflack Date: Thu, 8 Jan 2026 12:00:53 +0000 Subject: [PATCH 44/91] Adds tests for specific humidity from relative humidity --- tests/operators/test_humidity.py | 83 +++++++++++++++++++++++++++++++- 1 file changed, 82 insertions(+), 1 deletion(-) diff --git a/tests/operators/test_humidity.py b/tests/operators/test_humidity.py index f235f4320..d7b272941 100644 --- a/tests/operators/test_humidity.py +++ b/tests/operators/test_humidity.py @@ -356,7 +356,7 @@ def test_mixing_ratio_from_relative_humidity_cubelist( pressure_for_conversions_cube, relative_humidity_for_conversions_cube, ): - """Test calculation of mixing ration from relative humidity for a CubeList.""" + """Test calculation of mixing ratio from relative humidity for a CubeList.""" expected_data = humidity.saturation_mixing_ratio( temperature_for_conversions_cube, pressure_for_conversions_cube ) * misc.convert_units(relative_humidity_for_conversions_cube, "1") @@ -375,3 +375,84 @@ def test_mixing_ratio_from_relative_humidity_cubelist( ) for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) + + +def test_specific_humidity_from_relative_humidity( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, +): + """Test calculation of specific humidity from relative humidity.""" + expected_data = humidity.saturation_specific_humidity( + temperature_for_conversions_cube, pressure_for_conversions_cube + ) * misc.convert_units(relative_humidity_for_conversions_cube, "1") + assert np.allclose( + expected_data.data, + humidity.specific_humidity_from_relative_humidity( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_specific_humidity_from_relative_humidity_name( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, +): + """Test naming of specific humidity cube.""" + expected_name = "specific_humidity" + assert ( + expected_name + == humidity.specific_humidity_from_relative_humidity( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, + ).name() + ) + + +def test_specific_humidity_from_relative_humidity_units( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, +): + """Test units of specific humidity cube.""" + expected_units = cf_units.Unit("kg/kg") + assert ( + expected_units + == humidity.specific_humidity_from_relative_humidity( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, + ).units + ) + + +def test_specific_humidity_from_relative_humidity_cubelist( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, +): + """Test calculation of specific humidity from relative humidity for a CubeList.""" + expected_data = humidity.saturation_specific_humidity( + temperature_for_conversions_cube, pressure_for_conversions_cube + ) * misc.convert_units(relative_humidity_for_conversions_cube, "1") + expected_list = iris.cube.CubeList([expected_data, expected_data]) + temperature_list = iris.cube.CubeList( + [temperature_for_conversions_cube, temperature_for_conversions_cube] + ) + pressure_list = iris.cube.CubeList( + [pressure_for_conversions_cube, pressure_for_conversions_cube] + ) + rh_list = iris.cube.CubeList( + [relative_humidity_for_conversions_cube, relative_humidity_for_conversions_cube] + ) + actual_cubelist = humidity.specific_humidity_from_relative_humidity( + temperature_list, pressure_list, rh_list + ) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From 6e8ed4f0e6e53f5b2e5896bea93ff0f0143f7800 Mon Sep 17 00:00:00 2001 From: daflack Date: Thu, 8 Jan 2026 13:24:53 +0000 Subject: [PATCH 45/91] Adds tests for relative humidity from mixing ratio --- tests/operators/test_humidity.py | 70 ++++++++++++++++++++++++++++++++ 1 file changed, 70 insertions(+) diff --git a/tests/operators/test_humidity.py b/tests/operators/test_humidity.py index d7b272941..e41460ad3 100644 --- a/tests/operators/test_humidity.py +++ b/tests/operators/test_humidity.py @@ -377,6 +377,29 @@ def test_mixing_ratio_from_relative_humidity_cubelist( assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) +def test_mixing_ratio_from_relative_humidity_calculated( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + mixing_ratio_for_conversions_cube, +): + """Test mixing ratio from relative humidity using calculated relative humidity.""" + calculated_w = humidity.mixing_ratio_from_relative_humidity( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + humidity.relative_humidity_from_mixing_ratio( + mixing_ratio_for_conversions_cube, + temperature_for_conversions_cube, + pressure_for_conversions_cube, + ), + ) + assert np.allclose( + calculated_w.data, + mixing_ratio_for_conversions_cube.data, + rtol=1e-6, + atol=1e-2, + ) + + def test_specific_humidity_from_relative_humidity( temperature_for_conversions_cube, pressure_for_conversions_cube, @@ -456,3 +479,50 @@ def test_specific_humidity_from_relative_humidity_cubelist( ) for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) + + +def test_relative_humidity_from_mixing_ratio( + mixing_ratio_for_conversions_cube, + temperature_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test calculation of relative humidity from mixing ratio.""" + expected_data = ( + mixing_ratio_for_conversions_cube + / humidity.saturation_mixing_ratio( + temperature_for_conversions_cube, pressure_for_conversions_cube + ) + ) + assert np.allclose( + expected_data.data, + humidity.relative_humidity_from_mixing_ratio( + mixing_ratio_for_conversions_cube, + temperature_for_conversions_cube, + pressure_for_conversions_cube, + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_relative_humidity_from_mixing_ratio_calculated( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, +): + """Test relative humidity from mixing ratio using calculated mixing ratio.""" + calculated_rh = humidity.relative_humidity_from_mixing_ratio( + humidity.mixing_ratio_from_relative_humidity( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, + ), + temperature_for_conversions_cube, + pressure_for_conversions_cube, + ) + assert np.allclose( + calculated_rh.data, + misc.convert_units(relative_humidity_for_conversions_cube, "1").data, + rtol=1e-6, + atol=1e-2, + ) From 3a5ec8d2e15036006fd9f3c360d911d6e59b2f3d Mon Sep 17 00:00:00 2001 From: daflack Date: Thu, 8 Jan 2026 13:55:30 +0000 Subject: [PATCH 46/91] Update tests for relative humidity from mixing ratio and correct RH units --- src/CSET/operators/humidity.py | 2 + tests/operators/test_humidity.py | 67 +++++++++++++++++++++++++++++++- 2 files changed, 67 insertions(+), 2 deletions(-) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index a0305beaf..9ef248d9e 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -153,6 +153,7 @@ def relative_humidity_from_mixing_ratio( ): rel_h = W / saturation_mixing_ratio(T, P) rel_h.rename("relative_humidity") + rel_h = convert_units(rel_h, "%") RH.append(rel_h) if len(RH) == 1: return RH[0] @@ -175,6 +176,7 @@ def relative_humidity_from_specific_humidity( ): rel_h = Q / saturation_specific_humidity(T, P) rel_h.rename("relative_humidity") + rel_h = convert_units(rel_h, "%") RH.append(rel_h) if len(RH) == 1: return RH[0] diff --git a/tests/operators/test_humidity.py b/tests/operators/test_humidity.py index e41460ad3..16df29b6c 100644 --- a/tests/operators/test_humidity.py +++ b/tests/operators/test_humidity.py @@ -487,7 +487,7 @@ def test_relative_humidity_from_mixing_ratio( pressure_for_conversions_cube, ): """Test calculation of relative humidity from mixing ratio.""" - expected_data = ( + expected_data = 100.0 * ( mixing_ratio_for_conversions_cube / humidity.saturation_mixing_ratio( temperature_for_conversions_cube, pressure_for_conversions_cube @@ -522,7 +522,70 @@ def test_relative_humidity_from_mixing_ratio_calculated( ) assert np.allclose( calculated_rh.data, - misc.convert_units(relative_humidity_for_conversions_cube, "1").data, + relative_humidity_for_conversions_cube.data, rtol=1e-6, atol=1e-2, ) + + +def test_relative_humidity_from_mixing_ratio_name( + mixing_ratio_for_conversions_cube, + temperature_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test name of relative humidity cube.""" + expected_name = "relative_humidity" + assert ( + expected_name + == humidity.relative_humidity_from_mixing_ratio( + mixing_ratio_for_conversions_cube, + temperature_for_conversions_cube, + pressure_for_conversions_cube, + ).name() + ) + + +def test_relative_humidity_from_mixing_ratio_units( + mixing_ratio_for_conversions_cube, + temperature_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test units of relative humidity cube.""" + expected_units = cf_units.Unit("%") + assert ( + expected_units + == humidity.relative_humidity_from_mixing_ratio( + mixing_ratio_for_conversions_cube, + temperature_for_conversions_cube, + pressure_for_conversions_cube, + ).units + ) + + +def test_relative_humidity_from_mixing_ratio_cubelist( + mixing_ratio_for_conversions_cube, + temperature_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test calculation of relative humidity from mixing ratio for a CubeList.""" + expected_data = 100.0 * ( + mixing_ratio_for_conversions_cube + / humidity.saturation_mixing_ratio( + temperature_for_conversions_cube, pressure_for_conversions_cube + ) + ) + expected_list = iris.cube.CubeList([expected_data, expected_data]) + mixing_ratio_list = iris.cube.CubeList( + [mixing_ratio_for_conversions_cube, mixing_ratio_for_conversions_cube] + ) + temperature_list = iris.cube.CubeList( + [temperature_for_conversions_cube, temperature_for_conversions_cube] + ) + pressure_list = iris.cube.CubeList( + [pressure_for_conversions_cube, pressure_for_conversions_cube] + ) + actual_cubelist = humidity.relative_humidity_from_mixing_ratio( + mixing_ratio_list, temperature_list, pressure_list + ) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From 4f6522054050385a7448b505e5db7f32128e8991 Mon Sep 17 00:00:00 2001 From: daflack Date: Thu, 8 Jan 2026 14:05:28 +0000 Subject: [PATCH 47/91] Adds tests for relative humidity from specific humidity --- tests/operators/test_humidity.py | 133 +++++++++++++++++++++++++++++++ 1 file changed, 133 insertions(+) diff --git a/tests/operators/test_humidity.py b/tests/operators/test_humidity.py index 16df29b6c..9ff94d14d 100644 --- a/tests/operators/test_humidity.py +++ b/tests/operators/test_humidity.py @@ -481,6 +481,29 @@ def test_specific_humidity_from_relative_humidity_cubelist( assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) +def test_specific_humidity_from_relative_humidity_calculated( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + specific_humidity_for_conversions_cube, +): + """Test specific humidity from relative humidity using calculated relative humidity.""" + calculated_q = humidity.specific_humidity_from_relative_humidity( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + humidity.relative_humidity_from_specific_humidity( + specific_humidity_for_conversions_cube, + temperature_for_conversions_cube, + pressure_for_conversions_cube, + ), + ) + assert np.allclose( + calculated_q.data, + specific_humidity_for_conversions_cube.data, + rtol=1e-6, + atol=1e-2, + ) + + def test_relative_humidity_from_mixing_ratio( mixing_ratio_for_conversions_cube, temperature_for_conversions_cube, @@ -589,3 +612,113 @@ def test_relative_humidity_from_mixing_ratio_cubelist( ) for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) + + +def test_relative_humidity_from_specific_humidity( + specific_humidity_for_conversions_cube, + temperature_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test calculation of relative humidity from specific humidity.""" + expected_data = 100.0 * ( + specific_humidity_for_conversions_cube + / humidity.saturation_specific_humidity( + temperature_for_conversions_cube, pressure_for_conversions_cube + ) + ) + assert np.allclose( + expected_data.data, + humidity.relative_humidity_from_specific_humidity( + specific_humidity_for_conversions_cube, + temperature_for_conversions_cube, + pressure_for_conversions_cube, + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_relative_humidity_from_specific_humidity_calculated( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, +): + """Test relative humidity from specific humidity using calculated specific humidity.""" + calculated_rh = humidity.relative_humidity_from_specific_humidity( + humidity.specific_humidity_from_relative_humidity( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + relative_humidity_for_conversions_cube, + ), + temperature_for_conversions_cube, + pressure_for_conversions_cube, + ) + assert np.allclose( + calculated_rh.data, + relative_humidity_for_conversions_cube.data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_relative_humidity_from_specific_humidity_name( + specific_humidity_for_conversions_cube, + temperature_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test name of relative humidity cube.""" + expected_name = "relative_humidity" + assert ( + expected_name + == humidity.relative_humidity_from_specific_humidity( + specific_humidity_for_conversions_cube, + temperature_for_conversions_cube, + pressure_for_conversions_cube, + ).name() + ) + + +def test_relative_humidity_from_specific_humidity_units( + specific_humidity_for_conversions_cube, + temperature_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test units of relative humidity cube.""" + expected_units = cf_units.Unit("%") + assert ( + expected_units + == humidity.relative_humidity_from_specific_humidity( + specific_humidity_for_conversions_cube, + temperature_for_conversions_cube, + pressure_for_conversions_cube, + ).units + ) + + +def test_relative_humidity_from_specific_humidity_cubelist( + specific_humidity_for_conversions_cube, + temperature_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test calculation of relative humidity from specific humidity for a CubeList.""" + expected_data = 100.0 * ( + specific_humidity_for_conversions_cube + / humidity.saturation_specific_humidity( + temperature_for_conversions_cube, pressure_for_conversions_cube + ) + ) + expected_list = iris.cube.CubeList([expected_data, expected_data]) + specific_humidity_list = iris.cube.CubeList( + [specific_humidity_for_conversions_cube, specific_humidity_for_conversions_cube] + ) + temperature_list = iris.cube.CubeList( + [temperature_for_conversions_cube, temperature_for_conversions_cube] + ) + pressure_list = iris.cube.CubeList( + [pressure_for_conversions_cube, pressure_for_conversions_cube] + ) + actual_cubelist = humidity.relative_humidity_from_specific_humidity( + specific_humidity_list, temperature_list, pressure_list + ) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From 21671c881d68a4426223f315941b7218f7549c35 Mon Sep 17 00:00:00 2001 From: daflack Date: Thu, 8 Jan 2026 14:43:44 +0000 Subject: [PATCH 48/91] Adds tests and fixes calculation for dewpoint temperature --- src/CSET/operators/temperature.py | 18 +++-- tests/operators/test_temperature.py | 107 ++++++++++++++++++++++++++++ 2 files changed, 118 insertions(+), 7 deletions(-) create mode 100644 tests/operators/test_temperature.py diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index c1221bc5b..f432cff83 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -24,7 +24,10 @@ saturation_mixing_ratio, ) from CSET.operators.misc import convert_units -from CSET.operators.pressure import exner_pressure, vapour_pressure +from CSET.operators.pressure import ( + exner_pressure, + vapour_pressure_from_relative_humidity, +) def dewpoint_temperature( @@ -36,15 +39,16 @@ def dewpoint_temperature( for T, RH in zip( iter_maybe(temperature), iter_maybe(relative_humidity), strict=True ): - vp = vapour_pressure(T, RH) + vp = vapour_pressure_from_relative_humidity(T, RH) td = vp.copy() - td.data = (243.5 * np.log(vp.core_data()) - 440.8) / ( + td.data = ((243.5 * np.log(vp.core_data())) - 440.8) / ( 19.48 - np.log(vp.core_data()) ) - td.data[td.data - T0 < -35.0] = np.nan - td.data[td.data - T0 > 35.0] = np.nan - td.units = "K" - td.rename("dewpoint temperature") + td.data[T.data - T0 < -35.0] = np.nan + td.data[T.data - T0 > 35.0] = np.nan + td.units = "Celsius" + td.rename("dewpoint_temperature") + td = convert_units(td, "K") Td.append(td) if len(Td) == 1: return Td[0] diff --git a/tests/operators/test_temperature.py b/tests/operators/test_temperature.py new file mode 100644 index 000000000..e12b4f8d1 --- /dev/null +++ b/tests/operators/test_temperature.py @@ -0,0 +1,107 @@ +# © Crown copyright, Met Office (2022-2026) and CSET contributors. +# +# Licensed under the Apache License, Version 2.0 (the "License"); +# you may not use this file except in compliance with the License. +# You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software +# distributed under the License is distributed on an "AS IS" BASIS, +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +# See the License for the specific language governing permissions and +# limitations under the License. + +"""Test temperature operators.""" + +import cf_units +import iris.cube +import numpy as np + +from CSET.operators import _atmospheric_constants, misc, pressure, temperature + + +def test_dewpoint_temperature( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube +): + """Test calculation of dewpoint temperature.""" + vp = pressure.vapour_pressure_from_relative_humidity( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube + ) + expected_data = vp.copy() + expected_data.data = (243.5 * np.log(vp.core_data()) - 440.8) / ( + 19.48 - np.log(vp.core_data()) + ) + expected_data.data[ + temperature_for_conversions_cube.data - _atmospheric_constants.T0 < -35.0 + ] = np.nan + expected_data.data[ + temperature_for_conversions_cube.data - _atmospheric_constants.T0 > 35.0 + ] = np.nan + expected_data.units = "Celsius" + expected_data = misc.convert_units(expected_data, "K") + assert np.allclose( + expected_data.data, + temperature.dewpoint_temperature( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_dewpoint_temperature_name( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube +): + """Test name of dewpoint temperature cube.""" + expected_name = "dewpoint_temperature" + assert ( + expected_name + == temperature.dewpoint_temperature( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube + ).name() + ) + + +def test_dewpoint_temperature_unit( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube +): + """Test units for dewpoint temperature cube.""" + expected_units = cf_units.Unit("K") + assert ( + expected_units + == temperature.dewpoint_temperature( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube + ).units + ) + + +def test_dewpoint_temperature_cubelist( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube +): + """Test calculation of dewpoint temperature for a CubeList.""" + vp = pressure.vapour_pressure_from_relative_humidity( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube + ) + expected_data = vp.copy() + expected_data.data = (243.5 * np.log(vp.core_data()) - 440.8) / ( + 19.48 - np.log(vp.core_data()) + ) + expected_data.data[ + temperature_for_conversions_cube.data - _atmospheric_constants.T0 < -35.0 + ] = np.nan + expected_data.data[ + temperature_for_conversions_cube.data - _atmospheric_constants.T0 > 35.0 + ] = np.nan + expected_data.units = "Celsius" + expected_data = misc.convert_units(expected_data, "K") + expected_list = iris.cube.CubeList([expected_data, expected_data]) + temperature_list = iris.cube.CubeList( + [temperature_for_conversions_cube, temperature_for_conversions_cube] + ) + rh_list = iris.cube.CubeList( + [relative_humidity_for_conversions_cube, relative_humidity_for_conversions_cube] + ) + actual_cubelist = temperature.dewpoint_temperature(temperature_list, rh_list) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From 0bbed29db8862edf9851868173c175318b947358 Mon Sep 17 00:00:00 2001 From: daflack Date: Thu, 8 Jan 2026 16:47:28 +0000 Subject: [PATCH 49/91] Adds tests for virtual temperature --- tests/operators/test_temperature.py | 66 +++++++++++++++++++++++++++++ 1 file changed, 66 insertions(+) diff --git a/tests/operators/test_temperature.py b/tests/operators/test_temperature.py index e12b4f8d1..670c51185 100644 --- a/tests/operators/test_temperature.py +++ b/tests/operators/test_temperature.py @@ -105,3 +105,69 @@ def test_dewpoint_temperature_cubelist( actual_cubelist = temperature.dewpoint_temperature(temperature_list, rh_list) for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) + + +def test_virtual_temperature( + temperature_for_conversions_cube, mixing_ratio_for_conversions_cube +): + """Test to calculate virtual temperature.""" + expected_data = temperature_for_conversions_cube * ( + (mixing_ratio_for_conversions_cube + _atmospheric_constants.EPSILON) + / (_atmospheric_constants.EPSILON * (1 + mixing_ratio_for_conversions_cube)) + ) + assert np.allclose( + expected_data.data, + temperature.virtual_temperature( + temperature_for_conversions_cube, mixing_ratio_for_conversions_cube + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_virtual_temperature_name( + temperature_for_conversions_cube, mixing_ratio_for_conversions_cube +): + """Test name of virtual temperature cube.""" + expected_name = "virtual_temperature" + assert ( + expected_name + == temperature.virtual_temperature( + temperature_for_conversions_cube, mixing_ratio_for_conversions_cube + ).name() + ) + + +def test_virtual_temperature_units( + temperature_for_conversions_cube, mixing_ratio_for_conversions_cube +): + """Test units of virtual temperature cube.""" + expected_units = cf_units.Unit("K") + assert ( + expected_units + == temperature.virtual_temperature( + temperature_for_conversions_cube, mixing_ratio_for_conversions_cube + ).units + ) + + +def test_virtual_temperature_cubelist( + temperature_for_conversions_cube, mixing_ratio_for_conversions_cube +): + """Test to calculate virtual temperature for a CubeList.""" + expected_data = temperature_for_conversions_cube * ( + (mixing_ratio_for_conversions_cube + _atmospheric_constants.EPSILON) + / (_atmospheric_constants.EPSILON * (1 + mixing_ratio_for_conversions_cube)) + ) + expected_list = iris.cube.CubeList([expected_data, expected_data]) + temperature_list = iris.cube.CubeList( + [temperature_for_conversions_cube, temperature_for_conversions_cube] + ) + mixing_ratio_list = iris.cube.CubeList( + [mixing_ratio_for_conversions_cube, mixing_ratio_for_conversions_cube] + ) + actual_cubelist = temperature.virtual_temperature( + temperature_list, mixing_ratio_list + ) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From 201a1871c52a36d75c61fa726a1ec9b46b67612a Mon Sep 17 00:00:00 2001 From: daflack Date: Tue, 13 Jan 2026 13:59:01 +0000 Subject: [PATCH 50/91] Adds tests for wet-bulb temperature conversions --- tests/operators/test_temperature.py | 80 +++++++++++++++++++++++++++++ 1 file changed, 80 insertions(+) diff --git a/tests/operators/test_temperature.py b/tests/operators/test_temperature.py index 670c51185..c20ee5fe2 100644 --- a/tests/operators/test_temperature.py +++ b/tests/operators/test_temperature.py @@ -171,3 +171,83 @@ def test_virtual_temperature_cubelist( ) for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) + + +def test_wet_bulb_temperature( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube +): + """Test to calculate wet-bulb temperature.""" + T_units = misc.convert_units(temperature_for_conversions_cube, "Celsius") + RH_units = misc.convert_units(relative_humidity_for_conversions_cube, "%") + TW = ( + T_units * np.arctan(0.151977 * (RH_units.core_data() + 8.313659) ** 0.5) + + np.arctan(T_units.core_data() + RH_units.core_data()) + - np.arctan(RH_units.core_data() - 1.676331) + + 0.00391838 + * (RH_units.core_data()) ** (3.0 / 2.0) + * np.arctan(0.023101 * RH_units.core_data()) + - 4.686035 + ) + TW = misc.convert_units(TW, "K") + assert np.allclose( + TW.data, + temperature.wet_bulb_temperature( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_wet_bulb_temperature_units( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube +): + """Test units for wet bulb temperature cube.""" + expected_units = cf_units.Unit("K") + assert ( + expected_units + == temperature.wet_bulb_temperature( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube + ).units + ) + + +def test_wet_bulb_temperature_name( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube +): + """Test name for wet bulb temperature cube.""" + expected_name = "wet_bulb_temperature" + assert ( + expected_name + == temperature.wet_bulb_temperature( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube + ).name() + ) + + +def test_wet_bulb_temperature_cubelist( + temperature_for_conversions_cube, relative_humidity_for_conversions_cube +): + """Test to calculate wet-bulb temperature for a CubeList.""" + T_units = misc.convert_units(temperature_for_conversions_cube, "Celsius") + RH_units = misc.convert_units(relative_humidity_for_conversions_cube, "%") + TW = ( + T_units * np.arctan(0.151977 * (RH_units.core_data() + 8.313659) ** 0.5) + + np.arctan(T_units.core_data() + RH_units.core_data()) + - np.arctan(RH_units.core_data() - 1.676331) + + 0.00391838 + * (RH_units.core_data()) ** (3.0 / 2.0) + * np.arctan(0.023101 * RH_units.core_data()) + - 4.686035 + ) + TW = misc.convert_units(TW, "K") + expected_list = iris.cube.CubeList([TW, TW]) + temperature_list = iris.cube.CubeList( + [temperature_for_conversions_cube, temperature_for_conversions_cube] + ) + rh_list = iris.cube.CubeList( + [relative_humidity_for_conversions_cube, relative_humidity_for_conversions_cube] + ) + actual_cubelist = temperature.wet_bulb_temperature(temperature_list, rh_list) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From c4e383cf5a747a32fcddef613edf075c35e95a6b Mon Sep 17 00:00:00 2001 From: daflack Date: Wed, 14 Jan 2026 09:08:54 +0000 Subject: [PATCH 51/91] Adds tests for potenital temperature --- tests/operators/test_temperature.py | 62 +++++++++++++++++++++++++++++ 1 file changed, 62 insertions(+) diff --git a/tests/operators/test_temperature.py b/tests/operators/test_temperature.py index c20ee5fe2..d8eb851f9 100644 --- a/tests/operators/test_temperature.py +++ b/tests/operators/test_temperature.py @@ -251,3 +251,65 @@ def test_wet_bulb_temperature_cubelist( actual_cubelist = temperature.wet_bulb_temperature(temperature_list, rh_list) for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) + + +def test_potential_temperature( + temperature_for_conversions_cube, pressure_for_conversions_cube +): + """Test to calculate potenital temperature.""" + expected_data = temperature_for_conversions_cube / pressure.exner_pressure( + pressure_for_conversions_cube + ) + assert np.allclose( + expected_data.data, + temperature.potential_temperature( + temperature_for_conversions_cube, pressure_for_conversions_cube + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_potenital_temperature_name( + temperature_for_conversions_cube, pressure_for_conversions_cube +): + """Test name of potential temperature cube.""" + expected_name = "potential_temperature" + assert ( + expected_name + == temperature.potential_temperature( + temperature_for_conversions_cube, pressure_for_conversions_cube + ).name() + ) + + +def test_potential_temperature_units( + temperature_for_conversions_cube, pressure_for_conversions_cube +): + """Test units of potential temperature cueb.""" + expected_units = cf_units.Unit("K") + assert ( + expected_units + == temperature.potential_temperature( + temperature_for_conversions_cube, pressure_for_conversions_cube + ).units + ) + + +def test_potenital_temperature_cubelist( + temperature_for_conversions_cube, pressure_for_conversions_cube +): + """Test calculation of potenital temperature for a CubeList.""" + expected_data = temperature_for_conversions_cube / pressure.exner_pressure( + pressure_for_conversions_cube + ) + expected_list = iris.cube.CubeList([expected_data, expected_data]) + temperature_list = iris.cube.CubeList( + [temperature_for_conversions_cube, temperature_for_conversions_cube] + ) + pressure_list = iris.cube.CubeList( + [pressure_for_conversions_cube, pressure_for_conversions_cube] + ) + actual_cubelist = temperature.potential_temperature(temperature_list, pressure_list) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From 4c59f4118f70c203a0e9a67d98b4335462ed3eb0 Mon Sep 17 00:00:00 2001 From: daflack Date: Wed, 14 Jan 2026 09:33:04 +0000 Subject: [PATCH 52/91] Adds tests for virtual potential temperature --- tests/operators/test_temperature.py | 81 +++++++++++++++++++++++++++++ 1 file changed, 81 insertions(+) diff --git a/tests/operators/test_temperature.py b/tests/operators/test_temperature.py index d8eb851f9..aa49d45e4 100644 --- a/tests/operators/test_temperature.py +++ b/tests/operators/test_temperature.py @@ -313,3 +313,84 @@ def test_potenital_temperature_cubelist( actual_cubelist = temperature.potential_temperature(temperature_list, pressure_list) for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) + + +def test_virtual_potential_temperature( + temperature_for_conversions_cube, + mixing_ratio_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test calculation of virtual potential temperature.""" + expected_data = temperature.virtual_temperature( + temperature_for_conversions_cube, mixing_ratio_for_conversions_cube + ) / pressure.exner_pressure(pressure_for_conversions_cube) + assert np.allclose( + expected_data.data, + temperature.virtual_potential_temperature( + temperature_for_conversions_cube, + mixing_ratio_for_conversions_cube, + pressure_for_conversions_cube, + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_virtual_potential_temperature_name( + temperature_for_conversions_cube, + mixing_ratio_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test name of virtual potential temperature cube.""" + expected_name = "virtual_potential_temperature" + assert ( + expected_name + == temperature.virtual_potential_temperature( + temperature_for_conversions_cube, + mixing_ratio_for_conversions_cube, + pressure_for_conversions_cube, + ).name() + ) + + +def test_virtual_potential_temperature_units( + temperature_for_conversions_cube, + mixing_ratio_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test units of virtual potenital temperature cube.""" + expected_units = cf_units.Unit("K") + assert ( + expected_units + == temperature.virtual_potential_temperature( + temperature_for_conversions_cube, + mixing_ratio_for_conversions_cube, + pressure_for_conversions_cube, + ).units + ) + + +def test_virtual_potential_temperature_cubelist( + temperature_for_conversions_cube, + mixing_ratio_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test calculation of virtual potenital temperature for a CubeList.""" + expected_data = temperature.virtual_temperature( + temperature_for_conversions_cube, mixing_ratio_for_conversions_cube + ) / pressure.exner_pressure(pressure_for_conversions_cube) + expected_list = iris.cube.CubeList([expected_data, expected_data]) + temperature_list = iris.cube.CubeList( + [temperature_for_conversions_cube, temperature_for_conversions_cube] + ) + mixing_ratio_list = iris.cube.CubeList( + [mixing_ratio_for_conversions_cube, mixing_ratio_for_conversions_cube] + ) + pressure_list = iris.cube.CubeList( + [pressure_for_conversions_cube, pressure_for_conversions_cube] + ) + actual_cubelist = temperature.virtual_potential_temperature( + temperature_list, mixing_ratio_list, pressure_list + ) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From a0a3fe1b25958ff313a250220472667aefa3576d Mon Sep 17 00:00:00 2001 From: daflack Date: Thu, 15 Jan 2026 08:52:31 +0000 Subject: [PATCH 53/91] Adds tests for equivalent potential temperature --- tests/operators/test_temperature.py | 109 +++++++++++++++++++++++++++- 1 file changed, 108 insertions(+), 1 deletion(-) diff --git a/tests/operators/test_temperature.py b/tests/operators/test_temperature.py index aa49d45e4..a7afd0259 100644 --- a/tests/operators/test_temperature.py +++ b/tests/operators/test_temperature.py @@ -18,7 +18,7 @@ import iris.cube import numpy as np -from CSET.operators import _atmospheric_constants, misc, pressure, temperature +from CSET.operators import _atmospheric_constants, humidity, misc, pressure, temperature def test_dewpoint_temperature( @@ -394,3 +394,110 @@ def test_virtual_potential_temperature_cubelist( ) for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) + + +def test_equivalent_potential_temperature( + temperature_for_conversions_cube, + relative_humidity_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test calculation of equivalent potential temperature.""" + RH = misc.convert_units(relative_humidity_for_conversions_cube, "1") + theta = temperature.potential_temperature( + temperature_for_conversions_cube, pressure_for_conversions_cube + ) + w = humidity.mixing_ratio_from_relative_humidity( + temperature_for_conversions_cube, pressure_for_conversions_cube, RH + ) + second_term_power = -(w * _atmospheric_constants.RV) / _atmospheric_constants.CPD + second_term = RH.core_data() ** second_term_power.core_data() + third_term_power = ( + _atmospheric_constants.LV + * w + / (_atmospheric_constants.CPD * temperature_for_conversions_cube) + ) + third_term = np.exp(third_term_power.core_data()) + expected_data = theta * second_term * third_term + assert np.allclose( + expected_data.data, + temperature.equivalent_potential_temperature( + temperature_for_conversions_cube, + relative_humidity_for_conversions_cube, + pressure_for_conversions_cube, + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_equivalent_potenital_temperature_name( + temperature_for_conversions_cube, + relative_humidity_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test name of equivalent potential temperature cube.""" + expected_name = "equivalent_potential_temperature" + assert ( + expected_name + == temperature.equivalent_potential_temperature( + temperature_for_conversions_cube, + relative_humidity_for_conversions_cube, + pressure_for_conversions_cube, + ).name() + ) + + +def test_equivalent_potenital_temperature_units( + temperature_for_conversions_cube, + relative_humidity_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test units of equivalent potential temperature cube.""" + expected_units = cf_units.Unit("K") + assert ( + expected_units + == temperature.equivalent_potential_temperature( + temperature_for_conversions_cube, + relative_humidity_for_conversions_cube, + pressure_for_conversions_cube, + ).units + ) + + +def test_equivalent_potential_temperature_cubelist( + temperature_for_conversions_cube, + relative_humidity_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test calculation of equivalent potential temperature for a CubeList.""" + RH = misc.convert_units(relative_humidity_for_conversions_cube, "1") + theta = temperature.potential_temperature( + temperature_for_conversions_cube, pressure_for_conversions_cube + ) + w = humidity.mixing_ratio_from_relative_humidity( + temperature_for_conversions_cube, pressure_for_conversions_cube, RH + ) + second_term_power = -(w * _atmospheric_constants.RV) / _atmospheric_constants.CPD + second_term = RH.core_data() ** second_term_power.core_data() + third_term_power = ( + _atmospheric_constants.LV + * w + / (_atmospheric_constants.CPD * temperature_for_conversions_cube) + ) + third_term = np.exp(third_term_power.core_data()) + expected_data = theta * second_term * third_term + expected_list = iris.cube.CubeList([expected_data, expected_data]) + temperature_list = iris.cube.CubeList( + [temperature_for_conversions_cube, temperature_for_conversions_cube] + ) + rh_list = iris.cube.CubeList( + [relative_humidity_for_conversions_cube, relative_humidity_for_conversions_cube] + ) + pressure_list = iris.cube.CubeList( + [pressure_for_conversions_cube, pressure_for_conversions_cube] + ) + actual_cubelist = temperature.equivalent_potential_temperature( + temperature_list, rh_list, pressure_list + ) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From a6535a2f89dfa9153c126251df7f395e157332f4 Mon Sep 17 00:00:00 2001 From: daflack Date: Thu, 15 Jan 2026 15:30:39 +0000 Subject: [PATCH 54/91] Adds tests for saturation equivalent potenital temperature --- tests/operators/test_temperature.py | 91 +++++++++++++++++++++++++++++ 1 file changed, 91 insertions(+) diff --git a/tests/operators/test_temperature.py b/tests/operators/test_temperature.py index a7afd0259..15b431326 100644 --- a/tests/operators/test_temperature.py +++ b/tests/operators/test_temperature.py @@ -501,3 +501,94 @@ def test_equivalent_potential_temperature_cubelist( ) for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) + + +def test_saturation_equivalent_potential_temperature( + temperature_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test calculation of saturation equivalent potential temperature.""" + theta = temperature.potential_temperature( + temperature_for_conversions_cube, pressure_for_conversions_cube + ) + ws = humidity.saturation_mixing_ratio( + temperature_for_conversions_cube, pressure_for_conversions_cube + ) + third_term_power = ( + _atmospheric_constants.LV + * ws + / (_atmospheric_constants.CPD * temperature_for_conversions_cube) + ) + third_term = np.exp(third_term_power.core_data()) + expected_data = theta * third_term + assert np.allclose( + expected_data.data, + temperature.saturation_equivalent_potential_temperature( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_saturation_equivalent_potenital_temperature_name( + temperature_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test name of saturation potential temperature cube.""" + expected_name = "saturation_equivalent_potential_temperature" + assert ( + expected_name + == temperature.saturation_equivalent_potential_temperature( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + ).name() + ) + + +def test_saturation_equivalent_potenital_temperature_units( + temperature_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test units of saturation equivalent potential temperature cube.""" + expected_units = cf_units.Unit("K") + assert ( + expected_units + == temperature.saturation_equivalent_potential_temperature( + temperature_for_conversions_cube, + pressure_for_conversions_cube, + ).units + ) + + +def test_saturation_equivalent_potential_temperature_cubelist( + temperature_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test calculation of saturation equivalent potential temperature for a CubeList.""" + theta = temperature.potential_temperature( + temperature_for_conversions_cube, pressure_for_conversions_cube + ) + ws = humidity.saturation_mixing_ratio( + temperature_for_conversions_cube, pressure_for_conversions_cube + ) + third_term_power = ( + _atmospheric_constants.LV + * ws + / (_atmospheric_constants.CPD * temperature_for_conversions_cube) + ) + third_term = np.exp(third_term_power.core_data()) + expected_data = theta * third_term + expected_list = iris.cube.CubeList([expected_data, expected_data]) + temperature_list = iris.cube.CubeList( + [temperature_for_conversions_cube, temperature_for_conversions_cube] + ) + pressure_list = iris.cube.CubeList( + [pressure_for_conversions_cube, pressure_for_conversions_cube] + ) + actual_cubelist = temperature.saturation_equivalent_potential_temperature( + temperature_list, pressure_list + ) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From 9809d4662347eb0677d69d63b2314f28b8eba8c3 Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 19 Jan 2026 10:38:03 +0000 Subject: [PATCH 55/91] Adds tests for wet-bulb potential temperature and fixes its calculation --- src/CSET/operators/temperature.py | 2 +- tests/operators/test_temperature.py | 121 ++++++++++++++++++++++++++++ 2 files changed, 122 insertions(+), 1 deletion(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index f432cff83..f04db28a8 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -182,7 +182,7 @@ def wet_bulb_potential_temperature( iter_maybe(pressure), strict=True, ): - TH_E = equivalent_potential_temperature(T, P, RH) + TH_E = equivalent_potential_temperature(T, RH, P) X = TH_E / T0 X.units = "1" A0 = 7.101574 diff --git a/tests/operators/test_temperature.py b/tests/operators/test_temperature.py index 15b431326..fc4c076f5 100644 --- a/tests/operators/test_temperature.py +++ b/tests/operators/test_temperature.py @@ -592,3 +592,124 @@ def test_saturation_equivalent_potential_temperature_cubelist( ) for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) + + +def test_wet_bulb_potential_temperature( + temperature_for_conversions_cube, + relative_humidity_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test calculation for wet-bulb potential temperature.""" + TH_E = temperature.equivalent_potential_temperature( + temperature_for_conversions_cube, + relative_humidity_for_conversions_cube, + pressure_for_conversions_cube, + ) + X = TH_E / _atmospheric_constants.T0 + X.units = "1" + A0 = 7.101574 + A1 = -20.68208 + A2 = 16.11182 + A3 = 2.574631 + A4 = -5.205688 + B1 = -3.552497 + B2 = 3.781782 + B3 = -0.6899655 + B4 = -0.5929340 + exponent = (A0 + A1 * X + A2 * X**2 + A3 * X**3 + A4 * X**4) / ( + 1.0 + B1 * X + B2 * X**2 + B3 * X**3 + B4 * X**4 + ) + expected_data = TH_E.copy() + expected_data.data[:] = TH_E.core_data() - np.exp(exponent.core_data()) + expected_data.data[expected_data.data - _atmospheric_constants.T0 < -30.0] = np.nan + expected_data.data[expected_data.data - _atmospheric_constants.T0 > 50.0] = np.nan + assert np.allclose( + expected_data.data, + temperature.wet_bulb_potential_temperature( + temperature_for_conversions_cube, + relative_humidity_for_conversions_cube, + pressure_for_conversions_cube, + ).data, + rtol=1e-6, + atol=1e-2, + ) + + +def test_wet_bulb_potential_temperature_name( + temperature_for_conversions_cube, + relative_humidity_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test name of wet-bulb potenital temperature cube.""" + expected_name = "wet_bulb_potential_temperature" + assert ( + expected_name + == temperature.wet_bulb_potential_temperature( + temperature_for_conversions_cube, + relative_humidity_for_conversions_cube, + pressure_for_conversions_cube, + ).name() + ) + + +def test_wet_bulb_potential_temperature_units( + temperature_for_conversions_cube, + relative_humidity_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test units of wet-bulb potenital temperature cube.""" + expected_units = cf_units.Unit("K") + assert ( + expected_units + == temperature.wet_bulb_potential_temperature( + temperature_for_conversions_cube, + relative_humidity_for_conversions_cube, + pressure_for_conversions_cube, + ).units + ) + + +def test_wet_bulb_potential_temperature_cubelist( + temperature_for_conversions_cube, + relative_humidity_for_conversions_cube, + pressure_for_conversions_cube, +): + """Test calculation for wet-bulb potential temperature.""" + TH_E = temperature.equivalent_potential_temperature( + temperature_for_conversions_cube, + relative_humidity_for_conversions_cube, + pressure_for_conversions_cube, + ) + X = TH_E / _atmospheric_constants.T0 + X.units = "1" + A0 = 7.101574 + A1 = -20.68208 + A2 = 16.11182 + A3 = 2.574631 + A4 = -5.205688 + B1 = -3.552497 + B2 = 3.781782 + B3 = -0.6899655 + B4 = -0.5929340 + exponent = (A0 + A1 * X + A2 * X**2 + A3 * X**3 + A4 * X**4) / ( + 1.0 + B1 * X + B2 * X**2 + B3 * X**3 + B4 * X**4 + ) + expected_data = TH_E.copy() + expected_data.data[:] = TH_E.core_data() - np.exp(exponent.core_data()) + expected_data.data[expected_data.data - _atmospheric_constants.T0 < -30.0] = np.nan + expected_data.data[expected_data.data - _atmospheric_constants.T0 > 50.0] = np.nan + expected_list = iris.cube.CubeList([expected_data, expected_data]) + temperature_list = iris.cube.CubeList( + [temperature_for_conversions_cube, temperature_for_conversions_cube] + ) + rh_list = iris.cube.CubeList( + [relative_humidity_for_conversions_cube, relative_humidity_for_conversions_cube] + ) + pressure_list = iris.cube.CubeList( + [pressure_for_conversions_cube, pressure_for_conversions_cube] + ) + actual_cubelist = temperature.wet_bulb_potential_temperature( + temperature_list, rh_list, pressure_list + ) + for cube_a, cube_b in zip(expected_list, actual_cubelist, strict=True): + assert np.allclose(cube_a.data, cube_b.data, rtol=1e-6, atol=1e-2) From e1098fb117567029f6f65e035c6a321eb418d030 Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 19 Jan 2026 11:40:26 +0000 Subject: [PATCH 56/91] Add documentation to _atomspheric_constants.py --- src/CSET/operators/_atmospheric_constants.py | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/src/CSET/operators/_atmospheric_constants.py b/src/CSET/operators/_atmospheric_constants.py index 85390b730..41e55c358 100644 --- a/src/CSET/operators/_atmospheric_constants.py +++ b/src/CSET/operators/_atmospheric_constants.py @@ -14,25 +14,25 @@ """Constants for the atmosphere.""" # Reference pressure. -P0 = 1000.0 # hPa +P0 = 1000.0 # hPa. # Specific gas constant for dry air. -RD = 287.0 +RD = 287.0 # J/kg/K. # Specific gas constant for water vapour. -RV = 461.0 +RV = 461.0 # J/kg/K. # Specific heat capacity for dry air. -CPD = 1005.7 +CPD = 1005.7 # J/kg/K. # Latent heat of vaporization. -LV = 2.501e6 +LV = 2.501e6 # J/kg/K. # Reference vapour pressure. -E0 = 6.1078 # hPa +E0 = 6.1078 # hPa. # Reference temperature. -T0 = 273.15 # K +T0 = 273.15 # K. # Ratio between mixing ratio of dry and moist air. EPSILON = 0.622 From b8918a2cfd0dc474b8121f1299785e84bf3e83bf Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 19 Jan 2026 11:41:26 +0000 Subject: [PATCH 57/91] Update copyright years --- src/CSET/operators/__init__.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/operators/__init__.py b/src/CSET/operators/__init__.py index bb9247c75..cf583c30d 100644 --- a/src/CSET/operators/__init__.py +++ b/src/CSET/operators/__init__.py @@ -1,4 +1,4 @@ -# © Crown copyright, Met Office (2022-2024) and CSET contributors. +# © Crown copyright, Met Office (2022-2026) and CSET contributors. # # Licensed under the Apache License, Version 2.0 (the "License"); # you may not use this file except in compliance with the License. From 1370f60896a58e1a1babde6a3903c098978104b9 Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 19 Jan 2026 13:09:08 +0000 Subject: [PATCH 58/91] Add documentation to exner pressure and fix brackets in beaufort scale conversion --- docs/source/reference/operators.rst | 27 +++++++++++++++++++ src/CSET/operators/pressure.py | 40 ++++++++++++++++++++++++++++- src/CSET/operators/wind.py | 2 +- 3 files changed, 67 insertions(+), 2 deletions(-) diff --git a/docs/source/reference/operators.rst b/docs/source/reference/operators.rst index 259bf9815..db6ab6aa1 100644 --- a/docs/source/reference/operators.rst +++ b/docs/source/reference/operators.rst @@ -95,6 +95,15 @@ CSET.operators.ensembles .. automodule:: CSET.operators.ensembles :members: +Humidity Operators +~~~~~~~~~~~~~~~~~~ + +CSET.operators.humidity +----------------------- + +.. automodule:: CSET.operators.humidity + :members: + Image Processing Operators ~~~~~~~~~~~~~~~~~~~~~~~~~~ @@ -113,6 +122,24 @@ CSET.operators.mesoscale .. automodule:: CSET.operators.mesoscale :members: +Pressure Operators +~~~~~~~~~~~~~~~~~~ + +CSET.operators.pressure +----------------------- + +.. automodule:: CSET.operators.pressure + :members: + +Temperature Operators +~~~~~~~~~~~~~~~~~~~~~ + +CSET.operators.temperature +-------------------------- + +.. automodule:: CSET.operators.temperature + :members: + Wind Operators ~~~~~~~~~~~~~~~~~~~ diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index fe03c3b8c..f3c7d722d 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -63,7 +63,45 @@ def vapour_pressure_from_relative_humidity( def exner_pressure( pressure: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate the exner pressure.""" + r"""Calculate the exner pressure. + + Arguments + --------- + pressure: iris.cube.Cube | iris.cube.CubeList + Cubes of pressure to be converted. + + Returns + ------- + iris.cube.Cube | iris.cube.CubeList + Converted pressure. + + Notes + ----- + The Exner pressure is also referred to as the Exner function. It is a + dimensionless parameter that can be used either as a vertical coordinate or + more frequently as a means to simplifying conversions between different + parameters (e.g. Temperature and Potential Temperature; [Holton13]_). + It is calculated as + + .. math:: \Pi = \left(\frac{P}{P_0}\right)^{\frac{R_d}{c_p}} + + for :math:`\Pi` the Exner Pressure, P the pressure in hPa, :math:`P_0` a reference pressure + of 1000 hPa, :math:`R_d` the specific gas constant of dry air taken as 297 :math:`J kg^{-1} K^{-1}`, + :math:`c_p` the specific heat capacity at constant pressure taken as 1005.7 :math:`J kg^{-1} K^{-1}`. + + A value below one implies the pressure is higher than the reference pressure; + values above one implies the pressure is lower than the reference pressure; a + value of one implies the pressure is equal to the reference pressure. + + References + ---------- + .. [Holton13] Holton, J.R. and Hakim G.J. (2013) An Introduction to Dynamic + Meteorology. 5th Edition, Burlington, MA, Elsevier Academic Press, 532 pp. + + Examples + -------- + >>> Exner_pressure = pressure.exner_pressure(pressure) + """ pi = iris.cube.CubeList([]) for P in iter_maybe(pressure): PI = P.copy() diff --git a/src/CSET/operators/wind.py b/src/CSET/operators/wind.py index 4ea767f63..a9faf2fe6 100644 --- a/src/CSET/operators/wind.py +++ b/src/CSET/operators/wind.py @@ -42,7 +42,7 @@ def convert_to_beaufort_scale( The relationship used to convert the windspeed from m/s to the Beaufort Scale is an empirical relationship (e.g., [Beer96]_): - .. math:: F = (\frac{v}{0.836})^{2/3} + .. math:: F = \left(\frac{v}{0.836}\right)^{2/3} for F the Beaufort Force, and v the windspeed at 10 m in m/s. From fe67c082445574cb70ca1a0b14e9d176fb747b35 Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 19 Jan 2026 13:32:14 +0000 Subject: [PATCH 59/91] Adds documentation for vapour pressure conversion --- src/CSET/operators/pressure.py | 42 ++++++++++++++++++++++++++++++++-- 1 file changed, 40 insertions(+), 2 deletions(-) diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index f3c7d722d..72062f30b 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -25,7 +25,45 @@ def vapour_pressure( temperature: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate the vapour pressure of the atmosphere.""" + r"""Calculate the vapour pressure of the atmosphere. + + Arguments + --------- + temperature: iris.cube.Cube | iris.cube.CubeList + Cubes of temperature to be converted into vapour pressure. + + Returns + ------- + iris.cube.Cube | iris.cube.CubeList + Vapour pressure. + + Notes + ----- + The vapour pressure represents the pressure exerted by water vapour on the + atmosphere. It is related to atmospheric temperature through the + Calusius-Clapeyron equation. There are several different formulations based + on empirical relations that are used to calculate the vapour pressure. Here, + we calculate the vapour pressure based on [Tetens30]_ equation: + + .. math:: e = e_0 exp\left(17.27\frac{T - 273.16}{T - 35.86}\right) + + for e the vapour pressure, :math:`e_0` the saturation vapour pressure at + a reference temperature of 273.15 K with a value of 6.1078 hPa, and T the + temperature. + + If the (dry-bulb) temperature is used the value given by the vapour pressure + will be the saturation vapour pressure. On the other hand, if the dewpoint + temperature is used the vapour pressure will be calculated. + + References + ---------- + .. [Tetens30] Tetens, O. (1930) Über einige meteorologische Begriffe. + Z. Geophys 6: 297-309. + + Examples + -------- + >>> vapour_pressure = pressure.vapour_pressure(temperature) + """ v_pressure = iris.cube.CubeList([]) for T in iter_maybe(temperature): es = T.copy() @@ -73,7 +111,7 @@ def exner_pressure( Returns ------- iris.cube.Cube | iris.cube.CubeList - Converted pressure. + Exner pressure. Notes ----- From 8ef9d26fb871993a7472bd7eb8fcd295c86e2695 Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 19 Jan 2026 13:41:13 +0000 Subject: [PATCH 60/91] Adds documentation for vapour_pressure_from_relative_humidity --- src/CSET/operators/pressure.py | 42 +++++++++++++++++++++++++++++++--- 1 file changed, 39 insertions(+), 3 deletions(-) diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index 72062f30b..ee60c0afb 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -31,11 +31,12 @@ def vapour_pressure( --------- temperature: iris.cube.Cube | iris.cube.CubeList Cubes of temperature to be converted into vapour pressure. + Temperature must be provided in Kelvin. Returns ------- iris.cube.Cube | iris.cube.CubeList - Vapour pressure. + Vapour pressure in hPa. Notes ----- @@ -45,12 +46,15 @@ def vapour_pressure( on empirical relations that are used to calculate the vapour pressure. Here, we calculate the vapour pressure based on [Tetens30]_ equation: - .. math:: e = e_0 exp\left(17.27\frac{T - 273.16}{T - 35.86}\right) + .. math:: e = e_0 * exp\left(17.27\frac{T - 273.16}{T - 35.86}\right) for e the vapour pressure, :math:`e_0` the saturation vapour pressure at a reference temperature of 273.15 K with a value of 6.1078 hPa, and T the temperature. + The temperature is provided in Kelvin, and the vapour pressure is returned + in hPa. + If the (dry-bulb) temperature is used the value given by the vapour pressure will be the saturation vapour pressure. On the other hand, if the dewpoint temperature is used the vapour pressure will be calculated. @@ -82,7 +86,39 @@ def vapour_pressure_from_relative_humidity( temperature: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate the vapour pressure using RH.""" + r"""Calculate the vapour pressure using RH. + + Arguments + --------- + temperature: iris.cube.Cube | iris.cube.CubeList + Cubes of temperature to be converted into saturation vapour pressure. + Temperature must be provided in Kelvin. + relative_humidity: iris.cube.Cube | iris.cube.CubeList + Cubes of relative humidity to be converted into vapour pressure. + Relative humidity will be converted to a decimal by the operator. + + Returns + ------- + iris.cube.Cube | iris.cube.CubeList + Vapour pressure in hPa. + + Notes + ----- + The vapour pressure can be derived from the relative humidity and + temperature based on the following relation + + .. math:: e = RH * e_s + + for e the vapour pressure, :math:`e_s` the saturation vapour pressure, + and RH the relative humidity. + + The relative humidity is converted to a decimal. The saturation vapour + pressure is calculated using `pressure.vapour_pressure`. + + Examples + -------- + >>> vapour_pressure = pressure.vapour_pressure_from_relative_humidity(temperature, relative_humidity) + """ v_pressure = iris.cube.CubeList([]) for T, RH in zip( iter_maybe(temperature), iter_maybe(relative_humidity), strict=True From 84576377cc9e13214893ad65b544691d5a8aaf81 Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 19 Jan 2026 15:10:12 +0000 Subject: [PATCH 61/91] Adds documentation for specific humidity and mixing ratio conversions from each other Also adds placeholder for documentation for other fields. --- src/CSET/operators/humidity.py | 198 ++++++++++++++++++++++++++++-- src/CSET/operators/temperature.py | 184 +++++++++++++++++++++++++-- 2 files changed, 366 insertions(+), 16 deletions(-) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index 9ef248d9e..1d8c9e357 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -25,7 +25,36 @@ def mixing_ratio_from_specific_humidity( specific_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Convert specific humidity to mixing ratio.""" + r"""Convert specific humidity to mixing ratio. + + Arguments + --------- + specific_humidity: iris.cube.Cube | iris.cube.CubeList + Cubes of specific humidity to be converted to mixing ratio. + + Returns + ------- + iris.cube.Cube | iris.cube.CubeList + Converted mixing ratio. + + Notes + ----- + Atmospheric water vapour can be described by multiple quantities. Here, + we convert the specific humidity to the mixing ratio using the following + relation + + .. math:: w = \frac{q}{1 - q} + + for w the mixing ratio and q the specific humidity. + + Larger mixing ratios imply more moisture in the atmosphere. The mixing + ratio will have the same units as the specific humidity (kg/kg). + + + Examples + -------- + >>> w = humidity.mixing_ratio_from_specific_humidity(specific_humidity) + """ w = iris.cube.CubeList([]) for q in iter_maybe(specific_humidity): mr = q.copy() @@ -41,7 +70,34 @@ def mixing_ratio_from_specific_humidity( def specific_humidity_from_mixing_ratio( mixing_ratio: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Convert mixing ratio to specific humidity.""" + r"""Convert mixing ratio to specific humidity. + + Arguments + --------- + mixing_ratio: iris.cube.Cube | iris.Cube.CubeList + Cubes of mixing ratio to be converted to specific humidity. + + Returns + ------- + iris.cube.Cube | iris.cube.CubeList + Converted specific humidity. + + Notes + ----- + Here, we invert the relationship from `humidity.mixing_ratio_from_specific_humidity` + for the following relation + + .. math:: q = \frac{w}{1 + w} + + for q the specific humidity and w the mixing ratio. + + A larger specific humidity implies a more moist atmosphere. The specific + humidity will have the same units as the mixing ratio (kg/kg). + + Examples + -------- + >>> q = humidity.specific_humidity_from_mixing_ratio(mixing_ratio) + """ q = iris.cube.CubeList([]) for w in iter_maybe(mixing_ratio): sh = w.copy() @@ -58,7 +114,28 @@ def saturation_mixing_ratio( temperature: iris.cube.Cube | iris.cube.CubeList, pressure: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate saturation mixing ratio.""" + r"""Calculate saturation mixing ratio. + + Arguments + --------- + s + + Returns + ------- + s + + Notes + ----- + s + + References + ---------- + s + + Examples + -------- + >>> s + """ w = iris.cube.CubeList([]) for T, P in zip(iter_maybe(temperature), iter_maybe(pressure), strict=True): P = convert_units(P, "hPa") @@ -76,7 +153,28 @@ def saturation_specific_humidity( temperature: iris.cube.Cube | iris.cube.CubeList, pressure: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate saturation specific humidity.""" + r"""Calculate saturation specific humidity. + + Arguments + --------- + s + + Returns + ------- + s + + Notes + ----- + s + + References + ---------- + s + + Examples + -------- + >>> s + """ q = iris.cube.CubeList([]) for T, P in zip(iter_maybe(temperature), iter_maybe(pressure), strict=True): P = convert_units(P, "hPa") @@ -95,7 +193,28 @@ def mixing_ratio_from_relative_humidity( pressure: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate the mixing ratio from RH.""" + r"""Calculate the mixing ratio from RH. + + Arguments + --------- + s + + Returns + ------- + s + + Notes + ----- + s + + References + ---------- + s + + Examples + -------- + >>> s + """ w = iris.cube.CubeList([]) for T, P, RH in zip( iter_maybe(temperature), @@ -119,7 +238,28 @@ def specific_humidity_from_relative_humidity( pressure: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate the mixing ratio from RH.""" + r"""Calculate the mixing ratio from RH. + + Arguments + --------- + s + + Returns + ------- + s + + Notes + ----- + s + + References + ---------- + s + + Examples + -------- + >>> s + """ q = iris.cube.CubeList([]) for T, P, RH in zip( iter_maybe(temperature), @@ -143,7 +283,28 @@ def relative_humidity_from_mixing_ratio( temperature: iris.cube.Cube | iris.cube.CubeList, pressure: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Convert mixing ratio to relative humidity.""" + r"""Convert mixing ratio to relative humidity. + + Arguments + --------- + s + + Returns + ------- + s + + Notes + ----- + s + + References + ---------- + s + + Examples + -------- + >>> s + """ RH = iris.cube.CubeList([]) for W, T, P in zip( iter_maybe(mixing_ratio), @@ -166,7 +327,28 @@ def relative_humidity_from_specific_humidity( temperature: iris.cube.Cube | iris.cube.CubeList, pressure: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Convert specific humidity to relative humidity.""" + r"""Convert specific humidity to relative humidity. + + Arguments + --------- + s + + Returns + ------- + s + + Notes + ----- + s + + References + ---------- + s + + Examples + -------- + >>> s + """ RH = iris.cube.CubeList([]) for Q, T, P in zip( iter_maybe(specific_humidity), diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index f04db28a8..63e74378d 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -34,7 +34,28 @@ def dewpoint_temperature( temperature: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate the dewpoint temperature following Bolton (1980).""" + r"""Calculate the dewpoint temperature following Bolton (1980). + + Arguments + --------- + s + + Returns + ------- + s + + Notes + ----- + s + + References + ---------- + s + + Examples + -------- + >>> s + """ Td = iris.cube.CubeList([]) for T, RH in zip( iter_maybe(temperature), iter_maybe(relative_humidity), strict=True @@ -60,7 +81,28 @@ def virtual_temperature( temperature: iris.cube.Cube | iris.cube.CubeList, mixing_ratio: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate the virtual temperature.""" + r"""Calculate the virtual temperature. + + Arguments + --------- + s + + Returns + ------- + s + + Notes + ----- + s + + References + ---------- + s + + Examples + -------- + >>> s + """ Tv = iris.cube.CubeList([]) for T, W in zip(iter_maybe(temperature), iter_maybe(mixing_ratio), strict=True): virT = T * ((W + EPSILON) / (EPSILON * (1 + W))) @@ -76,7 +118,28 @@ def wet_bulb_temperature( temperature: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate the wet-bulb temperature following Stull (2011).""" + r"""Calculate the wet-bulb temperature following Stull (2011). + + Arguments + --------- + s + + Returns + ------- + s + + Notes + ----- + s + + References + ---------- + s + + Examples + -------- + >>> s + """ Tw = iris.cube.CubeList([]) for T, RH in zip( iter_maybe(temperature), iter_maybe(relative_humidity), strict=True @@ -105,7 +168,28 @@ def potential_temperature( temperature: iris.cube.Cube | iris.cube.CubeList, pressure: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate the potenital temperature.""" + r"""Calculate the potenital temperature. + + Arguments + --------- + s + + Returns + ------- + s + + Notes + ----- + s + + References + ---------- + s + + Examples + -------- + >>> s + """ theta = iris.cube.CubeList([]) for T, P in zip(iter_maybe(temperature), iter_maybe(pressure), strict=True): TH = T / exner_pressure(P) @@ -122,7 +206,28 @@ def virtual_potential_temperature( mixing_ratio: iris.cube.Cube | iris.cube.CubeList, pressure: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate the virtual potential temperature.""" + r"""Calculate the virtual potential temperature. + + Arguments + --------- + s + + Returns + ------- + s + + Notes + ----- + s + + References + ---------- + s + + Examples + -------- + >>> s + """ theta_v = iris.cube.CubeList([]) for T, W, P in zip( iter_maybe(temperature), @@ -144,7 +249,28 @@ def equivalent_potential_temperature( relative_humidity: iris.cube.Cube | iris.cube.CubeList, pressure: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate the equivalent potenital temperature.""" + r"""Calculate the equivalent potenital temperature. + + Arguments + --------- + s + + Returns + ------- + s + + Notes + ----- + s + + References + ---------- + s + + Examples + -------- + >>> s + """ theta_e = iris.cube.CubeList([]) for T, RH, P in zip( iter_maybe(temperature), @@ -174,7 +300,28 @@ def wet_bulb_potential_temperature( relative_humidity: iris.cube.Cube | iris.cube.CubeList, pressure: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate wet-bulb potential temperature after Davies-Jones (2008).""" + r"""Calculate wet-bulb potential temperature after Davies-Jones (2008). + + Arguments + --------- + s + + Returns + ------- + s + + Notes + ----- + s + + References + ---------- + s + + Examples + -------- + >>> s + """ theta_w = iris.cube.CubeList([]) for T, RH, P in zip( iter_maybe(temperature), @@ -213,7 +360,28 @@ def saturation_equivalent_potential_temperature( temperature: iris.cube.Cube | iris.cube.CubeList, pressure: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - """Calculate the saturation equivalent potenital temperature.""" + r"""Calculate the saturation equivalent potenital temperature. + + Arguments + --------- + s + + Returns + ------- + s + + Notes + ----- + s + + References + ---------- + s + + Examples + -------- + >>> s + """ theta_s = iris.cube.CubeList([]) for T, P in zip( iter_maybe(temperature), From 3f43d06a6899b7653ba4f90ea47f77b0dcee2b14 Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 19 Jan 2026 15:27:28 +0000 Subject: [PATCH 62/91] Adds documentation for specific humidity and mixing ratio from RH --- src/CSET/operators/humidity.py | 100 ++++++++++++++++++++++++--------- 1 file changed, 72 insertions(+), 28 deletions(-) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index 1d8c9e357..7297e407f 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -118,23 +118,33 @@ def saturation_mixing_ratio( Arguments --------- - s + temperature: iris.cube.Cube | iris.cube.CubeList + Cubes of temperature in Kelvin. + pressure: iris.cube.Cube | iris.cube.CubeList + Cubes of pressure. Returns ------- - s + iris.cube.Cube | iris.cube.CubeList + Saturation mixing ratio in kg/kg. Notes ----- - s + The saturation mixing ratio is required to help calculate the relative + humidity and other diagnostics with respect to the mixing ratio. It can + be calculated from - References - ---------- - s + .. math:: w = \epsilon \frac{e}{P - e} + + for w the mixing ratio, :math:`\epsilon` the ratio between the mixing ratio + of dry and moist air equating to 0.622, P the pressure and e the vapour + pressure. To ensure that the saturation mixing ratio (:math:`w_s`) is + calculated the vapour pressure calculated should be with the (dry-bulb) + temperature to ensure it is the saturated vapour pressure. Examples -------- - >>> s + >>> w_s = humidity.saturation_mixing_ratio(temperature, pressure) """ w = iris.cube.CubeList([]) for T, P in zip(iter_maybe(temperature), iter_maybe(pressure), strict=True): @@ -157,23 +167,33 @@ def saturation_specific_humidity( Arguments --------- - s + temperature: iris.cube.Cube | iris.cube.CubeList + Cubes of temperature in Kelvin. + pressure: iris.cube.Cube | iris.cube.CubeList + Cubes of pressure. Returns ------- - s + iris.cube.Cube | iris.cube.CubeList + Saturation specific humidity in kg/kg. Notes ----- - s + The saturation specific humidity is required to help calculate the relative + humidity and other diagnostics with respect to the mixing ratio. It can + be calculated from - References - ---------- - s + .. math:: q = \epsilon \frac{e}{P} + + for q the specific humidity, :math:`\epsilon` the ratio between the mixing ratio + of dry and moist air equating to 0.622, P the pressure and e the vapour + pressure. To ensure that the saturation specific humidity (:math:`q_{sat}`) is + calculated the vapour pressure calculated should be with the (dry-bulb) + temperature to ensure it is the saturated vapour pressure. Examples -------- - >>> s + >>> q_sat = humidity.saturation_specific_humidity(temperature, pressure) """ q = iris.cube.CubeList([]) for T, P in zip(iter_maybe(temperature), iter_maybe(pressure), strict=True): @@ -197,23 +217,35 @@ def mixing_ratio_from_relative_humidity( Arguments --------- - s + temperature: iris.cube.Cube | iris.cube.CubeList + Cubes of temperature in Kelvin. + pressure: iris.cube.Cube | iris.cube.CubeList + Cubes of pressure. + relative_humidity: iris.cube.Cube | iris.cube.CubeList + Cubes of relative humidity. Returns ------- - s + iris.cube.Cube | iris.cube.CubeList + Calculated mixing ratio from relative humidity in kg/kg. Notes ----- - s + The mixing ratio can be calculated from temperature, pressure, and + relative humidity using the following relation - References - ---------- - s + .. math:: w = RH * w_s + + for w the mixing ratio, :math:`w_s` the saturation mixing ratio, and + RH the relative humidity. + + The operator usese `humidity.saturation_mixing_ratio` to calculate the + saturation mixing ratio from the temperature and pressure. The relative + humidity is converted into a decimal before the multiplication occurs. Examples -------- - >>> s + >>> w = humidity.mixing_ratio_from_relative_humidity(T, P, RH) """ w = iris.cube.CubeList([]) for T, P, RH in zip( @@ -242,23 +274,35 @@ def specific_humidity_from_relative_humidity( Arguments --------- - s + temperature: iris.cube.Cube | iris.cube.CubeList + Cubes of temperature in Kelvin. + pressure: iris.cube.Cube | iris.cube.CubeList + Cubes of pressure. + relative_humidity: iris.cube.Cube | iris.cube.CubeList + Cubes of relative humidity. Returns ------- - s + iris.cube.Cube | iris.cube.CubeList + Calculated specific humidity from relative humidity in kg/kg. Notes ----- - s + The specific humidity can be calculated from temperature, pressure, and + relative humidity using the following relation - References - ---------- - s + .. math:: q = RH * q_{sat} + + for q the specific humidity, :math:`q_{sat}` the saturation specific + humidity, and RH the relative humidity. + + The operator usese `humidity.saturation_specific_humidity` to calculate the + saturation specific humidity from the temperature and pressure. The relative + humidity is converted into a decimal before the multiplication occurs. Examples -------- - >>> s + >>> q = humidity.specific_humidity_from_relative_humidity(T, P, RH) """ q = iris.cube.CubeList([]) for T, P, RH in zip( From 9b4a4cf4693680d7d76c21d40eaa61d6477232f3 Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 19 Jan 2026 15:45:02 +0000 Subject: [PATCH 63/91] Adds documentation for relative humidity conversions from different humidities --- src/CSET/operators/humidity.py | 62 ++++++++++++++++++++++++++-------- 1 file changed, 48 insertions(+), 14 deletions(-) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index 7297e407f..1633ed9ae 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -142,6 +142,8 @@ def saturation_mixing_ratio( calculated the vapour pressure calculated should be with the (dry-bulb) temperature to ensure it is the saturated vapour pressure. + All cubes need to be on the same grid. + Examples -------- >>> w_s = humidity.saturation_mixing_ratio(temperature, pressure) @@ -191,6 +193,8 @@ def saturation_specific_humidity( calculated the vapour pressure calculated should be with the (dry-bulb) temperature to ensure it is the saturated vapour pressure. + All cubes need to be on the same grid. + Examples -------- >>> q_sat = humidity.saturation_specific_humidity(temperature, pressure) @@ -243,6 +247,8 @@ def mixing_ratio_from_relative_humidity( saturation mixing ratio from the temperature and pressure. The relative humidity is converted into a decimal before the multiplication occurs. + All cubes need to be on the same grid. + Examples -------- >>> w = humidity.mixing_ratio_from_relative_humidity(T, P, RH) @@ -300,6 +306,8 @@ def specific_humidity_from_relative_humidity( saturation specific humidity from the temperature and pressure. The relative humidity is converted into a decimal before the multiplication occurs. + All cubes need to be on the same grid. + Examples -------- >>> q = humidity.specific_humidity_from_relative_humidity(T, P, RH) @@ -331,23 +339,36 @@ def relative_humidity_from_mixing_ratio( Arguments --------- - s + mixing_ratio: iris.cube.Cube | iris.cube.CubeList + Cubes of mixing ratio. + temperature: iris.cube.Cube | iris.cube.CuebList + Cubes of temperature in Kelvin. + pressure: iris.cube.Cube | iris.cube.CubeList + Cubes of pressure. Returns ------- - s + iris.cube.Cube | iris.cube.CubeList + Relative humidity calculated from mixing ratio. Notes ----- - s + The relative humidity can be calculated from the mixing ratio following + + .. math:: RH = \frac{w}{w_s} + + for RH the relative humidity, w the mixing ratio, and :math:`w_s` the + saturation mixing ratio. The saturation mixing ratio is calculated using + `humidity.saturation_mixing_ratio`. - References - ---------- - s + The RH varies predominatly between zero (completely dry) and one (saturated). + Values larger than one are possible and imply supersaturation. + + All cubes must be on the same grid. Examples -------- - >>> s + >>> RH = humidity.relative_humidity_from_mixing_ratio(w, T, P) """ RH = iris.cube.CubeList([]) for W, T, P in zip( @@ -375,23 +396,36 @@ def relative_humidity_from_specific_humidity( Arguments --------- - s + specific_humidity: iris.cube.Cube | iris.cube.CubeList + Cubes of specific humidity. + temperature: iris.cube.Cube | iris.cube.CuebList + Cubes of temperature in Kelvin. + pressure: iris.cube.Cube | iris.cube.CubeList + Cubes of pressure. Returns ------- - s + iris.cube.Cube | iris.cube.CubeList + Relative humidity calculated from specific humidity. Notes ----- - s + The relative humidity can be calculated from the specific humidity following + + .. math:: RH = \frac{q}{q_{sat}} + + for RH the relative humidity, q the specific humidity, and :math:`q_{sat}` the + saturation specific humidity. The saturation specific humidity is calculated using + `humidity.saturation_specific_humidity`. + + The RH varies predominatly between zero (completely dry) and one (saturated). + Values larger than one are possible and imply supersaturation. - References - ---------- - s + All cubes must be on the same grid. Examples -------- - >>> s + >>> RH = humidity.relative_humidity_from_specific_humidity(q, T, P) """ RH = iris.cube.CubeList([]) for Q, T, P in zip( From 1523c50a1bcb84add486933d7e97e87f1821cc9d Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 19 Jan 2026 16:14:03 +0000 Subject: [PATCH 64/91] Adds documentation for dewpoint temperature --- src/CSET/operators/pressure.py | 6 ++--- src/CSET/operators/temperature.py | 37 ++++++++++++++++++++++++++----- 2 files changed, 34 insertions(+), 9 deletions(-) diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index ee60c0afb..6fa057775 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -61,7 +61,7 @@ def vapour_pressure( References ---------- - .. [Tetens30] Tetens, O. (1930) Über einige meteorologische Begriffe. + .. [Tetens30] Tetens, O. (1930) "Über einige meteorologische Begriffe". Z. Geophys 6: 297-309. Examples @@ -169,8 +169,8 @@ def exner_pressure( References ---------- - .. [Holton13] Holton, J.R. and Hakim G.J. (2013) An Introduction to Dynamic - Meteorology. 5th Edition, Burlington, MA, Elsevier Academic Press, 532 pp. + .. [Holton13] Holton, J.R. and Hakim G.J. (2013) "An Introduction to Dynamic + Meteorology." 5th Edition, Burlington, MA, Elsevier Academic Press, 532 pp. Examples -------- diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 63e74378d..9d5c647e8 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -34,27 +34,52 @@ def dewpoint_temperature( temperature: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - r"""Calculate the dewpoint temperature following Bolton (1980). + r"""Calculate the dewpoint temperature. Arguments --------- - s + temperature: iris.cube.Cube | iris.cube.CubeList + Cubes of temperature in Kelvin. + relative_humidity: iris.cube.Cube | iris.cube.CubeList + Cubes of relative humidity. Returns ------- - s + iris.cube.Cube | iris.cueb.CubeList + Calculated dewpoint temperature in Kelvin. Notes ----- - s + The dewpoint temperature provides the temperature at which the air must be + cooled for dew to form (i.e. the water vapour condenses to liquid water). + It is calculated based upon [Bolton80] + + .. math:: T_d = \frac{243.5 * ln(e) - 440.8}{19.48 - ln(e)} + + for :math:`T_d` the dewpoint temperature and e the vapour pressure. Here the + vapour pressure is calculated from relative humidity using + `pressure.vapour_pressure_from_relative_humidity`. + + The dewpoint temperature is presented when -35.0 C < T < 35.0 C as this is + when the calculation is the most accurate [Bolton80] this roughly equates to + exclusion of dewpoints on pressures of 400 hPa and above (e.g. [Flack24]). + + When :math:`T_d` is equivalent to T the relative humidity will be 100 %. + + All cubes must be on the same grid. References ---------- - s + .. [Bolton80] Bolton. D. (1980) "The computation of equivalent potential + temperature". Monthly Weather Review, vol. 108, 1046-1053, + doi: 10.1175/1520-0493(1980)108<1046:TCOEPT>2.0.CO;2 + .. [Flack24] Flack, D.L.A. (2024) "Stratification of the vertical spread-skill + relation by radiosonde drift in a convective-scale ensemble." + Atmospheric Science Letters, vol. 25, e1194, doi: 10.1002/asl.1194 Examples -------- - >>> s + >>> Td = temperature.dewpoint_temperature(T, RH) """ Td = iris.cube.CubeList([]) for T, RH in zip( From 743f30be4df1df629b6d72990c15bd2092437d04 Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 19 Jan 2026 16:16:38 +0000 Subject: [PATCH 65/91] Fix documentaiton --- src/CSET/operators/temperature.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 9d5c647e8..bd1d44917 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -52,7 +52,7 @@ def dewpoint_temperature( ----- The dewpoint temperature provides the temperature at which the air must be cooled for dew to form (i.e. the water vapour condenses to liquid water). - It is calculated based upon [Bolton80] + It is calculated based upon [Bolton80]_ .. math:: T_d = \frac{243.5 * ln(e) - 440.8}{19.48 - ln(e)} @@ -61,8 +61,8 @@ def dewpoint_temperature( `pressure.vapour_pressure_from_relative_humidity`. The dewpoint temperature is presented when -35.0 C < T < 35.0 C as this is - when the calculation is the most accurate [Bolton80] this roughly equates to - exclusion of dewpoints on pressures of 400 hPa and above (e.g. [Flack24]). + when the calculation is the most accurate [Bolton80]_ this roughly equates to + exclusion of dewpoints on pressures of 400 hPa and above (e.g. [Flack24]_). When :math:`T_d` is equivalent to T the relative humidity will be 100 %. From 469eda7ad3cf99597dbe4b81366e13a2ff44ab04 Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 19 Jan 2026 16:21:51 +0000 Subject: [PATCH 66/91] Adds documentation for virtual temperature --- src/CSET/operators/temperature.py | 24 +++++++++++++++++------- 1 file changed, 17 insertions(+), 7 deletions(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index bd1d44917..25f8969a5 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -110,23 +110,33 @@ def virtual_temperature( Arguments --------- - s + temperature: iris.cube.Cube | iris.cube.CubeList + Cubes of temperature. + mixing_ratio: iris.cube.Cube | iris.cube.CubeList + Cubes of mixing ratio. Returns ------- - s + iris.cube.Cube | iris.cube.CubeList + Calculated virtual temperature. Notes ----- - s + The virtual temperature is that required for a dry parcel to have the + same density as a moist parcel at the same pressure. It is calculated + as - References - ---------- - s + .. math:: T_v = T * \frac{w + \epsilon}{\epsilon (1 + w)} + + for :math:`T_v` the virtual temperature, T the temperature, w the mixing + ratio, and :math:`\epsilon` the ratio between dry and moist air equating + to 0.622. + + All cubes must be on the same grid. Examples -------- - >>> s + >>> Tv = temperature.virtual_temperature(T, w) """ Tv = iris.cube.CubeList([]) for T, W in zip(iter_maybe(temperature), iter_maybe(mixing_ratio), strict=True): From e772aac92beea5294910da5cfbb3e35d18644ff2 Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 19 Jan 2026 16:42:28 +0000 Subject: [PATCH 67/91] Adds documentation for wet-bulb temperature --- src/CSET/operators/temperature.py | 40 +++++++++++++++++++++++------ tests/operators/test_temperature.py | 8 ++++++ 2 files changed, 40 insertions(+), 8 deletions(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 25f8969a5..84c0b5834 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -111,14 +111,14 @@ def virtual_temperature( Arguments --------- temperature: iris.cube.Cube | iris.cube.CubeList - Cubes of temperature. + Cubes of temperature in Kelvin. mixing_ratio: iris.cube.Cube | iris.cube.CubeList Cubes of mixing ratio. Returns ------- iris.cube.Cube | iris.cube.CubeList - Calculated virtual temperature. + Calculated virtual temperature in Kelvin. Notes ----- @@ -153,27 +153,47 @@ def wet_bulb_temperature( temperature: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - r"""Calculate the wet-bulb temperature following Stull (2011). + r"""Calculate the wet-bulb temperature. Arguments --------- - s + temperature: iris.cube.Cube | iris.cube.CubeList + Cubes of temperature. + relative_humidity: iris.cube.Cube | iris.cube.CubeList + Cubes of relative humidity. Returns ------- - s + iris.cube.Cube | iris.cube.CubeList + Calculated wet-bulb temperature. Notes ----- - s + The wet-bulb temperature can be calculated from temperature in Celsius + and relative humidity in percent following [Stull11]_ + + .. math:: T_w = T * arctan\left(0.151977*(RH + 8.313659)^{0.5}\right) + arctan(T + RH) - arctan(RH - 1.676331) + 0.00391838*RH^{\frac{3}{2}}*arctan(0.023101*RH) - 4.686035 + + for :math:`T_w` the wet-bulb temperature, T the temperature, and RH the relative + humidity. + + The equation is valid for 5% < RH < 99% and -20 C < T < 50 C, and results are + only presented for these values. + + The temperature and relative humidity unit conversions are applied in the + operator. + + All cubes should be on the same grid. References ---------- - s + .. [Stull11] Stull, R. (2011) "Wet-Bulb Temperature from Relative Humidity + and air temperature." Journal of Applied Meteorology and Climatology, vol. 50, + 2267-2269, doi: 10.1175/JAMC-D-11-0143.1 Examples -------- - >>> s + >>> Tw = temperature.wet_bulb_temperature(T, RH) """ Tw = iris.cube.CubeList([]) for T, RH in zip( @@ -192,6 +212,10 @@ def wet_bulb_temperature( ) wetT.rename("wet_bulb_temperature") wetT = convert_units(wetT, "K") + wetT.data[T.core_data() < -20.0] = np.nan + wetT.data[T.core_data() > 50.0] = np.nan + wetT.data[RH.core_data() < 5.0] = np.nan + wetT.data[RH.core_data() > 99.0] = np.nan Tw.append(wetT) if len(Tw) == 1: return Tw[0] diff --git a/tests/operators/test_temperature.py b/tests/operators/test_temperature.py index fc4c076f5..1af6adcae 100644 --- a/tests/operators/test_temperature.py +++ b/tests/operators/test_temperature.py @@ -189,6 +189,10 @@ def test_wet_bulb_temperature( - 4.686035 ) TW = misc.convert_units(TW, "K") + TW.data[T_units.data < -20.0] = np.nan + TW.data[T_units.data > 50.0] = np.nan + TW.data[RH_units.data < 5.0] = np.nan + TW.data[RH_units.data > 99.0] = np.nan assert np.allclose( TW.data, temperature.wet_bulb_temperature( @@ -241,6 +245,10 @@ def test_wet_bulb_temperature_cubelist( - 4.686035 ) TW = misc.convert_units(TW, "K") + TW.data[T_units.data < -20.0] = np.nan + TW.data[T_units.data > 50.0] = np.nan + TW.data[RH_units.data < 5.0] = np.nan + TW.data[RH_units.data > 99.0] = np.nan expected_list = iris.cube.CubeList([TW, TW]) temperature_list = iris.cube.CubeList( [temperature_for_conversions_cube, temperature_for_conversions_cube] From fb0f655a564b1afcc705a5c6fda05d29a79fa78a Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 19 Jan 2026 16:44:25 +0000 Subject: [PATCH 68/91] Adds further documentation to wet-bulb temperature --- src/CSET/operators/temperature.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 84c0b5834..6710cc4df 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -169,7 +169,8 @@ def wet_bulb_temperature( Notes ----- - The wet-bulb temperature can be calculated from temperature in Celsius + The wet-bulb temperature is the temperature the air reaches when all the + water has been evaporated out of it. It can be calculated from temperature in Celsius and relative humidity in percent following [Stull11]_ .. math:: T_w = T * arctan\left(0.151977*(RH + 8.313659)^{0.5}\right) + arctan(T + RH) - arctan(RH - 1.676331) + 0.00391838*RH^{\frac{3}{2}}*arctan(0.023101*RH) - 4.686035 From 47524919ba689c041c84073dccefb7b879996a44 Mon Sep 17 00:00:00 2001 From: daflack Date: Mon, 19 Jan 2026 16:51:01 +0000 Subject: [PATCH 69/91] Adds potential temperature documentation --- src/CSET/operators/temperature.py | 25 +++++++++++++++++-------- 1 file changed, 17 insertions(+), 8 deletions(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 6710cc4df..c65d28ae9 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -165,7 +165,7 @@ def wet_bulb_temperature( Returns ------- iris.cube.Cube | iris.cube.CubeList - Calculated wet-bulb temperature. + Calculated wet-bulb temperature in Kelvin. Notes ----- @@ -232,23 +232,32 @@ def potential_temperature( Arguments --------- - s + temperature: iris.cube.Cube | iris.cube.CubeList + Cubes of temperature. + pressure: iris.cube.Cube | iris.cube.CuebList + Cubes of pressure. Returns ------- - s + iris.cube.Cube | iris.cube.CubeList + Calculated potential temperature in Kelvin. Notes ----- - s + The potential temperature is the temperature at which an air parcel would reach + if it was moved adiabatically to a reference pressure. Here we use 1000 hPa + as the reference pressure. It is calculated from - References - ---------- - s + .. math:: \theta = \frac{T}{\Pi} + + for :math:`\theta` the potential temperature, T the temperature, and :math:`\Pi` + the exner pressure. The exner pressure is calculated using `pressure.exner_pressure`. + + All cubes must be on the same grid. Examples -------- - >>> s + >>> Theta = temperature.potential_temperature(T, P) """ theta = iris.cube.CubeList([]) for T, P in zip(iter_maybe(temperature), iter_maybe(pressure), strict=True): From d2499e107de185dfe944f64487abcdb5dd1dec74 Mon Sep 17 00:00:00 2001 From: daflack Date: Tue, 20 Jan 2026 08:55:56 +0000 Subject: [PATCH 70/91] Adds virtual potential temperature documentation --- src/CSET/operators/temperature.py | 27 ++++++++++++++++++++------- 1 file changed, 20 insertions(+), 7 deletions(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index c65d28ae9..6235d94c2 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -279,23 +279,36 @@ def virtual_potential_temperature( Arguments --------- - s + temperature: iris.cube.Cube | iris.cube.CubeList + Cubes of temperature. + mixing_ratio: iris.cube.Cube | iris.cube.CubeList + Cubes of mixing ratio. + pressure: iris.cube.Cube | iris.cube.CubeList + Cubes of pressure. Returns ------- - s + iris.cube.Cube | iris.cube.CubeList + Calculated virtual potential temperature in Kelvin. Notes ----- - s + The virtual potenital temperature is mechanistically equivalent to the + potential temperature, except rather than using the (dry-bulb) temperature + the virtual temperature is adiabatically moved to a reference pressure + (1000 hPa). It is calculated as - References - ---------- - s + .. math:: \theta_v = \frac{T_v}{\Pi} + + for :math:`\theta_v` the virtual potenital temperature, :math:`T_v` the + virtual temperature, and :math:`\Pi` the exner pressure. The exner pressure + is calculated using `pressure.exner_pressure`. + + All cubes must be on the same grid. Examples -------- - >>> s + >>> Theta_v = temperature.virtual_potenital_temperature(T, W, P) """ theta_v = iris.cube.CubeList([]) for T, W, P in zip( From 52b6e34fd0de44a466da414ce2f12f31cfe6d14a Mon Sep 17 00:00:00 2001 From: daflack Date: Tue, 20 Jan 2026 13:35:31 +0000 Subject: [PATCH 71/91] Adds documentation for wet-bulb potential temperature --- src/CSET/operators/temperature.py | 39 ++++++++++++++++++++++++++----- 1 file changed, 33 insertions(+), 6 deletions(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 6235d94c2..d163c0e5d 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -382,27 +382,54 @@ def wet_bulb_potential_temperature( relative_humidity: iris.cube.Cube | iris.cube.CubeList, pressure: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - r"""Calculate wet-bulb potential temperature after Davies-Jones (2008). + r"""Calculate wet-bulb potential temperature. Arguments --------- - s + temperature: iris.cube.Cube | iris.cube.CubeList + Cubes of temperature. + relative_humidity: iris.cube.Cube | iris.cube.CubeList + Cubes of relative humidity. + pressure: iris.cube.Cube | iris.cube.CubeList + Cubes of pressure. Returns ------- - s + iris.cube.Cube | iris.cube.CubeList + Calculated wet-bulb potential temperature in Kelvin. Notes ----- - s + The wet-bulb potential temperature represents the temperature an air parcel + would have if coold adiabatically until saturation and then taken down to + a reference pressure (1000 hPa) whilst conserving the moisture (i.e. along a + psuedoadiabat). It can be calculated either through a series of iterations, + or through empirical relations. Here, we use the [DaviesJones08]_ formulation: + + .. math:: \theta_w = \theta_e - exp\left(\frac{a_0 + a_1 X + a_2 X^2 + a_3 X^3 + a_4 X^4}{1.0 + b_1 X + b_2 X^2 + b_3 X^3 + b_4 X^4} \right) + + for :math:`\theta_w` the wet-bulb potential temperature, :math:`\theta_e` the + equivalent potenital temperature, X = :math:`\theta_e` / 273.15 K, and a + series of constants :math:`a_0` = 7.101574, :math:`a_1` = -20.68208, + :math:`a_2` = 16.11182, :math:`a_3` = 2.574631, :math:`a_4` = -5.205688, + :math:`b_1` = -3.552497, :math:`b_2` = 3.781782, :math:`b_3` = -0.6899655, and + :math:`b_4` = -0.5929340. + + The wet-bulb potential temperature calculated for this method is valid for + the range -30 C < :math:`\theta_w` < 50 C and is thus set to NaN outside of + this range. + + All cubes must be on the same grid. References ---------- - s + .. [DaviesJones08] Davies-Jones, R. (2008) "An Efficient and Accurate + Method for Computing the Wet-Bulb Temperature along Pseudoadiabats" + Monthly Weather Review, vol. 136, 2764-2785, doi: 10.1175/2007MWR2224.1 Examples -------- - >>> s + >>> Theta_w = temperature.wet_bulb_potential_temperature(T, RH, P) """ theta_w = iris.cube.CubeList([]) for T, RH, P in zip( From 3dccf83dbacc8e274c3247c3059cb601a4abf51c Mon Sep 17 00:00:00 2001 From: daflack Date: Tue, 20 Jan 2026 14:04:31 +0000 Subject: [PATCH 72/91] Adds equivalent potenital temperature documentation --- src/CSET/operators/temperature.py | 40 +++++++++++++++++++++++++++---- 1 file changed, 35 insertions(+), 5 deletions(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index d163c0e5d..dacc04346 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -335,23 +335,53 @@ def equivalent_potential_temperature( Arguments --------- - s + temperature: iris.cube.Cube | iris.cube.CubeList + Cubes of temperature. + relative_humidity: iris.cube.Cube | iris.cube.CubeList + Cubes of relative humidity. + pressure: iris.cube.Cube | iris.cube.CubeList + Cubes of pressure. Returns ------- - s + iris.cube.Cube | iris.cube.CubeList + Calculated equivalent potential temperature. Notes ----- - s + The equivalent potential temperature is the temperature an air parcel + would have if it was raised until all of the water vapour in it was + condensed out and then brought adiabatically to the reference pressure, + here taken as 1000 hPa. It is calculated as + + .. math:: \theta_e = \theta * RH^{- \frac{w R_v}{c_p}} * exp\left(\frac{L_v w}{c_p T} \right) + + for :math:`\theta_e` the equivalent potenital temperature, :math:`\theta` the + potential temperature, RH the relative humidity, w the mixing ratio, + :math:`R_v` the specific gas constant of water vapour (461 + :math:`J kg^{-1} K^{-1}`), :math:`c_p` the specific heat capacity of dry air + (1005.7 :math:`J kg^{-1} K^{-1}`), :math:`L_v` the latent heat of vapourization + (2.5 x :math:`10^6 J kg^{-1} K^{-1}`), and T the temperature. + + This is a simplification of [Paluch79]_ following [Bolton80b]_ and [Emanuel94]_, + whilst still holding reasonable accuracy. + + All cubes must be on the same grid. References ---------- - s + .. [Bolton80b] Bolton. D. (1980) "The computation of equivalent potential + temperature". Monthly Weather Review, vol. 108, 1046-1053, + doi: 10.1175/1520-0493(1980)108<1046:TCOEPT>2.0.CO;2 + .. [Emanuel94] Emanuel, K.A. (1994) "Atmospheric Convection" Oxford University + Press, 580 pp. + .. [Paluch79] Paluch, I.R., (1979) "The entrainment mechanism in Colorado + Cumuli" Journal of the Atmospheric Sciences, vol. 36, 2467-2478, + doi: 10.1175/1520-0469(1979)036<2467:TEMICC>2.0.CO;2 Examples -------- - >>> s + >>> Theta_e = temperature.equivalent_potential_temperature(T, RH, P) """ theta_e = iris.cube.CubeList([]) for T, RH, P in zip( From 472b2890db5e1c1d2a6a75e348c9db80dde899d5 Mon Sep 17 00:00:00 2001 From: daflack Date: Tue, 20 Jan 2026 14:07:25 +0000 Subject: [PATCH 73/91] Adds clarification to theta_e documentation --- src/CSET/operators/temperature.py | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index dacc04346..10f1a2142 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -363,6 +363,10 @@ def equivalent_potential_temperature( (1005.7 :math:`J kg^{-1} K^{-1}`), :math:`L_v` the latent heat of vapourization (2.5 x :math:`10^6 J kg^{-1} K^{-1}`), and T the temperature. + In this operator the mixing ratio is calculated from + `humidity.mixing_ratio_from_relative_humidity` and the potential temperature + is calculated from `temperature.potential_temperature`. + This is a simplification of [Paluch79]_ following [Bolton80b]_ and [Emanuel94]_, whilst still holding reasonable accuracy. From 24d6e6a0b9c2e5f7f0da62fc10508d8d46bea0ac Mon Sep 17 00:00:00 2001 From: daflack Date: Tue, 20 Jan 2026 14:25:52 +0000 Subject: [PATCH 74/91] Adds saturation equivalent potenital temperature documentation --- src/CSET/operators/temperature.py | 54 +++++++++++++++++++------------ 1 file changed, 34 insertions(+), 20 deletions(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 10f1a2142..7f9e29461 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -345,7 +345,7 @@ def equivalent_potential_temperature( Returns ------- iris.cube.Cube | iris.cube.CubeList - Calculated equivalent potential temperature. + Calculated equivalent potential temperature in Kelvin. Notes ----- @@ -367,16 +367,13 @@ def equivalent_potential_temperature( `humidity.mixing_ratio_from_relative_humidity` and the potential temperature is calculated from `temperature.potential_temperature`. - This is a simplification of [Paluch79]_ following [Bolton80b]_ and [Emanuel94]_, + This equation is a simplification of [Paluch79]_ following [Emanuel94]_, whilst still holding reasonable accuracy. All cubes must be on the same grid. References ---------- - .. [Bolton80b] Bolton. D. (1980) "The computation of equivalent potential - temperature". Monthly Weather Review, vol. 108, 1046-1053, - doi: 10.1175/1520-0493(1980)108<1046:TCOEPT>2.0.CO;2 .. [Emanuel94] Emanuel, K.A. (1994) "Atmospheric Convection" Oxford University Press, 580 pp. .. [Paluch79] Paluch, I.R., (1979) "The entrainment mechanism in Colorado @@ -507,25 +504,42 @@ def saturation_equivalent_potential_temperature( Arguments --------- - s + temperature: iris.cube.Cube | iris.cube.CubeList + Cubes of temperature. + pressure: iris.cube.Cube | iris.cube.CubeList + Cubes of pressure. Returns ------- - s + iris.cube.Cube | iris.cube.CubeList + Calculated saturation equivalent potenital temperature in Kelvin. Notes ----- - s + The saturation equivalent potential temperature, also referred to as the + saturation potential temperature is as the equivalent potenital temperature + following a saturated process throughout. It is calculated as - References - ---------- - s + .. math:: \theta_{es} = \theta * exp\left(\frac{L_v w}{c_p T} \right) + + for :math:`\theta_{es}` the saturation equivalent potenital temperature, + :math:`\theta` the potential temperature, w the mixing ratio, + :math:`c_p` the specific heat capacity of dry air (1005.7 + :math:`J kg^{-1} K^{-1}`), :math:`L_v` the latent heat of vapourization + (2.5 x :math:`10^6 J kg^{-1} K^{-1}`), and T the temperature. + + As a saturated process is assumed throughout the RH multiplier in the + equivalent potenital temperature will always be a value of one, and is thus + omitted. In this operator the potential temperature is calculated from + `temperature.potential_temperature`. + + All cubes must be on the same grid. Examples -------- - >>> s + >>> theta_es = temperature.saturation_equivalent_potential_temperature(T, P) """ - theta_s = iris.cube.CubeList([]) + theta_es = iris.cube.CubeList([]) for T, P in zip( iter_maybe(temperature), iter_maybe(pressure), @@ -535,11 +549,11 @@ def saturation_equivalent_potential_temperature( ws = saturation_mixing_ratio(T, P) second_term_power = LV * ws / (CPD * T) second_term = np.exp(second_term_power.core_data()) - TH_S = theta * second_term - TH_S.rename("saturation_equivalent_potential_temperature") - TH_S.units = "K" - theta_s.append(TH_S) - if len(theta_s) == 1: - return theta_s[0] + TH_ES = theta * second_term + TH_ES.rename("saturation_equivalent_potential_temperature") + TH_ES.units = "K" + theta_es.append(TH_ES) + if len(theta_es) == 1: + return theta_es[0] else: - return theta_s + return theta_es From 0f9eaaebeba0cb6e1f4e897cb500ef5d74f8e81c Mon Sep 17 00:00:00 2001 From: daflack Date: Tue, 20 Jan 2026 14:28:32 +0000 Subject: [PATCH 75/91] Update copyright years --- docs/source/conf.py | 2 +- docs/source/index.rst | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/docs/source/conf.py b/docs/source/conf.py index e104f4220..1834f0222 100644 --- a/docs/source/conf.py +++ b/docs/source/conf.py @@ -8,7 +8,7 @@ # https://www.sphinx-doc.org/en/master/usage/configuration.html#project-information project = "CSET" -copyright = "Crown copyright, Met Office (2022-2025) and CSET contributors." +copyright = "Crown copyright, Met Office (2022-2026) and CSET contributors." author = "Met Office and Partners" # -- General configuration --------------------------------------------------- diff --git a/docs/source/index.rst b/docs/source/index.rst index 1d103a465..53c35565c 100644 --- a/docs/source/index.rst +++ b/docs/source/index.rst @@ -104,7 +104,7 @@ Conduct`_. Licence ------- -© Crown copyright, Met Office (2022-2025) and CSET contributors. +© Crown copyright, Met Office (2022-2026) and CSET contributors. Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except in compliance with the License. You may obtain a From d19244892b8b5d94640afc31786221c30f127b77 Mon Sep 17 00:00:00 2001 From: Sylvia Bohnenstengel <62748926+Sylviabohnenstengel@users.noreply.github.com> Date: Mon, 26 Jan 2026 16:28:12 +0000 Subject: [PATCH 76/91] Update humidity.py --- src/CSET/operators/humidity.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index 1633ed9ae..12322ce10 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -45,7 +45,7 @@ def mixing_ratio_from_specific_humidity( .. math:: w = \frac{q}{1 - q} - for w the mixing ratio and q the specific humidity. + with w the mixing ratio and q the specific humidity. Larger mixing ratios imply more moisture in the atmosphere. The mixing ratio will have the same units as the specific humidity (kg/kg). From 7266f67c85f0ebca44ca6ab0f84b285425db60ad Mon Sep 17 00:00:00 2001 From: Sylvia Bohnenstengel <62748926+Sylviabohnenstengel@users.noreply.github.com> Date: Mon, 26 Jan 2026 16:29:09 +0000 Subject: [PATCH 77/91] Update humidity.py --- src/CSET/operators/humidity.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index 12322ce10..d7a47e7ce 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -89,7 +89,7 @@ def specific_humidity_from_mixing_ratio( .. math:: q = \frac{w}{1 + w} - for q the specific humidity and w the mixing ratio. + with q the specific humidity and w the mixing ratio. A larger specific humidity implies a more moist atmosphere. The specific humidity will have the same units as the mixing ratio (kg/kg). From 542103960d9b1ebb3217e76ad9e37d07b4e64d7f Mon Sep 17 00:00:00 2001 From: Sylvia Bohnenstengel <62748926+Sylviabohnenstengel@users.noreply.github.com> Date: Mon, 26 Jan 2026 16:52:40 +0000 Subject: [PATCH 78/91] Update humidity.py --- src/CSET/operators/humidity.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index d7a47e7ce..ed4b64b1f 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -241,7 +241,7 @@ def mixing_ratio_from_relative_humidity( .. math:: w = RH * w_s for w the mixing ratio, :math:`w_s` the saturation mixing ratio, and - RH the relative humidity. + RH the relative humidity. RH is converted to dimensionless fraction rather than percentage. The operator usese `humidity.saturation_mixing_ratio` to calculate the saturation mixing ratio from the temperature and pressure. The relative From 6de578c8bcd1af31630979847d3de74f97033e7a Mon Sep 17 00:00:00 2001 From: Sylvia Bohnenstengel <62748926+Sylviabohnenstengel@users.noreply.github.com> Date: Mon, 26 Jan 2026 16:54:58 +0000 Subject: [PATCH 79/91] Update humidity.py --- src/CSET/operators/humidity.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index ed4b64b1f..019b1cb04 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -241,9 +241,10 @@ def mixing_ratio_from_relative_humidity( .. math:: w = RH * w_s for w the mixing ratio, :math:`w_s` the saturation mixing ratio, and - RH the relative humidity. RH is converted to dimensionless fraction rather than percentage. + RH the relative humidity. RH is converted to dimensionless fraction rather + than percentage. - The operator usese `humidity.saturation_mixing_ratio` to calculate the + The operator uses `humidity.saturation_mixing_ratio` to calculate the saturation mixing ratio from the temperature and pressure. The relative humidity is converted into a decimal before the multiplication occurs. From bb9e6b5d1df0cbdee68034a9cf52ed1572bc6917 Mon Sep 17 00:00:00 2001 From: Sylvia Bohnenstengel <62748926+Sylviabohnenstengel@users.noreply.github.com> Date: Mon, 26 Jan 2026 17:28:37 +0000 Subject: [PATCH 80/91] Update humidity.py --- src/CSET/operators/humidity.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index 019b1cb04..f7d5a4d08 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -277,7 +277,7 @@ def specific_humidity_from_relative_humidity( pressure: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - r"""Calculate the mixing ratio from RH. + r"""Calculate the specific humidity from relative humidity. Arguments --------- From 081e872eac2fa014f9a4822e9b011949d51c6be1 Mon Sep 17 00:00:00 2001 From: Sylvia Bohnenstengel <62748926+Sylviabohnenstengel@users.noreply.github.com> Date: Mon, 26 Jan 2026 17:29:39 +0000 Subject: [PATCH 81/91] Update humidity.py --- src/CSET/operators/humidity.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/operators/humidity.py b/src/CSET/operators/humidity.py index f7d5a4d08..85856f14a 100644 --- a/src/CSET/operators/humidity.py +++ b/src/CSET/operators/humidity.py @@ -303,7 +303,7 @@ def specific_humidity_from_relative_humidity( for q the specific humidity, :math:`q_{sat}` the saturation specific humidity, and RH the relative humidity. - The operator usese `humidity.saturation_specific_humidity` to calculate the + The operator uses `humidity.saturation_specific_humidity` to calculate the saturation specific humidity from the temperature and pressure. The relative humidity is converted into a decimal before the multiplication occurs. From 01240c6f2807090b9109ddb25d01dcb01a87a0e5 Mon Sep 17 00:00:00 2001 From: Sylvia Bohnenstengel <62748926+Sylviabohnenstengel@users.noreply.github.com> Date: Tue, 27 Jan 2026 14:16:14 +0000 Subject: [PATCH 82/91] Update pressure.py --- src/CSET/operators/pressure.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index 6fa057775..2907a42ff 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -42,7 +42,7 @@ def vapour_pressure( ----- The vapour pressure represents the pressure exerted by water vapour on the atmosphere. It is related to atmospheric temperature through the - Calusius-Clapeyron equation. There are several different formulations based + Clausius-Clapeyron equation. There are several different formulations based on empirical relations that are used to calculate the vapour pressure. Here, we calculate the vapour pressure based on [Tetens30]_ equation: From dae1ab544fdf6da9b7b87a990c221224d8ffbfe7 Mon Sep 17 00:00:00 2001 From: Sylvia Bohnenstengel <62748926+Sylviabohnenstengel@users.noreply.github.com> Date: Tue, 27 Jan 2026 14:21:00 +0000 Subject: [PATCH 83/91] Update pressure.py --- src/CSET/operators/pressure.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index 2907a42ff..9c8a7b916 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -86,7 +86,7 @@ def vapour_pressure_from_relative_humidity( temperature: iris.cube.Cube | iris.cube.CubeList, relative_humidity: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - r"""Calculate the vapour pressure using RH. + r"""Calculate the vapour pressure using relative humidity RH. Arguments --------- From 4d595312ae0094f651a7259d9ee7a45d5bf82d9a Mon Sep 17 00:00:00 2001 From: Sylvia Bohnenstengel <62748926+Sylviabohnenstengel@users.noreply.github.com> Date: Tue, 27 Jan 2026 14:22:24 +0000 Subject: [PATCH 84/91] Update pressure.py --- src/CSET/operators/pressure.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/CSET/operators/pressure.py b/src/CSET/operators/pressure.py index 9c8a7b916..544ea19c0 100644 --- a/src/CSET/operators/pressure.py +++ b/src/CSET/operators/pressure.py @@ -95,7 +95,8 @@ def vapour_pressure_from_relative_humidity( Temperature must be provided in Kelvin. relative_humidity: iris.cube.Cube | iris.cube.CubeList Cubes of relative humidity to be converted into vapour pressure. - Relative humidity will be converted to a decimal by the operator. + Relative humidity must be provided in percentage and + will be converted to a decimal by the operator. Returns ------- From fe3aabab82615c27e304f2e96fee06505516af52 Mon Sep 17 00:00:00 2001 From: Sylvia Bohnenstengel <62748926+Sylviabohnenstengel@users.noreply.github.com> Date: Tue, 27 Jan 2026 14:33:15 +0000 Subject: [PATCH 85/91] Update temperature.py --- src/CSET/operators/temperature.py | 1 + 1 file changed, 1 insertion(+) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 7f9e29461..6604d649f 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -132,6 +132,7 @@ def virtual_temperature( ratio, and :math:`\epsilon` the ratio between dry and moist air equating to 0.622. + T is given in K and w is given in kg kg-1. All cubes must be on the same grid. Examples From afde16f91da2655ac0975665afee54b786901459 Mon Sep 17 00:00:00 2001 From: Sylvia Bohnenstengel <62748926+Sylviabohnenstengel@users.noreply.github.com> Date: Tue, 27 Jan 2026 14:37:49 +0000 Subject: [PATCH 86/91] Update temperature.py --- src/CSET/operators/temperature.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 6604d649f..7e240eefd 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -245,7 +245,7 @@ def potential_temperature( Notes ----- - The potential temperature is the temperature at which an air parcel would reach + The potential temperature is the temperature an air parcel would reach if it was moved adiabatically to a reference pressure. Here we use 1000 hPa as the reference pressure. It is calculated from From 7d439b1026bf5b0a41e40963c45e6d411ebcade8 Mon Sep 17 00:00:00 2001 From: Sylvia Bohnenstengel <62748926+Sylviabohnenstengel@users.noreply.github.com> Date: Tue, 27 Jan 2026 16:32:50 +0000 Subject: [PATCH 87/91] Update temperature.py --- src/CSET/operators/temperature.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 7e240eefd..5fd1529bb 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -254,6 +254,8 @@ def potential_temperature( for :math:`\theta` the potential temperature, T the temperature, and :math:`\Pi` the exner pressure. The exner pressure is calculated using `pressure.exner_pressure`. + Temperature must be in Kelvin. + All cubes must be on the same grid. Examples From d5032870709e20f020562ebf34c234b147cdc04d Mon Sep 17 00:00:00 2001 From: Sylvia Bohnenstengel <62748926+Sylviabohnenstengel@users.noreply.github.com> Date: Tue, 27 Jan 2026 16:39:12 +0000 Subject: [PATCH 88/91] Update temperature.py --- src/CSET/operators/temperature.py | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 5fd1529bb..76dec8130 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -296,10 +296,11 @@ def virtual_potential_temperature( Notes ----- - The virtual potenital temperature is mechanistically equivalent to the + The virtual potential temperature is mechanistically equivalent to the potential temperature, except rather than using the (dry-bulb) temperature - the virtual temperature is adiabatically moved to a reference pressure - (1000 hPa). It is calculated as + the virtual temperature used. The virtual potential temperature is the potential + temperature a parcel would have if its temperature was replaced by its virtual temperature + and then the parcel is brought adiabatically to 1000 hPa. It is calculated as .. math:: \theta_v = \frac{T_v}{\Pi} From a0e81cdb0499bfc6b2767068e07054fa0428ed82 Mon Sep 17 00:00:00 2001 From: Sylvia Bohnenstengel <62748926+Sylviabohnenstengel@users.noreply.github.com> Date: Tue, 27 Jan 2026 17:18:14 +0000 Subject: [PATCH 89/91] Update temperature.py --- src/CSET/operators/temperature.py | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 76dec8130..b449eacfa 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -367,6 +367,10 @@ def equivalent_potential_temperature( (1005.7 :math:`J kg^{-1} K^{-1}`), :math:`L_v` the latent heat of vapourization (2.5 x :math:`10^6 J kg^{-1} K^{-1}`), and T the temperature. + Potential temperature and temperature in K. + Relative humidity in percentage and will be converted to fraction. + Mixing ratio in kg kg-1 (dimensionless). + In this operator the mixing ratio is calculated from `humidity.mixing_ratio_from_relative_humidity` and the potential temperature is calculated from `temperature.potential_temperature`. From 1e3e9fc8d3757fb624ccc091d3d3dabc5e4c36d6 Mon Sep 17 00:00:00 2001 From: Sylvia Bohnenstengel <62748926+Sylviabohnenstengel@users.noreply.github.com> Date: Tue, 27 Jan 2026 17:22:22 +0000 Subject: [PATCH 90/91] Update temperature.py --- src/CSET/operators/temperature.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index b449eacfa..24caa3e37 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -365,7 +365,7 @@ def equivalent_potential_temperature( :math:`R_v` the specific gas constant of water vapour (461 :math:`J kg^{-1} K^{-1}`), :math:`c_p` the specific heat capacity of dry air (1005.7 :math:`J kg^{-1} K^{-1}`), :math:`L_v` the latent heat of vapourization - (2.5 x :math:`10^6 J kg^{-1} K^{-1}`), and T the temperature. + (2.5 x :math:`10^6 J kg^{-1}`), and T the temperature. Potential temperature and temperature in K. Relative humidity in percentage and will be converted to fraction. From fd2973645309d7873279c8ec0059229562a87812 Mon Sep 17 00:00:00 2001 From: Sylvia Bohnenstengel <62748926+Sylviabohnenstengel@users.noreply.github.com> Date: Tue, 27 Jan 2026 17:28:03 +0000 Subject: [PATCH 91/91] Update temperature.py changed units of L_v amongst other small changes. --- src/CSET/operators/temperature.py | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/src/CSET/operators/temperature.py b/src/CSET/operators/temperature.py index 24caa3e37..1ecbcadc7 100644 --- a/src/CSET/operators/temperature.py +++ b/src/CSET/operators/temperature.py @@ -440,7 +440,7 @@ def wet_bulb_potential_temperature( Notes ----- The wet-bulb potential temperature represents the temperature an air parcel - would have if coold adiabatically until saturation and then taken down to + would have if cooled adiabatically until saturation and then taken down to a reference pressure (1000 hPa) whilst conserving the moisture (i.e. along a psuedoadiabat). It can be calculated either through a series of iterations, or through empirical relations. Here, we use the [DaviesJones08]_ formulation: @@ -448,7 +448,7 @@ def wet_bulb_potential_temperature( .. math:: \theta_w = \theta_e - exp\left(\frac{a_0 + a_1 X + a_2 X^2 + a_3 X^3 + a_4 X^4}{1.0 + b_1 X + b_2 X^2 + b_3 X^3 + b_4 X^4} \right) for :math:`\theta_w` the wet-bulb potential temperature, :math:`\theta_e` the - equivalent potenital temperature, X = :math:`\theta_e` / 273.15 K, and a + equivalent potential temperature, X = :math:`\theta_e` / 273.15 K, and a series of constants :math:`a_0` = 7.101574, :math:`a_1` = -20.68208, :math:`a_2` = 16.11182, :math:`a_3` = 2.574631, :math:`a_4` = -5.205688, :math:`b_1` = -3.552497, :math:`b_2` = 3.781782, :math:`b_3` = -0.6899655, and @@ -508,7 +508,7 @@ def saturation_equivalent_potential_temperature( temperature: iris.cube.Cube | iris.cube.CubeList, pressure: iris.cube.Cube | iris.cube.CubeList, ) -> iris.cube.Cube | iris.cube.CubeList: - r"""Calculate the saturation equivalent potenital temperature. + r"""Calculate the saturation equivalent potential temperature. Arguments --------- @@ -525,16 +525,16 @@ def saturation_equivalent_potential_temperature( Notes ----- The saturation equivalent potential temperature, also referred to as the - saturation potential temperature is as the equivalent potenital temperature + saturation potential temperature is as the equivalent potential temperature following a saturated process throughout. It is calculated as .. math:: \theta_{es} = \theta * exp\left(\frac{L_v w}{c_p T} \right) - for :math:`\theta_{es}` the saturation equivalent potenital temperature, + for :math:`\theta_{es}` the saturation equivalent potential temperature, :math:`\theta` the potential temperature, w the mixing ratio, :math:`c_p` the specific heat capacity of dry air (1005.7 :math:`J kg^{-1} K^{-1}`), :math:`L_v` the latent heat of vapourization - (2.5 x :math:`10^6 J kg^{-1} K^{-1}`), and T the temperature. + (2.5 x :math:`10^6 J kg^{-1}`), and T the temperature. As a saturated process is assumed throughout the RH multiplier in the equivalent potenital temperature will always be a value of one, and is thus