Completed
Pull Request — main (#59)
by Angeline
01:52 queued 32s
created

test_Apex.igrf_file()   A

Complexity

Conditions 1

Size

Total Lines 12
Code Lines 9

Duplication

Lines 0
Ratio 0 %

Importance

Changes 0
Metric Value
eloc 9
dl 0
loc 12
rs 9.95
c 0
b 0
f 0
cc 1
nop 0
1
# -*- coding: utf-8 -*-
2
3
from __future__ import division, absolute_import, unicode_literals
4
5
import datetime as dt
6
import itertools
7
import numpy as np
8
from numpy.testing import assert_allclose
9
import os
10
import pytest
11
12
from apexpy import fortranapex as fa
13
from apexpy import Apex, ApexHeightError, helpers
14
15
16
# ----------------------------------------------------------------------------
17
# NOTE: whenever function outputs are tested against hard-coded numbers, the
18
# test results (numbers) were obtained by running the code that is tested.
19
# Therefore these tests below only check that nothing changes when refactoring,
20
# etc., and not if the results are actually correct.
21
# -----------------------------------------------------------------------------
22
23
24
# ============================================================================
25
#  Test initiating the Apex class
26
# ============================================================================
27
28
29
def test_init_defaults():
30
    Apex()
31
32
33
def test_init_date_int():
34
    apex_out = Apex(date=2015)
35
    assert apex_out.year == 2015
36
37
38
def test_init_date_float():
39
    apex_out = Apex(date=2015.5)
40
    assert apex_out.year == 2015.5
41
42
43
def test_init_date():
44
    date = dt.date(2015, 1, 1)
45
    apex_out = Apex(date=date)
46
    assert apex_out.year == helpers.toYearFraction(date)
47
48
49
def test_init_datetime():
50
    datetime = dt.datetime(2015, 6, 1, 18, 23, 45)
51
    apex_out = Apex(date=datetime)
52
    assert apex_out.year == helpers.toYearFraction(datetime)
53
54
55
def test_init_datafile_IOError():
56
    with pytest.raises(IOError):
57
        Apex(date=2015, datafile='foo/path/to/datafile.blah')
58
59
60
# ============================================================================
61
#  Test the low-level interfaces to the fortran wrappers
62
# ============================================================================
63
64
def test__geo2qd_scalar():
65
    apex_out = Apex(date=2000, refh=300)
66
    for lat in [0, 30, 60, 89]:
67
        for lon in [-179, -90, 0, 90, 180]:
68
            assert_allclose(apex_out._geo2qd(lat, lon, 100),
69
                            fa.apxg2q(lat, lon, 100, 0)[:2])
70
71
72
def test__geo2qd_array():
73
    apex_out = Apex(date=2000, refh=300)
74
    lats, lons = apex_out._geo2qd([[0, 30], [60, 90]], 15,
75
                                  [[100, 200], [300, 400]])
76
    lat1, lon1 = fa.apxg2q(0, 15, 100, 0)[:2]
77
    lat2, lon2 = fa.apxg2q(30, 15, 200, 0)[:2]
78
    lat3, lon3 = fa.apxg2q(60, 15, 300, 0)[:2]
79
    lat4, lon4 = fa.apxg2q(90, 15, 400, 0)[:2]
80
    assert_allclose(lats.astype(float), np.array([[lat1, lat2], [lat3, lat4]],
81
                                                 dtype=float))
82
    assert_allclose(lons.astype(float), np.array([[lon1, lon2], [lon3, lon4]],
83
                                                 dtype=float))
84
85
86
def test__geo2qd_longitude():
87
    apex_out = Apex(date=2000, refh=300)
88
    assert_allclose(apex_out._geo2qd(60, 180, 100),
89
                    fa.apxg2q(60, 180, 100, 0)[:2])
90
    assert_allclose(apex_out._geo2qd(60, -180, 100),
91
                    fa.apxg2q(60, -180, 100, 0)[:2])
92
    assert_allclose(apex_out._geo2qd(60, -180, 100),
93
                    apex_out._geo2qd(60, 180, 100))
94
    for i in range(-5, 5):
95
        for lat in [0, 30, 60, 90]:
96
            assert_allclose(apex_out._geo2qd(lat, 15 + i * 360, 100),
97
                            fa.apxg2q(lat, 15, 100, 0)[:2])
98
99
100
def test__geo2apex_scalar():
101
    apex_out = Apex(date=2000, refh=300)
102
    for lat in [0, 30, 60, 89]:
103
        for lon in [-179, -90, 0, 90, 180]:
104
            assert_allclose(apex_out._geo2apex(lat, lon, 100),
105
                            fa.apxg2all(lat, lon, 100, 300, 0)[2:4])
106
107
108
def test__geo2apex_array():
109
    apex_out = Apex(date=2000, refh=300)
110
    lats, lons = apex_out._geo2apex([[0, 30], [60, 90]], 15,
111
                                    [[100, 200], [300, 400]])
112
    lat1, lon1 = fa.apxg2all(0, 15, 100, 300, 0)[2:4]
113
    lat2, lon2 = fa.apxg2all(30, 15, 200, 300, 0)[2:4]
114
    lat3, lon3 = fa.apxg2all(60, 15, 300, 300, 0)[2:4]
115
    lat4, lon4 = fa.apxg2all(90, 15, 400, 300, 0)[2:4]
116
    assert_allclose(lats.astype(float), np.array([[lat1, lat2], [lat3, lat4]],
117
                                                 dtype=float))
118
    assert_allclose(lons.astype(float), np.array([[lon1, lon2], [lon3, lon4]],
119
                                                 dtype=float))
120
121
122
def test__geo2apex_longitude():
123
    apex_out = Apex(date=2000, refh=300)
124
    assert_allclose(apex_out._geo2apex(60, 180, 100),
125
                    fa.apxg2all(60, 180, 100, 300, 0)[2:4])
126
    assert_allclose(apex_out._geo2apex(60, -180, 100),
127
                    fa.apxg2all(60, -180, 100, 300, 0)[2:4])
128
    assert_allclose(apex_out._geo2apex(60, -180, 100),
129
                    apex_out._geo2apex(60, 180, 100))
130
    for i in range(-5, 5):
131
        for lat in [0, 30, 60, 90]:
132
            assert_allclose(apex_out._geo2apex(lat, 15 + i * 360, 100),
133
                            fa.apxg2all(lat, 15, 100, 300, 0)[2:4])
134
135
136
def test__geo2apexall_scalar():
137
    apex_out = Apex(date=2000, refh=300)
138
    for lat in [0, 30, 60, 89]:
139
        for lon in [-179, -90, 0, 90, 180]:
140
            ret1 = apex_out._geo2apexall(lat, lon, 100)
141
            ret2 = fa.apxg2all(lat, lon, 100, 300, 1)
142
            for r1, r2 in zip(ret1, ret2):
143
                assert_allclose(r1, r2)
144
145
146
def test__geo2apexall_array():
147
    apex_out = Apex(date=2000, refh=300)
148
    ret = apex_out._geo2apexall([[0, 30], [60, 90]], 15,
149
                                [[100, 200], [300, 400]])
150
    ret1 = fa.apxg2all(0, 15, 100, 300, 1)
151
    ret2 = fa.apxg2all(30, 15, 200, 300, 1)
152
    ret3 = fa.apxg2all(60, 15, 300, 300, 1)
153
    ret4 = fa.apxg2all(90, 15, 400, 300, 1)
154
    for i in range(len(ret)):
155
        try:
156
            # ret[i] is array of floats
157
            assert_allclose(ret[i].astype(float),
158
                            np.array([[ret1[i], ret2[i]], [ret3[i], ret4[i]]],
159
                                     dtype=float))
160
        except ValueError:
161
            # ret[i] is array of arrays
162
            assert_allclose(ret[i][0, 0], ret1[i])
163
            assert_allclose(ret[i][0, 1], ret2[i])
164
            assert_allclose(ret[i][1, 0], ret3[i])
165
            assert_allclose(ret[i][1, 1], ret4[i])
166
167
168
def test__qd2geo_scalar():
169
    apex_out = Apex(date=2000, refh=300)
170
    for lat in [0, 30, 60, 89]:
171
        for lon in [-179, -90, 0, 90, 180]:
172
            for prec in [-1, 1e-2, 1e-10]:
173
                assert_allclose(apex_out._qd2geo(lat, lon, 100, prec),
174
                                fa.apxq2g(lat, lon, 100, prec))
175
176
177
def test__qd2geo_array():
178
    apex_out = Apex(date=2000, refh=300)
179
    lats, lons, errs = apex_out._qd2geo([[0, 30], [60, 90]], 15,
180
                                        [[100, 200], [300, 400]], 1e-2)
181
    lat1, lon1, err1 = fa.apxq2g(0, 15, 100, 1e-2)
182
    lat2, lon2, err2 = fa.apxq2g(30, 15, 200, 1e-2)
183
    lat3, lon3, err3 = fa.apxq2g(60, 15, 300, 1e-2)
184
    lat4, lon4, err4 = fa.apxq2g(90, 15, 400, 1e-2)
185
    assert_allclose(lats.astype(float), np.array([[lat1, lat2], [lat3, lat4]],
186
                                                 dtype=float))
187
    assert_allclose(lons.astype(float), np.array([[lon1, lon2], [lon3, lon4]],
188
                                                 dtype=float))
189
    assert_allclose(errs.astype(float), np.array([[err1, err2], [err3, err4]],
190
                                                 dtype=float))
191
192
193
def test__qd2geo_longitude():
194
    apex_out = Apex(date=2000, refh=300)
195
    assert_allclose(apex_out._qd2geo(60, 180, 100, 1e-2),
196
                    fa.apxq2g(60, 180, 100, 1e-2))
197
    assert_allclose(apex_out._qd2geo(60, -180, 100, 1e-2),
198
                    fa.apxq2g(60, -180, 100, 1e-2))
199
    assert_allclose(apex_out._qd2geo(60, -180, 100, 1e-2),
200
                    apex_out._qd2geo(60, 180, 100, 1e-2))
201
    for i in range(-5, 5):
202
        for lat in [0, 30, 60, 90]:
203
            assert_allclose(apex_out._qd2geo(lat, 15 + i * 360, 100, 1e-2),
204
                            fa.apxq2g(lat, 15, 100, 1e-2))
205
206
207
def test__basevec_scalar():
208
    apex_out = Apex(date=2000, refh=300)
209
    for lat in [0, 30, 60, 89]:
210
        for lon in [-179, -90, 0, 90, 180]:
211
            assert_allclose(apex_out._basevec(lat, lon, 100),
212
                            fa.apxg2q(lat, lon, 100, 1)[2:4])
213
214
215
def test__basevec_array():
216
    apex_out = Apex(date=2000, refh=300)
217
    f1s, f2s = apex_out._basevec([[0, 30], [60, 90]], 15,
218
                                 [[100, 200], [300, 400]])
219
    f11, f21 = fa.apxg2q(0, 15, 100, 1)[2:4]
220
    f12, f22 = fa.apxg2q(30, 15, 200, 1)[2:4]
221
    f13, f23 = fa.apxg2q(60, 15, 300, 1)[2:4]
222
    f14, f24 = fa.apxg2q(90, 15, 400, 1)[2:4]
223
    assert_allclose(f1s[0, 0], f11)
224
    assert_allclose(f1s[0, 1], f12)
225
    assert_allclose(f1s[1, 0], f13)
226
    assert_allclose(f1s[1, 1], f14)
227
    assert_allclose(f2s[0, 0], f21)
228
    assert_allclose(f2s[0, 1], f22)
229
    assert_allclose(f2s[1, 0], f23)
230
    assert_allclose(f2s[1, 1], f24)
231
232
233
def test__basevec_longitude():
234
    apex_out = Apex(date=2000, refh=300)
235
    assert_allclose(apex_out._basevec(60, 180, 100),
236
                    fa.apxg2q(60, 180, 100, 1)[2:4])
237
    assert_allclose(apex_out._basevec(60, -180, 100),
238
                    fa.apxg2q(60, -180, 100, 1)[2:4])
239
    assert_allclose(apex_out._basevec(60, -180, 100),
240
                    apex_out._basevec(60, 180, 100))
241
    for i in range(-5, 5):
242
        for lat in [0, 30, 60, 90]:
243
            assert_allclose(apex_out._basevec(lat, 15 + i * 360, 100),
244
                            fa.apxg2q(lat, 15, 100, 1)[2:4])
245
246
247
# ============================================================================
248
#  Test the convert() method
249
# ============================================================================
250
251
252
def test_convert_geo2apex():
253
    apex_out = Apex(date=2000, refh=300)
254
    assert_allclose(apex_out.convert(60, 15, 'geo', 'apex', height=100),
255
                    apex_out.geo2apex(60, 15, 100))
256
257
258
def test_convert_geo2qd():
259
    apex_out = Apex(date=2000, refh=300)
260
    assert_allclose(apex_out.convert(60, 15, 'geo', 'qd', height=100),
261
                    apex_out.geo2qd(60, 15, 100))
262
263
264
def test_convert_geo2mlt_nodate():
265
    apex_out = Apex(date=2000, refh=300)
266
    with pytest.raises(ValueError):
267
        apex_out.convert(60, 15, 'geo', 'mlt')
268
269
270
def test_convert_geo2mlt():
271
    datetime = dt.datetime(2000, 3, 9, 14, 25, 58)
272
    apex_out = Apex(date=2000, refh=300)
273
    assert_allclose(apex_out.convert(60, 15, 'geo', 'mlt', height=100,
274
                                     ssheight=2e5, datetime=datetime)[1],
275
                    apex_out.mlon2mlt(apex_out.geo2apex(60, 15, 100)[1],
276
                                      datetime, ssheight=2e5))
277
278
279
def test_convert_apex2geo():
280
    apex_out = Apex(date=2000, refh=300)
281
    assert_allclose(apex_out.convert(60, 15, 'apex', 'geo', height=100,
282
                                     precision=1e-2),
283
                    apex_out.apex2geo(60, 15, 100, precision=1e-2)[:-1])
284
285
286
def test_convert_apex2qd():
287
    apex_out = Apex(date=2000, refh=300)
288
    assert_allclose(apex_out.convert(60, 15, 'apex', 'qd', height=100),
289
                    apex_out.apex2qd(60, 15, height=100))
290
291
292
def test_convert_apex2mlt():
293
    datetime = dt.datetime(2000, 3, 9, 14, 25, 58)
294
    apex_out = Apex(date=2000, refh=300)
295
    assert_allclose(apex_out.convert(60, 15, 'apex', 'mlt', height=100,
296
                                     datetime=datetime, ssheight=2e5)[1],
297
                    apex_out.mlon2mlt(15, datetime, ssheight=2e5))
298
299
300
def test_convert_qd2geo():
301
    apex_out = Apex(date=2000, refh=300)
302
    assert_allclose(apex_out.convert(60, 15, 'qd', 'geo', height=100,
303
                                     precision=1e-2),
304
                    apex_out.qd2geo(60, 15, 100, precision=1e-2)[:-1])
305
306
307
def test_convert_qd2apex():
308
    apex_out = Apex(date=2000, refh=300)
309
    assert_allclose(apex_out.convert(60, 15, 'qd', 'apex', height=100),
310
                    apex_out.qd2apex(60, 15, height=100))
311
312
313
def test_convert_qd2apex_at_equator():
314
    """Test the quasi-dipole to apex conversion at the magnetic equator."""
315
    apex_out = Apex(date=2000, refh=80)
316
    elat, elon = apex_out.convert(lat=0.0, lon=0, source='qd', dest='apex',
317
                                  height=120.0)
318
    clat, clon = apex_out.convert(lat=0.001, lon=0, source='qd', dest='apex',
319
                                  height=120.0)
320
    assert_allclose([elat, elon], [clat, clon], atol=1e-4)
321
322
323
def test_convert_qd2mlt():
324
    datetime = dt.datetime(2000, 3, 9, 14, 25, 58)
325
    apex_out = Apex(date=2000, refh=300)
326
    assert_allclose(apex_out.convert(60, 15, 'qd', 'mlt', height=100,
327
                                     datetime=datetime, ssheight=2e5)[1],
328
                    apex_out.mlon2mlt(15, datetime, ssheight=2e5))
329
330
331
def test_convert_mlt2geo():
332
    datetime = dt.datetime(2000, 3, 9, 14, 25, 58)
333
    apex_out = Apex(date=2000, refh=300)
334
    assert_allclose(apex_out.convert(60, 15, 'mlt', 'geo', height=100,
335
                                     datetime=datetime, precision=1e-2,
336
                                     ssheight=2e5),
337
                    apex_out.apex2geo(60, apex_out.mlt2mlon(15, datetime,
338
                                                            ssheight=2e5), 100,
339
                                      precision=1e-2)[:-1])
340
341
342
def test_convert_mlt2apex():
343
    datetime = dt.datetime(2000, 3, 9, 14, 25, 58)
344
    apex_out = Apex(date=2000, refh=300)
345
    assert_allclose(apex_out.convert(60, 15, 'mlt', 'apex', height=100,
346
                                     datetime=datetime, ssheight=2e5),
347
                    (60, apex_out.mlt2mlon(15, datetime, ssheight=2e5)))
348
349
350
def test_convert_mlt2qd():
351
    datetime = dt.datetime(2000, 3, 9, 14, 25, 58)
352
    apex_out = Apex(date=2000, refh=300)
353
    assert_allclose(apex_out.convert(60, 15, 'mlt', 'qd', height=100,
354
                                     datetime=datetime, ssheight=2e5),
355
                    apex_out.apex2qd(60, apex_out.mlt2mlon(15, datetime,
356
                                                           ssheight=2e5),
357
                                     height=100))
358
359
360
def test_convert_invalid_lat():
361
    apex_out = Apex(date=2000, refh=300)
362
    with pytest.raises(ValueError):
363
        apex_out.convert(91, 0, 'geo', 'geo')
364
    with pytest.raises(ValueError):
365
        apex_out.convert(-91, 0, 'geo', 'geo')
366
    apex_out.convert(90, 0, 'geo', 'geo')
367
    apex_out.convert(-90, 0, 'geo', 'geo')
368
369
    assert_allclose(apex_out.convert(90 + 1e-5, 0, 'geo', 'apex'),
370
                    apex_out.convert(90, 0, 'geo', 'apex'), rtol=0, atol=1e-8)
371
372
373
def test_convert_invalid_transformation():
374
    apex_out = Apex(date=2000, refh=300)
375
    with pytest.raises(NotImplementedError):
376
        apex_out.convert(0, 0, 'foobar', 'geo')
377
    with pytest.raises(NotImplementedError):
378
        apex_out.convert(0, 0, 'geo', 'foobar')
379
380
381
coord_names = ['geo', 'apex', 'qd']
382
383
384
@pytest.mark.parametrize('transform', itertools.product(coord_names,
385
                                                        coord_names))
386
def test_convert_withnan(transform):
387
    """Test Apex.convert success with NaN input."""
388
    num_nans = 5
389
    in_lat = np.arange(0, 10, dtype=float)
390
    in_lat[:num_nans] = np.nan
391
    in_lon = np.arange(0, 10, dtype=float)
392
    in_lon[:num_nans] = np.nan
393
    src, dest = transform
394
    apex_out = Apex(date=2000, refh=80)
395
    out_lat, out_lon = apex_out.convert(in_lat, in_lon, src, dest, height=120)
396
    assert np.all(np.isnan(out_lat[:num_nans]))
397
    assert np.all(np.isnan(out_lon[:num_nans]))
398
    assert np.all(np.isfinite(out_lat[num_nans:]))
399
    assert np.all(np.isfinite(out_lat[num_nans:]))
400
401
402
# ============================================================================
403
#  Test the geo2apex() method
404
# ============================================================================
405
406
407
def test_geo2apex():
408
    apex_out = Apex(date=2000, refh=300)
409
    lat, lon = apex_out.geo2apex(60, 15, 100)
410
    assert_allclose((lat, lon), apex_out._geo2apex(60, 15, 100))
411
    assert type(lat) != np.ndarray
412
    assert type(lon) != np.ndarray
413
414
415
def test_geo2apex_vectorization():
416
    apex_out = Apex(date=2000, refh=300)
417
    assert apex_out.geo2apex([60, 60], 15, 100)[0].shape == (2,)
418
    assert apex_out.geo2apex(60, [15, 15], 100)[0].shape == (2,)
419
    assert apex_out.geo2apex(60, 15, [100, 100])[0].shape == (2,)
420
421
422
def test_geo2apex_invalid_lat():
423
    apex_out = Apex(date=2000, refh=300)
424
    with pytest.raises(ValueError):
425
        apex_out.geo2apex(91, 0, 0)
426
    with pytest.raises(ValueError):
427
        apex_out.geo2apex(-91, 0, 0)
428
    apex_out.geo2apex(90, 0, 0)
429
    apex_out.geo2apex(-90, 0, 0)
430
431
    assert_allclose(apex_out.geo2apex(90 + 1e-5, 0, 0),
432
                    apex_out.geo2apex(90, 0, 0), rtol=0, atol=1e-8)
433
434
435
def test_geo2apex_undefined_warning(recwarn):
436
    """Test warning and fill values for an undefined location."""
437
    apex_out = Apex(date=2000, refh=10000)
438
    ret = apex_out.geo2apex(0, 0, 0)
439
440
    assert np.isnan(ret[0])
441
    assert len(recwarn) == 1
442
    assert issubclass(recwarn[-1].category, UserWarning)
443
    assert 'set to NaN where' in str(recwarn[-1].message)
444
445
446
# ============================================================================
447
#  Test the apex2geo() method
448
# ============================================================================
449
450
451
def test_apex2geo():
452
    apex_out = Apex(date=2000, refh=300)
453
    lat, lon, error = apex_out.apex2geo(60, 15, 100, precision=1e-2)
454
    assert_allclose((lat, lon, error),
455
                    apex_out.qd2geo(*apex_out.apex2qd(60, 15, 100), height=100,
456
                                    precision=1e-2))
457
458
    assert type(lat) != np.ndarray
459
    assert type(lon) != np.ndarray
460
    assert type(error) != np.ndarray
461
462
463
def test_apex2geo_vectorization():
464
    apex_out = Apex(date=2000, refh=300)
465
    assert apex_out.apex2geo([60, 60], 15, 100)[0].shape == (2,)
466
    assert apex_out.apex2geo(60, [15, 15], 100)[0].shape == (2,)
467
    assert apex_out.apex2geo(60, 15, [100, 100])[0].shape == (2,)
468
469
470
def test_apex2geo_invalid_lat():
471
    apex_out = Apex(date=2000, refh=300)
472
    with pytest.raises(ValueError):
473
        apex_out.apex2geo(91, 0, 0, 1e-2)
474
    with pytest.raises(ValueError):
475
        apex_out.apex2geo(-91, 0, 0, 1e-2)
476
    apex_out.apex2geo(90, 0, 0, 1e-2)
477
    apex_out.apex2geo(-90, 0, 0, 1e-2)
478
479
    assert_allclose(apex_out.apex2geo(90 + 1e-5, 0, 0, 1e-2),
480
                    apex_out.apex2geo(90, 0, 0, 1e-2), rtol=0, atol=1e-8)
481
482
483
# ============================================================================
484
#  Test the geo2qd() method
485
# ============================================================================
486
487
488
def test_geo2qd():
489
    apex_out = Apex(date=2000, refh=300)
490
    lat, lon = apex_out.geo2qd(60, 15, 100)
491
    assert_allclose((lat, lon), apex_out._geo2qd(60, 15, 100))
492
    assert type(lat) != np.ndarray
493
    assert type(lon) != np.ndarray
494
495
496
def test_geo2qd_vectorization():
497
    apex_out = Apex(date=2000, refh=300)
498
    assert apex_out.geo2qd([60, 60], 15, 100)[0].shape == (2,)
499
    assert apex_out.geo2qd(60, [15, 15], 100)[0].shape == (2,)
500
    assert apex_out.geo2qd(60, 15, [100, 100])[0].shape == (2,)
501
502
503
def test_geo2qd_invalid_lat():
504
    apex_out = Apex(date=2000, refh=300)
505
    with pytest.raises(ValueError):
506
        apex_out.geo2qd(91, 0, 0)
507
    with pytest.raises(ValueError):
508
        apex_out.geo2qd(-91, 0, 0)
509
    apex_out.geo2qd(90, 0, 0)
510
    apex_out.geo2qd(-90, 0, 0)
511
512
    assert_allclose(apex_out.geo2qd(90 + 1e-5, 0, 0), apex_out.geo2qd(90, 0, 0),
513
                    rtol=0, atol=1e-8)
514
515
516
# ============================================================================
517
#  Test the qd2geo() method
518
# ============================================================================
519
520
521
def test_qd2geo():
522
    apex_out = Apex(date=2000, refh=300)
523
    lat, lon, error = apex_out.qd2geo(60, 15, 100, precision=1e-2)
524
    assert_allclose((lat, lon, error), apex_out._qd2geo(60, 15, 100, 1e-2))
525
    assert type(lat) != np.ndarray
526
    assert type(lon) != np.ndarray
527
    assert type(error) != np.ndarray
528
529
530
def test_qd2geo_vectorization():
531
    apex_out = Apex(date=2000, refh=300)
532
    assert apex_out.qd2geo([60, 60], 15, 100)[0].shape == (2,)
533
    assert apex_out.qd2geo(60, [15, 15], 100)[0].shape == (2,)
534
    assert apex_out.qd2geo(60, 15, [100, 100])[0].shape == (2,)
535
536
537
def test_qd2geo_invalid_lat():
538
    apex_out = Apex(date=2000, refh=300)
539
    with pytest.raises(ValueError):
540
        apex_out.qd2geo(91, 0, 0, precision=1e-2)
541
    with pytest.raises(ValueError):
542
        apex_out.qd2geo(-91, 0, 0, precision=1e-2)
543
    apex_out.qd2geo(90, 0, 0, precision=1e-2)
544
    apex_out.qd2geo(-90, 0, 0, precision=1e-2)
545
546
    assert_allclose(apex_out.qd2geo(90 + 1e-5, 0, 0, 1e-2),
547
                    apex_out.qd2geo(90, 0, 0, 1e-2), rtol=0, atol=1e-8)
548
549
550
# ============================================================================
551
#  Test the apex2qd() method
552
# ============================================================================
553
554
555
def test_apex2qd():
556
    apex_out = Apex(date=2000, refh=300)
557
    lat, lon = apex_out.apex2qd(60, 15, 100)
558
    assert_allclose((lat, lon),
559
                    [60.498401, 15])
560
    assert type(lat) != np.ndarray
561
    assert type(lon) != np.ndarray
562
563
564
def test_apex2qd_vectorization():
565
    apex_out = Apex(date=2000, refh=300)
566
    assert apex_out.apex2qd([60, 60], 15, 100)[0].shape == (2,)
567
    assert apex_out.apex2qd(60, [15, 15], 100)[0].shape == (2,)
568
    assert apex_out.apex2qd(60, 15, [100, 100])[0].shape == (2,)
569
570
571
def test_apex2qd_invalid_lat():
572
    apex_out = Apex(date=2000, refh=300)
573
    with pytest.raises(ValueError):
574
        apex_out.apex2qd(91, 0, 0)
575
    with pytest.raises(ValueError):
576
        apex_out.apex2qd(-91, 0, 0)
577
    apex_out.apex2qd(90, 0, 0)
578
    apex_out.apex2qd(-90, 0, 0)
579
580
    assert_allclose(apex_out.apex2qd(90 + 1e-5, 0, 0),
581
                    apex_out.apex2qd(90, 0, 0), rtol=0, atol=1e-8)
582
583
584
def test_apex2qd_apexheight_close():
585
    apex_out = Apex(date=2000, refh=300)
586
    apex_out.apex2qd(0, 15, 300 + 1e-6)
587
588
589
def test_apex2qd_apexheight_over():
590
    apex_out = Apex(date=2000, refh=300)
591
    with pytest.raises(ApexHeightError):
592
        apex_out.apex2qd(0, 15, 301)
593
594
595
# ============================================================================
596
#  Test the qd2apex() method
597
# ============================================================================
598
599
600
def test_qd2apex():
601
    apex_out = Apex(date=2000, refh=300)
602
    lat, lon = apex_out.qd2apex(60, 15, 100)
603
    assert_allclose((lat, lon),
604
                    [59.491381, 15])
605
    assert type(lat) != np.ndarray
606
    assert type(lon) != np.ndarray
607
608
609
def test_qd2apex_vectorization():
610
    apex_out = Apex(date=2000, refh=300)
611
    assert apex_out.qd2apex([60, 60], 15, 100)[0].shape == (2,)
612
    assert apex_out.qd2apex(60, [15, 15], 100)[0].shape == (2,)
613
    assert apex_out.qd2apex(60, 15, [100, 100])[0].shape == (2,)
614
615
616
def test_qd2apex_invalid_lat():
617
    apex_out = Apex(date=2000, refh=300)
618
    with pytest.raises(ValueError):
619
        apex_out.qd2apex(91, 0, 0)
620
    with pytest.raises(ValueError):
621
        apex_out.qd2apex(-91, 0, 0)
622
    apex_out.qd2apex(90, 0, 0)
623
    apex_out.qd2apex(-90, 0, 0)
624
625
    assert_allclose(apex_out.qd2apex(90 + 1e-5, 0, 0),
626
                    apex_out.qd2apex(90, 0, 0), rtol=0, atol=1e-8)
627
628
629
def test_qd2apex_apexheight_close():
630
    apex_out = Apex(date=2000, refh=300)
631
    assert_allclose(apex_out.qd2apex(0, 15, 300 - 1e-5),
632
                    apex_out.qd2apex(0, 15, 300))
633
634
635
def test_qd2apex_apexheight_over():
636
    apex_out = Apex(date=2000, refh=300)
637
    with pytest.raises(ApexHeightError):
638
        apex_out.qd2apex(0, 15, 299)
639
640
641
# ============================================================================
642
#  Test mlon2mlt()
643
# ============================================================================
644
645
646
def test_mlon2mlt_scalar():
647
    apex_out = Apex(date=2000, refh=300)
648
    mlon = apex_out.mlon2mlt(0, dt.datetime(2000, 2, 3, 4, 5, 6))
649
    assert_allclose(mlon, 23.019629923502603)
650
    assert type(mlon) != np.ndarray
651
652
653
def test_mlon2mlt_ssheight():
654
    apex_out = Apex(date=2000, refh=300)
655
    mlt = apex_out.mlon2mlt(0, dt.datetime(2000, 2, 3, 4, 5, 6),
656
                            ssheight=50 * 2000)
657
    assert_allclose(mlt, 23.026712036132814)
658
659
660
def test_mlon2mlt_1Darray():
661
    apex_out = Apex(date=2000, refh=300)
662
    assert_allclose(apex_out.mlon2mlt([0, 180],
663
                                      dt.datetime(2000, 2, 3, 4, 5, 6)),
664
                    [23.019261, 11.019261], rtol=1e-4)
665
666
667
def test_mlon2mlt_2Darray():
668
    apex_out = Apex(date=2000, refh=300)
669
    assert_allclose(apex_out.mlon2mlt([[0, 180], [0, 180]],
670
                                      dt.datetime(2000, 2, 3, 4, 5, 6)),
671
                    [[23.019261, 11.019261], [23.019261, 11.019261]], rtol=1e-4)
672
673
674
def test_mlon2mlt_diffdates():
675
    apex_out = Apex(date=2000, refh=300)
676
    dtime1 = dt.datetime(2000, 2, 3, 4, 5, 6)
677
    dtime2 = dt.datetime(2000, 2, 3, 5, 5, 6)
678
    assert apex_out.mlon2mlt(0, dtime1) != apex_out.mlon2mlt(0, dtime2)
679
680
681
def test_mlon2mlt_offset():
682
    apex_out = Apex(date=2000, refh=300)
683
    date = dt.datetime(2000, 2, 3, 4, 5, 6)
684
    assert_allclose(apex_out.mlon2mlt(0, date),
685
                    apex_out.mlon2mlt(-15, date) + 1)
686
    assert_allclose(apex_out.mlon2mlt(0, date),
687
                    apex_out.mlon2mlt(-10 * 15, date) + 10)
688
689
690
def test_mlon2mlt_range():
691
    apex_out = Apex(date=2000, refh=300)
692
    assert_allclose(apex_out.mlon2mlt(range(0, 361, 30),
693
                                      dt.datetime(2000, 2, 3, 4, 5, 6)),
694
                    [23.01963, 1.01963, 3.01963, 5.01963, 7.01963,
695
                     9.01963, 11.01963, 13.01963, 15.01963, 17.01963,
696
                     19.01963, 21.01963, 23.01963],
697
                    rtol=1e-4)
698
699
700
# ============================================================================
701
#  Test mlt2mlon()
702
# ============================================================================
703
704
705
def test_mlt2mlon_scalar():
706
    apex_out = Apex(date=2000, refh=300)
707
    mlt = apex_out.mlt2mlon(0, dt.datetime(2000, 2, 3, 4, 5, 6))
708
    assert_allclose(mlt, 14.705551147460938)
709
    assert type(mlt) != np.ndarray
710
711
712
def test_mlt2mlon_ssheight():
713
    apex_out = Apex(date=2000, refh=300)
714
    mlt = apex_out.mlt2mlon(0, dt.datetime(2000, 2, 3, 4, 5, 6),
715
                            ssheight=50 * 2000)
716
    assert_allclose(mlt, 14.599319458007812)
717
718
719
def test_mlt2mlon_1Darray():
720
    apex_out = Apex(date=2000, refh=300)
721
    assert_allclose(apex_out.mlt2mlon([0, 12],
722
                                      dt.datetime(2000, 2, 3, 4, 5, 6)),
723
                    [14.705551, 194.705551], rtol=1e-4)
724
725
726
def test_mlt2mlon_2Darray():
727
    apex_out = Apex(date=2000, refh=300)
728
    assert_allclose(apex_out.mlt2mlon([[0, 12], [0, 12]],
729
                                      dt.datetime(2000, 2, 3, 4, 5, 6)),
730
                    [[14.705551, 194.705551], [14.705551, 194.705551]],
731
                    rtol=1e-4)
732
733
734
def test_mlt2mlon_diffdates():
735
    apex_out = Apex(date=2000, refh=300)
736
    dtime1 = dt.datetime(2000, 2, 3, 4, 5, 6)
737
    dtime2 = dt.datetime(2000, 2, 3, 5, 5, 6)
738
    assert apex_out.mlt2mlon(0, dtime1) != apex_out.mlt2mlon(0, dtime2)
739
740
741
def test_mlt2mlon_offset():
742
    apex_out = Apex(date=2000, refh=300)
743
    date = dt.datetime(2000, 2, 3, 4, 5, 6)
744
    assert_allclose(apex_out.mlt2mlon(0, date), apex_out.mlt2mlon(1, date) - 15)
745
    assert_allclose(apex_out.mlt2mlon(0, date),
746
                    apex_out.mlt2mlon(10, date) - 150)
747
748
749
def test_mlt2mlon_range():
750
    apex_out = Apex(date=2000, refh=300)
751
    assert_allclose(apex_out.mlt2mlon(range(0, 25, 2),
752
                                      dt.datetime(2000, 2, 3, 4, 5, 6)),
753
                    [14.705551, 44.705551, 74.705551, 104.705551, 134.705551,
754
                     164.705551, 194.705551, 224.705551, 254.705551, 284.705551,
755
                     314.705551, 344.705551, 14.705551],
756
                    rtol=1e-4)
757
758
759
# ============================================================================
760
#  Test mlt/mlon back and forth
761
# ============================================================================
762
763
764
def test_mlon2mlt2mlon():
765
    apex_out = Apex(date=2000, refh=300)
766
    date = dt.datetime(2000, 2, 3, 4, 5, 6)
767
    assert_allclose(apex_out.mlon2mlt(apex_out.mlt2mlon(0, date), date), 0)
768
    assert_allclose(apex_out.mlon2mlt(apex_out.mlt2mlon(6, date), date), 6)
769
    assert_allclose(apex_out.mlon2mlt(apex_out.mlt2mlon(12, date), date), 12)
770
    assert_allclose(apex_out.mlon2mlt(apex_out.mlt2mlon(18, date), date), 18)
771
    assert_allclose(apex_out.mlon2mlt(apex_out.mlt2mlon(24, date), date), 0)
772
773
774
def test_mlt2mlon2mlt():
775
    apex_out = Apex(date=2000, refh=300)
776
    date = dt.datetime(2000, 2, 3, 4, 5, 6)
777
    assert_allclose(apex_out.mlt2mlon(apex_out.mlon2mlt(0, date), date), 0)
778
    assert_allclose(apex_out.mlt2mlon(apex_out.mlon2mlt(90, date), date), 90)
779
    assert_allclose(apex_out.mlt2mlon(apex_out.mlon2mlt(180, date), date), 180)
780
    assert_allclose(apex_out.mlt2mlon(apex_out.mlon2mlt(270, date), date), 270)
781
    assert_allclose(apex_out.mlt2mlon(apex_out.mlon2mlt(360, date), date), 0)
782
783
784
# ============================================================================
785
#  Test the map_to_height() method
786
# ============================================================================
787
788
789
def test_map_to_height():
790
    apex_out = Apex(date=2000, refh=300)
791
    assert_allclose(apex_out.map_to_height(60, 15, 100, 10000, conjugate=False,
792
                                           precision=1e-10),
793
                    (31.841459274291992, 17.916629791259766, 0))
794
    assert_allclose(apex_out.map_to_height(30, 170, 100, 500, conjugate=False,
795
                                           precision=1e-2),
796
                    (25.727252960205078, 169.60546875, 0.00017655163537710905))
797
798
799
def test_map_to_height_same_height():
800
    apex_out = Apex(date=2000, refh=300)
801
    assert_allclose(apex_out.map_to_height(60, 15, 100, 100, conjugate=False,
802
                                           precision=1e-10),
803
                    (60, 15, 3.4150946248701075e-6), rtol=1e-5)
804
805
806
def test_map_to_height_conjugate():
807
    """Test results of map_to_height using conjugacy."""
808
    apex_out = Apex(date=2000, refh=300)
809
    assert_allclose(apex_out.map_to_height(60, 15, 100, 10000, conjugate=True,
810
                                           precision=1e-10),
811
                    (-25.424892425537109, 27.310417175292969,
812
                     1.2074182222931995e-6), atol=1e-6)
813
    assert_allclose(apex_out.map_to_height(30, 170, 100, 500, conjugate=True,
814
                                           precision=1e-2),
815
                    (-13.76642894744873, 164.24259948730469,
816
                     0.00056820799363777041), atol=1e-6)
817
818
819
def test_map_to_height_vectorization():
820
    apex_out = Apex(date=2000, refh=300)
821
    assert_allclose(apex_out.map_to_height([60, 60], 15, 100, 100),
822
                    ([60] * 2, [15] * 2, [3.4150946248701075e-6] * 2),
823
                    rtol=1e-5)
824
    assert_allclose(apex_out.map_to_height(60, [15, 15], 100, 100),
825
                    ([60] * 2, [15] * 2, [3.4150946248701075e-6] * 2),
826
                    rtol=1e-5)
827
    assert_allclose(apex_out.map_to_height(60, 15, [100, 100], 100),
828
                    ([60] * 2, [15] * 2, [3.4150946248701075e-6] * 2),
829
                    rtol=1e-5)
830
    assert_allclose(apex_out.map_to_height(60, 15, 100, [100, 100]),
831
                    ([60] * 2, [15] * 2, [3.4150946248701075e-6] * 2),
832
                    rtol=1e-5)
833
834
835
def test_map_to_height_ApexHeightError():
836
    apex_out = Apex(date=2000, refh=300)
837
    with pytest.raises(ApexHeightError):
838
        apex_out.map_to_height(0, 15, 100, 10000)
839
840
841
# ============================================================================
842
#  Test the map_E_to_height() method
843
# ============================================================================
844
845
846 View Code Duplication
def test_map_E_to_height():
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
847
    apex_out = Apex(date=2000, refh=300)
848
    out_60_15_100_500 = [0.7115211, 2.3562392, 0.57259707]
849
    out_60_15_100_500_234 = [1.560284, 3.439154, 0.782339]
850
    out_60_15_100_1000 = [0.677964, 2.089811, 0.558601]
851
    out_60_15_200_500 = [0.723773, 2.427366, 0.590826]
852
    out_60_30_100_500 = [0.686265, 2.375296, 0.600594]
853
    out_70_15_100_500 = [0.727605, 2.180817, 0.291414]
854
855
    # scalar
856
    assert_allclose(apex_out.map_E_to_height(60, 15, 100, 500, [1, 2, 3]),
857
                    out_60_15_100_500, rtol=1e-5)
858
    assert_allclose(apex_out.map_E_to_height(60, 15, 100, 500, [2, 3, 4]),
859
                    out_60_15_100_500_234, rtol=1e-5)
860
    assert_allclose(apex_out.map_E_to_height(60, 15, 100, 1000, [1, 2, 3]),
861
                    out_60_15_100_1000, rtol=1e-5)
862
    assert_allclose(apex_out.map_E_to_height(60, 15, 200, 500, [1, 2, 3]),
863
                    out_60_15_200_500, rtol=1e-5)
864
    assert_allclose(apex_out.map_E_to_height(60, 30, 100, 500, [1, 2, 3]),
865
                    out_60_30_100_500, rtol=1e-5)
866
    assert_allclose(apex_out.map_E_to_height(70, 15, 100, 500, [1, 2, 3]),
867
                    out_70_15_100_500, rtol=1e-5)
868
869
    # vectorize lat
870
    assert_allclose(apex_out.map_E_to_height([60, 70], 15, 100, 500,
871
                                             np.array([[1, 2, 3]] * 2).T),
872
                    np.array([out_60_15_100_500, out_70_15_100_500]).T,
873
                    rtol=1e-5)
874
875
    # vectorize lon
876
    assert_allclose(apex_out.map_E_to_height(60, [15, 30], 100, 500,
877
                                             np.array([[1, 2, 3]] * 2).T),
878
                    np.array([out_60_15_100_500, out_60_30_100_500]).T,
879
                    rtol=1e-5)
880
881
    # vectorize height
882
    assert_allclose(apex_out.map_E_to_height(60, 15, [100, 200], 500,
883
                                             np.array([[1, 2, 3]] * 2).T),
884
                    np.array([out_60_15_100_500, out_60_15_200_500]).T,
885
                    rtol=1e-5)
886
887
    # vectorize newheight
888
    assert_allclose(apex_out.map_E_to_height(60, 15, 100, [500, 1000],
889
                                             np.array([[1, 2, 3]] * 2).T),
890
                    np.array([out_60_15_100_500, out_60_15_100_1000]).T,
891
                    rtol=1e-5)
892
893
    # vectorize E
894
    assert_allclose(apex_out.map_E_to_height(60, 15, 100, 500,
895
                                             np.array([[1, 2, 3],
896
                                                       [2, 3, 4]]).T),
897
                    np.array([out_60_15_100_500, out_60_15_100_500_234]).T,
898
                    rtol=1e-5)
899
900
901
# ============================================================================
902
#  Test the map_V_to_height() method
903
# ============================================================================
904
905
906 View Code Duplication
def test_map_V_to_height():
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
907
    apex_out = Apex(date=2000, refh=300)
908
    out_60_15_100_500 = [0.819719, 2.845114, 0.695437]
909
    out_60_15_100_500_234 = [1.830277, 4.14345, 0.947624]
910
    out_60_15_100_1000 = [0.924577, 3.149964, 0.851343]
911
    out_60_15_200_500 = [0.803882, 2.793206, 0.682839]
912
    out_60_30_100_500 = [0.761412, 2.878837, 0.736549]
913
    out_70_15_100_500 = [0.846819, 2.592572, 0.347919]
914
915
    # scalar
916
    assert_allclose(apex_out.map_V_to_height(60, 15, 100, 500, [1, 2, 3]),
917
                    out_60_15_100_500, rtol=1e-5)
918
    assert_allclose(apex_out.map_V_to_height(60, 15, 100, 500, [2, 3, 4]),
919
                    out_60_15_100_500_234, rtol=1e-5)
920
    assert_allclose(apex_out.map_V_to_height(60, 15, 100, 1000, [1, 2, 3]),
921
                    out_60_15_100_1000, rtol=1e-5)
922
    assert_allclose(apex_out.map_V_to_height(60, 15, 200, 500, [1, 2, 3]),
923
                    out_60_15_200_500, rtol=1e-5)
924
    assert_allclose(apex_out.map_V_to_height(60, 30, 100, 500, [1, 2, 3]),
925
                    out_60_30_100_500, rtol=1e-5)
926
    assert_allclose(apex_out.map_V_to_height(70, 15, 100, 500, [1, 2, 3]),
927
                    out_70_15_100_500, rtol=1e-5)
928
929
    # vectorize lat
930
    assert_allclose(apex_out.map_V_to_height([60, 70], 15, 100, 500,
931
                                             np.array([[1, 2, 3]] * 2).T),
932
                    np.array([out_60_15_100_500, out_70_15_100_500]).T,
933
                    rtol=1e-5)
934
935
    # vectorize lon
936
    assert_allclose(apex_out.map_V_to_height(60, [15, 30], 100, 500,
937
                                             np.array([[1, 2, 3]] * 2).T),
938
                    np.array([out_60_15_100_500, out_60_30_100_500]).T,
939
                    rtol=1e-5)
940
941
    # vectorize height
942
    assert_allclose(apex_out.map_V_to_height(60, 15, [100, 200], 500,
943
                                             np.array([[1, 2, 3]] * 2).T),
944
                    np.array([out_60_15_100_500, out_60_15_200_500]).T,
945
                    rtol=1e-5)
946
947
    # vectorize newheight
948
    assert_allclose(apex_out.map_V_to_height(60, 15, 100, [500, 1000],
949
                                             np.array([[1, 2, 3]] * 2).T),
950
                    np.array([out_60_15_100_500, out_60_15_100_1000]).T,
951
                    rtol=1e-5)
952
953
    # vectorize E
954
    assert_allclose(apex_out.map_V_to_height(60, 15, 100, 500,
955
                                             np.array([[1, 2, 3],
956
                                                       [2, 3, 4]]).T),
957
                    np.array([out_60_15_100_500, out_60_15_100_500_234]).T,
958
                    rtol=1e-5)
959
960
961
# ============================================================================
962
#  Test basevectors_qd()
963
# ============================================================================
964
965
966
# test coords
967
968
def test_basevectors_qd_scalar_geo():
969
    apex_out = Apex(date=2000, refh=300)
970
    assert_allclose(apex_out.basevectors_qd(60, 15, 100, coords='geo'),
971
                    apex_out._basevec(60, 15, 100))
972
973
974
def test_basevectors_qd_scalar_apex():
975
    apex_out = Apex(date=2000, refh=300)
976
    glat, glon, _ = apex_out.apex2geo(60, 15, 100, precision=1e-2)
977
    assert_allclose(apex_out.basevectors_qd(60, 15, 100, coords='apex',
978
                                            precision=1e-2),
979
                    apex_out._basevec(glat, glon, 100))
980
981
982
def test_basevectors_qd_scalar_qd():
983
    apex_out = Apex(date=2000, refh=300)
984
    glat, glon, _ = apex_out.qd2geo(60, 15, 100, precision=1e-2)
985
    assert_allclose(apex_out.basevectors_qd(60, 15, 100, coords='qd',
986
                                            precision=1e-2),
987
                    apex_out._basevec(glat, glon, 100))
988
989
990
# test shapes and vectorization of arguments
991
def test_basevectors_qd_scalar_shape():
992
    apex_out = Apex(date=2000, refh=300)
993
    ret = apex_out.basevectors_qd(60, 15, 100)
994
    for r in ret:
995
        assert r.shape == (2,)
996
997
998
def test_basevectors_qd_vectorization():
999
    apex_out = Apex(date=2000, refh=300)
1000
    ret = apex_out.basevectors_qd([60, 60, 60, 60], 15, 100, coords='geo')
1001
    for r in ret:
1002
        assert r.shape == (2, 4)
1003
    ret = apex_out.basevectors_qd(60, [15, 15, 15, 15], 100, coords='geo')
1004
    for r in ret:
1005
        assert r.shape == (2, 4)
1006
    ret = apex_out.basevectors_qd(60, 15, [100, 100, 100, 100], coords='geo')
1007
    for r in ret:
1008
        assert r.shape == (2, 4)
1009
1010
1011
# test array return values
1012
1013
def test_basevectors_qd_array():
1014
    apex_out = Apex(date=2000, refh=300)
1015
    f1, f2 = apex_out.basevectors_qd([0, 30], 15, 100, coords='geo')
1016
    f1_lat0, f2_lat0 = apex_out._basevec(0, 15, 100)
1017
    f1_lat30, f2_lat30 = apex_out._basevec(30, 15, 100)
1018
    assert_allclose(f1[:, 0], f1_lat0)
1019
    assert_allclose(f2[:, 0], f2_lat0)
1020
    assert_allclose(f1[:, 1], f1_lat30)
1021
    assert_allclose(f2[:, 1], f2_lat30)
1022
1023
1024
# ============================================================================
1025
#  Test basevectors_apex()
1026
# ============================================================================
1027
1028
1029
# test against return from _geo2apexall for different coords
1030
1031
def test_basevectors_apex_scalar_geo():
1032
    apex_out = Apex(date=2000, refh=300)
1033
1034
    (f1, f2, f3, g1, g2, g3, d1, d2, d3, e1, e2,
1035
     e3) = apex_out.basevectors_apex(60, 15, 100, coords='geo')
1036
1037
    (_, _, _, _, f1_, f2_, _, d1_, d2_, d3_, _, e1_, e2_,
1038
     e3_) = apex_out._geo2apexall(60, 15, 100)
1039
1040
    assert_allclose(f1, f1_)
1041
    assert_allclose(f2, f2_)
1042
    assert_allclose(d1, d1_)
1043
    assert_allclose(d2, d2_)
1044
    assert_allclose(d3, d3_)
1045
    assert_allclose(e1, e1_)
1046
    assert_allclose(e2, e2_)
1047
    assert_allclose(e3, e3_)
1048
1049
1050 View Code Duplication
def test_basevectors_apex_scalar_apex():
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
1051
    apex_out = Apex(date=2000, refh=300)
1052
1053
    (f1, f2, f3, g1, g2, g3, d1, d2, d3, e1, e2,
1054
     e3) = apex_out.basevectors_apex(60, 15, 100, coords='apex', precision=1e-2)
1055
1056
    glat, glon, _ = apex_out.apex2geo(60, 15, 100, precision=1e-2)
1057
    (_, _, _, _, f1_, f2_, _, d1_, d2_, d3_, _, e1_, e2_,
1058
     e3_) = apex_out._geo2apexall(glat, glon, 100)
1059
1060
    assert_allclose(f1, f1_)
1061
    assert_allclose(f2, f2_)
1062
    assert_allclose(d1, d1_)
1063
    assert_allclose(d2, d2_)
1064
    assert_allclose(d3, d3_)
1065
    assert_allclose(e1, e1_)
1066
    assert_allclose(e2, e2_)
1067
    assert_allclose(e3, e3_)
1068
1069
1070 View Code Duplication
def test_basevectors_apex_scalar_qd():
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
1071
    apex_out = Apex(date=2000, refh=300)
1072
1073
    (f1, f2, f3, g1, g2, g3, d1, d2, d3, e1, e2,
1074
     e3) = apex_out.basevectors_apex(60, 15, 100, coords='qd', precision=1e-2)
1075
1076
    glat, glon, _ = apex_out.qd2geo(60, 15, 100, precision=1e-2)
1077
    (_, _, _, _, f1_, f2_, _, d1_, d2_, d3_, _, e1_, e2_,
1078
     e3_) = apex_out._geo2apexall(glat, glon, 100)
1079
1080
    assert_allclose(f1, f1_)
1081
    assert_allclose(f2, f2_)
1082
    assert_allclose(d1, d1_)
1083
    assert_allclose(d2, d2_)
1084
    assert_allclose(d3, d3_)
1085
    assert_allclose(e1, e1_)
1086
    assert_allclose(e2, e2_)
1087
    assert_allclose(e3, e3_)
1088
1089
1090
# test shapes and vectorization of arguments
1091
1092
def test_basevectors_apex_scalar_shape():
1093
    apex_out = Apex(date=2000, refh=300)
1094
    ret = apex_out.basevectors_apex(60, 15, 100, precision=1e-2)
1095
    for r in ret[:2]:
1096
        assert r.shape == (2,)
1097
    for r in ret[2:]:
1098
        assert r.shape == (3,)
1099
1100
1101
def test_basevectors_apex_vectorization():
1102
    apex_out = Apex(date=2000, refh=300)
1103
    ret = apex_out.basevectors_apex([60, 60, 60, 60], 15, 100)
1104
    for r in ret[:2]:
1105
        assert r.shape == (2, 4)
1106
    for r in ret[2:]:
1107
        assert r.shape == (3, 4)
1108
    ret = apex_out.basevectors_apex(60, [15, 15, 15, 15], 100)
1109
    for r in ret[:2]:
1110
        assert r.shape == (2, 4)
1111
    for r in ret[2:]:
1112
        assert r.shape == (3, 4)
1113
    ret = apex_out.basevectors_apex(60, 15, [100, 100, 100, 100])
1114
    for r in ret[:2]:
1115
        assert r.shape == (2, 4)
1116
    for r in ret[2:]:
1117
        assert r.shape == (3, 4)
1118
1119
1120
# test correct vectorization of height
1121 View Code Duplication
def test_basevectors_apex_vectorization_height():
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
1122
    apex_out = Apex(date=2000, refh=0)
1123
    (f1, f2, f3, g1, g2, g3, d1, d2, d3, e1, e2,
1124
     e3) = apex_out.basevectors_apex(60, 15, [200, 400], coords='geo')
1125
    (_, _, _, _, f1_1, f2_1, _, d1_1, d2_1, d3_1, _, e1_1, e2_1,
1126
     e3_1) = apex_out._geo2apexall(60, 15, 200)
1127
    (_, _, _, _, f1_2, f2_2, _, d1_2, d2_2, d3_2, _, e1_2, e2_2,
1128
     e3_2) = apex_out._geo2apexall(60, 15, 400)
1129
1130
    assert_allclose(f1[:, 0], f1_1)
1131
    assert_allclose(f2[:, 0], f2_1)
1132
    assert_allclose(d1[:, 0], d1_1)
1133
    assert_allclose(d2[:, 0], d2_1)
1134
    assert_allclose(d3[:, 0], d3_1)
1135
    assert_allclose(e1[:, 0], e1_1)
1136
    assert_allclose(e2[:, 0], e2_1)
1137
    assert_allclose(e3[:, 0], e3_1)
1138
1139
    assert_allclose(f3[:, 0], np.array([-0.088671, -0.018272, 0.993576]),
1140
                    rtol=1e-4)
1141
    assert_allclose(g1[:, 0], np.array([0.903098, 0.245273, 0.085107]),
1142
                    rtol=1e-4)
1143
    assert_allclose(g2[:, 0], np.array([-0.103495, 1.072078, 0.01048]),
1144
                    rtol=1e-4)
1145
    assert_allclose(g3[:, 0], np.array([0, 0, 1.006465]), rtol=1e-4)
1146
1147
    assert_allclose(f1[:, 1], f1_2)
1148
    assert_allclose(f2[:, 1], f2_2)
1149
    assert_allclose(d1[:, 1], d1_2)
1150
    assert_allclose(d2[:, 1], d2_2)
1151
    assert_allclose(d3[:, 1], d3_2)
1152
    assert_allclose(e1[:, 1], e1_2)
1153
    assert_allclose(e2[:, 1], e2_2)
1154
    assert_allclose(e3[:, 1], e3_2)
1155
1156
    assert_allclose(f3[:, 1], np.array([-0.085415, -0.021176, 0.989645]),
1157
                    rtol=1e-4)
1158
    assert_allclose(g1[:, 1], np.array([0.902695, 0.246919, 0.083194]),
1159
                    rtol=1e-4)
1160
    assert_allclose(g2[:, 1], np.array([-0.11051, 1.066094, 0.013274]),
1161
                    rtol=1e-4)
1162
    assert_allclose(g3[:, 1], np.array([0, 0, 1.010463]), rtol=1e-4)
1163
1164
1165
# test scalar return values
1166
1167
def test_basevectors_apex_scalar():
1168
    apex_out = Apex(date=2000, refh=300)
1169
1170
    (f1, f2, f3, g1, g2, g3, d1, d2, d3, e1, e2,
1171
     e3) = apex_out.basevectors_apex(0, 15, 100, coords='geo')
1172
    (_, _, _, _, f1_1, f2_1, _, d1_1, d2_1, d3_1, _, e1_1, e2_1,
1173
     e3_1) = apex_out._geo2apexall(0, 15, 100)
1174
1175
    assert_allclose(f1, f1_1)
1176
    assert_allclose(f2, f2_1)
1177
    assert_allclose(d1, d1_1)
1178
    assert_allclose(d2, d2_1)
1179
    assert_allclose(d3, d3_1)
1180
    assert_allclose(e1, e1_1)
1181
    assert_allclose(e2, e2_1)
1182
    assert_allclose(e3, e3_1)
1183
1184
    assert_allclose(f3, np.array([0.092637, -0.245951, 0.938848]), rtol=1e-4)
1185
    assert_allclose(g1, np.array([0.939012, 0.073416, -0.07342]), rtol=1e-4)
1186
    assert_allclose(g2, np.array([0.055389, 1.004155, 0.257594]), rtol=1e-4)
1187
    assert_allclose(g3, np.array([0, 0, 1.065135]), rtol=1e-4)
1188
1189
1190
# test 1D array return values
1191
1192 View Code Duplication
def test_basevectors_apex_array():
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
1193
    apex_out = Apex(date=2000, refh=300)
1194
    (f1, f2, f3, g1, g2, g3, d1, d2, d3, e1, e2,
1195
     e3) = apex_out.basevectors_apex([0, 30], 15, 100, coords='geo')
1196
    (_, _, _, _, f1_1, f2_1, _, d1_1, d2_1, d3_1, _, e1_1, e2_1,
1197
     e3_1) = apex_out._geo2apexall(0, 15, 100)
1198
    (_, _, _, _, f1_2, f2_2, _, d1_2, d2_2, d3_2, _, e1_2, e2_2,
1199
     e3_2) = apex_out._geo2apexall(30, 15, 100)
1200
1201
    assert_allclose(f1[:, 0], f1_1)
1202
    assert_allclose(f2[:, 0], f2_1)
1203
    assert_allclose(d1[:, 0], d1_1)
1204
    assert_allclose(d2[:, 0], d2_1)
1205
    assert_allclose(d3[:, 0], d3_1)
1206
    assert_allclose(e1[:, 0], e1_1)
1207
    assert_allclose(e2[:, 0], e2_1)
1208
    assert_allclose(e3[:, 0], e3_1)
1209
1210
    assert_allclose(f3[:, 0], np.array([0.092637, -0.245951, 0.938848]),
1211
                    rtol=1e-4)
1212
    assert_allclose(g1[:, 0], np.array([0.939012, 0.073416, -0.07342]),
1213
                    rtol=1e-4)
1214
    assert_allclose(g2[:, 0], np.array([0.055389, 1.004155, 0.257594]),
1215
                    rtol=1e-4)
1216
    assert_allclose(g3[:, 0], np.array([0, 0, 1.065135]), rtol=1e-4)
1217
1218
    assert_allclose(f1[:, 1], f1_2)
1219
    assert_allclose(f2[:, 1], f2_2)
1220
    assert_allclose(d1[:, 1], d1_2)
1221
    assert_allclose(d2[:, 1], d2_2)
1222
    assert_allclose(d3[:, 1], d3_2)
1223
    assert_allclose(e1[:, 1], e1_2)
1224
    assert_allclose(e2[:, 1], e2_2)
1225
    assert_allclose(e3[:, 1], e3_2)
1226
1227
    assert_allclose(f3[:, 1], np.array([-0.036618, -0.071019, 0.861604]),
1228
                    rtol=1e-4)
1229
    assert_allclose(g1[:, 1], np.array([0.844391, 0.015353, 0.037152]),
1230
                    rtol=1e-4)
1231
    assert_allclose(g2[:, 1], np.array([0.050808, 1.02131, 0.086342]),
1232
                    rtol=1e-4)
1233
    assert_allclose(g3[:, 1], np.array([0, 0, 1.160625]), rtol=1e-4)
1234
1235
1236
# test that vectors are calculated correctly
1237
1238
def test_basevectors_apex_delta():
1239
    apex_out = Apex(date=2000, refh=300)
1240
    for lat in range(0, 90, 10):
1241
        for lon in range(0, 360, 15):
1242
            (f1, f2, f3, g1, g2, g3, d1, d2, d3, e1, e2,
1243
             e3) = apex_out.basevectors_apex(lat, lon, 500)
1244
            f = [np.append(f1, 0), np.append(f2, 0), f3]
1245
            g = [g1, g2, g3]
1246
            d = [d1, d2, d3]
1247
            e = [e1, e2, e3]
1248
            for i, j in [(i, j) for i in range(3) for j in range(3)]:
1249
                delta = 1 if i == j else 0
1250
                assert_allclose(np.sum(f[i] * g[j]), delta, rtol=0, atol=1e-5)
1251
                assert_allclose(np.sum(d[i] * e[j]), delta, rtol=0, atol=1e-5)
1252
1253
1254
def test_basevectors_apex_invalid_scalar(recwarn):
1255
    """Test warning and fill values for calculating base vectors with bad value.
1256
    """
1257
    apex_out = Apex(date=2000, refh=10000)
1258
    base_vecs = apex_out.basevectors_apex(0, 0, 0)
1259
1260
    assert issubclass(recwarn[-1].category, UserWarning)
1261
    assert 'set to NaN where' in str(recwarn[-1].message)
1262
1263
    invalid = np.ones(3) * np.nan
1264
    for i, bvec in enumerate(base_vecs):
1265
        if i < 2:
1266
            assert not np.allclose(bvec, invalid[:2])
1267
        else:
1268
            assert_allclose(bvec, invalid)
1269
1270
1271
# ============================================================================
1272
#  Test the get_apex() method
1273
# ============================================================================
1274
1275
1276
def test_get_apex():
1277
    apex_out = Apex(date=2000, refh=300)
1278
    assert_allclose(apex_out.get_apex(10), 507.409702543805)
1279
    assert_allclose(apex_out.get_apex(60), 20313.026999999987)
1280
1281
1282
def test_get_apex_invalid_lat():
1283
    apex_out = Apex(date=2000, refh=300)
1284
    with pytest.raises(ValueError):
1285
        apex_out.get_apex(91)
1286
    with pytest.raises(ValueError):
1287
        apex_out.get_apex(-91)
1288
    apex_out.get_apex(90)
1289
    apex_out.get_apex(-90)
1290
1291
    assert_allclose(apex_out.get_apex(90 + 1e-5), apex_out.get_apex(90),
1292
                    rtol=0, atol=1e-8)
1293
1294
1295
# ============================================================================
1296
#  Test the set_epoch() method
1297
# ============================================================================
1298
1299
1300
def test_set_epoch():
1301
    """Test successful setting of Apex epoch."""
1302
    apex_out = Apex(date=2000.2, refh=300)
1303
    assert_allclose(apex_out.year, 2000.2)
1304
    ret_2000_2_py = apex_out._geo2apex(60, 15, 100)
1305
    apex_out.set_epoch(2000.8)
1306
    assert_allclose(apex_out.year, 2000.8)
1307
    ret_2000_8_py = apex_out._geo2apex(60, 15, 100)
1308
1309
    assert ret_2000_2_py != ret_2000_8_py
1310
1311
    fa.loadapxsh(apex_out.datafile, 2000.2)
1312
    ret_2000_2_apex = fa.apxg2all(60, 15, 100, 300, 0)[2:4]
1313
    fa.loadapxsh(apex_out.datafile, 2000.8)
1314
    ret_2000_8_apex = fa.apxg2all(60, 15, 100, 300, 0)[2:4]
1315
1316
    assert ret_2000_2_apex != ret_2000_8_apex
1317
1318
    assert_allclose(ret_2000_2_py, ret_2000_2_apex)
1319
    assert_allclose(ret_2000_8_py, ret_2000_8_apex)
1320
1321
1322
@pytest.fixture()
1323
def igrf_file():
1324
    # Ensure the coefficient file exists
1325
    original_file = os.path.join(os.path.dirname(helpers.__file__),
1326
                                 'igrf13coeffs.txt')
1327
    tmp_file = "temp_coeff.txt"
1328
    assert os.path.isfile(original_file)
1329
    # Move the coefficient file
1330
    os.rename(original_file, tmp_file)
1331
    yield original_file
1332
    # Move the coefficient file back
1333
    os.rename(tmp_file, original_file)
1334
1335
1336
def test_set_epoch_file_error(igrf_file):
1337
    """Test raises OSError when IGRF coefficient file is missing."""
1338
    # Test missing coefficient file failure
1339
    with pytest.raises(OSError) as oerr:
1340
        Apex(date=2000.2, refh=300)
1341
    error_string = "File {:} does not exist".format(igrf_file)
1342
    assert str(oerr.value).startswith(error_string)
1343
1344
1345
# ============================================================================
1346
#  Test the set_refh() method
1347
# ============================================================================
1348
1349
1350
def test_set_refh():
1351
    apex_out = Apex(date=2000, refh=300)
1352
    assert apex_out.refh, 300
1353
    ret_300 = apex_out._geo2apex(60, 15, 100)
1354
    apex_out.set_refh(500)
1355
    assert apex_out.refh == 500
1356
    ret_500 = apex_out._geo2apex(60, 15, 100)
1357
1358
    assert_allclose(ret_300, fa.apxg2all(60, 15, 100, 300, 0)[2:4])
1359
    assert_allclose(ret_500, fa.apxg2all(60, 15, 100, 500, 0)[2:4])
1360
1361
1362
# ============================================================================
1363
#  Test the get_babs() method
1364
# ============================================================================
1365
1366
1367
def test_get_babs():
1368
    inputs = [[[80], [100], [300]], [range(50, 90, 8), range(0, 360, 80),
1369
                                     [300] * 5], [90.0, 0, 1000]]
1370
    temp1 = np.array([4.22045410e-05, 5.15672743e-05, 4.98150200e-05,
1371
                     5.06769359e-05, 4.91028428e-05])
1372
    expected = [[5.1303124427795412e-05], temp1, [3.793962299823761e-05]]
1373
1374
    apex_out = Apex(date=2018.1, refh=0)
1375
    for i in range(len(inputs)):
1376
        outputs = apex_out.get_babs(*inputs[i])
1377
        if isinstance(outputs, np.float64):
1378
            outputs = [outputs]
1379
        for j, output in enumerate(outputs):
1380
            assert_allclose(output, expected[i][j], rtol=0, atol=1e-5)
1381
1382
1383
# ============================================================================
1384
#  Test the bvectors_apex() method
1385
# ============================================================================
1386
1387
1388
def test_bvectors_apex():
1389
    inputs = [[80, 81], [100, 120], [100, 200]]
1390
1391
    expected = (np.array([5.94623305e-05, 5.95450722e-05]),
1392
                np.array([[0.02008877, 0.00303204],
1393
                          [0.03571109, 0.03377986],
1394
                          [-0.94045794, -0.89848483]]),
1395
                np.array([5.26919505e-05, 4.81377429e-05]),
1396
                np.array([[0.02266997, 0.00375055],
1397
                          [0.04029961, 0.04178477],
1398
                          [-1.0612973, -1.1114012]]))
1399
1400
    apex_out = Apex(date=2018.1, refh=0)
1401
1402
    outputs = apex_out.bvectors_apex(*inputs, coords='geo', precision=1e-10)
1403
    for i, output in enumerate(outputs):
1404
        for j in range(output.size):
1405
            assert_allclose(output.ravel()[j], expected[i].ravel()[j], rtol=0,
1406
                            atol=1e-5)
1407
1408
1409
if __name__ == '__main__':
1410
    pytest.main()
1411