Completed
Push — master ( 5b2ef8...f509c8 )
by Ryan
01:01
created

test_no_el_multi_crossing()   A

Complexity

Conditions 3

Size

Total Lines 12

Duplication

Lines 0
Ratio 0 %

Importance

Changes 0
Metric Value
cc 3
dl 0
loc 12
rs 9.4285
c 0
b 0
f 0
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 (cape_cin, density, dewpoint, dewpoint_rh, dry_lapse, el,
10
                        equivalent_potential_temperature, lcl, lfc, mixing_ratio,
11
                        mixing_ratio_from_specific_humidity, moist_lapse,
12
                        most_unstable_parcel, parcel_profile, potential_temperature,
13
                        psychrometric_vapor_pressure_wet,
14
                        relative_humidity_from_mixing_ratio,
15
                        relative_humidity_from_specific_humidity,
16
                        relative_humidity_wet_psychrometric,
17
                        saturation_mixing_ratio,
18
                        saturation_vapor_pressure, vapor_pressure,
19
                        virtual_potential_temperature, virtual_temperature)
20
21
from metpy.calc.thermo import _find_append_zero_crossings
22
from metpy.testing import assert_almost_equal, assert_array_almost_equal, assert_nan
23
from metpy.units import units
24
25
26
def test_potential_temperature():
27
    """Test potential_temperature calculation."""
28
    temp = np.array([278., 283., 291., 298.]) * units.kelvin
29
    pres = np.array([900., 500., 300., 100.]) * units.mbar
30
    real_th = np.array([286.493, 344.961, 410.4335, 575.236]) * units.kelvin
31
    assert_array_almost_equal(potential_temperature(pres, temp), real_th, 3)
32
33
34
def test_scalar():
35
    """Test potential_temperature accepts scalar values."""
36
    assert_almost_equal(potential_temperature(1000. * units.mbar, 293. * units.kelvin),
37
                        293. * units.kelvin, 4)
38
    assert_almost_equal(potential_temperature(800. * units.mbar, 293. * units.kelvin),
39
                        312.2828 * units.kelvin, 4)
40
41
42
def test_fahrenheit():
43
    """Test that potential_temperature handles temperature values in Fahrenheit."""
44
    assert_almost_equal(potential_temperature(800. * units.mbar, 68. * units.degF),
45
                        (312.444 * units.kelvin).to(units.degF), 2)
46
47
48
def test_pot_temp_inhg():
49
    """Test that potential_temperature can handle pressure not in mb (issue #165)."""
50
    assert_almost_equal(potential_temperature(29.92 * units.inHg, 29 * units.degC),
51
                        301.019735 * units.kelvin, 4)
52
53
54
def test_dry_lapse():
55
    """Test dry_lapse calculation."""
56
    levels = np.array([1000, 900, 864.89]) * units.mbar
57
    temps = dry_lapse(levels, 303.15 * units.kelvin)
58
    assert_array_almost_equal(temps,
59
                              np.array([303.15, 294.16, 290.83]) * units.kelvin, 2)
60
61
62
def test_dry_lapse_2_levels():
63
    """Test dry_lapse calculation when given only two levels."""
64
    temps = dry_lapse(np.array([1000., 500.]) * units.mbar, 293. * units.kelvin)
65
    assert_array_almost_equal(temps, [293., 240.3723] * units.kelvin, 4)
66
67
68
def test_moist_lapse():
69
    """Test moist_lapse calculation."""
70
    temp = moist_lapse(np.array([1000., 800., 600., 500., 400.]) * units.mbar,
71
                       293. * units.kelvin)
72
    true_temp = np.array([293, 284.64, 272.81, 264.42, 252.91]) * units.kelvin
73
    assert_array_almost_equal(temp, true_temp, 2)
74
75
76
def test_moist_lapse_degc():
77
    """Test moist_lapse with Celsius temperatures."""
78
    temp = moist_lapse(np.array([1000., 800., 600., 500., 400.]) * units.mbar,
79
                       19.85 * units.degC)
80
    true_temp = np.array([293, 284.64, 272.81, 264.42, 252.91]) * units.kelvin
