Source code for irfpy.mels.scidata

""" MEX/ELS data accessor

.. note::

    As the VEX/ELS modules are highly dedicated,
    I gave up to spend time on refactoring and
    making coordinated class between MEX & VEX.

    Though, I am trying to make the interface as similar as possible.

*Preparation*

Download the tree structure of elsdata from spaghetti.irf.se.
Then, specify the downloaded path in ~/.irfpyrc or ./.irfpyrc
by the entry of elsmatbase under the section ``[mels]``.

*Access*

The following two functions are the easy accessor for users.

- :func:`getobstime` will give you the actual start time list
  of each ELS data (energy/direction spectra, usually 4s data)
  for the given time range.

- :func:`gettimeseries` will return the time series of the ELS count rate
  dataset in :class:`ElsCountTimeSeries` for the given time range.

*Other notes*

In this module, the simplest way of ELS data analysis is employed.
This means that we do not refer to the HV monitor level, but only the count data.
The energy and g-factor will be calculated from the static table.

As of the end of 2013, the energy table and g-factor is constant (presumably).
I referred to PDS, and this is the case.

The hyerarchy of dataset in this module is

* Single time data of counts. Full mode is 128 energy step & 16 directions every 4 s.
  This is implemented as :class:`ElsCountData` class.

* A series of single time data, implemented as :class:`ElsCountTimeSeries`.

* A file, in Matlab form, implemented as :class:`ElsCountFile`.  Both data and HV information is there.

* A collection of file (data file structure), implemented as :class:`ElsCountFileDB`.
"""
import os
import datetime
import scipy.io
import numpy as np

import irfpy.util.irfpyrc
import irfpy.util.tdict
import irfpy.util.ringcache
import irfpy.util.timeseriesdb
from irfpy.util import utc
from irfpy.util import exception as ex

