Completed
Pull Request — master (#468)
by
unknown
01:05
created

test_log_interp_2args()   A

Complexity

Conditions 1

Size

Total Lines 10

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