Skip to content

Commit 35a86d3

Browse files
committed
Fix methods names and qa in vertical
1 parent a593136 commit 35a86d3

File tree

5 files changed

+35
-14
lines changed

5 files changed

+35
-14
lines changed

src/earthkit/meteo/thermo/array/thermo.py

+2-2
Original file line numberDiff line numberDiff line change
@@ -1730,8 +1730,9 @@ def wet_bulb_potential_temperature_from_specific_humidity(t, q, p, ept_method="i
17301730
else:
17311731
return temperature_on_moist_adiabat(ept, constants.p0, ept_method=ept_method, t_method=t_method)
17321732

1733+
17331734
def specific_gas_constant(q):
1734-
"""Computes the specific gas constant of moist air
1735+
"""Computes the specific gas constant of moist air.
17351736
(specific content of cloud particles and hydrometeors are neglected).
17361737
17371738
Parameters
@@ -1747,4 +1748,3 @@ def specific_gas_constant(q):
17471748

17481749
R = constants.Rd + (constants.Rv - constants.Rd) * q
17491750
return R
1750-

src/earthkit/meteo/thermo/thermo.py

+3-2
Original file line numberDiff line numberDiff line change
@@ -161,5 +161,6 @@ def wet_bulb_potential_temperature_from_dewpoint(*args, **kwargs):
161161
def wet_bulb_potential_temperature_from_specific_humidity(*args, **kwargs):
162162
return array.wet_bulb_potential_temperature_from_specific_humidity(*args, **kwargs)
163163

164-
def specific_gas_constant(*arg, **kwargs):
165-
return array.specific_gas_constant(*args, **kwargs)
164+
165+
def specific_gas_constant(*args, **kwargs):
166+
return array.specific_gas_constant(*args, **kwargs)

src/earthkit/meteo/vertical/array/vertical.py

+15-9
Original file line numberDiff line numberDiff line change
@@ -7,17 +7,14 @@
77
# nor does it submit to any jurisdiction.
88

99
from typing import Any
10-
from typing import Dict
11-
from typing import List
12-
from typing import Optional
1310
from typing import Tuple
1411
from typing import Union
1512

1613
import numpy as np
1714
from numpy.typing import NDArray
1815

1916
from earthkit.meteo import constants
20-
from earthkit.meteo.thermo import specific_gas_consant
17+
2118

2219
def pressure_at_model_levels(
2320
A: NDArray[Any], B: NDArray[Any], sp: Union[float, NDArray[Any]]
@@ -82,22 +79,29 @@ def pressure_at_model_levels(
8279
# calculate alpha
8380
alpha = np.zeros(new_shape_full)
8481

85-
alpha[1:, ...] = 1.0 - p_half_level[1:-1, ...] / (p_half_level[2:, ...] - p_half_level[1:-1, ...]) * delta[1:, ...]
82+
alpha[1:, ...] = (
83+
1.0 - p_half_level[1:-1, ...] / (p_half_level[2:, ...] - p_half_level[1:-1, ...]) * delta[1:, ...]
84+
)
8685

8786
# pressure at highest half level <= 0.1
8887
if np.any(p_half_level[0, ...] <= PRESSURE_TOA):
8988
alpha[0, ...] = 1.0 # ARPEGE choice, ECMWF IFS uses log(2)
9089
# pressure at highest half level > 0.1
9190
else:
92-
alpha[0, ...] = 1.0 - p_half_level[0, ...] / (p_half_level[1, ...] - p_half_level[0, ...]) * delta[0, ...]
91+
alpha[0, ...] = (
92+
1.0 - p_half_level[0, ...] / (p_half_level[1, ...] - p_half_level[0, ...]) * delta[0, ...]
93+
)
9394

9495
# calculate pressure on model full levels
9596
# TODO: is there a faster way to calculate the averages?
9697
# TODO: introduce option to calculate full levels in more complicated way
97-
p_full_level = np.apply_along_axis(lambda m: np.convolve(m, np.ones(2) / 2, mode="valid"), axis=0, arr=p_half_level)
98+
p_full_level = np.apply_along_axis(
99+
lambda m: np.convolve(m, np.ones(2) / 2, mode="valid"), axis=0, arr=p_half_level
100+
)
98101

99102
return p_full_level, p_half_level, delta, alpha
100103

104+
101105
def relative_geopotential_thickness(alpha: NDArray[Any], q: NDArray[Any], T: NDArray[Any]) -> NDArray[Any]:
102106
"""Calculates the geopotential thickness w.r.t the surface on model full-levels.
103107
@@ -115,13 +119,15 @@ def relative_geopotential_thickness(alpha: NDArray[Any], q: NDArray[Any], T: NDA
115119
ndarray
116120
geopotential thickness of model full-levels w.r.t. the surface
117121
"""
122+
from earthkit.meteo import specific_gas_constant
118123

119-
R = calc_specific_gas_constant(q)
124+
R = specific_gas_constant(q)
120125
dphi = np.cumsum(np.flip(alpha * R * T, axis=0), axis=0)
121126
dphi = np.flip(dphi, axis=0)
122127

123128
return dphi
124129

130+
125131
def pressure_at_height_level(
126132
height: float, q: NDArray[Any], T: NDArray[Any], sp: NDArray[Any], A: NDArray[Any], B: NDArray[Any]
127133
) -> Union[float, NDArray[Any]]:
@@ -154,7 +160,7 @@ def pressure_at_height_level(
154160
tdphi = height * constants.g
155161

156162
# pressure(-related) variables
157-
p_full, p_half, _, alpha = model_level_pressure(A, B, sp)
163+
p_full, p_half, _, alpha = pressure_at_model_levels(A, B, sp)
158164

159165
# relative geopot. thickness of full levels
160166
dphi = relative_geopotential_thickness(alpha, q, T)

src/earthkit/meteo/vertical/vertical.py

+4-1
Original file line numberDiff line numberDiff line change
@@ -9,11 +9,14 @@
99

1010
from . import array
1111

12+
1213
def pressure_at_model_levels(*args, **kwargs):
1314
return array.pressure_at_model_levels(*args, **kwargs)
1415

16+
1517
def relative_geopotential_thickness(*arg, **kwargs):
1618
return array.relative_geopotential_thickness(*arg, **kwargs)
1719

20+
1821
def pressure_at_height_level(*args, **kwargs):
19-
return array.pressure_at_height_level(*args, **kwargs)
22+
return array.pressure_at_height_level(*args, **kwargs)

tests/vertical/test_vertical.py

+11
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,11 @@
1+
# (C) Copyright 2021 ECMWF.
2+
#
3+
# This software is licensed under the terms of the Apache Licence Version 2.0
4+
# which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
5+
# In applying this licence, ECMWF does not waive the privileges and immunities
6+
# granted to it by virtue of its status as an intergovernmental organisation
7+
# nor does it submit to any jurisdiction.
8+
#
9+
10+
11+
# TODO: add tests

0 commit comments

Comments
 (0)