as3_delaylook.mainΒΆ


import datetime
import matplotlib

import irfpy.imacommon.ql.panels

matplotlib.rcParams['font.size'] = 7
import matplotlib.pyplot as plt
import numpy as np
import matplotlib.colors as mplcolors
import matplotlib.dates as mpldates

from irfpy.util import timeseries

import irfpy.mima.scidata_util as msu
from irfpy.imacommon import imascipac

from irfpy.mima import massring as mmassring
from irfpy.mima import energy as menergy
from irfpy.mima import fov as mfov
import irfpy.mexpvat.mexspice as ms
ms.init()

from irfpy.mima.ql import panels as qlpanels


class Dataset:
    """ Represent the dataset at the given ``t``.
    """

    def __init__(self, t):
        """ Data for given *t*. """

        self.longdt = datetime.timedelta(hours=12)
        self.middt = datetime.timedelta(hours=3)
        self.shortdt = datetime.timedelta(hours=1)


        ### Count data
        self.longterm2d = msu.getarray2d(t - self.longdt, t + self.longdt, emulate_full=True)
        self.midterm2d = self.longterm2d.clipped(t - self.middt, t + self.middt)
        self.shortterm2d = self.midterm2d.clipped(t - self.shortdt, t + self.shortdt)
        self.shortterm3d = imascipac.convert2to3(self.shortterm2d, emulate_full=True)


        ### At the exact moment.  count3d is CntMatrix (or None)
        t3d = self.shortterm3d.getobstime()
        import bisect
        idx = bisect.bisect(t3d, t)
        if idx == 0:
            self.count3d = None    # No data for the instance
        elif idx == len(t3d) and t - t3d[-1] > 193:
            self.count3d = None
        else:
            self.count3d = self.shortterm3d[idx - 1][1]

        ### Position, refer to SPICE, via global variable ``ms``.

        # PACC=4 pre-assumed here.
        self.pacc = 4
        self.energytable = menergy.get_default_table_v5_late()
        self.mlines = {1: mmassring.massline(1, 4, enestep=self.energytable)[1],
                       2: mmassring.massline(2, 4, enestep=self.energytable)[1],
                       4: mmassring.massline(4, 4, enestep=self.energytable)[1],
                      16: mmassring.massline(16, 4, enestep=self.energytable)[1],
                      32: mmassring.massline(32, 4, enestep=self.energytable)[1],
                      44: mmassring.massline(44, 4, enestep=self.energytable)[1],}
        self.fov = mfov.TableFOV()


def _vec2deg(vec):
    ''' Vector to angle converter.

    :param vec: Vector with (3, N) shape
    :returns: (R, theta, phi) where R is the length, theta is the polar angle
        in degrees (x-y plane to be 0) and phi is the azimuth angle in degrees.

    >>> vec = [[0, 0, 1], [0, 2, 0], [4, 0, 0], [1, 1, 0], [0, 1, 1]]
    >>> vec = np.array(vec).T
    >>> print(vec.shape)
    (3, 5)
    >>> r, t, p = _vec2deg(vec)
    >>> print('{v[0]:.3f} {v[1]:.3f} {v[2]:.3f} {v[3]:.3f} {v[4]:.3f}'.format(v=r))
    1.000 2.000 4.000 1.414 1.414
    >>> print('{v[0]:.3f} {v[1]:.3f} {v[2]:.3f} {v[3]:.3f} {v[4]:.3f}'.format(v=t))
    90.000 0.000 0.000 0.000 45.000
    >>> print('{v[0]:.3f} {v[1]:.3f} {v[2]:.3f} {v[3]:.3f} {v[4]:.3f}'.format(v=p))
    0.000 90.000 0.000 45.000 90.000
    '''
    r = np.sqrt(vec[0] ** 2 + vec[1] ** 2 + vec[2] ** 2)   # (N,)
    vec2 = vec / r    # (3, N) shape
    t = 90 -np.rad2deg(np.arccos(vec2[2]))
    p = np.rad2deg(np.arctan2(vec[1], vec[0]))
#    p = _normalize_phi(p)

    return r, t, p