[docs]class ElsCountData: """ Container of ELS count data. Base class. """ def __init__(self, t, counts): self.counts = np.ma.array(counts, dtype=np.int_) # As the dtype could be np.uint8, it is good to clearly specify the type. self.t = t def __str__(self): s = "<{clsnmn:s}: {t:%FT%T.%f} Max={max:d} Min={min:d}>".format( clsnmn=self.__class__.__name__, t=self.t, max=self.counts.max(), min=self.counts.min()) return s def __getitem__(self, item): return self.counts.__getitem__(item)
[docs]class ElsCountDataBurst(ElsCountData): """ Container of ELS count data. (128, 16) array. ELS count data container. The object contains only two attributes. - :attr t:`The time` - :attr counts: Count matrix (128, 16) shaped np.ma.array. Usually the users do not need to create this object. >>> # Prepare the matrix and time for dummy data >>> t = datetime.datetime(2014, 10, 18, 13, 15, 27, 123456) >>> matx = np.arange(128 * 16).reshape([128, 16]) >>> els_data = ElsCountDataBurst(t, matx) >>> print(els_data) <ElsCountDataBurst: 2014-10-18T13:15:27.123456 Max=2047 Min=0> >>> print(els_data.counts.shape) (128, 16) >>> print(els_data[3, 12]) 60 >>> print(els_data[3, 12:16]) [60 61 62 63] """ def __init__(self, t, counts): if counts.shape != (128, 16): raise ex.IrfpyException("Format of count array is different (given: {}}".format( counts.shape )) ElsCountData.__init__(self, t, counts)
[docs]class ElsCountTimeSeries(irfpy.util.tdict.TimeDict): """ Time series of :class:`ElsCountData` (and its derivable) object. The object stores the ELS count data. Usually produced by :func:`gettimeseries`, and to make the array form, use :meth:`toarray128`. >>> t0 = datetime.datetime(2014, 10, 19, 18, 30) >>> t1 = datetime.datetime(2014, 10, 19, 18, 32) # 2 minutes >>> timeser = gettimeseries(t0, t1) # You get ElsCountTimeSeries object >>> print(len(timeser)) 27 >>> tlist, arr = timeser.toarray128() # You get time information together with matrix form of counts. >>> print(len(tlist)) 27 >>> print(arr.shape) (128, 16, 27) """ def __init__(self): irfpy.util.tdict.TimeDict.__init__(self)
[docs] def toarray128(self): """ Return the np.array expression of count data """ n = len(self) val = np.zeros([128, 16, n]) # Shape (E128, D16, Tn) tlist = [] for i, (t, v) in enumerate(self): # t is datetime, v is ElsCountData128 tlist.append(t) val[:, :, i] = v.counts return tlist, val
[docs]class ElsCountFile: """ A container of ELS matlab file. Represents a single file contents. Usually, the object is made via calling, :func:`els_count_file_factory`. >>> els_count_file = els_count_file_factory("/Volumes/scidata/data/mars/elsdata/201410/20141019180000.mat") # doctest: +SKIP """ def __init__(self, filename): self.filename = filename self.mat = None def _load(self): if self.mat is None: self.mat = scipy.io.loadmat(self.filename) tlist = self.mat['elstimes'].flatten() self.tlist = [utc.mat2dt(i) for i in tlist] self.cnts = self.mat['elsmatrix'] # (128, 16, Nt) shape
[docs] def getobstime(self): self._load() return self.tlist
[docs] def get_timeseries(self): """ Return the time series of the count data. :returns: Time series of the count data :rtype: :class:`ElsCountTimeSeries` """ self._load() timeser = ElsCountTimeSeries() tlist = self.getobstime() ndat = len(tlist) for idat, t in enumerate(tlist): els_count_data = ElsCountData(t, self.cnts[:, :, idat]) timeser.append(t, els_count_data, dup='ignore') return timeser
class _ElsCountFileFactory: """ :class:`ElsCountFile` is cached. """ def __init__(self): self.cache = irfpy.util.ringcache.RingCache(500) def get(self, filename): if not self.cache.hasKey(filename): els_count_file = ElsCountFile(filename) self.cache.add(filename, els_count_file) return self.cache.get(filename) def __call__(self, filename): return self.get(filename) els_count_file_factory = _ElsCountFileFactory() def _fn2t0(filename): """ Get the start time from the file name. """ els_count_file = els_count_file_factory.get(filename) return els_count_file.getobstime()[0]
[docs]class ElsCountFileDB: """ Represent the database of ELS counts file in Matlab form. One have to download the Matlab file from spaghetti.irf.se to start with. Then, one has to setup ~/.irfpyrc or ./.irfpyrc as follows. :: [mels] elsmatbase = /path/to/els/matlab/file For usual use, you may want to get a file name where the specific observation is obtained. >>> elsdb = ElsCountFileDB() >>> t = datetime.datetime(2014, 10, 19, 18, 30) # Specify the time >>> datafile = elsdb.get(t) # Obtain the filename that included the data at t. >>> print(os.path.basename(datafile)) 20141019180000.mat >>> t1 = datetime.datetime(2014, 10, 19, 20, 25) >>> files = elsdb.getfiles(t, t1) # Obtain the list of filenames >>> print(len(files)) 3 >>> print(os.path.basename(files[1])) 20141019190000.mat """ def __init__(self): rc = irfpy.util.irfpyrc.Rc() rootdir = rc.get('mels', 'elsmatbase') if rootdir is None: raise ex.IrfpyError('Error in .irfpyrc. [mels] elsmatbase should be defined') if not os.path.isdir(rootdir): raise IOError('Data directory {} not found. Check .irfpyrc [mels] elsmatbase'.format(rootdir)) guess_db = irfpy.util.timeseriesdb.DB() ### First obtain a database from time stamp in the file anem for pp, dd, ff in os.walk(rootdir): ppbase = os.path.basename(pp) # print(ppbase) for ffn in ff: if len(ffn) != 18: continue if not ffn.endswith('.mat'): continue # print(ffn) try: yr = int(ffn[0:4]) mo = int(ffn[4:6]) dy = int(ffn[6:8]) hr = int(ffn[8:10]) # print(ffn, yr, mo, dy, hr) guess_db.append(os.path.join(pp, ffn), datetime.datetime(yr, mo, dy, hr, 0, 0)) except: continue self.db = irfpy.util.timeseriesdb.FazzyDB(guess_db, _fn2t0)
[docs] def get(self, t): """ Return the filename. :param t: Time of observation :return: The file name where the data is included. """ return self.db.get(t)
[docs] def getfiles(self, t0, t1): """ Return the filename list for the given time range. :param t0: Start of inquiry :param t1: End of inquiry :return: The file names. """ return self.db.getfiles(t0, t1)
_elsdb = None
[docs]def gettimeseries(t0, t1): """ Return the time series for each energy spectra :param t0: Start of inquiry :param t1: Stop of inquiry :return: Time series of ELS count data :rtype: :class:`ElsCountTimeSeries` >>> t0 = datetime.datetime(2014, 10, 19, 17, 55) >>> t1 = datetime.datetime(2014, 10, 19, 19, 8) >>> tseries = gettimeseries(t0, t1) >>> print(len(tseries)) 954 >>> t0r, dat0 = tseries[0] # The first data >>> print(t0r, dat0) 2014-10-19 17:54:55.271000 <ElsCountData: 2014-10-19T17:54:55.271000 Max=20 Min=0> >>> t1r, dat1 = tseries[-1] # The last data >>> print(t1r, dat1) 2014-10-19 19:07:56.202000 <ElsCountData: 2014-10-19T19:07:56.202000 Max=76 Min=0> """ global _elsdb if _elsdb is None: _elsdb = ElsCountFileDB() files = _elsdb.getfiles(t0, t1) # Obtain the file names elstimeseries = ElsCountTimeSeries() # ELS time series, zero size. for file in files: # Walk through the file # Instance the ElsCountFile object els_count_file = els_count_file_factory(file) elstimeseries = elstimeseries + els_count_file.get_timeseries() return elstimeseries.clipped(t0, t1)
[docs]def getobstime(t0, t1): """ Return the list of observation time. :param t0: Start of inquiry :param t1: Stop of inquiry :return: List of observation time >>> t0 = datetime.datetime(2014, 10, 19, 17, 55) >>> t1 = datetime.datetime(2014, 10, 19, 19, 8) >>> obstime = getobstime(t0, t1) >>> print(len(obstime)) 954 >>> print(obstime[0]) 2014-10-19 17:54:55.271000 >>> print(obstime[-1]) 2014-10-19 19:07:56.202000 """ timeser = gettimeseries(t0, t1) return timeser.getobstime()
[docs]def get_counts(t0, t1): """ Return the numpy array of the ELS counts. :param t0: Start of inquiry :param t1: Stop of inquiry :return: A tuple with length of two. ``(obstime, counts)`` where ``obstime`` is the list of the observation time and ``counts`` is the array of the counts. >>> t0 = datetime.datetime(2014, 10, 19, 17, 55) >>> t1 = datetime.datetime(2014, 10, 19, 19, 8) >>> obstime, counts = get_counts(t0, t1) >>> print(len(obstime)) 954 >>> print(obstime[0]) 2014-10-19 17:54:55.271000 >>> print(obstime[-1]) 2014-10-19 19:07:56.202000 >>> print(counts.shape) (128, 16, 954) >>> print(counts.max()) 432.0 """ timeser = gettimeseries(t0, t1) return timeser.toarray128()
[docs]def get_cps(t0, t1): raise NotImplementedError('')
[docs]def get_dnf(t0, t1): raise NotImplementedError('')
[docs]def get_def(t0, t1): raise NotImplementedError('')