cena_tser_hflux

Cena line plot. Hydrogen ENA flux.

#!/usr/bin/env python
''' Cena line plot. Hydrogen ENA flux.'''
import matplotlib
matplotlib.use('Agg')
import sys
from numpy import *
from pylab import *
import datetime
import dateutil.parser
from irfpy.util.utc import *

from irfpy.util.julday import *
import irfpy.cy1orb.pvat2 as pvat


if __name__ == '__main__':

    if len(sys.argv) != 3:
        print("USAGE:   %s   t0   t1"%sys.argv[0], file=sys.stderr)
        sys.exit(-1)

    t0 = dateutil.parser.parse(sys.argv[1])
    t1 = dateutil.parser.parse(sys.argv[2])
    print('T = %s til %s' % (t0, t1), file=sys.stderr)

    f=figure()

    ### Orbit data

    dt = datetime.timedelta(seconds=4)  # 4 sec.

    #a2=f.add_subplot(212, sharex=a)
    a2=f.add_subplot(212)

    lon_l = []
    lat_l = []
    hei_l = []
    sza_l=[]

    # Make tlist
    tlist = []
    t = t0
    while t <= t1:
        lon, lat, hei = pvat.getmepos(t)
        sza = pvat.getsza(t)
        tlist.append(t)
        lon_l.append(lon)
        lat_l.append(lat)
        hei_l.append(hei)
        sza_l.append(sza)

        t = t + dt
    

    a2.plot(tlist, lon_l, label='LON')
    a2.plot(tlist, lat_l, label='LAT')
    a2.set_ylim(-90, 360)

    a2.plot(tlist, sza_l, label='SZA')
    a2.axhline(90, color='r', linestyle='--')

    a2.legend()

    # Read CENA data
    from irfpy.cena.cena_mass_time import *

    dat = getfullHdeflux(t0, t1)

    obst, flux = dat  # obst is (N,) shape.  flux is (16, 7, N) shape

    flux = flux[:, 2:5, :].sum(1)  # Integrate CH2 to 4.

    off = 10

    a=f.add_subplot(211, sharex=a2)
    a.plot(obst, flux[11] * off ** 0, 'b', label='E11')
    a.plot(obst, flux[10] * off ** 1, 'g', label='E10')
    a.plot(obst, flux[9] * off ** 2, 'r', label='E9')
    a.plot(obst, flux[8] * off ** 3, 'c', label='E8')
    a.plot(obst, flux[7] * off ** 4, 'm', label='E7')
    a.plot(obst, flux[6] * off ** 5, 'y', label='E6')
    a.plot(obst, flux[5] * off ** 6, 'k', label='E5')
    a.plot(obst, flux[4] * off ** 7, 'b', label='E4')

    a.text(obst[0], 1e3 * off ** 0, 'E11', rotation=10, color='b')
    a.text(obst[0], 1e3 * off**1, 'E10 (x%.2e)'%(off**1), rotation=10, color='g')
    a.text(obst[0], 1e3 * off**2, 'E9 (x%.2e)'%(off**2), rotation=10, color='r')
    a.text(obst[0], 1e3 * off**3, 'E8 (x%.2e)'%(off**3), rotation=10, color='c')
    a.text(obst[0], 1e3 * off**4, 'E7 (x%.2e)'%(off**4), rotation=10, color='m')
    a.text(obst[0], 1e3 * off**5, 'E6 (x%.2e)'%(off**5), rotation=10, color='y')
    a.text(obst[0], 1e3 * off**6, 'E5 (x%.2e)'%(off**6), rotation=10, color='k')
    a.text(obst[0], 1e3 * off**7, 'E4 (x%.2e)'%(off**7), rotation=10, color='b')

    a.set_xlim(t0, t1)
#    a.set_ylim(0,off*8+5)

    a.set_yscale('log')
    a.set_title('Diff. Flux: DIR=2,3,4 / MASS=0-63 ')
    a.set_ylabel('Diff. Flux [# / cm2 sr eV s]')

    f.autofmt_xdate()

    savefig('cena_%s-%s.png' % (t0.strftime('%FT%T'), t1.strftime('%T')))