How to handle IMA science data

Due to highly complex data structure of the IMA data (too many modes, time series of 4-dimensional, many processing levels, different operations, similar sensors on MEX & VEX but not fully compatible operations and processing, …) the access to the IMA science data is highly complicated.

The irfpy.aspera library give a quick access to such dataset.

It is noted that the “final form of the matrix” depends highly on the users. Thus, the interface of the data access should be as simple and flexible as possible.

Getting the matrix at the specific time

To get the full matrix (192s data), you may just use the following syntax.

>>> import datetime
>>> from irfpy.vima import scidata_util      # Use from irfpy.mima import scidata_util for MEX.
>>> tstart, data = scidata_util.get3d(datetime.datetime(2008, 1, 9, 14, 30))

The returned tstart is the time of the start of the data. The data includes the VEX/IMA data.

The data.matrix is the matrix that you want.

>>> print(data.matrix)

This is (32, 16, 96, 8) shaped array.

Getting the time series matrix

To get the time series of the matrix, you may use the following syntax.

>>> import datetime
>>> from irfpy.vima import scidata_util
>>> timeser = scidata_util.getarray2d(datetime.datetime(2008, 1, 9, 14, 30), (datetime.datetime(2008, 1, 9, 15, 30)))
>>> print(timeser)
<TimeSeriesCntMatrix2D:len=105 from 2008-01-09 14:29:40.863880 to 2008-01-09 15:11:16.908056>

A bit more in details

ASPERA/IMA data structure

The ASPERA/IMA sensor can measure ions with mass, energy and direction discrimilated.

At a certain time, IMA measures ions at a specific energy and polar angle. Simultaneously, the mass and azimuth angle can be discrimimated. It means that the data in IMA at a certain time has (M=32, A=16) shaped array. The accumulation time is 125 ms.

Then, the energy is swept. Usually, energy has 96 steps with 12 sec. Considering this, the 12-sec resolution data with (M=32, A=16, E=96) can be obtained.

The energy sweep is repeated with different polar angle scan. 16 steps (=192 s) are sweeped. The resulting dataset is the 192-sec resolution data with (M=32, A=16, E=96, P=16).

This 192-sec resolution data is transferred as is if the mode is 24. In this mode, the data is 192-sec resolution data with (M=32, A=16, E=96, P=16). For MEX, mode=24 is used near the pericenter.

If the mode is 25, the neighboring polar angle is added up on board, and sent to the ground. In this mode, the data is 192-sec resolution data with (M=32, A=16, E=96, P=8). For VEX, mode=25 has been used for the beginning of the mission. Later the mode=24 is used, but with 1 count reduction, for VEX.

In summary, the data at a certain time can be considered as (depending on the way of data analysis)

1. 3D (azimuth, energy, polar + mass) data with time resolution of 192s

2. 2D (azimuth, energy + mass) data with time resolution of 12s (mode 24)

3. 2D (azimuth, energy + mass) data with time resolution of 24s (mode 25)

4. 1D (azimuth) data with resolution of 125 ms at specific energy and (mode 24)

5. 1D (azimuth) data with resolution of 250 ms at specific energy and (mode 25)

Each of them is indeed a time series.

More advanced dataset

For IMA, there are other “operation”, such as 32 second mode, and fast-sweep mode (24-s resolution). 32 second mode has not used so that IRFPY will not support this. Fast sweep mode was prepared for Phobos flyby for MEX. This is intended to be implemented into IRPFY.

Object oriented dataset

For irfpy.aspera, each of these data is represented by a specific class in an object oriented way. Below we review the class used for each dataset.

Data	At specific time	Time series
3D (Az, En, Po; Ma)	irfpy.imacommon.imascipac.CntMatrix	irfpy.imacommon.imascipac.TimeSeriesCntMatrix
2D (Az, En; Ma)	irfpy.imacommon.imascipac.CntMatrix2D	irfpy.imacommon.imascipac.TimeSeriesCntMatrix2D
1D (Az; Ma)	N/A (To be developed)	N/A (To be developed)

3D (azimuth, energy, polar + mass)