81
    assert_array_almost_equal(temp, true_temp, 2)
82
83
84
def test_parcel_profile():
85
    """Test parcel profile calculation."""
86
    levels = np.array([1000., 900., 800., 700., 600., 500., 400.]) * units.mbar
87
    true_prof = np.array([303.15, 294.16, 288.026, 283.073, 277.058, 269.402,
88
                          258.966]) * units.kelvin
89
90
    prof = parcel_profile(levels, 30. * units.degC, 20. * units.degC)
91
    assert_array_almost_equal(prof, true_prof, 2)
92
93
94
def test_parcel_profile_saturated():
95
    """Test parcel_profile works when LCL in levels (issue #232)."""
96
    levels = np.array([1000., 700., 500.]) * units.mbar
97
    true_prof = np.array([296.95, 284.381, 271.123]) * units.kelvin
98
99
    prof = parcel_profile(levels, 23.8 * units.degC, 23.8 * units.degC)
100
    assert_array_almost_equal(prof, true_prof, 2)
101
102
103
def test_sat_vapor_pressure():
104
    """Test saturation_vapor_pressure calculation."""
105
    temp = np.array([5., 10., 18., 25.]) * units.degC
106
    real_es = np.array([8.72, 12.27, 20.63, 31.67]) * units.mbar
107
    assert_array_almost_equal(saturation_vapor_pressure(temp), real_es, 2)
108
109
110
def test_sat_vapor_pressure_scalar():
111
    """Test saturation_vapor_pressure handles scalar values."""
112
    es = saturation_vapor_pressure(0 * units.degC)
113
    assert_almost_equal(es, 6.112 * units.mbar, 3)
114
115
116
def test_sat_vapor_pressure_fahrenheit():
117
    """Test saturation_vapor_pressure handles temperature in Fahrenheit."""
118
    temp = np.array([50., 68.]) * units.degF
119
    real_es = np.array([12.2717, 23.3695]) * units.mbar
120
    assert_array_almost_equal(saturation_vapor_pressure(temp), real_es, 4)
121
122
123
def test_basic_dewpoint_rh():
124
    """Test dewpoint_rh function."""
125
    temp = np.array([30., 25., 10., 20., 25.]) * units.degC
126
    rh = np.array([30., 45., 55., 80., 85.]) / 100.
127
128
    real_td = np.array([11, 12, 1, 16, 22]) * units.degC
129
    assert_array_almost_equal(real_td, dewpoint_rh(temp, rh), 0)
130
131
132
def test_scalar_dewpoint_rh():
133
    """Test dewpoint_rh with scalar values."""
134
    td = dewpoint_rh(10.6 * units.degC, 0.37)
135
    assert_almost_equal(td, 26. * units.degF, 0)
136
137
138
def test_dewpoint():
139
    """Test dewpoint calculation."""
140
    assert_almost_equal(dewpoint(6.112 * units.mbar), 0. * units.degC, 2)
141
142
143
def test_dewpoint_weird_units():
144
    """Test dewpoint using non-standard units.
145
146
    Revealed from odd dimensionless units and ending up using numpy.ma math
147
    functions instead of numpy ones.
148
    """
149
    assert_almost_equal(dewpoint(15825.6 * units('g * mbar / kg')),
150
                        13.8564 * units.degC, 4)
151
152
153
def test_mixing_ratio():
154
    """Test mixing ratio calculation."""
155
    p = 998. * units.mbar
156
    e = 73.75 * units.mbar
157
    assert_almost_equal(mixing_ratio(e, p), 0.04963, 2)
158
159
160
def test_vapor_pressure():
161
    """Test vapor pressure calculation."""
162
    assert_almost_equal(vapor_pressure(998. * units.mbar, 0.04963),
163
                        73.74925 * units.mbar, 5)
164
165
166
def test_lcl():
167
    """Test LCL calculation."""
168
    lcl_pressure, lcl_temperature = lcl(1000. * units.mbar, 30. * units.degC, 20. * units.degC)
169
    assert_almost_equal(lcl_pressure, 864.761 * units.mbar, 2)