def _normalize_phi(phi):
    """ Phi range to be -180 and 180.

    >>> phi = [-1000, -600, -200, -180, 0, 45, 180, 200, 350]
    >>> print(_normalize_phi(phi))
    [  80  120  160 -180    0   45 -180 -160  -10]
    """
    p = np.array(phi)
    while p.min() < -180:
        p = np.where(p < -180, p + 360, p)
    while p.max() >= 180:
        p = np.where(p >= 180, p - 360, p)

    return p

def _subplot_pa_matrix(dataset, mass=None, ene=None, ax=None):
    count3d = dataset.count3d

    if ax is None:
        fig = plt.figure()
        ax = fig.add_subplot(111)

    if mass is None:
        mass = [i for i in range(32)]
    if ene is None:
        ene = [i for i in range(96)]

    if count3d is None:
        ax.plot([0])
    else:
        ### count3d.matrix is [M32, A16, E96, P16] matrix
        _matx = count3d.matrix[:, :, ene, :].sum(2)   # (M32, A16, P16)
        az_pol = _matx[mass, :, :].sum(0)   # (A16, P16)
        plt.pcolormesh(np.arange(17) - 0.5, np.arange(17) - 0.5, az_pol.T, vmin=0.1, vmax=100, norm=mplcolors.LogNorm())


    ax.set_xlim(-0.5, 15.5)
    ax.set_ylim(-0.5, 15.5)

    return ax

def _subplot_em_matrix(dataset, polar=None, azim=None, ax=None):
    ### E-M matrix
    if ax is None:
        fig = plt.figure()
        ax = fig.add_subplot(111)

    if polar is None:
        polar = [i for i in range(16)]
    if azim is None:
        azim = [i for i in range(16)]

    count3d = dataset.count3d

    for i in dataset.mlines.keys():
        ax.plot(dataset.mlines[i], dataset.energytable, 'k')

    if count3d is None:   # No instant data
        pass
    else:
        ### count3d.matrix is [M32, A16, E96, P16] matrix.
        _matx = count3d.matrix[:, :, :, polar].sum(3)   # (M32, A16, E96)
        mass_ene = _matx[:, azim, :].sum(1)   # (M32, E96)
        plt.pcolormesh(np.arange(33), dataset.energytable, mass_ene.T, vmin=0.1, vmax=100, norm=mplcolors.LogNorm())

    ax.set_xlim(0, 32)
    ax.set_ylim(1, 16000)
    ax.set_yscale('log')

    return ax

