#!/usr/local/bin/python -i
import sys, itertools, aa, collections
from ROOT import *

"""
  This example shows how to read time-slices in aanet. 
  
  It plots a delta-t plot for all combinations of PMTs on 
  a certain DOM.

  timeslices are treated as Evt's in aanet. They can be read by 
  setting the 'read_timeslices' static member of EventFile to true. 

  The script also shows the ability to compile some c++ code on-the-fly
  and to use this immediately from within the python program. This is 
  usefull to speed up the time-critical parts.

"""

compile = False

# We define the function fill_delta_t, either in c++ or in python

if not compile :

    def fill_delta_t( hist, hits, dom_id = 808966287 ) :
        "Fill hist will delta-t for all hits on the dom with id = dom_id."
        
        selhits = [ h for h in hits if h.dom_id == 808966287 ]
        
        for pair in itertools.combinations( selhits , 2 ):            
            hdt.Fill( pair[0].t- pair[1].t )

else :

    cxx ("""
           // Fill hist will delta-t for all hits on the dom with id = dom_id.

           void fill_delta_t( TH1& hist, vector<Hit>& hits, int dom_id = 808966287 )
           {
              vector<Hit*> selhits;
              foreach( h, hits ) 
                 {
                    if (h.dom_id == dom_id ) selhits.push_back( &h );
                 }

              for (unsigned i = 0; i<selhits.size(); i++) 
                for (unsigned j = i+1 ; j<selhits.size(); j++)
                  {
                     hist.Fill( selhits[i]->t - selhits[j]->t );
                  }
           }
     """)




# open det and event files

det = Det("/sps/km3net/users/heijboer/example_files/KM3NeT_00000007_02122015_zest_DR_PMT.detx")

EventFile.read_timeslices = True;  # Timeslices mode!
f = EventFile("/in2p3/km3net/data/raw/sea/KM3NeT_00000007/1/KM3NeT_00000007_00001921.root") #irods staged automatically

hdt = TH1D("hdt","time difference;delta-t(ns);entries",200,-100,100)

# --- loop over slices

for evt in f:
    
    print f.index
    if f.index > 99 : break
    det.apply( evt ); # applies timing offsets from detx file
    fill_delta_t ( hdt, evt.hits,  808966287);


hdt.Draw()
hdt.Fit('gaus','','',-7,7)