cena_massL_cent

Cena line plot.

''' Cena line plot.'''
import matplotlib
matplotlib.use('Agg')
import sys
from numpy import *
from pylab import *
import datetime
from irfpy.util.utc import *

from Sadppac import *
from irfpy.cy1orb.Cy1Orbit import Cy1Orbit
from irfpy.cy1orb.Cy1OrbitNr import Cy1OrbitNr
from irfpy.util.julday import *

if __name__ == '__main__':

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

    orb=int(sys.argv[1])
    print('ORB=%d'%orb, file=sys.stderr)
    f=figure()

    ### Orbit data

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

    o=Cy1Orbit()
    pos=o.getSelenographicOfOrb(orb)
    sza=o.getSzaOfOrb(orb)

    utc2_l=[]
    lon_l=[]
    lat_l=[]
    hei_l=[]

    print(len(pos[0,:]))

    for i in range(0, len(pos[0,:]), 20): 
        utc2 = convert(Julday(pos[0,i]), datetime.datetime)

        utc2_l.append(utc2)
        lon_l.append(pos[1,i])
        lat_l.append(pos[2,i])
        hei_l.append(pos[3,i])

    utc3_l=[]
    sza_l=[]

    for i in range(0, len(sza[0,:]), 20):
        utc3 = convert(Julday(sza[0,i]), datetime.datetime)
        utc3_l.append(utc3)
        sza_l.append(sza[1,i])

    a2.plot(utc2_l, lon_l, label='LON')
    a2.axhline(321, color='b', linestyle='--')
    a2.plot(utc2_l, lat_l, label='LAT')
    a2.axhline(-82, color='g', linestyle='--')
    a2.set_ylim(-90, 360)

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

    a2.legend()

    # Read CENA data


    from irfpy.cena.cena_mass import *

    orbnr=Cy1OrbitNr()
    orbnr.setFromDefaultUri()
    t0=convert(orbnr.getStartTime(orb), datetime.datetime)
    t1=convert(orbnr.getStopTime(orb), datetime.datetime)

    obst=getobstime(timerange=[t0,t1])

    utc_list=[]
    cnt_list=[]

    for t in obst:
        print(t)
        try:
            dat=getdata(t)
        except RuntimeError as e:
            print('Data missing for %s\n***%s'%(t,e))
            continue

        jd=dat.getJd()
        utc=jd2dt(jd)
        cnt_8_7_128 = dat.getData()
        cnt_light_center = cnt_8_7_128[:, 2:5, :64]
        cnt_impo = cnt_light_center.sum(2).sum(1)

        utc_list.append(utc)
        cnt_list.append(cnt_impo)
    

    cnts=array(cnt_list)
    print(cnts.shape)

    off=5

    a=f.add_subplot(211, sharex=a2)
    a.plot(utc_list, cnts[:,7])
    a.plot(utc_list, cnts[:,6]+off*1)
    a.plot(utc_list, cnts[:,5]+off*2)
    a.plot(utc_list, cnts[:,4]+off*3)
    a.plot(utc_list, cnts[:,3]+off*4)
    a.plot(utc_list, cnts[:,2]+off*5)
    a.plot(utc_list, cnts[:,1]+off*6)
    a.plot(utc_list, cnts[:,0]+off*7)

    a.text(utc_list[0], 0, 'E7 (+0)', rotation=10)
    a.text(utc_list[0], off*1, 'E6 (+%d)'%(off*1), rotation=10)
    a.text(utc_list[0], off*2, 'E5 (+%d)'%(off*2), rotation=10)
    a.text(utc_list[0], off*3, 'E4 (+%d)'%(off*3), rotation=10)
    a.text(utc_list[0], off*4, 'E3 (+%d)'%(off*4), rotation=10)
    a.text(utc_list[0], off*5, 'E2 (+%d)'%(off*5), rotation=10)
    a.text(utc_list[0], off*6, 'E1 (+%d)'%(off*6), rotation=10)
    a.text(utc_list[0], off*7, 'E0 (+%d)'%(off*7), rotation=10)

    ceb_time=datetime.datetime(2009,4,18,1,30,0)

    a.axvline(ceb_time, linestyle='--', color='k')
    a2.axvline(ceb_time, linestyle='--', color='k')
    a.text(ceb_time, off*8, 'Cebeus crossing', rotation=0, horizontalalignment='right')
    a.set_xlim(t0, t1)
    a.set_ylim(0,off*8+5)

    a.set_title('Counts: DIR=2,3,4 / MASS=0-63 ')

    f.autofmt_xdate()

    savefig('cena_%04d.png'%orb)

    #try:
    #     __IPYTHON__
    #except:
    #    show()