To represent the 3D distribution data at a certain time, irfpy.imacommon.imascipac.CntMatrix class is used. The object contains the count rate matrix, mode, and observation time. Count rate matrix is a masked array with the shape corresponding to the mode (for Mode 24, the shape is (32, 16, 96, 16)).

2D (azimuth, energy + mass) data

To represent the 2D distribution at a ceratin time, irfpy.imacommon.imascipac.CntMatrix2D class is used. The object contains the count rate matrix, mode, polar index, and observation time. Count rate matrix is a masked array with the shape corresponding to the mode (but both mode 24 and 25, the shape is (32, 16, 96)). The polar index is a pair of (start, stop), both inclusive.

Time series of the data

Time series of 3D data is irfpy.imacommon.imascipac.TimeSeriesCntMatrix. Use irfpy.vima.scidata_util.getarray3d to get the time series data.

Time series of 2D data is irfpy.imacommon.imascipac.TimeSeriesCntMatrix2D. Use irfpy.vima.scidata_util.getarray2d to get the time series data.

To read (VEX/IMA)

The first step to read the VEX/IMA is getting time series of 2D data.

Reading time series of 2D data

Use irfpy.vima.scidata_util.getarray2d() function. With time intervals as argument, the function will return the time series of 2D data (irfpy.imacommon.imascipac.TimeSeriesCntMatrix2D.

>>> import datetime
>>> from irfpy.vima import scidata_util
>>> t0 = datetime.datetime(2012, 10, 5, 0, 1, 23)
>>> t1 = datetime.datetime(2012, 10, 6, 4, 59, 43)
>>> dataset = scidata_util.getarray2d(t0, t1)

The dataset represent the time series of VEX/IMA 2D data. Let’s look at the dataset.

>>> print(dataset)
<TimeSeriesCntMatrix2D:len=3362 from 2012-10-05 00:01:11.198276 to 2012-10-06 03:44:03.177666>

You may find that there are 3362 2-D data.

You can use the getobstime() method to get the list of observation time.

>>> obstime = dataset.getobstime()
>>> print(obstime)
[datetime.datetime(2012, 10, 5, 0, 1, 11, 198276),
 datetime.datetime(2012, 10, 5, 0, 1, 23, 198281),
 datetime.datetime(2012, 10, 5, 0, 1, 36, 167030),
 datetime.datetime(2012, 10, 5, 0, 1, 48, 167035),
 datetime.datetime(2012, 10, 5, 0, 2, 0, 167040),
 ...]

You can get the array expression of the dataset by getarray() method.

>>> dataarray = dataset.getarray()
>>> print(dataarray.shape)
(32, 16, 96, 3362)

Note on the start time

It could be slightly minor issue, but the first data in the example above is slightly earlier than the given time. Given time is t0=2012-10-05T00:01:23, but the first data is observed by obstime[0]=2012-10-05T00:01:11. This is because the data at t0 might be included in the data at index 0.

This could produce a problem if t0 falls at a large data gap. For example,

>>> tg0 = datetime.datetime(2006, 12, 5, 0, 0, 0)
>>> tg1 = datetime.datetime(2006, 12, 6, 0, 0, 0)
>>> dataset2 = scidata_util.getarray2d(tg0, tg1)
>>> obstime2 = dataset2.getobstime()
>>> print(obstime2)
[datetime.datetime(2006, 12, 4, 9, 41, 0, 48676),
 datetime.datetime(2006, 12, 5, 5, 21, 31, 41914),
 datetime.datetime(2006, 12, 5, 5, 21, 55, 41914),
 datetime.datetime(2006, 12, 5, 5, 22, 19, 41925),
 datetime.datetime(2006, 12, 5, 5, 22, 43, 41926),
 ...]

As you see, the first data is observed 14 hour before the given time. Thus, this is users’ responsibility to check if the first data (index=0) is valid for the data analysis or not.

One idea is to use the following syntax.

>>> if t0 - obstime2[0] > datetime.timedelta(seconds=12):
...     dataset2 = dataset2.clipped(obstime2[1], obstime2[-1])
>>> obstime2 = dataset2.getobstime()
>>> print(obstime2)
[datetime.datetime(2006, 12, 5, 5, 21, 31, 41914),
 datetime.datetime(2006, 12, 5, 5, 21, 55, 41914),
 datetime.datetime(2006, 12, 5, 5, 22, 19, 41925),
 datetime.datetime(2006, 12, 5, 5, 22, 43, 41926),
 ...]

Reading 2D data at a certain time

You can get the single 2-D data using the index of the time series data, iteration, and get() method. In any cases, the obtained value is a tuple of (datetime, irfpy.imacommon.imascipac.CntMatrix2D).

• Using index

>>> data0 = dataset[0]    # Return the first data.
>>> print(data0)
(datetime.datetime(2012, 10, 5, 0, 1, 11, 198276), <irfpy.imacommon.imascipac.CntMatrix2D object at 0x10f8b16d0>)

• Using iteration

>>> for obstime, data in dataset:
...     print(obstime, data)
2012-10-05 00:01:11.198276 <<class 'irfpy.imacommon.imascipac.CntMatrix2D'>@2012-10-05T00:01:11.198276:MOD=24 >>24<<:POL=14/14:CNTmax=13>
2012-10-05 00:01:23.198281 <<class 'irfpy.imacommon.imascipac.CntMatrix2D'>@2012-10-05T00:01:23.198281:MOD=24 >>24<<:POL=15/15:CNTmax=5>
2012-10-05 00:01:36.167030 <<class 'irfpy.imacommon.imascipac.CntMatrix2D'>@2012-10-05T00:01:36.167030:MOD=24 >>24<<:POL=00/00:CNTmax=2>
2012-10-05 00:01:48.167035 <<class 'irfpy.imacommon.imascipac.CntMatrix2D'>@2012-10-05T00:01:48.167035:MOD=24 >>24<<:POL=01/01:CNTmax=4>
...

• Using get method.

>>> data_at_0100 = dataset.get(datetime.datetime(2012, 10, 5, 1))  # Data at 1 oclock.
>>> print(data_at_0100)
(datetime.datetime(2012, 10, 5, 0, 59, 48, 230124), <irfpy.imacommon.imascipac.CntMatrix2D object at 0x10fad8190>)

Now let’s look at the irfpy.imacommon.imascipac.CntMatrix2D object more carefully.

>>> dat2d = data0[1]
>>> print(dat2d)
<<class 'irfpy.imacommon.imascipac.CntMatrix2D'>@2012-10-05T00:01:11.198276:MOD=24:POL=14/14:CNTmax=13>

You may see that the CntMatrix2D contains the information on the counts and some other misc data. The attributes below includes them.

• t

  Time in datetime object

• matrix

  Matrix of count rate. Shape depends on the mode, but for nominal modes it is (M32,A16,E96).

• polar

  Tuple providing the start polar index and the stop polar index (inclusive)

• mode

  Mode object (irfpy.imacommon.imamode.Mode)

See example below.

>>> print(dat.t)
2012-10-05 00:01:11.198276
>>> print(dat.polar)
[14, 14]
>>> print(dat.mode)
<IMA:Mode-24 (Exm-0) M32/A16/E96/P16>
>>> print(dat.matrix)
[[[0 0 0 ..., 0 0 0]
  [0 0 0 ..., 0 0 0]
  [0 0 0 ..., 0 0 0]
  ...,
  [0 0 0 ..., 0 0 0]
  [0 0 0 ..., 0 0 0]
  [0 0 0 ..., 0 0 0]]
  ...
>>> print(dat.matrix.shape)
(32, 16, 96)

Time series of 3D data

Time series of 3-D data can be obtained by using irfpy.vima.scidata_util.getarray3d() function.

>>> import datetime
>>> from irfpy.vima import scidata_util
>>> t0 = datetime.datetime(2012, 10, 5, 0, 1, 23)
>>> t1 = datetime.datetime(2012, 10, 6, 4, 59, 43)
>>> dataset3d = scidata_util.getarray3d(t0, t1)

The dataset3d is an object of irfpy.imacommon.imascipac.TimeSeriesCntMatrix, as seen below.

>>> print(dataset3d)
<TimeSeriesCntMatrix:len=211 from 2012-10-05 00:01:11.198276 to 2012-10-06 03:41:04.177631>

For the given time range, you can realize 211 3-D data is available.

Of course the functionality works for the mode 25 data.

>>> t0 = datetime.datetime(2008, 10, 5, 2, 1, 23)
>>> t1 = datetime.datetime(2008, 10, 5, 8, 59, 43)
.. >>> dataset = scidata_util.getarray2d(t0, t1)
.. >>> dataset3d = imascipac.convert2to3(dataset)
>>> dataset3d = scidata_util.getarray3d(t0, t1)
>>> print(dataset3d)
<TimeSeriesCntMatrix:len=54 from 2008-10-05 00:38:38.833861 to 2008-10-05 08:56:40.876417>

You can look at the observation time as below; the time resolution is 24 s for this case.

>>> print(dataset3d.getobstime())
[datetime.datetime(2008, 10, 5, 0, 38, 38, 833861),
 datetime.datetime(2008, 10, 5, 6, 7, 3, 868441),
 datetime.datetime(2008, 10, 5, 6, 10, 15, 868476),
 datetime.datetime(2008, 10, 5, 6, 13, 27, 837259),
 datetime.datetime(2008, 10, 5, 6, 16, 39, 831038),
 datetime.datetime(2008, 10, 5, 6, 19, 51, 799822),
 datetime.datetime(2008, 10, 5, 6, 23, 3, 793610),
 ...

You can get the array expression of the dataset by getarray() method.

>>> dataarray = dataset3d.getarray()
>>> print(dataarray.shape)
(32, 16, 96, 8, 54)

Specific 3D data

To get the 3D data of the specific time, you can use the index, or get() method.

>>> print(dataset3d[15])
(datetime.datetime(2008, 10, 5, 6, 51, 51, 881401), <irfpy.imacommon.imascipac.CntMatrix object at ...>)
>>> dat3d15 = dataset3d[15][1]    # Index 0 is time, and index 1 is the data (CntMatrix object)
>>> print(dat3d15.mode)
<IMA:Mode-25 (Exm-1) M32/A16/E96/P 8>

The matrix contains the count rate data.

>>> print(dat3d15.matrix)
[[[[ 0.  0.  0. ...,  0.  0.  0.]
   [ 0.  0.  0. ...,  0.  0.  0.]
   [ 0.  0.  0. ...,  0.  0.  0.]
   ...,
   [ 0.  0.  0. ...,  0.  0.  0.]
   [ 0.  0.  0. ...,  0.  0.  0.]
   [ 0.  1.  0. ...,  0.  0.  0.]]

  [[ 0.  0.  0. ...,  0.  0.  0.]
   [ 0.  0.  0. ...,  0.  0.  0.]
   [ 0.  0.  0. ...,  0.  0.  0.]
   ...,
   [ 0.  0.  0. ...,  0.  0.  0.]
   [ 0.  0.  0. ...,  0.  0.  0.]
   ...

>>> print(dat3d15.matrix.shape)
(32, 16, 96, 8)

Processing data

There are some functions that can convert the object-oriented dataset to a simple numpy array.

irfpy.imacommon.imascipac.timeseries2d_to_energytime_simple() will convert the irfpy.imacommon.imascipac.TimeSeriesCntMatrix2D to (96, N) numpy array. This is useful to make the E-t diagram.

>>> import datetime
>>> from irfpy.vima import scidata_util
>>> t0 = datetime.datetime(2012, 10, 5, 0, 1, 23)
>>> t1 = datetime.datetime(2012, 10, 6, 4, 59, 43)
>>> dataset = scidata_util.getarray2d(t0, t1)
>>> print(dataset)
<TimeSeriesCntMatrix2D:len=3362 from 2012-10-05 00:01:11.198276 to 2012-10-06 03:44:03.177666>
>>> from irfpy.imacommon import imascipac
>>> etmatrix = imascipac.timeseries2d_to_energytime_simple(dataset)     # etmatrix is numpy array.
>>> print(etmatrix)
[[ 0.  0.  0. ...,  0.  0.  0.]
 [ 0.  0.  0. ...,  0.  0.  0.]
 ...,
 [ 0.  1.  0. ...,  0.  0.  0.]
 [ 0.  0.  0. ...,  0.  0.  0.]]
>>> print(etmatrix.shape)
(96, 3362)