def _mima_quicklook_ima(t, dataset, fig, axes):
    fig.text(0.5, 0.95, t.strftime('%FT%T.%f'), ha='center', fontsize=12)

    ### E-t diagram

    ax = fig.add_axes([0.1, 0.8, 0.3, 0.15])
    axes['etlong'] = ax
    irfpy.imacommon.ql.panels.axis_etdiagram(dataset.energytable, dataset.longterm2d, ax=ax)
    ax.axvline(t, linestyle='--', color='m')
    longformat = mpldates.DateFormatter('%y%m%dT%H%M')
    ax.set_ylabel('Energy')
    ax.set_title('all mass / all azim')
    ax.set_xlim(t - dataset.longdt, t + dataset.longdt)
    ax.xaxis.set_major_formatter(longformat)

    ax = fig.add_axes([0.1, 0.62, 0.3, 0.15])
    axes['etmid'] = ax
    irfpy.imacommon.ql.panels.axis_etdiagram(dataset.energytable, dataset.midterm2d, ax=ax)
    ax.axvline(t, linestyle='--', color='m')
    midformat = mpldates.DateFormatter('%H:%M:%S')
    ax.set_xlim(t - dataset.middt, t + dataset.middt)
    ax.xaxis.set_major_formatter(midformat)

    ax = fig.add_axes([0.1, 0.44, 0.3, 0.15])
    axes['etshort'] = ax
    irfpy.imacommon.ql.panels.axis_etdiagram(dataset.energytable, dataset.shortterm2d, ax=ax)
    ax.axvline(t, linestyle='--', color='m')
    shortformat = mpldates.DateFormatter('%H:%M:%S')
    ax.set_xlim(t - dataset.shortdt, t + dataset.shortdt)
    ax.xaxis.set_major_formatter(shortformat)

    ax = fig.add_axes([0.1, 0.35, 0.3, 0.06])
    axes['etshort1'] = ax
    irfpy.imacommon.ql.panels.axis_etdiagram(dataset.energytable, dataset.shortterm2d, ax=ax, azim=[0, 1, 2, 3])
    ax.axvline(t, linestyle='--', color='m')
    ax.set_ylabel('(Az:0-3)')
    ax.set_xlim(t - dataset.shortdt, t + dataset.shortdt)
    plt.setp(ax.get_xticklabels(), visible=False)

    ax = fig.add_axes([0.1, 0.29, 0.3, 0.06])
    axes['etshort2'] = ax
    irfpy.imacommon.ql.panels.axis_etdiagram(dataset.energytable, dataset.shortterm2d, ax=ax, azim=[4, 5, 6, 7])
    ax.axvline(t, linestyle='--', color='m')
    ax.set_ylabel('(Az:4-7)')
    ax.set_xlim(t - dataset.shortdt, t + dataset.shortdt)
    plt.setp(ax.get_xticklabels(), visible=False)

    ax = fig.add_axes([0.1, 0.23, 0.3, 0.06])
    axes['etshort3'] = ax
    irfpy.imacommon.ql.panels.axis_etdiagram(dataset.energytable, dataset.shortterm2d, ax=ax, azim=[8, 9, 10, 11])
    ax.axvline(t, linestyle='--', color='m')
    ax.set_ylabel('(Az:8-11)')
    ax.set_xlim(t - dataset.shortdt, t + dataset.shortdt)
    plt.setp(ax.get_xticklabels(), visible=False)

    ax = fig.add_axes([0.1, 0.17, 0.3, 0.06])
    axes['etshort4'] = ax
    irfpy.imacommon.ql.panels.axis_etdiagram(dataset.energytable, dataset.shortterm2d, ax=ax, azim=[12, 13, 14, 15])
    ax.axvline(t, linestyle='--', color='m')
    ax.set_ylabel('(Az:12-15)')
    ax.set_xlim(t - dataset.shortdt, t + dataset.shortdt)
    ax.xaxis.set_major_formatter(shortformat)


    ### E-M matrix (for all over Azim & Polar)
    ax = fig.add_axes([0.45, 0.7, 0.15, 0.15])
    axes['emall'] = ax
    _subplot_em_matrix(dataset, ax=ax)

    for mlim in (8, 16, 24):
        ax.axvline(mlim, linestyle='--', color='0.5')
    for elim in (22, 44, 66):
        ax.axhline(dataset.energytable[elim], linestyle='--', color='0.5')

    ax.set_xlabel('Mass channel')
    ax.set_ylabel('Energy [eV]')
    ax.set_title('All az & All pol')

    ### E-M matrix (for selected polar and azim bin)
    ax = fig.add_axes([0.66, 0.44, 0.06, 0.06])
    axes['em11'] = ax
    _subplot_em_matrix(dataset, polar=np.arange(12, 16), azim=np.arange(0, 4), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    ax.set_ylabel('Energy (P12-15)')

    ax = fig.add_axes([0.72, 0.44, 0.06, 0.06])
    axes['em12'] = ax
    _subplot_em_matrix(dataset, polar=np.arange(12, 16), azim=np.arange(4, 8), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.78, 0.44, 0.06, 0.06])
    axes['em13'] = ax
    _subplot_em_matrix(dataset, polar=np.arange(12, 16), azim=np.arange(8, 12), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.84, 0.44, 0.06, 0.06])
    axes['em14'] = ax
    _subplot_em_matrix(dataset, polar=np.arange(12, 16), azim=np.arange(12, 16), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.66, 0.38, 0.06, 0.06])
    axes['em21'] = ax
    _subplot_em_matrix(dataset, polar=np.arange(8, 12), azim=np.arange(0, 4), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    ax.set_ylabel('Energy (P8-11)')

    ax = fig.add_axes([0.72, 0.38, 0.06, 0.06])
    axes['em22'] = ax
    _subplot_em_matrix(dataset, polar=np.arange(8, 12), azim=np.arange(4, 8), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.78, 0.38, 0.06, 0.06])
    axes['em23'] = ax
    _subplot_em_matrix(dataset, polar=np.arange(8, 12), azim=np.arange(8, 12), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.84, 0.38, 0.06, 0.06])
    axes['em24'] = ax
    _subplot_em_matrix(dataset, polar=np.arange(8, 12), azim=np.arange(12, 16), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.66, 0.32, 0.06, 0.06])
    axes['em31'] = ax
    _subplot_em_matrix(dataset, polar=np.arange(4, 8), azim=np.arange(0, 4), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    ax.set_ylabel('Energy (P4-7)')

    ax = fig.add_axes([0.72, 0.32, 0.06, 0.06])
    axes['em32'] = ax
    _subplot_em_matrix(dataset, polar=np.arange(4, 8), azim=np.arange(4, 8), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.78, 0.32, 0.06, 0.06])
    axes['em33'] = ax
    _subplot_em_matrix(dataset, polar=np.arange(4, 8), azim=np.arange(8, 12), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.84, 0.32, 0.06, 0.06])
    axes['em34'] = ax
    _subplot_em_matrix(dataset, polar=np.arange(4, 8), azim=np.arange(12, 16), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.66, 0.26, 0.06, 0.06])
    axes['em41'] = ax
    _subplot_em_matrix(dataset, polar=np.arange(0, 4), azim=np.arange(0, 4), ax=ax)
    ax.set_xlabel('Mass (A0-3)')
    ax.set_ylabel('Energy (P0-3)')

    ax = fig.add_axes([0.72, 0.26, 0.06, 0.06])
    axes['em42'] = ax
    _subplot_em_matrix(dataset, polar=np.arange(0, 4), azim=np.arange(4, 8), ax=ax)
    ax.set_xlabel('Mass (A4-7)')
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.78, 0.26, 0.06, 0.06])
    axes['em43'] = ax
    _subplot_em_matrix(dataset, polar=np.arange(0, 4), azim=np.arange(8, 12), ax=ax)
    ax.set_xlabel('Mass (A8-11)')
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.84, 0.26, 0.06, 0.06])
    axes['em44'] = ax
    _subplot_em_matrix(dataset, polar=np.arange(0, 4), azim=np.arange(12, 16), ax=ax)
    ax.set_xlabel('Mass (A12-15)')
    plt.setp(ax.get_yticklabels(), visible=False)


    ### P-A matrix (for all mass & energy step 0-65)

    ax = fig.add_axes([0.45, 0.35, 0.15, 0.15])
    axes['paall'] = ax

    _subplot_pa_matrix(dataset, mass=None, ene=[i for i in range(66)], ax=ax)

    for lim in (4, 8, 12):
        ax.axhline(lim, color='0.5', linestyle='--')
        ax.axvline(lim, color='0.5', linestyle='--')

    ax.set_xlabel('Azimuth (index)')
    ax.set_ylabel('Polar (index)')
    ax.set_title('All mass & Ene 0-65')

    ### P-A matrix for selected mass & energy bin

    ax = fig.add_axes([0.66, 0.79, 0.06, 0.06])
    axes['pa11'] = ax
    _subplot_pa_matrix(dataset, mass=np.arange(0, 8), ene=np.arange(0, 22), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    ax.set_ylabel('Polar (E0-21)')

    ax = fig.add_axes([0.72, 0.79, 0.06, 0.06])
    axes['pa12'] = ax
    _subplot_pa_matrix(dataset, mass=np.arange(8, 16), ene=np.arange(0, 22), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.78, 0.79, 0.06, 0.06])
    axes['pa13'] = ax
    _subplot_pa_matrix(dataset, mass=np.arange(16, 24), ene=np.arange(0, 22), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.84, 0.79, 0.06, 0.06])
    axes['pa14'] = ax
    _subplot_pa_matrix(dataset, mass=np.arange(24, 32), ene=np.arange(0, 22), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.66, 0.73, 0.06, 0.06])
    axes['pa21'] = ax
    _subplot_pa_matrix(dataset, mass=np.arange(0, 8), ene=np.arange(22, 44), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    ax.set_ylabel('Polar (E22-43)')

    ax = fig.add_axes([0.72, 0.73, 0.06, 0.06])
    axes['pa22'] = ax
    _subplot_pa_matrix(dataset, mass=np.arange(8, 16), ene=np.arange(22, 44), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.78, 0.73, 0.06, 0.06])
    axes['pa23'] = ax
    _subplot_pa_matrix(dataset, mass=np.arange(16, 24), ene=np.arange(22, 44), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.84, 0.73, 0.06, 0.06])
    axes['pa24'] = ax
    _subplot_pa_matrix(dataset, mass=np.arange(24, 32), ene=np.arange(22, 44), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.66, 0.67, 0.06, 0.06])
    axes['pa31'] = ax
    _subplot_pa_matrix(dataset, mass=np.arange(0, 8), ene=np.arange(44, 65), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    ax.set_ylabel('Polar (E44-65)')

    ax = fig.add_axes([0.72, 0.67, 0.06, 0.06])
    axes['pa32'] = ax
    _subplot_pa_matrix(dataset, mass=np.arange(8, 16), ene=np.arange(44, 65), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.78, 0.67, 0.06, 0.06])
    axes['pa33'] = ax
    _subplot_pa_matrix(dataset, mass=np.arange(16, 24), ene=np.arange(44, 65), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.84, 0.67, 0.06, 0.06])
    axes['pa34'] = ax
    _subplot_pa_matrix(dataset, mass=np.arange(24, 32), ene=np.arange(44, 65), ax=ax)
    plt.setp(ax.get_xticklabels(), visible=False)
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.66, 0.60, 0.06, 0.06])
    axes['pa41'] = ax
    _subplot_pa_matrix(dataset, mass=np.arange(0, 8), ene=np.arange(66, 96), ax=ax)
    ax.set_xlabel('Azimuth (M0-7)')
    ax.set_ylabel('Polar* (E66-95)')

    ax = fig.add_axes([0.72, 0.60, 0.06, 0.06])
    axes['pa42'] = ax
    _subplot_pa_matrix(dataset, mass=np.arange(8, 16), ene=np.arange(66, 96), ax=ax)
    ax.set_xlabel('Azimuth (M8-15)')
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.78, 0.60, 0.06, 0.06])
    axes['pa43'] = ax
    _subplot_pa_matrix(dataset, mass=np.arange(16, 24), ene=np.arange(66, 96), ax=ax)
    ax.set_xlabel('Azimuth (M16-23)')
    plt.setp(ax.get_yticklabels(), visible=False)

    ax = fig.add_axes([0.84, 0.60, 0.06, 0.06])
    axes['pa44'] = ax
    _subplot_pa_matrix(dataset, mass=np.arange(24, 32), ene=np.arange(66, 96), ax=ax)
    ax.set_xlabel('Azimuth (M24-31)')
    plt.setp(ax.get_yticklabels(), visible=False)




