swim_mass_spectra

A sample code for SWIM mass spectra

''' A sample code for SWIM mass spectra
'''

import os
import sys
import logging
logging.basicConfig()
import datetime
import math

import matplotlib.pyplot as plt
import numpy as np
import scipy as sp
import matplotlib.dates as mdate

from irfpy.swim import swim_mass, swim_start


def main():
    ''' Main code.
    '''
    t0 = datetime.datetime(2009, 6, 24, 3)
    t1 = datetime.datetime(2009, 6, 24, 4)

    tlist = swim_start.getobstime(timerange=[t0, t1])

    ### NO filter
#    f = []

    ### With L2 filter
#    l2 = swim_mass.filter.remove_one_count_event_l2
#    f = [l2,]

    ### With L4 filter
    l4 = swim_mass.filter.remove_one_count_event_l4
    f = [l4,]

    # Start count rate (total)
    cnt0 = np.array([swim_start.getdata(t).getData() for t in tlist]) # Nx16x16
    cntm = np.array([swim_mass.getdata(t, filter=f).getData() for t in tlist]) # Nx16x16x32 

    tlistd = [j.getDatetime() for j in tlist]
    tlistn = [mdate.date2num(j) for j in tlistd]

    fig = plt.figure(figsize=(8, 10))
    ax0 = fig.add_subplot(6, 1, 1)
    ax0.plot_date(tlistd, cnt0.sum(2).sum(1))
    ax0.set_yscale('log')

    ax1 = fig.add_subplot(6, 1, 2)
    x, y = np.meshgrid(tlistn, np.arange(17))
    ax1.pcolor(x, y, np.log10(cnt0.sum(2).T+.1))
    ax1.set_xlim(tlistn[0], tlistn[-1])

    ax3 = fig.add_subplot(6, 1, 3)
    x, y = np.meshgrid(tlistn, np.arange(33))
    ax3.pcolor(x, y, np.log10(cntm.sum(2).sum(1).T+.1))
    ax3.set_xlim(tlistn[0], tlistn[-1])


    # Mass spectra
    t0sw = datetime.datetime(2009, 6, 24, 3, 17)
    t1sw = datetime.datetime(2009, 6, 24, 3, 22)
    tlistsw = swim_mass.getobstime(timerange=[t0sw, t1sw])
    cntmsw = np.array([swim_mass.getdata(t, filter=f).getData() for t in tlistsw])

    # N x E x D x M

    # Solar wind, full
    ax2 = fig.add_subplot(2, 1, 2)
    ax2.plot(np.arange(32), cntmsw.sum(2).sum(1).sum(0), 'b-o', label='3:17-3:22 (SW)')

    # Solar wind proton (E 6,7) and alphs (9, 10)
    ax2.plot(np.arange(16), cntmsw.sum(2).sum(0)[6:8, :16].sum(0), 'b--x', label='SW-H+/SW-al')
    ax2.plot(np.arange(16), cntmsw.sum(2).sum(0)[9:11, :16].sum(0), 'b--x')

    t0ev = datetime.datetime(2009, 6, 24, 3, 28)
    t1ev = datetime.datetime(2009, 6, 24, 3, 40)
    tlistev = swim_mass.getobstime(timerange=[t0ev, t1ev])
    cntmev = np.array([swim_mass.getdata(t, filter=f).getData() for t in tlistev])
    ax2.plot(np.arange(32), cntmev.sum(2).sum(1).sum(0), 'r-o', label='3:28-3:40 (Unknown pop.)')

    ax2.plot(np.arange(16), cntmev.sum(2).sum(0)[7:9, :16].sum(0), 'r--v', label='E7,8 (Lower)')
    ax2.plot(np.arange(16), cntmev.sum(2).sum(0)[10, :16], 'r-.^', label='E10 (Upper)')

    t0bg = datetime.datetime(2009, 6, 24, 3, 43)
    t1bg = datetime.datetime(2009, 6, 24, 3, 55)
    tlistbg = swim_mass.getobstime(timerange=[t0bg, t1bg])
    cntmbg = np.array([swim_mass.getdata(t, filter=f).getData() for t in tlistbg])
    ax2.plot(np.arange(32), cntmbg.sum(2).sum(1).sum(0), 'g-o', label='3:43-3:55 (Background)')

    
    ax2.set_yscale('log')
    ax2.legend(loc='best')

    fig.savefig('swim_mass_spectra.png')

if __name__ == '__main__':
    main()