source: fact/tools/pyscripts/sandbox/kraehenb/CalFitsPerformance.py@ 13447

#!/usr/bin/python -tt
# ********************************
# Test script for the CalFits class
# 
# written by Thomas Kraehenbuehl, ETH Zurich
# tpk@phys.ethz.ch, +41 44 633 3973
# April 2012
# ********************************

datafilename = '/fact/raw/2012/04/17/20120417_004.fits.gz'
calibfilename = '/fact/raw/2012/04/17/20120417_003.drs.fits.gz'

import numpy as np
import itertools

from ROOT import gSystem
gSystem.Load("calfits_h.so")
from ROOT import *
print "Testing object creation: "
caltest = CalFits(datafilename,calibfilename)
npcalevent  = np.empty( caltest.npix * caltest.nroi, np.float64) #.reshape(caltest.npix ,caltest.nroi)
caltest.SetNpcaldataPtr(npcalevent)

print "Common variables:"
print "ROI: ", caltest.nroi
print "#Pix: ", caltest.npix
print "Number of events: ", caltest.nevents
print

print "Information per Event:"
caltest.GetCalEvent()
print "Event ID: ", caltest.event_id
print "Trigger type: ", caltest.event_triggertype
print "Uncalibrated data: ", caltest.event_data
print "Calibrated data: ", caltest.npcaldata
print "Board times: ", caltest.event_boardtimes
print "Trigger offsets: ", caltest.event_offset
print

print "Examples of other information"
print "Calibfile ROI: ", caltest.calib_nroi
print "Column size BaselineMean: ", caltest.calibfile.GetN("BaselineMean")
print "Datafile ROI: ", caltest.data_nroi
print "Data: ", caltest.datafile.GetN("Data")
print "StartCellData: ", caltest.datafile.GetN("StartCellData")
print "Direct datafile access: ", caltest.datafile.GetN("StartCellData")
print
print "Columns of the datafile: "
caltest.datafile.PrintColumns()

def offsum(x,y):
    return np.sum(npcalevent[x,y+5:y+15])
vecoffsum = np.vectorize(offsum) #specifying otypes=[np.float/double] does not improve speed

while caltest.GetCalEvent():
    if caltest.event_id>30:
        break
    npcalevent = npcalevent.reshape((caltest.npix, caltest.nroi))
    npthr = npcalevent>2.5 #true if above threshold
    npthr = npthr[:,1:-14] * np.logical_not(npthr[:,:-15]) #only true for the zero crossings, shape (1400,285) [smaller due to necessary integration range]
#    print [(x,y) for x,y in zip(npthr.nonzero()[0],npthr.nonzero()[1])] #print the coordinates, range 0-1399,0-284
    
    #Various versions to get the integral, all with approximately the same miserable speed (3 Hz), except for the last one (vectorized function) with ~4.3 Hz
    #Missing: add deadtime after an integration, remove start & end of the array
#    integrals = [np.sum(npcalevent[x,y+5:y+15]) for x,y in np.transpose(npthr.nonzero())]
#    integrals = [np.sum(npcalevent[x,y+5:y+15]) for x,y in zip(npthr.nonzero()[0],npthr.nonzero()[1])]
#    integrals = [np.sum(npcalevent[x,y+5:y+15]) for x,y in itertools.izip(npthr.nonzero()[0],npthr.nonzero()[1])]
#    integrals = map((lambda index_a,index_b: np.sum(npcalevent[index_a,index_b+5:index_b+15])),npthr.nonzero()[0],npthr.nonzero()[1])
    integrals = vecoffsum(npthr.nonzero()[0],npthr.nonzero()[1])
    print len(integrals)
    
#    print caltest.event_id, caltest.event_triggertype, caltest.event_caldata[10]
#    pass

del caltest