def _mima_quicklook_spice(t, dataset, fig, axes):

    fov = mfov.SimpleFOV()

    #### SPICE. Overplot.

    et = ms.spice.str2et(t.strftime('%FT%T'))

    # Sun-direction in IMA frame.
    sun, lt = ms.spice.spkezr('SUN', et, 'MEX_ASPERA_IMAS', 'NONE', 'MARS_EXPRESS')
    xyzsun = sun[:3]
    velsun = sun[3:]

    print(xyzsun)
    print(_vec2deg(xyzsun))
    azsun, elsun = fov.numbering_imasframe(xyzsun)
    print(azsun, elsun)


    # Mars direction in IMA frame.
    mars, lt = ms.spice.spkezr('MARS', et, 'MEX_ASPERA_IMAS', 'NONE', 'MARS_EXPRESS')
    xyzmars = mars[:3]
    velmars = mars[3:]

    print(xyzmars)
    print(_vec2deg(xyzmars))
    azmars, elmars = fov.numbering_imasframe(xyzmars)
    print(azmars, elmars)

    for axlbl in ('paall', 'pa11', 'pa12', 'pa13', 'pa14',
                    'pa21', 'pa22', 'pa23', 'pa24',
                    'pa31', 'pa32', 'pa33', 'pa34',):
        _subplot_direction(azsun, elsun, 'orange', axes[axlbl])
        _subplot_direction(azmars, elmars, 'r', axes[axlbl])

    ### Spice, addition (orbit plots)
    trange = datetime.timedelta(hours=3)
    tc = datetime.datetime(t.year, t.month, t.day, t.hour, (t.minute // 10) * 10, 0)

    t0 = tc - trange
    t1 = tc + trange

    ttraj = timeseries.dtrange(t0, t1, datetime.timedelta(minutes=1))  # Trajectory 1 min
    ttick = timeseries.dtrange(t0, t1, datetime.timedelta(minutes=10))   # Tick every 10 min
    tlabel = [tt for tt in ttick if tt.minute == 0]   # Only for the 0 minutes.

    # Get the position in MSO frame.
    msotick = ms.get_positions(ttick) / 3394    # (N, 3) shape, in Rm
    msotraj = ms.get_positions(ttraj) / 3394    # (N', 3) shape, in Rm
    msotc = ms.get_position(t) / 3394   # (3, ) shape in Rm
    msotlabel = ms.get_positions(tlabel) / 3394

    ax = fig.add_axes([0.45, 0.08, 0.15, 0.15])
    axes['pos_xr'] = ax

    ax.plot(msotraj[:, 0], np.sqrt(msotraj[:, 1] ** 2 + msotraj[:, 2] ** 2), 'k-')
    ax.plot(msotick[:, 0], np.sqrt(msotick[:, 1] ** 2 + msotick[:, 2] ** 2), 'k.')
    ax.plot(msotc[0], np.sqrt(msotc[1] ** 2 + msotc[2] ** 2), 'ko')
    ax.text(msotc[0], np.sqrt(msotc[1] ** 2 + msotc[2] ** 2), '{:%T}'.format(t), color="k")
    for tt_i in range(len(tlabel)):
        ax.text(msotlabel[tt_i, 0], np.sqrt(msotlabel[tt_i, 1] ** 2 + msotlabel[tt_i, 2] ** 2),  '{:%T}'.format(tlabel[tt_i]), color="0.5")

    ### Mars.
    ax.plot(np.cos(np.linspace(0, np.pi, 100)), np.sin(np.linspace(0, np.pi, 100)), color='r')

    ax.set_aspect('equal')
    ax.set_xlim(-4, 4)
    ax.set_ylim(0, 4)
    ax.set_xlabel('x [MSO; Rm]')
    ax.set_ylabel('r [MSO; Rm]')


def mima_quicklook(tn):
    fig = plt.figure(figsize=(18, 12))

    axes = {}

    dataset = Dataset(tn)    # dataset.(long|mid|short)term(2|3)d will show the data.

    count3d = dataset.count3d

    if count3d is None:
        t = tn
    else:
        t = dataset.count3d.t

    _mima_quicklook_ima(t, dataset, fig, axes)
    _mima_quicklook_spice(t, dataset, fig, axes)
    return dataset, fig, axes



def main(t0, t1, formats=None):

    if formats is None: formats = []

#    data2d = msu.getarray2d(t0, t1, emulate_full=True)
#    data3d = imascipac.convert2to3(data2d, emulate_full=True)
    data3d = msu.getarray3d(t0, t1)

    tlist = data3d.getobstime()
    print(len(tlist))

    for t in tlist:
        try:
            dataset, fig, axes = mima_quicklook(t)
            for fmt in formats:
                outfile = 'as3_dl_{:%y%m%d-%H%M%S}.{}'.format(t, fmt)
                print(outfile)
                fig.savefig(outfile)
            plt.close(fig)
        except Exception as e:
            print('!' * 80)
            print('ERROR for {}'.format(t))
            print(e)
            print('!' * 80)
        

def _subplot_direction(azsun, elsun, color, ax):

    if elsun >= 15.5:
        ax.plot([azsun], [15.5], color=color, marker='^')
    elif elsun <= -0.5:
        ax.plot([azsun], [-0.5], color=color, marker='v')
    else:
        ax.plot([azsun], [elsun], color=color, marker='o')


def mainexec():
    import argparse
    parser = argparse.ArgumentParser(description='USAGE: %prog t0 t1')

    parser.add_argument('t0', help="Start time")
    parser.add_argument('t1', help="End time")
    parser.add_argument('-F', '--format', default=['png'],
                    help="Specify the output graphics format. Default is png",
                    action='append',
                    dest='format')

    args = parser.parse_args()

    import dateutil.parser
    t0 = dateutil.parser.parse(args.t0)
    t1 = dateutil.parser.parse(args.t1)

    main(t0, t1, formats=list(set(args.format)))

if __name__ == '__main__':
    mainexec()