|
1
|
|
|
# -*- coding: utf-8 -*- |
|
2
|
|
|
# vim:fileencoding=utf-8 |
|
3
|
|
|
# |
|
4
|
|
|
# Copyright (c) 2018 Stefan Bender |
|
5
|
|
|
# |
|
6
|
|
|
# This file is part of sciapy. |
|
7
|
|
|
# sciapy is free software: you can redistribute it or modify it |
|
8
|
|
|
# under the terms of the GNU General Public License as published by |
|
9
|
|
|
# the Free Software Foundation, version 2. |
|
10
|
|
|
# See accompanying LICENSE file or http://www.gnu.org/licenses/gpl-2.0.html. |
|
11
|
1 |
|
"""SCIAMACHY level 2 post-processed number densities interface |
|
12
|
|
|
|
|
13
|
|
|
Interface classes for the level 2 post-processed retrieval results. |
|
14
|
|
|
""" |
|
15
|
|
|
|
|
16
|
1 |
|
from __future__ import absolute_import, division, print_function |
|
17
|
|
|
|
|
18
|
1 |
|
import datetime as dt |
|
19
|
|
|
|
|
20
|
1 |
|
import numpy as np |
|
21
|
1 |
|
try: |
|
22
|
1 |
|
from netCDF4 import Dataset as netcdf_file |
|
23
|
1 |
|
fmtargs = {"format": "NETCDF4"} |
|
24
|
|
|
except ImportError: |
|
25
|
|
|
try: |
|
26
|
|
|
from scipy.io.netcdf import netcdf_file |
|
27
|
|
|
fmtargs = {"version": 1} |
|
28
|
|
|
except ImportError: |
|
29
|
|
|
from pupynere import netcdf_file |
|
30
|
|
|
fmtargs = {"version": 1} |
|
31
|
|
|
|
|
32
|
1 |
|
from .density import scia_densities, _UTC |
|
33
|
1 |
|
from .. import __version__ |
|
34
|
|
|
|
|
35
|
1 |
|
__all__ = ["scia_densities_pp", "scia_density_day"] |
|
36
|
|
|
|
|
37
|
|
|
|
|
38
|
1 |
|
class scia_densities_pp(scia_densities): |
|
39
|
|
|
"""Post-processed SCIAMACHY number densities |
|
40
|
|
|
|
|
41
|
|
|
Extends `scia_densities` with additional post-processing |
|
42
|
|
|
attributes such as (MSIS) temperature and density, local |
|
43
|
|
|
solar time, solar zenith angle, and geomagnetic latitudes |
|
44
|
|
|
and longitudes. |
|
45
|
|
|
|
|
46
|
|
|
This class only supports writing ascii files but reading to |
|
47
|
|
|
and writing from netcdf. |
|
48
|
|
|
|
|
49
|
|
|
Attributes |
|
50
|
|
|
---------- |
|
51
|
|
|
temperature |
|
52
|
|
|
NRLMSISE-00 temperatures |
|
53
|
|
|
noem_no |
|
54
|
|
|
NOEM NO nuimber densities |
|
55
|
|
|
vmr |
|
56
|
|
|
NO vmr using the NRLMSISE-00 total air number densities |
|
57
|
|
|
lst |
|
58
|
|
|
Apparent local solar times |
|
59
|
|
|
mst |
|
60
|
|
|
Mean local solar times |
|
61
|
|
|
sza |
|
62
|
|
|
Solar zenith angles |
|
63
|
|
|
utchour |
|
64
|
|
|
UTC hours into measurement day |
|
65
|
|
|
utcdays |
|
66
|
|
|
Number of days since reference date |
|
67
|
|
|
gmlats |
|
68
|
|
|
IGRF-12 geomagentic latitudes |
|
69
|
|
|
gmlons |
|
70
|
|
|
IGRF-12 geomagentic longitudes |
|
71
|
|
|
aacgmgmlats |
|
72
|
|
|
AACGM geomagentic latitudes |
|
73
|
|
|
aacgmgmlons |
|
74
|
|
|
AACGM geomagentic longitudes |
|
75
|
|
|
|
|
76
|
|
|
Methods |
|
77
|
|
|
------- |
|
78
|
|
|
write_to_textfile |
|
79
|
|
|
write_to_netcdf |
|
80
|
|
|
read_from_netcdf |
|
81
|
|
|
to_xarray |
|
82
|
|
|
""" |
|
83
|
1 |
|
def __init__(self, ref_date="2000-01-01", |
|
84
|
|
|
ver=None, data_ver=None): |
|
85
|
1 |
|
self.filename = None |
|
86
|
1 |
|
self.temperature = None |
|
87
|
1 |
|
self.noem_no = None |
|
88
|
1 |
|
self.vmr = None |
|
89
|
1 |
|
self.lst = None |
|
90
|
1 |
|
self.mst = None |
|
91
|
1 |
|
self.sza = None |
|
92
|
1 |
|
self.utchour = None |
|
93
|
1 |
|
self.utcdays = None |
|
94
|
1 |
|
self.gmlats = None |
|
95
|
1 |
|
self.gmlons = None |
|
96
|
1 |
|
self.aacgmgmlats = None |
|
97
|
1 |
|
self.aacgmgmlons = None |
|
98
|
1 |
|
super(scia_densities_pp, self).__init__( |
|
99
|
|
|
ref_date=ref_date, ver=ver, data_ver=data_ver) |
|
100
|
|
|
|
|
101
|
1 |
|
def write_to_textfile(self, filename): |
|
102
|
|
|
"""Write the variables to ascii table files |
|
103
|
|
|
|
|
104
|
|
|
Parameters |
|
105
|
|
|
---------- |
|
106
|
|
|
filename: str or file object or io.TextIOBase.buffer |
|
107
|
|
|
The filename or stream to write the data to. For writing to |
|
108
|
|
|
stdout in python 3, pass `sys.stdout.buffer`. |
|
109
|
|
|
""" |
|
110
|
|
|
if hasattr(filename, 'seek'): |
|
111
|
|
|
f = filename |
|
112
|
|
|
else: |
|
113
|
|
|
f = open(filename, 'w') |
|
114
|
|
|
|
|
115
|
|
|
header = "%5s %13s %12s %13s %14s %12s %14s %12s %12s %12s %12s %12s" % ("GP_ID", |
|
116
|
|
|
"Max_Hoehe[km]", "Hoehe[km]", "Min_Hoehe[km]", |
|
117
|
|
|
"Max_Breite[°]", "Breite[°]", "Min_Breite[°]", |
|
118
|
|
|
"Laenge[°]", |
|
119
|
|
|
"Dichte[cm^-3]", "Fehler Mess[cm^-3]", |
|
120
|
|
|
"Fehler tot[cm^-3]", "Gesamtdichte[cm^-3]") |
|
121
|
|
|
if self.apriori is not None: |
|
122
|
|
|
header = header + " %12s" % ("apriori[cm^-3]",) |
|
123
|
|
|
if self.temperature is not None: |
|
124
|
|
|
header = header + " %12s" % ("T[K]",) |
|
125
|
|
|
if self.noem_no is not None: |
|
126
|
|
|
header = header + " %12s" % ("NOEM_NO[cm^-3]",) |
|
127
|
|
|
if self.akdiag is not None: |
|
128
|
|
|
header = header + " %12s" % ("AKdiag",) |
|
129
|
|
|
if self.lst is not None: |
|
130
|
|
|
header = header + " %12s" % ("LST",) |
|
131
|
|
|
if self.mst is not None: |
|
132
|
|
|
header = header + " %12s" % ("MST",) |
|
133
|
|
|
if self.sza is not None: |
|
134
|
|
|
header = header + " %12s" % ("SZA",) |
|
135
|
|
|
if self.utchour is not None: |
|
136
|
|
|
header = header + " %12s" % ("Hour",) |
|
137
|
|
|
if self.utcdays is not None: |
|
138
|
|
|
header = header + " %12s" % ("Days",) |
|
139
|
|
|
if self.gmlats is not None: |
|
140
|
|
|
header = header + " %12s" % ("GeomagLat",) |
|
141
|
|
|
if self.gmlons is not None: |
|
142
|
|
|
header = header + " %12s" % ("GeomagLon",) |
|
143
|
|
|
if self.aacgmgmlats is not None: |
|
144
|
|
|
header = header + " %12s" % ("AACGMGeomagLat",) |
|
145
|
|
|
if self.aacgmgmlons is not None: |
|
146
|
|
|
header = header + " %12s" % ("AACGMGeomagLon",) |
|
147
|
|
|
print(header, file=f) |
|
148
|
|
|
|
|
149
|
|
|
oformat = "%5i %+1.5E %+1.5E %+1.5E %+1.5E %+1.5E %+1.5E %+1.5E %+1.5E %+1.5E %+1.5E %+1.5E" |
|
150
|
|
|
oformata = " %+1.5E" |
|
151
|
|
|
|
|
152
|
|
|
for i, a in enumerate(self.alts): |
|
153
|
|
|
for j, b in enumerate(self.lats): |
|
154
|
|
|
print(oformat % (i * self.nlat + j, |
|
155
|
|
|
self.alts_max[i], a, self.alts_min[i], |
|
156
|
|
|
self.lats_max[j], b, self.lats_min[j], self.lons[j], |
|
157
|
|
|
self.densities[j, i], self.dens_err_meas[j, i], |
|
158
|
|
|
self.dens_err_tot[j, i], self.dens_tot[j, i]), |
|
159
|
|
|
end="", file=f) |
|
160
|
|
|
if self.apriori is not None: |
|
161
|
|
|
print(" " + oformata % self.apriori[j, i], end="", file=f) |
|
162
|
|
|
if self.temperature is not None: |
|
163
|
|
|
print(" " + oformata % self.temperature[j, i], end="", file=f) |
|
164
|
|
|
if self.noem_no is not None: |
|
165
|
|
|
print(" " + oformata % self.noem_no[j, i], end="", file=f) |
|
166
|
|
|
if self.akdiag is not None: |
|
167
|
|
|
print(" " + oformata % self.akdiag[j, i], end="", file=f) |
|
168
|
|
|
if self.lst is not None: |
|
169
|
|
|
print(" " + oformata % self.lst[j], end="", file=f) |
|
170
|
|
|
if self.mst is not None: |
|
171
|
|
|
print(" " + oformata % self.mst[j], end="", file=f) |
|
172
|
|
|
if self.sza is not None: |
|
173
|
|
|
print(" " + oformata % self.sza[j], end="", file=f) |
|
174
|
|
|
if self.utchour is not None: |
|
175
|
|
|
print(" " + oformata % self.utchour[j], end="", file=f) |
|
176
|
|
|
if self.utcdays is not None: |
|
177
|
|
|
print(" " + oformata % self.utcdays[j], end="", file=f) |
|
178
|
|
|
if self.gmlats is not None: |
|
179
|
|
|
print(" " + oformata % self.gmlats[j], end="", file=f) |
|
180
|
|
|
if self.gmlons is not None: |
|
181
|
|
|
print(" " + oformata % self.gmlons[j], end="", file=f) |
|
182
|
|
|
if self.aacgmgmlats is not None: |
|
183
|
|
|
print(" " + oformata % self.aacgmgmlats[j], end="", file=f) |
|
184
|
|
|
if self.aacgmgmlons is not None: |
|
185
|
|
|
print(" " + oformata % self.aacgmgmlons[j], end="", file=f) |
|
186
|
|
|
print("", file=f) |
|
187
|
|
|
|
|
188
|
1 |
|
def write_to_netcdf(self, filename, close=True): |
|
189
|
|
|
"""Write variables to netcdf files |
|
190
|
|
|
|
|
191
|
|
|
This function has no stream, i.e. file object, support. |
|
192
|
|
|
|
|
193
|
|
|
Parameters |
|
194
|
|
|
---------- |
|
195
|
|
|
filename: str |
|
196
|
|
|
The name of the file to write the data to. |
|
197
|
|
|
""" |
|
198
|
|
|
# write the base variables first and keep the file open for appending |
|
199
|
1 |
|
ncf = scia_densities.write_to_netcdf(self, filename, close=False) |
|
200
|
|
|
|
|
201
|
1 |
|
if self.temperature is not None: |
|
202
|
1 |
|
ftemp = ncf.createVariable('temperature', |
|
203
|
|
|
np.dtype('float64').char, ('time', 'latitude', 'altitude')) |
|
204
|
1 |
|
ftemp.units = 'K' |
|
205
|
1 |
|
ftemp.long_name = 'temperature' |
|
206
|
1 |
|
ftemp.model = 'NRLMSIS-00' |
|
207
|
1 |
|
ftemp[0, :] = self.temperature |
|
208
|
|
|
|
|
209
|
1 |
|
if self.noem_no is not None: |
|
210
|
1 |
|
fnoem_no = ncf.createVariable('NOEM_density', |
|
211
|
|
|
np.dtype('float64').char, ('time', 'latitude', 'altitude')) |
|
212
|
1 |
|
fnoem_no.units = 'cm^{-3}' |
|
213
|
1 |
|
fnoem_no.long_name = 'NOEM NO number density' |
|
214
|
1 |
|
fnoem_no[0, :] = self.noem_no |
|
215
|
|
|
|
|
216
|
1 |
|
if self.vmr is not None: |
|
217
|
1 |
|
fvmr = ncf.createVariable('VMR', |
|
218
|
|
|
np.dtype('float64').char, ('time', 'latitude', 'altitude')) |
|
219
|
1 |
|
fvmr.units = 'ppb' |
|
220
|
1 |
|
fvmr.long_name = 'volume mixing ratio' |
|
221
|
1 |
|
fvmr[0, :] = self.vmr |
|
222
|
|
|
|
|
223
|
1 |
|
if self.lst is not None: |
|
224
|
1 |
|
flst = ncf.createVariable('app_lst', |
|
225
|
|
|
np.dtype('float64').char, ('time', 'latitude',)) |
|
226
|
1 |
|
flst.units = 'hours' |
|
227
|
1 |
|
flst.long_name = 'apparent local solar time' |
|
228
|
1 |
|
flst[0, :] = self.lst |
|
229
|
|
|
|
|
230
|
1 |
|
if self.mst is not None: |
|
231
|
1 |
|
fmst = ncf.createVariable('mean_lst', |
|
232
|
|
|
np.dtype('float64').char, ('time', 'latitude',)) |
|
233
|
1 |
|
fmst.units = 'hours' |
|
234
|
1 |
|
fmst.long_name = 'mean local solar time' |
|
235
|
1 |
|
fmst[0, :] = self.mst |
|
236
|
|
|
|
|
237
|
1 |
|
if self.sza is not None: |
|
238
|
1 |
|
fsza = ncf.createVariable('mean_sza', |
|
239
|
|
|
np.dtype('float64').char, ('time', 'latitude',)) |
|
240
|
1 |
|
fsza.units = 'degrees' |
|
241
|
1 |
|
fsza.long_name = 'mean solar zenith angle' |
|
242
|
1 |
|
fsza[0, :] = self.sza |
|
243
|
|
|
|
|
244
|
1 |
|
if self.utchour is not None: |
|
245
|
1 |
|
futc = ncf.createVariable('utc_hour', |
|
246
|
|
|
np.dtype('float64').char, ('time', 'latitude',)) |
|
247
|
1 |
|
futc.units = 'hours' |
|
248
|
1 |
|
futc.long_name = 'measurement utc time' |
|
249
|
1 |
|
futc[0, :] = self.utchour |
|
250
|
|
|
|
|
251
|
1 |
|
if self.utcdays is not None: |
|
252
|
1 |
|
futcd = ncf.createVariable('utc_days', |
|
253
|
|
|
np.dtype('float64').char, ('time', 'latitude',)) |
|
254
|
1 |
|
futcd.long_name = 'measurement day' |
|
255
|
1 |
|
futcd.units = 'days since {0}'.format(self.date0.isoformat(sep=' ')) |
|
256
|
1 |
|
futcd[0, :] = self.utcdays |
|
257
|
|
|
|
|
258
|
1 |
|
if self.gmlats is not None: |
|
259
|
1 |
|
fgmlats = ncf.createVariable('gm_lats', |
|
260
|
|
|
np.dtype('float64').char, ('time', 'latitude',)) |
|
261
|
1 |
|
fgmlats.long_name = 'geomagnetic_latitude' |
|
262
|
1 |
|
fgmlats.model = 'IGRF' |
|
263
|
1 |
|
fgmlats.units = 'degrees_north' |
|
264
|
1 |
|
fgmlats[0, :] = self.gmlats |
|
265
|
|
|
|
|
266
|
1 |
|
if self.gmlons is not None: |
|
267
|
1 |
|
fgmlons = ncf.createVariable('gm_lons', |
|
268
|
|
|
np.dtype('float64').char, ('time', 'latitude',)) |
|
269
|
1 |
|
fgmlons.long_name = 'geomagnetic_longitude' |
|
270
|
1 |
|
fgmlons.model = 'IGRF' |
|
271
|
1 |
|
fgmlons.units = 'degrees_east' |
|
272
|
1 |
|
fgmlons[0, :] = self.gmlons |
|
273
|
|
|
|
|
274
|
1 |
|
if self.aacgmgmlats is not None: |
|
275
|
1 |
|
faacgmgmlats = ncf.createVariable('aacgm_gm_lats', |
|
276
|
|
|
np.dtype('float64').char, ('time', 'latitude',)) |
|
277
|
1 |
|
faacgmgmlats.long_name = 'geomagnetic_latitude' |
|
278
|
1 |
|
faacgmgmlats.model = 'AACGM' |
|
279
|
1 |
|
faacgmgmlats.units = 'degrees_north' |
|
280
|
1 |
|
faacgmgmlats[0, :] = self.aacgmgmlats |
|
281
|
|
|
|
|
282
|
1 |
|
if self.aacgmgmlons is not None: |
|
283
|
1 |
|
faacgmgmlons = ncf.createVariable('aacgm_gm_lons', |
|
284
|
|
|
np.dtype('float64').char, ('time', 'latitude',)) |
|
285
|
1 |
|
faacgmgmlons.long_name = 'geomagnetic_longitude' |
|
286
|
1 |
|
faacgmgmlons.model = 'AACGM' |
|
287
|
1 |
|
faacgmgmlons.units = 'degrees_east' |
|
288
|
1 |
|
faacgmgmlons[0, :] = self.aacgmgmlons |
|
289
|
|
|
|
|
290
|
1 |
|
if close: |
|
291
|
1 |
|
ncf.close() |
|
292
|
|
|
else: |
|
293
|
|
|
return ncf |
|
294
|
|
|
|
|
295
|
1 |
|
def read_from_netcdf(self, filename, close=True): |
|
296
|
|
|
"""Read post-processed level 2 orbit files |
|
297
|
|
|
|
|
298
|
|
|
Parameters |
|
299
|
|
|
---------- |
|
300
|
|
|
filename: str |
|
301
|
|
|
The name of the netcdf file. |
|
302
|
|
|
""" |
|
303
|
|
|
# read the base variables first and keep the file open for reading |
|
304
|
1 |
|
ncf = scia_densities.read_from_netcdf(self, filename, close=False) |
|
305
|
|
|
|
|
306
|
|
|
# additional data... |
|
307
|
|
|
# MSIS temperature |
|
308
|
|
|
try: |
|
309
|
|
|
self.temperature = ncf.variables['temperature'][:] |
|
310
|
|
|
except: |
|
311
|
|
|
pass |
|
312
|
|
|
# NOEM density |
|
313
|
|
|
try: |
|
314
|
|
|
self.noem_no = ncf.variables['NOEM_density'][:] |
|
315
|
|
|
except: |
|
316
|
|
|
pass |
|
317
|
|
|
# calculated vmr |
|
318
|
|
|
try: |
|
319
|
|
|
self.vmr = ncf.variables['VMR'][:] |
|
320
|
|
|
except: |
|
321
|
|
|
pass |
|
322
|
|
|
# apparent local solar time |
|
323
|
|
|
try: |
|
324
|
|
|
self.lst = ncf.variables['app_lst'][:] |
|
325
|
|
|
except: |
|
326
|
|
|
pass |
|
327
|
|
|
# mean local solar time |
|
328
|
|
|
try: |
|
329
|
|
|
self.mst = ncf.variables['mean_lst'][:] |
|
330
|
|
|
except: |
|
331
|
|
|
pass |
|
332
|
|
|
# mean solar zenith angle |
|
333
|
|
|
try: |
|
334
|
|
|
self.sza = ncf.variables['mean_sza'][:] |
|
335
|
|
|
except: |
|
336
|
|
|
pass |
|
337
|
|
|
# utc hours |
|
338
|
|
|
try: |
|
339
|
|
|
self.utchour = ncf.variables['utc_hour'][:] |
|
340
|
|
|
except: |
|
341
|
|
|
pass |
|
342
|
|
|
# utc days |
|
343
|
|
|
try: |
|
344
|
|
|
self.utcdays = ncf.variables['utc_days'][:] |
|
345
|
|
|
except: |
|
346
|
|
|
pass |
|
347
|
|
|
|
|
348
|
|
|
if close: |
|
349
|
|
|
ncf.close() |
|
350
|
|
|
else: |
|
351
|
|
|
return ncf |
|
352
|
|
|
|
|
353
|
1 |
|
def to_xarray(self, dateo, orbit): |
|
354
|
|
|
"""Convert the data to :class:`xarray.Dataset` |
|
355
|
|
|
|
|
356
|
|
|
This is a very simple approach, it dumps the data to a temporary |
|
357
|
|
|
netcdf file and reads that using :func:`xarray.open_dataset()`. |
|
358
|
|
|
|
|
359
|
|
|
Parameters |
|
360
|
|
|
---------- |
|
361
|
|
|
dateo: float |
|
362
|
|
|
The days since the reference date at the equator |
|
363
|
|
|
crossing of the orbit. Used to set the `time` |
|
364
|
|
|
dimension of the dataset. |
|
365
|
|
|
orbit: int |
|
366
|
|
|
The SCIAMACHY/Envisat orbit number of the retrieved data. |
|
367
|
|
|
|
|
368
|
|
|
Returns |
|
369
|
|
|
------- |
|
370
|
|
|
dataset: xarray.Dataset |
|
371
|
|
|
""" |
|
372
|
1 |
|
import tempfile |
|
373
|
1 |
|
try: |
|
374
|
1 |
|
import xarray as xr |
|
375
|
|
|
except ImportError: |
|
376
|
|
|
print("Error: xarray not available!") |
|
377
|
|
|
return None |
|
378
|
1 |
|
with tempfile.NamedTemporaryFile() as tf: |
|
379
|
1 |
|
self.write_to_netcdf(tf.name) |
|
380
|
1 |
|
with xr.open_dataset(tf.name, decode_cf=False) as sdorb: |
|
381
|
1 |
|
sdorb = sdorb.drop(["alt_min", "alt_max", "lat_min", "lat_max"]) |
|
382
|
1 |
|
sdorb["time"] = np.array([dateo], dtype=np.float64) |
|
383
|
1 |
|
sdorb["orbit"] = orbit |
|
384
|
1 |
|
sdorb.load() |
|
385
|
1 |
|
return sdorb |
|
386
|
|
|
|
|
387
|
|
|
|
|
388
|
1 |
|
class scia_density_day(object): |
|
389
|
|
|
"""SCIAMACHY daily number densities combined |
|
390
|
|
|
|
|
391
|
|
|
Contains a stacked version of the post-processed orbit data |
|
392
|
|
|
for multiple orbits on a day. Used to combine the results. |
|
393
|
|
|
|
|
394
|
|
|
Parameters |
|
395
|
|
|
---------- |
|
396
|
|
|
name: str |
|
397
|
|
|
Name of the retrieved species, default: "NO". |
|
398
|
|
|
Used to name the netcdf variables appropriately. |
|
399
|
|
|
ref_date: str, optional |
|
400
|
|
|
The reference date on which to base the date calculations on. |
|
401
|
|
|
Default: "2000-01-01" |
|
402
|
|
|
|
|
403
|
|
|
Attributes |
|
404
|
|
|
---------- |
|
405
|
|
|
alts |
|
406
|
|
|
Retrieval grid altitudes |
|
407
|
|
|
lats |
|
408
|
|
|
Retrieval grid latitudes |
|
409
|
|
|
no_dens |
|
410
|
|
|
Retrieved number densities |
|
411
|
|
|
no_errs |
|
412
|
|
|
Measurement uncertainty |
|
413
|
|
|
no_etot |
|
414
|
|
|
Total uncertainty |
|
415
|
|
|
no_rstd |
|
416
|
|
|
Relative measurement uncertainty |
|
417
|
|
|
no_akd |
|
418
|
|
|
Averaging kernel diagonal elements |
|
419
|
|
|
no_apri |
|
420
|
|
|
Prior number density |
|
421
|
|
|
temperature |
|
422
|
|
|
NRLMSISE-00 temperatures |
|
423
|
|
|
noem_no |
|
424
|
|
|
NOEM NO nuimber densities |
|
425
|
|
|
vmr |
|
426
|
|
|
NO vmr using the NRLMSISE-00 total air number densities |
|
427
|
|
|
lst |
|
428
|
|
|
Apparent local solar times |
|
429
|
|
|
mst |
|
430
|
|
|
Mean local solar times |
|
431
|
|
|
sza |
|
432
|
|
|
Solar zenith angles |
|
433
|
|
|
utchour |
|
434
|
|
|
UTC hours into measurement day |
|
435
|
|
|
utcdays |
|
436
|
|
|
Number of days since reference date |
|
437
|
|
|
gmlats |
|
438
|
|
|
IGRF-12 geomagentic latitudes |
|
439
|
|
|
gmlons |
|
440
|
|
|
IGRF-12 geomagentic longitudes |
|
441
|
|
|
aacgmgmlats |
|
442
|
|
|
AACGM geomagentic latitudes |
|
443
|
|
|
aacgmgmlons |
|
444
|
|
|
AACGM geomagentic longitudes |
|
445
|
|
|
""" |
|
446
|
1 |
|
def __init__(self, name="NO", ref_date="2000-01-01", author="unknown"): |
|
447
|
1 |
|
self.date0 = (dt.datetime.strptime(ref_date, "%Y-%m-%d") |
|
448
|
|
|
.replace(tzinfo=_UTC)) |
|
449
|
1 |
|
self.alts = None |
|
450
|
1 |
|
self.lats = None |
|
451
|
1 |
|
self.version = None |
|
452
|
1 |
|
self.data_version = None |
|
453
|
1 |
|
self.name = name |
|
454
|
1 |
|
self.author = author |
|
455
|
1 |
|
self.date = [] |
|
456
|
1 |
|
self.time = [] |
|
457
|
1 |
|
self.orbit = [] |
|
458
|
1 |
|
self.no_dens = None |
|
459
|
1 |
|
self.no_errs = None |
|
460
|
1 |
|
self.no_etot = None |
|
461
|
1 |
|
self.no_rstd = None |
|
462
|
1 |
|
self.no_akd = None |
|
463
|
1 |
|
self.no_apri = None |
|
464
|
1 |
|
self.no_noem = None |
|
465
|
1 |
|
self.temperature = None |
|
466
|
1 |
|
self.tot_dens = None |
|
467
|
1 |
|
self.no_vmr = None |
|
468
|
1 |
|
self.lons = None |
|
469
|
1 |
|
self.lst = None |
|
470
|
1 |
|
self.mst = None |
|
471
|
1 |
|
self.sza = None |
|
472
|
1 |
|
self.utchour = None |
|
473
|
1 |
|
self.utcdays = None |
|
474
|
1 |
|
self.gmlats = None |
|
475
|
1 |
|
self.gmlons = None |
|
476
|
1 |
|
self.aacgmgmlats = None |
|
477
|
1 |
|
self.aacgmgmlons = None |
|
478
|
|
|
|
|
479
|
1 |
|
def append(self, cdata): |
|
480
|
|
|
"""Append (stack) the data from one orbit |
|
481
|
|
|
|
|
482
|
|
|
Parameters |
|
483
|
|
|
---------- |
|
484
|
|
|
cdata: :class:`scia_densities_pp` instance |
|
485
|
|
|
Post-processed level 2 orbital data. |
|
486
|
|
|
""" |
|
487
|
|
|
self.time.extend(cdata.time) |
|
488
|
|
|
self.date.extend(cdata.date) |
|
489
|
|
|
self.no_dens = np.ma.dstack((self.no_dens, cdata.no_dens)) |
|
490
|
|
|
self.no_errs = np.ma.dstack((self.no_errs, cdata.no_errs)) |
|
491
|
|
|
self.no_etot = np.ma.dstack((self.no_etot, cdata.no_etot)) |
|
492
|
|
|
self.no_rstd = np.ma.dstack((self.no_rstd, cdata.no_rstd)) |
|
493
|
|
|
self.no_akd = np.ma.dstack((self.no_akd, cdata.no_akd)) |
|
494
|
|
|
self.no_apri = np.ma.dstack((self.no_apri, cdata.no_apri)) |
|
495
|
|
|
self.no_noem = np.ma.dstack((self.no_noem, cdata.no_noem)) |
|
496
|
|
|
self.tot_dens = np.ma.dstack((self.tot_dens, cdata.tot_dens)) |
|
497
|
|
|
self.no_vmr = np.ma.dstack((self.no_vmr, cdata.no_vmr)) |
|
498
|
|
|
self.lons = np.ma.dstack((self.lons, cdata.lons)) |
|
499
|
|
|
self.lst = np.ma.dstack((self.lst, cdata.lst)) |
|
500
|
|
|
self.mst = np.ma.dstack((self.mst, cdata.mst)) |
|
501
|
|
|
self.sza = np.ma.dstack((self.sza, cdata.sza)) |
|
502
|
|
|
self.utchour = np.ma.dstack((self.utchour, cdata.utchour)) |
|
503
|
|
|
self.utcdays = np.ma.dstack((self.utcdays, cdata.utcdays)) |
|
504
|
|
|
self.gmlats = np.ma.dstack((self.gmlats, cdata.gmlats)) |
|
505
|
|
|
self.gmlons = np.ma.dstack((self.gmlons, cdata.gmlons)) |
|
506
|
|
|
self.aacgmgmlats = np.ma.dstack((self.aacgmgmlats, cdata.aacgmgmlats)) |
|
507
|
|
|
self.aacgmgmlons = np.ma.dstack((self.aacgmgmlons, cdata.aacgmgmlons)) |
|
508
|
|
|
|
|
509
|
1 |
|
def append_data(self, date, orbit, equtime, scia_dens): |
|
510
|
|
|
"""Append (stack) the data from one orbit |
|
511
|
|
|
|
|
512
|
|
|
Updates the data in place. |
|
513
|
|
|
|
|
514
|
|
|
Parameters |
|
515
|
|
|
---------- |
|
516
|
|
|
date: float |
|
517
|
|
|
Days since `ref_date` for the time coordinate |
|
518
|
|
|
orbit: int |
|
519
|
|
|
SCIAMACHY/Envisat orbit number |
|
520
|
|
|
equtime: float |
|
521
|
|
|
UTC hour into the day at the equator |
|
522
|
|
|
scia_dens: :class:`scia_densities_pp` instance |
|
523
|
|
|
The post-processed orbit data set |
|
524
|
|
|
""" |
|
525
|
1 |
|
def _vstack_or_new(a, b): |
|
526
|
|
|
# Check if we 'stack' for the first time (a is None), |
|
527
|
|
|
# in that case we assign first. |
|
528
|
1 |
|
if a is None: |
|
529
|
1 |
|
return b[None] |
|
530
|
1 |
|
return np.ma.vstack((a, b[None])) |
|
531
|
|
|
|
|
532
|
1 |
|
self.version = scia_dens.version |
|
533
|
1 |
|
self.data_version = scia_dens.data_version |
|
534
|
1 |
|
self.date.append(date) |
|
535
|
1 |
|
self.orbit.append(orbit) |
|
536
|
1 |
|
self.time.append(equtime) |
|
537
|
1 |
|
_dens = scia_dens.densities |
|
538
|
1 |
|
_errs = scia_dens.dens_err_meas |
|
539
|
1 |
|
_etot = scia_dens.dens_err_tot |
|
540
|
1 |
|
_r_std = np.abs(_errs / _dens) * 100.0 |
|
541
|
1 |
|
if self.alts is None: |
|
542
|
|
|
# we need altitudes and latitudes only once |
|
543
|
1 |
|
self.alts = scia_dens.alts |
|
544
|
1 |
|
self.lats = scia_dens.lats |
|
545
|
1 |
|
self.no_dens = _vstack_or_new(self.no_dens, _dens) |
|
546
|
1 |
|
self.no_errs = _vstack_or_new(self.no_errs, _errs) |
|
547
|
1 |
|
self.no_etot = _vstack_or_new(self.no_etot, _etot) |
|
548
|
1 |
|
self.no_rstd = _vstack_or_new(self.no_rstd, _r_std) |
|
549
|
1 |
|
self.no_akd = _vstack_or_new(self.no_akd, scia_dens.akdiag) |
|
550
|
1 |
|
self.no_apri = _vstack_or_new(self.no_apri, scia_dens.apriori) |
|
551
|
1 |
|
self.temperature = _vstack_or_new(self.temperature, scia_dens.temperature) |
|
552
|
1 |
|
self.no_noem = _vstack_or_new(self.no_noem, scia_dens.noem_no) |
|
553
|
1 |
|
self.tot_dens = _vstack_or_new(self.tot_dens, scia_dens.dens_tot) |
|
554
|
1 |
|
self.no_vmr = _vstack_or_new(self.no_vmr, scia_dens.vmr) |
|
555
|
1 |
|
self.lons = _vstack_or_new(self.lons, scia_dens.lons) |
|
556
|
1 |
|
self.lst = _vstack_or_new(self.lst, scia_dens.lst) |
|
557
|
1 |
|
self.mst = _vstack_or_new(self.mst, scia_dens.mst) |
|
558
|
1 |
|
self.sza = _vstack_or_new(self.sza, scia_dens.sza) |
|
559
|
1 |
|
self.utchour = _vstack_or_new(self.utchour, scia_dens.utchour) |
|
560
|
1 |
|
self.utcdays = _vstack_or_new(self.utcdays, scia_dens.utcdays) |
|
561
|
1 |
|
self.gmlats = _vstack_or_new(self.gmlats, scia_dens.gmlats) |
|
562
|
1 |
|
self.gmlons = _vstack_or_new(self.gmlons, scia_dens.gmlons) |
|
563
|
1 |
|
self.aacgmgmlats = _vstack_or_new(self.aacgmgmlats, scia_dens.aacgmgmlats) |
|
564
|
1 |
|
self.aacgmgmlons = _vstack_or_new(self.aacgmgmlons, scia_dens.aacgmgmlons) |
|
565
|
|
|
|
|
566
|
1 |
|
def write_to_netcdf(self, filename): |
|
567
|
|
|
"""Write variables to netcdf files |
|
568
|
|
|
|
|
569
|
|
|
Parameters |
|
570
|
|
|
---------- |
|
571
|
|
|
filename: str |
|
572
|
|
|
The name of the file to write the data to. |
|
573
|
|
|
""" |
|
574
|
|
|
_var_dicts = { |
|
575
|
|
|
"2.1": { |
|
576
|
|
|
"dens_tot": { |
|
577
|
|
|
"name": "TOT_DENS", |
|
578
|
|
|
"long_name": "total number density (NRLMSIS-00)", |
|
579
|
|
|
"model": None, |
|
580
|
|
|
}, |
|
581
|
|
|
"temperature": { |
|
582
|
|
|
"name": "temperature", |
|
583
|
|
|
"long_name": "temperature", |
|
584
|
|
|
"model": "NRLMSIS-00", |
|
585
|
|
|
}, |
|
586
|
|
|
}, |
|
587
|
|
|
"2.2": { |
|
588
|
|
|
"dens_tot": { |
|
589
|
|
|
"name": "MSIS_Dens", |
|
590
|
|
|
"long_name": "MSIS total number density", |
|
591
|
|
|
"model": "NRLMSIS-00", |
|
592
|
|
|
}, |
|
593
|
|
|
"temperature": { |
|
594
|
|
|
"name": "MSIS_Temp", |
|
595
|
|
|
"long_name": "MSIS temperature", |
|
596
|
|
|
"model": "NRLMSIS-00", |
|
597
|
|
|
}, |
|
598
|
|
|
}, |
|
599
|
|
|
} |
|
600
|
|
|
|
|
601
|
|
|
ncf = netcdf_file(filename, 'w', **fmtargs) |
|
602
|
|
|
o = np.asarray(self.orbit) |
|
603
|
|
|
d = np.asarray(self.date) |
|
604
|
|
|
t = np.asarray(self.time) |
|
605
|
|
|
|
|
606
|
|
|
if self.version is not None: |
|
607
|
|
|
ncf.version = self.version |
|
608
|
|
|
if self.data_version is not None: |
|
609
|
|
|
ncf.L2_data_version = self.data_version |
|
610
|
|
|
ncf.software = "sciapy {0}".format(__version__) |
|
611
|
|
|
ncf.creation_time = dt.datetime.utcnow().strftime("%a %b %d %Y %H:%M:%S +00:00 (UTC)") |
|
612
|
|
|
ncf.author = self.author |
|
613
|
|
|
|
|
614
|
|
|
ncf.createDimension('time', None) |
|
615
|
|
|
ncf.createDimension('altitude', self.alts.size) |
|
616
|
|
|
ncf.createDimension('latitude', self.lats.size) |
|
617
|
|
|
forbit = ncf.createVariable('orbit', np.dtype('int32').char, ('time',)) |
|
618
|
|
|
forbit.axis = 'T' |
|
619
|
|
|
forbit.calendar = 'standard' |
|
620
|
|
|
forbit.long_name = 'orbit' |
|
621
|
|
|
forbit.standard_name = 'orbit' |
|
622
|
|
|
forbit.units = 'orbit number' |
|
623
|
|
|
# the time coordinate is actually called "date" here |
|
624
|
|
|
ftime = ncf.createVariable('time', 'f8', ('time',)) |
|
625
|
|
|
ftime.axis = 'T' |
|
626
|
|
|
ftime.calendar = 'standard' |
|
627
|
|
|
ftime.long_name = 'equatorial crossing time' |
|
628
|
|
|
ftime.standard_name = 'time' |
|
629
|
|
|
ftime.units = 'days since {0}'.format(self.date0.isoformat(sep=' ')) |
|
630
|
|
|
#ftime.units = 'days since {0}'.format(self.date0.strftime('%Y-%m-%d %H:%M:%S%z (%Z)')) |
|
631
|
|
|
#fdate = ncf.createVariable('date', np.dtype('float64').char, ('time',)) |
|
632
|
|
|
#fdate.axis = 'T' |
|
633
|
|
|
#fdate.calendar = 'standard' |
|
634
|
|
|
#fdate.long_name = 'date' |
|
635
|
|
|
#fdate.standard_name = 'date' |
|
636
|
|
|
#fdate.units = 'days since 1950-01-01 00:00:00' |
|
637
|
|
|
falts = ncf.createVariable('altitude', 'f8', ('altitude',)) |
|
638
|
|
|
falts.axis = 'Z' |
|
639
|
|
|
falts.long_name = 'altitude' |
|
640
|
|
|
falts.standard_name = 'altitude' |
|
641
|
|
|
falts.units = 'km' |
|
642
|
|
|
falts.positive = 'up' |
|
643
|
|
|
flats = ncf.createVariable('latitude', 'f8', ('latitude',)) |
|
644
|
|
|
flats.axis = 'Y' |
|
645
|
|
|
flats.long_name = 'latitude' |
|
646
|
|
|
flats.standard_name = 'latitude' |
|
647
|
|
|
flats.units = 'degrees_north' |
|
648
|
|
|
flons = ncf.createVariable('longitude', 'f8', ('time', 'latitude',)) |
|
649
|
|
|
flons.long_name = 'longitude' |
|
650
|
|
|
flons.standard_name = 'longitude' |
|
651
|
|
|
flons.units = 'degrees_east' |
|
652
|
|
|
fdens = ncf.createVariable('%s_DENS' % self.name, 'f8', ('time', 'latitude', 'altitude')) |
|
653
|
|
|
fdens.units = 'cm^{-3}' |
|
654
|
|
|
fdens.long_name = '%s number density' % self.name |
|
655
|
|
|
ferrs = ncf.createVariable('%s_ERR' % self.name, 'f8', ('time', 'latitude', 'altitude')) |
|
656
|
|
|
ferrs.units = 'cm^{-3}' |
|
657
|
|
|
ferrs.long_name = '%s density measurement error' % self.name |
|
658
|
|
|
fetot = ncf.createVariable('%s_ETOT' % self.name, 'f8', ('time', 'latitude', 'altitude')) |
|
659
|
|
|
fetot.units = 'cm^{-3}' |
|
660
|
|
|
fetot.long_name = '%s density total error' % self.name |
|
661
|
|
|
frstd = ncf.createVariable('%s_RSTD' % self.name, 'f8', ('time', 'latitude', 'altitude')) |
|
662
|
|
|
frstd.units = '%' |
|
663
|
|
|
frstd.long_name = '%s relative standard deviation' % self.name |
|
664
|
|
|
fakd = ncf.createVariable('%s_AKDIAG' % self.name, 'f8', ('time', 'latitude', 'altitude')) |
|
665
|
|
|
fakd.units = '1' |
|
666
|
|
|
fakd.long_name = '%s averaging kernel diagonal element' % self.name |
|
667
|
|
|
fapri = ncf.createVariable('%s_APRIORI' % self.name, 'f8', ('time', 'latitude', 'altitude')) |
|
668
|
|
|
fapri.units = 'cm^{-3}' |
|
669
|
|
|
fapri.long_name = '%s apriori density' % self.name |
|
670
|
|
|
ftemp = ncf.createVariable(_var_dicts[self.version]["temperature"]["name"], |
|
671
|
|
|
'f8', ('time', 'latitude', 'altitude')) |
|
672
|
|
|
ftemp.long_name = _var_dicts[self.version]["temperature"]["long_name"] |
|
673
|
|
|
ftemp.model = 'NRLMSIS-00' |
|
674
|
|
|
ftemp.units = 'K' |
|
675
|
|
|
fdens_tot = ncf.createVariable(_var_dicts[self.version]["dens_tot"]["name"], |
|
676
|
|
|
'f8', ('time', 'latitude', 'altitude')) |
|
677
|
|
|
fdens_tot.long_name = _var_dicts[self.version]["dens_tot"]["long_name"] |
|
678
|
|
|
fdens_tot.units = 'cm^{-3}' |
|
679
|
|
|
if _var_dicts[self.version]["dens_tot"]["model"] is not None: |
|
680
|
|
|
fdens_tot.model = _var_dicts[self.version]["dens_tot"]["model"] |
|
681
|
|
|
fnoem = ncf.createVariable('%s_NOEM' % self.name, 'f8', ('time', 'latitude', 'altitude')) |
|
682
|
|
|
fnoem.units = 'cm^{-3}' |
|
683
|
|
|
fnoem.long_name = 'NOEM %s number density' % self.name |
|
684
|
|
|
fvmr = ncf.createVariable('%s_VMR' % self.name, 'f8', ('time', 'latitude', 'altitude')) |
|
685
|
|
|
fvmr.units = 'ppb' |
|
686
|
|
|
fvmr.long_name = '%s volume mixing ratio' % self.name |
|
687
|
|
|
flst = ncf.createVariable('app_LST', 'f8', ('time', 'latitude')) |
|
688
|
|
|
flst.units = 'hours' |
|
689
|
|
|
flst.long_name = 'apparent local solar time' |
|
690
|
|
|
fmst = ncf.createVariable('mean_LST', 'f8', ('time', 'latitude')) |
|
691
|
|
|
fmst.units = 'hours' |
|
692
|
|
|
fmst.long_name = 'mean local solar time' |
|
693
|
|
|
fsza = ncf.createVariable('mean_SZA', 'f8', ('time', 'latitude')) |
|
694
|
|
|
fsza.units = 'degrees' |
|
695
|
|
|
fsza.long_name = 'solar zenith angle at mean altitude' |
|
696
|
|
|
futc = ncf.createVariable('UTC', 'f8', ('time', 'latitude')) |
|
697
|
|
|
futc.units = 'hours' |
|
698
|
|
|
futc.long_name = 'measurement utc time' |
|
699
|
|
|
futcd = ncf.createVariable('utc_days', 'f8', ('time', 'latitude')) |
|
700
|
|
|
futcd.long_name = 'measurement utc day' |
|
701
|
|
|
futcd.units = 'days since {0}'.format(self.date0.isoformat(sep=' ')) |
|
702
|
|
|
|
|
703
|
|
|
fgmlats = ncf.createVariable('gm_lats', 'f8', ('time', 'latitude',)) |
|
704
|
|
|
fgmlats.long_name = 'geomagnetic_latitude' |
|
705
|
|
|
fgmlats.model = 'IGRF' |
|
706
|
|
|
fgmlats.units = 'degrees_north' |
|
707
|
|
|
|
|
708
|
|
|
fgmlons = ncf.createVariable('gm_lons', 'f8', ('time', 'latitude',)) |
|
709
|
|
|
fgmlons.long_name = 'geomagnetic_longitude' |
|
710
|
|
|
fgmlons.model = 'IGRF' |
|
711
|
|
|
fgmlons.units = 'degrees_east' |
|
712
|
|
|
|
|
713
|
|
|
faacgmgmlats = ncf.createVariable('aacgm_gm_lats', 'f8', ('time', 'latitude',)) |
|
714
|
|
|
faacgmgmlats.long_name = 'geomagnetic_latitude' |
|
715
|
|
|
faacgmgmlats.model = 'AACGM' |
|
716
|
|
|
faacgmgmlats.units = 'degrees_north' |
|
717
|
|
|
|
|
718
|
|
|
faacgmgmlons = ncf.createVariable('aacgm_gm_lons', 'f8', ('time', 'latitude',)) |
|
719
|
|
|
faacgmgmlons.long_name = 'geomagnetic_longitude' |
|
720
|
|
|
faacgmgmlons.model = 'AACGM' |
|
721
|
|
|
faacgmgmlons.units = 'degrees_east' |
|
722
|
|
|
|
|
723
|
|
|
forbit[:] = o |
|
724
|
|
|
ftime[:] = d |
|
725
|
|
|
falts[:] = self.alts |
|
726
|
|
|
flats[:] = self.lats |
|
727
|
|
|
flons[:] = self.lons |
|
728
|
|
|
fdens[:] = np.ma.atleast_3d(self.no_dens[:]) |
|
729
|
|
|
ferrs[:] = np.ma.atleast_3d(self.no_errs[:]) |
|
730
|
|
|
fetot[:] = np.ma.atleast_3d(self.no_etot[:]) |
|
731
|
|
|
frstd[:] = np.ma.atleast_3d(self.no_rstd[:]) |
|
732
|
|
|
fakd[:] = np.ma.atleast_3d(self.no_akd[:]) |
|
733
|
|
|
fapri[:] = np.ma.atleast_3d(self.no_apri[:]) |
|
734
|
|
|
ftemp[:] = np.ma.atleast_3d(self.temperature[:]) |
|
735
|
|
|
fnoem[:] = np.ma.atleast_3d(self.no_noem[:]) |
|
736
|
|
|
fdens_tot[:] = np.ma.atleast_3d(self.tot_dens[:]) |
|
737
|
|
|
fvmr[:] = np.ma.atleast_3d(self.no_vmr[:]) |
|
738
|
|
|
flst[:] = np.ma.atleast_2d(self.lst[:]) |
|
739
|
|
|
fmst[:] = np.ma.atleast_2d(self.mst[:]) |
|
740
|
|
|
fsza[:] = np.ma.atleast_2d(self.sza[:]) |
|
741
|
|
|
futc[:] = np.ma.atleast_2d(self.utchour[:]) |
|
742
|
|
|
futcd[:] = np.ma.atleast_2d(self.utcdays[:]) |
|
743
|
|
|
fgmlats[:] = np.ma.atleast_2d(self.gmlats[:]) |
|
744
|
|
|
fgmlons[:] = np.ma.atleast_2d(self.gmlons[:]) |
|
745
|
|
|
faacgmgmlats[:] = np.ma.atleast_2d(self.aacgmgmlats[:]) |
|
746
|
|
|
faacgmgmlons[:] = np.ma.atleast_2d(self.aacgmgmlons[:]) |
|
747
|
|
|
ncf.close() |
|
748
|
|
|
|
|
749
|
1 |
|
def to_xarray(self): |
|
750
|
|
|
"""Convert the combined orbit data to :class:`xarray.Dataset` |
|
751
|
|
|
|
|
752
|
|
|
Exports the data using the same data variables as |
|
753
|
|
|
when writing to netcdf. |
|
754
|
|
|
|
|
755
|
|
|
Returns |
|
756
|
|
|
------- |
|
757
|
|
|
dataset: xarray.Dataset |
|
758
|
|
|
""" |
|
759
|
1 |
|
try: |
|
760
|
1 |
|
import xarray as xr |
|
761
|
|
|
except ImportError: |
|
762
|
|
|
print("Error: xarray not available!") |
|
763
|
|
|
return None |
|
764
|
1 |
|
o = np.asarray(self.orbit) |
|
765
|
1 |
|
d = np.asarray(self.date) |
|
766
|
|
|
|
|
767
|
1 |
|
xr_dens = xr.DataArray(self.no_dens, coords=[d, self.lats, self.alts], |
|
768
|
|
|
dims=["time", "latitude", "altitude"], |
|
769
|
|
|
attrs={"units": "cm^{-3}", |
|
770
|
|
|
"long_name": "{0} number density".format(self.name)}, |
|
771
|
|
|
name="{0}_DENS".format(self.name)) |
|
772
|
|
|
|
|
773
|
1 |
|
xr_errs = xr.DataArray(self.no_errs, coords=[d, self.lats, self.alts], |
|
774
|
|
|
dims=["time", "latitude", "altitude"], |
|
775
|
|
|
attrs={"units": "cm^{-3}", |
|
776
|
|
|
"long_name": "{0} density measurement error".format(self.name)}, |
|
777
|
|
|
name="{0}_ERR".format(self.name)) |
|
778
|
|
|
|
|
779
|
1 |
|
xr_etot = xr.DataArray(self.no_etot, coords=[d, self.lats, self.alts], |
|
780
|
|
|
dims=["time", "latitude", "altitude"], |
|
781
|
|
|
attrs={"units": "cm^{-3}", |
|
782
|
|
|
"long_name": "{0} density total error".format(self.name)}, |
|
783
|
|
|
name="{0}_ETOT".format(self.name)) |
|
784
|
|
|
|
|
785
|
1 |
|
xr_rstd = xr.DataArray(self.no_rstd, coords=[d, self.lats, self.alts], |
|
786
|
|
|
dims=["time", "latitude", "altitude"], |
|
787
|
|
|
attrs=dict(units='%', |
|
788
|
|
|
long_name='{0} relative standard deviation'.format(self.name)), |
|
789
|
|
|
name="{0}_RSTD".format(self.name)) |
|
790
|
|
|
|
|
791
|
1 |
|
xr_akd = xr.DataArray(self.no_akd, coords=[d, self.lats, self.alts], |
|
792
|
|
|
dims=["time", "latitude", "altitude"], |
|
793
|
|
|
attrs=dict(units='1', |
|
794
|
|
|
long_name='{0} averaging kernel diagonal element'.format(self.name)), |
|
795
|
|
|
name="{0}_AKDIAG".format(self.name)) |
|
796
|
|
|
|
|
797
|
1 |
|
xr_apri = xr.DataArray(self.no_apri, coords=[d, self.lats, self.alts], |
|
798
|
|
|
dims=["time", "latitude", "altitude"], |
|
799
|
|
|
attrs=dict(units='cm^{-3}', long_name='{0} apriori density'.format(self.name)), |
|
800
|
|
|
name="{0}_APRIORI".format(self.name)) |
|
801
|
|
|
|
|
802
|
1 |
|
xr_noem = xr.DataArray(self.no_noem, coords=[d, self.lats, self.alts], |
|
803
|
|
|
dims=["time", "latitude", "altitude"], |
|
804
|
|
|
attrs=dict(units='cm^{-3}', long_name='NOEM {0} number density'.format(self.name)), |
|
805
|
|
|
name="{0}_NOEM".format(self.name)) |
|
806
|
|
|
|
|
807
|
1 |
|
xr_vmr = xr.DataArray(self.no_vmr, coords=[d, self.lats, self.alts], |
|
808
|
|
|
dims=["time", "latitude", "altitude"], |
|
809
|
|
|
attrs=dict(units='ppb', long_name='{0} volume mixing ratio'.format(self.name)), |
|
810
|
|
|
name="{0}_VMR".format(self.name)) |
|
811
|
|
|
|
|
812
|
1 |
|
xr_dtot = xr.DataArray(self.tot_dens, coords=[d, self.lats, self.alts], |
|
813
|
|
|
dims=["time", "latitude", "altitude"], |
|
814
|
|
|
attrs=dict(units='cm^{-3}', long_name='MSIS total number density', |
|
815
|
|
|
model='NRLMSIS-00'), |
|
816
|
|
|
name="MSIS_Dens") |
|
817
|
|
|
|
|
818
|
1 |
|
xr_temp = xr.DataArray(self.temperature, coords=[d, self.lats, self.alts], |
|
819
|
|
|
dims=["time", "latitude", "altitude"], |
|
820
|
|
|
attrs=dict(units='K', long_name='MSIS temperature', |
|
821
|
|
|
model="NRLMSIS-00"), |
|
822
|
|
|
name="MSIS_Temp") |
|
823
|
|
|
|
|
824
|
1 |
|
xr_lons = xr.DataArray(self.lons, coords=[d, self.lats], |
|
825
|
|
|
dims=["time", "latitude"], |
|
826
|
|
|
attrs=dict(long_name='longitude', standard_name='longitude', |
|
827
|
|
|
units='degrees_east'), |
|
828
|
|
|
name='longitude') |
|
829
|
|
|
|
|
830
|
1 |
|
xr_lst = xr.DataArray(self.lst, coords=[d, self.lats], |
|
831
|
|
|
dims=["time", "latitude"], |
|
832
|
|
|
attrs=dict(units='hours', long_name='apparent local solar time'), |
|
833
|
|
|
name="app_LST") |
|
834
|
|
|
|
|
835
|
1 |
|
xr_mst = xr.DataArray(self.mst, coords=[d, self.lats], |
|
836
|
|
|
dims=["time", "latitude"], |
|
837
|
|
|
attrs=dict(units='hours', long_name='mean local solar time'), |
|
838
|
|
|
name="mean_LST") |
|
839
|
|
|
|
|
840
|
1 |
|
xr_sza = xr.DataArray(self.sza, coords=[d, self.lats], |
|
841
|
|
|
dims=["time", "latitude"], |
|
842
|
|
|
attrs=dict(units='degrees', |
|
843
|
|
|
long_name='solar zenith angle at mean altitude'), |
|
844
|
|
|
name="mean_SZA") |
|
845
|
|
|
|
|
846
|
1 |
|
xr_utch = xr.DataArray(self.utchour, coords=[d, self.lats], |
|
847
|
|
|
dims=["time", "latitude"], |
|
848
|
|
|
attrs=dict(units='hours', |
|
849
|
|
|
long_name='measurement utc time'), |
|
850
|
|
|
name="UTC") |
|
851
|
|
|
|
|
852
|
1 |
|
xr_utcd = xr.DataArray(self.utcdays, coords=[d, self.lats], |
|
853
|
|
|
dims=["time", "latitude"], |
|
854
|
|
|
attrs=dict(units='days since {0}'.format(self.date0.isoformat(sep=' ')), |
|
855
|
|
|
long_name='measurement utc day'), |
|
856
|
|
|
name="utc_days") |
|
857
|
|
|
|
|
858
|
1 |
|
xr_gmlats = xr.DataArray(self.gmlats, coords=[d, self.lats], |
|
859
|
|
|
dims=["time", "latitude"], |
|
860
|
|
|
attrs=dict(long_name='geomagnetic_latitude', |
|
861
|
|
|
model='IGRF', units='degrees_north'), |
|
862
|
|
|
name="gm_lats") |
|
863
|
|
|
|
|
864
|
1 |
|
xr_gmlons = xr.DataArray(self.gmlons, coords=[d, self.lats], |
|
865
|
|
|
dims=["time", "latitude"], |
|
866
|
|
|
attrs=dict(long_name='geomagnetic_longitude', |
|
867
|
|
|
model='IGRF', units='degrees_east'), |
|
868
|
|
|
name="gm_lons") |
|
869
|
|
|
|
|
870
|
1 |
|
xr_aacgmgmlats = xr.DataArray(self.aacgmgmlats, coords=[d, self.lats], |
|
871
|
|
|
dims=["time", "latitude"], |
|
872
|
|
|
attrs=dict(long_name='geomagnetic_latitude', |
|
873
|
|
|
model='AACGM', units='degrees_north'), |
|
874
|
|
|
name="aacgm_gm_lats") |
|
875
|
|
|
|
|
876
|
1 |
|
xr_aacgmgmlons = xr.DataArray(self.aacgmgmlons, coords=[d, self.lats], |
|
877
|
|
|
dims=["time", "latitude"], |
|
878
|
|
|
attrs=dict(long_name='geomagnetic_longitude', |
|
879
|
|
|
model='AACGM', units='degrees_east'), |
|
880
|
|
|
name="aacgm_gm_lons") |
|
881
|
|
|
|
|
882
|
1 |
|
xr_orbit = xr.DataArray(o, coords=[d], dims=["time"], |
|
883
|
|
|
attrs=dict(axis='T', calendar='standard', long_name='orbit', |
|
884
|
|
|
standard_name='orbit', units='orbit number'), |
|
885
|
|
|
name="orbit") |
|
886
|
|
|
|
|
887
|
1 |
|
xr_ds = xr.Dataset({da.name: da for da in |
|
888
|
|
|
[xr_dens, xr_errs, xr_etot, xr_rstd, xr_akd, xr_apri, xr_noem, |
|
889
|
|
|
xr_vmr, xr_dtot, xr_temp, xr_lons, xr_lst, xr_mst, xr_sza, |
|
890
|
|
|
xr_utch, xr_utcd, xr_gmlats, xr_gmlons, xr_aacgmgmlats, |
|
891
|
|
|
xr_aacgmgmlons, xr_orbit]}) |
|
892
|
|
|
|
|
893
|
1 |
|
xr_ds["time"].attrs = dict(axis='T', calendar='standard', |
|
894
|
|
|
long_name='equatorial crossing time', |
|
895
|
|
|
standard_name='time', |
|
896
|
|
|
units='days since {0}'.format(self.date0.isoformat(sep=' '))) |
|
897
|
|
|
|
|
898
|
1 |
|
xr_ds["altitude"].attrs = dict(axis='Z', long_name='altitude', |
|
899
|
|
|
standard_name='altitude', units='km', positive='up') |
|
900
|
|
|
|
|
901
|
1 |
|
xr_ds["latitude"].attrs = dict(axis='Y', long_name='latitude', |
|
902
|
|
|
standard_name='latitude', units='degrees_north') |
|
903
|
|
|
|
|
904
|
1 |
|
if self.version is not None: |
|
905
|
1 |
|
xr_ds.attrs["version"] = self.version |
|
906
|
1 |
|
if self.data_version is not None: |
|
907
|
1 |
|
xr_ds.attrs["L2_data_version"] = self.data_version |
|
908
|
1 |
|
xr_ds.attrs["software"] = "sciapy {0}".format(__version__) |
|
909
|
1 |
|
xr_ds.attrs["creation_time"] = dt.datetime.utcnow().strftime("%a %b %d %Y %H:%M:%S +00:00 (UTC)") |
|
910
|
1 |
|
xr_ds.attrs["author"] = self.author |
|
911
|
|
|
|
|
912
|
|
|
return xr_ds |
|
913
|
|
|
|
st-bender /
sciapy
