Completed
Pull Request — master (#494)
by
unknown
58s
created

test_mucape_mucin_surface()   A

Complexity

Conditions 1

Size

Total Lines 8

Duplication

Lines 0
Ratio 0 %

Importance

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