170
    assert_almost_equal(lcl_temperature, 17.676 * units.degC, 2)
171
172
173
def test_lcl_convergence():
174
    """Test LCL calculation convergence failure."""
175
    with pytest.raises(RuntimeError):
176
        lcl(1000. * units.mbar, 30. * units.degC, 20. * units.degC, max_iters=2)
177
178
179 View Code Duplication
def test_lfc_basic():
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
180
    """Test LFC calculation."""
181
    levels = np.array([959., 779.2, 751.3, 724.3, 700., 269.]) * units.mbar
182
    temperatures = np.array([22.2, 14.6, 12., 9.4, 7., -49.]) * units.celsius
183
    dewpoints = np.array([19., -11.2, -10.8, -10.4, -10., -53.2]) * units.celsius
184
    l = lfc(levels, temperatures, dewpoints)
185
    assert_almost_equal(l[0], 727.468 * units.mbar, 2)
186
    assert_almost_equal(l[1], 9.705 * units.celsius, 2)
187
188
189
def test_no_lfc():
190
    """Test LFC calculation when there is no LFC in the data."""
191
    levels = np.array([959., 867.9, 779.2, 647.5, 472.5, 321.9, 251.]) * units.mbar
192
    temperatures = np.array([22.2, 17.4, 14.6, 1.4, -17.6, -39.4, -52.5]) * units.celsius
193
    dewpoints = np.array([9., 4.3, -21.2, -26.7, -31., -53.3, -66.7]) * units.celsius
194
    lfc_pressure, lfc_temperature = lfc(levels, temperatures, dewpoints)
195
    assert assert_nan(lfc_pressure, levels.units)
196
    assert assert_nan(lfc_temperature, temperatures.units)
197
198
199
def test_lfc_inversion():
200
    """Test LFC when there is an inversion to be sure we don't pick that."""
201
    levels = np.array([963., 789., 782.3, 754.8, 728.1, 727., 700.,
202
                       571., 450., 300., 248.]) * units.mbar
203
    temperatures = np.array([25.4, 18.4, 17.8, 15.4, 12.9, 12.8,
204
                             10., -3.9, -16.3, -41.1, -51.5]) * units.celsius
205
    dewpoints = np.array([20.4, 0.4, -0.5, -4.3, -8., -8.2, -9.,
206
                          -23.9, -33.3, -54.1, -63.5]) * units.celsius
207
    l = lfc(levels, temperatures, dewpoints)
208
    assert_almost_equal(l[0], 706.0103 * units.mbar, 2)
209
    assert_almost_equal(l[1], 10.6232 * units.celsius, 2)
210
211
212
def test_lfc_equals_lcl():
213
    """Test LFC when there is no cap and the lfc is equal to the lcl."""
214
    levels = np.array([912., 905.3, 874.4, 850., 815.1, 786.6, 759.1,
215
                       748., 732.2, 700., 654.8]) * units.mbar
216
    temperatures = np.array([29.4, 28.7, 25.2, 22.4, 19.4, 16.8,
217
                             14.3, 13.2, 12.6, 11.4, 7.1]) * units.celsius
218
    dewpoints = np.array([18.4, 18.1, 16.6, 15.4, 13.2, 11.4, 9.6,
219
                          8.8, 0., -18.6, -22.9]) * units.celsius
220
    l = lfc(levels, temperatures, dewpoints)
221
    assert_almost_equal(l[0], 777.0333 * units.mbar, 2)
222
    assert_almost_equal(l[1], 15.8714 * units.celsius, 2)
223
224
225
def test_saturation_mixing_ratio():
226
    """Test saturation mixing ratio calculation."""
227
    p = 999. * units.mbar
228
    t = 288. * units.kelvin
229
    assert_almost_equal(saturation_mixing_ratio(p, t), .01068, 3)
230
231
232
def test_equivalent_potential_temperature():
233
    """Test equivalent potential temperature calculation."""
234
    p = 999. * units.mbar
235
    t = 288. * units.kelvin
236
    ept = equivalent_potential_temperature(p, t)
