Skip to content
Snippets Groups Projects
REINFORCE_Timeslices.py 4.14 KiB
Newer Older
#!/usr/bin/env python
# coding: utf-8

# import matplotlib.pyplot as plt   # our plotting module
import matplotlib as mpl
import pandas as pd               # the main HDF5 reader
import numpy as np                # must have
import km3pipe as kp              # some KM3NeT related helper functions
import seaborn as sns              # beautiful statistical plots!
import sys, itertools, aa, collections
from ROOT import EventFile, Det, Timer
import matplotlib.pylab as pylab
pylab.rcParams['figure.figsize'] = 10, 6
import glob
import os, sys, glob, numpy, matplotlib, ROOT
from matplotlib import pyplot
import pylab as P
import csv, math
import matplotlib as plt
get_ipython().magic(u'pylab inline')
import matplotlib.pylab as pylab
pylab.rcParams['figure.figsize'] = 14, 10
import numpy as np
import pylab as P
from numpy import sin, cos, pi, linspace
from numpy.random import randn
from scipy.signal import lfilter, lfilter_zi, filtfilt, butter

from matplotlib.pyplot import plot, legend, show, grid, figure, savefig
pylab.rcParams['figure.figsize'] = 14, 10
from datetime import datetime, timedelta, time


from ROOT import EventFile, Det
import aa
EventFile.read_timeslices = True
EventFile.timeslice_treename = "KM3NET_TIMESLICE_SN"


det = Det("/pbs/throng/km3net/detectors/KM3NeT_00000029_20170920.detx")
f = EventFile("/sps/km3net/users/gdewasse/L0-files/KM3NeT_00000029_00002984.root") #3094 - 2984 #Feel free to add more runs!


R = []
R2 = []
R3 = []
Rsum = [] 
TIME = [] 
RR = []
RR2 = []
RR3 = []
RSUM = []
T = [] 
DATE = []
i = 0 
for evt in f:
    if i == 0 :
        #print('time', evt.t.AsString())
        time_first_event = evt.t.GetSec() + 1e-9 * evt.t.GetNanoSec()
        
    i+=1
    if f.index > 50000: break
    f.get_rates( det ) 
    event_time = evt.t.GetSec() + 1e-9 * evt.t.GetNanoSec()
    # now, det knows the rate on each PMT
    for dom_id, dom in det.doms :
        rates = [ pmt.rate for pmt in dom.pmts ]
        #print( dom_id , rates[:5] )
        if dom_id == 808982077 : #808432835:809544061:808982077:808997793  #some DOM IDs
                R.append(rates[0])
                R2.append(rates[3])
                R3.append(rates[6])
                #print(event_time)
                TIME.append(event_time - time_first_event)
                Rsum.append(sum(rates))
    #print(max(TIME) - min(TIME)) 
    if max(TIME) - min(TIME) > 120: 
        DATE.append(evt.t.AsString())
        T.append(TIME)
        RSUM.append(Rsum)
        RR.append(R)
        RR2.append(R2)
        RR3.append(R3)
        R = []
        R2 = []
        R3 = []
        Rsum = []
        TIME = []
        time_first_event = evt.t.GetSec() + 1e-9 * evt.t.GetNanoSec()
    #print("\n\n\n")


i = 0 
while i < len(RR):
    
    fig = pylab.figure()
    
    plt.plot(T[i], RSUM[i], linestyle = 'None', marker = 'o')
    v = [min(TIME), max(TIME), 190e3, 350e3]#220000, 310000]
    axis(v)
    pylab.xlabel(r'Time [s]', size = 20)
    pylab.ylabel(r'Rate [kHZ]', size = 20)
    pylab.xticks([0,10,20,30, 40,50,60, 70,80, 90,100,110, 120], size=20)
    pylab.yticks([190e3,200e3,210e3, 220e3,230e3,240e3,250e3,260e3,270e3,280e3,290e3,300e3,310e3,320e3,330e3,340e3,350e3],['190', '200', '210','220','230','240','250','260','270','280','290','300','310', '320', '330', '340','350'],size=20)
    plt.grid()
    plt.savefig((str(DATE[i])))
    #fig.save()
    plt.show()
    
    i+=1
    #plt.legend(bbox_to_anchor=(0., 1.02, 1., .102), loc=3,
    #           ncol=2,fontsize = 15, mode="expand", borderaxespad=0.)


INFO_T = []
i = 0 
while i < len(T):
    j = 0 
    info_T = []
    while j < len(T[i]):
        info_T.append([T[i][j], RSUM[i][j], DATE[i]])
        j+=1
    INFO_T.append(info_T)
    i+=1


i = 0 
info_sorted = []
while i < len(INFO_T):
    info_sorted.append(sorted(INFO_T[i], key=lambda u: u[0], reverse=False))
    i+=1


i = 0
f = open('REINFORCE_DATA.txt', 'w')

while i < 2:#len(info_sorted):
    j = 0
    f.write('New event \n')
    f.write(info_sorted[i][0][2])
    while j < len(info_sorted[i]):
        print(info_sorted[i][j][0],info_sorted[i][j][1])
        L = [str(info_sorted[i][j][0]), ' ', str(info_sorted[i][j][1]), '\n']
        f.writelines(L) 
        j+=1
    i+=1
f.close()