6
6
import numpy as np
7
7
import xarray
8
8
from typing import Sequence
9
+
9
10
#
10
11
import iri90fort
11
12
12
13
rdir = Path (__file__ ).parent
13
- Ts = ['Tn' , 'Ti' , 'Te' ]
14
- simout = ['ne' , 'Tn' , 'Ti' , 'Te' , ' nO+' , ' nH+' , ' nHe+' , ' nO2+' , ' nNO+' ]
14
+ Ts = ["Tn" , "Ti" , "Te" ]
15
+ simout = ["ne" , "Tn" , "Ti" , "Te" , " nO+" , " nH+" , " nHe+" , " nO2+" , " nNO+" ]
15
16
16
17
17
18
def datetimerange (start : datetime , end : datetime , step : timedelta ) -> list :
@@ -20,72 +21,93 @@ def datetimerange(start: datetime, end: datetime, step: timedelta) -> list:
20
21
assert isinstance (end , datetime )
21
22
assert isinstance (step , timedelta )
22
23
23
- return [start + i * step for i in range ((end - start ) // step )]
24
-
24
+ return [start + i * step for i in range ((end - start ) // step )]
25
25
26
- def runiri (time : datetime , altkm : float , glatlon : tuple , f107 : float , f107a : float , ap : int ) -> xarray .DataArray :
27
26
27
+ def runiri (
28
+ time : datetime , altkm : np .ndarray , glatlon : tuple , f107 : float , f107a : float , ap : int
29
+ ) -> xarray .DataArray :
28
30
def _collect_output () -> xarray .DataArray :
29
- """ collect IRI90 output into xarray.DataArray with metadata"""
30
- iono = xarray .DataArray (outf [:9 , :].T ,
31
- coords = {'alt_km' : altkm , 'sim' : simout },
32
- dims = ['alt_km' , 'sim' ],
33
- attrs = {'f107' : f107 , 'f107a' : f107a , 'ap' : ap , 'glatlon' : glatlon , 'time' : time })
34
-
35
- # i=(iono['Ti']<iono['Tn']).values
36
- # iono.ix[i,'Ti'] = iono.ix[i,'Tn']
37
-
38
- # i=(iono['Te']<iono['Tn']).values
39
- # iono.ix[i,'Te'] = iono.ix[i,'Tn']
40
-
41
- # %% iri90 outputs percentage of Ne
42
- iono .loc [:, ['nO+' , 'nH+' , 'nHe+' , 'nO2+' , 'nNO+' ]] *= iono .loc [:, 'ne' ]/ 100.
43
- # %% These two parameters only output if JF(6)=False, otherwise bogus values
31
+ """collect IRI90 output into xarray.DataArray with metadata"""
32
+ iono = xarray .DataArray (
33
+ outf [:9 , :].T ,
34
+ coords = {"alt_km" : altkm , "sim" : simout },
35
+ dims = ["alt_km" , "sim" ],
36
+ attrs = {
37
+ "f107" : f107 ,
38
+ "f107a" : f107a ,
39
+ "ap" : ap ,
40
+ "glatlon" : glatlon ,
41
+ "time" : time ,
42
+ },
43
+ )
44
+
45
+ # i=(iono['Ti']<iono['Tn']).values
46
+ # iono.ix[i,'Ti'] = iono.ix[i,'Tn']
47
+
48
+ # i=(iono['Te']<iono['Tn']).values
49
+ # iono.ix[i,'Te'] = iono.ix[i,'Tn']
50
+
51
+ # %% iri90 outputs percentage of Ne
52
+ iono .loc [:, ["nO+" , "nH+" , "nHe+" , "nO2+" , "nNO+" ]] *= iono .loc [:, "ne" ] / 100.0
53
+ # %% These two parameters only output if JF(6)=False, otherwise bogus values
44
54
# iono['nClusterIons'] = iono['ne'] * outf[9,:]/100.
45
55
# iono['nN+'] = iono['ne'] * outf[10,:]/100.
46
- # %% negative indicates undefined
56
+ # %% negative indicates undefined
47
57
for c in iono .sim :
48
- iono .loc [iono .loc [:, c ] <= 0. , c ] = np .nan
58
+ iono .loc [iono .loc [:, c ] <= 0.0 , c ] = np .nan
49
59
50
60
return iono
51
61
52
62
glat , glon = glatlon
53
63
jmag = 0 # coordinates are: 0:geographic 1: magnetic
54
64
55
- JF = np .array ((1 , 1 , 1 , 1 , 0 , 1 ,
56
- 1 , 1 , 1 , 1 , 1 , 1 ), bool ) # Solomon 1993 version of IRI
65
+ JF = np .array ((1 , 1 , 1 , 1 , 0 , 1 , 1 , 1 , 1 , 1 , 1 , 1 ), bool ) # Solomon 1993 version of IRI
57
66
# JF = (1,1,1) + (0,0,0) +(1,)*14 + (0,1,0,1,1,1,1,0,0,0,1,1,0,1,0,1,1,1) #for 2013 version of IRI
58
67
59
- monthday = time .month * 100 + time .day # yep, that's how the IRI code wants it, NOT as character.
68
+ monthday = (
69
+ time .month * 100 + time .day
70
+ ) # yep, that's how the IRI code wants it, NOT as character.
60
71
# + 25 hours for UTC time
61
- hourfrac = (time .hour + 25 ) + time .minute / 60 + time .second / 3600
62
- datadir = str (rdir / 'data' )+ '/'
63
- # %% call IRI
64
- outf , oarr = iri90fort .iri90 (JF , jmag ,
65
- glat , glon % 360. ,
66
- - f107 ,
67
- monthday , # integer
68
- hourfrac ,
69
- altkm ,
70
- datadir )
72
+ hourfrac = (time .hour + 25 ) + time .minute / 60 + time .second / 3600
73
+ datadir = str (rdir / "data" ) + "/"
74
+ # %% call IRI
75
+ outf , oarr = iri90fort .iri90 (
76
+ JF ,
77
+ jmag ,
78
+ glat ,
79
+ glon % 360.0 ,
80
+ - f107 ,
81
+ monthday , # integer
82
+ hourfrac ,
83
+ altkm ,
84
+ datadir ,
85
+ )
71
86
72
87
iono = _collect_output ()
73
88
74
89
return iono
75
90
76
91
77
- def timeprofile (tlim : Sequence [datetime ], dt : timedelta ,
78
- altkm : np .ndarray , glatlon : tuple ,
79
- f107 : float , f107a : float , ap : int ) -> xarray .DataArray :
92
+ def timeprofile (
93
+ tlim : Sequence [datetime ],
94
+ dt : timedelta ,
95
+ altkm : np .ndarray ,
96
+ glatlon : tuple ,
97
+ f107 : float ,
98
+ f107a : float ,
99
+ ap : int ,
100
+ ) -> xarray .DataArray :
80
101
"""compute IRI90 at a single altiude, over time range"""
81
102
82
103
T = datetimerange (tlim [0 ], tlim [1 ], dt )
83
104
84
- iono = xarray .DataArray (np .empty ((len (T ), altkm .size , 9 )),
85
- coords = {'time' : T , 'alt_km' : altkm , 'sim' : simout },
86
- dims = ['time' , 'alt_km' , 'sim' ],
87
- attrs = {'f107' : f107 , 'f107a' : f107a , 'ap' : ap , 'glatlon' : glatlon }
88
- )
105
+ iono = xarray .DataArray (
106
+ np .empty ((len (T ), altkm .size , 9 )),
107
+ coords = {"time" : T , "alt_km" : altkm , "sim" : simout },
108
+ dims = ["time" , "alt_km" , "sim" ],
109
+ attrs = {"f107" : f107 , "f107a" : f107a , "ap" : ap , "glatlon" : glatlon },
110
+ )
89
111
90
112
for t in T :
91
113
iono .loc [t , ...] = runiri (t , altkm , glatlon , f107 , f107a , ap )
0 commit comments