237
    assert_almost_equal(ept, 315.9548 * units.kelvin, 3)
238
239
240
def test_virtual_temperature():
241
    """Test virtual temperature calculation."""
242
    t = 288. * units.kelvin
243
    qv = .0016  # kg/kg
244
    tv = virtual_temperature(t, qv)
245
    assert_almost_equal(tv, 288.2796 * units.kelvin, 3)
246
247
248
def test_virtual_potential_temperature():
249
    """Test virtual potential temperature calculation."""
250
    p = 999. * units.mbar
251
    t = 288. * units.kelvin
252
    qv = .0016  # kg/kg
253
    theta_v = virtual_potential_temperature(p, t, qv)
254
    assert_almost_equal(theta_v, 288.3620 * units.kelvin, 3)
255
256
257
def test_density():
258
    """Test density calculation."""
259
    p = 999. * units.mbar
260
    t = 288. * units.kelvin
261
    qv = .0016  # kg/kg
262
    rho = density(p, t, qv).to(units.kilogram / units.meter ** 3)
263
    assert_almost_equal(rho, 1.2072 * (units.kilogram / units.meter ** 3), 3)
264
265
266
def test_el():
267
    """Test equilibrium layer calculation."""
268
    levels = np.array([959., 779.2, 751.3, 724.3, 700., 269.]) * units.mbar
269
    temperatures = np.array([22.2, 14.6, 12., 9.4, 7., -38.]) * units.celsius
270
    dewpoints = np.array([19., -11.2, -10.8, -10.4, -10., -53.2]) * units.celsius
271
    el_pressure, el_temperature = el(levels, temperatures, dewpoints)
272
    assert_almost_equal(el_pressure, 520.8700 * units.mbar, 3)
273
    assert_almost_equal(el_temperature, -11.7027 * units.degC, 3)
274
275
276
def test_no_el():
277
    """Test equilibrium layer calculation when there is no EL in the data."""
278
    levels = np.array([959., 867.9, 779.2, 647.5, 472.5, 321.9, 251.]) * units.mbar
279
    temperatures = np.array([22.2, 17.4, 14.6, 1.4, -17.6, -39.4, -52.5]) * units.celsius
280
    dewpoints = np.array([19., 14.3, -11.2, -16.7, -21., -43.3, -56.7]) * units.celsius
281
    el_pressure, el_temperature = el(levels, temperatures, dewpoints)
282
    assert assert_nan(el_pressure, levels.units)
283
    assert assert_nan(el_temperature, temperatures.units)
284
285
286
def test_no_el_multi_crossing():
287
    """Test el calculation with no el and severel parcel path-profile crossings."""
288
    levels = np.array([918., 911., 880., 873.9, 850., 848., 843.5, 818., 813.8, 785.,
289
                       773., 763., 757.5, 730.5, 700., 679., 654.4, 645.,
290
                       643.9]) * units.mbar
291
    temperatures = np.array([24.2, 22.8, 19.6, 19.1, 17., 16.8, 16.5, 15., 14.9, 14.4, 16.4,
292
                             16.2, 15.7, 13.4, 10.6, 8.4, 5.7, 4.6, 4.5]) * units.celsius
293
    dewpoints = np.array([19.5, 17.8, 16.7, 16.5, 15.8, 15.7, 15.3, 13.1, 12.9, 11.9, 6.4,
294
                          3.2, 2.6, -0.6, -4.4, -6.6, -9.3, -10.4, -10.5]) * units.celsius
295
    el_pressure, el_temperature = el(levels, temperatures, dewpoints)
296
    assert assert_nan(el_pressure, levels.units)
297
    assert assert_nan(el_temperature, temperatures.units)
298
299
300
def test_el_lfc_equals_lcl():
301
    """Test equilibrium layer calculation when the lfc equals the lcl."""
302
    levels = np.array([912., 905.3, 874.4, 850., 815.1, 786.6, 759.1, 748.,
303
                       732.3, 700., 654.8, 606.8, 562.4, 501.8, 500., 482.,
304
                       400., 393.3, 317.1, 307., 300., 252.7, 250., 200.,
305
                       199.3, 197., 190., 172., 156.6, 150., 122.9, 112.,
306
                       106.2, 100.]) * units.mbar
