source: fact/tools/pyscripts/sandbox/vogler/CalFitsPerformanceWeave_4.py@ 16858

#!/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
# ********************************
#
# modified and adapted py Patrick Vogler
#
# ################################
from ROOT import gSystem
gSystem.Load("calfactfits_h.so")   # according to new naming scheme
from ROOT import *

#define filenames

# 2012 04 17
#calibfilename = '/fact/raw/2012/04/17/20120417_003.drs.fits.gz'    # NROI 300, pedestal
calibfilename = '/fact/raw/2012/04/17/20120417_033.drs.fits.gz'    # NROI 300, pedestal

#datafilename = '/fact/raw/2012/04/17/20120417_015.fits.gz'         # NROI 300, LP_ext
#datafilename = '/fact/raw/2012/04/17/20120417_021.fits.gz'         # NROI 300, LP_ext
#datafilename = '/fact/raw/2012/04/17/20120417_036.fits.gz'         # NROI 300, LP_ext
datafilename = '/fact/raw/2012/04/17/20120417_042.fits.gz'         # NROI 300, LP_ext




# 2012 01 25
#datafilename = '/fact/raw/2012/01/25/20120125_042.fits.gz'     # light-pulser-ext
#  = '/fact/raw/2012/01/25/20120125_095.fits.gz'     # pedestal
#  = '/fact/raw/2012/01/25/20120125_094.fits.gz'     # pedestal
#  = '/fact/raw/2012/01/25/20120125_093.fits.gz'     # pedestal
#datafilename = '/fact/raw/2012/01/25/20120125_075.fits.gz'
#calibfilename = '/fact/raw/2012/01/25/20120125_088.drs.fits.gz'


import numpy as np
from scipy import weave
from scipy.weave import converters

from plotters import Plotter         # ADC display
from plotters import CamPlotter      # event display
#from drs_spikes import DRSSpikes


print "Testing object creation: "
caltest = CalFactFits(datafilename, calibfilename)
npcalevent  = np.empty( caltest.npix * caltest.nroi, np.float64) #.reshape(caltest.npix ,caltest.nroi)
caltest.SetNpcaldataPtr(npcalevent)


numroi = np.int64(caltest.nroi)
numpix = np.int64(caltest.npix)
numevents = np.int64(caltest.nevents)



print "Common variables run information: "
print "ROI: ", numroi
print "#Pix: ", numpix
print "Number of events: ", numevents
print


event = np.zeros((numpix, numroi))  # create an array to store an event in the format numpix * numroi (2-dim array)
print "Array event created"
run_average = np.zeros((numpix, numroi))   # create an array to store the "average event" of a run
print "Array run_average created"
extracted_average = np.zeros(numpix)
print "Array extracted_average created "
data_average_run = np.zeros(numroi)
print "Arrays created "



#####################################################################################################################
#print "creating plotters ..."

#make a Plotter class ... this is an easy way for plotting ... but there are 
# many was to plot data ... without this class ... it was written for convenience, but there is no strong reason to use it...
#myplotter = Plotter('titel of the plot', xlabel='time in slices', ylabel='amplitude calibrated data ... in mV')

#myplotter2 = Plotter('titel of the plot', xlabel='time in slices', ylabel='amplitude calibrated data ... in mV')

# make a CamPlotter class
#mycamplotter = CamPlotter('titel of the plot',  map_file_path = '../map_dn.txt', vmin=0, vmax=350)  # for Lightpulser data
#mycamplotter = CamPlotter('titel of the plot',  map_file_path = '../map_dn.txt', vmin=-1.2, vmax=1.2)    # for pedestal data

#print "Plotters created "
###################################################################################################################

print "... looping..."

while caltest.GetCalEvent():    # Loop  ueber alle events
    #print npcalevent  ## Daten, Array   1 dim Laenge caltest.npix * caltest.nroi   1 Event
    event = np.reshape(npcalevent, (numpix, -1))     # bring the event the shape numpix * numroi (2-dim array)

    run_average += event

    # A first test of Adrians idea of primitive signal extractor
    # signal and background
    extracted_signal = np.zeros(numpix)
    for signalslice in range(90, 100):
        extracted_signal += event[0:(numpix), signalslice]
    # background subtraction
    for backgroundslice in range(15, 25):
        extracted_signal -= event[0:(numpix), backgroundslice]
    extracted_average += extracted_signal

print "Looped ... "

for i in range (0, (numpix - 1) ):
    data_average_run += run_average[i]


print "Looped second loop "

del caltest
print "caltest deleted "



print "Common variables run information: "
print "ROI: ", numroi
print "#Pix: ", numpix
print "Number of events: ", numevents
print



#####################################################################################################################
print "creating plotters ..."

#make a Plotter class ... this is an easy way for plotting ... but there are 
# many was to plot data ... without this class ... it was written for convenience, but there is no strong reason to use it...
myplotter = Plotter('titel of the plot', xlabel='time in slices', ylabel='amplitude calibrated data ... in mV')

myplotter2 = Plotter('titel of the plot', xlabel='time in slices', ylabel='amplitude calibrated data ... in mV')

# make a CamPlotter class
mycamplotter = CamPlotter('titel of the plot',  map_file_path = 'map_dn.txt', vmin=0, vmax=350)  # for Lightpulser data
#mycamplotter = CamPlotter('titel of the plot',  map_file_path = '../map_dn.txt', vmin=0, vmax=350)  # for Lightpulser data
#mycamplotter = CamPlotter('titel of the plot',  map_file_path = '../map_dn.txt', vmin=-1.2, vmax=1.2)    # for pedestal data
#mycamplotter = CamPlotter('titel of the plot',  map_file_path = 'map_dn.txt', vmin=-2, vmax=2)    # for pedestal data


print "Plotters created "
###################################################################################################################






myplotter((run_average[22]) / numevents, 'pix 22, average' )

myplotter2( data_average_run / (numpix * numevents) , 'average signal (all pixel) over the whole run' )


#mycamplotter( data[0:1443, 100] )    # plot slice 100 of all pixel
mycamplotter( extracted_average / (10 * numevents) )      # plot the extracted signal of all pixel

answer = raw_input('type "quit" to quit ')
while not 'quit' in answer:
    answer = raw_input('type "quit" to quit ')




del caltest
print "caltest deleted "