| 1 | #!/usr/bin/python -tt
|
|---|
| 2 | # ********************************
|
|---|
| 3 | # Test script for the CalFits class
|
|---|
| 4 | #
|
|---|
| 5 | # written by Thomas Kraehenbuehl, ETH Zurich
|
|---|
| 6 | # tpk@phys.ethz.ch, +41 44 633 3973
|
|---|
| 7 | # April 2012
|
|---|
| 8 | # ********************************
|
|---|
| 9 | #
|
|---|
| 10 | # modified and adapted py Patrick Vogler
|
|---|
| 11 | #
|
|---|
| 12 | # ################################
|
|---|
| 13 | from ROOT import gSystem
|
|---|
| 14 | gSystem.Load("calfactfits_h.so") # according to new naming scheme
|
|---|
| 15 | from ROOT import *
|
|---|
| 16 |
|
|---|
| 17 |
|
|---|
| 18 |
|
|---|
| 19 |
|
|---|
| 20 |
|
|---|
| 21 |
|
|---|
| 22 |
|
|---|
| 23 | #define filenames
|
|---|
| 24 |
|
|---|
| 25 | # 2011 12 05
|
|---|
| 26 | calibfilename = '/fact/raw/2011/12/05/20111205_011.drs.fits.gz' # NROI 300, pedestal
|
|---|
| 27 |
|
|---|
| 28 | datafilename = '/fact/raw/2011/12/05/20111205_022.fits.gz' # NROI 300, 5 minutes physics data with
|
|---|
| 29 |
|
|---|
| 30 |
|
|---|
| 31 | import numpy as np
|
|---|
| 32 | from scipy import weave
|
|---|
| 33 | from scipy.weave import converters
|
|---|
| 34 |
|
|---|
| 35 | from plotters import Plotter # ADC display
|
|---|
| 36 | # from plotters import CamPlotter # event display
|
|---|
| 37 | #from drs_spikes import DRSSpikes
|
|---|
| 38 |
|
|---|
| 39 |
|
|---|
| 40 | print "Testing object creation: "
|
|---|
| 41 | caltest = CalFactFits(datafilename, calibfilename)
|
|---|
| 42 | npcalevent = np.empty( caltest.npix * caltest.nroi, np.float64) #.reshape(caltest.npix ,caltest.nroi)
|
|---|
| 43 | caltest.SetNpcaldataPtr(npcalevent)
|
|---|
| 44 |
|
|---|
| 45 |
|
|---|
| 46 | numroi = np.int64(caltest.nroi)
|
|---|
| 47 | numpix = np.int64(caltest.npix)
|
|---|
| 48 | numevents = np.int64(caltest.nevents)
|
|---|
| 49 |
|
|---|
| 50 | num_LP_ev = 0 # number of ext Lightpulser events
|
|---|
| 51 | num_ped = 0 # number of pedestal events
|
|---|
| 52 |
|
|---|
| 53 |
|
|---|
| 54 | print "Common variables run information: "
|
|---|
| 55 | print "ROI: ", numroi
|
|---|
| 56 | print "#Pix: ", numpix
|
|---|
| 57 | print "Number of events: ", numevents
|
|---|
| 58 | print
|
|---|
| 59 |
|
|---|
| 60 |
|
|---|
| 61 | event = np.zeros((numpix, numroi)) # create an array to store an event in the format numpix * numroi (2-dim array)
|
|---|
| 62 | print "Array event created"
|
|---|
| 63 | run_average = np.zeros((numpix, numroi)) # create an array to store the "average event" of a run
|
|---|
| 64 | data_average_run = np.zeros(numroi)
|
|---|
| 65 | print "Arrays created "
|
|---|
| 66 |
|
|---|
| 67 |
|
|---|
| 68 | print "... looping..."
|
|---|
| 69 |
|
|---|
| 70 | while caltest.GetCalEvent(): # Loop ueber alle events
|
|---|
| 71 | # print npcalevent ## Daten, Array 1 dim Laenge caltest.npix * caltest.nroi 1 Event
|
|---|
| 72 |
|
|---|
| 73 |
|
|---|
| 74 | if (caltest.event_triggertype == 4): #ext LP event
|
|---|
| 75 | print 'event id:', caltest.event_id, ' Trigger Type:', caltest.event_triggertype
|
|---|
| 76 | num_LP_ev = num_LP_ev + 1
|
|---|
| 77 | event = np.reshape(npcalevent, (numpix, -1)) # bring the event the shape numpix * numroi (2-dim array)
|
|---|
| 78 | run_average += event
|
|---|
| 79 |
|
|---|
| 80 | print "Looped ... "
|
|---|
| 81 |
|
|---|
| 82 |
|
|---|
| 83 |
|
|---|
| 84 |
|
|---|
| 85 | for i in range (0, (numpix - 1) ):
|
|---|
| 86 | data_average_run += run_average[i]**2
|
|---|
| 87 | print "Looped second loop "
|
|---|
| 88 | data_average_run = data_average_run/(numpix * num_LP_ev)
|
|---|
| 89 | data_average_run = np.sqrt(data_average_run)
|
|---|
| 90 |
|
|---|
| 91 |
|
|---|
| 92 |
|
|---|
| 93 | del caltest
|
|---|
| 94 | print "caltest deleted "
|
|---|
| 95 |
|
|---|
| 96 |
|
|---|
| 97 |
|
|---|
| 98 | print "Common variables run information: "
|
|---|
| 99 | print "ROI: ", numroi
|
|---|
| 100 | print "#Pix: ", numpix
|
|---|
| 101 | print "Number of events: ", numevents
|
|---|
| 102 | print "Number of external Lightpulser events: ", num_LP_ev
|
|---|
| 103 | print "Number of pedestal: ", num_ped
|
|---|
| 104 | print
|
|---|
| 105 |
|
|---|
| 106 |
|
|---|
| 107 |
|
|---|
| 108 | #####################################################################################################################
|
|---|
| 109 | print "creating plotters ..."
|
|---|
| 110 |
|
|---|
| 111 | #make a Plotter class ... this is an easy way for plotting ... but there are
|
|---|
| 112 | # many was to plot data ... without this class ... it was written for convenience, but there is no strong reason to use it...
|
|---|
| 113 | #myplotter = Plotter('titel of the plot', xlabel='time in slices', ylabel='amplitude calibrated data ... in mV')
|
|---|
| 114 |
|
|---|
| 115 | myplotter = Plotter('titel of the plot', xlabel='time in slices', ylabel='amplitude calibrated data ... in mV')
|
|---|
| 116 | myplotter2 = Plotter('titel of the plot', xlabel='time in slices', ylabel='amplitude calibrated data ... in mV')
|
|---|
| 117 | #myplotter3 = Plotter('titel of the plot', xlabel='time in slices', ylabel='amplitude calibrated data ... in mV')
|
|---|
| 118 |
|
|---|
| 119 | myplotter4 = Plotter('titel of the plot', xlabel='time in slices', ylabel='amplitude calibrated data ... in mV')
|
|---|
| 120 |
|
|---|
| 121 | print "Plotters created "
|
|---|
| 122 | ###################################################################################################################
|
|---|
| 123 |
|
|---|
| 124 | myplotter((run_average[20]) / num_LP_ev, 'pix 20, average' )
|
|---|
| 125 |
|
|---|
| 126 | myplotter2((run_average[21]) / num_LP_ev, 'pix 21, average' )
|
|---|
| 127 |
|
|---|
| 128 | #myplotter3((run_average[22]) / num_LP_ev, 'pix 22, average' )
|
|---|
| 129 |
|
|---|
| 130 | #myplotter3(data_average_run / (numpix * num_LP_ev), 'average signal (all pixel) over the whole run' )
|
|---|
| 131 |
|
|---|
| 132 | myplotter4(data_average_run , 'RMS (all pixel) over the whole run' )
|
|---|
| 133 |
|
|---|
| 134 |
|
|---|
| 135 | answer = raw_input('type "quit" to quit ')
|
|---|
| 136 | while not 'quit' in answer:
|
|---|
| 137 | answer = raw_input('type "quit" to quit ')
|
|---|
| 138 |
|
|---|