307
    temperatures = np.array([29.4, 28.7, 25.2, 22.4, 19.4, 16.8, 14.3,
308
                             13.2, 12.6, 11.4, 7.1, 2.2, -2.7, -10.1,
309
                             -10.3, -12.4, -23.3, -24.4, -38., -40.1, -41.1,
310
                             -49.8, -50.3, -59.1, -59.1, -59.3, -59.7, -56.3,
311
                             -56.9, -57.1, -59.1, -60.1, -58.6, -56.9]) * units.celsius
312
    dewpoints = np.array([18.4, 18.1, 16.6, 15.4, 13.2, 11.4, 9.6, 8.8, 0.,
313
                          -18.6, -22.9, -27.8, -32.7, -40.1, -40.3, -42.4, -53.3,
314
                          -54.4, -68., -70.1, -70., -70., -70., -70., -70., -70.,
315
                          -70., -70., -70., -70., -70., -70., -70., -70.]) * units.celsius
316
    el_pressure, el_temperature = el(levels, temperatures, dewpoints)
317
    assert_almost_equal(el_pressure, 175.8684 * units.mbar, 3)
318
    assert_almost_equal(el_temperature, -57.0307 * units.degC, 3)
319
320
321
def test_wet_psychrometric_vapor_pressure():
322
    """Test calculation of vapor pressure from wet and dry bulb temperatures."""
323
    p = 1013.25 * units.mbar
324
    dry_bulb_temperature = 20. * units.degC
325
    wet_bulb_temperature = 18. * units.degC
326
    psychrometric_vapor_pressure = psychrometric_vapor_pressure_wet(dry_bulb_temperature,
327
                                                                    wet_bulb_temperature, p)
328
    assert_almost_equal(psychrometric_vapor_pressure, 19.3673 * units.mbar, 3)
329
330
331
def test_wet_psychrometric_rh():
332
    """Test calculation of relative humidity from wet and dry bulb temperatures."""
333
    p = 1013.25 * units.mbar
334
    dry_bulb_temperature = 20. * units.degC
335
    wet_bulb_temperature = 18. * units.degC
336
    psychrometric_rh = relative_humidity_wet_psychrometric(dry_bulb_temperature,
337
                                                           wet_bulb_temperature, p)
338
    assert_almost_equal(psychrometric_rh, 82.8747 * units.percent, 3)
339
340
341
def test_wet_psychrometric_rh_kwargs():
342
    """Test calculation of relative humidity from wet and dry bulb temperatures."""
343
    p = 1013.25 * units.mbar
344
    dry_bulb_temperature = 20. * units.degC
345
    wet_bulb_temperature = 18. * units.degC
346
    coeff = 6.1e-4 / units.kelvin
347
    psychrometric_rh = relative_humidity_wet_psychrometric(dry_bulb_temperature,
348
                                                           wet_bulb_temperature, p,
349
                                                           psychrometer_coefficient=coeff)
350
    assert_almost_equal(psychrometric_rh, 82.9701 * units.percent, 3)
351
352
353
def test_rh_mixing_ratio():
354
    """Tests relative humidity from mixing ratio."""
355
    p = 1013.25 * units.mbar
356
    temperature = 20. * units.degC
357
    w = 0.012
358
    rh = relative_humidity_from_mixing_ratio(w, temperature, p)
359
    assert_almost_equal(rh, 81.7219 * units.percent, 3)
360
361
362
def test_mixing_ratio_from_specific_humidity():
363
    """Tests mixing ratio from specific humidity."""
364 View Code Duplication
    q = 0.012
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
365
    w = mixing_ratio_from_specific_humidity(q)
366
    assert_almost_equal(w, 0.01215, 3)
367
368
369
def test_rh_specific_humidity():
370
    """Tests relative humidity from specific humidity."""
371
    p = 1013.25 * units.mbar
372
    temperature = 20. * units.degC
373
    q = 0.012
