1
|
|
|
# Copyright (c) 2008-2015 MetPy Developers. |
2
|
|
|
# Distributed under the terms of the BSD 3-Clause License. |
3
|
|
|
# SPDX-License-Identifier: BSD-3-Clause |
4
|
|
|
"""Test the `thermo` module.""" |
5
|
|
|
|
6
|
|
|
import numpy as np |
7
|
|
|
import pytest |
8
|
|
|
|
9
|
|
|
from metpy.calc import (density, dewpoint, dewpoint_rh, dry_lapse, el, |
10
|
|
|
equivalent_potential_temperature, lcl, lfc, mixing_ratio, moist_lapse, |
11
|
|
|
parcel_profile, potential_temperature, |
12
|
|
|
psychrometric_vapor_pressure_wet, relative_humidity_wet_psychrometric, |
13
|
|
|
saturation_mixing_ratio, saturation_vapor_pressure, vapor_pressure, |
14
|
|
|
virtual_potential_temperature, virtual_temperature) |
15
|
|
|
|
16
|
|
|
from metpy.testing import assert_almost_equal, assert_array_almost_equal |
17
|
|
|
from metpy.units import units |
18
|
|
|
|
19
|
|
|
|
20
|
|
|
def test_potential_temperature(): |
21
|
|
|
"""Test potential_temperature calculation.""" |
22
|
|
|
temp = np.array([278., 283., 291., 298.]) * units.kelvin |
23
|
|
|
pres = np.array([900., 500., 300., 100.]) * units.mbar |
24
|
|
|
real_th = np.array([286.493, 344.961, 410.4335, 575.236]) * units.kelvin |
25
|
|
|
assert_array_almost_equal(potential_temperature(pres, temp), real_th, 3) |
26
|
|
|
|
27
|
|
|
|
28
|
|
|
def test_scalar(): |
29
|
|
|
"""Test potential_temperature accepts scalar values.""" |
30
|
|
|
assert_almost_equal(potential_temperature(1000. * units.mbar, 293. * units.kelvin), |
31
|
|
|
293. * units.kelvin, 4) |
32
|
|
|
assert_almost_equal(potential_temperature(800. * units.mbar, 293. * units.kelvin), |
33
|
|
|
312.2828 * units.kelvin, 4) |
34
|
|
|
|
35
|
|
|
|
36
|
|
|
def test_fahrenheit(): |
37
|
|
|
"""Test that potential_temperature handles temperature values in Fahrenheit.""" |
38
|
|
|
assert_almost_equal(potential_temperature(800. * units.mbar, 68. * units.degF), |
39
|
|
|
(312.444 * units.kelvin).to(units.degF), 2) |
40
|
|
|
|
41
|
|
|
|
42
|
|
|
def test_pot_temp_inhg(): |
43
|
|
|
"""Test that potential_temperature can handle pressure not in mb (issue #165).""" |
44
|
|
|
assert_almost_equal(potential_temperature(29.92 * units.inHg, 29 * units.degC), |
45
|
|
|
301.019735 * units.kelvin, 4) |
46
|
|
|
|
47
|
|
|
|
48
|
|
|
def test_dry_lapse(): |
49
|
|
|
"""Test dry_lapse calculation.""" |
50
|
|
|
levels = np.array([1000, 900, 864.89]) * units.mbar |
51
|
|
|
temps = dry_lapse(levels, 303.15 * units.kelvin) |
52
|
|
|
assert_array_almost_equal(temps, |
53
|
|
|
np.array([303.15, 294.16, 290.83]) * units.kelvin, 2) |
54
|
|
|
|
55
|
|
|
|
56
|
|
|
def test_dry_lapse_2_levels(): |
57
|
|
|
"""Test dry_lapse calculation when given only two levels.""" |
58
|
|
|
temps = dry_lapse(np.array([1000., 500.]) * units.mbar, 293. * units.kelvin) |
59
|
|
|
assert_array_almost_equal(temps, [293., 240.3723] * units.kelvin, 4) |
60
|
|
|
|
61
|
|
|
|
62
|
|
|
def test_moist_lapse(): |
63
|
|
|
"""Test moist_lapse calculation.""" |
64
|
|
|
temp = moist_lapse(np.array([1000., 800., 600., 500., 400.]) * units.mbar, |
65
|
|
|
293. * units.kelvin) |
66
|
|
|
true_temp = np.array([293, 284.64, 272.81, 264.42, 252.91]) * units.kelvin |
67
|
|
|
assert_array_almost_equal(temp, true_temp, 2) |
68
|
|
|
|
69
|
|
|
|
70
|
|
|
def test_moist_lapse_degc(): |
71
|
|
|
"""Test moist_lapse with Celsius temperatures.""" |
72
|
|
|
temp = moist_lapse(np.array([1000., 800., 600., 500., 400.]) * units.mbar, |
73
|
|
|
19.85 * units.degC) |
74
|
|
|
true_temp = np.array([293, 284.64, 272.81, 264.42, 252.91]) * units.kelvin |
75
|
|
|
assert_array_almost_equal(temp, true_temp, 2) |
76
|
|
|
|
77
|
|
|
|
78
|
|
|
def test_parcel_profile(): |
79
|
|
|
"""Test parcel profile calculation.""" |
80
|
|
|
levels = np.array([1000., 900., 800., 700., 600., 500., 400.]) * units.mbar |
81
|
|
|
true_prof = np.array([303.15, 294.16, 288.026, 283.073, 277.058, 269.402, |
82
|
|
|
258.966]) * units.kelvin |
83
|
|
|
|
84
|
|
|
prof = parcel_profile(levels, 30. * units.degC, 20. * units.degC) |
85
|
|
|
assert_array_almost_equal(prof, true_prof, 2) |
86
|
|
|
|
87
|
|
|
|
88
|
|
|
def test_parcel_profile_saturated(): |
89
|
|
|
"""Test parcel_profile works when LCL in levels (issue #232).""" |
90
|
|
|
levels = np.array([1000., 700., 500.]) * units.mbar |
91
|
|
|
true_prof = np.array([296.95, 284.381, 271.123]) * units.kelvin |
92
|
|
|
|
93
|
|
|
prof = parcel_profile(levels, 23.8 * units.degC, 23.8 * units.degC) |
94
|
|
|
assert_array_almost_equal(prof, true_prof, 2) |
95
|
|
|
|
96
|
|
|
|
97
|
|
|
def test_sat_vapor_pressure(): |
98
|
|
|
"""Test saturation_vapor_pressure calculation.""" |
99
|
|
|
temp = np.array([5., 10., 18., 25.]) * units.degC |
100
|
|
|
real_es = np.array([8.72, 12.27, 20.63, 31.67]) * units.mbar |
101
|
|
|
assert_array_almost_equal(saturation_vapor_pressure(temp), real_es, 2) |
102
|
|
|
|
103
|
|
|
|
104
|
|
|
def test_sat_vapor_pressure_scalar(): |
105
|
|
|
"""Test saturation_vapor_pressure handles scalar values.""" |
106
|
|
|
es = saturation_vapor_pressure(0 * units.degC) |
107
|
|
|
assert_almost_equal(es, 6.112 * units.mbar, 3) |
108
|
|
|
|
109
|
|
|
|
110
|
|
|
def test_sat_vapor_pressure_fahrenheit(): |
111
|
|
|
"""Test saturation_vapor_pressure handles temperature in Fahrenheit.""" |
112
|
|
|
temp = np.array([50., 68.]) * units.degF |
113
|
|
|
real_es = np.array([12.2717, 23.3695]) * units.mbar |
114
|
|
|
assert_array_almost_equal(saturation_vapor_pressure(temp), real_es, 4) |
115
|
|
|
|
116
|
|
|
|
117
|
|
|
def test_basic_dewpoint_rh(): |
118
|
|
|
"""Test dewpoint_rh function.""" |
119
|
|
|
temp = np.array([30., 25., 10., 20., 25.]) * units.degC |
120
|
|
|
rh = np.array([30., 45., 55., 80., 85.]) / 100. |
121
|
|
|
|
122
|
|
|
real_td = np.array([11, 12, 1, 16, 22]) * units.degC |
123
|
|
|
assert_array_almost_equal(real_td, dewpoint_rh(temp, rh), 0) |
124
|
|
|
|
125
|
|
|
|
126
|
|
|
def test_scalar_dewpoint_rh(): |
127
|
|
|
"""Test dewpoint_rh with scalar values.""" |
128
|
|
|
td = dewpoint_rh(10.6 * units.degC, 0.37) |
129
|
|
|
assert_almost_equal(td, 26. * units.degF, 0) |
130
|
|
|
|
131
|
|
|
|
132
|
|
|
def test_dewpoint(): |
133
|
|
|
"""Test dewpoint calculation.""" |
134
|
|
|
assert_almost_equal(dewpoint(6.112 * units.mbar), 0. * units.degC, 2) |
135
|
|
|
|
136
|
|
|
|
137
|
|
|
def test_dewpoint_weird_units(): |
138
|
|
|
"""Test dewpoint using non-standard units. |
139
|
|
|
|
140
|
|
|
Revealed from odd dimensionless units and ending up using numpy.ma math |
141
|
|
|
functions instead of numpy ones. |
142
|
|
|
""" |
143
|
|
|
assert_almost_equal(dewpoint(15825.6 * units('g * mbar / kg')), |
144
|
|
|
13.8564 * units.degC, 4) |
145
|
|
|
|
146
|
|
|
|
147
|
|
|
def test_mixing_ratio(): |
148
|
|
|
"""Test mixing ratio calculation.""" |
149
|
|
|
p = 998. * units.mbar |
150
|
|
|
e = 73.75 * units.mbar |
151
|
|
|
assert_almost_equal(mixing_ratio(e, p), 0.04963, 2) |
152
|
|
|
|
153
|
|
|
|
154
|
|
|
def test_vapor_pressure(): |
155
|
|
|
"""Test vapor pressure calculation.""" |
156
|
|
|
assert_almost_equal(vapor_pressure(998. * units.mbar, 0.04963), |
157
|
|
|
73.74925 * units.mbar, 5) |
158
|
|
|
|
159
|
|
|
|
160
|
|
|
def test_lcl(): |
161
|
|
|
"""Test LCL calculation.""" |
162
|
|
|
lcl_pressure, lcl_temperature = lcl(1000. * units.mbar, 30. * units.degC, 20. * units.degC) |
163
|
|
|
assert_almost_equal(lcl_pressure, 864.761 * units.mbar, 2) |
164
|
|
|
assert_almost_equal(lcl_temperature, 17.676 * units.degC, 2) |
165
|
|
|
|
166
|
|
|
|
167
|
|
|
def test_lcl_convergence(): |
168
|
|
|
"""Test LCL calculation convergence failure.""" |
169
|
|
|
with pytest.raises(RuntimeError): |
170
|
|
|
lcl(1000. * units.mbar, 30. * units.degC, 20. * units.degC, max_iters=2) |
171
|
|
|
|
172
|
|
|
|
173
|
|
View Code Duplication |
def test_lfc_basic(): |
|
|
|
|
174
|
|
|
"""Test LFC calculation.""" |
175
|
|
|
levels = np.array([959., 779.2, 751.3, 724.3, 700., 269.]) * units.mbar |
176
|
|
|
temperatures = np.array([22.2, 14.6, 12., 9.4, 7., -49.]) * units.celsius |
177
|
|
|
dewpoints = np.array([19., -11.2, -10.8, -10.4, -10., -53.2]) * units.celsius |
178
|
|
|
l = lfc(levels, temperatures, dewpoints) |
179
|
|
|
assert_almost_equal(l[0], 727.468 * units.mbar, 2) |
180
|
|
|
assert_almost_equal(l[1], 9.705 * units.celsius, 2) |
181
|
|
|
|
182
|
|
|
|
183
|
|
|
def test_no_lfc(): |
184
|
|
|
"""Test LFC calculation when there is no LFC in the data.""" |
185
|
|
|
levels = np.array([959., 867.9, 779.2, 647.5, 472.5, 321.9, 251.]) * units.mbar |
186
|
|
|
temperatures = np.array([22.2, 17.4, 14.6, 1.4, -17.6, -39.4, -52.5]) * units.celsius |
187
|
|
|
dewpoints = np.array([9., 4.3, -21.2, -26.7, -31., -53.3, -66.7]) * units.celsius |
188
|
|
|
lfc_pressure, lfc_temperature = lfc(levels, temperatures, dewpoints) |
189
|
|
|
assert lfc_pressure is None |
190
|
|
|
assert lfc_temperature is None |
191
|
|
|
|
192
|
|
|
|
193
|
|
View Code Duplication |
def test_lfc_inversion(): |
|
|
|
|
194
|
|
|
"""Test LFC when there is an inversion to be sure we don't pick that.""" |
195
|
|
|
levels = np.array([963., 789., 782.3, 754.8, 728.1, 727., 700., |
196
|
|
|
571., 450., 300., 248.]) * units.mbar |
197
|
|
|
temperatures = np.array([25.4, 18.4, 17.8, 15.4, 12.9, 12.8, |
198
|
|
|
10., -3.9, -16.3, -41.1, -51.5]) * units.celsius |
199
|
|
|
dewpoints = np.array([20.4, 0.4, -0.5, -4.3, -8., -8.2, -9., |
200
|
|
|
-23.9, -33.3, -54.1, -63.5]) * units.celsius |
201
|
|
|
l = lfc(levels, temperatures, dewpoints) |
202
|
|
|
assert_almost_equal(l[0], 706.0103 * units.mbar, 2) |
203
|
|
|
assert_almost_equal(l[1], 10.6232 * units.celsius, 2) |
204
|
|
|
|
205
|
|
|
|
206
|
|
|
def test_saturation_mixing_ratio(): |
207
|
|
|
"""Test saturation mixing ratio calculation.""" |
208
|
|
|
p = 999. * units.mbar |
209
|
|
|
t = 288. * units.kelvin |
210
|
|
|
assert_almost_equal(saturation_mixing_ratio(p, t), .01068, 3) |
211
|
|
|
|
212
|
|
|
|
213
|
|
|
def test_equivalent_potential_temperature(): |
214
|
|
|
"""Test equivalent potential temperature calculation.""" |
215
|
|
|
p = 999. * units.mbar |
216
|
|
|
t = 288. * units.kelvin |
217
|
|
|
ept = equivalent_potential_temperature(p, t) |
218
|
|
|
assert_almost_equal(ept, 315.9548 * units.kelvin, 3) |
219
|
|
|
|
220
|
|
|
|
221
|
|
|
def test_virtual_temperature(): |
222
|
|
|
"""Test virtual temperature calculation.""" |
223
|
|
|
t = 288. * units.kelvin |
224
|
|
|
qv = .0016 # kg/kg |
225
|
|
|
tv = virtual_temperature(t, qv) |
226
|
|
|
assert_almost_equal(tv, 288.2796 * units.kelvin, 3) |
227
|
|
|
|
228
|
|
|
|
229
|
|
|
def test_virtual_potential_temperature(): |
230
|
|
|
"""Test virtual potential temperature calculation.""" |
231
|
|
|
p = 999. * units.mbar |
232
|
|
|
t = 288. * units.kelvin |
233
|
|
|
qv = .0016 # kg/kg |
234
|
|
|
theta_v = virtual_potential_temperature(p, t, qv) |
235
|
|
|
assert_almost_equal(theta_v, 288.3620 * units.kelvin, 3) |
236
|
|
|
|
237
|
|
|
|
238
|
|
|
def test_density(): |
239
|
|
|
"""Test density calculation.""" |
240
|
|
|
p = 999. * units.mbar |
241
|
|
|
t = 288. * units.kelvin |
242
|
|
|
qv = .0016 # kg/kg |
243
|
|
|
rho = density(p, t, qv).to(units.kilogram / units.meter ** 3) |
244
|
|
|
assert_almost_equal(rho, 1.2072 * (units.kilogram / units.meter ** 3), 3) |
245
|
|
|
|
246
|
|
|
|
247
|
|
|
def test_el(): |
248
|
|
|
"""Test equilibrium layer calculation.""" |
249
|
|
|
levels = np.array([959., 779.2, 751.3, 724.3, 700., 269.]) * units.mbar |
250
|
|
|
temperatures = np.array([22.2, 14.6, 12., 9.4, 7., -38.]) * units.celsius |
251
|
|
|
dewpoints = np.array([19., -11.2, -10.8, -10.4, -10., -53.2]) * units.celsius |
252
|
|
|
el_pressure, el_temperature = el(levels, temperatures, dewpoints) |
253
|
|
|
assert_almost_equal(el_pressure, 520.8420 * units.mbar, 3) |
254
|
|
|
assert_almost_equal(el_temperature, -11.7055 * units.degC, 3) |
255
|
|
|
|
256
|
|
|
|
257
|
|
|
def test_no_el(): |
258
|
|
|
"""Test equilibrium layer calculation when there is no EL in the data.""" |
259
|
|
|
levels = np.array([959., 867.9, 779.2, 647.5, 472.5, 321.9, 251.]) * units.mbar |
260
|
|
|
temperatures = np.array([22.2, 17.4, 14.6, 1.4, -17.6, -39.4, -52.5]) * units.celsius |
261
|
|
|
dewpoints = np.array([19., 14.3, -11.2, -16.7, -21., -43.3, -56.7]) * units.celsius |
262
|
|
|
el_pressure, el_temperature = el(levels, temperatures, dewpoints) |
263
|
|
|
assert el_pressure is None |
264
|
|
|
assert el_temperature is None |
265
|
|
|
|
266
|
|
|
|
267
|
|
|
def test_wet_psychrometric_vapor_pressure(): |
268
|
|
|
"""Test calculation of vapor pressure from wet and dry bulb temperatures.""" |
269
|
|
|
p = 1013.25 * units.mbar |
270
|
|
|
dry_bulb_temperature = 20. * units.degC |
271
|
|
|
wet_bulb_temperature = 18. * units.degC |
272
|
|
|
psychrometric_vapor_pressure = psychrometric_vapor_pressure_wet(dry_bulb_temperature, |
273
|
|
|
wet_bulb_temperature, p) |
274
|
|
|
assert_almost_equal(psychrometric_vapor_pressure, 19.3673 * units.mbar, 3) |
275
|
|
|
|
276
|
|
|
|
277
|
|
|
def test_wet_psychrometric_rh(): |
278
|
|
|
"""Test calculation of relative humidity from wet and dry bulb temperatures.""" |
279
|
|
|
p = 1013.25 * units.mbar |
280
|
|
|
dry_bulb_temperature = 20. * units.degC |
281
|
|
|
wet_bulb_temperature = 18. * units.degC |
282
|
|
|
psychrometric_rh = relative_humidity_wet_psychrometric(dry_bulb_temperature, |
283
|
|
|
wet_bulb_temperature, p) |
284
|
|
|
assert_almost_equal(psychrometric_rh, 82.8747 * units.percent, 3) |
285
|
|
|
|
286
|
|
|
|
287
|
|
|
def test_wet_psychrometric_rh_kwargs(): |
288
|
|
|
"""Test calculation of relative humidity from wet and dry bulb temperatures.""" |
289
|
|
|
p = 1013.25 * units.mbar |
290
|
|
|
dry_bulb_temperature = 20. * units.degC |
291
|
|
|
wet_bulb_temperature = 18. * units.degC |
292
|
|
|
coeff = 6.1e-4 / units.kelvin |
293
|
|
|
psychrometric_rh = relative_humidity_wet_psychrometric(dry_bulb_temperature, |
294
|
|
|
wet_bulb_temperature, p, |
295
|
|
|
psychrometer_coefficient=coeff) |
296
|
|
|
assert_almost_equal(psychrometric_rh, 82.9701 * units.percent, 3) |
297
|
|
|
|