Completed
Pull Request — master (#485)
by Ryan
54s
created

test_get_layer_float32()   A

Complexity

Conditions 1

Size

Total Lines 15

Duplication

Lines 0
Ratio 0 %

Importance

Changes 0
Metric Value
cc 1
c 0
b 0
f 0
dl 0
loc 15
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
"""Tests for `calc.tools` module."""
5
6
import numpy as np
7
import numpy.ma as ma
8
import pytest
9
10
from metpy.calc import (find_intersections, get_layer, interp, interpolate_nans, log_interp,
11
                        nearest_intersection_idx, pressure_to_height_std,
12
                        reduce_point_density, resample_nn_1d)
13
from metpy.calc.tools import (_get_bound_pressure_height, _next_non_masked_element,
14
                              delete_masked_points)
15
from metpy.testing import assert_array_almost_equal, assert_array_equal
16
from metpy.units import units
17
18
19
def test_resample_nn():
20
    """Test 1d nearest neighbor functionality."""
21
    a = np.arange(5.)
22
    b = np.array([2, 3.8])
23
    truth = np.array([2, 4])
24
25
    assert_array_equal(truth, resample_nn_1d(a, b))
26
27
28
def test_nearest_intersection_idx():
29
    """Test nearest index to intersection functionality."""
30
    x = np.linspace(5, 30, 17)
31
    y1 = 3 * x**2
32
    y2 = 100 * x - 650
33
    truth = np.array([2, 12])
34
35
    assert_array_equal(truth, nearest_intersection_idx(y1, y2))
36
37
38
@pytest.mark.parametrize('direction, expected', [
39
    ('all', np.array([[8.88, 24.44], [238.84, 1794.53]])),
40
    ('increasing', np.array([[24.44], [1794.53]])),
41
    ('decreasing', np.array([[8.88], [238.84]]))
42
])
43
def test_find_intersections(direction, expected):
44
    """Test finding the intersection of two curves functionality."""
45
    x = np.linspace(5, 30, 17)
46
    y1 = 3 * x**2
47
    y2 = 100 * x - 650
48
    # Note: Truth is what we will get with this sampling, not the mathematical intersection
49
    assert_array_almost_equal(expected, find_intersections(x, y1, y2, direction=direction), 2)
50
51
52
def test_find_intersections_no_intersections():
53
    """Test finding the intersection of two curves with no intersections."""
54
    x = np.linspace(5, 30, 17)
55
    y1 = 3 * x + 0
56
    y2 = 5 * x + 5
57
    # Note: Truth is what we will get with this sampling, not the mathematical intersection
58
    truth = np.array([[],
59
                      []])
60
    assert_array_equal(truth, find_intersections(x, y1, y2))
61
62
63
def test_find_intersections_invalid_direction():
64
    """Test exception if an invalid direction is given."""
65
    x = np.linspace(5, 30, 17)
66
    y1 = 3 * x ** 2
67
    y2 = 100 * x - 650
68
    with pytest.raises(ValueError):
69
        find_intersections(x, y1, y2, direction='increaing')
70
71
72
def test_interpolate_nan_linear():
73
    """Test linear interpolation of arrays with NaNs in the y-coordinate."""
74
    x = np.linspace(0, 20, 15)
75
    y = 5 * x + 3
76
    nan_indexes = [1, 5, 11, 12]
77
    y_with_nan = y.copy()
78
    y_with_nan[nan_indexes] = np.nan
79
    assert_array_almost_equal(y, interpolate_nans(x, y_with_nan), 2)
80
81
82
def test_interpolate_nan_log():
83
    """Test log interpolation of arrays with NaNs in the y-coordinate."""
84
    x = np.logspace(1, 5, 15)
85
    y = 5 * np.log(x) + 3
86
    nan_indexes = [1, 5, 11, 12]
87
    y_with_nan = y.copy()
88
    y_with_nan[nan_indexes] = np.nan
89
    assert_array_almost_equal(y, interpolate_nans(x, y_with_nan, kind='log'), 2)
90
91
92
def test_interpolate_nan_invalid():
93
    """Test log interpolation with invalid parameter."""
94
    x = np.logspace(1, 5, 15)
95
    y = 5 * np.log(x) + 3
96
    with pytest.raises(ValueError):
97
        interpolate_nans(x, y, kind='loog')
98
99
100
@pytest.mark.parametrize('mask, expected_idx, expected_element', [
101
    ([False, False, False, False, False], 1, 1),
102
    ([False, True, True, False, False], 3, 3),
103
    ([False, True, True, True, True], None, None)
104
])
105
def test_non_masked_elements(mask, expected_idx, expected_element):
106
    """Test with a valid element."""
107
    a = ma.masked_array(np.arange(5), mask=mask)
108
    idx, element = _next_non_masked_element(a, 1)
109
    assert idx == expected_idx
110
    assert element == expected_element
111
112
113
@pytest.fixture
114
def thin_point_data():
115
    r"""Provide scattered points for testing."""
116
    xy = np.array([[0.8793620, 0.9005706], [0.5382446, 0.8766988], [0.6361267, 0.1198620],
117
                   [0.4127191, 0.0270573], [0.1486231, 0.3121822], [0.2607670, 0.4886657],
118
                   [0.7132257, 0.2827587], [0.4371954, 0.5660840], [0.1318544, 0.6468250],
119
                   [0.6230519, 0.0682618], [0.5069460, 0.2326285], [0.1324301, 0.5609478],
120
                   [0.7975495, 0.2109974], [0.7513574, 0.9870045], [0.9305814, 0.0685815],
121
                   [0.5271641, 0.7276889], [0.8116574, 0.4795037], [0.7017868, 0.5875983],
122
                   [0.5591604, 0.5579290], [0.1284860, 0.0968003], [0.2857064, 0.3862123]])
123 View Code Duplication
    return xy
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
124
125
126
@pytest.mark.parametrize('radius, truth',
127
                         [(2.0, np.array([1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
128
                                          0, 0, 0, 0, 0, 0, 0, 0, 0, 0], dtype=np.bool)),
129
                          (1.0, np.array([1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
130
                                          0, 0, 0, 0, 0, 0, 0, 0, 1, 0], dtype=np.bool)),
131
                          (0.3, np.array([1, 1, 1, 0, 1, 0, 0, 1, 1, 0, 0,
132
                                          0, 0, 0, 0, 0, 1, 0, 0, 0, 0], dtype=np.bool)),
133
                          (0.1, np.array([1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1,
134
                                          0, 1, 1, 1, 1, 1, 1, 1, 1, 1], dtype=np.bool))
135
                          ])
136
def test_reduce_point_density(thin_point_data, radius, truth):
137
    r"""Test that reduce_point_density works."""
138 View Code Duplication
    assert_array_equal(reduce_point_density(thin_point_data, radius=radius), truth)
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
139
140
141
@pytest.mark.parametrize('radius, truth',
142
                         [(2.0, np.array([0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
143
                                          0, 0, 0, 0, 0, 0, 0, 0, 0, 1], dtype=np.bool)),
144
                          (0.7, np.array([1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
145
                                          0, 0, 0, 1, 0, 0, 0, 0, 0, 1], dtype=np.bool)),
146
                          (0.3, np.array([1, 1, 0, 1, 0, 0, 1, 0, 1, 0, 0,
147
                                          0, 0, 0, 1, 0, 0, 0, 1, 0, 1], dtype=np.bool)),
148
                          (0.1, np.array([1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1,
149
                                          0, 1, 1, 1, 1, 1, 1, 1, 1, 1], dtype=np.bool))
150
                          ])
151
def test_reduce_point_density_priority(thin_point_data, radius, truth):
152
    r"""Test that reduce_point_density works properly with priority."""
153
    key = np.array([8, 6, 2, 8, 6, 4, 4, 8, 8, 6, 3, 4, 3, 0, 7, 4, 3, 2, 3, 3, 9])
154
    assert_array_equal(reduce_point_density(thin_point_data, radius, key), truth)
155
156
157
def test_reduce_point_density_1d():
158
    r"""Test that reduce_point_density works with 1D points."""
159
    x = np.array([1, 3, 4, 8, 9, 10])
160
    assert_array_equal(reduce_point_density(x, 2.5),
161
                       np.array([1, 0, 1, 1, 0, 0], dtype=np.bool))
162
163
164
def test_delete_masked_points():
165
    """Test deleting masked points."""
166
    a = ma.masked_array(np.arange(5), mask=[False, True, False, False, False])
167
    b = ma.masked_array(np.arange(5), mask=[False, False, False, True, False])
168
    expected = np.array([0, 2, 4])
169
    a, b = delete_masked_points(a, b)
170
    assert_array_equal(a, expected)
171
    assert_array_equal(b, expected)
172
173
174
def test_log_interp():
175
    """Test interpolating with log x-scale."""
176
    x_log = np.array([1e3, 1e4, 1e5, 1e6])
177
    y_log = np.log(x_log) * 2 + 3
178
    x_interp = np.array([5e3, 5e4, 5e5])
179
    y_interp_truth = np.array([20.0343863828, 24.6395565688, 29.2447267548])
180
    y_interp = log_interp(x_interp, x_log, y_log)
181
    assert_array_almost_equal(y_interp, y_interp_truth, 7)
182
183
184
def test_log_interp_units():
185
    """Test interpolating with log x-scale with units."""
186
    x_log = np.array([1e3, 1e4, 1e5, 1e6]) * units.hPa
187
    y_log = (np.log(x_log.m) * 2 + 3) * units.degC
188
    x_interp = np.array([5e5, 5e6, 5e7]) * units.Pa
189
    y_interp_truth = np.array([20.0343863828, 24.6395565688, 29.2447267548]) * units.degC
190
    y_interp = log_interp(x_interp, x_log, y_log)
191
    assert_array_almost_equal(y_interp, y_interp_truth, 7)
192
193
194
@pytest.fixture
195
def get_bounds_data():
196
    """Provide pressure and height data for testing layer bounds calculation."""
197
    pressures = np.linspace(1000, 100, 10) * units.hPa
198
    heights = pressure_to_height_std(pressures)
199
    return pressures, heights
200
201
202
@pytest.mark.parametrize('pressure, bound, hgts, interp, expected', [
203
    (get_bounds_data()[0], 900 * units.hPa, None, True,
204
     (900 * units.hPa, 0.9880028 * units.kilometer)),
205
    (get_bounds_data()[0], 900 * units.hPa, None, False,
206
     (900 * units.hPa, 0.9880028 * units.kilometer)),
207
    (get_bounds_data()[0], 870 * units.hPa, None, True,
208
     (870 * units.hPa, 1.2665298 * units.kilometer)),
209
    (get_bounds_data()[0], 870 * units.hPa, None, False,
210
     (900 * units.hPa, 0.9880028 * units.kilometer)),
211
    (get_bounds_data()[0], 0.9880028 * units.kilometer, None, True,
212
     (900 * units.hPa, 0.9880028 * units.kilometer)),
213
    (get_bounds_data()[0], 0.9880028 * units.kilometer, None, False,
214
     (900 * units.hPa, 0.9880028 * units.kilometer)),
215
    (get_bounds_data()[0], 1.2665298 * units.kilometer, None, True,
216
     (870 * units.hPa, 1.2665298 * units.kilometer)),
217
    (get_bounds_data()[0], 1.2665298 * units.kilometer, None, False,
218
     (900 * units.hPa, 0.9880028 * units.kilometer)),
219
    (get_bounds_data()[0], 900 * units.hPa, get_bounds_data()[1], True,
220
     (900 * units.hPa, 0.9880028 * units.kilometer)),
221
    (get_bounds_data()[0], 900 * units.hPa, get_bounds_data()[1], False,
222
     (900 * units.hPa, 0.9880028 * units.kilometer)),
223
    (get_bounds_data()[0], 870 * units.hPa, get_bounds_data()[1], True,
224
     (870 * units.hPa, 1.2643214 * units.kilometer)),
225
    (get_bounds_data()[0], 870 * units.hPa, get_bounds_data()[1], False,
226
     (900 * units.hPa, 0.9880028 * units.kilometer)),
227
    (get_bounds_data()[0], 0.9880028 * units.kilometer, get_bounds_data()[1], True,
228
     (900 * units.hPa, 0.9880028 * units.kilometer)),
229
    (get_bounds_data()[0], 0.9880028 * units.kilometer, get_bounds_data()[1], False,
230
     (900 * units.hPa, 0.9880028 * units.kilometer)),
231
    (get_bounds_data()[0], 1.2665298 * units.kilometer, get_bounds_data()[1], True,
232
     (870.9869087 * units.hPa, 1.2665298 * units.kilometer)),
233
    (get_bounds_data()[0], 1.2665298 * units.kilometer, get_bounds_data()[1], False,
234
     (900 * units.hPa, 0.9880028 * units.kilometer)),
235
    (get_bounds_data()[0], 0.98800289 * units.kilometer, get_bounds_data()[1], True,
236
     (900 * units.hPa, 0.9880028 * units.kilometer))
237
])
238
def test_get_bound_pressure_height(pressure, bound, hgts, interp, expected):
239
    """Test getting bounds in layers with various parameter combinations."""
240
    bounds = _get_bound_pressure_height(pressure, bound, heights=hgts, interpolate=interp)
241
    assert_array_almost_equal(bounds[0], expected[0], 5)
242
    assert_array_almost_equal(bounds[1], expected[1], 5)
243
244
245
def test_get_bound_invalid_bound_units():
246
    """Test that value error is raised with invalid bound units."""
247
    p = np.arange(900, 300, -100) * units.hPa
248
    with pytest.raises(ValueError):
249
        _get_bound_pressure_height(p, 100 * units.degC)
250
251
252
def test_get_bound_pressure_out_of_range():
253
    """Test when bound is out of data range in pressure."""
254
    p = np.arange(900, 300, -100) * units.hPa
255
    with pytest.raises(ValueError):
256
        _get_bound_pressure_height(p, 100 * units.hPa)
257
    with pytest.raises(ValueError):
258
        _get_bound_pressure_height(p, 1000 * units.hPa)
259
260
261
def test_get_bound_height_out_of_range():
262
    """Test when bound is out of data range in height."""
263
    p = np.arange(900, 300, -100) * units.hPa
264
    h = np.arange(1, 7) * units.kilometer
265
    with pytest.raises(ValueError):
266
        _get_bound_pressure_height(p, 8 * units.kilometer, heights=h)
267
    with pytest.raises(ValueError):
268
        _get_bound_pressure_height(p, 100 * units.meter, heights=h)
269
270
271
def test_get_layer_float32():
272
    """Test that get_layer works properly with float32 data."""
273
    p = np.asarray([940.85083008, 923.78851318, 911.42022705, 896.07220459,
274
                    876.89404297, 781.63330078], np.float32) * units('hPa')
275
    hgt = np.asarray([563.671875, 700.93817139, 806.88098145, 938.51745605,
276
                      1105.25854492, 2075.04443359], dtype=np.float32) * units.meter
277
278
    true_p_layer = np.asarray([940.85083008, 923.78851318, 911.42022705, 896.07220459,
279
                               876.89404297, 831.86472819], np.float32) * units('hPa')
280
    true_hgt_layer = np.asarray([563.671875, 700.93817139, 806.88098145, 938.51745605,
281
                                 1105.25854492, 1549.8079], dtype=np.float32) * units.meter
282
283
    p_layer, hgt_layer = get_layer(p, hgt, heights=hgt, depth=1000. * units.meter)
284
    assert_array_almost_equal(p_layer, true_p_layer, 4)
285
    assert_array_almost_equal(hgt_layer, true_hgt_layer, 4)
286
287
288
def test_get_layer_ragged_data():
289
    """Tests that error is raised for unequal length pressure and data arrays."""
290
    p = np.arange(10) * units.hPa
291
    y = np.arange(9) * units.degC
292
    with pytest.raises(ValueError):
293
        get_layer(p, y)
294
295
296
def test_get_layer_invalid_depth_units():
297
    """Tests that error is raised when depth has invalid units."""
298
    p = np.arange(10) * units.hPa
299
    y = np.arange(9) * units.degC
300
    with pytest.raises(ValueError):
301
        get_layer(p, y, depth=400 * units.degC)
302
303
304
@pytest.fixture
305
def layer_test_data():
306
    """Provide test data for testing of layer bounds."""
307
    pressure = np.arange(1000, 10, -100) * units.hPa
308
    temperature = np.linspace(25, -50, len(pressure)) * units.degC
309
    return pressure, temperature
310
311
312
@pytest.mark.parametrize('pressure, variable, heights, bottom, depth, interp, expected', [
313
    (layer_test_data()[0], layer_test_data()[1], None, None, 150 * units.hPa, True,
314
     (np.array([1000, 900, 850]) * units.hPa,
315
      np.array([25.0, 16.666666, 12.62262]) * units.degC)),
316
    (layer_test_data()[0], layer_test_data()[1], None, None, 150 * units.hPa, False,
317
     (np.array([1000, 900]) * units.hPa, np.array([25.0, 16.666666]) * units.degC)),
318
    (layer_test_data()[0], layer_test_data()[1], None, 2 * units.km, 3 * units.km, True,
319
     (np.array([794.85264282, 700., 600., 540.01696548]) * units.hPa,
320
      np.array([7.93049516, 0., -8.33333333, -13.14758845]) * units.degC))
321
])
322
def test_get_layer(pressure, variable, heights, bottom, depth, interp, expected):
323 View Code Duplication
    """Tests get_layer functionality."""
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
324
    p_layer, y_layer = get_layer(pressure, variable, heights=heights, bottom=bottom,
325
                                 depth=depth, interpolate=interp)
326
    assert_array_almost_equal(p_layer, expected[0], 5)
327
    assert_array_almost_equal(y_layer, expected[1], 5)
328
329
330
def test_log_interp_2d():
331
    """Test interpolating with log x-scale in 2 dimensions."""
332
    x_log = np.array([[1e3, 1e4, 1e5, 1e6], [1e3, 1e4, 1e5, 1e6]])
333 View Code Duplication
    y_log = np.log(x_log) * 2 + 3
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
334
    x_interp = np.array([5e3, 5e4, 5e5])
335
    y_interp_truth = np.array([20.0343863828, 24.6395565688, 29.2447267548])
336
    y_interp = log_interp(x_interp, x_log, y_log, axis=1)
337
    assert_array_almost_equal(y_interp[1], y_interp_truth, 7)
338
339
340
def test_log_interp_3d():
341
    """Test interpolating with log x-scale 3 dimensions along second axis."""
342
    x_log = np.ones((3, 4, 3)) * np.array([1e3, 1e4, 1e5, 1e6]).reshape(-1, 1)
343 View Code Duplication
    y_log = np.log(x_log) * 2 + 3
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
344
    x_interp = np.array([5e3, 5e4, 5e5])
345
    y_interp_truth = np.array([20.0343863828, 24.6395565688, 29.2447267548])
346
    y_interp = log_interp(x_interp, x_log, y_log, axis=1)
347
    assert_array_almost_equal(y_interp[0, :, 0], y_interp_truth, 7)
348
349
350
def test_log_interp_4d():
351
    """Test interpolating with log x-scale 4 dimensions."""
352
    x_log = np.ones((2, 2, 3, 4)) * np.array([1e3, 1e4, 1e5, 1e6])
353
    y_log = np.log(x_log) * 2 + 3
354
    x_interp = np.array([5e3, 5e4, 5e5])
355
    y_interp_truth = np.array([20.0343863828, 24.6395565688, 29.2447267548])
356
    y_interp = log_interp(x_interp, x_log, y_log, axis=3)
357
    assert_array_almost_equal(y_interp[0, 0, 0, :], y_interp_truth, 7)
358
359
360
def test_log_interp_2args():
361
    """Test interpolating with log x-scale with 2 arguments."""
362
    x_log = np.array([1e3, 1e4, 1e5, 1e6])
363
    y_log = np.log(x_log) * 2 + 3
364
    y_log2 = np.log(x_log) * 2 + 3
365
    x_interp = np.array([5e3, 5e4, 5e5])
366
    y_interp_truth = np.array([20.0343863828, 24.6395565688, 29.2447267548])
367
    y_interp = log_interp(x_interp, x_log, y_log, y_log2)
368
    assert_array_almost_equal(y_interp[1], y_interp_truth, 7)
369
    assert_array_almost_equal(y_interp[0], y_interp_truth, 7)
370
371
372
def test_log_interp_set_nan_above():
373
    """Test interpolating with log x-scale setting out of bounds above data to nan."""
374
    x_log = np.array([1e3, 1e4, 1e5, 1e6])
375
    y_log = np.log(x_log) * 2 + 3
376
    x_interp = np.array([1e7])
377
    y_interp_truth = np.nan
378
    with pytest.warns(Warning):
379
        y_interp = log_interp(x_interp, x_log, y_log)
380
    assert_array_almost_equal(y_interp, y_interp_truth, 7)
381
382
383
def test_log_interp_no_extrap():
384
    """Test interpolating with log x-scale setting out of bounds value error."""
385
    x_log = np.array([1e3, 1e4, 1e5, 1e6])
386
    y_log = np.log(x_log) * 2 + 3
387
    x_interp = np.array([1e7])
388
    with pytest.raises(ValueError):
389
        log_interp(x_interp, x_log, y_log, fill_value=None)
390
391
392
def test_log_interp_set_nan_below():
393
    """Test interpolating with log x-scale setting out of bounds below data to nan."""
394
    x_log = np.array([1e3, 1e4, 1e5, 1e6])
395
    y_log = np.log(x_log) * 2 + 3
396
    x_interp = 1e2
397
    y_interp_truth = np.nan
398
    with pytest.warns(Warning):
399
        y_interp = log_interp(x_interp, x_log, y_log)
400
    assert_array_almost_equal(y_interp, y_interp_truth, 7)
401
402
403
def test_interp_2args():
404
    """Test interpolation with 2 arguments."""
405
    x = np.array([1., 2., 3., 4.])
406
    y = x
407
    y2 = x
408
    x_interp = np.array([2.5000000, 3.5000000])
409
    y_interp_truth = np.array([2.5000000, 3.5000000])
410
    y_interp = interp(x_interp, x, y, y2)
411
    assert_array_almost_equal(y_interp[0], y_interp_truth, 7)
412
    assert_array_almost_equal(y_interp[1], y_interp_truth, 7)
413
414
415
def test_interp_decrease():
416
    """Test interpolation with decreasing interpolation points."""
417
    x = np.array([1., 2., 3., 4.])
418
    y = x
419
    x_interp = np.array([3.5000000, 2.5000000])
420
    y_interp_truth = np.array([3.5000000, 2.5000000])
421
    y_interp = interp(x_interp, x, y)
422
    assert_array_almost_equal(y_interp, y_interp_truth, 7)
423
424
425
def test_interp_decrease_xp():
426
    """Test interpolation with decreasing order."""
427
    x = np.array([4., 3., 2., 1.])
428
    y = x
429
    x_interp = np.array([3.5000000, 2.5000000])
430
    y_interp_truth = np.array([3.5000000, 2.5000000])
431
    y_interp = interp(x_interp, x, y)
432
    assert_array_almost_equal(y_interp, y_interp_truth, 7)
433
434
435
def test_interp_end_point():
436
    """Test interpolation with point at data endpoints."""
437
    x = np.array([1., 2., 3., 4.])
438
    y = x
439
    x_interp = np.array([1.0, 4.0])
440
    y_interp_truth = np.array([1.0, 4.0])
441
    y_interp = interp(x_interp, x, y)
442
    assert_array_almost_equal(y_interp, y_interp_truth, 7)
443