374
    rh = relative_humidity_from_specific_humidity(q, temperature, p)
375
    assert_almost_equal(rh, 82.7145 * units.percent, 3)
376
377
378
def test_cape_cin():
379
    """Tests the basic CAPE and CIN calculation."""
380
    p = np.array([959., 779.2, 751.3, 724.3, 700., 269.]) * units.mbar
381
    temperature = np.array([22.2, 14.6, 12., 9.4, 7., -38.]) * units.celsius
382
    dewpoint = np.array([19., -11.2, -10.8, -10.4, -10., -53.2]) * units.celsius
383
    parcel_prof = parcel_profile(p, temperature[0], dewpoint[0])
384
    cape, cin = cape_cin(p, temperature, dewpoint, parcel_prof)
385
    assert_almost_equal(cape, 58.0368212 * units('joule / kilogram'), 6)
386 View Code Duplication
    assert_almost_equal(cin, -89.8073512 * units('joule / kilogram'), 6)
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
387
388
389
def test_cape_cin_no_el():
390
    """Tests that CAPE works with no EL."""
391
    p = np.array([959., 779.2, 751.3, 724.3]) * units.mbar
392
    temperature = np.array([22.2, 14.6, 12., 9.4]) * units.celsius
393
    dewpoint = np.array([19., -11.2, -10.8, -10.4]) * units.celsius
394
    parcel_prof = parcel_profile(p, temperature[0], dewpoint[0]).to('degC')
395
    cape, cin = cape_cin(p, temperature, dewpoint, parcel_prof)
396
    assert_almost_equal(cape, 0.08750805 * units('joule / kilogram'), 6)
397
    assert_almost_equal(cin, -89.8073512 * units('joule / kilogram'), 6)
398
399
400
def test_cape_cin_no_lfc():
401
    """Tests that CAPE is zero with no LFC."""
402
    p = np.array([959., 779.2, 751.3, 724.3, 700., 269.]) * units.mbar
403
    temperature = np.array([22.2, 24.6, 22., 20.4, 18., -10.]) * units.celsius
404
    dewpoint = np.array([19., -11.2, -10.8, -10.4, -10., -53.2]) * units.celsius
405
    parcel_prof = parcel_profile(p, temperature[0], dewpoint[0]).to('degC')
406
    cape, cin = cape_cin(p, temperature, dewpoint, parcel_prof)
407
    assert_almost_equal(cape, 0.0 * units('joule / kilogram'), 6)
408
    assert_almost_equal(cin, 0.0 * units('joule / kilogram'), 6)
409
410 View Code Duplication
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
411
def test_find_append_zero_crossings():
412
    """Tests finding and appending zero crossings of an x, y series."""
413
    x = np.arange(11) * units.hPa
414
    y = np.array([3, 2, 1, -1, 2, 2, 0, 1, 0, -1, 2]) * units.degC
415
    x2, y2 = _find_append_zero_crossings(x, y)
416
417
    x_truth = np.array([0., 1., 2., 2.5, 3., 3.33333333, 4., 5.,
418
                        6., 7., 8., 9., 9.33333333, 10.]) * units.hPa
419
    y_truth = np.array([3, 2, 1, 0, -1, 0, 2, 2, 0, 1, 0, -1, 0, 2]) * units.degC
420
    assert_array_almost_equal(x2, x_truth, 6)
421
    assert_almost_equal(y2, y_truth, 6)
422
423
424
def test_most_unstable_parcel():
425
    """Tests calculating the most unstable parcel."""
426
    levels = np.array([1000., 959., 867.9]) * units.mbar
427
    temperatures = np.array([18.2, 22.2, 17.4]) * units.celsius
428
    dewpoints = np.array([19., 19., 14.3]) * units.celsius
429
    ret = most_unstable_parcel(levels, temperatures, dewpoints, depth=100 * units.hPa)
430
    assert_almost_equal(ret[0], 959.0 * units.hPa, 6)
431
    assert_almost_equal(ret[1], 22.2 * units.degC, 6)
432
    assert_almost_equal(ret[2], 19.0 * units.degC, 6)
433