Completed
Push — master ( 10e33f...5d2a6a )
by Ryan
17s
created

test_interp_end_point()   A

Complexity

Conditions 1

Size

Total Lines 8

Duplication

Lines 0
Ratio 0 %

Importance

Changes 1
Bugs 0 Features 1
Metric Value
cc 1
c 1
b 0
f 1
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
"""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_ragged_data():
272
    """Tests that error is raised for unequal length pressure and data arrays."""
273
    p = np.arange(10) * units.hPa
274
    y = np.arange(9) * units.degC
275
    with pytest.raises(ValueError):
276
        get_layer(p, y)
277
278
279
def test_get_layer_invalid_depth_units():
280
    """Tests that error is raised when depth has invalid units."""
281
    p = np.arange(10) * units.hPa
282
    y = np.arange(9) * units.degC
283
    with pytest.raises(ValueError):
284
        get_layer(p, y, depth=400 * units.degC)
285
286
287
@pytest.fixture
288
def layer_test_data():
289
    """Provide test data for testing of layer bounds."""
290
    pressure = np.arange(1000, 10, -100) * units.hPa
291
    temperature = np.linspace(25, -50, len(pressure)) * units.degC
292
    return pressure, temperature
293
294
295
@pytest.mark.parametrize('pressure, variable, heights, bottom, depth, interp, expected', [
296
    (layer_test_data()[0], layer_test_data()[1], None, None, 150 * units.hPa, True,
297
     (np.array([1000, 900, 850]) * units.hPa,
298
      np.array([25.0, 16.666666, 12.62262]) * units.degC)),
299
    (layer_test_data()[0], layer_test_data()[1], None, None, 150 * units.hPa, False,
300
     (np.array([1000, 900]) * units.hPa, np.array([25.0, 16.666666]) * units.degC)),
301
    (layer_test_data()[0], layer_test_data()[1], None, 2 * units.km, 3 * units.km, True,
302
     (np.array([794.85264282, 700., 600., 540.01696548]) * units.hPa,
303
      np.array([7.93049516, 0., -8.33333333, -13.14758845]) * units.degC))
304
])
305
def test_get_layer(pressure, variable, heights, bottom, depth, interp, expected):
306
    """Tests get_layer functionality."""
307
    p_layer, y_layer = get_layer(pressure, variable, heights=heights, bottom=bottom,
308
                                 depth=depth, interpolate=interp)
309
    assert_array_almost_equal(p_layer, expected[0], 5)
310
    assert_array_almost_equal(y_layer, expected[1], 5)
311
312
313
def test_log_interp_2d():
314
    """Test interpolating with log x-scale in 2 dimensions."""
315
    x_log = np.array([[1e3, 1e4, 1e5, 1e6], [1e3, 1e4, 1e5, 1e6]])
316
    y_log = np.log(x_log) * 2 + 3
317
    x_interp = np.array([5e3, 5e4, 5e5])
318
    y_interp_truth = np.array([20.0343863828, 24.6395565688, 29.2447267548])
319
    y_interp = log_interp(x_interp, x_log, y_log, axis=1)
320
    assert_array_almost_equal(y_interp[1], y_interp_truth, 7)
321
322
323 View Code Duplication
def test_log_interp_3d():
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
324
    """Test interpolating with log x-scale 3 dimensions along second axis."""
325
    x_log = np.ones((3, 4, 3)) * np.array([1e3, 1e4, 1e5, 1e6]).reshape(-1, 1)
326
    y_log = np.log(x_log) * 2 + 3
327
    x_interp = np.array([5e3, 5e4, 5e5])
328
    y_interp_truth = np.array([20.0343863828, 24.6395565688, 29.2447267548])
329
    y_interp = log_interp(x_interp, x_log, y_log, axis=1)
330
    assert_array_almost_equal(y_interp[0, :, 0], y_interp_truth, 7)
331
332
333 View Code Duplication
def test_log_interp_4d():
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
334
    """Test interpolating with log x-scale 4 dimensions."""
335
    x_log = np.ones((2, 2, 3, 4)) * np.array([1e3, 1e4, 1e5, 1e6])
336
    y_log = np.log(x_log) * 2 + 3
337
    x_interp = np.array([5e3, 5e4, 5e5])
338
    y_interp_truth = np.array([20.0343863828, 24.6395565688, 29.2447267548])
339
    y_interp = log_interp(x_interp, x_log, y_log, axis=3)
340
    assert_array_almost_equal(y_interp[0, 0, 0, :], y_interp_truth, 7)
341
342
343 View Code Duplication
def test_log_interp_2args():
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
344
    """Test interpolating with log x-scale with 2 arguments."""
345
    x_log = np.array([1e3, 1e4, 1e5, 1e6])
346
    y_log = np.log(x_log) * 2 + 3
347
    y_log2 = np.log(x_log) * 2 + 3
348
    x_interp = np.array([5e3, 5e4, 5e5])
349
    y_interp_truth = np.array([20.0343863828, 24.6395565688, 29.2447267548])
350
    y_interp = log_interp(x_interp, x_log, y_log, y_log2)
351
    assert_array_almost_equal(y_interp[1], y_interp_truth, 7)
352
    assert_array_almost_equal(y_interp[0], y_interp_truth, 7)
353
354
355
def test_log_interp_set_nan_above():
356
    """Test interpolating with log x-scale setting out of bounds above data to nan."""
357
    x_log = np.array([1e3, 1e4, 1e5, 1e6])
358
    y_log = np.log(x_log) * 2 + 3
359
    x_interp = np.array([1e7])
360
    y_interp_truth = np.nan
361
    with pytest.warns(Warning):
362
        y_interp = log_interp(x_interp, x_log, y_log)
363
    assert_array_almost_equal(y_interp, y_interp_truth, 7)
364
365
366
def test_log_interp_no_extrap():
367
    """Test interpolating with log x-scale setting out of bounds value error."""
368
    x_log = np.array([1e3, 1e4, 1e5, 1e6])
369
    y_log = np.log(x_log) * 2 + 3
370
    x_interp = np.array([1e7])
371
    with pytest.raises(ValueError):
372
        log_interp(x_interp, x_log, y_log, fill_value=None)
373
374
375
def test_log_interp_set_nan_below():
376
    """Test interpolating with log x-scale setting out of bounds below data to nan."""
377
    x_log = np.array([1e3, 1e4, 1e5, 1e6])
378
    y_log = np.log(x_log) * 2 + 3
379
    x_interp = 1e2
380
    y_interp_truth = np.nan
381
    with pytest.warns(Warning):
382
        y_interp = log_interp(x_interp, x_log, y_log)
383
    assert_array_almost_equal(y_interp, y_interp_truth, 7)
384
385
386
def test_interp_2args():
387
    """Test interpolation with 2 arguments."""
388
    x = np.array([1., 2., 3., 4.])
389
    y = x
390
    y2 = x
391
    x_interp = np.array([2.5000000, 3.5000000])
392
    y_interp_truth = np.array([2.5000000, 3.5000000])
393
    y_interp = interp(x_interp, x, y, y2)
394
    assert_array_almost_equal(y_interp[0], y_interp_truth, 7)
395
    assert_array_almost_equal(y_interp[1], y_interp_truth, 7)
396
397
398
def test_interp_decrease():
399
    """Test interpolation with decreasing interpolation points."""
400
    x = np.array([1., 2., 3., 4.])
401
    y = x
402
    x_interp = np.array([3.5000000, 2.5000000])
403
    y_interp_truth = np.array([3.5000000, 2.5000000])
404
    y_interp = interp(x_interp, x, y)
405
    assert_array_almost_equal(y_interp, y_interp_truth, 7)
406
407
408
def test_interp_decrease_xp():
409
    """Test interpolation with decreasing order."""
410
    x = np.array([4., 3., 2., 1.])
411
    y = x
412
    x_interp = np.array([3.5000000, 2.5000000])
413
    y_interp_truth = np.array([3.5000000, 2.5000000])
414
    y_interp = interp(x_interp, x, y)
415
    assert_array_almost_equal(y_interp, y_interp_truth, 7)
416
417
418
def test_interp_end_point():
419
    """Test interpolation with point at data endpoints."""
420
    x = np.array([1., 2., 3., 4.])
421
    y = x
422
    x_interp = np.array([1.0, 4.0])
423
    y_interp_truth = np.array([1.0, 4.0])
424
    y_interp = interp(x_interp, x, y)
425
    assert_array_almost_equal(y_interp, y_interp_truth, 7)
426