'''The module provide access data through the official ELS data set for VEX.
See step by step :ref:`tutorial <stepbystep>` for data analysis.
.. todo::
More user-oriented data analysis explanation should be here.
.. todo::
Class relation figure to be added in the document.
This is a model of implementation of data accessor.
.. todo::
Consider the mixature of ELS 128 energy mode and 32 energy mode.
This may influence the basic implementaion, thus, this change might be
implemented as scidata2 module.
- Data getter:
Simple methods to get data.
:meth:`get_countdata`, :meth:`get_count`, :meth:`get_def`, :meth:`get_dnf` and :meth:`get_psd` is implemented.
Partly supports 32 energy step mode. But all the returned values are matrix
of 128 elements, which means that the user must reshape the data manually
into 32 element array.
.. note::
It is important to know what time of the data is returned.
If user specified the time of *t*, then the returned data is at the time *t0*, where
*t0* is the latest time satisfying *t0 < t*.
Thus, if you give time in the data gap between observation sessions,
then the returned time is the last data of the previous observation session.
User should be careful on this behaviour.
Introduced from v4.2.3a6.
- Observation time getter:
Obtain the time of the observation. Each 4s time returned. :meth:`get_obstime`.
- Single matlab file:
A single file of the ELS matlab-formatted file is expressed by
:class:`ElsCountFile`.
Each single file contains 1 hour data of ELS.
- Matlab data base:
An aggregation of the ELS matlab-formatted file is called data base.
The data base is implemented in :class:`ElsCountFileDB`.
The method :meth:`ElsCountFileDB.get_filename` returns the filename
that includes the data of the specified time.
.. note::
This is not for user functionality. Document to be removed.
- Default Matlab data base:
The method :meth:`get_defaultdb` will return the default data base
specifieid in ``.irfpyrc``.
.. note::
This is not for user functionality. Document revised to adopt to the user oriented way.
- ELS data at a time:
:class:`ElsCountData` is a *record* of the data. It is produced
from single matlab file with specification of the time using
the method :meth:`ElsCountFile.get_elscountdata`
Simple accessor is :meth:`get_countdata` function.
'''
import os
import datetime
import logging
import warnings
import bisect
import numpy as np
import scipy.io
from irfpy.asperacommon.excepts import DbPathNotFoundError, DatafileNotFoundError
from irfpy.vels import mode
from irfpy.vels import energy
from irfpy.vels.bible import energy as bible_energy
from irfpy.vels.bible import flux as bible_flux
from irfpy.vels.bible import kfactor as bible_kfactor
import irfpy.util.timeseriesdb
import irfpy.util.ringcache
import irfpy.util.irfpyrc
from irfpy.util import utc
[docs]class ElsCountData:
''' ELS count data at a specified time (4-s data with full energy scan). No interpretation.
ELS data consists of the following information.
*Start time* (``starttime``)
The start time of the observation. ``datetime.datetime`` instance.
*Energy spectra* (``count``)
The energy spectra as (nE, nA) array, where nE is the energy step and nA is the anode size.
*HV levels* (``level``)
HV levels to calculate the energy more precisely than the prediction.
``None`` is also allowed. ``None`` may happen at the boundary of the data files.
This class is for raw, non-interpreted data. Non-interpreted data means that the mode is
unknown. Therefore, this class can represent both E128 mode and E32 mode, but the
four E32 mode data is stacked into one :class:`ElsCountData` object.
This calss does not support the energy table.
It is recommended to use the interpreted data (:class:`ElsCountDataE128` and :class:`ElsCountDataE32`)
for data analysis
'''
def __init__(self, t, count, level, level_time=None):
if not isinstance(t, datetime.datetime):
raise TypeError('Wrong argument (0), should be datetime instance')
try:
if count.shape[1] != 16:
raise ValueError('Wrong shape: {}'.format(str(count.shape)))
except AttributeError:
raise TypeError('Wrong argument. Count should be numpy array')
except IndexError:
raise ValueError('Wrong dimension. The cound data should be (Ne, 16) shape.')
if level is not None:
try:
if level.shape[0] != count.shape[0]:
raise ValueError('Inconsistent shape in level and in energy')
except AttributeError:
raise TypeError('Wrong argument. Level should be numpy array')
self.starttime = t
self.count = count.astype(None)
self.level = level
self.leveltime = level_time
if self.level is not None:
self.level = self.level.astype(None)
[docs] def hasLevel(self):
return (self.level is not None)
def __str__(self):
if self.hasLevel():
if self.leveltime is None:
islevel = 'Defined (@???)'
else:
islevel = 'Defined (@{:%FT%T})'.format(self.leveltime)
else:
islevel = 'Not defined'
s = '<%s at %s, count_max=%.1f, level=(%s)>' % (self.__class__.__name__, self.starttime.strftime("%FT%T"),
self.count.max(), islevel)
return s
[docs] def get_count(self):
''' Return the count rate data.
:returns: The count rate of ELS
:rtype: ``numpy.array`` with (128, 16) shape.
'''
return self.count.copy()
[docs] def get_energy(self, elsmode=None):
''' Return the energy table that is calculated from TM level.
:param elsmode: ELS mode to interpret this data.
See :meth:`ElsCountData.get_tm_level` for details.
:returns: Energy table in electron volts.
:rtype: ``numpy.array`` of (128, 16) shape
.. warning::
It is not recommended to use this method for data analysis.
Use :meth:`irfpy.vels.energy.get_default_table_128` and :meth:`irfpy.vels.energy.get_default_table_32`.
'''
kfact = bible_kfactor.kfactors()
tm = self.get_tm_level()
hv = [bible_energy.tm2ref(val) for val in tm]
etbl = np.zeros((128, 16))
for ie in range(128):
for id in range(16):
etbl[ie, id] = hv[ie] * kfact[id]
return etbl
[docs] def get_tm_level(self, elsmode=None):
''' Return the TM level in the data.
:param elsmode: ELS mode to interpret this data.
If ``None`` is specified, the ELS mode is guessed by guess_mode.
When guessing failed, default mode (:attr:`irfpy.vels.mode.MODE_E128`) is used.
If you want to explicitly specify the mode, you can use either :attr:`irfpy.vels.mode.MODE_E128`
of :attr:`irfpy.vels.mode.MODE_E32`. :attr:`irfpy.vels.mode.MODE_UNKNOWN` will be
interpreted as the default, i.e. :attr:`irfpy.vels.mode.MODE_E128`.
:returns: TM level. Low power supply corresponding to negative values as well as the trick
implied in :meth:`irfpy.vels.bible.energy.tm2ref` method.
'''
if elsmode is None:
elsmode =self.guess_mode()
if elsmode == mode.MODE_UNKNOWN:
elsmode = mode.MODE_E128
if elsmode == mode.MODE_E128:
tm = self.level.copy()
tm[60:] = -tm[60:] # First 60 steps is for high power supply, and the rest is for low.
elif elsmode == mode.MODE_E32:
tm = -self.level.copy() # All the steps are low power suppy.
else:
raise ValueError('Unknown mode (%s). Check the mode or specify with ``elsmode`` keyword.' % str(elsmode))
return tm
[docs] def guess_mode(self):
''' Guess the mode.
:returns: Mode instance
:rtype: An ELS mode constant defined in :mod:`irfpy.vels.mode`.
Guess the mode from the HV level.
See :func:`irfpy.vels.mode.guess_mode_from_level`.
'''
if self.level is None:
return mode.MODE_UNKNOWN
return mode.guess_mode_from_level(self.level)
[docs]def interpret_rawdata(elsdata, elsmode=None):
""" Interpret the data, and return the interpreted data.
The data is interpreted
:param elsmode: ELS mode to interpret this data.
See :meth:`ElsCountData.get_tm_level` for details.
:param elsdata: An object of :class:`ElsCountData`.
:return: A tuple. Either a tuple of 1-size :class:`irfpy.vels.ElsCountDataE128` object
or a tuple of 4-size :class:`irfpy.vels.ElsCountDataE32` objects.
"""
from irfpy.vels import mode
if elsmode is None:
elsmode = mode.MODE_UNKNOWN
if elsmode == mode.MODE_UNKNOWN:
elsmode = mode.guess_mode_from_level(elsdata.level)
logger = logging.getLogger(__name__)
logger.debug('MODE={}'.format(str(elsmode)))
if elsmode == mode.MODE_UNKNOWN:
elsmode = mode.MODE_E128 # Unknwon case is E128.
if elsmode == mode.MODE_E128:
return (ElsCountDataE128(elsdata.starttime, elsdata.count, elsdata.level), )
elif elsmode == mode.MODE_E32:
return _interpret_as_E32(elsdata)
else:
raise ValueError('Either MODE_E128 or MODE_E32 is supported. (ELSMODE={})'.format(elsmode))
def _interpret_as_E32(elsdata):
""" The ELS data is interpreted ad E32 sequence.
>>> elsdata = ElsCountData(datetime.datetime(2001, 3, 14), np.arange(128*16).reshape(128, 16), np.arange(128)) # Emulate ELS data
>>> print(elsdata)
<ElsCountData at 2001-03-14T00:00:00, count_max=2047.0, level=(Defined (@???))>
>>> elsdata32 = _interpret_as_E32(elsdata)
>>> print(elsdata32[0])
<ElsCountDataE32 at 2001-03-14T00:00:00, count_max=511.0, level=(Defined)>
>>> print(elsdata32[1])
<ElsCountDataE32 at 2001-03-14T00:00:01, count_max=1023.0, level=(Defined)>
>>> print(elsdata32[2])
<ElsCountDataE32 at 2001-03-14T00:00:02, count_max=1535.0, level=(Defined)>
>>> print(elsdata32[3])
<ElsCountDataE32 at 2001-03-14T00:00:03, count_max=2047.0, level=(Defined)>
"""
# Time should be separated into four.
t0, t1, t2, t3 = [elsdata.starttime + i * datetime.timedelta(seconds=1) for i in range(4)]
# Data should be separated into four. This could cause error for non-E128 data... but
c0 = elsdata.count[0:32]
c1 = elsdata.count[32:64]
c2 = elsdata.count[64:96]
c3 = elsdata.count[96:128]
# Level as well
if elsdata.level is None:
l0 = l1 = l2 = l3 = None
else:
l0 = elsdata.level[0:32]
l1 = elsdata.level[32:64]
l2 = elsdata.level[64:96]
l3 = elsdata.level[96:128]
return (ElsCountDataE32(t0, c0, l0),
ElsCountDataE32(t1, c1, l1),
ElsCountDataE32(t2, c2, l2),
ElsCountDataE32(t3, c3, l3),
)
# I could not find document but the following code does not provide convincing data... (2009-06-07 evening data)
# return (ElsCountDataE32(t0, c3, l3),
# ElsCountDataE32(t1, c2, l2),
# ElsCountDataE32(t2, c1, l1),
# ElsCountDataE32(t3, c0, l0),
# )
[docs]class ElsCountDataE128(ElsCountData):
''' ELS count data at a specific time (4-s resolution) of E128 mode.
'''
def __init__(self, t, count, level):
''' As the entry of matlab file has not analyzed, this is still under notimplemented
ELS data consists of the following information.
*Start time* (``starttime``)
The start time of the observation. ``datetime.datetime`` instance.
*Energy spectra* (``count``)
The energy spectra as (nE, nA) array, where nE is the energy step and nA is the anode size.
*HV levels* (``level``)
HV levels to calculate the energy more precisely than the prediction.
``None`` is also allowed. ``None`` may happen at the boundary of the data files.
This class is for E128 mode data.
'''
ElsCountData.__init__(self, t, count, level)
def __str__(self):
if self.hasLevel():
islevel = 'Defined'
else:
islevel = 'Not defined'
s = '<%s at %s, count_max=%.1f, level=(%s)>' % (self.__class__.__name__, self.starttime.strftime("%FT%T"),
self.count.max(), islevel)
return s
[docs] def guess_mode(self):
""" Return the mode of the object. It is always :attr:`irfpy.vels.mode.MODE_E128`.
:return: :irfpy.vels.mode.MODE_E128` always
"""
return mode.MODE_E128
[docs]class ElsCountDataE32(ElsCountData):
''' The class is ELS count data for 32 energy step mode.
'''
def __init__(self, t, count, level):
ElsCountData.__init__(self, t, count, level)
def __str__(self):
if self.hasLevel():
isLevel = 'Defined'
else:
isLevel = 'Not Defined'
s = '<%s at %s, count_max=%.1f, level=(%s)>' % (self.__class__.__name__, self.starttime.strftime("%FT%T"),
self.count.max(), isLevel)
return s
[docs] def guess_mode(self):
return mode.MODE_E32
[docs]class ElsCountFile:
''' Class for a single ELS count rate file.
:param filename: File name of ELS count rate file.
ELS count rate file is in matlab format.
Normally this is to reading the file, the instance should be
regarded as *immutable*.
>>> filename = ElsCountFile.get_sample_filename()
>>> elsfile = ElsCountFile(filename)
'''
def __init__(self, filename):
self.logger = logging.getLogger(self.__class__.__name__)
self.filename = filename
self.obstime_mpl = None
self.lvltime_mpl = None
try:
self.mat = scipy.io.loadmat(self.filename)
except IOError as e:
emsg = 'Failed to open the matlab file: %s' % self.filename
self.logger.error(emsg)
raise DatafileNotFoundError(emsg)
def __str__(self):
s = '<irfpy.vels.ElsCountFile: {}>'.format(self.filename)
return s
[docs] def get_elscountdata(self, t, level_file=None):
''' Pick the data of specified time
:param t: Time.
:type t: ``datetime.datetime``
:keyword level_file: It specifies an instance of :class:`ElsCountFile`
from which the level data will be read.
As default (*None*), the level data is read from *self*, however,
there may have a case level data is stored different file
(normally the file before).
:returns: The corresponding :class:`ElsCountData`.
:rtype: :class:`ElsCountData`
>>> countfile = ElsCountFile(ElsCountFile.get_sample_filename())
>>> els = countfile.get_elscountdata(datetime.datetime(2006, 12, 6, 8, 0, 0))
>>> print(els.starttime)
2006-12-06 07:59:57.845000
'''
elstimes = self.get_obstime()
idx = bisect.bisect(elstimes, t)
if idx <= 0:
raise ValueError('No data for %s in this datafile (too early)' % t.strftime('%FT%T'))
# if idx >= len(elstimes):
# if t - elstimes[-1] > datetime.timedelta(seconds=4):
# raise ValueError('No data for %s in this datafile (too late)' % t.strftime('%FT%T'))
idx = idx - 1
t0 = elstimes[idx]
if self.mat['elsmatrix'].ndim == 2:
count = self.mat['elsmatrix']
else:
count = self.mat['elsmatrix'][:, :, idx]
if level_file:
lvlfile = level_file
else:
lvlfile = self
elsleveltimes = lvlfile.get_leveltime()
elslevels = lvlfile.mat['elslevels']
if len(elsleveltimes) != elslevels.shape[1]:
self.logger.warning('Data file {} could be damaged. Continue without reading ELS level data.'.format(self))
level = None
else:
idx_lvl = bisect.bisect(elsleveltimes, t)
if idx_lvl <= 0:
level = None # In this case, level can not be defined.
else: # In this case, we know already data should exist as count data exist.
level = elslevels[:, idx_lvl-1]
return ElsCountData(t0, count, level)
[docs] def get_leveltime(self):
''' Return the level time list.
:returns: Level time list of ``datetime.datetime`` instance
:rtype: ``list`` of ``datetime.datetime``.
>>> elsfile = ElsCountFile(ElsCountFile.get_sample_filename())
>>> timelist = elsfile.get_leveltime()
>>> print(len(timelist))
112
>>> print(timelist[0].strftime('%FT%T'))
2006-12-06T07:00:11
>>> print(timelist[-1].strftime('%FT%T'))
2006-12-06T07:59:41
'''
if self.lvltime_mpl is not None:
return self.lvltime_mpl
lvltime_mat = self.mat.get('elsleveltimes').flatten()
self.lvltime_mpl = [utc.mat2dt(i) for i in lvltime_mat]
return self.lvltime_mpl
[docs] def get_obstime(self):
'''Return the observation time list.
:returns: Observation time list of ``datetime.datetime`` instance
:rtype: ``list`` of ``datetime.datetime``
>>> elsfile = ElsCountFile(ElsCountFile.get_sample_filename())
>>> timelist = elsfile.get_obstime()
>>> print(len(timelist))
783
>>> print(timelist[0].strftime('%FT%T'))
2006-12-06T07:00:03
>>> print(timelist[-1].strftime('%FT%T'))
2006-12-06T07:59:57
'''
if self.obstime_mpl is not None:
return self.obstime_mpl
obstime_mat = self.mat.get('elstimes').flatten()
self.obstime_mpl = [utc.mat2dt(i) for i in obstime_mat]
return self.obstime_mpl
[docs] def get_timerange(self, start_offset=0, stop_offset=0):
''' Return the time range of the observation.
:returns: Range of data file time as a list of
``datetime.datetime`` instances.
[*t0* + ``start_offset``, *t1* + ``stop_offset``] is returned,
where *t0* is the start time (the time of the first packet)
and *t1* is the end time (the time of the last packet).
This means that the end time *t1* does not agree with
the end of the observation time.
The end of the observation time is the end time plus
data acquisition time (4 sec, typically).
Specifying ``stop_offset`` will compensate this issue.
:keyword start_offset: Offset in seconds to be added to the start time.
:keyword stop_offset: Offset in seconds to be added to the stop time.
>>> elscnt = ElsCountFile(ElsCountFile.get_sample_filename())
>>> print(elscnt.get_timerange())
[datetime.datetime(2006, 12, 6, 7, 0, 3, 612000), datetime.datetime(2006, 12, 6, 7, 59, 57, 845000)]
>>> print(elscnt.get_timerange(stop_offset=4))
[datetime.datetime(2006, 12, 6, 7, 0, 3, 612000), datetime.datetime(2006, 12, 6, 8, 0, 1, 845000)]
'''
tlist = self.get_obstime()
off0 = datetime.timedelta(seconds=start_offset)
off1 = datetime.timedelta(seconds=stop_offset)
return [min(tlist) + off0, max(tlist) + off1]
[docs] @classmethod
def get_sample_filename(cls):
''' Returns sample file name.
>>> fnam = ElsCountFile.get_sample_filename()
>>> fnam.endswith('20061206070000.mat')
True
'''
from pkg_resources import resource_filename
# fn = resource_filename(__name__, os.path.join(
# 'sample', 'elsdata', '200612', '20061206070000.mat'))
rc = irfpy.util.irfpyrc.Rc()
p = rc.get('vels', 'elsmatbase')
fn = os.path.join(p, '200612', '20061206070000.mat')
logging.debug(fn)
return fn
[docs]class ElsCountFileDB:
''' Class for ELS count data base.
:param rootfolder: Root folder.
Data base is assumed to have a strucuture like
::
<rootfolder> -+- 200512 -+- 20051209130000.mat
| +- 20051209140000.mat
| +- ...
|
+- 200601 -+- 20060108100000.mat
| +- 20060108110000.mat
| +- ...
|
+- ...
The data base is regularly produced in spaghetti, and in this case
the ``<rootfolder>`` is ``/venus/elsdata/``.
For quick search, each file is assumed to contain the data
from the exact seconds to exactly one hour.
This means that ``20060108110000.mat`` is assumed contain data
between 2006-01-08T11:00:00 and 2006-01-08T12:00:00.
'''
def __init__(self, rootfolder, enable_cache=True):
'''Initiate ElsCountFileDB.
:param rootfolder: The root folder.
:param enable_cache: Enable cache. ~50% performance up according to ``scidata_test2.py``.
'''
self.logger = logging.getLogger(self.__class__.__name__)
# self.logger.setLevel(logging.WARN)
self.db = irfpy.util.timeseriesdb.DB()
self.enable_cache = enable_cache
if enable_cache:
self.elsfilecache = irfpy.util.ringcache.RingCache(10)
if os.path.isdir(rootfolder):
list = os.listdir(rootfolder) # <rootfolder>/######
for l in list:
if len(l) != 6: # Only 6 chars
continue
try: # Only integer
int(l)
except ValueError:
continue
full_l = os.path.join(rootfolder, l)
if not os.path.isdir(full_l): # Only directory
continue
self.logger.debug("# Examine filder: ``%s``" % full_l)
filelist = os.listdir(full_l) # <full_l>/**********0000.mat
for f in filelist:
if len(f) != 18:
continue
if not f.startswith(full_l[-6:]):
continue
if not f.endswith('.mat'):
continue
try:
int(f[:14])
except ValueError:
continue
full_f = os.path.join(full_l, f)
yr, mo, dy, hr = int(f[0:4]), int(f[4:6]), int(f[6:8]), int(f[8:10])
starttime = datetime.datetime(yr, mo, dy, hr, 0, 0)
self.logger.debug('### %s - %s' % (full_f, starttime))
self.db.append(full_f, starttime)
self.logger.debug("# Match: ``%s`` @ %s" % (full_f, starttime.strftime('%FT%T')))
else:
raise IOError(
'``%s`` does not exist or is not directory.' % rootfolder)
self.logger.info("# Data base loaded. Total %d files." % len(self.db))
[docs] def get_elscountfile(self, filename):
if self.enable_cache:
if self.elsfilecache.hasKey(filename):
els = self.elsfilecache.get(filename)
else:
els = ElsCountFile(filename)
self.elsfilecache.add(filename, els)
return els
else:
return ElsCountFile(filename)
[docs] def get_filename_fromdb(self, t):
''' Return the file name registered in the data base.
No clever check will be done as :meth:`get_filename`.
Just return the file name registered in the data base.
'''
return self.db.get(t)
[docs] def get_filename(self, t, offset=4):
''' Return the filename of the Matlab data file.
:param t: Time
:type t: ``datetime.datetime``
:keyword offset: It is deprecated, and it does not offer anything.
From version 4.2.3a6.
If the given time is not in the data base, i.e. no observation
considering offset, DataNotInDbError is raised.
'''
# self.logger.setLevel(logging.DEBUG)
fnam = self.db.get(t)
self.logger.debug('Sepcified time %s' % t)
self.logger.debug('Filename %s' % fnam)
self.logger.debug(' - Time (registered): %s' % self.db.gettime(fnam))
elsfile = self.get_elscountfile(fnam)
trange = elsfile.get_timerange()
self.logger.debug(' - Time (file) %s-%s' % (trange[0], trange[1]))
if t < trange[0]:
self.logger.debug('*** Range error. Too early. Check the previous data file.')
fnam = self.db.previousof(fnam)
self.logger.debug('Filename %s' % fnam)
self.logger.debug(' - Time (registered): %s' % self.db.gettime(fnam))
elsfile = self.get_elscountfile(fnam)
trange = elsfile.get_timerange()
self.logger.debug(' - Time (file) %s-%s' % (trange[0], trange[1]))
# t1 = trange[1] + datetime.timedelta(seconds=offset)
# if t > t1:
# raise irfpy.util.timeseriesdb.DataNotInDbError(
# 'Range error. No data corresponding to %s.\nThe last data is for %s' % (t, t1))
return fnam
elif t > trange[1]:
# In this case, there are three possibilities. Return the original file,
# return the next file, or raise exception.
self.logger.debug('*** Range error. Too late. Check the previous data file.')
try:
fnam2 = self.db.nextof(fnam)
except irfpy.util.timeseriesdb.DataNotInDbError:
return fnam
self.logger.debug('Filename %s' % fnam2)
self.logger.debug(' - Time (registered): %s' % self.db.gettime(fnam2))
elsfile = self.get_elscountfile(fnam2)
trange2 = elsfile.get_timerange()
self.logger.debug(' - Time (file): %s-%s' % (trange2[0], trange2[1]))
if t >= trange2[0]: # The specified data is in the range of the second file.
return fnam2
# else: # The specified data is not in the range of the second file.
# # Here we check the first file time + offset included the specified data.
# t1 = trange[1] + datetime.timedelta(seconds=offset) # Refer to the previous data.
# if t <= t1:
# return fnam
# else:
# raise irfpy.util.timeseriesdb.DataNotInDbError(
# 'Range error. No data corresponding to %s.\nThe last data is %s' % (t, t1))
else:
return fnam
else:
return fnam
# self.logger.setLevel(logging.WARN) # Set back to default.
[docs] def previousof(self, filename):
return self.db.previousof(filename)
[docs] def nextof(self, filename):
return self.db.nextof(filename)
[docs] @classmethod
def get_defaultdb(cls, enable_cache=True):
''' Return the instance of ElsCountFileDB2 from default rootfolder.
The default rootfolder can be defined by ``[vels]-elsmatbase``
in ``.irfpyrc``.
'''
logger = logging.getLogger(cls.__name__)
# logger.setLevel(logging.INFO)
rc = irfpy.util.irfpyrc.Rc()
rootfolder = rc.get('vels', 'elsmatbase')
logger.info('Default DB loading from ``%s``.' % rootfolder)
if rootfolder is None:
emesg = "No entry of [vels]-elsmatbase in .irfpyrc."
raise RuntimeError(emesg)
# logger.setLevel(logging.WARN)
return ElsCountFileDB(rootfolder, enable_cache=enable_cache)
[docs]def get_defaultdb(enable_cache=True):
''' This is an easy-accessor for :meth:`ElsCountFileDB.get_defaultdb`.
'''
return ElsCountFileDB.get_defaultdb(enable_cache=enable_cache)
### Below is for easy access of the dataset to return the data.
class __DataFactory:
''' Factory of the ELS data access. Used for package methods having same names.
The factory would consist only of class methods.
The class should be considered as a singleton by the developer of the moduel :mod:`scidata`.
The only class member is the default instance of :class:`ElsCountFileDB`.
The instance is generated by :meth:`get_defaultdb`.
'''
__database = None # ElsCountFileDB, this should not be referred explictily.
# ``db = cls._mkdb()`` will return this database (and create it if necessary).
@classmethod
def _isdb(cls):
try:
cls._mkdb()
return True
except Exception as e:
logging.warn(e)
return False
@classmethod
def _mkdb(cls):
if cls.__database is None:
cls.__database = get_defaultdb()
return cls.__database
@classmethod
def get_countdata(cls, t):
''' Implementation of :meth:`get_countdata`.
'''
logger = logging.getLogger('getcount')
# logger.setLevel(logging.INFO)
db = cls._mkdb() # Should always called.
fn = db.get_filename(t)
elscountfile = db.get_elscountfile(fn)
els = elscountfile.get_elscountdata(t)
if not els.hasLevel():
logger.info('No level data found. Try previous one.')
try:
fn2 = db.previousof(fn)
elslevelfile = db.get_elscountfile(fn2)
els = elscountfile.get_elscountdata(t, level_file=elslevelfile)
except irfpy.util.timeseriesdb.DataNotInDbError:
logger.info('No file in DB.')
t0 = els.starttime
# logger.setLevel(logging.WARN)
return (t0, els)
@classmethod
def get_count(cls, t):
''' Implementation of :meth:`get_count`.
'''
logger = logging.getLogger('get_count')
# logger.setLevel(logging.INFO)
db = cls._mkdb() # Should always called.
t0, els = cls.get_countdata(t)
elscnt = els.get_count()
# logger.setLevel(logging.WARN)
return (t0, elscnt)
@classmethod
def get_cps(cls, t):
''' Implementation of :meth:`get_cps`.
'''
db = cls._mkdb()
t0, els = cls.get_countdata(t)
elscnt = els.get_count() # Return the raw count. 128x16 array
elscps = bible_flux.raw2cps_np(elscnt)
return (t0, elscps)
@classmethod
def __energy_table_handling(cls, energy_table, elscountdata):
''' get_??? functions have energy_table keyword. Handles it.
This function returns the energy table as numpy.array according to the
keyword energy_table.
Details of the energy_table keyword can be found :meth:`get_dnf` function.
'''
logger = logging.getLogger()
cnts = elscountdata.get_count()
ne, nd = cnts.shape
try: # Check if energy_table is numpy.array or not.
ets = energy_table.shape
if len(ets) == 1:
if ets[0] != ne:
raise ValueError('Given energy table (%d) do not match the data (%d).' % (ets[0], ne))
etbl = energy_table
elif len(ets) == 2:
if ets != (ne, nd):
raise ValueError('Given energy table (%d, %d) do not match the data (%d, %d).' % (ets[0], ets[1], ne, nd))
etbl = energy_table
else:
raise ValueError('Given energy table (dim=%d) should have dimension 1 or 2.' % (len(ets)))
except AttributeError as e:
if energy_table in (None, "default", "defaultE128"):
logger.debug('Energy table default E128')
etbl = energy.get_default_table_128()
logger.debug('etbl: %.1f %.1f %.1f ... %.1f' % (etbl[0], etbl[1], etbl[2], etbl[-1]))
elif energy_table in ("defaultE32"):
logger.debug('Energy table default E32')
etbl = energy.get_default_table_32_as_128()
logger.debug('etbl: %.1f %.1f %.1f ... %.1f' % (etbl[0], etbl[1], etbl[2], etbl[-1]))
elif energy_table in ("data"):
logger.debug('Energy table guessing from data.')
etbl = elscountdata.get_energy()
else:
raise ValueError('Unknown energy_table %s.' % energy_table)
return etbl
@classmethod
def get_dnf(cls, t, energy_table=None):
''' Implementation of :meth:`get_dnf`.
'''
logger = logging.getLogger('get_dnf')
# logger.setLevel(logging.DEBUG)
warnings.warn("The irfpy.vels.scidata.get_dnf is yet experimental. Use with the users' responsibility."
"It is highly recommended to contact to Futaana for details")
t0, elscountdata = cls.get_countdata(t)
cnts = elscountdata.get_count()
ne, nd = cnts.shape
logger.debug('Ne=%d, Nd=%d' % (ne, nd))
etbl = cls.__energy_table_handling(energy_table, elscountdata)
dlist = np.arange(nd)
dnf = bible_flux.raw2dnf_np(cnts, etbl, dlist)
return t0, etbl, dnf
@classmethod
def get_def(cls, t, energy_table=None):
''' Implementation of :meth:`get_def`.
'''
warnings.warn('Use get_def is deprecated. Use equivalent method of ``get_def2``.', DeprecationWarning)
logger = logging.getLogger('get_def')
# logger.setLevel(logging.DEBUG)
t0, elscountdata = cls.get_countdata(t)
cnts = elscountdata.get_count()
ne, nd = cnts.shape
logger.debug('Ne=%d, Nd=%d' % (ne, nd))
etbl = cls.__energy_table_handling(energy_table, elscountdata)
dlist = np.arange(nd)
deflux = bible_flux.raw2def_np(cnts, etbl, dlist)
return t0, etbl, deflux
@classmethod
def get_def2(cls, t, energy_table=None, constant_background=0.0, percentile_background=None):
''' Implementation of :meth:`get_def2`.
:param t: Time.
:param energy_table: Energy table to specify.
:param constant_background: Subtract background before calculating the flux.
:param percentile_background: Specify the background level from the percentile of the data counts.
Useful for high background count events. It overwrite the value in ``constant_background``.
For example, specifying 10 will produce the constant background by calculating the
10% percentile as a constant background. Specifying 50 will lead the median of the counts
to be a constant background.
:type percentile_background: ``int``
:returns: A list of the following. ``time``, energy table, differential energy flux, and background.
Background is also a list including the background count and the number of removal in counts.
'''
logger = logging.getLogger('get_def')
# logger.setLevel(logging.INFO)
warnings.warn("The irfpy.vels.scidata.get_def2 is yet experimental. Use with the users' responsibility."
"It is highly recommended to contact to Futaana for details")
t0, elscountdata = cls.get_countdata(t)
cnts = elscountdata.get_count()
ne, nd = cnts.shape
logger.debug('Ne=%d, Nd=%d' % (ne, nd))
etbl = cls.__energy_table_handling(energy_table, elscountdata)
dlist = np.arange(nd)
totcnts = cnts.sum()
if percentile_background is not None:
constant_background = np.percentile(cnts, percentile_background)
logger.debug('Percentile background (pc=%f): BG=%f' % (percentile_background, constant_background))
cnts -= constant_background
cnts = np.where(cnts > 0, cnts, np.zeros_like(cnts))
totvalidcnts = cnts.sum()
deflux = bible_flux.raw2def_np(cnts, etbl, dlist)
return t0, etbl, deflux, (constant_background, totcnts - totvalidcnts)
@classmethod
def get_psd(cls, t, energy_table=None):
''' Implementation of :meth:`get_psd`.
'''
logger = logging.getLogger('get_psd')
# logger.setLevel(logging.DEBUG)
warnings.warn("The irfpy.vels.scidata.get_psd is yet experimental. Use with the users' responsibility."
"It is highly recommended to contact to Futaana for details")
t0, elscountdata = cls.get_countdata(t)
cnts = elscountdata.get_count()
ne, nd = cnts.shape
logger.debug('Ne=%d, Nd=%d' % (ne, nd))
etbl = cls.__energy_table_handling(energy_table, elscountdata)
dlist = np.arange(nd)
psd = bible_flux.raw2psd_np(cnts, etbl, dlist)
return t0, etbl, psd
@classmethod
def get_obstime(cls, t0, t1):
''' Implementation of :meth:`get_obstime`.
'''
logger = logging.getLogger('get_obstime')
db = cls._mkdb()
obstime = []
t_file = db.get_filename_fromdb(t0)
try:
t_file = db.previousof(t_file)
except irfpy.util.timeseriesdb.DataNotInDbError:
t_file = t_file
els_file = db.get_elscountfile(t_file)
t_range = els_file.get_timerange()
while t_range[0] <= t1:
obstime = obstime + els_file.get_obstime()
try:
t_file = db.nextof(t_file)
except irfpy.util.timeseriesdb.DataNotInDbError:
logger.info('No further file.')
break
els_file = db.get_elscountfile(t_file)
t_range = els_file.get_timerange()
otime = []
for t in obstime:
if t >= t0 and t<= t1:
otime.append(t)
otime.sort()
return otime
@classmethod
def get_counts(cls, t0, t1):
""" Obtain an array of ELS counts between the given time period.
:returns: (tlist, cnts) where tlist is the list of obstime and cnts as arrawy (128, 16, N) shape
"""
tlist = cls.get_obstime(t0, t1)
tobslist = np.array(tlist, copy=True)
ntlist = len(tlist)
arr = np.zeros([128, 16, ntlist])
for i in range(ntlist):
tobs, val = cls.get_count(tlist[i])
tobslist[i] = tobs
arr[:, :, i] = val
return tobslist, arr
[docs]def get_countdata(t):
''' Return the count data as an array of (``datetime``, :class:`ElsCountData` instance).
Returns count rate data as an array of the real observation time and :class:`ElsCountData` instance.
The returned data time, *t0*, should satisfy *t0* <= *t* where *t* is the specified time is the specified time.
The *t0* is the real observation data time that contains the data obtained at *t*.
:param t: Time
:type t: ``datetime.datetime``
:returns: An array of observation time and corresponding :class:`ElsCountData` instance.
If no data found, DataNotInDbError is raised.
'''
return __DataFactory.get_countdata(t)
[docs]def get_count(t):
''' Return the ELS count data at the specified time, *t*.
:param t: Time
:type t: ``datetime.datetime``.
:returns: An array of the observation time and
corresponding numpy array with (128, 16) shape.
Return the ELS count rate date as numpy.array type with (128, 16) shape.
>>> t0 = datetime.datetime(2006, 12, 5, 8)
>>> t, els = get_count(t0)
>>> print(t) # This may fail if data base changed.
2006-12-05 07:59:59.821000
>>> print(els.shape) # This may fail if data base changed.
(128, 16)
>>> print(els.max()) # This may fail if data base changed.
25.0
'''
return __DataFactory.get_count(t)
[docs]def get_cps(t):
''' Return the ELS count rate data at the time *t*.
:param t: Time
:type t: ``datetime.datetime``.
:returns: An array of the observation time and
corresponding numpy array with (128, 16) shape.
'''
return __DataFactory.get_cps(t)
[docs]def get_dnf(t, energy_table=None):
''' Return the ELS differential number flux.
:param t: Time
:type t: ``datetime.datetime``
:keyword energy_table: Energy table used for the calculation.
Allowed are ``None``, ``"default"``, ``"defaultE128"``, ``"defaultE32"``, ``"data"``,
or ``numpy.array`` with (*nE*, ) or (*nE*, *nD*) shape.
:returns: List of the observation time *t0* (including data on *t*), energy table, and its DNF.
:rtypes: (t0, etbl, dnf) where t0 is ``datetime.datetime``, etbl is ``numpy.array`` with (nE, nD) shape
in the unit of eV, and dnf is ``numpy.array`` with (nE, nD) shape in the unit of /eV cm2 sr s.
Energy table is not trivial for ELS, so the function support several ways of
using energy table.
======================================== =============================================================================
*energy_table* keyword Table used
======================================== =============================================================================
``None`` "default" or "defaultE128" Default table by :meth:`irfpy.vels.energy.get_default_table_128`.
"defaultE32" Default table by :meth:`irfpy.vels.energy.get_default_table_32_as_128`.
"data" Guessed from the data obtained. See :meth:`ElsCountData.guess_mode`
``numpy.array`` with (*nE*,) shape Used the given table itself. No anode dependence (repeated *nD* times)
``numpy.array`` with (*nE*, *nD*) shape Used the given table itself.
======================================== =============================================================================
'''
return __DataFactory.get_dnf(t, energy_table=energy_table)
[docs]def get_def(t, **args):
''' Return the ELS differential energy flux in the unit of eV/eV cm2 sr s.
For details see :meth:`get_dnf`.
'''
warnings.warn('Use get_def is deprecated. Use equivalent method of ``get_def2``.', DeprecationWarning)
return __DataFactory.get_def(t, **args)
[docs]def get_def2(t, **args):
''' Return the ELS differential energy flux in the unit of eV/eV cm2 sr s.
See details also in :meth:`get_dnf`.
'''
return __DataFactory.get_def2(t, **args)
[docs]def get_psd(t, energy_table=None):
''' Return the ELS velocity distribution functions, in the unit of s3/cm6.
For details see :meth:`get_dnf`.
'''
return __DataFactory.get_psd(t, energy_table=energy_table)
[docs]def get_obstime(t0, t1):
''' Get the obseration times (4-s resolution).
:keyword t0: Start time
:type t0: ``datetime.datetime``
:keyword t1: Start time
:type t1: ``datetime.datetime``
:returns: List of ``datetime.datetime``.
The following tests shows simple samples of using
:meth:`get_obstime`.
This is based on database on 110704, so tests may be failed
according to the change in the database.
>>> t0 = datetime.datetime(2006, 12, 5)
>>> t1 = datetime.datetime(2006, 12, 6)
>>> t = get_obstime(t0, t1)
>>> print(len(t)) # Database change may cause fails of the test.
3619
>>> print(t[0]) # Database change may cause fails of the test.
2006-12-05 05:17:27.936000
>>> print(t[-1]) # Database change may cause fails of the test.
2006-12-05 09:54:25.980000
'''
return __DataFactory.get_obstime(t0, t1)
[docs]def get_counts(t0, t1):
""" Return the counts between the given time period.
:param t0: Start time.
:param t1: End time
:returns: An array with a shape of (128, 16, N), where N is the number of data (time).
>>> t0 = datetime.datetime(2006, 12, 5, 6, 10)
>>> t1 = datetime.datetime(2006, 12, 5, 6, 15)
>>> tlist, elsdat = get_counts(t0, t1)
>>> print(len(tlist))
65
>>> print(elsdat.shape)
(128, 16, 65)
>>> print(elsdat.max())
44.0
"""
return __DataFactory.get_counts(t0, t1)
[docs]def isdb():
''' Return if the default database is available or not.
'''
return __DataFactory._isdb()